Skip to content

Instantly share code, notes, and snippets.

@marcandre
Last active December 5, 2020 10:57
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save marcandre/b263bdae488e76dabdda84daf73733b9 to your computer and use it in GitHub Desktop.
Save marcandre/b263bdae488e76dabdda84daf73733b9 to your computer and use it in GitHub Desktop.
Miller Rabin prime test
<<~RESULT
Testing 10..100
current: 389.3 i/s
faster_prime: 269.8 i/s - 1.44x (± 0.00) slower
miller_rabin: 226.5 i/s - 1.72x (± 0.00) slower
Testing 100..1000
current: 323.4 i/s
faster_prime: 254.4 i/s - 1.27x (± 0.00) slower
miller_rabin: 213.0 i/s - 1.52x (± 0.00) slower
Testing 1000..10000
current: 252.3 i/s
faster_prime: 250.5 i/s - same-ish: difference falls within error
miller_rabin: 178.1 i/s - 1.42x (± 0.00) slower
Testing 10000..100000
faster_prime: 242.2 i/s
current: 176.0 i/s - 1.38x (± 0.00) slower
miller_rabin: 168.7 i/s - 1.44x (± 0.00) slower
Testing 100000..1000000
miller_rabin: 153.6 i/s
current: 97.6 i/s - 1.57x (± 0.00) slower
faster_prime: 63.7 i/s - 2.41x (± 0.00) slower
Testing 1000000..10000000
miller_rabin: 146.7 i/s
faster_prime: 58.4 i/s - 2.51x (± 0.00) slower
current: 45.0 i/s - 3.26x (± 0.00) slower
Testing 10000000..100000000
miller_rabin: 131.4 i/s
faster_prime: 56.5 i/s - 2.33x (± 0.00) slower
current: 20.1 i/s - 6.53x (± 0.00) slower
Testing 100000000..1000000000
miller_rabin: 125.0 i/s
faster_prime: 47.7 i/s - 2.62x (± 0.00) slower
current: 7.4 i/s - 16.88x (± 0.00) slower
Testing 1000000000..10000000000
miller_rabin: 54.1 i/s
faster_prime: 26.0 i/s - 2.08x (± 0.00) slower
current: 2.6 i/s - 20.47x (± 0.00) slower
Testing 10000000000..100000000000
miller_rabin: 24.0 i/s
faster_prime: 17.1 i/s - 1.41x (± 0.00) slower
current: 1.0 i/s - 23.86x (± 0.00) slower
Testing 100000000000..1000000000000
miller_rabin: 20.2 i/s
faster_prime: 15.1 i/s - 1.34x (± 0.00) slower
current: 0.4 i/s - 56.60x (± 0.00) slower
Testing 1000000000000..10000000000000
miller_rabin: 14.5 i/s
faster_prime: 10.8 i/s - 1.35x (± 0.00) slower
current: 0.1 i/s - 132.11x (± 0.00) slower
Testing 10000000000000..100000000000000
miller_rabin: 11.9 i/s
faster_prime: 9.6 i/s - same-ish: difference falls within error
Testing 100000000000000..1000000000000000
miller_rabin: 10.6 i/s
faster_prime: 0.4 i/s - 25.93x (± 0.00) slower
Testing 1000000000000000..10000000000000000
miller_rabin: 8.9 i/s
faster_prime: 0.3 i/s - 26.44x (± 0.00) slower
Testing 10000000000000000..100000000000000000
miller_rabin: 8.3 i/s
faster_prime: 0.3 i/s - 25.18x (± 0.00) slower
Testing 100000000000000000..1000000000000000000
miller_rabin: 7.6 i/s
faster_prime: 0.4 i/s - 20.22x (± 0.00) slower
Testing 1000000000000000000..10000000000000000000
miller_rabin: 13.7 i/s
faster_prime: 0.3 i/s - 41.51x (± 0.00) slower
Testing 10000000000000000000..100000000000000000000
miller_rabin: 27.5 i/s
faster_prime: 0.3 i/s - 109.57x (± 0.00) slower
Testing 100000000000000000000..1000000000000000000000
miller_rabin: 28.0 i/s
faster_prime: 0.1 i/s - 344.33x (± 0.00) slower
Testing 1000000000000000000000..10000000000000000000000
miller_rabin: 23.9 i/s
faster_prime: 0.3 i/s - 85.94x (± 0.00) slower
Testing 10000000000000000000000..100000000000000000000000
miller_rabin: 27.7 i/s
faster_prime: 0.3 i/s - 98.29x (± 0.00) slower
Testing 100000000000000000000000..1000000000000000000000000
miller_rabin: 26.6 i/s
faster_prime: 0.2 i/s - 144.06x (± 0.00) slower
RESULT
class Integer
def prime2?
return self >= 2 if self <= 3
return false if even?
r = 0
d = self >> 1
while d.even?
d >>= 1
r += 1
end
self_minus_1 = self-1
miller_rabin_bases.each do |a|
x = a.pow(d, self)
next if x == 1 || x == self_minus_1 || a == self
return false if r.times do
x = x.pow(2, self)
break if x == self_minus_1
end
end
true
end
MILLER_RABIN_BASES = [
[2],
[2,3],
[31,73],
[2,3,5],
[2,3,5,7],
[2,7,61],
[2,13,23,1662803],
[2,3,5,7,11],
[2,3,5,7,11,13],
[2,3,5,7,11,13,17],
[2,3,5,7,11,13,17,19,23],
[2,3,5,7,11,13,17,19,23,29,31,37],
[2,3,5,7,11,13,17,19,23,29,31,37,41],
].map!(&:freeze).freeze
private def miller_rabin_bases
# Miller-Rabin's complexity is O(k log^3n).
# So we can reduce the complexity by reducing the number of bases tested.
# Using values from https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test
i = case
when self < 2_047 then 0
when self < 1_373_653 then 1
when self < 9_080_191 then 2
when self < 25_326_001 then 3
when self < 3_215_031_751 then 4
when self < 4_759_123_141 then 5
when self < 1_122_004_669_633 then 6
when self < 2_152_302_898_747 then 7
when self < 3_474_749_660_383 then 8
when self < 341_550_071_728_321 then 9
when self < 3_825_123_056_546_413_051 then 10
when self < 318_665_857_834_031_151_167_461 then 11
when self < 3_317_044_064_679_887_385_961_981 then 12
else raise
end
MILLER_RABIN_BASES[i]
end
end
require 'bundler/inline'
gemfile do
source 'https://rubygems.org'
gem 'benchmark-ips'
gem 'faster_prime', git: 'https://github.com/mame/faster_prime.git', require: false
end
require 'prime'
require "faster_prime/base"
srand(42)
2.upto(24) do |exp|
r = 10 ** (exp - 1) .. 10 ** exp
samples = 10_000.times.map { rand(r) }
puts "Testing #{r}"
Benchmark.ips do |x|
x.config(:time => 1, :warmup => 0.1)
x.report('current') { samples.each { |i| i.prime? } } if exp < 14
x.report('miller_rabin') { samples.each { |i| i.prime2? } }
x.report('faster_prime') { samples.each { |i| FasterPrime.prime?(i) } }
x.compare!
end
end
@mame
Copy link

mame commented Dec 5, 2020

Looks like samples is not defined?

@marcandre
Copy link
Author

Oops, sorry, I deleted a line by mistake. Fixed.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment