分解因子算法——Pollard ρ \rho ρ 算法
攻击RSA密码体制最明显的方法就是分解大整数。关于大整数分解的算法有许多,常见的有试除法、 p ± 1 p\pm1 p±1算法、Pollard ρ \rho ρ算法、数域筛法、二次筛法等等。这篇分享的就是其中之一的Pollard ρ \rho ρ算法。
我们说因子分解,很多时候并不一定要把大整数 n n n彻底分解成素数乘积,而是求出 n n n的某个非平凡因子即可。求大整数的一个因子固然困难,但是Euclidean算法告诉我们求两个数的最大公因子是可以在 O ( log n ) O(\log n) O(logn)的时间得到的。Pollard ρ \rho ρ算法利用的就是这个思想。
问题描述
设大整数为 n n n,我们要找到 n n n的一个因子 p p p。
Pollard ρ \rho ρ 算法
为了利用Euclidean算法找到 p p p,我们设想存在两个整数 x , x ′ ∈ Z n x,x'\in \mathbb{Z}_{n} x,x′∈Zn,满足 x ≡ x ′ ( m o d p ) x\equiv x'\pmod{p} x≡x′(modp),这样的话就有 p ∣ ( x − x ′ ) p|(x-x') p∣(x−x′),又因为 p ∣ n p|n p∣n,所以有 p ≤ g c d ( x − x ′ , n ) < n p\le gcd(x-x',n)<n p≤gcd(x−x′,n)<n——这句话的意思是: g c d ( x − x ′ , n ) gcd(x-x',n) gcd(x−x′,n)一定是 n n n的非平凡因子。这就完成了因子分解。问题归结为:如何找到满足条件的 x x x和 x ′ x' x′。
一个朴素的想法就是:随机选择一个子集 X ⊂ Z n X\subset \mathbb{Z}_{n} X⊂Zn,然后对所有不同的 x , x ′ ∈ X x,x'\in X x,x′∈X,判断 g c d ( x − x ′ , n ) gcd(x-x',n) gcd(x−x′,n)是否为1。这个方法能够成功当且仅当映射 x → x ( m o d p ) x \rarr x\pmod{p} x→x(modp)在 X X X中至少存在一个碰撞。根据生日攻击可以分析:如果 ∣ X ∣ ≈ 1.71 p |X| \approx 1.71\sqrt{p} ∣X∣≈1.71p时,至少存在一个碰撞的概率是50%,为了找到 X X X中存在的这个碰撞,我们需要求 ( ∣ X ∣ 2 ) \begin{pmatrix}|X|\\2 \end{pmatrix} (∣X∣2)次gcd。
Pollard
ρ
\rho
ρ算法减少了求解gcd的次数。它并不是随机产生一个子集
X
X
X,而是事先确定一个整系数多项式
f
f
f,并取
x
1
∈
Z
n
x_{1}\in \mathbb{Z}_{n}
x1∈Zn,令
x
i
+
1
=
f
(
x
i
)
m
o
d
n
,
i
=
1
,
2
,
3
⋯
x_{i+1}=f(x_{i})\mod{n},i=1,2,3\cdots
xi+1=f(xi)modn,i=1,2,3⋯
这就产生了集合
X
=
{
x
1
,
x
2
,
⋯
,
x
m
}
X=\lbrace x_{1},x_{2},\cdots,x_{m} \rbrace
X={x1,x2,⋯,xm}。下面描述这个集合是如何减少gcd计算次数的。
核心思想:如果 x i ≡ x j ( m o d p ) x_{i}\equiv x_{j} \pmod{p} xi≡xj(modp),由于 f f f是整系数多项式,所以 f ( x i ) ≡ f ( x j ) ( m o d p ) f(x_{i})\equiv f(x_{j}) \pmod{p} f(xi)≡f(xj)(modp) 。又因为 p ∣ n p|n p∣n,所以
x i + 1 m o d p = ( f ( x i ) m o d n ) m o d p = f ( x i ) m o d p x_{i+1}\mod{p}=(f(x_{i})\mod{n})\mod{p}=f(x_{i})\mod{p} xi+1modp=(f(xi)modn)modp=f(xi)modp
同理有 x j + 1 m o d p = f ( x j ) m o d p x_{j+1}\mod{p}=f(x_{j})\mod{p} xj+1modp=f(xj)modp。因此 x i + 1 ≡ x j + 1 ( m o d p ) x_{i+1} \equiv x_{j+1}\pmod{p} xi+1≡xj+1(modp)。进一步有:如果 x i ≡ x j ( m o d p ) x_{i}\equiv x_{j} \pmod{p} xi≡xj(modp),则对 ∀ δ ≥ 0 \forall \delta \ge0 ∀δ≥0 ,有 x i + δ ≡ x j + δ ( m o d p ) x_{i+\delta}\equiv x_{j+\delta} \pmod{p} xi+δ≡xj+δ(modp)。
记 l = j − i l=j-i l=j−i,可知 x i ′ ≡ x j ′ ( m o d p ) x_{i'} \equiv x_{j'} \pmod{p} xi′≡xj′(modp)如果 j ′ ≥ i ′ ≥ i j'\ge i'\ge i j′≥i′≥i且 j ′ − i ′ ≡ 0 ( m o d l ) j'-i'\equiv0\pmod{l} j′−i′≡0(modl)。
这样序列X就可以看成带有一个“尾巴”
x 1 → x 2 → ⋯ → x i ( m o d p ) x_{1}\rarr x_{2} \rarr \cdots \rarr x_{i} \pmod{p} x1→x2→⋯→xi(modp)
和一个“环”
x i → x i + 1 → ⋯ → x j = x i ( m o d p ) x_{i} \rarr x_{i+1} \rarr \cdots \rarr x_{j}=x_{i} \pmod{p} xi→xi+1→⋯→xj=xi(modp)
看成图的话就像希腊字母 ρ \rho ρ,这正是算法名字的由来。对这个算法,我们需要计算 ( j − 1 2 ) + i \begin{pmatrix}j-1\\2 \end{pmatrix}+i (j−12)+i次gcd才能得到第一个碰撞 x i , x j x_{i},x_{j} xi,xj。
进一步改进:我们只要求得一个碰撞即可,不一定是 x i x_{i} xi和 x j x_{j} xj,也可以是 x i + δ x_{i+\delta} xi+δ和 x j + δ x_{j+\delta} xj+δ。所以在找碰撞的时候,不必逐个逐个求gcd,可以跳跃着求。我们断言:在第一个“环” x i → x i + 1 → ⋯ → x j = x i ( m o d p ) x_{i} \rarr x_{i+1} \rarr \cdots \rarr x_{j}=x_{i} \pmod{p} xi→xi+1→⋯→xj=xi(modp)中,一定存在 i ′ i' i′满足 x 2 i ′ ≡ x i ′ ( m o d p ) x_{2i'} \equiv x_{i'}\pmod{p} x2i′≡xi′(modp)。事实上,只要 2 i ′ − i ′ = i ′ = k l 2i'-i'=i'=kl 2i′−i′=i′=kl,就有 x 2 i ′ ≡ x i ′ ( m o d p ) x_{2i'} \equiv x_{i'}\pmod{p} x2i′≡xi′(modp)成立;而 i ≤ i ′ = k l ≤ j − 1 i\le i'=kl\le j-1 i≤i′=kl≤j−1,所以最多找 j j j次就可以得到 i ′ i' i′。
算法流程
Pollard ρ \rho ρ (n, x 1 x_{1} x1)
x = x 1 , x ′ = f ( x ) m o d n x=x_{1},x'=f(x)\mod{n} x=x1,x′=f(x)modn
p = g c d ( x − x ′ , n ) p=gcd(x-x',n) p=gcd(x−x′,n)
while p = = 1 p==1 p==1
//第i次迭代
x = x i , x ′ = x 2 i x=x_{i},x'=x_{2i} x=xi,x′=x2i
x = f ( x ) m o d p x=f(x)\mod{p} x=f(x)modp
x ′ = f ( x ′ ) m o d p x'=f(x')\mod{p} x′=f(x′)modp
x ′ = f ( x ′ ) m o d p x'=f(x')\mod{p} x′=f(x′)modp
p = g c d ( x − x ′ , n ) p=gcd(x-x',n) p=gcd(x−x′,n)
end
if p==n
return “failure”
else
return p
算法性能分析
假设 f = x 2 + 1 f=x^{2}+1 f=x2+1,则 j j j的最大值为 p \sqrt{p} p(因为 x j x_{j} xj是第一个环中的数),所以最多需要 p \sqrt{p} p次gcd计算。而 p < n p<\sqrt{n} p<n,所以算法的期望复杂度是 O ( n 1 4 ) O(n^{\frac{1}{4}}) O(n41)。需要强调的是,这是一个启发式算法,不是严格的数学证明。算法可能会失败,此时是因为子集 X X X中没有碰撞,需要重新选择初始值 x 1 x_{1} x1或者选择不同的函数 f f f产生新的集合 X X X。
例子
设 n = 7171 n=7171 n=7171, f ( x ) = x 2 + 1 , x 1 = 1 f(x) = x^{2}+1,x_{1} = 1 f(x)=x2+1,x1=1。
计算集合
X
X
X部分元素为:
1
,
2
,
5
,
26
,
677
,
6557
,
4105
,
6347
,
4903
,
2218
,
219
,
4936
,
4210
,
4560
,
4872
,
375
,
4377
,
4389
,
2016
,
5471
,
88
⋯
1,2,5,26,677,6557,4105,6347,4903,2218,219,4936,4210,4560,4872,375,4377,4389,2016,5471,88\cdots
1,2,5,26,677,6557,4105,6347,4903,2218,219,4936,4210,4560,4872,375,4377,4389,2016,5471,88⋯
计算
g
c
d
(
x
1
−
x
2
,
n
)
,
g
c
d
(
x
2
−
x
4
,
n
)
,
g
c
d
(
x
3
−
x
6
,
n
)
⋯
gcd(x_{1}-x_{2},n),gcd(x_{2}-x_{4},n),gcd(x_{3}-x_{6},n)\cdots
gcd(x1−x2,n),gcd(x2−x4,n),gcd(x3−x6,n)⋯,发现
g
c
d
(
x
11
−
x
22
,
n
)
=
71
gcd(x_{11}-x_{22}, n)=71
gcd(x11−x22,n)=71,这就找到了
n
n
n的一个因子71。
参考书籍:Stinson D , 斯廷森, 冯登国. 密码学原理与实践[M]. 电子工业出版社, 2009.