# Project Euler, problem # 9, Pythagorean triplet

A Pythagorean triplet is a set of three natural numbers $$\a < b < c\$$ for which $$\a^2 + b^2 = c^2\$$.

For example, $$\3^2 + 4^2 = 9 + 16 = 25 = 5^2\$$.

There exists exactly one Pythagorean triplet for which $$\a + b + c = 1000\$$. Find the product $$\a b c\$$.

Here is my implementation in Python, awaiting your feedback.

def get_triplet():
for c in range(2, 1000):
for b in range(2, c):
a = 1000 - c - b
if a ** 2 + b ** 2 == c ** 2:
return 'n1 = %s\nn2 = ' \
'%s\nn3 = %s\nproduct = %s' \
% (a, b, c, a * b * c)

if __name__ == '__main__':
print(get_triplet())

• Please ensure that your code is posted with the intended formatting. One way to post code is to paste it into the question editor, highlight it, and press Ctrl-K to mark it as a code block. Jul 14, 2019 at 5:53
• There is an answer invalidation issue that has happened because of edits to the indenting of the return line... but I believe the indent issue resulted from bad question formatting... more than bad code formatting. part of one answer has been invalidated, but the overall effect of keeping the edit is better than rolling it back. Jul 16, 2019 at 19:16
• @rolfl I apologize for that; I'm the one who proposed the indent formatting fix; I hadn't read any of the answers when I suggested that, in part b/c code indent issues have usually been the more "obvious" fixes. I'm hoping the answers that are referencing the indents can be updated, or if it's too much, then I suppose my fix may need a rollback and the question marked with something for others to avoid fixing the same issue over and over. Jul 16, 2019 at 19:26

## Search Space

Your search space is too large. Since 0 < a < b < c and a + b + c = 1000, you can put hard limits on both a and c. The maximum a can be is 332, since 333 + 334 + 335 > 1000. Similarly, the minimum value c can be is 335, since 332 + 333 + 334 < 1000. As well, the maximum value c can have is 997 (1 + 2 + 997 = 1000).

Based on the above, if you wish your outside loop to loop over c values, your first loop should be:

for c in range(335, 997 + 1):
# ...


And the inner loop should be over a values:

    for a in range(1, 332 + 1):
b = 1000 - c - a
# ...


