Project Euler 39

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.