Compares the ratios of ending digits of number to next number, for random sets of numbers, in comparison with the primes. See https://www.quantamagazine.org/20160313-mathematicians-discover-prime-conspiracy/

prime_property_sim

from random import random
from math import log
from gmpy2 import next_prime

limit = 1000000

rands = []
primes = []

cur_prime = 11

for i in range(11, limit)[::2]:
    if random() < (1 / log(i) * 10 / 4) and i % 5 != 0: # 10/4 because we only look at 1, 3, 7, 9 out of 10 possible digits
        rands.append(i)
        primes.append(cur_prime)
        cur_prime = next_prime(cur_prime)

data_rands = {}
data_primes = {}

for i in range(len(rands) - 1):
    if i != 0:
        if (rands[i - 1] % 10, rands[i] % 10) in data_rands:
            data_rands[(rands[i - 1] % 10, rands[i] % 10)] += 1
        else:
            data_rands[(primes[i - 1] % 10, primes[i] % 10)] = 0
        if (primes[i - 1] % 10, primes[i] % 10) in data_primes:
            data_primes[(primes[i - 1] % 10, primes[i] % 10)] += 1
        else:
            data_primes[(primes[i - 1] % 10, primes[i] % 10)] = 0



for i in [1, 3, 7, 9]:
    total_rands = 0
    total_primes = 0
    print("cumulative:")
    for j in [1, 3, 7, 9]:
        print(str(i) + " to " + str(j) + ": " + str(data_rands[(i, j)]) + " rand, " + str(data_primes[(i, j)]) + " prime")
        total_rands += data_rands[(i, j)]
        total_primes += data_primes[(i, j)]
    print("ratios:")
    for j in [1, 3, 7, 9]:
        print(str(i) + " to " + str(j) + ": " + str(data_rands[(i, j)] / total_rands) + " rand, " + str(data_primes[(i, j)] / total_primes) + " prime")
    print("")