Java Square Root Integer Operations Without Casting?

33,958

Solution 1

The cast technique is much better than it may seem to be at first sight.

Conversion from int to double is exact.

Math.sqrt is specified, for normal positive numbers, to return "the double value closest to the true mathematical square root of the argument value". If the input was a perfect square int, the result will be an integer-valued double that is in the int range.

Similarly, conversion of that integer-valued double back to an int is also exact.

The net effect is that the cast technique does not introduce any rounding error in your situation.

If your program would otherwise benefit from use of the Guava intMath API, you can use that to do the square root, but I would not add dependency on an API just to avoid the cast.

Solution 2

You can always use the Google Guava IntMath API.

final int sqrtX = IntMath.sqrt(x, RoundingMode.HALF_EVEN);

Solution 3

maybe someone interested, found a link:

http://atoms.alife.co.uk/sqrt/index.html

/*
ATOMS
Fast integer square roots in Java

Here are several fast integer square root methods, written in Java, including:
sqrt(x) agrees completely with (int)(java.lang.Math.sqrt(x)) for x < 2147483648 (i.e. 2^31), while executing about three times faster than it1.

fast_sqrt(x) agrees completely with (int)(java.lang.Math.sqrt(x)) for x < 289 (and provides a reasonable approximation thereafter), while executing about five times faster than it1.

SquareRoot.java - fast integer square root class;
SquareRootTest.java - basic speed and accuracy test class.
Thanks to Paul Hsieh's square root page for the accurate algorithm.
This code has been placed in the public domain.

Can you improve on the speed or accuracy of these methods (without chewing an "unreasonable" quantity of cache with a huge look-up table)?

If so, drop us a line, and hopefully we'll put your name up in lights.

[1: your mileage may vary]



[email protected] | http://atoms.org.uk/
*/

/*
 * Integer Square Root function
 * Contributors include Arne Steinarson for the basic approximation idea, Dann 
 * Corbit and Mathew Hendry for the first cut at the algorithm, Lawrence Kirby 
 * for the rearrangement, improvments and range optimization, Paul Hsieh 
 * for the round-then-adjust idea, Tim Tyler, for the Java port
 * and Jeff Lawson for a bug-fix and some code to improve accuracy.
 * 
 * 
 * v0.02 - 2003/09/07
 */

/**
 * Faster replacements for (int)(java.lang.Math.sqrt(integer))
 */
public class SquareRoot {
  final static int[] table = {
     0,    16,  22,  27,  32,  35,  39,  42,  45,  48,  50,  53,  55,  57,
     59,   61,  64,  65,  67,  69,  71,  73,  75,  76,  78,  80,  81,  83,
     84,   86,  87,  89,  90,  91,  93,  94,  96,  97,  98,  99, 101, 102,
     103, 104, 106, 107, 108, 109, 110, 112, 113, 114, 115, 116, 117, 118,
     119, 120, 121, 122, 123, 124, 125, 126, 128, 128, 129, 130, 131, 132,
     133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 144, 145,
     146, 147, 148, 149, 150, 150, 151, 152, 153, 154, 155, 155, 156, 157,
     158, 159, 160, 160, 161, 162, 163, 163, 164, 165, 166, 167, 167, 168,
     169, 170, 170, 171, 172, 173, 173, 174, 175, 176, 176, 177, 178, 178,
     179, 180, 181, 181, 182, 183, 183, 184, 185, 185, 186, 187, 187, 188,
     189, 189, 190, 191, 192, 192, 193, 193, 194, 195, 195, 196, 197, 197,
     198, 199, 199, 200, 201, 201, 202, 203, 203, 204, 204, 205, 206, 206,
     207, 208, 208, 209, 209, 210, 211, 211, 212, 212, 213, 214, 214, 215,
     215, 216, 217, 217, 218, 218, 219, 219, 220, 221, 221, 222, 222, 223,
     224, 224, 225, 225, 226, 226, 227, 227, 228, 229, 229, 230, 230, 231,
     231, 232, 232, 233, 234, 234, 235, 235, 236, 236, 237, 237, 238, 238,
     239, 240, 240, 241, 241, 242, 242, 243, 243, 244, 244, 245, 245, 246,
     246, 247, 247, 248, 248, 249, 249, 250, 250, 251, 251, 252, 252, 253,
     253, 254, 254, 255
  };

