# Prime Sieve and brute force

Got caught up with this morning's C++ question.
So I thought I would have a go.

SO here is my attempt at printing the first 1001 primes.

#include <iostream>
#include <vector>
#include <cmath>

bool isPrime(std::size_t num, std::vector<std::size_t> const& primes)
{
std::size_t max = std::sqrt(num);
for(auto const& prime: primes)
{
if (prime > max) {
return true;
}
if (num % prime == 0)
{
return false;
}
}
return true;
}

std::size_t sieve(std::size_t size, std::vector<std::size_t>& results)
{
// primeCandidates holds 0/1 for each potential prime candidate.
// The index represents the potential prime (index * 2 + 1)
// This allows us to ignore all factors of 2
// max is one past the highest prime candidate we can test for and store in primeCandidates
std::size_t const           max = size * 2 + 1;
std::vector<std::size_t>    primeCandidates(size, 1);

// Add some default well know primes.
results.push_back(2);
results.push_back(3);

// We will use the technique of incrementing by 2 then 4 then 2 then 4
// This means skip all factors of 2 and 3 automatically.
std::size_t       inc = 2;
std::size_t loop = 5;
for(; loop < max; loop += inc, inc = 6 - inc) {
std::size_t index = loop/2;

// If we find a candidate that is valid then add it to results.
if (primeCandidates[index]) {
results.push_back(loop);

// Now strain out all factors of the prime we just found.
for(; index < primeCandidates.size(); index += loop) {
primeCandidates[index] = 0;
}
}
}
return loop;
}

int main()
{
std::size_t constexpr   primeRequired   = 1001;
std::size_t constexpr   siveStartSize   = 1000;

std::vector<std::size_t>    result;
result.reserve(primeRequired);

// Use the Sieve of Eratosthenes to get a basic set of primes.
std::size_t next = sieve(siveStartSize, result);

// Want to re-use the 2/4/2/4 increment pattern
// So work out the correct start point and increment value.
std::size_t inc  = next % 6 == 5 ? 2 : 4;

// We now use brute force checking each number against all the primes
// that we have already found.
for(; result.size() <= primeRequired; next += inc, inc  = 6 - inc) {
if (isPrime(next, result)) {
result.emplace_back(next);
}
}

// Print out the primes we have found
for(auto val: result) {
std::cout << val << " ";
}
std::cout << "\n";
}

• Missing <cstddef> for std::size_t. – yuri Jun 4 '19 at 19:38
• Yeah ... now that's a mixture of the sieve and brute force. From the first sentence here I thought it's possible to implement solely with the sieve part. – Thomas Weller Jun 4 '19 at 20:16
• @ThomasWeller The "Sieve of Eratosthenes" can be implemented in many ways. I just took a simplistic approach here. There is a way of estimating the number of primes lower than a particular value. See Prime Counting but there is also the factor of how much space your computer has. So for the sieve you don't want to allocate more space than your computer has and you don't want to overestimate (too much) how many primes there are below the value you want. This usually means a guess that is low and then brute forcing the rest. – Martin York Jun 4 '19 at 20:22
• You can make a good guess by for the size of the array by solving x/ln(x) = <No of Primes> X= <Size of Array>. But is that the optimal size of the array? At some point you are doing a lot of work maintaining the "potential candidates" when it could be simpler to brute force? (not sure about that). – Martin York Jun 4 '19 at 20:34
• How is the performance? – pacmaninbw Jun 13 '19 at 14:30

When you find a prime number, say 7, you begin crossing out all odd multiples of 7 from the primeCandidates vector. You do this by incrementing index += loop, where loop is the prime number, and but primeCandidates only holds odd candidates, so the step size in natural numbers would be 2*loop.

The issue is you start by removing 7 from primeCandidates, then 21, then 35. At this point in the sieve, you've already eliminated multiples of 3 and 5, so marking off 3*7=21 and 5*7=35 is just busy work. And marking off 1*7=7 is just pointless. What you should be doing is starting to eliminate at 7*7=49. It saves only three eliminations during the 7 loop, but the saving get more substantial as the prime gets larger. For instance, when you get to 199, you wouldn't eliminate 1*199, 3*199, 5*199, 7*199, ... 191*199, 193*199, 195*199, 197*199. You would just start at 199*199 and go up from there.

Code change:

// Now strain out all factors of the prime we just found, starting with prime^2
index = loop * loop / 2


std::vector<std::size_t> is overkill for primeCandidates. You only ever store 1 or 0, so allocating 4 or 8 bytes for each candidate is wasting memory.

std::vector<char> would reduce the memory load to 1 byte per candidate.

