Work on Idea #164

A Script for Generating and Testing Prime Number Formulas

In Idea #164, we proposed generating expressions at random and testing if they were capable of rendering the prime numbers. Here I present one of such scripts made up in Python. I also produced a handy executable. You can download them at the bottom of the page.

The algebraic expressions it generates are made up of integer polynomials terms of n, that are joined by operations such as sum, subtraction, multiplication or division. Brackets, even nested brackets, all generated randomly, add to the possibilities. Positive and negative signs are also randomly selected.

Naturally, testing all the prime numbers would be unfeasible, so the program tests only the first 9592 prime numbers. It is designed to stop testing an expression if that has already failed to generate all the primes, and move on to the next, what actually makes testing expressions a very fast process. 

Here are the standard running parameters of the program, but it is all up for customization on the Python script:

Note that the possibilities are vast. I can't even begin to calculate the odds. It could find the right expression in the next 5 minutes, or it may never find it, and there is no way of telling. That is why this is no different than a lottery. Except perhaps that it is lottery for the advancement of Science.  The idea is just that, you sit to work and meanwhile, you put the program to run in the background. At least it feels nice to know that the computer is working for you trying to solve this ancient puzzle. 

Every 5000 expressions tested it bumps the screen just to let you know the program is still alive. 

If it does find an expression, it will write "Eureka" on the screen and write the expression, as well as registering on the output file Prime_Function.txt for further investigation. That means that from time to time, check that file to see if anything came up. An expression registered on the file means it is a Prime Function candidate, it does not guarantee that it can generate all the primes, merely the first 9592, that is already on itself quite an achievement. Further testing is necessary for a larger set of primes in such case.

I certainly hope that the registering on the file is working, I never did find a function, for obvious reasons, to test that part of the code. The program relies on Eval() which some people say it is untrustworthy feature of Python, but honestly I have not had many issues when it is carefully structured in a try-except block. 

The prime number table, primes.txt, must be always on the same place as the script ou executável or else an error will be thrown off. 

Performance Issues

Upon running some very simple tests on the performance of the script, I determined that running it does not affect significantly the performance of the computer and it can be safely used in the background. You may notice performance issues if you are running activities that really require a lot of CPU Power, such as video editing, rendering images, data processing, but that would be true for any kind of multi-tasking. On my computer, an old Intel Core i7-4710HQ CPU @ 2,50GHz, the program alone consumed 15% of total CPU Power, with every day tasks not raising things beyond 65% of CPU power.

The Complete Code

# This script generates a randomized algebraic expression and evaluates it, checking if

# the first 9592 prime numbers were produced by the expression.

# The expression is rational, formed by integer polynomials, with random brackets.

# If the expression did indeed produced the prime numbers, it is registered for further

# studies in the output file Prime_Function.txt.

# Note that being able to reproduce the first 9592 prime numbers makes it a Prime Number

# Function candidate, nothing more. There are plenty other prime numbers out there.

# Every 5000 attempts it outputs in the screen the current screen, otherwise it is a very

# silent program, designed to be run in the background. It also indicates in the screen

# if it by chance ran into a Prime Number expression.

# Note that the possibities are vast. It could find an expression in  the next 5 min,

# it could take the age of the Universe, or it even may never find it because it is

# impossible.

import math

import random

# Execution Parameters

max_number_of_terms = +25

max_numeric_term = +1000000000

min_literal_exponent = -50

max_literal_exponent = +50

max_literal_coefficient = +1000000000

number_of_attempts = +1000000000


# Seeds the pseudo-random number generator:


# This function sorts out if the term is going to be numeric (a constant) or literal (a polynomial of n).

def integer_or_numeric_term():

    # type_term = 0 --> Numeric

    # type_term = 1 --> Literal

    random_float = random.random()

    if random_float <= 0.5:

        type_term = 0

    if random_float > 0.5:

        type_term = 1

    return type_term

# This function sorts if a sign will be positive or negative.

