findprimes.c

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <stdint.h>
#include <string.h>

#define ARGLEN  10
#define PRIMESIZE sizeof(uint64_t)

static uint64_t * primes_found = NULL;
static uint64_t   prime_count  = 0;

/*
 *   findPrimes(maxPrime, maxCount)
 *
 *  Finds all primes up to either maxPrime, or finds as many primes as
 *  maxCount, whichever is reached first. 0 in one means "not a factor".
 *  0 in both is disallowed, and returns -1.
 *
 *  Upon normal completion, it will return the amount of primes that
 *  have been found in all executions.
 *
 */

uint64_t findPrimes(uint64_t maxPrime, uint64_t maxCount)
{
    if (maxPrime == 0 && maxCount == 0) { return -1; }

    char primeCond = maxPrime != 0;
    char countCond = maxCount != 0;

    if (primes_found == NULL)
    {
        primes_found = malloc(PRIMESIZE);
        primes_found[0] = 2;
        prime_count     = 1;
    }

    // If we've already calculated past maxPrime, we're done here

    uint64_t prime = primes_found[prime_count-1] + 1;

    long double maxTestNum_f;
    uint32_t    maxTestNum_i;

    uint64_t    curIndex;
    uint64_t    curPrime;
    uint8_t     isPrime;

    long double divResult_f;
    uint64_t    divResult_i;

    while (1)
    {
        if (primeCond && maxPrime <= prime)       { break; }
        if (countCond && maxCount <= prime_count) { break; }

        // No number has a prime factor larger than the square root of it
        maxTestNum_f = sqrtl((long double)prime);
        maxTestNum_i = (uint32_t)maxTestNum_f;

        isPrime = 1;

        // Test against all previous primes

        for (curIndex = 0; curIndex < prime_count; curIndex++)
        {
            curPrime = primes_found[curIndex];
            if (curPrime > maxTestNum_i) { break; }

            divResult_f = (long double)prime / (long double)curPrime;
            divResult_i = prime / curPrime;

            if (divResult_f == (long double)divResult_i)
            {
                isPrime = 0;
                break;
            }
        }

        if (isPrime)
        {
            // Extend by one prime
            primes_found = realloc(primes_found, PRIMESIZE * (++prime_count));
            primes_found[prime_count - 1] = prime;
        }

        prime++;
    }

    return prime_count;
}

/*   nthPrime(index)
 *
 *  Finds the nth prime, assuming it's been calculated. If it hasn't
 *  been calculated, return 0.
 *
 *  For reference: nthPrime(1) = 2, nthPrime(2) = 3, nthPrime(3) = 5
 *
 */

uint64_t nthPrime(uint64_t index)
{
    if (--index >= prime_count) { return 0; }
    return primes_found[index];
}


int main(void)
{
    uint64_t gotoPrime = 10001;

    printf("Found %lu primes\n", findPrimes(0, gotoPrime));
    printf("Prime %lu is %lu\n", gotoPrime, nthPrime(gotoPrime));

    return 0;
}