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

Sieve of Sundaram

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

Problem

Every integer can be written in the form \$2^i*v\$ where \$v\$ is odd. This algorithm generates a list of numbers that are only odd. We want to ensure that all the composite odd numbers are marked out.

In the final list, all numbers will be of the form \$2^0 v\$. So, we only need to worry about odd factors. An odd number will be composite if it is of the form \$(2i + 1)(2j + 1)\$, where \$i\$ and \$j\$ are strictly positive integers. So, a composite odd number \$q = 2k + 1 = 4ij + 2i + 2j + 1\$. This means \$k = i + j + 2ij\$. Whenever \$k\$ is of the form, \$i + j + 2ij, 2k+ 1\$ is composite and whenever an odd integer \$q\$ is composite, the condition holds for \$k, q = 2k + 1\$.

All we seek to do is mark all numbers that can be written as \$i +j + 2ij\$. And for any unmarked integer \$k\$, we know that \$2k + 1\$ will be prime. Note - We need to hardcode 2 since this algorithm only generates the odd primes.

Please provide feedback.

```
#include
#include
#define target 50000

//Function prototypes
void initial_marking(long int[]);
void marking(long int[]);
void make_prime_list(long int[],long int[]);
void print_prime(long int[]);

int main()
{
//Half of the numbers from 1 to target will be even so the number of primes will never exceed target/2.
//The algorithm only worries about the first half of the range. That's why list also goes till target/2.
long list[target/2], primes[target/2];

//Sieve of Sundaram
initial_marking(list);
//Crossing out all numbers of the form i + j + 2ij
marking(list);

//Making the list of primes
make_prime_list(primes, list);

//Displaying the list
print_primes(primes);

return 0;
}

//Initially, everything is marked 1
void initial_marking(long int list[])
{
int i;

for(i = 1; i <= target/2; i++)
{
list[i] = 1;
}
}

//Crossing out numbers of the form i + j + 2ij
void marking(long int list[])
{
int i, j, crossed_out_num;

//All numbers of the form i +

Solution

Let \$f(i,j) = i + j + 2ij\$

Then \$f(i,j+1) = i + j + 1 + 2i(j+1) = i+j+2ij+2i+1 = f(i,j)+2i+1\$

Using this, the marking() function (which is the highest complexity function) can be rewritten as such:

void marking(long int list[])
{
    long long int i, crossed_out_num, i_limit, cross_limit, increment;
    i_limit = target/6; // i > i_limit will never give f(i,j) in range

    //All numbers of the form i + j + 2ij need to be marked out
    for(i = 1; i  2*i*(i+1)) cross_limit = 2*i*(i+1);

        // Find crossed_out_num directly by starting with j=1 and using the above result
        for(crossed_out_num = 3*i + 1; crossed_out_num <= cross_limit; crossed_out_num += increment)
        {
            list[crossed_out_num] = 0;
        }
    }
}


This should speed up the code considerably and avoid overflow when i becomes large.

Code Snippets

void marking(long int list[])
{
    long long int i, crossed_out_num, i_limit, cross_limit, increment;
    i_limit = target/6; // i > i_limit will never give f(i,j) in range

    //All numbers of the form i + j + 2ij need to be marked out
    for(i = 1; i <= i_limit; i++)
    {
        increment = 2*i + 1;

        // f(i,j) <= target/2 and j <= i
        cross_limit = target/2;
        if(cross_limit > 2*i*(i+1)) cross_limit = 2*i*(i+1);

        // Find crossed_out_num directly by starting with j=1 and using the above result
        for(crossed_out_num = 3*i + 1; crossed_out_num <= cross_limit; crossed_out_num += increment)
        {
            list[crossed_out_num] = 0;
        }
    }
}

Context

StackExchange Code Review Q#154328, answer score: 4

Revisions (0)

No revisions yet.