Never been to DZone Snippets before?

Snippets is a public source code repository. Easily build up your personal collection of code snippets, categorize them with tags / keywords, and share them with the world

« Newer Snippets
Older Snippets »
Showing 1-3 of 3 total  RSS 

Sieve of Zakiya -- Number Theory Sieves (NTS)

// description of your code here
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]<