  /**
   * A faster replacement for (int)(java.lang.Math.sqrt(x)).  Completely accurate for x < 2147483648 (i.e. 2^31)...
   */
  static int sqrt(int x) {
    int xn;

    if (x >= 0x10000) {
      if (x >= 0x1000000) {
        if (x >= 0x10000000) {
          if (x >= 0x40000000) {
            xn = table[x >> 24] << 8;
          } else {
            xn = table[x >> 22] << 7;
          }
        } else {
          if (x >= 0x4000000) {
            xn = table[x >> 20] << 6;
          } else {
            xn = table[x >> 18] << 5;
          }
        }

        xn = (xn + 1 + (x / xn)) >> 1;
        xn = (xn + 1 + (x / xn)) >> 1;
        return ((xn * xn) > x) ? --xn : xn;
      } else {
        if (x >= 0x100000) {
          if (x >= 0x400000) {
            xn = table[x >> 16] << 4;
          } else {
            xn = table[x >> 14] << 3;
          }
        } else {
          if (x >= 0x40000) {
            xn = table[x >> 12] << 2;
          } else {
            xn = table[x >> 10] << 1;
          }
        }

        xn = (xn + 1 + (x / xn)) >> 1;

        return ((xn * xn) > x) ? --xn : xn;
      }
    } else {
      if (x >= 0x100) {
        if (x >= 0x1000) {
          if (x >= 0x4000) {
            xn = (table[x >> 8]) + 1;
          } else {
            xn = (table[x >> 6] >> 1) + 1;
          }
        } else {
          if (x >= 0x400) {
            xn = (table[x >> 4] >> 2) + 1;
          } else {
            xn = (table[x >> 2] >> 3) + 1;
          }
        }

        return ((xn * xn) > x) ? --xn : xn;
      } else {
        if (x >= 0) {
          return table[x] >> 4;
        }
      }
    }

    illegalArgument();
    return -1;
  }

  /**
   * A faster replacement for (int)(java.lang.Math.sqrt(x)).  Completely accurate for x < 2147483648 (i.e. 2^31)...
   * Adjusted to more closely approximate 
   * "(int)(java.lang.Math.sqrt(x) + 0.5)"
   * by Jeff Lawson.
   */
  static int accurateSqrt(int x) {
    int xn;

    if (x >= 0x10000) {
      if (x >= 0x1000000) {
        if (x >= 0x10000000) {
          if (x >= 0x40000000) {
            xn = table[x >> 24] << 8;
          } else {
            xn = table[x >> 22] << 7;
          }
        } else {
          if (x >= 0x4000000) {
            xn = table[x >> 20] << 6;
          } else {
            xn = table[x >> 18] << 5;
          }
        }

        xn = (xn + 1 + (x / xn)) >> 1;
        xn = (xn + 1 + (x / xn)) >> 1;
        return adjustment(x, xn);
      } else {
        if (x >= 0x100000) {
          if (x >= 0x400000) {
            xn = table[x >> 16] << 4;
          } else {
            xn = table[x >> 14] << 3;
          }
        } else {
          if (x >= 0x40000) {
            xn = table[x >> 12] << 2;
          } else {
            xn = table[x >> 10] << 1;
          }
        }

        xn = (xn + 1 + (x / xn)) >> 1;

         return adjustment(x, xn);
      }
    } else {
      if (x >= 0x100) {
        if (x >= 0x1000) {
          if (x >= 0x4000) {
            xn = (table[x >> 8]) + 1;
          } else {
            xn = (table[x >> 6] >> 1) + 1;
          }
        } else {
          if (x >= 0x400) {
            xn = (table[x >> 4] >> 2) + 1;
          } else {
            xn = (table[x >> 2] >> 3) + 1;
          }
        }

        return adjustment(x, xn);
      } else {
        if (x >= 0) {
          return adjustment(x, table[x] >> 4);
        }
      }
    }

    illegalArgument();
    return -1;
  }

  private static int adjustment(int x, int xn) {
    // Added by Jeff Lawson:
    // need to test:
    //   if  |xn * xn - x|  >  |x - (xn-1) * (xn-1)|  then xn-1 is more accurate
    //   if  |xn * xn - x|  >  |(xn+1) * (xn+1) - x|  then xn+1 is more accurate
    // or, for all cases except x == 0:
    //    if  |xn * xn - x|  >  x - xn * xn + 2 * xn - 1 then xn-1 is more accurate
    //    if  |xn * xn - x|  >  xn * xn + 2 * xn + 1 - x then xn+1 is more accurate
    int xn2 = xn * xn;

    // |xn * xn - x|
    int comparitor0 = xn2 - x;
    if (comparitor0 < 0) {
      comparitor0 = -comparitor0;
    }

    int twice_xn = xn << 1;

    // |x - (xn-1) * (xn-1)|
    int comparitor1 = x - xn2 + twice_xn - 1;
    if (comparitor1 < 0) { // need to correct for x == 0 case?
      comparitor1 = -comparitor1; // only gets here when x == 0
    }

    // |(xn+1) * (xn+1) - x|
    int comparitor2 = xn2 + twice_xn + 1 - x;

    if (comparitor0 > comparitor1) {
      return (comparitor1 > comparitor2) ? ++xn : --xn;
    }

    return (comparitor0 > comparitor2) ? ++xn : xn;
  }

