浅谈扩展欧几里得算法

什么是拓展欧几里得？简单的说，就是求关于x,y的方程 ax + by = gcd(a,b) 的所有整数解

现在我们来解决四个问题

什么是裴属定理，如何证明裴属定理？
怎么用扩展欧几里得来求ax + by = gcd(a,b) 的特解？
怎么求由特解推出其他的所有解？
扩展欧几里得的三大应用

一、裴蜀定理内容即证明

内容

设a, b是不全为零的整数，则存在整数x, y 使得 ax + by = gcd(a, b).

证明：

a, b中假如有一个为0，则上述定理完全吻合，如a = 0, 则gcd(a, b) = b , 显然成立
若a, b 都不等于0

设gcd(a, b) = d,则d != 0

等式两边同时除以d，得：

a₁x + b₁y = 1

再利用辗转相除法：gcd(a, b) --> gcd(b, a % b) --> ……

设模出来的余数叫做x, 则有：

gcd(a₁, b₁) = gcd(b₁, r₁) = gcd(r₁, r₂) = …… = gcd(r_n-1, r_n) = 1

把辗转相除法中的运算展开，做成带余数的除法，得
- a₁ = x₁b + r₁
- b₁ = x₂r₁ + r₂
- r₁ = x₃r₂ + r₃
  
  ……
- r_n-3 = x_n-1r_n-2 + r_n-1
- r_n-2 = x_nr_n-1 + r_n
- r_n-1 = x_n+1r_n
不妨令辗转相除法在除到互质的时候退出则 r_n = 1 , 由倒数第二行得
- r_n-2 = x_nr_n-1 + 1 —> 1 = r_n-2 - x_nr_n-1
由倒数第三行得：**r_n-1 = r_n-3 - x_n-1r_n-2**带入上式

得：1 = (1 + x_nx_n-1)r_n-2 - x_nr_n-3

采取同样的方法，逐个消去r_n-2,……, r₁

最终会得到1 = f(x)a₁ + g(x)b₁

其中f(x)和g(x) 是关于x的函数

设f(x) = x, g(x) = y，则可以得到1 = a₁x + b₁y, 得证!

正常人都不喜欢看上面写的奇奇怪怪的推导，那就看下面我这个具体的推导叭，这里是假设到 r4 = 1，然后往前推导滴

二、求特解

我们都知道，欧几里得公式可以由这个式子表达

gcd(a, b) = gcd(b, a % b)

通过这个式子，我们可以不断递推到b = 0，此时a即为a和b的最大公约数

现在我们对这个式子进行展开，看看有什么好玩的嘛

gcd(a, b) = a * x1 + b * x2

gcd(b, a % b) = b * x2 + (a % b) * y2

其中a % b = a - (a / b) * b

故由第一行的欧几里得公式得：

a * x1 + b * y1 = b * x2 + {a - (a / b) * b}y2

化简得a * x1 + b * y1 = a * y2 + b * {x2 - (a / b) * y2}

根据待定系数法得到

x1 = y2
y1 = x2 - (a / b) * y2

也就是说，如果有了gcd(b, a % b)的解x2, y2就可以推出gcd(a, b)的解x1, y1

那么这就和求gcd的过程一样，一直推到最后x = 1， y = 0，就可以回溯回去，得到gcd(a, b)的特解x1, y1

特解其实是最小的一个解

二、怎么利用特解推出其他所有整数解

先说结论：

x=x0+k∗bgcd(a,b)x = x0 + k * \frac{b}{gcd(a,b)}x=x0+k∗gcd(a,b)b

y=y0−k∗agcd(a,b)y = y0 - k * \frac{a}{gcd(a,b)}y=y0−k∗gcd(a,b)a

任意的x + y = x0 + y0

x0, y0就是方程的特解

对于a * x + b * y = g来说，让x增加b/gcd,让y增加a / gcd

这样就相当于加a * b / gcd，减去a * b / gcd,一加一减，最终不变

