Over a million developers have joined DZone.
{{announcement.body}}
{{announcement.title}}

DZone's Guide to

# First N Primes In Ruby

·
Free Resource

Comment (0)

Save
{{ articles[0].views | formatCount}} Views
``````

#!/usr/local/bin/ruby -w

require 'benchmark'

class Integer
def prime?         # cf. http://snippets.dzone.com/posts/show/4636
n = self.abs()
return true if n == 2
return false if n == 1 || n & 1 == 0
d = n-1
d >>= 1 while d & 1 == 0
20.times do                               # 20 = k from above
a = rand(n-2) + 1
t = d
y = ModMath.pow(a,t,n)                  # implemented below
while t != n-1 && y != 1 && y != n-1
y = (y * y) % n
t <<= 1
end
return false if y != n-1 && t & 1 == 0
end
return true
end
end

module ModMath
def ModMath.pow(base, power, mod)
result = 1
while power > 0
result = (result * base) % mod if power & 1 == 1
base = (base * base) % mod
power >>= 1;
end
result
end
end

class Integer

def primes   # cf. http://snippets.dzone.com/posts/show/3734

sieve = []
3.step(self, 2) { |i| sieve[i] = i }
sieve[1] = nil
sieve[2] = 2

3.step(Math.sqrt(self).floor, 2) do |i|
next unless sieve[i]
(i*i).step(self, i) do |j|
sieve[j] = nil
end
end

sieve.compact!

end

def primes2       # cf. http://snippets.dzone.com/posts/show/3734

primes = [nil, nil].concat((2..self).to_a)

(2 .. Math.sqrt(self)).each do |i|
next unless primes[i]
(i*i).step(self, i) do |j|
primes[j] = nil
end
end

primes.compact!

end

end

class Integer

@primes_cache ||= 100.primes
#@primes_cache ||= [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31]
class << self; attr_accessor :primes_cache; end

@nthprime_limit ||= 5_000_000
class << self; attr_reader :nthprime_limit; end
#class << self; attr_accessor :nthprime_limit; end

def sieve_size
num = self.to_i.abs
logn = Math.log(num)
if num < 15985     # cf. http://primes.utm.edu/howmany.shtml
( num * (logn + Math.log(logn) - 1 + 1.8 * Math.log(logn) / logn) ).floor
else
( num * (logn + Math.log(logn) - 0.9427) ).floor
end
end

def nthprime

num = self.to_i.abs

# set a limit; cf. http://primes.utm.edu/nthprime/: The 5,000,000th prime is 86,028,121.
raise "#{num}: number too big for Integer#nthprime" if num > Integer.nthprime_limit

primes_cache_size = Integer.primes_cache.size

if num > primes_cache_size

if num >= 50_000 && num < 100_000 && (num - primes_cache_size) < 5000 then return num.nthprime_add end
if num >= 100_000 && num < 200_000 && (num - primes_cache_size) < 15000 then return num.nthprime_add end
if num >= 200_000 && (num - primes_cache_size) < 35000 then return num.nthprime_add end

logn = Math.log(num)

if num < 15985       # cf. http://primes.utm.edu/howmany.shtml
limit = ( num * (logn + Math.log(logn) - 1 + 1.8 * Math.log(logn) / logn) ).floor
Integer.primes_cache = limit.primes
else
limit = ( num * (logn + Math.log(logn) - 0.9427) ).floor
Integer.primes_cache = limit.primes
end

=begin
elsif num >= 15985 && num <= 1_000_000
limit = ( num * (logn + Math.log(logn) - 0.9427) ).floor
Integer.primes_cache = limit.primes
elsif num > 1_000_000
Integer.primes_cache = 20_000_000.primes
end
=end

else
return Integer.primes_cache.at(num-1)
end

if num <= Integer.primes_cache.size
return Integer.primes_cache.at(num-1)
end

if Integer.primes_cache.size > 2_500_000
else
num.nthprime_add   # faster for a (relatively) small prime cache
end

end

def nthprime_add       #  add primes to Integer.primes_cache up to the nth prime

num = self.to_i.abs

if num <= Integer.primes_cache.size
return Integer.primes_cache.at(num-1)
end

last_prime = Integer.primes_cache.last
last_prime_divmod = last_prime.divmod(6)

if last_prime_divmod.last == 1
i = last_prime_divmod.first
Integer.primes_cache.pop      # avoid a duplicate prime
else
i = last_prime_divmod.first + 1
end

while Integer.primes_cache.size < num

n1 = 6*i+1       # cf. http://betterexplained.com/articles/another-look-at-prime-numbers/ and
n2 = n1+4        # http://everything2.com/index.pl?node_id=1176369
i += 1

[n1, n2].each do |p|

next if p % 5 == 0 || p % 7 == 0

next_num_sqrt = Math.sqrt(p).floor

Integer.primes_cache.each do |d|
if d > next_num_sqrt
Integer.primes_cache << p
#print "\r\e[0K#{Integer.primes_cache.size}"
#\$stdout.flush
break
elsif p % d == 0
break
end
end
end   # each

end  # while

return Integer.primes_cache.at(num-1)

end

def nthprime_add_mr       #  add Miller-Rabin primes to Integer.primes_cache up to the nth prime

num = self.to_i.abs

if num <= Integer.primes_cache.size
return Integer.primes_cache.at(num-1)
end