def sort_sign():

    random_float = random.random()

    if random_float <= 0.5:

        sign_str = '-'

    if random_float > 0.5:

        sign_str = '+'

    return sign_str

# This function generates a randomized expression that could be the Prime Function Generator:

def genExpression():

    number_of_terms = random.randint(2, max_number_of_terms)

    collected_terms = []

    for item in range(number_of_terms):

        type_term = integer_or_numeric_term()

        # Numeric Term

        if type_term == 0:

            sign_str = sort_sign()

            coefficient = random.randint(0, max_numeric_term+1)

            coefficient_str = '{:20d}'.format(coefficient)

            #coefficient_str = '{:d}'.format(coefficient)

            term_str = '(' + sign_str + coefficient_str + '*n**( +0))'

        # Literal Term

        if type_term == 1:

            sign_str = sort_sign()

            coefficient = random.randint(0, max_literal_coefficient+1)

            coefficient_str = '{:20d}'.format(coefficient)

            #coefficient_str = '{:d}'.format(coefficient)

            exponent = random.randint(min_literal_exponent, max_literal_exponent+1)

            exponent_str = '{:+3d}'.format(exponent)

            term_str = '(' + sign_str + coefficient_str + '*n**(' + exponent_str + '))'


    number_of_operations = number_of_terms-1

    # Collects the operations (+,-,*,/) betweeen terms, all generated randomly.

    operations = []

    for operation in range(number_of_operations):

        random_float = random.random()

        if random_float < 0.25:

            operation_str = '+'

        if random_float >= 0.25 and random_float < 0.5:

            operation_str = '-'

        if random_float >= 0.5 and random_float < 0.75:

            operation_str = '*'

        if random_float > 0.75:

            operation_str = '/'


    # Assembles the expression:

    expression = ''

    brackets_opened = 0

    for item in range(number_of_operations):

        open_bracket = ''

        closed_bracket = ''

        # The possibility of a bracket to open:

        random_float = random.random()

        if random_float <= (0.5):

            open_bracket = '('

            brackets_opened = brackets_opened + 1


            bracket = ''

        # The possibility of a bracket to close:

        random_float = random.random()

        if random_float <= (0.5) and (brackets_opened > 0):

            closed_bracket = ')'

            brackets_opened = brackets_opened - 1


            bracket = ''

        expression = expression + open_bracket + collected_terms[item] + closed_bracket + operations[item]

    # Adds the final term:

    expression = expression + collected_terms[item+1]

    # Closes any un-closed bracket:

    while brackets_opened > 0:

        expression = expression + ')'

        brackets_opened = brackets_opened - 1

    return expression

# Evaluate expression:

# Opens output file:

resultsFile = open('Prime_Function.txt', 'a')

# Loads up prime numbers table:

with open('primes.txt', 'r') as primesFile:

    content = primesFile.readlines()

# Clears input of garbage:

primes = []

for item in content:

    prime = item.strip('\n')


# Attempts to find an expression for the prime numbers:

for run in range(number_of_attempts):

    expression = genExpression()

    #print(run ,  ':',  expression)

    if (run % 5000 == 0) : print("RUN: ", run)

    # Tries, for every prime, to see if there is a match, and stops trying that expression in negative case:

    is_possible_p_function = True

    for integer in range(len(primes)):

        n = integer+1

        if is_possible_p_function == True:


                expression_value = eval(expression)

                if expression_value != primes[integer]:

                    is_possible_p_function = False

                    #print('Not prime function.')


                #print('Error: EVAL() failed.')

                is_possible_p_function = False

                #print('Not prime function')

    if is_possible_p_function == True:



        #Writes to file:




print("END PROGRAM..")

# Note:

# The genExpression() function sometimes generate expressions with 0 on the denominator.

# This is catched by the try-except block of Eval. They cannot be evaluated.

# It was easier to simply ignore the problem.

BANNER IMAGE CREDITS: ESA/Hubble & NASA, A. Filippenko, R. Jansen 

Want to know more about this image? Follow this external link.