Eratosthenes sieve in Java: a puzzle and some optimizations

Question

I quickly implemented the SoE algo implementation in Java (code at the end). The output on my dual core AMD processor:

  Allocation: 31
 Meat: 10140
 Listing: 10171
 Preparing end: 10187

• The Meat section consumes the maximum amount of time as expected.

• One of my observations was that using Math.pow(variable, 2) slower than (variable * variable) . I think that besides the jump function, there may be other overhead.

• Does Math.pow (x, 2) have optimization for degrees 2, 3, etc.? I ask because there are some Java custom libraries where there are much faster multiplication algorithms than native Java.

Here are my questions:

• What arithmetic optimizations can you offer in the "Meat" section? Is it possible to avoid the module operator at all?

• The function does not work when start == end. If I make a sieve (4, 4), the returned array has a length of 1: [4]. What am I doing wrong? It should return [] (basically a new int (0)).

• What are some of the famous fast number / math java libraries?

Thank you for reading. Finally, here is the code I wrote: "Not GangOfFour / TopCoder, but not too pathetic (hopefully formatting the code in SO is ... weird?):

 import java.util.ArrayList; import java.util.Iterator; import java.util.List; public class Sieve { public static void main(String[] args) { /* small test */ int[] primes = sieve(1, 1000000); } /** * returns an array of prime integers * in the given range * * @param start range start * @param end range end * @return */ private static int[] sieve(int start, int end) { long startTime = System.currentTimeMillis(); /* some basic range checks */ if(end < start || start < 1 || end < 1) { throw new ArithmeticException("Messed up input"); } /* generate ints within range */ int[] naturals = new int[end-start+1]; for (int j = 0; j < end - start + 1; j++) { naturals[j] = start + j; } System.out.println("Allocation: \t" + (System.currentTimeMillis() - startTime)); /* init running prime to start, and increment until * running prime squared is greater than the end */ for (int runningPrime = (start == 1 ? 2: start); end > runningPrime*runningPrime; runningPrime++) { for (int i = runningPrime; i < naturals.length; i++) { if(-1 != naturals[i]) { if(naturals[i] % runningPrime == 0) { naturals[i] = -1; } } } } System.out.println("Meat: \t\t" + (System.currentTimeMillis() - startTime)); if(naturals[0] == 1) { naturals[0] = -1; } /* list primes */ List list = new ArrayList(); for (int i = 0; i < naturals.length; i++) { if(-1 != naturals[i]) list.add(naturals[i]); } System.out.println("Listing: \t" + (System.currentTimeMillis() - startTime)); /* create the return int array */ int[] primes = new int[list.size()]; int k = 0; for (Iterator iterator = list.iterator(); iterator.hasNext();) { primes[k++] = ((Integer) iterator.next()).intValue(); } System.out.println("Preparing end: \t" + (System.currentTimeMillis() - startTime)); return primes; } } 

Thanks for the feedback. This is the fixed version below (until someone can break it down again :)

 import java.util.ArrayList; import java.util.Iterator; import java.util.List; public class Sieve { public static void main(String[] args) { /* small test */ int[] primes = sieve(2, 5); System.out.println("Number of primes: " + primes.length); for (int i : primes) { System.out.println(i); } } /** * returns an array of prime integers * in the given range * * @param start range start * @param end range end * @return */ private static int[] sieve(int start, int end) { long startTime = System.currentTimeMillis(); /* some basic range checks */ if(end < start || start < 1 || end < 1) { throw new ArithmeticException("Messed up input"); } /* generate ints within range */ int[] naturals = new int[(int)Math.floor((end-start+1) / 2) + 1]; int allocator = 0; for (int j = 0; j < end - start + 1; j++) { if(!((start + j) % 2 == 0)) { naturals[allocator++] = start + j; } } System.out.println("Allocation: \t" + (System.currentTimeMillis() - startTime)); /* init running prime to 2, and increment until * running prime squared is greater than the end */ for (int runningPrime = 2; end >= runningPrime*runningPrime; runningPrime++) { for (int i = 0; i < naturals.length; i++) { if(-1 != naturals[i]) { if(naturals[i] != runningPrime && naturals[i] % runningPrime == 0) { naturals[i] = -1; } } } } System.out.println("Meat: \t\t" + (System.currentTimeMillis() - startTime)); if(naturals[0] == 1) { naturals[0] = -1; } /* list primes */ List list = new ArrayList(); for (int i = 0; i < naturals.length; i++) { if(-1 != naturals[i]) list.add(naturals[i]); } System.out.println("Listing: \t" + (System.currentTimeMillis() - startTime)); /* create the return int array */ int size = list.size(); int k = 0; /* tricky tricky :) */ if(start <= 2) { size += 1; k = 1; } int[] primes = new int[size]; if(start <= 2) { primes[0] = 2; } for (Iterator iterator = list.iterator(); iterator.hasNext();) { primes[k++] = ((Integer) iterator.next()).intValue(); } System.out.println("Preparing end: \t" + (System.currentTimeMillis() - startTime)); return primes; } } 

