# A fast approach to prime number sieving (non-threaded array)

I am posting this question for input on comparisons between other sieving methods and opinions on what I have discovered. I believe this might be the fastest (singular, non-threaded) array implementation of prime sieving to date, correct me if I am wrong.

My rudimentary benchmark tests stand as follows:

1. My optimized code: 4 seconds N = 1,000,000,000

2. My original code: 7 seconds N = 1,000,000,000

3. Standard Sieve-of-Eratosthenes: 16 seconds N = 1,000,000,000

I have seen PrimeSieve and as I understand it utilizes better means to compute faster, but as it seems - it must dump the data somewhere and retrieve it to utilize it again, because it cycles the computations through the cache only.

My original thread with same pseudo-code is here, however, my original Java code did not match the pseudo-code I developed (but essentially performed similar functions at a lesser speed and easier complexity to code).

Pseudo-code as follows:

1. Initialize list with initial values $(2, 3)$.
2. Use the formula $6k ± 1$ where $k = 1$ and must increment by $1$. Store values to list.
3. Compute step two up to a product $n$, where $n > 7$ if $k = 1$.
4. Starting with $x = 5$, remove all products of $x \times y$ (up to $n$), where $y$ starts at $x$ and
increments to the next number in the list as $x$ remains the same.
5. Repeat step four by assigning $x$ to the next number in the list. Stop when $x \times y > n$.
6. Final step: We must remove all factors of the square of every prime number $(25, 49)$

//-----------------------------------------------------------------
//  Class:   SixesSieve.java
//  Author:  Alex Lieberman
//  Purpose: Finds all prime numbers up to a specified
//           maxNumber (exclusive) and marks them as true in the boolean array.
//-----------------------------------------------------------------

public class SixesSieve {
public static void main (String[] args) {

int maxNumber = 1000000000;
boolean[] isPrime = new boolean [maxNumber];
isPrime = true;
isPrime = true;

for (int i = 5, j = 7; i < maxNumber; i+=6, j+=6) {
isPrime[i] = true;
if (j < maxNumber)
isPrime[j] = true;
}

for (int a = 5; a*a <= maxNumber; a+=2) {
if (isPrime[a]) {

for (int b = a, z = a*a; b*a <= maxNumber; ) {
if(b*a == z && z < maxNumber/2) {
while (true){
if (z <= maxNumber) {
isPrime[z] = false;
z = z+(a*a*2);
}
else
{b+=2; break;}
}
}
else {
isPrime[b*a] = false;
b+=2;
while (!isPrime[b])
{b+=2;}
}
}
}
}

for (int i = 0; i < maxNumber; i++) {
if (isPrime[i])
System.out.print(i + " ");
}
}
}

• What exactly is your question? Are you looking for a code review of the implementation? Or a peer review of your math? – Emily L. Jun 19 '14 at 20:10
• Peer review of the math more than anything and interest in comparisons between other implementations. – Alex Lieberman Jun 19 '14 at 20:18
• For benchmarking you should use some test harness (why. However, for something taking that long it might be OK. I just hope you didn't include the printing (as it mostly depends on the console and nothing else). – maaartinus Jun 19 '14 at 20:46
• Didn't use the printing. Figured that out recently too. – Alex Lieberman Jun 19 '14 at 21:00
• I took the liberty of reverting the question. Please do not edit the code in place when you receive answers, it invalidates these answers and is confusing ("why does that guy say that the variable naming is bad, I think these are good?!"). You can add additional iteration, though I'd advice against it as it always adds the following of a "never ending" question. The most appropriate way would be to revise the code with all suggestions you've got, and come back here with a follow-up question (that is if a significant portion of the code changed and you did not only rename variables). – Bobby Jun 20 '14 at 19:28

Just some minor comments to the code style. Let's start with your fileheader...in these times and days it should be unnecessary. First it doesn't make sense to have the file/classname in the comment, it's only something that will get out of date because refactoring can't reach it. Second, the other parts should be in the JavaDoc of the class.

Variable names. I know it is hard to find good names for such...mathematical stuff, but try to find them. 'a', 'b', 'z' are meaningless. It is hard to read the code if variables have meaningless names.

The code formatting is partly a bloody mess. If you use an IDE, make use of these formatting features. It is so much easier to read code if it is well formatted.

First, don't omit braces, even if it is just a one line 'if' or 'for' block. Second, be consistent about your brace style, and be consequent.

System.out.print(i + " ");


I still don't understand why the Java compiler does allow this. This is not a "cast", this is a hack in the best of cases.

System.out.print(Integer.toString(i) + " ");

• I think the i + " " is there to put a space between each number. – David Harkness Jun 19 '14 at 20:26
• @DavidHarkness: Aye, good point. – Bobby Jun 19 '14 at 20:29
• @Bobby It's a concatenation enforcing an implicit cast, see here and it's damn common. – maaartinus Jun 19 '14 at 20:43
• @maaartinus: Yes, I still don't like it. I have seen too many horrors because of implicit casting. – Bobby Jun 20 '14 at 19:23

Don't write a big main, this makes the code unusable for anything but the current demo. Whenever main gets more than a few lines, extract a method. Btw., Java loves short methods and optimizes them way better than long ones.

Move constants out of methods. Actually, int maxNumber = 1000; should be an argument. Then a comment like

// maxNumber must be an even number


is obviously not enough. You need a Javadoc. And then don't be afraid of

if (maxNumber % 2 != 0) {
throw new IllegalArgumentException("maxNumber must be an even");
}


as trying to process wrong input is waste of time at best.

Try to separate the output from the computation, even for small programs. This is sometimes more work, but you get a reusable result ad it helps you write cleaner code.

If you're interested in optimization, be assured that it's hard. And also fun. In Java it's especially complicated and usually you need some benchmarking framework (they're free).

