Tonelli–Shanks Algorithm

参考：http://en.wikipedia.org/wiki/Tonelli%E2%80%93Shanks_algorithm

Tonelli–Shanks算法是一个求解二次平方根的算法

$$x^2 \equiv n \pmod p$$

其中n是p的二次剩余，p为奇素数

算法

输入:

• p,一个素数
• n,一个模p的二次剩余

输出:

• r,

算法:

1. 从p-1分解出2的幂次，即有:p-1=Q*2^S，其中Q是奇素数。如果S=1，即，然后直接返回。
2. 计算z，使得其满足L(z,p)=-1，令
3. 令
4. 循环：
   • 如果，返回r
   • 否则，找到最小的一个i(0<i<m)，且
   • 令 ，再令 ,,,

如果得到一个解r，另一个解就是p-r

证明

首先有p-1=Q*2^S，令，注意到有，这一同余式在每次循环中都保持正确；如果在某一时间点，，则有，于是就找到n的二次平方根

如果，那就考虑二次非剩余z；令，然后就有，并且，这意味着c的阶是

类似地，，故t的阶能整除；假设t的阶是，由于n是模p的二次剩余，

现在令,,,,和之前一样，任然成立；然而现在的t和c’的阶数都是，这意味着t’的阶数满足

如果，则，算法终止，返回，否则重新执行循环，定义……直到停止；由于S序列严格递减，算法一定会结束

实现(python)

def legendre(a,p):
    symbol = pow(a, (p - 1) // 2, p)
    if symbol == p - 1:
        return -1
    return symbol

def tonelli(a,p):
    if a == 0 or legendre(a,p) != 1:
        return 0
    q = p - 1
    s = 0
    while q % 2 == 0:
        q //= 2
        s += 1
    if s == 1:
        return pow(n, (p + 1) // 4, p)
    
    z = 2
    while legendre(z, p) != -1:
        z += 1
        
    m = s
    c = pow(z, q, p)
    t = pow(a, q, p)
    r = pow(a, (q +1) // 2, p)
    
    while t!= 1:
        t2 = t
        i = 0
        while t2 != 1 and i < m:
            t2 = pow(t2, 2, p)
            i += 1
        b = pow(c, 2 ** (m - i - 1), p)
        m = i
        c = (b * b) % p
        t = (t * c) % p
        r = (r * b) % p
        
    return r