It was proposed by Christian Goldbach that every odd composite number can be written as the sum of a prime and twice a square.
$$ 9 = 7 + 2 \cdot 1^2 \\ 15 = 7 + 2 \cdot 2^2 \\ 21 = 3 + 2 \cdot 3^2 \\ 25 = 7 + 2 \cdot 3^2 \\ 27 = 19 + 2 \cdot 2^2 \\ 33 = 31 + 2 \cdot 1^2 $$It turns out that the conjecture was false.
What is the smallest odd composite that cannot be written as the sum of a prime and twice a square?
Formally stated, Goldbach's Other Conjecture says that all odd composite numbers $n$ can be expressed as
$$ n = 2k^2 + p $$for some integer $k$ and prime number $p$. Let
$$ S_n = \{ n - 2k^2 : k = 1, 2, \dotsc, \lfloor \sqrt{\frac{n}{2}} \rfloor \} $$If any element $n - 2k^2$ of $S_n$ is prime, then we call $k$ a witness to Goldbach's Other Conjecture for $n$. We are asked to find the smallest $n$ that has no witness, i.e. such that all elements of $S_n$ are composite.
Let $P_n$ be the set of all prime numbers stricly smaller than $n$, then our algorithm searches for the smallest $n$ such that $P_n \cap S_n = \emptyset$, providing a counterexample to Goldbach's Other Conjecture.
In [1]:
from IPython.display import Math, display
from collections import defaultdict
from itertools import count
We modify and augment a very space-efficient implementation (due to David Eppstein, UC Irvine) of the Sieve of Erastothenes to generate composite numbers and keep track of all prime numbers below it. In the outermost-loop below, primes will always be the list of all prime number strictly below $n$, i.e. it is equivalent to $P_n$. Note that if $n$ is odd and composite, then the largest prime below it is no greater than $n-2$, since $n-1$ is even. Since $\max S_n = n - 2$ searching through $P_n$ is sufficient (i.e. if $n-2$ is prime, then $n-2 \in P_n$.) The loop terminates when we encounter an $n$ with no witnesses.
In [2]:
factors = defaultdict(list)
witness = {}
primes = set()
for n in count(2):
if factors[n]:
# n is composite
for m in factors.pop(n):
factors[n+m].append(m)
if n % 2:
# n is odd and composite
for k in range(1, int((n/2)**.5)+1):
p = n - 2*k**2
if p in primes:
# TODO: `in` is O(len(primes))
# could optimize by using set
witness[n] = k
break
if not n in witness:
break
else:
factors[n*n].append(n)
primes.add(n)
In [3]:
n
Out[3]:
The answer is 5777.
Note that not only is this implementation space-efficient, it is also very time efficient, since we incrementally build our list of primes, composites and witnesses incrementally from bottom-up. If we hadn't augmented the implementation of the prime sieve, we would have had to use the prime sieve to obtain all odd composites, and perform a primality test on all elements of $S_n$, which would (usually) require a prime factorization algorithm, which in turn (usually) requires a prime sieve, not to mention the fact that there would be many overlapping subproblems, i.e. many $n < m$ such that $S_n \cap S_m \ne \emptyset$, so we'd have to use dynamic programming, by, for example memoizing the primality testing function or some other optimization - just a whole mess of redundancies and inefficiencies that we happily avoided with this method :)
Let's list out the witnesses to the first 100 numbers.
In [4]:
lines = [r'{0} &= {1} + 2 \cdot {2}^2 \\'.format(n, n - 2*witness[n]**2, witness[n]) for n in sorted(witness)]
In [5]:
Math(r"""
\begin{{align}}
{body}
\end{{align}}
""".format(body='\n'.join(lines[:100])))
Out[5]: