The prime sieve of Atkin (http://en.wikipedia.org/wiki/Sieve_of_Atkin)

is, if done right, faster than the sieve of Eratosthenes. I’ve

included a reimplementation of the mathn Prime class using this sieve,

but it’s slower than the one in ruby 1.9.

If anyone wants to try their hand at improving this (or starting anew

if this is a bad start), feel free. Perhaps we can beat the one in

ruby 1.9 so we can replace it before 1.9 (1.10 ?) becomes the stable

version.

class Prime_Atkin

include Enumerable

# @@next_to_check is a multiple of 12.

@@next_to_check = 240

# @@primes should contain all primes in 1…@@next_to_check

@@primes = [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,

73,79,83,89,97,101,103,107,109,113,127,131,137,139,149,151,157,

163,167,173,179,181,191,193,197,199,211,223,227,229,233,239]

# @@primes_squared = @@primes.map { |prime| prime**2 }[2…-1]

@@primes_squared = [25,49,121,169,289,361,529,841,961,1369,1681,

1849,2209,2809,3481,3721,4489,5041,5329,6241,6889,7921,9409,

10201,10609,11449,11881,12769,16129,17161,18769,19321,22201,

22801,24649,26569,27889,29929,32041,32761,36481,37249,38809,

39601,44521,49729,51529,52441,54289,57121]

# SieveWidth is a multiple of 12.

# Ensure that @@next_to_check + SieveWidth <= @@primes_squared.last

SieveWidth = 56880

@@new_primes = Array.new(SieveWidth, false)

def initialize

@index = -1

end

def succ

@index += 1

while @index >= @@primes.length

range_end = @@next_to_check + SieveWidth

high_x = Math.sqrt(range_end).floor

```
1.upto(high_x) do |x|
x_sq = x**2
low_y = @@next_to_check - 4*x_sq
if low_y < 1
low_y = 1
else
low_y = Math.sqrt(low_y).ceil
end
high_y = [Math.sqrt(3*x_sq - @@next_to_check),
```

Math.sqrt(range_end - 3*x_sq)].max.floor

```
high_y.downto(low_y) do |y|
y_sq = y**2
n = 4*x_sq + y_sq - @@next_to_check
@@new_primes[n] = (not @@new_primes.at(n)) if (n >= 0 and n <
```

SieveWidth and (n % 12 == 1 or n % 12 == 5))

n -= x_sq

@@new_primes[n] = (not @@new_primes.at(n)) if (n >= 0 and n <

SieveWidth and n % 12 == 7)

n -= 2*y_sq

@@new_primes[n] = (not @@new_primes.at(n)) if (n >= 0 and n <

SieveWidth and x > y and n % 12 == 11)

end

end

```
@@primes_squared.each do |prime_squared|
multiple_index = prime_squared *
```

(@@next_to_check/prime_squared).ceil - @@next_to_check

while multiple_index < SieveWidth

@@new_primes[multiple_index] = false

multiple_index += prime_squared

end

end

```
@@new_primes.each_with_index do |prime_test, index|
if prime_test
prime = @@next_to_check + index
@@primes << prime
@@primes_squared << prime**2
end
end
@@next_to_check = range_end
@@new_primes.fill false
end
@@primes.at @index
```

end

alias next succ

def each

loop do

yield succ

end

end

end