std::bitset<N> will reduce this to 1 bit per candidate. Alternately, the vector specialization std::vector<bool> will also give 1 bit per candidate, with a non-compile-time fixed number of bits.

• Because the primeCandidates only holds odd values incrementing by 7 is actually equivalent of incrementing by 14 each time. Since we start at 7 then increment by 14 each time we always hit every odd value (there is no point in hitting the even ones). Note we get the index first by starting at index = 7/2 then add 2*7 to the index each loop through. That is actually a nice optimization that fell out of writting this. – Martin York Jun 5 '19 at 16:08
• Nice optimization to start at loop * loop in the sieve – Martin York Jun 5 '19 at 16:10
• @MartinYork My first paragraph wasn’t a review comment on the code; it was just an acknowledgement that I grokked that primeCandidates held just the odd candidates, so you can increment by 7 to eliminate just the odd multiples of 7 (which occur every 14 natural numbers). Perhaps I need to reword it. The second paragraph starts “The issue is...”, which is what I was pointing out in the review. – AJNeufeld Jun 5 '19 at 17:03

# Commenting

You could improve the functions with some introductory comments. In particular, the isPrime() predicate has an extra argument compared with the conceptual version - we should be clear what that's for (i.e. it's an ordered set of all primes up to √num). Similarly, sieve()'s results argument is assumed to be empty, but that's not clearly communicated. For a program this small, that's probably a very minor consideration, but still a good habit to be in.

# Presentation

Choose one brace style and stick with it. Here, we have a mix of open brace styles:

    if (prime > max) {
return true;
}
if (num % prime == 0)
{
return false;
}


I'm also not very keen on the trailing whitespace, though that's easily fixed.

Some spellings are, let's say, unconventional: siveStartSize really ought to be spelt sieveStartSize, and well know primes probably means well-known primes. I think that ignore all factors of 2 actually meant all multiples of 2; the same substitution is present in skip all factors of 2 and 3.

# Storage

Is there a good reason why primeCandidates stores std::size_t values holding either 0 or 1? std::vector<char> works just as well for me.

Following this, I'd consider replacing

results.push_back(2);
results.push_back(3);


with

results = {2, 3};


so we don't have to assume that results is initially empty.

# Structure

It seems strange that half of the work is in main(), where we re-compute the value of inc we had inside sieve in order to continue with the second half of the algorithm. I'd be inclined to keep that loop in sieve(), separated by a comment. Then main() simply do its job of choosing parameters, calling the function and printing results.

# Flexibility

Why are primeRequired and siveStartSize constants? A useful program allows the user to obtain their choice of result.

First of all: great code. I'd love to read code like this in our implementations. For me that's a nice mixture of useful comments when needed with self-explaining code.

I also like the modern way of writing using constexpr instead of #define, which I still see a lot.

I'm not a C++ pro, rather coming from the C# side, so I notice the following:

• std::size_t is IMHO thought for the sizeof. You use it almost everywhere. I'd prefer to read int, long or even using ll = long long;. Using size_t for me adds semantic: this is of type size, so I e.g. use it as the end condition for a loop, use it for memory allocation or similar. That's not the case in the prime program.

• isPrime() takes a number and a list of primes, but it's not documented what needs to be in that list in order to make the function work. I could potentially call it with a large number but an empty list.

• I dislike crippled for loops. Even worse with two statements in the increment part. What's wrong with a while loop?

while(result.size() <= primeRequired) {
if (isPrime(next, result)) {
result.emplace_back(next);
}
next += inc;
inc  = 6 - inc;
}

• typo: primeRequired should IMHO be primesRequired, because you don't want to go up to a number but up to a count.

• you could split the main method in two methods, one for calculating (testable by unit tests) and one for printing

Just a though: instead of using math for doing the inc magic, would C++ support something like

int nextinc() {
while(true) {
yield 2;
yield 4;
}
}

• int nextInc() {static inc = 4; inc = 6 - inc; return inc;} – Martin York Jun 4 '19 at 20:42
• Can't argue with any of that. :-) – Martin York Jun 4 '19 at 20:47
• size_t is preferable to a signed type, since there's nothing here that needs negative numbers. It will also (generally) be a size that is easy/natural for the CPU to access and use. And I much prefer the for loop to the while. Fewer problems when you later add a continue. – 1201ProgramAlarm Jun 5 '19 at 4:01
• In response to the question about yield, C++ 20 has co_yield although it requires "a little" more code to work in practice – sudo rm -rf slash Jun 5 '19 at 6:09
• @1201ProgramAlarm: There is still unsigned int if you want positive numbers only. – Thomas Weller Jun 5 '19 at 7:07