Project Euler 138

http://projecteuler.net/index.php?section=problems&id=138


(b / 2, h, L)でピタゴラス数で、hが偶数だとhbの差が偶数になってしまうので、

h = m2 - n2 b = 4mn L = m2 + n2
m2 - n2 - 4mn = ±1
(m - 2n)2 - 5n2 = ±1

これもペル方程式です。

from itertools import *

def unfold(f, s):
    while True:
        x, s = f(s)
        yield x

def gen_Pells():
    def mul(p, q):
        x1, y1 = p
        x2, y2 = q
        return (x1 * x2 + 5 * y1 * y2, x1 * y2 + y1 * x2)
    
    a = (2, 1)
    return unfold(lambda s: (s, mul(s, a)), a)

N = 12
print sum((x + 4 * y) * x + 5 * y * y
                    for x, y in islice(gen_Pells(), N))