fermat_little_theorem_test.py

import numpy as np

def get_primes(n):
    """Returns  a list of primes < n
    https://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n
    """
    sieve = [True] * n
    for i in range(3, int(n ** 0.5) + 1, 2):
        if sieve[i]:
            sieve[i * i :: 2 * i] = [False] * ((n - i * i - 1) // (2 * i) + 1)
    return [2] + [i for i in range(3, n, 2) if sieve[i]]

def test_property(candidate: int, max_progression: int=100, sample_size: int=100):
    """Returns True if Fermat's little theorem's property is verified for a
    limited number of progressions

    Params:
    - candidate: the value whose power - 1 may be divisible by itself
    - max_progression: the highest progression ratio that will be tested
    - sample_size: if max_progression is high, a fixed amount of random
      progressions to test the property on
    """
    test_progressions = range(2, max_progression)
    if sample_size < len(test_progressions):
        test_progressions = np.random.choice(test_progressions, sample_size)

    for progression in range(2, max_progression):
        if progression % candidate == 0:
            # if the progression's ratio is proportional to the power, move on
            continue
        elif (progression ** (candidate - 1)) % candidate != 1:
            # the property isn't verified if any of the n-th-1 power mod p is
            # a number other than 1
            return False

    return True

def test_first_n_integers(n):
    """Test the first n integers and record if they were prime."""
    primes_to_n = get_primes(n + 1)

    # Variable to store stats
    count_primes = 0
    count_passes_and_prime = 0
    count_passes_and_not_prime = 0
    count_prime_but_does_not_pass = 0

    # We iterate over integers starting from 2
    for i in range(2, n + 1):
        print('Testing {}'.format(i))
        if i in primes_to_n:
            count_primes += 1
        if test_property(i, max_progression=i, sample_size=500):
            if i in primes_to_n:
                count_passes_and_prime += 1
            else:
                count_passes_and_not_prime += 1
        else:
            if i in primes_to_n:
                count_prime_but_does_not_pass += 1

    return count_primes, count_passes_and_prime, \
           count_passes_and_not_prime, count_prime_but_does_not_pass


if __name__ == '__main__':
    N = 100000
    a, b, c, d = test_first_n_integers(N)
    print(('{:,.0f} integers tested. {:,.0f} primes detected. {:,.0f} integers '
           'verified Fermat\'s property, of which {:,.0f} are primes. {:,.0f} '
           'primes did not verify the property').format(
               N, a, b + c, b, d))