离散对数的波拉德ρ算法

离散对数的波拉德ρ算法是约翰·波拉德（英语：John Pollard (mathematician)）1978年所发明解决离散对数问题的算法。

算法的目标是求 $\gamma$ 使得 $\alpha ^{\gamma }=\beta$ ，其中 $\beta$ 属于一个由 $\alpha$ 生成的 循环群 $G$ 。该算法寻找 $a$ , $b$ , $A$ , $B$ 使得 $\alpha ^{a}\beta ^{b}=\alpha ^{A}\beta ^{B}$ 。若他们基于的群是一个 $n$ 阶的循环群，则 $\gamma$ 是方程 $(B-b)\gamma =(a-A){\pmod {n}}$ 的其中一个解。

为求得 $a$ , $b$ , $A$ , $B$ ，该算法使用 Floyd判圈算法 在数列 $x_{i}=\alpha ^{a_{i}}\beta ^{b_{i}}$ 中寻找一个环。假设映射 $f:x_{i}\mapsto x_{i+1}$ 是近似于随机的，则有可能在约 ${\sqrt {\frac {\pi n}{2}}}$ 步后发现一个环。可使用一下规则来生成一个此类映射：将 $G$ 分割为三个不相交的子集 $S_{0}$ ， $S_{1}$ ， $S_{2}$ ，且其所含元素数量大致相等，如果 $x_{i}\in S_{0}$ 则将 $a$ 和 $b$ 加倍；如果 $x_{i}\in S_{1}$ 则将 $a$ 自增；如果 $x_{i}\in S_{2}$ 则将 $b$ 自增。

算法

使 $G$ 是一个 $p$ 阶的 循环群, 且有 $\alpha ,\beta \in G$ , 以及一个分割 $G=S_{0}\cup S_{1}\cup S_{2}$ , 定义映射 $f:G\to G$

f(x)={\begin{cases}\beta x&x\in S_{0}\\x^{2}&x\in S_{1}\\\alpha x&x\in S_{2}\end{cases}}

并据以下方式定义映射 $g:G\times \mathbb {Z} \to \mathbb {Z}$ 和 $h:G\times \mathbb {Z} \to \mathbb {Z}$

{\begin{aligned}g(x,n)&={\begin{cases}n&x\in S_{0}\\2n{\pmod {p}}&x\in S_{1}\\n+1{\pmod {p}}&x\in S_{2}\end{cases}}\\h(x,n)&={\begin{cases}n+1{\pmod {p}}&x\in S_{0}\\2n{\pmod {p}}&x\in S_{1}\\n&x\in S_{2}\end{cases}}\end{aligned}}

 输入 a: a 是 G 的生成元, b: G 的一个元素
 输出 整数 x 使得 a^x = b, 或者失败

 初始化 a₀ ← 0, b₀ ← 0, x₀ ← 1 ∈ G, 

 i ← 1
 loop
     x_i ← f(x_i-1), 
     a_i ← g(x_i-1, a_i-1), 
     b_i ← h(x_i-1, b_i-1)

     x_2i ← f(f(x_2i-2)), 
     a_2i ← g(f(x_2i-2), g(x_2i-2, a_2i-2)), 
     b_2i ← h(f(x_2i-2), h(x_2i-2, b_2i-2))

     if x_i = x_2i then
         r ← b_i - b_2i
         if r = 0 return failure
         x ← r⁻¹(a_2i - a_i) mod p
         return x
     else # x_i ≠ x_2i
         i ← i+1, 
         break loop
     end if
  end loop

举例

例如一个由 2 模 $N=1019$ 生成的群（群的阶是 $n=1018$ ，2是生成元，生成群的元素模1019同余)。这个算法可由以下 C++ 程序实现。

 #include <stdio.h>
 
 const int n = 1018, N = n + 1;  /* N = 1019 -- 素数       */
 const int alpha = 2;            /* 生成元                 */
 const int beta = 5;             /* 2^{10} = 1024 = 5 (N) */
 
 void new_xab( int& x, int& a, int& b ) {
   switch( x%3 ) {
   case 0: x = x*x     % N;  a =  a*2  % n;  b =  b*2  % n;  break;
   case 1: x = x*alpha % N;  a = (a+1) % n;                  break;
   case 2: x = x*beta  % N;                  b = (b+1) % n;  break;
   }
 }
 
 int main(void) {
   int x=1, a=0, b=0;
   int X=x, A=a, B=b;
   for(int i = 1; i < n; ++i ) {
     new_xab( x, a, b );
     new_xab( X, A, B );
     new_xab( X, A, B );
     printf( "%3d  %4d %3d %3d  %4d %3d %3d\n", i, x, a, b, X, A, B );
     if( x == X ) break;
   }
   return 0;
 }

结果如下 (已截断):

 i     x   a   b     X   A   B
------------------------------
 1     2   1   0    10   1   1
 2    10   1   1   100   2   2
 3    20   2   1  1000   3   3
 4   100   2   2   425   8   6
 5   200   3   2   436  16  14
 6  1000   3   3   284  17  15
 7   981   4   3   986  17  17
 8   425   8   6   194  17  19
..............................
48   224 680 376    86 299 412
49   101 680 377   860 300 413
50   505 680 378   101 300 415
51  1010 681 378  1010 301 416

可见 $2^{681}5^{378}=1010=2^{301}5^{416}{\pmod {1019}}$ 以及 $(416-378)\gamma =681-301{\pmod {1018}}$ 。正如预期，其中 $\gamma _{1}=10$ 是一个解。由于 $n=1018$ 不是素数，因此存在另一个解 $\gamma _{2}=519$ ，使得 $2^{519}=1014=-5{\pmod {1019}}$ 成立。

复杂度

时间复杂度近似于 ${\mathcal {O}}({\sqrt {n}})$ 。如果配合使用 Pohlig-Hellman 算法（英语：Pohlig-Hellman algorithm），则整体时间复杂度近似于 ${\mathcal {O}}({\sqrt {p}})$ ，其中 $p$ 是 $n$ 的最大质因数。

参考文献

Pollard, J. M. Monte Carlo methods for index computation (mod p). Mathematics of Computation. 1978, 32 (143): 918–924. doi:10.2307/2006496.
Menezes, Alfred J.; van Oorschot, Paul C.; Vanstone, Scott A. Chapter 3 (PDF). Handbook of Applied Cryptography. 2001 [2018-04-16]. （原始内容 (PDF)存档于2012-02-08）.