Issue #13219 has been updated by Nathan Zook.


Following the discussion from the second answer at http://cs.stackexchange.com/questions/37596/arbitrary-precision-integer-square-root-algorithm, I implemented an integer Newton - Raphson on the inverse square root.   This involved correcting an error and clarifying some things.  The results were disappointing, as the asymptotic behavior is clearly beaten by the transcribed Zimmerman algorithm.  The crossover was somewhat beyond 10 ** 1000.  For reference, I made a faster Newton method by using Math.sqrt for the initial guess.

I also adjusted my transcription of the Zimmerman algorithm to reflect the fact that Matz implemented a simple , portable Math.sqrt.  Before we use it, we would need to do serious testing, of course.

One more round of benchmarks.  I have removed the warmup & comparison sections.  I also implemented a newton_faster method, which uses Math.sqrt for the initial guess.

~~~
integer squareroot tests for n = 10**50
Calculating -------------------------------------
     newtons_fast(n)    176.714k ( 1.3%) i/s -    891.495k in   5.045686s
    newton_faster(n)    389.709k ( 0.8%) i/s -      1.952M in   5.009508s
           sqrt_z(n)    166.894k ( 1.4%) i/s -    847.715k in   5.080312s
   inverse Newton(n)    250.008k ( 0.4%) i/s -      1.271M in   5.083042s

integer squareroot tests for n = 10**500
Calculating -------------------------------------
     newtons_fast(n)     31.530k ( 3.4%) i/s -    158.496k in   5.032977s
    newton_faster(n)     61.087k ( 3.1%) i/s -    306.020k in   5.014480s
           sqrt_z(n)     33.167k ( 0.5%) i/s -    167.128k in   5.039057s
   inverse Newton(n)     51.252k ( 3.6%) i/s -    257.972k in   5.041610s

integer squareroot tests for n = 10**1000
Calculating -------------------------------------
     newtons_fast(n)     14.085k ( 0.8%) i/s -     71.656k in   5.087601s
    newton_faster(n)     24.078k ( 3.3%) i/s -    121.074k in   5.034512s
           sqrt_z(n)     25.410k ( 2.3%) i/s -    127.500k in   5.020581s
   inverse Newton(n)     28.921k ( 2.7%) i/s -    145.350k in   5.030223s

integer squareroot tests for n = 10**2000
Calculating -------------------------------------
     newtons_fast(n)      3.994k ( 0.7%) i/s -     20.247k in   5.069925s
    newton_faster(n)      6.698k ( 0.7%) i/s -     34.017k in   5.079015s
           sqrt_z(n)     19.176k ( 2.3%) i/s -     96.951k in   5.058731s
   inverse Newton(n)  13.724k ( 1.7%) i/s -     68.952k in   5.025785s

integer squareroot tests for n = 10**4000
Calculating -------------------------------------
     newtons_fast(n)      1.035k ( 0.4%) i/s -      5.253k in   5.074204s
    newton_faster(n)      1.600k ( 0.9%) i/s -      8.000k in   5.001751s
           sqrt_z(n)     12.478k ( 0.7%) i/s -     62.475k in   5.006905s
   inverse Newton(n)      6.094k ( 0.7%) i/s -     30.957k in   5.079890s

integer squareroot tests for n = 10**5000
Calculating -------------------------------------
     newtons_fast(n)    628.308  ( 0.6%) i/s -      3.162k in   5.032783s
    newton_faster(n)    915.131  ( 0.4%) i/s -      4.641k in   5.071494s
           sqrt_z(n)     10.107k ( 0.6%) i/s -     50.796k in   5.026155s
   inverse Newton(n)      4.539k ( 3.0%) i/s -     22.850k in   5.038822s

~~~

It is possible that "inverse Newton" could be modified to advantage the particulars of our implementation of bignum, or GMP, should it be linked.  I don't want spam any more performances after this until it is determined if we want to drop down to C for this.  

One last thing to play with might be Halley's method.  This converges cubically, but requires 4 multiplies.  This should compute faster.  But don't expect to see it very soon.

----------------------------------------
Feature #13219: bug in Math.sqrt(n).to_i, to compute integer squareroot,  new word to accurately fix it
https://bugs.ruby-lang.org/issues/13219#change-63155

* Author: Jabari Zakiya
* Status: Open
* Priority: Normal
* Assignee: 
* Target version: 
----------------------------------------
In doing a math application using **Math.sqrt(n).to_i** to compute integer squareroots 
of integers I started noticing errors for numbers > 10**28.


I coded an algorithm that accurately computes the integer squareroot for arbirary sized numbers
but its significantly slower than the math library written in C.

Thus, I request a new method **Math.intsqrt(n)** be created, that is coded in C and part of the
core library, that will compute the integer squareroots of integers accurately and fast.

Here is working highlevel code to accurately compute the integer squareroot.

