Generating prime spirals as described on numberphile https://www.youtube.com/watch?v=iFuR97YcSLM

prime_spirals.py

from matplotlib import pyplot as plt
from functools import lru_cache

SIEVE = 1000

@lru_cache()
def sieve(limit):
    factors = [2] * limit
    for base in range(2, int(limit**0.5 + 1)):
        # if factors[base] == 2:
        for a in range(base+base, limit, base):
            factors[a] += 1
    factors[0] = factors[1] = 2
    return factors

def is_prime(n):
    return sieve(SIEVE)[n] == 2

def num_factors(f):
    return sieve(SIEVE)[f]


def grid(n):
    g = []
    for i in range(n):
        l = []
        for j in range(n):
            l.append(0)
        g.append(l)

    m = (n-1)//2
    c = 1
    g[m][m] = 1

    for i in range(m-1, -1, -1):
        for j in range(i+1, n-i):
            c += 1
            g[i][j] = c
        for j in range(i+1, n-i):
            c += 1
            g[j][n-i-1] = c
        for j in range(i+1, n-i):
            c += 1
            g[n-i-1][n-j-1] = c
        for j in range(i+1, n-i):
            c += 1
            g[n-j-1][i] = c

    return g

def print_grid(g):
    for l in g:
        for i in l:
            print("%3d" % i, end=' ')
        print()

def plot_grid(g):
    for l in g:
        for i in range(len(l)):
            # l[i] = 2.0 / num_factors(l[i]) 
            l[i] = 1.0 if is_prime(l[i]) else 0.0
    plt.imshow(g, cmap='plasma')
    plt.show()

def main():
    # print_grid(grid(31))
    # plot_grid(grid(501))
    # for i in range(100):
    #     print(i, num_factors(i))
    global SIEVE
    SIEVE = 200*200
    plot_grid(grid(199))

if __name__ == "__main__":
    main()