HiveBrain v1.2.0
Get Started
← Back to all entries
patterncppMinor

Non-Sieve Prime Number generator

Submitted by: @import:stackexchange-codereview··
0
Viewed 0 times
numbernonsievegeneratorprime

Problem

#include 
#include 
#include 
using namespace std;
int main()
{
int l=0;
bool isprime;
vectorprimes;
primes.push_back(2);
for (int nur = 3;true;nur+=2)
{

     isprime=true;
     double ilimit= sqrt(nur);
    for (int primecount=0 ;primes[primecount]<=ilimit;primecount++)
    {
        if (nur % primes[primecount] == 0)
        {
            isprime = false;
            break;
        }

    }

     if (isprime)
     {
             cout<<nur<<endl;
          primes.push_back(nur);
             // l++;
     }

     if (nur<0)
     break;

}
cout<<"Finished"<<endl;
//cout<<primes[l];

  return 0;
 }


Is there anyway i can optimize this further without using a sieve?

Solution

Your code includes the optimization of evaluating only odd-valued candidates (by iterating with a step size of 2) because the even-valued candidates are known to contain a factor of 2, and are therefore always composite.

This candidate-skipping optimization can be generalized to include prime factors greater than 2, by varying the step size in a certain way.

I made a version of your code which skips candidates (values of nur) containing factors of 2 and/or 3. This is done fairly cheaply by alternating step size between 2 and 4.

Extending this optimization to even greater primes gives diminishing returns for two reasons: The step size adjustment calculation becomes more complex (and slower), eventually costing more time than what is saved by skipping candidates; and higher prime factors exist more sparsely than low ones in the candidate population. Skipping 2 alone catches 1/2 of the composite candidates, and skipping 3 alone would catch 1/3. But if used together, half of the composites 3 would skip were already skipped by the 2 skip, so skipping 3 just catches an additional 1/6 of the composite candidates. Including a 5 skip would only eliminate an additional 1/30 of the composite candidates, and a 7 skip 1/210.

I experimented with replacing the floating-point sqrt call with an integer approximation, but it turns out that sqrt is very fast (at least on my machine). And then I gave up entirely after seeing Vedran Šego's optimization, which is brilliant, especially because as nur gets bigger, the gaps between calls to sqrt become larger, to the point where it doesn't really matter much how slow your sqrt calculation is. So, I had to add this to my example, thank you Vedran Šego.

One other thing, the nur loop must terminate. If nur is signed (like int), it will overflow and wrap around to a negative value. This is a completely pathological situation -- calling sqrt with a negative number, plus other things, depedning on optimizations employed.

Even for unsigned nur, duplicate primes will be added to the primes vector until the process crashes when it runs out of memory.

So my example terminates before nur overflows. I also defined a prime_t type, to make it easy to switch the integer type used for calculating and storing primes. I use unsigned int in my example, because it gives results a larger maximum value at no extra cost compared to int. Try unsigned long long if you really want to see the program continue for much longer than you will be willing to wait.

#include 
#include 
#include 

typedef unsigned int prime_t;
// max_prime_t is the maximum value which may be held in a variable of type prime_t
// this assumes prime_t is an integer type (signed or unsigned)
// It dupicates information available from std::numeric_limits
const prime_t max_prime_t =
     prime_t(0)-prime_t(1) > prime_t(0) ?
         ~prime_t(0) :
         ~prime_t(0) ^ (prime_t(1)  primes = { 2, 3 };
    prime_t ilimit = 0;
    prime_t ilimit_limit = 0;

    static const prime_t max_nur = max_prime_t - 4;

    // "step" alternates between 2 and 4,
    // to skip values containing factors of 2 and/or 3
    for(prime_t nur=5, step=4; nur= ilimit_limit) {
            // Vedran Sego's optimization
            ilimit = static_cast(std::sqrt(nur));
            ilimit_limit = (ilimit + 1) * (ilimit + 1);  
        }
        for(auto one_prime : primes) {
            if(one_prime > ilimit) {
                std::cout << nur << '\n';
                primes.push_back(nur);
                break;
            }
            if(nur % one_prime == 0) {
                break;
            }
        }
    }
}

Code Snippets

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

typedef unsigned int prime_t;
// max_prime_t is the maximum value which may be held in a variable of type prime_t
// this assumes prime_t is an integer type (signed or unsigned)
// It dupicates information available from std::numeric_limits
const prime_t max_prime_t =
     prime_t(0)-prime_t(1) > prime_t(0) ?
         ~prime_t(0) :
         ~prime_t(0) ^ (prime_t(1) << (sizeof(prime_t)*8-1));

int main() {
    std::vector<prime_t> primes = { 2, 3 };
    prime_t ilimit = 0;
    prime_t ilimit_limit = 0;

    static const prime_t max_nur = max_prime_t - 4;

    // "step" alternates between 2 and 4,
    // to skip values containing factors of 2 and/or 3
    for(prime_t nur=5, step=4; nur<=max_nur; nur += (step ^= 6)) {
        if(nur >= ilimit_limit) {
            // Vedran Sego's optimization
            ilimit = static_cast<prime_t>(std::sqrt(nur));
            ilimit_limit = (ilimit + 1) * (ilimit + 1);  
        }
        for(auto one_prime : primes) {
            if(one_prime > ilimit) {
                std::cout << nur << '\n';
                primes.push_back(nur);
                break;
            }
            if(nur % one_prime == 0) {
                break;
            }
        }
    }
}

Context

StackExchange Code Review Q#131904, answer score: 4

Revisions (0)

No revisions yet.