Sieve of Zakiya -- Number Theory Sieves (NTS)
I present with this Ruby code a new class of Number Theory Sieves (NTS)
to generate primes numbers upto some number N, which I believe now
represent the fastest way to generate prime numbers. The general number
theory can utilize a class of prime generator functions, which have the
same form, and which can be implemented in the same manner. These NTS
can be inherently implemented in parallel with multiple processors.
The classical brute-force method to generate primes upto some N is the
Sieve of Eratothenes (SoE). In 2002/3 Atkin & Bernstein proposed what has
become known as the Sieve of Atkin (SoA), which was known to be the fastest
method to generate primes. My general number theory based Sieve of Zakiya (SoZ)
represents a multiple times increase in performance over the SoA, while being
easier to understand, code, and extend with newer prime generator functions.
I also used the number theory to created a primality tester, whose code
is short and elegant, which I also believe is the fastest method to test
for primality while being deterministic. Thus, it should replace other
non-deterministic or probabilistic methods, such as Miller-Rabin, et al.
I started this work this first week in May 2008, and present this code
herein on Friday, June 6, 2008.
The development of this work is presented in the paper:
"Ultimate Prime Sieve -- Sieve of Zakiya"
the pdf can be download from here:
http://www.4shared.com/dir/7467736/97bd7b71/sharing.html
#!/usr/local/bin/ruby -w require 'benchmark' class Integer def primes1 # classical brute-force Sieve of Eratosthenes (SoE) # initialze sieve array with all integers starting at i=2 sieve = [nil, nil].concat((2..self).to_a) (2 .. Math.sqrt(self).to_i).each do |i| next unless sieve[i] (i*i).step(self, i) { |j| sieve[j] = nil } end sieve.compact! end def primes2 # classical brute-force Sieve of Eratosthenes (SoE) # initialize sieve array with odd integers for odd indexes, then 2 sieve = []; 3.step(self, 2) { |i| sieve[i] = i }; sieve[2] = 2 3.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[i] (i*i).step(self, i<<1) { |j| sieve[j] = nil } end sieve.compact! end def primes3 # http://everything2.com/index.pl?node_id=1176369 # all prime candidates > 3 are of form 6*k+1 and 6*k+5 # initialize sieve array with only these candidate values n1, n2 = 1, 5; sieve = [] while n2 < self n1 += 6; n2 += 6; sieve[n1] = n1; sieve[n2] = n2 end # now initialize sieve array with first 3 primes, resize array sieve[2]=2; sieve[3]=3; sieve[5]=5; sieve=sieve[0..self] 5.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[i] (i*i).step(self, i<<1) {|j| sieve[j] = nil } end sieve.compact! end def primesP3 # all prime candidates > 3 are of form 6*k+1 and 6*k+5 # initialize sieve array with only these candidate values n1, n2 = 1, 5; sieve = [] while n2 < self n1 += 6; n2 += 6; sieve[n1] = n1; sieve[n2] = n2 end # now initialize sieve array with primes < 6, resize array sieve[2]=2; sieve[3]=3; sieve[5]=5; sieve=sieve[0..self] 5.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[i] # p5= 5*i, k = 6*i, p7 = 7*i p1 = 5*i; k = p1+i; p2 = k+i while p1 <= self sieve[p1] = nil; sieve[p2] = nil; p1 += k; p2 += k end end sieve.compact! end def primesP3a # all prime candidates > 3 are of form 6*k+1 and 6*k+5 # initialize sieve array with only these candidate values # where sieve contains the odd integers representatives # convert integers to array indices/vals by i = (n-3)>>1 = (n>>1)-1 n1, n2 = -1, 1; lndx= (self-1) >>1; sieve = [] while n2 < lndx n1 +=3; n2 += 3; sieve[n1] = n1; sieve[n2] = n2 end #now initialize sieve array with (odd) primes < 6, resize array sieve[0] =0; sieve[1]=1; sieve=sieve[0..lndx-1] 5.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[(i>>1) - 1] # p5= 5*i, k = 6*i, p7 = 7*i # p1 = (5*i-3)>>1; p2 = (7*i-3)>>1; k = (6*i)>>1 i6 = 6*i; p1 = (i6-i-3)>>1; p2 = (i6+i-3)>>1; k = i6>>1 while p1 < lndx sieve[p1] = nil; sieve[p2] = nil; p1 += k; p2 += k end end return [2] if self < 3 [2]+([nil]+sieve).compact!.map {|i| (i<<1) +3 } end def primesP5 # all prime candidates > 5 are of form 30*k+(1,7,11,13,17,19,23,29) # initialize sieve array with only these candidate values n1, n2, n3, n4, n5, n6, n7, n8 = 1, 7, 11, 13, 17, 19, 23, 29; sieve = [] while n8 < self n1 +=30; n2 += 30; n3 += 30; n4 += 30 n5 +=30; n6 += 30; n7 += 30; n8 += 30 sieve[n1] = n1; sieve[n2] = n2; sieve[n3] = n3; sieve[n4] = n4 sieve[n5] = n5; sieve[n6] = n6; sieve[n7] = n7; sieve[n8] = n8 end # now initialize sieve with the primes < 30, resize array sieve[2]=2; sieve[3]=3; sieve[5]=5; sieve[7]=7; sieve[11]=11; sieve[13]=13 sieve[17]=17; sieve[19]=19; sieve[23]=23; sieve[29]=29; sieve=sieve[0..self] 7.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[i] # p1=7*i, p2=11*i, p3=13*i ---- p6=23*i, p7=29*i, p8=31*i, k=30*i n6 = 6*i; p1 = i+n6; p2 = n6+n6-i; p3 = p1+n6; p4 = p2+n6 p5 = p3+n6; p6 = p4+n6; p7 = p6+n6; k = p7+i; p8 = k+i while p1 <= self sieve[p1] = nil; sieve[p2] = nil; sieve[p3] = nil; sieve[p4] = nil sieve[p5] = nil; sieve[p6] = nil; sieve[p7] = nil; sieve[p8] = nil p1 += k; p2 += k; p3 += k; p4 += k; p5 += k; p6 += k; p7 += k; p8 += k end end sieve.compact! end def primesP5a # all prime candidates > 5 are of form 30*k+(1,7,11,13,17,19,23,29) # initialize sieve array with only these candidate values # where sieve contains the odd integers representatives # convert integers to array indices/vals by i = (n-3)>>1 = n>>1 - 1 n1, n2, n3, n4, n5, n6, n7, n8 = -1, 2, 4, 5, 7, 8, 10, 13 lndx= (self-1) >>1; sieve = [] while n8 < lndx n1 +=15; n2 += 15; n3 += 15; n4 += 15 n5 +=15; n6 += 15; n7 += 15; n8 += 15 sieve[n1] = n1; sieve[n2] = n2; sieve[n3] = n3; sieve[n4] = n4 sieve[n5] = n5; sieve[n6] = n6; sieve[n7] = n7; sieve[n8] = n8 end # now initialize sieve with the (odd) primes < 30, resize array sieve[0]=0; sieve[1]=1; sieve[2]=2; sieve[4]=4; sieve[5]=5 sieve[7]=7; sieve[8]=8; sieve[10]=10; sieve[13]=13 sieve = sieve[0..lndx-1] 7.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[(i>>1) - 1] # p1=7*i, p2=11*i, p3=13*i ---- p6=23*i, p7=29*i, p8=31*i, k=30*i # maps to: p1= (7*i-3)>>1, --- p8=(31*i-3)>>1, k=30*i>>1 n2=i<<1; n4=i<<2; n8=i<<3; n16=i<<4; n32=i<<5; n12=n8+n4 p1 = (n8-i-3)>>1; p2 = (n12-i-3)>>1; p3 = (n12+i-3)>>1; p4 = (n16+i-3)>>1 p5 = (n16+n4-i-3)>>1; p6 = (n16+n8-i-3)>>1; p7 = (n32-n4+i-3)>>1 p8 = (n32-i-3)>>1; k = (n32-n2)>>1 while p1 < lndx sieve[p1] = nil; sieve[p2] = nil; sieve[p3] = nil; sieve[p4] = nil sieve[p5] = nil; sieve[p6] = nil; sieve[p7] = nil; sieve[p8] = nil p1 += k; p2 += k; p3 += k; p4 += k; p5 += k; p6 += k; p7 += k; p8 += k end end return [2] if self < 3 [2]+([nil]+sieve).compact!.map {|i| (i<<1) +3 } end def primesP7 # all prime candidates > 7 are of form 210*k+(1,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,121,127,131 # 137,139,143,149,151,157,163,167,169,173,179,181,187,191,193,197,199,209) # initialize sieve array with only these candidate values n1, n2, n3, n4, n5, n6, n7, n8, n9, n10 = 1, 11, 13, 17, 19, 23, 29, 31, 37, 41 n11, n12, n13, n14, n15, n16, n17, n18 = 43, 47, 53, 59, 61, 67, 71, 73 n19, n20, n21, n22, n23, n24, n25, n26 = 79, 83, 89, 97, 101, 103, 107, 109 n27, n28, n29, n30, n31, n32, n33, n34 = 113, 127, 131, 137, 139, 149, 151, 157 n35, n36, n37, n38, n39, n40, n41, n42, n43 = 163, 167, 173, 179, 181, 191, 193, 197, 199 n44, n45, n46, n47, n48 = 121, 143, 169, 187, 209 num = self - self[0]^1; sieve = [] while n48 < num n1 +=210; n2 += 210; n3 += 210; n4 += 210; n5 += 210 n6 +=210; n7 += 210; n8 += 210; n9 += 210; n10 += 210 n11 +=210; n12 += 210; n13 += 210; n14 += 210; n15 += 210 n16 +=210; n17 += 210; n18 += 210; n19 += 210; n20 += 210 n21 +=210; n22 += 210; n23 += 210; n24 += 210; n25 += 210 n26 +=210; n27 += 210; n28 += 210; n29 += 210; n30 += 210 n31 +=210; n32 += 210; n33 += 210; n34 += 210; n35 += 210 n36 +=210; n37 += 210; n38 += 210; n39 += 210; n40 += 210 n41 +=210; n42 += 210; n43 += 210; n44 += 210; n45 += 210 n46 += 210; n47 += 210; n48 += 210 sieve[n1] = n1; sieve[n2] = n2; sieve[n3] = n3; sieve[n4] = n4 sieve[n5] = n5; sieve[n6] = n6; sieve[n7] = n7; sieve[n8] = n8 sieve[n9] = n9; sieve[n10] = n10; sieve[n11] = n11; sieve[n12] = n12 sieve[n13] = n13; sieve[n14] = n14; sieve[n15] = n15; sieve[n16] = n16 sieve[n17] = n17; sieve[n18] = n18; sieve[n19] = n19; sieve[n20] = n20 sieve[n21] = n21; sieve[n22] = n22; sieve[n23] = n23; sieve[n24] = n24 sieve[n25] = n25; sieve[n26] = n26; sieve[n27] = n27; sieve[n28] = n28 sieve[n29] = n29; sieve[n30] = n30; sieve[n31] = n31; sieve[n32] = n32 sieve[n33] = n33; sieve[n34] = n34; sieve[n35] = n35; sieve[n36] = n36 sieve[n37] = n37; sieve[n38] = n38; sieve[n39] = n39; sieve[n40] = n40 sieve[n41] = n41; sieve[n42] = n42; sieve[n43] = n43; sieve[n44] = n44 sieve[n45] = n44; sieve[n46] = n46; sieve[n47] = n47; sieve[n48] = n48 end # now initialize sieve with the (odd) primes < 210, resize array sieve[2]=2; sieve[3]=3; sieve[5]=5; sieve[7]=7; sieve[11]=11; sieve[13]=13 sieve[17]=17; sieve[19]=19; sieve[23]=23; sieve[29]=29; sieve[31]=31 sieve[37]=37; sieve[41]=41; sieve[43]=43; sieve[47]=47; sieve[53]=53 sieve[59]=59; sieve[61]=61; sieve[67]=67; sieve[71]=71; sieve[73]=73 sieve[79]=79; sieve[83]=83; sieve[89]=89; sieve[97]=97; sieve[101]=101 sieve[103]=103; sieve[107]=107; sieve[109]=109; sieve[113]=113; sieve[127]=127 sieve[131]=131; sieve[137]=137; sieve[139]=139; sieve[149]=149; sieve[151]=151 sieve[157]=157; sieve[163]=163; sieve[167]=167; sieve[173]=173; sieve[179]=179 sieve[181]=181; sieve[191]=191; sieve[193]=193; sieve[197]=197; sieve[199]=199 sieve = sieve[0..self] 11.step(Math.sqrt(self).to_i, 2) do |i| next unless sieve[i] # p1=7*i, p2=11*i, p3=13*i ---- p6=23*i, p7=29*i, p8=31*i, k=30*i # maps to: p1= (7*i-3)>>1, --- p8=(31*i)>>1, k=30*i>>1 #n6=6*i; n7=n6+i; n11=n6+n6-i; n13=n6+n7; n17=n11+n6 #n19=n13+n6; n23=n17+n6; n29=n23+n6; n31=n29+(i<<1) #n4=i<<2; n8=i<<3; n16=i<<4; n7=n8-i; n11=n7+n4; n13=n11+(i<<1) #n17=n16+i; n19=n11+n8; n23=n16+n7; n29=n16+n13; n31=(i<<5)-i p1 = (11*i); p2 = (13*i); p3 = (17*i); p4 = (19*i) p5 = (23*i); p6 = (29*i); p7 = (31*i); p8 = (37*i) p9 = (41*i); p10 = (43*i); p11 = (47*i); p12 = (53*i) p13 = (59*i); p14 = (61*i); p15 = (67*i); p16 = (71*i) p17 = (73*i); p18 = (79*i); p19 = (83*i); p20 = (89*i) p21 = (97*i); p22 = (101*i); p23 = (103*i); p24 = (107*i) p25 = (109*i); p26 = (113*i); p27 = (127*i); p28 = (131*i) p29 = (137*i); p30 = (139*i); p31 = (149*i); p32 = (151*i) p33 = (157*i); p34 = (163*i); p35 = (167*i); p36 = (173*i) p37 = (179*i); p38 = (181*i); p39 = (191*i); p40 = (193*i) p41 = (197*i); p42 = (199*i); p43 = (211*i) p44 = (121*i); p45 = 143*i; p46 = 169*i; p47 = 187*i; p48 = 209*i k = (210*i) while p1 <= num sieve[p1] = nil; sieve[p2] = nil; sieve[p3] = nil; sieve[p4] = nil sieve[p5] = nil; sieve[p6] = nil; sieve[p7] = nil; sieve[p8] = nil sieve[p9] = nil; sieve[p10] = nil; sieve[p11] = nil; sieve[p12] = nil sieve[p13] = nil; sieve[p14] = nil; sieve[p15] = nil; sieve[p16] = nil sieve[p17] = nil; sieve[p18] = nil; sieve[p19] = nil; sieve[p20] = nil sieve[p21] = nil; sieve[p22] = nil; sieve[p23] = nil; sieve[p24] = nil sieve[p25] = nil; sieve[p26] = nil; sieve[p27] = nil; sieve[p28] = nil sieve[p29] = nil; sieve[p30] = nil; sieve[p31] = nil; sieve[p32] = nil sieve[p33] = nil; sieve[p34] = nil; sieve[p35] = nil; sieve[p36] = nil sieve[p37] = nil; sieve[p38] = nil; sieve[p39] = nil; sieve[p40] = nil sieve[p41] = nil; sieve[p42] = nil; sieve[p43] = nil; sieve[p44] = nil sieve[p45] = nil; sieve[p46] = nil; sieve[p47] = nil; sieve[p48] = nil p1 += k; p2 += k; p3 += k; p4 += k; p5 +=k; p6 += k; p7 += k p8 += k; p9 += k; p10 +=k; p11 += k; p12 +=k; p13 +=k; p14 +=k p15 +=k; p16 +=k; p17 +=k; p18 += k; p19 +=k; p20 +=k; p21 +=k p22 +=k; p23 +=k; p24 +=k; p25 += k; p26 +=k; p27 +=k; p28 +=k p29 +=k; p30 +=k; p31 +=k; p32 += k; p33 +=k; p34 +=k; p35 +=k p36 +=k; p37 +=k; p38 +=k; p39 += k; p40 +=k; p41 +=k; p42 +=k p43 +=k; p44 +=k; p45 +=k; p46 +=k; p47 +=k; p48 +=k end end sieve.compact! end def primesP7a # all prime candidates > 7 are of form 210*k+(1,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,121,127,131 # 137,139,143,149,151,157,163,167,169173,179,181,187,191,193,197,199,209) # initialize sieve array with only these candidate values # where sieve contains the odd integers representatives # convert integers to array indices/vals by i = (n )>>1 = (n>>1)-1 n1, n2, n3, n4, n5, n6, n7, n8, n9, n10 = -1, 4, 5, 7, 8, 10, 13, 14, 17, 19 n11, n12, n13, n14, n15, n16, n17, n18 = 20, 22, 25, 28, 29, 32, 34, 35 n19, n20, n21, n22, n23, n24, n25, n26 = 38, 40, 43, 47, 49, 50, 52, 53 n27, n28, n29, n30, n31, n32, n33, n34 = 55, 62, 64, 67, 68, 73, 74, 77 n35, n36, n37, n38, n39, n40, n41, n42, n43 = 80, 82, 85, 88, 89, 94, 95, 97, 98 n44, n45, n46, n47, n48 = 59, 70, 83, 92, 103 lndx= (self-1)>>1; sieve = [] while n48 < lndx n1 +=105; n2 += 105; n3 += 105; n4 += 105; n5 += 105 n6 +=105; n7 += 105; n8 += 105; n9 += 105; n10 += 105 n11 +=105; n12 += 105; n13 += 105; n14 += 105; n15 += 105 n16 +=105; n17 += 105; n18 += 105; n19 += 105; n20 += 105 n21 +=105; n22 += 105; n23 += 105; n24 += 105; n25 += 105 n26 +=105; n27 += 105; n28 += 105; n29 += 105; n30 += 105 n31 +=105; n32 += 105; n33 += 105; n34 += 105; n35 += 105 n36 +=105; n37 += 105; n38 += 105; n39 += 105; n40 += 105 n41 +=105; n42 += 105; n43 += 105; n44 += 105; n45 += 105 n46 +=105; n47 += 105; n48 += 105 sieve[n1] = n1; sieve[n2] = n2; sieve[n3] = n3; sieve[n4] = n4 sieve[n5] = n5; sieve[n6] = n6; sieve[n7] = n7; sieve[n8] = n8 sieve[n9] = n9; sieve[n10] = n10; sieve[n11] = n11; sieve[n12] = n12 sieve[n13] = n13; sieve[n14] = n14; sieve[n15] = n15; sieve[n16] = n16 sieve[n17] = n17; sieve[n18]<