The Ulam Spiral

Posted on 02 September 2016

The Ulam Spiral is a simple way to visualize the prime numbers which reveals unexpected, irregular structure. Starting at the centre of a rectangular grid of numbers and spiralling outwards, mark the primes as 1 and the composite numbers as 0.

The following code generates an image of the Ulam Spiral. The hard part is generating a "spiralled"-version of a two-dimensional array, which is the job of the make_spiral. This function takes a two-dimensional array of indices into the array, idx, and successively appends the top row to a new list of indices, spiral_idx. After each append, the appended row is removed from idx and the remainder of the array rotated by $90^\circ$. The concatenated array of spiral_idx is then used to index into the array to be "spiralled".

This code is also available on my GitHub page.

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm

def make_spiral(arr):
    nrows, ncols= arr.shape
    idx = np.arange(nrows*ncols).reshape(nrows,ncols)[::-1]
    spiral_idx = []
    while idx.size:
        spiral_idx.append(idx[0])
        # Remove the first row (the one we've just appended to spiral).
        idx = idx[1:]
        # Rotate the rest of the array anticlockwise
        idx = idx.T[::-1]
    # Make a flat array of indices spiralling into the array.
    spiral_idx = np.hstack(spiral_idx)
    # Index into a flattened version of our target array with spiral indices.
    spiral = np.empty_like(arr)
    spiral.flat[spiral_idx] = arr.flat[::-1]
    return spiral

# edge size of the square array.
w = 251
# Prime numbers up to and including w**2.
primes = np.array([n for n in range(2,w**2+1) if all(
                        (n % m) != 0 for m in range(2,int(np.sqrt(n))+1))])
# Create an array of boolean values: 1 for prime, 0 for composite
arr = np.zeros(w**2, dtype='u1')
arr[primes-1] = 1
# Spiral the values clockwise out from the centre
arr = make_spiral(arr.reshape((w,w)))

plt.matshow(arr, cmap=cm.binary)
plt.axis('off')
plt.show()

enter image description here

The prominent diagonal structures indicate that there are unexpectedly many primes of the form $f(n) = 4n^2 + bn + c$ for integer constants $b$ and $c$ and $n=1,2,3,\cdots$.