Issue #13219 has been updated by Akira Tanaka. Jabari Zakiya wrote: > In the library file **prime.rb** is the method **Integer#prime?**. > > ``` > class Integer > ... > ... > def prime? > return self >= 2 if self <= 3 > return true if self == 5 > return false unless 30.gcd(self) == 1 > (7..Math.sqrt(self).to_i).step(30) do |p| > return false if > self%(p) == 0 || self%(p+4) == 0 || self%(p+6) == 0 || self%(p+10) == 0 || > self%(p+12) == 0 || self%(p+16) == 0 || self%(p+22) == 0 || self%(p+24) == 0 > end > true > end > ... > ... > ``` > > As in typical prime algorithms, you just check for the primes <= the number's squareroot. > Here **Math.sqrt(self).to_i** determines the integer squareroot, but as we know (and as I > found out, which initiated all this), as the numbers gets ``> ~10**28 ``, the approximation for > the squareroot becomes more and more larger than the actual value. This will cause **prime?** > to do increasingly more work than necessary as the numbers get bigger. > (Completely as an aside, for numbers this size you should probably look to use: n.to_bn.prime?) Math.sqrt(n) can be smaller than mathematical sqrt(n). ``` % ruby -ve '1.upto(50) {|i| j = 10**i; j2 = Math.sqrt(j*j).to_i; p [i, j2 - j] if j != j2 }' ruby 2.5.0dev (2017-02-21 trunk 57675) [x86_64-linux] [23, -8388608] [24, -16777216] [25, 905969664] [26, 4764729344] [27, 13287555072] [28, -416880263168] [29, -8566849142784] [30, 19884624838656] [31, -364103705034752] [32, 5366162204393472] [33, -54424769012957184] [34, -544247690129571840] [35, -3136633892082024448] [36, 42420637374017961984] [37, -461237341797878857728] [38, -2251190176543965970432] [39, -60290833628396821413888] [40, 303786028427003666890752] [41, 620008645040778319495168] [42, 44885712678075916785549312] [43, 139372116959414099130712064] [44, -10985679223259661757684121600] [45, -70242710975464448780069240832] [46, -68601809640529787052340805632] [47, 4384584304507619735463404765184] [48, 43845843045076197354634047651840] [49, -535097230524518206803127585210368] [50, 7629769841091887003294964970946560] ``` So, Integer#prime? can underestimate the limit to search factors. Assume i is 1000000000000000031. I think (i*i).prime? returns true. 1000000000000000031 is not searched because Math.sqrt(i*i).to_i is 1000000000000000000. (I couldn't complete (i*i).prime? because it takes too long time.) ``` % ruby -e 'i=1000000000000000031; p i, Math.sqrt(i**2).to_i' 1000000000000000031 1000000000000000000 ``` I feel this is a good evidence that we need integer sqrt. > It would be interesting to see what an audit of the entire Ruby code base would reveal where > this is used to provide the integer squareroot, versus the floating point one. Simple search result of 'sqrt\(.*to_i' using gem-codesearch. https://github.com/akr/gem-codesearch ``` Oompa-euler-1.0.3/lib/euler.rb: max = Math.sqrt(to_factor).to_i + 1 Oompa-euler-1.0.3/lib/euler.rb: 5.step(Math.sqrt(self).to_i, 2) do |x| Oompa-euler-1.0.3/lib/euler.rb: 3.step(Math.sqrt(max).to_i, 2) { |i| algebra-0.2.3/lib/algebra/polynomial-factor-int.rb: (Math.log(c * Math.sqrt(s))/Math.log(u)).to_i + 1 blackwinter-gsl-1.15.3.2/examples/random/randomwalk.rb:sigma = Math::sqrt(N).to_i boggle-0.0.2/bin/boggle:opts[:size] = Math.sqrt(opts[:board].size).to_i if !opts[:board].empty? boggler-0.0.1/lib/boggler/grid.rb: @size = Math.sqrt(@dice.length).to_i boqwij-0.1.0/lib/boqwij/integer_extensions.rb: 5.step(Math.sqrt(self).to_i, 2) do |x| boqwij-0.1.0/lib/boqwij/integer_extensions.rb: max = Math.sqrt(to_factor).to_i + 1 clusterer-0.1.9/lib/clusterer/clustering.rb: options[:no_of_clusters] ||= Math.sqrt(objects.size).to_i collective-0.2.1/lib/collective/checker.rb: @estimate > 0 ? Math.sqrt(@estimate).to_i : 1 cw-0.3.3/lib/cw/read.rb: # @magnitude = Math.sqrt(magnitude_squared).to_i euler-1.1.0/lib/euler.rb: max = Math.sqrt(to_factor).to_i + 1 euler-1.1.0/lib/euler.rb: 5.step(Math.sqrt(self).to_i, 2) do |x| euler-1.1.0/lib/euler.rb: 3.step(Math.sqrt(max).to_i, 2) { |i| euler-1.1.0/lib/rudoku.rb: @sqrt_n = Math.sqrt(@n).to_i gecoder-1.1.1/lib/gecoder/interface/constraints/int/arithmetic.rb: Gecode::Raw::sqrt(@model.active_space, int_op.to_int_var.bind, gecoder-with-gecode-1.1.1.1/lib/gecoder/interface/constraints/int/arithmetic.rb: Gecode::Raw::sqrt(@model.active_space, int_op.to_int_var.bind, gphys-1.5.1/lib/numru/gphys/interpolate.rb: ndiv[i] = Math.sqrt( (kx[i+1]-kx[i])**2 + (ky[i+1]-ky[i])**2).to_i green_shoes-1.1.374/samples/sample56.rb: .text(lambda {|d| d.class_name[0, Math.sqrt(d.node_value).to_i / 8]}) gsl-2.1.0.2/examples/random/randomwalk.rb:sigma = Math::sqrt(N).to_i gsl-nmatrix-1.15.3.2/examples/random/randomwalk.rb:sigma = Math::sqrt(N).to_i hotcocoa-0.6.3/lib/hotcocoa/graphics/elements/sandpainter.rb: #@grains = (sqrt((ox-x)*(ox-x)+(oy-y)*(oy-y))).to_i huge_enumerable-0.0.2/lib/huge_enumerable.rb: increment = next_prime(( 2 * domain_size / (1 + Math.sqrt(5)) ).to_i) ifmapper-2.0.4/lib/IFMapper/FXSpline.rb: num = Math.sqrt( x2 + y2 ).to_i / 16 jwthumbs-0.1.0/lib/shutter.rb: gridsize = Math.sqrt(@images.length).ceil.to_i lazylist-0.3.2/examples/examples.rb: not (2..Math.sqrt(x).to_i).find { |d| x % d == 0 } magic_cloud-0.0.4/lib/magic_cloud/layouter/place.rb: case (Math.sqrt(1 + 4 * sign * t) - sign).to_i & 3 mapbaidu-0.0.7/lib/map_baidu.rb: x = Math.sqrt(radius**2 - y**2).ceil.to_i mapcache-0.1.1/lib/mgmaps_export.rb: tpf_x = Math.sqrt(tpf).to_i mapcache-0.1.1/lib/mgmaps_export.rb: tpf_x = Math.sqrt(tpf * 2).to_i matrix_disp-0.0.3/lib/matrix_disp.rb: for i in (0...Math.sqrt(other[0].size).to_i) do matrix_disp-0.0.3/lib/matrix_disp.rb: for x in (0...Math.sqrt(other[0].size).to_i) do matrix_disp-0.0.3/lib/matrix_disp.rb: for y in (0...Math.sqrt(other[0].size).to_i) do matrix_disp-0.0.3/lib/matrix_disp.rb: adj[z] = other[0][x+(y*Math.sqrt(other[0].size).to_i)] matrix_disp-0.0.3/lib/matrix_disp.rb: resultado = resultado + tmp[i+(j*Math.sqrt(other[0].size).to_i)]*det(adj) matrix_disp-0.0.3/lib/matrix_disp.rb: resultado = resultado - tmp[i+(j*Math.sqrt(other[0].size).to_i)]*det(adj) minesweeper-console-1.0.1/test/minesweeper/console/prettyprinter/minefield_pretty_printer_test.rb: Math.sqrt(@string_representation.length).to_i mlanett-hive-0.4.0/lib/hive/checker.rb: @estimate > 0 ? Math.sqrt(@estimate).to_i : 1 mtah-ruby-treemap-0.0.3.2/lib/treemap/node.rb: (base_size * Math.sqrt(self.bounds.width * self.bounds.height) / 125).to_i naga_perfectsquare-1.0.1/lib/naga_perfectsquare.rb: Math.sqrt(self) == Math.sqrt(self).to_i narray-0.6.1.2/src/lib/narray/narray_ext.rb: m = ((n+Math::sqrt(n))*1.27).to_i narray-bigmem-0.0.0/lib/narray_ext.rb: m = ((n+Math::sqrt(n))*1.27).to_i nsudoku-0.2.0/lib/nsudoku/checker.rb: @width = Math.sqrt(@sudoku.length).to_i nsudoku-0.2.0/lib/nsudoku/checker.rb: @block_width = Math.sqrt(@width).to_i nsudoku-0.2.0/lib/nsudoku/solver.rb: @width = Math.sqrt(sudoku.length).to_i nsudoku-0.2.0/lib/nsudoku/solver.rb: @block_width = Math.sqrt(@width).to_i nsudoku-0.2.0/lib/nsudoku.rb: @width = Math.sqrt(@sudoku.length).to_i nsudoku-0.2.0/lib/nsudoku.rb: @block_width = Math.sqrt(@width).to_i numb-0.186.0/lib/numb/sum_of_squares.rb: s = Math.sqrt(i).to_i numeric_memoist-0.0.2/examples/isqrt_example.rb: $stderr.puts "calculate sqrt(#{self}).to_i (<= MaxIntFloat)" numeric_memoist-0.0.2/examples/isqrt_example.rb: return Math.sqrt(self).to_i numeric_memoist-0.0.2/examples/isqrt_example.rb: $stderr.puts "calculate sqrt(#{self}).to_i" numeric_memoist-0.0.2/examples/isqrt_example.rb: # STDERR> calculate sqrt(123456789).to_i (<= MaxIntFloat) numeric_memoist-0.0.2/examples/isqrt_example.rb: # STDERR> calculate sqrt(986960440108935861883449099987613301336842162813430234017512025).to_i numru-narray-1.0.3/lib/numru/narray_ext.rb: m = ((n+Math::sqrt(n))*1.27).to_i pixeldistance-0.2.1/lib/pixeldistance.rb: Math.sqrt((x1-x2)**2 + (y1-y2)**2).to_i >> (21 - zoom) prime-multiplication-1.0/lib/PrimeMultiplication.rb: (2..(Math.sqrt(number).to_i)).each do |i| prime_numbers-0.0.4/lib/prime_numbers/prime_int.rb: 2.upto(Math.sqrt(num).to_i) do |x| prime_table-0.0.2/lib/prime_table/prime.rb: sqrt = Math.sqrt(self).to_i primes-utils-2.7.0/lib/primes/utils.rb: sqrtN = Math.sqrt(n).to_i primes-utils-2.7.0/lib/primes/utils.rb: sqrtN = Math.sqrt(n).to_i primes-utils-2.7.0/lib/primes/utils.rb: factors << p; r -=1; n /= p; sqrtN = Math.sqrt(n).to_i primes-utils-2.7.0/lib/primes/utils.rb: sqrtN = Math.sqrt(num).to_i # sqrt of end_num (end of range) primes-utils-2.7.0/lib/primes/utils.rb: if start_num <= Math.sqrt(num).to_i # for one array of primes upto N print_square-0.0.2/lib/print_square/command_runner.rb: size = Math.sqrt(number).to_i quasi-0.0.6/lib/halton.rb: m=Math.sqrt(q).to_int rlsm-1.8.1/lib/rlsm/monoid.rb: @order = Math::sqrt(@table.size).to_i royw-drbman-0.0.6/examples/primes/lib/primes/sieve_of_eratosthenes.rb: sqr_primes = primes(Math.sqrt(maximum).to_i, drbman) ruby-compiler-0.1.1/vendor/ruby/lib/prime.rb: (7..Math.sqrt(self).to_i).step(30) do |p| ruby-spark-1.2.1/benchmark/comparison/ruby.rb: upper = Math.sqrt(x.to_f).to_i ruby-spark-1.2.1/lib/spark/rdd.rb: max_sample_size = Integer::MAX - (num_st_dev * Math.sqrt(Integer::MAX)).to_i ruby-treemap-0.0.4/lib/treemap/node.rb: (base_size * Math.sqrt(self.bounds.width * self.bounds.height) / 125).to_i ruby-treemap-fork-0.0.4/lib/treemap/node.rb: (base_size * Math.sqrt(self.bounds.width * self.bounds.height) / 125).to_i rubyvis-0.6.1/examples/5_pv_hierarchies/bubble_charts.rb: .text(lambda {|d| d.class_name[0, Math.sqrt(d.node_value).to_i / 8]}) rusby-0.1.2/lib/rusby/profiler.rb: SQRT_ITERATIONS = Math.sqrt(1000).to_i simstring_pure-1.1.1/lib/simstring_pure.rb: (alpha * Math.sqrt(query_size * y_size)).ceil.to_i sudoku_maker-0.0.2/lib/sudoku_maker/maker.rb: @decision_num = Math.sqrt(num).to_i sudoku_solver-1.0.0/lib/sudoku_solver/solver.rb: @dimension = Math.sqrt(@width).to_i tactical_tic_tac_toe-1.0.0/lib/board.rb: Math.sqrt(marked_spaces.count).to_i taskjuggler-3.6.0/lib/taskjuggler/reports/ChartPlotter.rb: rr = (r / Math.sqrt(2.0)).to_i tdiary-contrib-5.0.3/plugin/category_to_tagcloud.rb: level = ((Math.sqrt(count) - min) * factor).to_i tic_tac_toes-0.1.8/lib/tic_tac_toes/ui/serializer.rb: row_size = Math.sqrt(board_structure.count).to_i ul-wukong-4.1.1/lib/wukong/widget/reducers/bin.rb: self.num_bins = Math.sqrt(total_count).to_i whitestone-1.0.2/etc/extra_tests/output_examples_code.rb: (2..Math.sqrt(n).to_i).map { |i| n % i }.all? { |rem| rem != 0 } wukong-4.0.0/lib/wukong/widget/reducers/bin.rb: self.num_bins = Math.sqrt(total_count).to_i ``` ---------------------------------------- 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-63135 * 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>