Infinite Sequence Part IV - Prime

Countinue from last post, Fibonacci, which is considerable a simple sequence, because it uses only 2 preceding numbers to generate next number. This post, we move on to harder sequence which next number has to generate from all preceding numbers. The sequence is Prime.

Prime is natural number which there is only 2 natural number divisors. The numbers are 1 and itself. Prime below 100 are 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97.

The technique to store prime that I used is Sieve of Eratosthenes. I had optimized it by skip all even numbers. Memory usage for storing prime is adaptable. Initial capacity is 1,024 odd numbers and it can be extended to 1,073,741,824 odd numbers which cover all interger values.

Here is the code for class to store prime:

public class PrimeStore {
    const int defaultCapacity = 1024;
    BitArray _store;
    int _lastIndex;

    public PrimeStore() {
        _store = new BitArray(defaultCapacity, true);
    }
    public PrimeStore(int capacity) {
        _store = new BitArray(capacity, true);
    }

    public bool IsPrime(int num) {
        if (num < 2)
            return false;
        if (num == 2)
            return true;
        if ((num & 1) == 0)
            return false;

        int index = (num - 3) >> 1;
        testPrime(index);
        return _store.Get(index);
    }
    void testPrime(int index) {
        if (index <= _lastIndex)
            return;

        ensureCapacity(index);

        for (int prime; _lastIndex < index; _lastIndex++) {
            if (_store.Get(_lastIndex)) {
                prime = (_lastIndex << 1) + 3;
                applyPrime(prime, _lastIndex);
            }
        }
    }
    void ensureCapacity(int index) {
        int lastIndex = _store.Count - 1;
        if (index <= lastIndex)
            return;

        int newCapacity = _store.Count;
        while (index >= newCapacity)
            newCapacity = newCapacity * 2;

        int newSize = newCapacity / 32;
        int[] array = new int[newSize];
        _store.CopyTo(array, 0);
        for (int i = _store.Count / 32; i < newSize; i++)
            array[i] = -1;
        _store = new BitArray(array);

        for (int i = 0, prime, num; i < _lastIndex; i++) {
            if (_store.Get(i)) {
                prime = (i << 1) + 3;
                num = lastIndex - i;
                applyPrime(prime, num - (num % prime) + i);
            }
        }
    }
    void applyPrime(int prime, int index) {
        index += prime;
        int capacity = _store.Count;
        for (; index < capacity && index > 0; index += prime)
            _store.Set(index, false);
    }
}

Usage of PrimeStore class is easy. Just pass a number to IsPrime method.

bool isPrime = new PrimeStore().IsPrime(34567);

But our objective is not just for test a prime number. We will generate prime sequence. Our prime sequence will use PrimeStore to store prime numbers. Any prime numbers beyond PrimeStore's capacity will be calculated by testing there is no prime divisor for that numbers.

To limit growth of PrimeStore, I also defined MaxStoredValue. Default value is 65,536. Generating prime below MaxStoredValue is fast. Generating prime beyond MaxStoredValue but below MaxStoredValue^2 is acceptable speed. Generating prime beyond MaxStoredValue^2 is considerable slow. Default value is suited with general prime generation.

Here is the prime sequence code:

public class Prime : BaseSequence<long> {
    PrimeStore _primeStore;
    const long maxStoredValue = 2147483647L;

    public Prime() {
        _primeStore = new PrimeStore();
    }
    public Prime(int capacity) {
        _primeStore = new PrimeStore(capacity);
    }

    protected override IEnumerable<long> Enumerate() {
        yield return 2L;

        long value = 3L;
        yield return value;

        while (true) {
            do
                value += 2L;
            while (!Contains(value));
            yield return value;
        }
    }
    public bool Contains(long value) {
        if (value < maxStoredValue)
            return _primeStore.IsPrime((int)value);

        return Factors(value).First() == value;
    }
    public IEnumerable<long> Factors(long value) {
        if (value < 2L)
            yield break;

        using (IEnumerator<long> enumerator = Enumerate().GetEnumerator()) {
            enumerator.MoveNext();

            long factor = enumerator.Current;
            long pow = factor * factor;
            long num = value;
            long mod, num2;

            while (num >= pow) {
                num2 = Math.DivRem(num, factor, out mod);
                if (mod == 0L) {
                    num = num2;
                    yield return factor;
                } else {
                    enumerator.MoveNext();
                    factor = enumerator.Current;
                    pow = factor * factor;
                }
            }
            yield return num;
        }
    }
}

Other methods beside generating prime sequence are Contains method and Factors method. Contains method is to test whether a number is prime. And Factors method is to enumerate prime factors for a number.

As usual, test prime sequence with projecteuler.net.
Find the largest prime factor of 317584931803.

Console.WriteLine(new Prime().Factors(317584931803).Max());
//3919

Find the 10001^st prime.

Console.WriteLine(new Prime().ElementAt(10000));
//(10000 is index of position 10001st)
//104743

Even creating a prime sequence has many details, but using of prime sequence is simple and fast.

Next we will create sequence which is not the number. Continue next post :)