# 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 - Prime Number Sieving Algorithm - 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)$

Code:

//-----------------------------------------------------------------
//  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[2] = true;
isPrime[3] = 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 at 20:10
Peer review of the math more than anything and interest in comparisons between other implementations. –  Alex Lieberman Jun 19 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 at 20:46
Didn't use the printing. Figured that out recently too. –  Alex Lieberman Jun 19 at 21:00
@AlexLieberman This site is about code review. For math questions/peer review, please use math.stackexchange.com. –  Emily L. Jun 19 at 21:15

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 at 20:26
@DavidHarkness: Aye, good point. –  Bobby Jun 19 at 20:29
@Bobby It's a concatenation enforcing an implicit cast, see here and it's damn common. –  maaartinus Jun 19 at 20:43
@maaartinus: Yes, I still don't like it. I have seen too many horrors because of implicit casting. –  Bobby Jun 20 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 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 at 1:13
True. I appreciate your comments and find them useful :) –  Alex Lieberman Jun 20 at 1:24
@AlexLieberman Glad to hear it... now I was being rather rude. Did you invent the algorithm? –  maaartinus Jun 20 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 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 at 22:26

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[2] = true;
isPrime[3] = 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 at 5:47