How to prevent overflow in sqrt?

Question

I have the code

unsigned long long f(unsigned long long x, unsigned long long y){
   return sqrt( ( (5*x*x + 1) * (5*y*y +1) - 1 ) / 4 );         
}

but if x or y is too big, the thing overflows even though the output is supposed to be small. is there a way around this?

Answer 1

Split the argument of the square root up algebraically, maybe?:

return sqrt((3*x*x+1)/2) * sqrt((6*y*y-5)/2);

Or split it up further based on your needs.

If x is big enough you can ignore the +1 and make the first term:

sqrt((3*x*x)/2) = fabs(x) * sqrt(3.0/2.0)

and similarly with y for the second term, making it

sqrt((6*y*y)/2) = fabs(y) * sqrt(3.0);

EDIT: After OP edited his question to be:

return sqrt(((3*x*x+1)*(6*y*y-5)-1)/4);  

You can in fact split things up. You just have to be a little bit more careful. Bottom line is that if x is really big, then the +1 can be ignored. If y is really big then the -5 can be ignored. If both (3*x*x+1) and (6*y*y-5) are positive and either is really big, then the -1 can be ignored. You can use these tips and some additional surrounding logic to break this down a bit further. Like this:

 if(fabs(x) > BIGNUMBER && fabs(y) > BIGNUMBER)
 {
     return fabs(x) * fabs(y) * sqrt(18.0/4.0);
 }
 if(fabs(x) > BIGNUMBER && fabs(y) > 1.0)  // x big and y term positive
 {
     return fabs(x) * sqrt(6*y*y-5) * sqrt(3.0/2.0);
 }
 if(fabs(y) > BIGNUMBER)  // x term positive and y big
 {
     return sqrt(3*x*x+1) * fabs(y) * sqrt(6.0/2.0);
 }
 return sqrt(((3*x*x+1)*(6*y*y-5)-1)/4); 

You can optimize this, but this is just meant to show the point.

Answer 2

I believe that at large x or y, the -1 can be neglected w.r.t (x*x)*(y*y) ... and since your function returns long long, the floating point precision won't matter.

You can check whether x or y are big or not, and accordingly, you can choose to neglect the -1 and do as Chris or Rob say.