Project Euler 100

Problem 100
もし箱の中に21個の色がついたディスクがあり、それが15個の青いディスクと6個の赤いディスクからなっており、2つのディスクをランダムに取ったとすると、2つの青いディスクを取る確率は、P(BB) = (15/21)×(14/20) = 1/2であることがわかる。
ランダムで2つの青いディスクを取る確率がちょうど50%であるようなこの次の構成は、85個の青いディスクと35個の赤いディスクである。
あわせて1012を超えるディスクを含む最初の構成を見つけて、青いディスク数を決めよ。
http://projecteuler.net/index.php?section=problems&id=100

こういう解が飛び飛びで巨大な数を求める問題は、Pell方程式であると想像されます。
さて、青いディスクの数をm、合計のディスクの数をnとすると、

(m/n)×(m-1)/(n-1) = 1/2

となります。

2m(m - 1) = n(n - 1)
8(m2 - m) = 4(n2 - n)
2(2m - 1)2 - 2 = (2n - 1)2 - 1
(2n - 1)2 - 2(2m - 1)2 = -1

つまり、

x2 - 2y2 = -1
x, y : 奇数

を解けばいいだけです。

from itertools import ifilter

def head(a):
    for e in a:
        return e

def gen_solutions():
    a, b = 1, 1
    while True:
        a, b = a * 3 + b * 4, a * 2 + b * 3
        if a % 2 == 1 and b % 2 == 1:
            m, n = (b + 1) / 2, (a + 1) / 2
            yield m, n

N = 10 ** 12
print head(x[0] for x in ifilter(lambda x: x[1] > N, gen_solutions()))

この企画はこれで終わりです。