Prime Sieves

Math 242 Modern Computational Math

Sieve of Eratosthenes

Implementation 1

Here is a simple implementation of the Sieve of Eratosthenes.

sieveErat1[n_Integer]:=Module[
{nums=Range[2,n], (* we will remove composite numbers from this list *)
i=1, (* index of d, a number whose multiples we will remove *)
j            (* index of a number being tested to see if it is a multiple of d *)
},

While[nums[[i]]≤Sqrt[n], (* we will remove multiples of nums[[i]] *)
j=i+1;                                  (* start at the number that follows nums[[i]] *)
While[j≤Length[nums],  (* loop until j reaches the end of the list *)
If[Divisible[nums[[j]],nums[[i]]],  (* if nums[[i]] divides nums[[j]] ... *)
nums=Delete[nums,j],                     (* ... then delete nums[[j]] ... *)
j++]                                                              (* ... otherwise, increment j *)
];

(* TESTING: print the list *)
(*Print["removed multiples of ",nums[[i]],": ",nums];*)

i++   (* increment i *)
];
nums   (* return the list, which now only contains primes *)
]

Try it out:

sieveErat1[100]

This works, but it’s awkward to have to keep track of two indexes, i and j. Also, the divisibility testing should be unnecessary: when doing the Sieve of Eratosthenes by hand, we never test for divisiblity; we just cross off every dth number in the list.

sieveErat1[1000]

How fast is this?

Timing[sieveErat1[1000]][[1]]

Timing[sieveErat1[100000]][[1]]

Implementation 2

Here is another implementation idea: Start with a list of positive integers from 1 to n. Instead of removing non-primes, set them to zero. When finished, select all nonzero elements. Write this implementation, and see how its runtime compares to your other implementations.

sieveErat2[n_Integer]:=Module[{nums=Range[n],m=Floor[Sqrt[n]]},
Do[
(* if nums[[k]] is not zero, then set all multiples to zero *)
If[ nums[[k]]≠0,nums[[k^2;;n;;k]]=0 ]
,{k,2,m}];

(* select all indexes greater than 1 whose values are nonzero *)
Select[ Range[2,n], nums[[#]]≠0 &]
]

Try it out:

sieveErat2[100]

Timing[sieveErat2[100000]][[1]]

Timing[sieveErat2[1000000]][[1]]

Timing[sieveErat2[10000000]][[1]]

Sieve of Sundaram

Start with the list of integers 1,2,...,n. For every pair of integers (i,j) satisfying 1≤i≤j and i+j+2ij≤n, remove the number i+j+2ij from the list. Take all the numbers still in the list: double them and add 1 to each. This results in a list of all the odd primes below 2n+2.

Implement the Sieve of Sundaram.

sieveSund:
	Argument: an integer n
	Return: a list of primes less than or equal to 2n+2

For example:  sieveSund[20] ⟶ {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41}

Runtime

How does the runtime of the Sieve of Sundaram compare with the runtime of the Sieve of Eratosthenes?

Bonus: Fermat’s Conjecture

Fermat conjectured that all numbers of the form (2^(2^n))+1 are prime. For example, (2^(2^1))+1=4+1=5 is prime, and so is (2^(2^2))+1=(2^4)+1=17; these are called Fermat primes. Investigate Fermat’s conjecture for other values of n. Do you think the conjecture is true? How many Fermat primes can you find?