#!/usr/bin/env python
# coding: utf-8

# # Euler's Sum of Powers Conjecture
# 
# In 1769, Leonhard Euler [conjectured that](https://en.wikipedia.org/wiki/Euler%27s_sum_of_powers_conjecture) for all integers *n* and *k* greater than 1, if the sum of *n* *k*th powers of positive integers is itself a *k*th power, then *n* is greater than or equal to *k*. For example, this would mean that no sum of a pair of cubes (a<sup>3</sup> + b<sup>3</sup>) can be equal to another cube (c<sup>3</sup>), but a sum of three cubes can, as in 3<sup>3</sup> + 4<sup>3</sup> + 5<sup>3</sup> = 6<sup>3</sup>. 
# 
# It took 200 years to disprove the conjecture: in 1966 L. J. Lander and T. R. Parkin published a refreshingly short [article](https://projecteuclid.org/download/pdf_1/euclid.bams/1183528522) giving a counterexample of four fifth-powers that summed to another fifth power. They found it via a program that did an exhaustive search.  Can we duplicate their work and find integers greater than 1 such that 
# *a*<sup>5</sup> + *b*<sup>5</sup> + *c*<sup>5</sup> + *d*<sup>5</sup> = *e*<sup>5</sup> ?
# 
# ## Algorithm
# 
# An exhaustive *O*(*m*<sup>4</sup>) algorithm woud be to look at all values of *a, b, c, d* < *m* and check if *a*<sup>5</sup> + *b*<sup>5</sup> + *c*<sup>5</sup> + *d*<sup>5</sup>  is a fifth power. But we can do better: a sum of four numbers is a sum of two pairs of numbers, so we
# are looking for
# 
# &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;*pair*<sub>1</sub> + *pair*<sub>2</sub> = *e*<sup>5</sup> &nbsp;&nbsp; **where** &nbsp; *pair*<sub>1</sub> =  *a*<sup>5</sup> + *b*<sup>5</sup> **and**  *pair*<sub>2</sub> = *c*<sup>5</sup> + *d*<sup>5</sup>.
# 
# We will define *pairs* be a dict of `{`*a*<sup>5</sup> + *b*<sup>5</sup>`: (`*a*<sup>5</sup>`, ` *b*<sup>5</sup>`)}` entries for all *a*  ≤ *b* < *m*; for example, for *a*=2 and *b*=10, the entry is  `{100032: (32, 100000)}`.
# Then we can ask for each *pair*<sub>1</sub>, and for each *e*, whether there is a *pair*<sub>2</sub> in the `dict` that makes the equation work. There are *O*(*m*<sup>2</sup>) pairs and *O*(*m*) values of *e*, and `dict` lookup is *O*(1), so the whole algorithm is *O*(*m*<sup>3</sup>):

# In[1]:


import itertools

def euler(m):
    """Yield tuples (a, b, c, d, e) such that a^5 + b^5 + c^5 + d^5 = e^5,
    where all are integers, and 1 < a ≤ b ≤ c ≤ d < e < m."""
    powers = [e**5 for e in range(2, m)] 
    pairs  = {sum(pair): pair 
              for pair in itertools.combinations_with_replacement(powers, 2)}
    for pair1 in pairs:
        for e5 in powers:
            pair2 = e5 - pair1
            if pair2 in pairs:
                yield fifthroots(pairs[pair1] + pairs[pair2] + (e5,))
    
def fifthroots(nums): 
    "Sorted integer fifth roots of a collection of numbers." 
    return tuple(sorted(int(round(x ** (1/5))) for x in nums))


# Let's look for a solution (arbitrarily choosing *m*=500):

# In[2]:


get_ipython().run_line_magic('time', 'next(euler(500))')


# That was easy, and it turns out this is the same answer that Lander and Parkin got: 27<sup>5</sup> + 84<sup>5</sup> +  110<sup>5</sup> +  133<sup>5</sup> = 144<sup>5</sup>.
# 
# We can keep going, collecting all the solutions up to `*m*=1000`:

# In[3]:


get_ipython().run_line_magic('time', 'set(euler(1000))')


# All the answers are multiples of the first one (this is easiest to see in the middle column: 110, 220, 330, ...).
# Since 1966 other mathematicians have found [other solutions](https://en.wikipedia.org/wiki/Euler%27s_sum_of_powers_conjecture), but all we need is one to disprove Euler's conjecture.