The code is multiplying possibly big numbers, but I can't see any guarantee that they don't overflow. When I increased the limit to 1.3e9, problems came. That's not that bad as it's close enough to the maximum possible input (so adding a check would do).

The code is badly missing a test. The idea is good (though I still don't understand it fully), but it may produce wrong results. Testing each prime below Integer.MAX_VALUE via BigInteger.isProbablePrime should be doable.

Some code pieces should be improved, e.g.

secondPrime += 2;
while (!isPrime[secondPrime]) {
secondPrime += 2;
}


is a do-while loop. This piece

while (true) {
if (square <= maxNumber) {
isPrime[square] = false;
square += (prime * prime * 2);
}
else {
secondPrime += 2;
break;
}
}


is nothing but

while (square <= maxNumber) {
isPrime[square] = false;
square += (prime * prime * 2);
}
secondPrime += 2;


After a few similar changes and some debugging, I finally understood what it's all about.

It's pretty fast and clever. The speed comes from testing only

• products p*q with p being a prime, isPrime[q] being true1, and p<q
• multiples of p**2 for a prime p

Finding this out from the code alone is a like the twelve labours of Hercules. No offense, this happens quite often with optimized code. That's what CR is for.

1 Because of q being a either prime or not yet detected composite number.

• I ran a test of BigInteger.isProbablePrime to print out if false on every index of the boolean array that held true up to N=1,000,000,000. It hung for quite some time (about a minute or two), but didn't print a thing. – Alex Lieberman Jun 20 '14 at 1:03
• @AlexLieberman There could an error in the opposite direction, too. The full test would take quite long... in the meantime I've verified it against some old code of mine (which is surely correct, as I used it a lot in project Euler). So it's OK, but when I say test, I mean JUnit or whatever running without someone watching the output (it must pass or throw). – maaartinus Jun 20 '14 at 1:13
• True. I appreciate your comments and find them useful :) – Alex Lieberman Jun 20 '14 at 1:24
• @AlexLieberman Glad to hear it... now I was being rather rude. Did you invent the algorithm? – maaartinus Jun 20 '14 at 1:30
• Yes I did. I didn't do a single bit of research actually until it was developed for the most part. I just grabbed a list of numbers and hacked away in Excel for many VERY late nights. I've found some similarities and minor mentions in very few places to some of the steps of the algorithm, but nothing pieced together like this. Thanks for your coding improvement suggestions. – Alex Lieberman Jun 20 '14 at 1:44

The revisions you've already made (variable names, and especially consistent indentation of {, }, and the things between them, helps a lot.

It's still better to organize the code in several functions instead of one main() function. Among other things, this helps reduce the required indentation of some of the code. But mainly it makes functionality a lot easier to see. At a minimum, the code that does the sieving and the code that prints the results should each be in independent functions.

Consider this code fragment:

                    secondPrime += 2;
while (!isPrime[secondPrime]) {
secondPrime += 2;
}


This appears to serve the function of changing secondPrime to the next larger prime. You could put something like this in the body of a function which you would then call like this:

                    secondPrime = getPrimeAfter(secondPrime);


