There are a few things you could do to improve your code.
Use braces
Constructions like this are a bug waiting to happen:
for (unsigned long int k = i + j; k < MAX; k += j)
if (a[k])
//Open the following comment very carefully
//std::cout << "removes " << 2 * k + 1 << "\n";
a[k] = false;
If anyone ever uncomments that line, the a[k] = false will then be outside the for loop and the program will not work correctly. Such errors are easily avoided by the use of braces:
for (unsigned long int k = i + j; k < MAX; k += j) {
if (a[k]) {
//no particular care is now required
//std::cout << "removes " << 2 * k + 1 << "\n";
a[k] = false;
}
}
Check for allocation failure
The code doesn't currently check for an allocation failure, but it should. Calls to malloc or calloc can fail, and your program should be robust enough not to fail if that happens.
Reduce the range of search for primes
We already know that the largest usable prime is \$<\sqrt{\text{MAX}}\$, but then there is this line:
while (((++i)<MAX) && (!a[i]));
We don't really need to search all the way to MAX. Also, for clarity, I'd rewrite that as a for loop:
// skip over composites to next prime
for (++i; i<sqrtMAX && a[i]; ++i)
{}
Skip more composite numbers
In the actual sieve, the code skips forward and then marks all odd multiples of the most recently discovered prime. However, the code currently doesn't skip as many as it could.
For instance, if the most recently discovered prime is 7, the code currently skips to 21, but we can do better. We can skip to 7*7 = 49 because any other smaller multiple of 7 will already have been eliminated. That transforms the inner loop from this:
for (unsigned long int k = i + j; k < MAX; k += j) {
To this:
for (unsigned long int k = 2*i*(i+1); k < MAX; k += j) {
Print as you go
Rather than searching the entire array at the end, only the part of the array after \$\sqrt{\text{MAX}}\$ really needs to be searched since all of the other primes were identified during sieveing. We can print them as they are discovered and save some time searching.
Increment rather than multiply
In the various loops, the expression 2*i+1 (stored in j in the first loop) is used to convert from an index to the actual numerical value it represents. It speeds things up a bit if the code simply increments j each loop increment.
Clean up the loop
Instead of using nested while loops, it seems to me to be cleaner if you use a for loop and inner if to determine whether the number is composite and should simply be skipped, or prime and should be sieved.
A new sieve:
Using all of these suggestions, I get this, which is about 8% faster on my machine than the original.
void sieve2()
{
bool *a = (bool*)malloc(MAX * sizeof(bool));
if (a == NULL) {
return; // bad alloc
}
memset(a, true, MAX);
const unsigned long sqrtMAX = ceil(sqrt(MAX));
unsigned long i;
unsigned long j;
std::cout << "done\n";
for (i=1, j=2*i+1; i < sqrtMAX; ++i, j+=2) {
if (a[i]) { // it's a prime
std::cout << j << '\n';
for (unsigned long int k = 2*i*(i+1); k < MAX; k += j) {
if (a[k]) {
a[k] = false;
}
}
}
}
for ( ; i < MAX; ++i, j+=2) {
if (a[i]) { // it's a prime
std::cout << j << '\n';
}
}
free(a);
}
Further improvements
There are many other things you might explore to further speed the code:
- Use bits rather than bytes to store values.
- Reimplement with pointers rather than array indexing.
- Create a multithreaded version
Use classes
Since this is supposed to be C++, why not use a class? Here's a very simple class template that implements something like a bitset.
template <typename T>
class Bits
{
public:
Bits(size_t sz) :
bitsperunit{8*sizeof(T)},
units{sz/bitsperunit + 1},
bits{new T[units]}
{
for (unsigned i=0; i < units; ++i) {
bits[i] = ~0;
}
}
void clear(unsigned index) {
unsigned n=index/bitsperunit;
T mask=1u<<(index%bitsperunit);
bits[n] &= ~mask;
}
bool operator[](unsigned index) const {
unsigned n=index/bitsperunit;
T mask=1u<<(index%bitsperunit);
return bits[n] & mask;
}
virtual ~Bits() {
delete[] bits;
}
private:
const unsigned bitsperunit;
const size_t units;
T *bits;
};
Now let's reimplement using that class:
void sieve3()
{
Bits<uint8_t> a(MAX);
const unsigned long sqrtMAX = std::ceil(std::sqrt(MAX));
unsigned long i;
unsigned long j;
for (i=1, j=2*i+1; i < sqrtMAX; ++i, j+=2) {
if (a[i]) { // it's a prime
for (unsigned long int k = 2*i*(i+1); k < MAX; k += j) {
a.clear(k);
}
}
}
std::cout << "done\n";
for (i=1, j=3 ; i < MAX; ++i, j+=2) {
if (a[i]) { // it's a prime
std::cout << j << '\n';
}
}
}
On my machine this is still faster than either of the other two implementations (original or the one above).
std::vector<bool>? \$\endgroup\$reserve ()orresize()once to prevent you from recurring reallocations. \$\endgroup\$std::vector<bool>has to use bit references. This mandatory inefficiency pretty much relegates its use to homework assignments and toy problems, since the overhead is difficult to eliminate even for current optimizing compilers.std::bitarray<>is better but fixed in size. Hence we have to roll our own or use third party libraries like boost. A simpler and more flexible alternative is templated bit accessor functions plusvector<unsigned>, as shown in Sieve of Eratosthenes - segmented with pre-sieving #2 \$\endgroup\$std::bitset<>, notbitarray<>\$\endgroup\$