Answer 1

You can avoid modulo by rewriting the inner loop:

  for (int i = runningPrime; i < naturals.length; i++) { if(-1 != naturals[i]) { if(naturals[i] % runningPrime == 0) { naturals[i] = -1; } } } 

as

  for (int i = runningPrime; i < naturals.length; i+=runningPrime) { naturals[i] = -1; } 

I'm also slightly worried that turning on the start parameter complicates the situation (given the case of sieve(4, 10) ).

Answer 2

Assuming that I did not miss something, I will write:

  for(int runningPrime = (start == 1 ? 2: start); end > runningPrime*runningPrime; runningPrime++) 

as

 int limit = Math.sqrt(end); for(int runningPrime = (start == 1 ? 2: start); runningPrime < limit; runningPrime++) 

to prevent unnecessary multiplication by each iteration. Also, I would just fill the array with odd numbers, effectively halving its length.

Answer 3

Your decision is not a sieve of Eratosthenes. This is obvious because you are using the modulo operator in your code; the correct sieve of Eratosthenes uses only the addition in the internal circuit, and not division or modulo. Here is a simple version of the Eratosthenes sieve that imports the BitSet and LinkedList from java.util and returns a LinkedList of primes less than n:

 public static LinkedList sieve(int n) { BitSet b = new BitSet(n); LinkedList ps = new LinkedList(); b.set(0,n); for (int p=2; p<n; p++) { if (b.get(p)) { ps.add(p); for (int i=p+p; i<n; i+=p) { b.clear(i); } } } return ps; } 

The basic idea is to create a sieve ( BitSet b) with each element originally set to Prime (represented as a set bit), iterating through a sieve that searches and reports each next prime sheet, and when one is hit, all of it multiple of the sieve, marking it Composite (represented as a cleared bit). Multipliers are detected by adding, not dividing, and the inner loop consists only of addition, the operation of clearing bits, comparing to find the end of the sieve and go to the beginning of the loop, so it is very fast.

Optimization options are available that make the Eratosthenes sieve run even faster, but that should be enough to get you started. When you're ready for more, I modestly recommend this essay on my blog.

If you need prime numbers in a range that doesn't start from zero, you need a segmented strainer from Eratosthenes. I discussed the Eratosthenes segmented sieve earlier on Stack Overflow, and also discuss it on my blog.

Answer 4

Filling only with coefficients (except 2), increasing with the runningPrime indicator and losing the divisibility checks that have already been proposed are probably the most important optimizations.

Java Math.pow for paired! It has no optimizations for squaring, basically bc, it immediately returns 2 as double.

Answer 5

In my opinion, before you start optimization, you have to fix two serious errors.

I compiled your code as a Java program and then tried to compute

 sieve(1, 9) 

and

 sieve(4,10); 

The first case worked correctly, except that 9 was considered simple. The square root of 9 is a prime number, but your loop condition stops the screening just before you get there.

In the second case, the estimated primes are 4, 5, 6, 7, 8, 9, 10. This is because you missed the screening of any of the primes below the beginning of the range. This, I'm afraid, is too optimistic :-)