last_prime = Integer.primes_cache.last
last_prime_divmod = last_prime.divmod(6)

if last_prime_divmod.last == 1
i = last_prime_divmod.first
Integer.primes_cache.pop
else
i = last_prime_divmod.first + 1
end

while Integer.primes_cache.size < num

search_next_prime = true

while search_next_prime

n1 = 6*i+1
n2 = n1+4
i += 1

[n1, n2].each do |p|

next if p % 5 == 0 || p % 7 == 0

if p.prime?
Integer.primes_cache << p
search_next_prime = false
#print "\r\e[0K#{Integer.primes_cache.size}"
#\$stdout.flush
end
end
end

end  #  first while loop

return Integer.primes_cache.at(num-1)

end

def next_primes_in_cache    # next primes in Integer.primes_cache

search_next_primes = self.to_i.abs

last_prime = Integer.primes_cache.last
last_prime_divmod = last_prime.divmod(6)

if last_prime_divmod.last == 1
i = last_prime_divmod.first
Integer.primes_cache.pop
else
i = last_prime_divmod.first + 1
end

while search_next_primes > 0

n1 = 6*i+1
n2 = n1+4
i += 1

[n1, n2].each do |p|

next if p % 5 == 0 || p % 7 == 0
next_num_sqrt = Math.sqrt(p).floor

Integer.primes_cache.each do |d|
if d > next_num_sqrt
Integer.primes_cache << p
search_next_primes -= 1
#print "\r\e[0K#{Integer.primes_cache.size}"
#\$stdout.flush
break
elsif p % d == 0
break
end
end
end
end   # while

Integer.primes_cache.last(self.to_i.abs)

end

def next_mr_prime     # next Miller-Rabin prime

last_num = self.to_i.abs
last_num_divmod = last_num.divmod(6)

if last_num_divmod.last == 1
i = last_num_divmod.first
else
i = last_num_divmod.first + 1
end

next_prime = nil
search_next_prime = true

while search_next_prime

n1 = 6*i+1
n2 = n1+4
i += 1

[n1, n2].each do |p|

next if p % 5 == 0 || p % 7 == 0

if p > last_num && p.prime?
next_prime = p
search_next_prime = false
break
end
end
end

next_prime

end

def next_mr_primes(n)     # next Miller-Rabin primes

last_num = self.to_i.abs
last_num_divmod = last_num.divmod(6)

if last_num_divmod.last == 1
i = last_num_divmod.first
else
i = last_num_divmod.first + 1
end

next_primes = []
search_next_primes = n.to_i.abs

while search_next_primes > 0

n1 = 6*i+1
n2 = n1+4
i += 1

[n1, n2].each do |p|

next if p % 5 == 0 || p % 7 == 0

if p > last_num && p.prime?
next_primes << p
search_next_primes -= 1
end
end
end

next_primes

end

end

num1 = 10_000
num1 = 1_000
num1 = 5_000

num2 = 210_349
num2 = 100_125
num2 = 55_000

ret1 = nil
ret2 = nil
ret3 = nil
ret4 = nil

Benchmark.bm(16) do |t|

t.report("first #{num1} primes: ") do
end

Integer.primes_cache.clear
Integer.primes_cache.concat(100.primes)

t.report("first #{num1} primes: ") do
end

Integer.primes_cache.clear
Integer.primes_cache.concat(100.primes)

t.report("first #{num1} primes: ") do
ret3 = num1.nthprime
end

Integer.primes_cache.clear
Integer.primes_cache.concat(100.primes)

t.report("first #{num2} primes: ") do
ret4 = num2.nthprime
end

end

puts
puts "the #{num1}th prime number: #{ret1}"
puts "the #{num1}th prime number: #{ret2}"
puts "the #{num1}th prime number: #{ret3}"
puts "the #{num2}th prime number: #{ret4}"
puts
puts "the #{num1}th prime: #{Integer.primes_cache.at(num1-1)}"
puts "the #{num2}th prime: #{Integer.primes_cache.at(num2-1)}"
puts
p Integer.primes_cache.first(10)
p Integer.primes_cache.last(10)
p Integer.primes_cache.size
puts

p 15.next_primes_in_cache

puts 594_213.next_mr_prime

p 149_137.next_mr_primes(5)

puts
p 1_000.sieve_size
p 1_000_000.sieve_size
p 2_500_000.sieve_size
p 5_000_000.sieve_size
p 10_000_000.sieve_size
p 100_000_000.sieve_size

=begin

# check with primegen.c, http://cr.yp.to/primegen.html
primes 2 48611 | nl | tail -n 1
primes 2 104729 | nl | tail -n 1
primes 2 679277 | nl | tail -n 1
primes 2 1301423 | nl | tail -n 1
primes 2 2903603 | nl | tail -n 1
primes 1303129 1303283 | nl
primes 594213 \$((594213 + 500)) | head -n 1
primes 149137 \$((149137 + 1000)) | head -n 5

# further prime tools from http://libtom.org
- LibTomMath
bn_mp_prime_is_prime.c
bn_mp_prime_miller_rabin.c
- TomsFastMath
fp_isprime.c
fp_prime_miller_rabin.c

=end
``````
Topics:

Comment (0)

Save
{{ articles[0].views | formatCount}} Views

Opinions expressed by DZone contributors are their own.