  /**
  * A *much* faster replacement for (int)(java.lang.Math.sqrt(x)).  Completely accurate for x < 289...
  */
  static int fastSqrt(int x) {
    if (x >= 0x10000) {
      if (x >= 0x1000000) {
        if (x >= 0x10000000) {
          if (x >= 0x40000000) {
            return (table[x >> 24] << 8);
          } else {
            return (table[x >> 22] << 7);
          }
        } else if (x >= 0x4000000) {
          return (table[x >> 20] << 6);
        } else {
          return (table[x >> 18] << 5);
        }
      } else if (x >= 0x100000) {
        if (x >= 0x400000) {
          return (table[x >> 16] << 4);
        } else {
          return (table[x >> 14] << 3);
        }
      } else if (x >= 0x40000) {
        return (table[x >> 12] << 2);
      } else {
        return (table[x >> 10] << 1);
      }
    } else if (x >= 0x100) {
      if (x >= 0x1000) {
        if (x >= 0x4000) {
          return (table[x >> 8]);
        } else {
          return (table[x >> 6] >> 1);
        }
      } else if (x >= 0x400) {
        return (table[x >> 4] >> 2);
      } else {
        return (table[x >> 2] >> 3);
      }
    } else if (x >= 0) {
      return table[x] >> 4;
    }
    illegalArgument();
    return -1;
  }

  private static void illegalArgument() {
    throw new IllegalArgumentException("Attemt to take the square root of negative number");
  }

  /** From http://research.microsoft.com/~hollasch/cgindex/math/introot.html
     * where it is presented by Ben Discoe ([email protected])
     * Not terribly speedy...
     */

  /*
     static int unrolled_sqrt(int x) {
        int v;
        int t = 1<<30;
        int r = 0;
        int s;

        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;} t >>= 2;
        s = t + r; r>>= 1; 
        if (s <= x) { x -= s; r |= t;}

        return r;
     }
  */

  /**
   * Mark Borgerding's algorithm...
   * Not terribly speedy...
   */

  /*
     static int mborg_sqrt(int val) {
        int guess=0;
        int bit = 1 << 15;
        do {
           guess ^= bit;  
           // check to see if we can set this bit without going over sqrt(val)...
           if (guess * guess > val )
              guess ^= bit;  // it was too much, unset the bit...
        } while ((bit >>= 1) != 0);

        return guess;
     }
    */

  /** 
   * Taken from http://www.jjj.de/isqrt.cc
   * Code not tested well...
   * Attributed to: http://www.tu-chemnitz.de/~arndt/joerg.html / email: [email protected]
   * Slow.
   */

  /*
     final static int BITS = 32;
     final static int NN = 0;  // range: 0...BITSPERLONG/2

     final static int test_sqrt(int x) {
        int i;
        int a = 0;                   // accumulator...
        int e = 0;                   // trial product...
        int r;

        r=0;                         // remainder...

        for (i=0; i < (BITS/2) + NN; i++)
        {
           r <<= 2;
           r +=  (x >> (BITS - 2));
           x <<= 2;

           a <<= 1;
           e = (a << 1)+1;

           if(r >= e)
           {
              r -= e;
              a++;
           }
        }

        return a;
     }
  */

  /*
  // Totally hopeless performance...
     static int test_sqrt(int n) {
        float r = 2.0F;
        float s = 0.0F;
        for(; r < (float)n / r; r *= 2.0F);
        for(s = (r + (float)n / r) / 2.0F; r - s > 1.0F; s = (r + (float)n / r) / 2.0F) {
           r = s;
        }

        return (int)s;
     }
    */
}


/*
 * Test Integer Square Root function
 */

public class SquareRootTest {
  public static void main(String args[]) {
    int x = -2;
    debug("" + (x == -2L));

    // perform timings...
    //testSpeed();

    // accuracy test...
    testAccuracy();

    // hang...
    do {} while (true);
  }

