# The sieve of Eratosthenes

### Question

• Hi all,

I have some code here that I found online for the sieve but it seems limited to in the range of 1000 or so.  Do you know which part of it requires change to let there be no limit to the values for the sieve?

Or do you know what needs changing to allow for a great value to be the limit of the sieve?  I'm thinking of numbers approaching 1000 digits. I realise that would take an awfully long time to output.

```/* Sieve Of Erathosthenes by Denis Sureau */
#include <stdlib.h>
#include <stdio.h>
#include <iostream>
#include <vector>
#include "BigIntegerLibrary.hh"
using namespace std;
void eratosthenes(int top)
{

vector<int> all(top);
int idx = 0;

std::cout << "1 ";

int prime = 3;

while(prime <= top)
{
for(int x = 0; x < top; x++)
{
if(all[x] == prime) goto skip;
}

std::cout << prime << " ";
int j = prime;
while(j <= (top / prime))
{
all[idx++] = prime * j;
j += 1;
}

skip:
prime+=2;
}
std::cout << std::endl;
return;
}

int main(int argc, char **argv)
{
if(argc == 2) eratosthenes(atoi(argv[1]));
else eratosthenes(1000);
cin.get();
return 0;
}
```

Thursday, January 13, 2011 9:56 PM

• On 1/13/2011 6:11 PM, Stephan T. Lavavej - MSFT wrote:

I can imagine the Mersenne prime 2^43,112,609 with 12,978,189 decimal
digits just fine.

Can you really? How thick a stack of paper would it take to print it out?

And its primality was verified by a computer
(though with the Lucas-Lehmer test, not the Sieve of Eratosthenes).

It's much easier to verify whether a given number is prime, than it is to find a prime in the first place. That's why people concentrate on numbers of a special form, which are analytically shown to be much more likely to be prime than a number of comparable magnitude just picked at random.

Sieve of Eratosthenes is an algorithm for finding primes, and is largely a toy, impractical beyond some, rather low, threshold.

Igor Tandetnik

Thursday, January 13, 2011 11:32 PM
• > Can you really? How thick a stack of paper would it take to print it out?

Using a reasonable font size, and estimates of 80 columns and 56 rows, and two-sided printing, 1449 pages.

> It's much easier to verify whether a given number is prime, than it is to find a prime in the first place.

Actually, finding primes isn't too hard.

