Ulam’s spiral in Mathematica

Today we’re having a little fun with Mathematica. We will generate a visualization of the Ulam spiral. To this end the natural numbers are arranged in a spiral pattern,

37 36 35 34 33 32 31
38 17 16 15 14 13 30
39 18  5  4  3 12 29 .
40 19  6  1  2 11 28 .
41 20  7  8  9 10 27 .
42 21 22 23 24 25 26 51
43 44 45 46 47 48 49 50

Now all non-prime numbers are removed:

37                31
   17          13   
       5     3    29  
   19        2 11     
41     7              
         23            
43          47         

Finally, one represents each prime number by a black square and the rest remains white.

That’s a perfect job for Mathematica! We generate a square array of dimension n, where n is odd, and fill it according to the above pattern. We choose n to be odd since we want the number “1” to be in the center. After that we use Mathematica’s built-in functions for prime testing and array plotting. Getting the array right is the only slightly non-trivial part. There are probably, no, certainly a lot of ways to achieve the goal. I’ve chosen to do it with two nested loops, basically in the spirit I would use in usual programming languages:

n = 101;
arr = ConstantArray[0, {n, n}];
arr[[n, n]] = n^2;
For[k = 1, k < (n + 1)/2, k++,
  arr[[(n - 1)/2 + k, (n - 1)/2 + k]] = (2 k - 1)^2;
  For[i = 0, i <= 2 k - 1, i++,
    arr[[(n - 1)/2 + k - i, (n + 1)/2 + k]] = (2 k - 1)^2 + 1 + i;
    arr[[(n + 1)/2 - k, (n + 1)/2 + k - i]] = (2 k - 1)^2 + 2 k + i;
    arr[[(n + 1)/2 - k + i, (n + 1)/2 - k]] = (2 k - 1)^2 + 4 k + i;
    arr[[(n + 1)/2 + k, (n + 1)/2 - k + i]] = (2 k - 1)^2 + 6 k + i;
  ]
]
ArrayPlot[Boole@PrimeQ@arr, PixelConstrained -> {2, 2}]

This code produces the picture
ulam
Looks funny, right? Each black square represents a prime number. It looks as if there are certain lines in this picture, where unusually many primes accumulate. Is that coincidence? Not quite. If one produces a larger picture with n=201
ulam200
or with n=1001
ulam1000
the existence of of these lines seem to persist even at large scales. I invite you to search through the literature and find out, what features of the Ulam spiral are actually understood on a mathematically rigorous level. Feel free to share your findings with me and the other readers through the comment section.

Now that we have an array filled with Ulam’s prescription, we can continue playing. For example we can plot the number of divisors (reverting back to the 201 square array for convenience)

Map [Length, Map [Divisors, arr, {2}], {2}]
ArrayPlot[%, PixelConstrained -> {2, 2}]

numberdivisors
or the sum of divisors

Map [Total, Map [Divisors, arr, {2}], {2}]
ArrayPlot[%, PixelConstrained -> {2, 2}]

sumofdivisors
Let’s stop here before we fill an entire photo album. You can of course continue doodling around. Try the Euler totient function or the Moebius function. Or whatever comes to your mind. Maybe something interesting pops up. Supposedly Ulam found the spiral by fooling around, and the outcome became famous.

Advertisements

About goobypl5

pizza baker, autodidact, particle physicist
This entry was posted in Mathematica and tagged , , , . Bookmark the permalink.

Share your thoughts

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s