Next level Improved Sieve of Eratosthenes

Question

This is the next development for the Sieve of Eratosthenes algorithm, where all the multiples of 2,3 and 5 are eliminated, which will not take any memory or time.

As my last question Improved Sieve of Eratosthenes I will follow your instructions and study it deeply as soon as I finish my exams.

I wanted to know if there are any specific improvements for this one!

This algorithm generates prime numbers up to a given limit N where it skips the first 10 primes when N>530.

N.B: I can make an infinite number of these developments, as I go higher it gets more and more complicated.

#include <iostream>

int D3SOE(unsigned int n1_m) {
    unsigned int n1 = 0, g = 0, z = 0, f1 = 0, f2 = 0, f3 = 0, f4 = 0, f5 = 0, f6 = 0, f7 = 0, f8 = 0, f9 = 0, l = 0, f10 = 0;
    int p = 0;

    unsigned char* prmLst = new unsigned char[n1_m + 1];

    // Initialising the D3SOE array with false values 
    for (int i = 0; i < n1_m; i++)
        prmLst[i] = false;

    // The main elimination theorem 
    for (n1 = 1; n1 <= ceil((8 * (ceil((9 * ceil((2 * (sqrt(((2 * (floor((3 * floor((10 * floor((9 * (n1_m)+1) / 8.0) + 8) / 9.0) + 1) / 2.0)) + 1) / 4.0)) - 0.5) - 1) / 3) - 8) / 10.0)) - 1) / 9); n1++)
    {
        if (prmLst[n1] != 0)
            continue;
        z = floor(((3 * floor((10 * floor((9 * (n1)+1) / 8.0) + 8) / 9.0) + 1) / 2.0));
        p = ((2 * floor(((3 * floor((10 * floor((9 * (n1)+1) / 8.0) + 8) / 9.0) + 1) / 2.0))) + 1);
        f1 = ceil((8 * (ceil((9 * ceil((p - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f2 = ceil((8 * (ceil((9 * ceil(((7 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f3 = ceil((8 * (ceil((9 * ceil(((11 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f4 = ceil((8 * (ceil((9 * ceil(((13 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f5 = ceil((8 * (ceil((9 * ceil(((17 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f6 = ceil((8 * (ceil((9 * ceil(((19 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f7 = ceil((8 * (ceil((9 * ceil(((23 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f8 = ceil((8 * (ceil((9 * ceil(((29 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
        f9 = ceil((8 * (ceil((9 * ceil(((31 * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9.0);

        l = ceil((8 * (ceil((9 * ceil(((4 * ((z * z) + z)) - 1) / 3.0) - 8) / 10.0)) - 1) / 9);
        
        for (g = l; g < n1_m; g += (f9 - f1))
        {
            prmLst[g] = true;
        }
        if ((p - 1) % 30 == 0)
        {
            for (g = l + (f2 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f3 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f4 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f5 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f6 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f7 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 7) % 30 == 0)
        {
            for (g = l - (f2 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f3 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f4 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f5 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f6 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f7 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 11) % 30 == 0)
        {
            for (g = l - (f3 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f3 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f4 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f5 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f6 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f7 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 13) % 30 == 0)
        {
            for (g = l - (f4 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f4 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f4 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f5 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f6 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f7 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 17) % 30 == 0)
        {
            for (g = l - (f5 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f5 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f5 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f5 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f6 - f5); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f7 - f5); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f5); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 19) % 30 == 0)
        {
            for (g = l - (f6 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f6 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f6 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f6 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f6 - f5); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f7 - f6); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f6); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 23) % 30 == 0)
        {
            for (g = l - (f7 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f7 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f7 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f7 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f7 - f5); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f7 - f6); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l + (f8 - f7); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
        else if ((p - 29) % 30 == 0)
        {
            for (g = l - (f8 - f1); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f8 - f2); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f8 - f3); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f8 - f4); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f8 - f5); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f8 - f6); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
            for (g = l - (f8 - f7); g < n1_m; g += (f9 - f1))
            {
                prmLst[g] = true;
            }
        }
    }

    // printing for loop
    for (int n1 = 1; n1 < n1_m; n1++)
        if (!prmLst[n1]) {
            z = (2 * floor(((3 * floor((10 * floor((9 * (n1)+1) / 8.0) + 8) / 9.0) + 1) / 2.0)) + 1);
            std :: cout << z << " ";
        }
    return 0;
}

// derivation function 
int main() {
    unsigned int n1_m, N = 0;
    std::cout << "\n Enter limit : ";
    std::cin >> N;
    n1_m = ceil((8 * (ceil((9 * ceil(((N)-2) / 3.0) - 8) / 10.0)) - 1) / 9);
    D3SOE(n1_m);
    return 0;
}

Answer 1

Naming things

Your code is hard to read, in a large part because your function and variable names don't convey any meaning. What does D3SOE mean? Why is the parameter of that function named n1_m? Try to use concise but meaningful names.

The function should be named after what it does. How it does it is less important. So instead of sieve_of_eratosthenes(), call it print_primes(). The parameter is the highest number to check, so call it highest_number, although max_n or n_max are fine as well and more concise, as n is a very commonly used abbreviation for number.

Don't abbreviate unnecessarily. Instead of prmLst, write primeList or primes.

There are lots of complex-looking expressions, like ceil((8 * (ceil((9 * ceil(((N)-2) / 3.0) - 8) / 10.0)) - 1) / 9). You can give those expressions a name by creating a new function. This makes the code more readable and more self-documenting.

Don't repeat yourself

There is a lot of repetition in your code. Apart from being more work to write, it also increases the likelihood of errors being introduced. Try to find some ways to avoid code duplication.

Consider the variables f1...f9 (f10 is declared but never used). You could make an array f[] instead, and intialize them like so:

static const int multipliers[9] = {1, 7, 11, 13, 17, 19, 23, 29, 31};
...
int f[9];
for (std::size_t i = 0; i < 9; ++i)
{
    f[i] = ceil((8 * (ceil((9 * ceil(((multipliers[i] * p) - 2) / 3.0) - 8) / 10.0)) - 1) / 9);
}

Once you have an array of fs, you can also get rid of the duplicated for-loops, for example:

if ((p - 1) % 30 == 0)
{
    for (std::size_t i = 1; i < 8; ++i)
    {
        for (g = l + (f[i] - f[i - 1]); g < n1_m; g += f[8] - f[0])
        {
            prmLst[g] = true;
        }
    }
}

But even the chain of if-else-statements can be avoided, by adding another loop that checks if (p - multipliers[i]) % 30 == 0, and calculates which values in f the inner loop should subtract from each other.

The above code still hardcodes the number of fs to 9, it would be even better to make it so generic that if you wanted to add another multiplier to check, you only would have to add another prime to multipliers[], and the rest of the code can be left unchanged.

Memory leak

You allocate memory for prmLst[], but you never delete[] it. This means your program has a memory leak. Maybe it's not important here, but in larger projects these things quickly become problematic. Either make sure you clean up memory properly, or avoid having to do manual memory allocations in the first place by using STL containers, such as std::vector:

std::vector<bool> prmLst(n1_m + 1);

Missing #includes and namespace issues

If you call math functions like ceil() and sqrt(), you should #include <cmath>. In C++, the functions in the standard library are all declared inside the namespace std, and there is no guarantee that the math functions will also be available in the global namespace. So you should write std::ceil(), std::sqrt() and so on.

Return value

Since D3SOE() doesn't return any useful value, just make it return void. Alternatively, instead of having D3SOE() both generate and print the primes, have it only generate a list of primes, and return it. This allows the caller to decide what to do with those primes. This reduces the responsibilities of D3SOE(), thus simplifying that function, while at the same time making it more versatile.

Answer 2

This is similar to a wheel sieve with basis {2,3,5}.

You can simplify everything and write it compactly like this:

#include <iostream>
#include <cmath>
#include <algorithm>
#include <vector>
#include <cstdlib>
#include <stdint.h>


const  int64_t C1[64]={-22,  -32,  -28,  -32,  -28,  -8,   -8,   -22,
                       -30,  -18,  -30,  -30,  -12,  -30,  -12,  -18,
                       -34,  -26,  -16,  -14,  -34,  -26,  -14,  -16,
                       -42,  -42,  -18,  -12,  -18,  -18,  -18,  -12,
                       -46,  -20,  -34,  -26,  -10,  -14,  -20,  -10,
                       -24,  -36,  -36,  -24,  -24,  -6,   -24,  -6,
                       -36,  -30,  -24,  -36,  -30,  -24,   0,    0,
                       -40,  -38,  -40,  -20,  -22,  -20,  -2,    2};

const int64_t C2[64]={ 0,  9,  7,  9,  7, 1, 1, 0,
                       6,  0,  8,  8,  1, 6, 1, 0,
                       9,  5,  0,  1,  9, 5, 1, 0,
                      15, 15,  1,  0,  3, 3, 1, 0,
                      18,  1, 10,  6,  0, 2, 1, 0,
                       1, 11, 11,  5,  5, 0, 1, 0,
                      10,  7,  4, 10,  7, 4, 0, 0,
                      13, 12, 13,  3,  4, 3, 0, 0};

const int64_t RW[8]={-23, -19, -17, -13,  -11,-7 , -1, 1};

const int64_t bW=30;

const int64_t nR=8;

void Sieve_Wheel_30(int64_t  k_end, std::vector<char> &Primes)
{
    
    int64_t  m,mmin,i,j;
    int64_t  kmax=(int64_t)std::sqrt(k_end/bW)+2;
    for (int64_t  k = 1; k  <= kmax; k++)
    {
        for (j = 0; j  < nR; j++)
        {
            if(Primes[k+j*(k_end+1)])
            {
                for (i = 0; i  < nR; i++)
                {
                    // Cancelling out all the multiples of RW[j]+bW*k   
                    mmin=bW*k*k + k*C1[i*nR+j] + C2[i*nR+j];
                    for ( m =mmin; m <= k_end; m += bW*k+RW[j])
                        Primes[m+i*(k_end+1)] = false;
                }
            }
        }
    }
}

int main()
{
    int64_t N = 0;
    std::cout << "\n Enter limit : ";
    std::cin >> N;

    int64_t  k_end=(N < bW) ? 2 : N/bW+1;

    std::vector<char> Primes(nR*(k_end + 1), true);
    Sieve_Wheel_30(k_end, Primes);

    std::cout << 2 << " " << 3 <<" " << 5 <<" ";
    for (int64_t k = 1; k < k_end; k++)
        for (int64_t i = 0; i  < nR; i++)
            if(Primes[k+i*(k_end+1)])
                std::cout << RW[i]+k*bW << " ";
    for (int64_t i = 0; i<nR; i++)
        if(Primes[k_end+i*(k_end+1)] and RW[i]<(int64_t)(N%bW-bW))
            std::cout << RW[i]+k_end*bW  << " ";

    return 0;
}

EDIT

As you wrote, the strengths of this sieve wheel is the use of a smaller memory equal to phi(bW)*N/bW with bW the wheel modulus and phi() is Euler's totient function.

This makes it possible to create a segmented version using a bit space phi(bW)*sqrt(N)/bW and you can also add the use of numbers smaller than N/bW.

I have not done a benchmark but this is the segmented version of the same sieve with the possible choice of increasing the value of the wheel modulus.

The sieves I have seen use sqrt(N) as memory for segmentation instead this wheel sieve uses the product of the prime numbers following the basis so that a presieving can be done. This way the memory used is slightly higher than phi(bW)*sqrt(N)/bW but is always less than sqrt(N).

I think you can also create m different blocks for multithreading in the second part of the code, using m different Segment_t vectors and increasing for each block k_low by m*segment_size.