ρ法(モンテカルロ法)

素因数分解の方法の一つ。根本的な原理は、ランダムに選ばれた２つの整数a, b(＜N)で
GCD(a - b, N)≠1, N
となればいいなぁ、ってことで、それ故に「モンテカルロ法」の別名もあるわけです。

その「ランダムに選ぶ」方法をアルゴリズムとして定義しましょう。よく使われるのが

です。本当は擬似乱数列を生成するものなら何でも良いわけですが、変動の大きさと計算量の小ささから上記の式が好まれる(？)わけです。

さて、上記の式で数列 x_i が定義されます。この数列で GCD(x_i-x_j,N) (i≠j) を取りまくるのが原理なのですが、これをまんまやるのは（過去の計算結果を残すため）メモリ的にも（同じ計算をするため）計算量的にも無駄が多くなります。
ρ法に使う疑似数列は二項漸化式で決まるので，有限体では適当な周期をもってループします．ここで mod p での疑似乱数が周期 t で周回するとすると j=2i だけ計算すれば i=t になった時に結果が出るわけです．（若干嘘言ってますけど．）
なので

def rho(n) :
    x = y = 1
    while true : # 無限ループ…にするのはあまりよろしくない
        x = (x * x + 1) % n
        y = (y * y + 1) % n
        y = (y * y + 1) % n
        p = gcd(abs(x - y), n)
        if 1 < p and p < n :
            return p

という風になりますが，計算に無駄が多いのです．具体的には j=2→4 の移動の際に y は x₃ を経由する訳ですが，i=3 になったときにまた x が x₃ を再計算するわけです．

で，さらに Brent により，i=2^k-1に対して2^k+1−2^k-1≦j＜2^k+1−1を計算すると，そういう無駄を削減し素因数が出る確率は変わらない方法を提案されました．
結構適当なコードにすると

def rho2(n):
    x = y = 1
    k = 2
    while true :
        y = x
        for i in range(1, k / 2) :
            x = (x * x + 1) % n
        for i in range(1, k / 2) :
            x = (x * x + 1) % n
            p = gcd(abs(x - y), n)
            if 1 < p and p < n :
                return p
        x = (x * x + 1) % n
        k *= 2

こんな感じ？若干嘘コードかも．