```
def intsqrt(n)
  bits_shift = (n.to_s(2).size)/2 + 1
  bitn_mask = root = 1 << bits_shift
  while true
    root ^= bitn_mask if (root * root) > n
    bitn_mask >>= 1
    return root if bitn_mask == 0
    root |= bitn_mask
  end
end

def intsqrt1(n)
  return n if n | 1 == 1   # if n is 0 or 1
  bits_shift = (Math.log2(n).ceil)/2 + 1
  bitn_mask = root = 1 << bits_shift
  while true
    root ^= bitn_mask if (root * root) > n
    bitn_mask >>= 1
    return root if bitn_mask == 0
    root |= bitn_mask
  end
end

require "benchmark/ips"

Benchmark.ips do |x|
  n = 10**40
  puts "integer squareroot tests for n = #{n}"
  x.report("intsqrt(n)"       ) { intsqrt(n)  }
  x.report("intsqrt1(n)"      ) { intsqrt1(n) }
  x.report("Math.sqrt(n).to_i") { Math.sqrt(n).to_i }
  x.compare!
end
```
Here's why it needs to be done in C.

```
2.4.0 :178 > load 'intsqrttest.rb'
integer squareroot tests for n = 10000000000000000000000000000000000000000
Warming up --------------------------------------
          intsqrt(n)     5.318k i/100ms
         intsqrt1(n)     5.445k i/100ms
   Math.sqrt(n).to_i   268.281k i/100ms
Calculating -------------------------------------
          intsqrt(n)     54.219k ( 5.5%) i/s -    271.218k in   5.017552s
         intsqrt1(n)     55.872k ( 5.4%) i/s -    283.140k in   5.082953s
   Math.sqrt(n).to_i      5.278M ( 6.1%) i/s -     26.560M in   5.050707s

Comparison:
   Math.sqrt(n).to_i:  5278477.8 i/s
         intsqrt1(n):    55871.7 i/s - 94.47x  slower
          intsqrt(n):    54219.4 i/s - 97.35x  slower

 => true 
2.4.0 :179 > 

```
Here are examples of math errors using **Math.sqrt(n).to_i** run on Ruby-2.4.0.

```
2.4.0 :101 > n = 10**27; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c   
1000000000000000000000000000
31622776601683
999999999999949826038432489
31622776601683
999999999999949826038432489
31622776601683
999999999999949826038432489
 => nil 
2.4.0 :102 > n = 10**28; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c 
10000000000000000000000000000
100000000000000
10000000000000000000000000000
100000000000000
10000000000000000000000000000
100000000000000
10000000000000000000000000000
 => nil 
2.4.0 :103 > n = 10**29; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c 
100000000000000000000000000000
316227766016837
99999999999999409792567484569
316227766016837
99999999999999409792567484569
316227766016837
99999999999999409792567484569
 => nil 
2.4.0 :104 > n = 10**30; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c  
1000000000000000000000000000000
1000000000000000
1000000000000000000000000000000
1000000000000000
1000000000000000000000000000000
1000000000000000
1000000000000000000000000000000
 => nil 
2.4.0 :105 > n = 10**31; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c 
10000000000000000000000000000000
3162277660168379
9999999999999997900254631487641
3162277660168379
9999999999999997900254631487641
3162277660168379
9999999999999997900254631487641
 => nil 
2.4.0 :106 > n = 10**32; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c 
100000000000000000000000000000000
10000000000000000
100000000000000000000000000000000
10000000000000000
100000000000000000000000000000000
10000000000000000
100000000000000000000000000000000
 => nil 
2.4.0 :107 > n = 10**33; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c
1000000000000000000000000000000000
31622776601683793
999999999999999979762122758866849
31622776601683793
999999999999999979762122758866849
31622776601683792
999999999999999916516569555499264
 => nil 
2.4.0 :108 > n = 10**34; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c
10000000000000000000000000000000000
100000000000000000
10000000000000000000000000000000000
100000000000000000
10000000000000000000000000000000000
100000000000000000
10000000000000000000000000000000000
 => nil 
2.4.0 :109 > n = 10**35; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c
100000000000000000000000000000000000
316227766016837933
99999999999999999873578871987712489
316227766016837933
99999999999999999873578871987712489
316227766016837952
100000000000000011890233980627554304
 => nil 
2.4.0 :110 > n = 10**36; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c
1000000000000000000000000000000000000
1000000000000000000
1000000000000000000000000000000000000
1000000000000000000
1000000000000000000000000000000000000
1000000000000000000
1000000000000000000000000000000000000
 => nil 
2.4.0 :111 > n = 10**37; puts n, (a = intsqrt(n)), a*a, (b = intsqrt1(n)), b*b, (c = Math.sqrt(n).to_i), c*c
10000000000000000000000000000000000000
3162277660168379331
9999999999999999993682442519108007561
3162277660168379331
9999999999999999993682442519108007561
3162277660168379392
10000000000000000379480317059650289664
 => nil 
2.4.0 :112 > 
```



-- 
https://bugs.ruby-lang.org/

Unsubscribe: <mailto:ruby-core-request / ruby-lang.org?subject=unsubscribe>
<http://lists.ruby-lang.org/cgi-bin/mailman/options/ruby-core>