For small primes - i.e. up to 2^32 and maybe up to 2^40..2^50ish regions - you can improve speed by several orders of magnitude if you use a Sieve of Eratosthenes.
In order to get decent speed the factor primes - i.e. those up to
sqrt(limit), whatever your limit might be - should be sieved up front and stored in an array, alongside with the current working offset for each prime (initialised to the square of the prime). The sieving of the small factor primes can also be done with a - simpler - Sieve of Eratosthenes.
The size of the sieve should be equal to the L1 cache size of the expected target machine (usually 32 KByte), and for maximum efficiency you should use a packed odds-only bitmap.
BitArray is slow as molasses, so either use
byte[1<<15] with manual bit indexing (i.e.
sieve[i>>3] |= 1<<(i&7) to set bit i). Using wheels like mod30 can improve efficiency further, but only in connection with major unrolling of loops and lots of other complications.
Another performance boost is offered by presieving. Instead of clearing the sieve prior to sieving a new window's worth of primes, you simply blast over the sieve a pre-computed bit pattern that corresponds to sieving by a certain number of (odd) small primes. Presieving up to 13 hits a sweet spot in a situation like this, but going up to 17 or 19 can realise further (tiny) speed gains at the cost of much increased size of the precomputed pattern.
In order to illustrate what the code would look like and how the various techniques connect, here's the actual sieving function from my SPOJ PRIME1 submission:
internal static byte presieve_pattern;
internal static int presieve_modulus;
internal static int first_sieve_prime_index;
internal static byte sieve;
internal static int result;
internal static int result_size;
internal static void sieve_primes_between_ (int m, int n)
if (m > n)
Trace.Assert(first_sieve_prime_index == 0 || m > small_odd_primes[first_sieve_prime_index - 1]);
Trace.Assert((m & 1) != 0 && (n & 1) != 0);
Trace.Assert(m > 1);
int sqrt_n = (int)Math.Sqrt(n);
int max_bit = (n - m) >> 1; // (the +1 is redundant here since both are odd)
Trace.Assert(max_bit < sieve.Length * 8);
if (first_sieve_prime_index == 0) // without presieving
Array.Clear(sieve, 0, (max_bit >> 3) + 1);
int base_bit = (m >> 1) % presieve_modulus;
while ((base_bit & 7) != 0)
base_bit += presieve_modulus;
Buffer.BlockCopy(presieve_pattern, base_bit >> 3, sieve, 0, (max_bit >> 3) + 1);
for (int i = first_sieve_prime_index, e = small_odd_primes.Length; i < e; ++i)
int prime = small_odd_primes[i];
if (prime > sqrt_n)
int start = prime * prime, stride = prime << 1;
if (start < m)
start = (stride - 1) - (m - start - 1) % stride;
start -= m;
for (int j = start >> 1; j <= max_bit; j += prime)
sieve[j >> 3] |= (byte)(1 << (j & 7));
Trace.Assert(m > 3 && m % 3 != 0);
for (int i = 0, d = 3 - m % 3, e = (n - m) >> 1; i <= e; i += d, d ^= 3)
if ((sieve[i >> 3] & (1 << (i & 7))) == 0)
result[result_size++] = m + (i << 1); // m is odd already
The code should answer quite a few questions that I simply glossed over in the text. This version is geared towards experimentation with various aspects (like number of presieve primes and so on), and thus it lacks a lot of - uglifying and code-bloating - optimisations.
The wrapper function needs to take care of things like non-sieve primes and normalising arguments, which is why I'm showing it here as well:
internal static void sieve_primes_between (int m, int n)
result_size = 0;
if (m < small_odd_primes[first_sieve_prime_index])
for (int i = -1; i < first_sieve_prime_index; ++i)
int prime = i < 0 ? 2 : small_odd_primes[i];
if (m <= prime && prime <= n)
result[result_size++] = prime;
m = small_odd_primes[first_sieve_prime_index];
n -= 1 - (n & 1); // make it odd
m |= 1; // make it odd
m += m % 3 == 0 ? 2 : 0; // skip if multiple of 3
This illustrates the kind of complications that you have to deal with when certain optimisations are applied, like special handling for the wheel primes (in this case only the number 2) and the presieve primes.
What is not shown here - because it does not fit the profile of the PRIME1 task - is the remembering of offsets. I don't have working C# code that uses the technique at the moment, because I just started learning the language. Remembering the working offsets shouldn't be too difficult to implement, despite the peculiarities of C# like lack of reference types.
The code takes 3..4 milliseconds to sieve 10 x 100,000 numbers just below the PRIME1 limit (10^9). A version using
bool could be marginally faster for the low PRIME1 limit but certain further optimisations require the use of
int, especially as C# has major performance problems with indexing.
This type of sieve performs decently in the 32-bit range. It can be modified to operate beyond 32 bits but in higher regions of the 64-bit range it gets bogged down by the sheer number of potential factor primes (hundreds of millions, at the end of the range) each of which needs to be considered for sieving each window unless a technique like bucketing is used.