Extended Euclidean Algorithm
一、欧几里得算法
在整数域中,存在如下事实:
( a , b ) = ( a − k b , b ) , k ∈ Z (a, b) = (a - kb, b), k \in \Z (a,b)=(a−kb,b),k∈Z
因此在求解两整数 a , b ( a > b ) a, b(a > b) a,b(a>b)的最大公约数时,令上式中 k = ⌊ a b ⌋ k = \lfloor \frac{a}{b} \rfloor k=⌊ba⌋,有:
( a , b ) = ( a − ⌊ a b ⌋ b , b ) = ( a m o d b , b ) (a, b) = (a - {\lfloor \frac{a}{b} \rfloor}b, b) = (a \thinspace mod \thinspace b, b) (a,b)=(a−⌊ba⌋b,b)=(amodb,b)
借此便缩小了所求问题的规模.分别用 b , a m o d b b, a \thinspace mod \thinspace b b,amodb替代 a , b a, b a,b,将这一操作重复地进行下去,直到 b = 0 b = 0 b=0,此时 ( a , 0 ) = a (a, 0) = a (a,0)=a,其中的 a a a便是最初两整数的最大公约数.
二、扩展欧几里得算法
扩展欧几里得算法用于求解下述问题:
对于给定的 a , b ∈ Z a, b \in \Z a,b∈Z,找到一组对应的 x , y ∈ Z x, y \in \Z x,y∈Z,使得 a x + b y = ( a , b ) ax + by = (a, b) ax+by=(a,b)
1.递归实现
(1)递归关系式
根据裴蜀定理 (Bézout’s lemma) ,存在整数 x , y , x ′ , y ′ x, y, x', y' x,y,x′,y′,使得:
{ a x + b y = ( a , b ) b x ′ + ( a − ⌊ a b ⌋ b ) y ′ = ( b , a − ⌊ a b ⌋ b ) = ( b , a m o d b ) \begin{cases} ax + by = (a, b) \\ bx' + (a - {\lfloor \frac{a}{b} \rfloor}b)y' = (b, a - {\lfloor \frac{a}{b} \rfloor}b) = (b, a \thinspace mod \thinspace b) \end{cases} {ax+by=(a,b)bx′+(a−⌊ba⌋b)y′=(b,a−⌊ba⌋b)=(b,amodb)
又由于 ( b , a m o d b ) = ( a , b ) (b, a \thinspace mod \thinspace b) = (a, b) (b,amodb)=(a,b),即两式右侧相等.得到:
a x + b y = b x ′ + ( a − ⌊ a b ⌋ b ) y ′ ax + by = bx' + (a - {\lfloor \frac{a}{b} \rfloor}b)y' ax+by=bx′+(a−⌊ba⌋b)y′
整理得到:
a ( x − y ′ ) + b ( y − ( x ′ − ⌊ a b ⌋ y ′ ) ) = 0 a(x - y') + b(y - (x'- {\lfloor \frac{a}{b} \rfloor}y')) = 0 a(x−y′)+b(y−(x′−⌊ba⌋y′))=0
上式对于任意 a , b a, b a,b均成立,因此:
{ x = y ′ y = x ′ − ⌊ a b ⌋ y ′ \begin{cases} x = y'\\ y = x' - {\lfloor \frac{a}{b} \rfloor}y' \end{cases} {x=y′y=x′−⌊ba⌋y′
只要知道了 x ′ , y ′ x', y' x′,y′,就能得到 x , y x, y x,y,借此便缩小了求解问题的规模.
(2)终止条件
每次递归时,分别用 b , a m o d b b, a \thinspace mod \thinspace b b,amodb替代 a , b a, b a,b,直到 b = 0 b = 0 b=0.此时 x = 1 , y = 0 x = 1, y = 0 x=1,y=0.
代码详见Extended Euclidean Algorithm.py
.
2.迭代实现
递归关系用矩阵表示:
( x y ) = ( 0 1 1 − d 1 ) ( x ′ y ′ ) \begin{pmatrix} x \\ y \end{pmatrix} = \begin{pmatrix} 0 & 1 \\ 1 & -d_1 \end{pmatrix} \begin{pmatrix} x' \\ y' \end{pmatrix} (xy)=(011−d1)(x′y′)
因此:
( x y ) = ( 0 1 1 − d 1 ) ( 0 1 1 − d 2 ) . . . ( 0 1 1 − d n ) ( 1 0 ) \begin{pmatrix} x \\ y \end{pmatrix} = { \begin{pmatrix} 0 & 1 \\ 1 & -d_1 \end{pmatrix} } { \begin{pmatrix} 0 & 1 \\ 1 & -d_2 \end{pmatrix} } ... \begin{pmatrix} 0 & 1 \\ 1 & -d_n \end{pmatrix} { \begin{pmatrix} 1 \\ 0 \end{pmatrix} } (xy)=(011−d1)(011−d2)...(011−dn)(10)
d i d_i di是第 i i i轮迭代时的 ⌊ a b ⌋ {\lfloor \frac{a}{b} \rfloor} ⌊ba⌋.
代码详见Extended Euclidean Algorithm.py
.
3.手工计算
记 a x + b y = d ax + by = d ax+by=d,分别令 x = 1 , y = 0 ; x = 0 y = 1 x = 1, y = 0; x = 0 y = 1 x=1,y=0;x=0y=1:
1 ⋅ a + 0 ⋅ b = a ( 1 ) 0 ⋅ a + 1 ⋅ b = b ( 2 ) 1·a + 0·b = a \quad (1) \\ 0·a + 1·b = b \quad (2) 1⋅a+0⋅b=a(1)0⋅a+1⋅b=b(2)
利用 ( 1 ) ( 2 ) (1)(2) (1)(2)两式,计算 ( 1 ) − ⌊ a b ⌋ ( 2 ) (1) - {\lfloor \frac{a}{b} \rfloor}(2) (1)−⌊ba⌋(2):
a − ⌊ a b ⌋ b = a m o d b ( 3 ) a - {\lfloor \frac{a}{b} \rfloor}b = a \thinspace mod \thinspace b \quad (3) a−⌊ba⌋b=amodb(3)
类似地,利用 ( 2 ) ( 3 ) (2)(3) (2)(3)两式,计算 ( 2 ) − ⌊ b a m o d b ⌋ ( 3 ) (2) - {\lfloor \frac{b}{a \thinspace mod \thinspace b} \rfloor}(3) (2)−⌊amodbb⌋(3)
反复执行这一操作,直到等式右侧为 ( a , b ) (a, b) (a,b)
三、在GF( p )中求乘法逆元
在 G F ( p ) GF(p) GF(p)中求 b b b的乘法逆元 y = b − 1 m o d p y = b^{-1} \thinspace mod \thinspace p y=b−1modp,即是求解满足 p x + b y = 1 m o d p px + by = 1 \thinspace mod \thinspace p px+by=1modp的 y y y.
利用上述扩展欧几里得算法求解即得.
def recursive_ext_euclid(a: int, b: int):
# recursive
if a < b:
return recursive_ext_euclid(b, a)
if b == 0:
return 1, 0
else:
x, y = recursive_ext_euclid(b, a % b)
q = a // b
return y, (x - q * y)
def iterative_ext_euclid(a: int, b: int):
# iterative
if a < b:
return iterative_ext_euclid(b, a)
# diag(2)
m11, m12, m21, m22 = 1, 0, 0, 1
while True:
if b == 0:
return m11, m21
q = a // b
a, b = b, a % b
m11, m12, m21, m22 = m12, m11 - q*m12, m22, m21 - q*m22