def solutions(p): solution = 0 for a in xrange(1, p/2): for b in xrange(a, p-a): c = p - a - b if a*a + b*b == c*c: solution += 1 return solution max_solutions = 0 max_p = 0 for p in xrange(3, 1000): n = solutions(p) if n > max_solutions: max_solutions = n max_p = p print max_p