# Conway's Subprime Fib Sequence

### A developer gives a quick tutorial on how to use Python to perform high-level mathematics, using the subprime fib sequence as an example.

Join the DZone community and get the full member experience.

Join For FreeHere's how to construct what John Horton Conway calls his "subprime fib" sequence. Seed the sequence with any two positive integers. Then at each step, sum the last two elements of the sequence. If the result is prime, append it to the sequence. Otherwise, divide it by its smallest prime factor and append that.

If we start the sequence with (1, 1) we get the sequence

1, 1, 2, 3, 5, 4, 3, 7, 5, 6, 11, 17, 14, 31, 15, 23, 19, 21, 20, 41, 61, 51, 56, 107, 163, 135, 149, 142, 97, 239, 168, 37, 41, 39, 40, 79, 17, 48, 13, 61, 37, 49, 43, 46, 89, 45, 67, 56, 41, 97, 69, 83, 76, 53, 43, 48, 13, ...

We know it will keep repeating because 48 followed by 13 has occurred before, and the "memory" of the sequence only stretches back two steps.

If we start with (6, 3) we get

6, 3, 3, 3, ...

It's easy to see that if the sequence ever repeats a value *n* then it will produce only that value from then on, because *n* + *n* = 2 *n* is obviously not prime, and its smallest prime factor is 2.

Conway believes that this sequence will always get stuck in a cycle but this has not been proven.

Here's Python code to try out Conway's conjecture:

```
from sympy import isprime, primefactors
def subprime(a, b):
seq = [a, b]
repeated = False
while not repeated:
s = seq[-1] + seq[-2]
if not isprime(s):
s //= primefactors(s)[0]
repeated = check_repeat(seq, s)
seq.append(s)
return seq
def check_repeat(seq, s):
if not s in seq:
return False
if seq[-1] == s:
return True
# all previous occurances of the last element in seq
indices = [i for i, x in enumerate(seq) if x == seq[-1]]
for i in indices[:-1]:
if seq[i+1] == s:
return True
return False
```

I've verified by brute force that the sequence repeats for all starting values less than 1000.

Published at DZone with permission of John Cook, DZone MVB. See the original article here.

Opinions expressed by DZone contributors are their own.

Comments