But now that we’ve selected a c value, can tighten the limits on a even further. When c is 600, then a+b must be 400, and since a < b must hold, the maximum a can have is 400 / 2. So:

    for a in range(1, min(332, (1000 - c) // 2) + 1):
b = 1000 - c - a
# ...


Moreover, since we’ve selected a c value, we have established a maximum value for b as well, since b < c must hold. If c is 400, then b can’t be larger than 399, so a must be at least 1000 - 400 - 399 = 201, establishing a lower limit for our inner loop:

    for a in range(max(1, 1000 - c - (c - 1)), min(332, (1000 - c) // 2) + 1):
b = 1000 - c - a
# ...


## Geometry

As pointed out in a comment by @Peter Taylor below, we can tighten the limits on c much further.

First, triangles must be triangular. In every triangle, the length of any side must be less than the sum of the lengths of the other two sides. In particular, $$\c < a + b\$$ must be true. Since $$\a + b + c = 1000\$$, we can conclude $$\c < 1000 - c\$$, or more simply, $$\c < 500\$$, establishing a much smaller upper limit.

The smallest value c could obtain would be when a and b are at their largest. If we relax the conditions slightly, allowing $$\a \leq b\$$ and $$\a, b, c, \in \Re\$$; we have an isosceles right-angle triangle, $$\a = b = c / \sqrt 2\$$. Then:

$$\a + b + c = 1000\$$

$$\c = 1000 \sqrt 2 / (2 + \sqrt 2) \approx 414.2\$$

giving us a larger lower limit.

Thus, the range of the outer loop can be reduced to:

for c in range(413, 500):
# ...


## Repeated Calculations

When searching for the Pythagorean triplet, you are calculating c ** 2 each and every iteration of the inner loop. The value of c is constant during execution of the inner loop, so c ** 2 is also constant. Since exponentiation is a relatively expensive operation, you should get a speed improvement by moving this calculation out of the inner loop:

for c in range(413, 500):
c2 = c ** 2
for a in range(max(1, 1000 - c - (c - 1)), min(332, (1000 - c) // 2) + 1):
b = 1000 - c - a
if a ** 2 + b ** 2 == c2:
return a, b, c


1000 - c is also repeated in the inner loop, so you could compute this outside the loop also; it should provide a slight performance improvement.

See other answers for important PEP8, f-string, docstring, and return value improvements not duplicated here.

• Since "exponentiation is a relatively expensive operation", wouldn't it be faster to also cache the squares in a first loop? There are a lot of repetitions at runtime of a**2 and b**2. Jul 15, 2019 at 9:15
• @OlivierGrégoire Don’t cache the squares in the first loop! That way leads to conditionals (if not cached[n]: ...) which will slow your whole code down more than you gain, and convolute your code. Cache the squares before the first loop. And don’t store the cached values in a dict; store them in a list (preferably an array.array()). Jul 15, 2019 at 15:16
• If you're going to analyse the ranges, you could also consider that at one extreme $a=1, b \approx c \approx 500$ and at the other extreme $a \approx b \approx c / \sqrt 2$, so $c \approx 1000 \sqrt 2 / (2 + \sqrt 2) \approx 414$. Tidy that up a tiny bit and the range is greatly reduced. Jul 15, 2019 at 15:26
• @AJNeufeld Actually, I wrote "in a first loop", not "in the first loop". ;-) Jul 15, 2019 at 16:16
• 3Blue1Brown has a great video that includes techniques to reduce the search space to basically zero by directly generating the next triples from the previous one. I did something similar in this answer from 2011. Jul 15, 2019 at 19:38

Note that your code actually gives the wrong answer, it returns a > b. Since a and b are the two smaller variables, you should actually start looping over them in order, so it looks more like:

for a in range(1, 1000):  # Don't know why you assumed a > 1 here
for b in range(a + 1, 1000):
c = 1000 - a - b


Also, it doesn't matter for this problem, but for future ones it's important to notice when you can break early to speed up loops:

c2 = a**2 + b**2
if c2 == c ** 2:
return a, b, c
elif c2 > c ** 2:
break

• Uhm, why is it better to loop starting from the smallest variable first? It seems to me that your version tests many more triples than OP's one. Jul 15, 2019 at 12:30
• @FedericoPoloni As I stated, OP's code actually gives the wrong result. The right triples, but in the wrong order. Jul 16, 2019 at 0:35
• @FedericoPoloni not to mention once you add the shortcut break in, my code does less checks than OPs, even if you include the shortcut with OPs method. Jul 16, 2019 at 2:21
• The wrong order seems irrelevant --- just rename a to b and vice versa. Also, you still haven't explained why in general it's a better idea to start looping from the smallest variable. Jul 16, 2019 at 6:25
• @FedericoPoloni The smallest variable will go through the smallest number of iterations. If a < c, you'd rather go through your outer loop a times than c times. Jul 16, 2019 at 22:37

Instead of returning a custom string, just return the values and leave the printing to the caller. This way it is at least feasible for this function to be used elsewhere.

def get_triplet():
...
return a, b, c

if __name__ == '__main__':
a, b, c = get_triplet()
print(f"a = {a}, b = {b}, c = {c}, product = {a*b*c}")


Here I also used an f-string (Python 3.6+) to make the formatting easier.

You should also have a look at Python's official style-guide, PEP8, which recommends using a consistent number of spaces (preferably 4) as indentation. Your return line does not follow this practice. In addition, functions (and classes) should have only two newlines separating them.

You should also think about adding a docstring to describe what the function does. The name get_triplet does not convey all the information (e.g. that they need to sum to a target value).

If you want a solution that doesn't involve any coding, you can use the fact that Pythagorean triples $$\a < b < c\$$ are of the form

$$\a = 2pqr\$$

$$\b = p(q^2 - r^2)\$$

$$\c = p(q^2 + r^2)\$$

or the same equations with $$\a\$$ and $$\b\$$ switched.

Here $$\p, q, r\$$ are positive integers, which are uniquely determined by the condition $$\p = \gcd(a, b, c)\$$.

This is known as Euclid's formula, a proof of which can be found in the Wikipedia article on Pythagorean triples. However, it may be more enlightening to prove this yourself, since all it requires is an understanding of the unique factorization of integers into products of powers of primes. The trick is to rewrite $$\a^2+b^2=c^2\$$ as

$$\a^2 = c^2-b^2 = (c+b)(c-b)\$$

First suppose that $$\\gcd(a,b,c)=1\$$. What are the possible common factors of $$\c+b\$$ and $$\c-b\$$? Given that $$\(c+b)(c-b) = a^2\$$ is a square, what does this imply?

• I might try doing another version that uses GCD, sounds cool.
– user203258
Jul 15, 2019 at 16:21
• Hi @Pieter, can you link to a source of these statements? Or show the derivation? Jul 16, 2019 at 9:17
• @RobAu I've edited my answer. Jul 16, 2019 at 18:25
• Thanks, this is really helpful Jul 16, 2019 at 19:03

# Solution 1

You have two equations and three unknowns. Since three minus two is one, you should only require one search loop.

a + b + c = 1000
=> c = 1000 - a - b

a2 + b2 = c2
=> a2 + b2 = (1000 - a - b)2 => a2 + b2 = 10002 - 2000a - 2000b + 2ab + a2 + b2
=> 0 = 1000000 - 2000a - 2000b + 2ab
=> 0 = 500000 - 1000a - 1000b + ab
=> 1000b - ab = 500000 - 1000a
=> (1000 - a)b = 1000(500 - a)
=> b = 1000(500 - a)/(1000 - a)

We can now try various values of a in a loop and solve for b and c relative to a. These solutions will all fulfill a + b + c = 1000 and a2 + b2 = c2. We just need to limit the actual solution to integral values of b (which will guarantee integral values of a and c because of the first equation), and stop the search once a >= b (which will guarantee a < b < c).

S = 1000

a = 1
while True:
b = S * (S - 2 * a) / (2 * S - 2 * a)
if a >= b:
break
elif b.is_integer():
b = int(b)
c = S - a - b

print 'a = %s, b = %s, c = %s, a * b * c = %s' % (a, b, c, a * b * c)

return a * b * c
a += 1


# Solution 2

Using Euclid's Formula from Wikipedia, we can drastically reduce this search. Euclid's Formula tells us that:

a = (m2 - n2)/2
b = mn
c = (m2 + n2)/2

First we can constrain these to solutions that meet our requirements.

a + b + c = 1000 => (m2 - n2)/2 + mn + (m2 + n2)/2 = 1000
=> m2 + mn = 1000
=> mn = 1000 - m2
=> n = 1000/m - m

Plugging that back into the equations above, we get the following:

a = (m2 - n2)/2
=> a = (m2 - (10002/m2 - 1000 - 1000 + m2))/2
=> a = (2000 - 1000000/m2)/2
=> a = 1000 - 500000/m2

b = mn
=> b = m(1000/m - m)
=> b = 1000 - m2

c = (m2 + n2)/2
=> c = (m2 + (10002/m2 - 1000 - 1000 + m2))/2
=> c = (2m2 + 1000000/m2 - 2000)/2
=> c = m2 + 500000/m2 - 1000

Further, we can constrain the search to 0 < a < c.

0 < a
=> 0 < 1000 - 500000/m2
=> 500000/m2 < 1000
=> 500000 < 1000m2
=> 500 < m2
=> sqrt(500) < m
=> 22.36068 < m

a < c
=> 1000 - 500000/m2 < m2 + 500000/m2 - 1000
=> 2000 < m2 + 1000000/m2
=> 2000m2 < m4 + 1000000
=> 0 < m4 - 2000m2 + 1000000
=> m2 < (-(-2000) +/- sqrt((-2000)2 - 4 * 1 * 1000000)) / (2 * 1)
=> m2 < (2000 +/- sqrt(4000000 - 4000000)) / 2
=> m2 < 2000 / 2
=> m2 < 1000
=> m < sqrt(1000)
=> m < 31.62278

So that means we only need to test integer values of m 23 to 31 (up to nine cases).

S = 1000

for m in range(math.ceil(math.sqrt(S / 2)),
math.floor(math.sqrt(S) + 1)):
a = S - S ** 2 / (2 * m ** 2)
if a.is_integer():
a = int(a)
b = S - m ** 2
c = S - a - b

print('a = %s, b = %s, c = %s, a * b * c = %s' %
(min(a, b), max(a, b), c, a * b * c))

return a * b * c

• Great algorithmic answer. Consider changing the print to f strings or printing from the caller. What is the computational complexity of both algorithms you showed us? Jul 16, 2019 at 17:09
• Thanks! If we say we're looking for a sum of S=1000, the first solution is O(S) and the second solution is O(sqrt(S)). And yeah, I'm not a fan of the print either and mostly just had it for debug purposes. I think the original question really only requires returning the product, right? Jul 16, 2019 at 17:42
• m^2+mn = 1000 can be factored to m(m+n) = 1000, so m is a factor of 1000. So now we're looking for factors of 1000 between 23 and 26. Jul 16, 2019 at 17:47
• That's an awesome point! Is there an algorithm for finding a factor within a range that doesn't require trying all values in that range? (Also note that the value of 26 was incorrect, and I updated my answer with the real upper limit of 31. Obviously it works for 1000, where the answer occurs at m=25, but for other values of S, like 90, it would not have always worked out.) Jul 16, 2019 at 17:58
• Euclid's formula actually includes a common factor, which is important. In this case, you can find values of m and n for which the common factor is 1 (although this isn't the least common factor), but for a Pythagorean triple like 9, 12, 15, there is no solution in integers of the form (m^2-n^2)/2, m*n, (m^2+n^2)/2, since 2*15 = 30 can't be written as the sum of two squares. Jul 18, 2019 at 14:15

# simplify a and fix search space

Given the nested for loop, your variable a can be decremented simpler:

Moreover, your search space is far too large. Other authors have addressed his in more detail, but if you want to check every possible triplet from (1-1000), then you need to change your second for loop to:

def get_triplet():
for c in range(2, 1000):
a = 1000 - c - 1
for b in range(1, c):
a -= 1


# cut down on for loop size

You obviously can’t have a Pythagorean triple such that $$\b=c\$$, so you may simplify the first for loop to:

for c in range(2, 1000):


# return

First, you may consider printing directly from the function get_triplets cut down on one call of return. This is not good practice unless you are trying to optimize at all costs.

However, your return may be better suited using f strings, which are faster and newer (f for fast). Perhaps you could also not suddenly rename a,b,c to n1,n2,n3:

return  f'a = {a}, b = {b}, c = {c}, product = {a*b*c}'

• I would also suggest not to suddenly rename a,b and c to n1,n2 and n3. Jul 14, 2019 at 8:44
• The way you modify the use of a does not seem simpler but rather more complicated. Now I have to look at two lines to understand what a is Jul 14, 2019 at 14:54
• Furthermore, I don't see how printing from within get_triplets cuts down on a function call, and printing from within a function that finds this kind of data is generally not recommend and at the very least would require renaming the function Jul 14, 2019 at 14:55
• @DreamConspiracy when we simplify a as shown above, the computer performs far fewer subtraction operations. Moreover, the author does not have to print from the function. It would be more reasonable for readability to return and then print, but if they print directly from the function, it eliminates one call of return Jul 14, 2019 at 15:23
• Your calculations seem off. If c is 2, a becomes 998, and then b starts at 2. Now a+b+c = 1002. As b increase, c decreases, maintaining 1002. Jul 14, 2019 at 18:21