  static void testSpeed() {
    int i;
    long newtime;
    long oldtime;
    int N = 10000000;
    int mask = (1 << 16) - 1;

    // SquareRoot.fast_sqrt()...
    oldtime = System.currentTimeMillis();

    for (i = 0; i < N; i++) {
      int temp = SquareRoot.fastSqrt(i & mask);
    }

    newtime = System.currentTimeMillis();
    debug("SquareRoot.fast_sqrt:" + (newtime - oldtime));

    // SquareRoot.sqrt()...
    oldtime = System.currentTimeMillis();

    for (i = 0; i < N; i++) {
     int temp = SquareRoot.sqrt(i & mask);
    }

    newtime = System.currentTimeMillis();
    debug("SquareRoot.sqrt:" + (newtime - oldtime));

    /*
       // SquareRoot.mborg_sqrt()...
       oldtime = System.currentTimeMillis();

       for (i = 0; i < N; i++) {
          temp = SquareRoot.mborg_sqrt(i & mask);
       }

       newtime = System.currentTimeMillis();
       debug("SquareRoot.mborg_sqrt:" + (newtime - oldtime));


       // SquareRoot.test_sqrt()...
       oldtime = System.currentTimeMillis();

       for (i = 0; i < N; i++) {
          temp = SquareRoot.test_sqrt(i & mask);
       }

       newtime = System.currentTimeMillis();
       debug("SquareRoot.test_sqrt:" + (newtime - oldtime));
    */

    // java.lang.Math.sqrt()...
    oldtime = System.currentTimeMillis();

    for (i = 0; i < N; i++) {
     int temp = (int) (java.lang.Math.sqrt(i & mask));
    }

    newtime = System.currentTimeMillis();
    debug("java.lang.Math.sqrt:" + (newtime - oldtime));
  }

  static void testAccuracy() {
    int i;
    int a;
    int b;
    int c;
    int d;
    int e;
    int start = 0;
    int last_wrong_value = 0;

    for (i = start; ; i++) {
      a = (int) (java.lang.Math.sqrt(i));
      b = (int) (java.lang.Math.sqrt(i) + 0.5);
      c = SquareRoot.sqrt(i);
      d = SquareRoot.fastSqrt(i);
      e = SquareRoot.accurateSqrt(i);
      // e = SquareRoot.mborg_sqrt(i);
      // f = SquareRoot.test_sqrt(i);

      if (b != e) {
        //if (a > (last_wrong_value * 1.05)) { // don't print too many wrong values - just a sample...
          last_wrong_value = a;
          debug("N:" + i + " - Math.sqrt:" + a + " - Math.sqrt+:" + b + " - accurateSqrt:" + e + " - sqrt:" + c);
        }
      //}
    }
  }

  final static void debug(String o) {
    System.out.println(o);
  }
}

; ;

Solution 4

If you will only ever compute the square root of a perfect square have you considered the option of pre-calculating them into a table? There are a limited number of integer perfect squares and most JVMs could hold all of them at once without much effort. If space is at a premium I am sure there would be some other option in this direction.

Although there are 65535 of them, if that is just too many you could store one in perhaps four and calculate the ones in between. After all (x+1)^2 = (x+1)(x+1) = x^2 + 2x + 1.

Share:
33,958
user1427661
Author by

user1427661

Updated on February 14, 2020

Comments

  • user1427661
    user1427661 over 4 years

    The Math.sqrt() function takes a double as an argument and returns a double. I'm working with a program that uses a perfect square grid in all circumstances, and I need to obtain the square root in integer form. Is the only way to do this to cast the argument as an int and then return an int-casted double?

  • Alex
    Alex about 11 years
    I'd hesitate to do this for anything other than very small values - I wonder at what size the benefit of not computing the sqrt is offset by the cost of searching in the table
  • Louis Wasserman
    Louis Wasserman about 11 years
    If the values are always perfect squares, you could even use RoundingMode.UNNECESSARY to enforce that.
  • Alex
    Alex about 11 years
    You can always write your own intSqrt method to hide the casting if you need to do it in multiple places, if you're like me and dislike casting for aesthetic reasons :)
  • Patricia Shanahan
    Patricia Shanahan about 11 years
    @Alex An intSqrt method would also be a good place to put a comment on why the Math.sqrt and cast technique is exact.
  • OldCurmudgeon
    OldCurmudgeon about 11 years
    @Alex - I'd hesitate too. I'd certainly gauge the cost of a Map lookup against the calculation. However, perhaps in an embedded environment with no FPU but plenty of memory this would be a good fit.
  • usr
    usr almost 10 years
    This is true and a good analysis. If you're not convinced of this analysis for some reason you can put in an assert ((int)sqrt(x) * (int)sqrt(x) == (int)x) to make 100% sure that all results are correct.
  • Vitruvie
    Vitruvie almost 9 years
    Guava just uses (int)Math.sqrt(x) and then fixes the rounding if you're not rounding down. (source)
  • mrog
    mrog almost 8 years
    If memory allows, an array of integers would provide an extremely fast way to look up square roots.