Andrew / Joe, Thank you for the much appreciated comments.
On Oct 3, 2015, at 11:45 AM, joe darcy <joe.da...@oracle.com> wrote: > For an initial implementation, I think it is acceptable to use a simple > algorithm that is clearly correct. That algorithm can then be replaced with > faster ones once adequate tests are built up (with the original > implementation perhaps retiring to the regression tests area where it can be > used for comparison purposes). This was in fact my intention for the initial implementation: something simple that is numerically accurate. On Oct 3, 2015, at 2:38 AM, Andrew Haley <a...@redhat.com> wrote: > On 02/10/15 21:41, Brian Burkhalter wrote: > >> Any suggestions as to the improvement of the approach concerning >> memory use or any other aspect of the algorithm would be >> appreciated, as would be suggestions regarding the test. > > This algorithm doesn't look like the best to me because it's got this > slow division in the middle. If we calculate 1/sqrt(x) instead we can > use a version of Newton's iteration which requires no division. > > Starting with an approximation y = 1/sqrt(x) using the first few digits, > > y = y * (3 - x*y^2) > --------------- > 2 > > ... and fix it all up with a single iteration of Heron's algorithm at the > end. I will investigate the foregoing and Louis suggestion for this pass of the implementation. Thanks for the most welcome ideas. > But even better is Karatsuba Square Root: > https://hal.inria.fr/inria-00072854/document > > I guess it all depends on how much work we think this deserves. In > general the core algorithms in BigInteger are among the best, and it > would be nice to continue this tradition. The Karatsuba Square Root approach by Zimmerman was in fact my original intention and you will find a link to the same work already in the issue on file. On Oct 3, 2015, at 11:58 AM, joe darcy <joe.da...@oracle.com> wrote: > Here is a trick to consider for future performance tuning: all large > floating-point numbers are integers. Once the size of the positive exponent > exceeds the number of bits of precision, the value must be an integer. For > double, that means values greater than 2^54 are integers and the full double > exponent range goes out to 2^1023. > > Consequently, if a BigInteger is less than 2^1023 and the spread of 1 bits in > the BigInteger is contained within a double , the result could be directly > calculated as > > BigInteger(Math.floor(Math.sqrt(bi.doubleValue()))) > > without needing any iterations on the BigInteger side. The bit spread could > be able to be calculated from something like > > bi.bitLenth() - bi.getLowestSetBit() Good observation: I’ll look into this as well. Thanks for all for the ideas. Brian