The probability of a number around 10^1000 (~= 2^3322, so a very realistic size) being prime is around 1 / ln(10^1000) == 1 / (1000 * ln(10)) == 1 / 2302.585... (see http://en.wikipedia.org/wiki/Prime_number_theorem ) and that's counting numbers with tiny divisors like 2, 3, 5. It gets better, because it takes only one failed probabilistic primality test to exclude a number from consideration (see http://en.wikipedia.org/wiki/Primality_test ). After narrowing down the candidates, then you can throw a bunch of probabilistic primality tests at a survivor (this is "good enough" for encryption), or a deterministic test if you want that (the probabilistic tests beforehand are cheap in comparison).

> That's why people concentrate on numbers of a special form, which are analytically
> shown to be much more likely to be prime than a number of comparable magnitude just picked at random.

Actually, Mersenne numbers are special because they have an exceedingly fast deterministic primality test (the fastest one known).

I don't know if Mersenne numbers can be said to be "much more likely" to be prime than ordinary numbers (this is beyond my amateur mathematical ability). Let's play around with C++:

```C:\Temp>type primes.cpp
#include <math.h>
#include <iostream>
#include <ostream>
using namespace std;

bool prime(const int n) {
if (n < 2) {
return false;
}

for (int i = 2; i <= n / i; ++i) {
if (n % i == 0) {
return false;
}
}

return true;
};

int main() {
double total_prime_exponents = 0.0;
double total_all_exponents = 0.0;

for (int i = 2; i <= 43112609; ++i) {
const double d = i < 31 ? log((1 << i) - 1.0) : i * log(2.0);

if (prime(i)) {
total_prime_exponents += 1 / d;
}

total_all_exponents += 1 / d;
}

cout << "total_prime_exponents: " << total_prime_exponents << endl;
cout << " total_all_exponents: " << total_all_exponents << endl;
}

C:\Temp>cl /EHsc /nologo /W4 /MT /O2 /GL primes.cpp
primes.cpp
Generating code
Finished generating code

C:\Temp>primes
total_prime_exponents: 4.73798
total_all_exponents: 24.9863
```

Looking at only those Mersenne numbers with prime exponents, the Prime Number Theorem says we should expect about 5 of them to be prime, when in fact 47 of them in this range are known to be prime (there may be more). It certainly appears that by exploiting an easy-to-prove fact about Mersenne numbers (for 2^N - 1 to be prime, N must be prime) and by performing a tiny amount of work (testing N for primality), the remaining Mersenne numbers with prime exponents are about 10 times more likely to be prime than numbers chosen at random.

Without exploiting that fact, Mersenne numbers in general don't appear to be "especially prime". The ultra-trivial fact about them is that they're all odd, and the chance of a random odd number being prime is double the chance of any random number being prime, so the Prime Number Theorem tells us to expect about 50 primes in this range, and we've found 47.

(Ooh, conjecture! It would be awesome if there were exactly 3 currently-unknown Mersenne primes with exponents between 20,996,011 and 43,112,609.)

That was fun.

> Sieve of Eratosthenes is an algorithm for finding primes,
> and is largely a toy, impractical beyond some, rather low, threshold.

Of course.

Saturday, January 15, 2011 12:59 AM
• I have a small sieve program, I might post it on my forum if I am bored. I use it in my Fibonacci program that is free on my developer site.

My MVP is for Windows XP, Vista and Windows 7 IT, and I am getting increasingly good with Visual Studio.

Saturday, January 15, 2011 2:47 AM

### All replies

• On 1/13/2011 4:56 PM, Michael.1975 wrote:

I have some code here that I found online for the sieve but it seems
limited to in the range of 1000 or so.

That doesn't seem to be the case. The program accepts the range on the command line, and only defaults to 1000 when none is given. Have you tried running it with a larger value as a command line parameter?

Or do you know what needs changing to allow for a great value to be
the limit of the sieve?  I'm thinking of numbers approaching 1000
digits.

I doubt that. Very few people (if any) can imagine numbers appoaching 10^1000. For comparison, the observable universe is thought to contain approximately 10^80 atoms.

I realise that would take an awfully long time to output.

In fact, probably longer than the age of the universe. In other words, completely impractical.

Igor Tandetnik

Thursday, January 13, 2011 10:54 PM
• I don't know if this is right or not, I'm just skimming and I found this:

```int main(int argc, char **argv)
{
if(argc == 2) eratosthenes(atoi(argv[1]));
else eratosthenes(1000);
cin.get();
return 0;
}

```

if you change the "else eratosthenes(1000)" To whatever you want number to change it to.

I don't know though, just skimming

Thursday, January 13, 2011 10:57 PM
• I can imagine the Mersenne prime 2^43,112,609 with 12,978,189 decimal digits just fine. And its primality was verified by a computer (though with the Lucas-Lehmer test, not the Sieve of Eratosthenes).
Thursday, January 13, 2011 11:11 PM
• On 1/13/2011 6:11 PM, Stephan T. Lavavej - MSFT wrote:

I can imagine the Mersenne prime 2^43,112,609 with 12,978,189 decimal
digits just fine.

Can you really? How thick a stack of paper would it take to print it out?

And its primality was verified by a computer
(though with the Lucas-Lehmer test, not the Sieve of Eratosthenes).

It's much easier to verify whether a given number is prime, than it is to find a prime in the first place. That's why people concentrate on numbers of a special form, which are analytically shown to be much more likely to be prime than a number of comparable magnitude just picked at random.

Sieve of Eratosthenes is an algorithm for finding primes, and is largely a toy, impractical beyond some, rather low, threshold.

Igor Tandetnik

Thursday, January 13, 2011 11:32 PM
• > Can you really? How thick a stack of paper would it take to print it out?

Using a reasonable font size, and estimates of 80 columns and 56 rows, and two-sided printing, 1449 pages.

> It's much easier to verify whether a given number is prime, than it is to find a prime in the first place.

Actually, finding primes isn't too hard.

The probability of a number around 10^1000 (~= 2^3322, so a very realistic size) being prime is around 1 / ln(10^1000) == 1 / (1000 * ln(10)) == 1 / 2302.585... (see http://en.wikipedia.org/wiki/Prime_number_theorem ) and that's counting numbers with tiny divisors like 2, 3, 5. It gets better, because it takes only one failed probabilistic primality test to exclude a number from consideration (see http://en.wikipedia.org/wiki/Primality_test ). After narrowing down the candidates, then you can throw a bunch of probabilistic primality tests at a survivor (this is "good enough" for encryption), or a deterministic test if you want that (the probabilistic tests beforehand are cheap in comparison).

> That's why people concentrate on numbers of a special form, which are analytically
> shown to be much more likely to be prime than a number of comparable magnitude just picked at random.

Actually, Mersenne numbers are special because they have an exceedingly fast deterministic primality test (the fastest one known).

I don't know if Mersenne numbers can be said to be "much more likely" to be prime than ordinary numbers (this is beyond my amateur mathematical ability). Let's play around with C++:

```C:\Temp>type primes.cpp
#include <math.h>
#include <iostream>
#include <ostream>
using namespace std;

bool prime(const int n) {
if (n < 2) {
return false;
}

for (int i = 2; i <= n / i; ++i) {
if (n % i == 0) {
return false;
}
}

return true;
};

int main() {
double total_prime_exponents = 0.0;
double total_all_exponents = 0.0;

for (int i = 2; i <= 43112609; ++i) {
const double d = i < 31 ? log((1 << i) - 1.0) : i * log(2.0);

if (prime(i)) {
total_prime_exponents += 1 / d;
}

total_all_exponents += 1 / d;
}

cout << "total_prime_exponents: " << total_prime_exponents << endl;
cout << " total_all_exponents: " << total_all_exponents << endl;
}

C:\Temp>cl /EHsc /nologo /W4 /MT /O2 /GL primes.cpp
primes.cpp
Generating code
Finished generating code

C:\Temp>primes
total_prime_exponents: 4.73798
total_all_exponents: 24.9863
```

Looking at only those Mersenne numbers with prime exponents, the Prime Number Theorem says we should expect about 5 of them to be prime, when in fact 47 of them in this range are known to be prime (there may be more). It certainly appears that by exploiting an easy-to-prove fact about Mersenne numbers (for 2^N - 1 to be prime, N must be prime) and by performing a tiny amount of work (testing N for primality), the remaining Mersenne numbers with prime exponents are about 10 times more likely to be prime than numbers chosen at random.

Without exploiting that fact, Mersenne numbers in general don't appear to be "especially prime". The ultra-trivial fact about them is that they're all odd, and the chance of a random odd number being prime is double the chance of any random number being prime, so the Prime Number Theorem tells us to expect about 50 primes in this range, and we've found 47.

(Ooh, conjecture! It would be awesome if there were exactly 3 currently-unknown Mersenne primes with exponents between 20,996,011 and 43,112,609.)

That was fun.

> Sieve of Eratosthenes is an algorithm for finding primes,
> and is largely a toy, impractical beyond some, rather low, threshold.

Of course.

Saturday, January 15, 2011 12:59 AM