indy256's blog

Accuracy of Math.sqrt(long)
By indy256, 10 years ago, In English

Last SRM's div1-250-task raised two questions:

  • how to calculate square root of a perfect square?
  • how to do a perfect square test?

The following two methods solve these problems correctly in Java:

  • long rootOfPerfectSquare(long xx) { return (long) Math.sqrt(xx); }

  • boolean isPerfectSquare(long xx) { long x = (long) Math.sqrt(xx); return x * x == xx; }

What makes correctness non-obvious is implicit conversion of Math.sqrt() argument from long to double, where precision can be lost.

Correctness can be verified using the following code:

for (long x = 0; x <= Math.sqrt(Long.MAX_VALUE); x++)
    if (x != (long) Math.sqrt(x * x)) System.out.println(x);

Is there any explanation of why it works?
Does this method work in other languages, e.g. C++?

Can you prove or disprove the following statement: for any long x >= 0: (long) Math.sqrt(x) ?

  • Vote: I like it
  • +21
  • Vote: I do not like it

Comments Comments (7)
| Write comment?
»
10 years ago, # |
  Vote: I like it +8 Vote: I do not like it

I've tested it in C++.

for(long long x=0; x<=sqrt(LLONG_MAX);x++){
	if(x!=(long long)sqrt(x*x)) printf("%lld\n",x);
}

It works well. But, this is only true when you are calculating square root of a perfect square. In other word, you can't test if number x is perfect square by sqrt(x)*sqrt(x)==x. It's because, in IEEE754 double is defined 52-bit accuracy and there is no way to distinguish x2 + 1 from x2 if it's larger than 252. Both Math.sqrt() in Java and sqrt() in C++ receive double as argument.

→ Reply
»
10 years ago, # |
  Vote: I like it +27 Vote: I do not like it 
boolean isPerfectSquare(long x) {
    long s = (long) Math.sqrt(x);
    while (s * s < x) s++;
    while (s * s > x) s--;
    return s * s == x;
}
→ Reply
  • »
    »
    7 years ago, # ^ |
      Vote: I like it -43 Vote: I do not like it

    I don't think this is necessary. Just consider the following cases.

    Case 1: x is a perfect square If x is a perfect square, then there's a y such that y*y = x. As loujunjie has pointed out, if x is a perfect square, then (long) sqrt(x) will return the exact square root of x, which in this case is y. So in this case, it's correct to say that if x is a perfect square, then (long) sqrt(x) will give us it's exact square root.

    Case 2: x is not a perfect square If x is not a perfect square, then there is no such y such that y*y = x. So it doesn't matter what (long) sqrt(x) gives us, (long)sqrt(x) * (long)sqrt(x) will never be equal to x.

    That is why it is sufficient to do just the following:

    boolean isPerfectSquare(long x) {
    long y = (long) Math.sqrt(x);
    return y * y == x;
}
    → Reply
    • »
      »
      »
      7 years ago, # ^ |
        Vote: I like it +20 Vote: I do not like it

      I don't know if it's necessary or not, and I don't want to know. I just write code that is 100% working and doesn't care about implementation of floating point types.

      → Reply
»
10 years ago, # |
  Vote: I like it +13 Vote: I do not like it

It's hardcoded in the IEEE 754 standard, the result of sqrt(x) will be the closest representable number to the square root of x. If the square root is representable then sqrt(x) will be exact.

→ Reply
  • »
    »
    10 years ago, # ^ |
      Vote: I like it +13 Vote: I do not like it

    In Java if we pass 64-bit long value x to Math.sqrt(x), then x will be implicitly converted to 64-bit double and certainly some precision can be lost. In spite of this, (long) Math.sqrt(x) gives exact answer for perfect squares. How can it be explained?

    I'm not sure how C++ converts long long argument passed to sqrt(): is it converted to double or long double and can there be any loss of precision?

    → Reply
    • »
      »
      »
      10 years ago, # ^ |
      ← Rev. 8   Vote: I like it 0 Vote: I do not like it

      In C, sqrt is defined 

      double sqrt(double);

      In C++98 sqrt is defined

      long double sqrt (long double x);

      So it may work. But long double is machine-dependent. We should avoid it.

      However, in C++11, one overloaded sqrt is defined

      double sqrt (T x);           // additional overloads for integral types

      I think it will work well for long long argument. I'll do some tests.

      UPD1: Test of C

      long long a=2000000000;
long long b=a*a-1;
printf("%lld\n",(long long)sqrt(b));

      result is 2000000000. It's not 1999999999, which is what we supposed.

      UPD2: Test of C++98 and C++11 also give same result, that is 2000000000.

      The C++ reference said:

      "Additional overloads are provided in this header () for the integral types: These overloads effectively cast x to a double before calculations (defined for T being any integral type)."

      → Reply