In fact, maybe you could use this to write

            for (int secondPrime = prime, square = prime * prime;
secondPrime * prime <= maxNumber;
secondPrime = getPrimeAfter(secondPrime)) {


which would make the for statement a bit better organized. (I am always uncomfortable when the "advance to the next iteration value" code is not between the second ; and the ) of the for loop control structure. The only question about this is that in the other main branch of logic within that same for loop, you wrote

                            secondPrime += 2;


... that is, in this case you don't increase secondPrime to the next prime, but rather just to the next odd number. The fact that you chose secondPrime as a name of this variable, however, suggests that you didn't actually mean to stop at the next odd number when that number isn't prime, so maybe you do really want to increase secondPrime in the same way on every iteration of the for (int secondPrime loop. Defining and using a function such as getPrimeAfter(int) will help avoid this error (if it is an error) and make it clearer what is actually supposed to happen.

It looks like the variable named square is actually multipleOfSquare; that is, it takes on the values prime*prime, 3*prime*prime, 5*prime*prime, and so forth.

Consider this line:

                if (secondPrime * prime == square && square < maxNumber / 2) {


If secondPrime is, in fact, a prime, then this condition can be true only on the first iteration, when secondPrime == prime. That would suggest that it might be better to execute the code of this if branch just once before the loop and only do the else code during the loop. (On the other hand, since the code as it stands at this moment does not appear to guarantee that secondPrime is prime, it's hard to determine exactly what would happen if you reorganized that loop as I just suggested.)

These three comments have me a bit confused:

                // if the product of 'prime' and 'secondPrime' equals
// 'square'

// searches for all multiples of square that may exist in
// array

// if 'prime' and 'secondPrime' are not equal to square
// then search for all multiples of 'prime' by 'secondPrime'


The first comment is pretty useless, since it just says exactly what the first part of the if condition below says. In the second comment, we know that secondPrime * prime == square, so actually this is searching for all multiples of secondPrime*prime. But that's just what the third comment says. So do the if and else branches of if (secondPrime * prime == square (the branch with the second comment and the branch with the third comment) actually do anything different from each other? Do you really need two branches? If so, that's what you should explain in the comments.

You might be able to save some computation here:

    for (int prime = 5; prime * prime <= maxNumber; prime += 2) {
// the if-statement below skips the testing of all non-prime numbers
if (isPrime[prime]) {


If you compute the square root of maxNumber just once at the very start of the sieve, you could write prime <= maxSquareRoot instead of prime * prime <= maxNumber, and save a multiplication on every iteration. Another question is whether you could use prime = getPrimeAfter(prime) instead of prime += 2. Before doing that, you should show that getPrimeAfter(prime) will always terminate before its internal iterator runs off the end of the isPrime array. But come to think of it, you really should also establish that fact in the place where you write

                    while (!isPrime[secondPrime]) {


deep within the bowels of your sieve. On the other hand, it might slow your code down a tiny bit if you change that line to

                    while (secondPrime <= maxNumber && !isPrime[secondPrime]) {


but it would guarantee no running off the end of the array, and it would take far less development time to achieve that guarantee. Alternatively, you could allocate an additional element in the isPrime array that would act as a "guard byte": that is, you allocate space for isPrime[maxNumber + 1] and set isPrime[maxNumber + 1] = true, and rely on your other logic to ensure that you will never set isPrime[maxNumber + 1] = false (and also never actually use maxNumber + 1 as a prime number).

• Thank you for your time and assistance. I have to take a lot of it into account and try it out. – Alex Lieberman Jun 20 '14 at 22:26

Regarding "I believe this might be the fastest (singular, non-threaded) array implementation of prime sieving to date", I disagree, at least for C implementations.

First, never use printing for benchmarking unless that is specifically what you want to test. For good sieves, the I/O will take more time than the sieving. A custom written buffered syswrite can be far faster than printf or cout, and I imagine there is something analogous in Java. I can beat primesieve by over 2x if we include printing, merely because that isn't something primesieve has optimized at all (it's irrelevant for reasonable amounts of output so there doesn't seem like there is a need).

Counting can be also be done faster than you're doing it. Clearly using advanced prime counting methods are orders of magnitude faster than any sieve, but I mean that you can use methods to count bytes or words which can be faster, and most sieve examples trying to show fast performance do this.

What other sieves did you test? How did you test them to reach this conclusion?

It's possible I messed something up in the translation to C, but here are my results for counting to 1011, using the June 26th code, all single threaded:

Monolithic:

423s   Lieberman
416s   Trivial odd SoE
234s   Mathisen 1998
147s   MPU 0.36 monolithic SoE
96s   Pari 2.6.2 monolithic


Segmented:

208s   Walisch byte segment example
104s   Oliveira e Silva v1 example
92s   Walisch bit segment example
81s   tverniquet
79s   Flammenkamp
78s   Mario Roy algorithm 3
55s   Oliveira e Silva v2 example
43s   Math::Prime::Util
34s   primegen 0.97 (tuned cache size)
23s   primesieve 5.0
20s   yafu


Note that this is just at one size, and starting from 2. Doing ranges will be vastly more efficient with the segmented sieves, and doing counts to larger values starts to choke the monolithic sieves. The segmented sieves of course use less memory. At 1012 yafu and primesieve are at 5 minutes; primegen, MPU, and OeS's v2 are all at 10 minutes; while the others are at 16+ minutes or worse.

• Thank you for your insight. I suppose I never found the other sieves as well as you have. I mostly found what was readily available on the surface. Basic SoE and Atkins. I haven't implemented my code in any other language and I wonder if it was implemented correctly - but yes I understand there are faster methods that I have more recently found. One question is - can these other methods also keep the data stored in an array at the same time? That is what I was referring to - something that finds the primes and stores them as well in some way in that time allotted. – Alex Lieberman Aug 27 '14 at 21:47
• And I have been unable to test up to 10^12 because of the limit of my computing power - I only tested up to 1,000,000,000. – Alex Lieberman Aug 27 '14 at 21:48
• I'm interested in how you implemented it in C, because I have looked up trivial odd SoE and I don't understand the logic of why mine would be slower than this sieve. – Alex Lieberman Aug 27 '14 at 21:52
• It seems all the good sieves on the web are hidden. Oh well. – Alex Lieberman Aug 27 '14 at 22:50
• can you share results up to 1 Billion? Does your formula store all the data for later use or just find the number of primes? – Alex Lieberman Aug 28 '14 at 2:02

Updated the code with as many suggestions as I felt optimized.

I tested it a lot of different ways and some ways it was slower - even with more methods.

import java.util.*;
import java.lang.*;
import java.math.*;

public class SixSieveFinal {

static int maxNumber = 1000000000;
static boolean[] isPrime = new boolean[maxNumber];
static int limit = (int) Math.sqrt(maxNumber);
static long timeStart, timeEnd, timeDifference;
static double timeElapsed;
static int totalPrimes;

public static void main(String[] args) {

isPrime = true;
isPrime = true;

System.out.println("Starting Computation.");
timeStart = System.currentTimeMillis();

getPotentialPrimes();
eliminateProducts();

timeEnd = System.currentTimeMillis();
timeDifference = timeEnd - timeStart;
timeElapsed = timeDifference / 1000.000000;

System.out.println("Finished.");
System.out.println("Time elapsed: " + timeElapsed);

totalPrimes = countPrimes(totalPrimes);
System.out.println("Total primes: " + totalPrimes);

queryPrint();

}

static void getPotentialPrimes() {

for (int countBySix = 6, countToSkip = 2; countBySix < maxNumber;) {

if (countToSkip == 5) {

isPrime[(countBySix) - 1] = true;
countToSkip++;

} else if (countToSkip == 7) {

isPrime[(countBySix) + 1] = true;
countToSkip = 3;

} else {

isPrime[(countBySix) - 1] = true;
isPrime[(countBySix) + 1] = true;
countToSkip++;

}

countBySix += 6;

}

}

static void eliminateProducts() {

int productOfSquare;
int possiblePrime;

for (int primeOne = 7; primeOne * primeOne <= maxNumber; primeOne += 2) {

if (isPrime[primeOne]) {

if (primeOne <= limit) {

productOfSquare = primeOne * primeOne;

while (productOfSquare < maxNumber && productOfSquare > 0) {

isPrime[productOfSquare] = false;
productOfSquare += (primeOne * primeOne * 2);

}

}

possiblePrime = primeOne;

do {

possiblePrime += 2;

} while (!isPrime[possiblePrime]);

while (possiblePrime * primeOne <= maxNumber) {

isPrime[possiblePrime * primeOne] = false;

do {

possiblePrime += 2;

} while (!isPrime[possiblePrime]);
}

}
}

}

static void queryPrint() {

System.out.println("Would you like to print the primes (Y/N)?");

Scanner inputQuery = new Scanner(System.in);
char printResponse;

do {

printResponse = inputQuery.next().charAt(0);
if (printResponse == 'Y') {

printPrimes();

} else {

System.out.println("Exiting program.");
return;

}

} while (printResponse != 'N' || printResponse != 'N');

}

static void printPrimes() {

for (int primeIndex = 0; primeIndex < maxNumber; primeIndex++) {

if (isPrime[primeIndex]) {

System.out.print(primeIndex + " ");

}

}

return;

}

static int countPrimes(int countPrimes) {

for (int primeIndex = 0; primeIndex < maxNumber; primeIndex++) {

if (isPrime[primeIndex]) {

countPrimes++;

}

}

return countPrimes;

}

}

• Any comments on new code? – Alex Lieberman Jul 12 '14 at 5:47
• We cannot review an answer. If you'd like further review, please post a new follow-up question. – Jamal Aug 27 '14 at 5:24
• Yes, repost it as a new question, I'm looking forward to tearing it down :D:D (no, seriously, I liked this question and the special handling of squares was an excellent idea). – maaartinus Oct 17 '14 at 18:01