三、扩展欧几里得三大应用

用扩展欧几里得算法解不定方程ax+by=c

如果 c % gcd(a, b) == 0，即c 是 gcd(a, b)的整数倍，则方程有解，且解就在上面的基础上乘以cgcd(a,b)\frac{c}{gcd(a, b)}gcd(a,b)c即可

void e_gcd(int a, int b, int &gcd, int &x, int &y)
{if (b == 0){x = 1;y = 0;gcd = a;}else{e_gcd(b, a % b, gcd, y, x);y -= x * (a / b);}
}int main(){int a, b, x, y, gcd;cin>>a>>b;e_gcd(a, b, gcd, x, y);cout<<x<<' '<<y<<endl;cout<<gcd<<endl;return 0;
}

用扩展欧几里得算法求解模线性方程ax≡b (mod n)

对于模线性方程ax≡b (mod n)可以化简为ax + ny = b，设d = gcd(a, n) 当且仅当b % d == 0时有解，且有d个解

当b % d == 0时，设ax + ny = d,的特解为x,y

等式两边同时乘以bd\frac{b}{d}db，则方程变成axbd+nybd=bax\frac{b}{d} + ny\frac{b}{d} = baxdb+nydb=b

则原方程ax + ny = b的解为**x′=x∗bd,y′=y∗bdx' = x*\frac{b}{d},y' = y*\frac{b}{d}x′=x∗db,y′=y∗db**

故**ax≡b (mod n)**的特解为 x0 = x * (b / d) % n

d个解为xi= x0 + i*(n/d) mod n，{i = 0……n-1}

解的间隔为n/d

void exgcd(int a, int b, int &gcd, int &x, int &y)
{if (b == 0){x = 1;y = 0;gcd = a;}else{exgcd(b, a % b, gcd, y, x);y -= x * (a / b);}
}bool modular_linear_equation(int a, int b, int n){int x, y, x0, gcd;exgcd(a, n, gcd, x, y);int d = gcd;//d为gcd(a, n)if(b % d)return false;//无解x0 = x * (b/d) % n;//特解for(int i = 0; i < d; ++i)cout<<x0 + i * (n / d)<<endl;//所有的解return true;
}

用欧几里德算法求乘法逆元：

什么叫乘法逆元？

ax ≡ 1 (mod p)

这里我们称x为a关于p的逆元

上述式子可以还原为：ax + py = 1,根据上面讲的扩展欧几里得算法，我们知道gcd(a, p) = 1时才有解，通过扩展欧几里得算法得到的解x₀,利用x₀%p就可以得到最小解

为什么呢？

因为通解为x = x₀ + p * k

那么，也就是说， a 关于 p 的逆元是一个关于 p 同余的，那么根据最小整数原理，一定存在一个最小的正整数，它是 a 关于p 的逆元，而最小的肯定是在（0 , p）之间的，而且只有一个，这就好解释了。

但是，由于问题的特殊性，有时候我们得到的特解 x₀ 是一个负数，还有的时候我们的 p 也是一个负数这怎么办？

当 p 是负数的时候，我们取 p 的绝对值就行了，当 x₀ 是负数的时候，x₀% p 的结果仍然是负数（在计算机计算的结果上是这样的，虽然定义的时候不是这样的），这时候，我们仍然让 x₀ 对abs§ 取模，然后结果再加上abs§ 就行了，于是，我们不难写出下面的代码求解一个数 a 对于另一个数 p 的乘法逆元：

int a, b, x, y, gcd, ans, p;void exgcd(int a, int b, int &gcd, int &x, int &y)
{if (b == 0){x = 1;y = 0;gcd = a;}else{exgcd(b, a % b, gcd, y, x);y -= x * (a / b);}
}inline int cal(int a, int p){exgcd(a, p, gcd, x, y);if(1 % gcd != 0)return -1;x = x * 1 / gcd;p = abs(p);ans = x % p;if(ans <= 0)ans += p;return ans;
}