If p is the perimeter of a right angle triangle with integral length
sides, {a,b,c}, there are exactly three solutions for p = 120.

{20,48,52}, {24,45,51}, {30,40,50}

For which value of p < 1000, is the number of solutions maximised?


You may remember from school the difference of 2 squares

$$a^2+b^2=c^2$$

$$a^2 = c^2 - b^2$$

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

let $$m=c+b$$ and

$$n=c-b$$ so $$m+n=2c$$

$$m-n=2b$$

and substituting these into

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

gives

$$a^2=mn$$

$$a=\sqrt{mn}$$

to get rid of that root, set

$$2p^2=m$$

and

$$2q^2=n$$

and you get

$$a=2pq$$

$$b=p^2-q^2$$

$$c=p^2+q^2$$

Primitive triples are ones in lowest terms (if you know (3,4,5) you can get (6,8,10) etc by multiplying each side by some integer)

If p and q are coprime and different parity you just get primitive ones (i did not bother and just used set() to remove dupes from the list)

from __future__ import division
from collections import defaultdict
from time import time
start = time()

def triples():
for q in range(1, 35):
#35 is greater than sqrt(1000)
for p in range(q+1, 35):
a = p**2 - q**2
b = 2 * p * q
c = p**2 + q**2
for m in range(1,int(1000/c)):
#generate multiples till c is above 1000
#triples with perimeter > 1000 are filtered later
yield (m*min(a,b),m*max(a,b), m*c)

counter = defaultdict(list)

for triple in triples():
perimeter = sum(triple)
if perimeter < 1001:
counter[perimeter] += [triple]

print set((counter[120]))

maxlen = 1
for i in counter:
length = len(set(counter[i]))
if length > maxlen:
maxlen = length
print i
end = time()
print "took ", end-start


Which outputs:

set([(24, 45, 51), (20, 48, 52), (30, 40, 50)])
520
528
630
660
720
840
took  0.0140960216522


I tried brute forcing it today and knocked up the following

from collections import defaultdict
counter = defaultdict(int)
from time import time</p>
def run():
start = time()
for perimeter in range(1,1001):
#print perimeter
for a in range(1,perimeter):
for b in range(1,perimeter-a):
c = perimeter - a - b
if a**2 + b**2 == c**2:
counter[a+b+c] += 1
mosttriples = 0
for i in counter:
if counter[i] > mosttriples:
mosttriples = counter[i]
print i
end = time()
print "TOOK", end-start</p>

print "without psyco"
run()

print "with"
import psyco
psyco.full()
run()


Which gives

without psyco
516
520
528
630
660
720
840
TOOK 134.663383007

with
516
520
528
630
660
720
840
TOOK 21.0988409519


So psyco makes it 6 times faster but a smarter algorithm is 2000 times faster still.