Algorithm for finding the N first prime numbers - Atkin Sieve
In the process of implementing one program, I was faced with the task of finding primes up to the number N of order . On a habr already repeatedly wrote about various ways and methods, but did not mention about the main method of the decision that I will try to correct today.
The algorithm has asymptotic complexity and requires a bit of memory. At the input values of the order, the Atkin sieve is 9.2 times faster than the Erastofen sieve. Here is a graph of the growth of the superiority of the Atkin algorithm on numbers from 2 to :
As a result, you can observe the following execution speed:
For convenience and brevity, we will create an abstract class from which we will inherit the Atkin sieve implementations. In the future, we will only create various class variants.
Let's create a basic implementation of the algorithm:
The execution time of the basic (classic) version of the algorithm:
Recall that in C #, operations with the ulong type take about 3 times longer than with the int type . In the program, when calculating large numbers, the values still exceed int.MaxValue, so we will find an empirical usage bar. Here it is 858599509. Insert the condition and get the final acceleration:
That's all. Thank you for your attention - I tried not to be very dry.
The algorithm has asymptotic complexity and requires a bit of memory. At the input values of the order, the Atkin sieve is 9.2 times faster than the Erastofen sieve. Here is a graph of the growth of the superiority of the Atkin algorithm on numbers from 2 to :
As a result, you can observe the following execution speed:
10,000,000 | 0.15 sec |
100,000,000 | 2.16 sec |
1,000,000,000 | 48.76 sec |
For convenience and brevity, we will create an abstract class from which we will inherit the Atkin sieve implementations. In the future, we will only create various class variants.
public abstract class AAtkin
{
internal readonly int _limit;
public bool[] IsPrimes;
protected AAtkin(int limit)
{
_limit = limit;
FindPrimes();
}
public abstract void FindPrimes();
}
Let's create a basic implementation of the algorithm:
public override void FindPrimes()
{
IsPrimes = new bool[_limit + 1];
double sqrt = Math.Sqrt(_limit);
var limit = (ulong) _limit;
for (ulong x = 1; x <= sqrt; x++)
for (ulong y = 1; y <= sqrt; y++)
{
ulong x2 = x*x;
ulong y2 = y*y;
ulong n = 4*x2 + y2;
if (n <= limit && (n%12 == 1 || n%12 == 5))
IsPrimes[n] ^= true;
n -= x2;
if (n <= limit && n%12 == 7)
IsPrimes[n] ^= true;
n -= 2 * y2;
if (x > y && n <= limit && n%12 == 11)
IsPrimes[n] ^= true;
}
for (ulong n = 5; n <= sqrt; n += 2)
if (IsPrimes[n])
{
ulong s = n*n;
for (ulong k = s; k <= limit; k += s)
IsPrimes[k] = false;
}
IsPrimes[2] = true;
IsPrimes[3] = true;
}
The execution time of the basic (classic) version of the algorithm:
1,000,000 | 0.03 sec |
10,000,000 | 0.39 sec |
100,000,000 | 4.6 sec |
1,000,000,000 | 48.92 seconds |
Recall that in C #, operations with the ulong type take about 3 times longer than with the int type . In the program, when calculating large numbers, the values still exceed int.MaxValue, so we will find an empirical usage bar. Here it is 858599509. Insert the condition and get the final acceleration:
IsPrimes = new bool[_limit + 1];
double sqrt = Math.Sqrt(_limit);
if (sqrt < 29301)
{
for (int x = 1; x <= sqrt; x++)
for (int y = 1; y <= sqrt; y++)
{
int x2 = x * x;
int y2 = y * y;
int n = 4 * x2 + y2;
if (n <= _limit && (n % 12 == 1 || n % 12 == 5))
IsPrimes[n] ^= true;
n -= x2;
if (n <= _limit && n % 12 == 7)
IsPrimes[n] ^= true;
n -= 2 * y2;
if (x > y && n <= _limit && n % 12 == 11)
IsPrimes[n] ^= true;
}
for (int n = 5; n <= sqrt; n += 2)
if (IsPrimes[n])
{
int s = n * n;
for (int k = s; k <= _limit; k += s)
IsPrimes[k] = false;
}
}
else
{
var limit = (ulong) _limit;
for (ulong x = 1; x <= sqrt; x++)
for (ulong y = 1; y <= sqrt; y++)
{
ulong x2 = x*x;
ulong y2 = y*y;
ulong n = 4*x2 + y2;
if (n <= limit && (n%12 == 1 || n%12 == 5))
IsPrimes[n] ^= true;
n -= x2;
if (n <= limit && n%12 == 7)
IsPrimes[n] ^= true;
n -= 2 * y2;
if (x > y && n <= limit && n%12 == 11)
IsPrimes[n] ^= true;
}
for (ulong n = 5; n <= sqrt; n += 2)
if (IsPrimes[n])
{
ulong s = n*n;
for (ulong k = s; k <= limit; k += s)
IsPrimes[k] = false;
}
}
IsPrimes[2] = true;
IsPrimes[3] = true;
}
That's all. Thank you for your attention - I tried not to be very dry.