扩展欧几里得算法（详细推导+代码实现+应用）

1.扩展欧几里得算法

贝祖定理：若a，b是整数，且gcd(a,b)=d，那么对于任意的整数x，y，ax+by=m中的m一定是d的倍数。（特别地，如果a、b是整数，那么一定存在整数x、y使得ax+by=gcd(a,b)。）

那么贝祖定理的一个直接应用就是：如果ax+by=1有解，那么gcd(a,b)=1（将原公式两边同时除以gcd(a,b)）。

扩展欧几里得算法用来解决这样一个问题：给定两个非零的整数a和b，求一组整数解（x，y）使得ax+by = gcd(a,b)成立，其中gcd(a,b)表示a和b的最大公约数。（其中我们通过前面的贝祖定理可知解一定存在）

回忆我们知道的欧几里得算法，它总是把gcd(a,b)转化为求解gcd(b,a%b)，而当b变为0时返回a，此时的a就等于gcd。也就是说，欧几里得算法结束时变量a中存放的是gcd，变量b中存放的是0，因此此时显然有$a\times 1+b\times 0=gcd$成立，此时有x=1、y=0成立。

我们不妨利用上面的欧几里得算法的过程来计算x和y。目前已知的是递归边界成立时为x=1、y=0，需要想办法反推出最初始的x和y。

当计算gcd(a,b)时，有$a{x}_1+b{x}_1=gcd$成立；

而在下一步计算gcd(b,a%b)时，又有$b{x}_2+(a\%b)y_2=gcd$成立。

因此$a{x}_1+b{y}_1=b{x}_2+(a\%b)y_2$成立。

又考虑到有关系$a\%b=a-(a/b)\times b$成立（此处除法为整除），因此$a{x}_1+b{y}_1=a{y}_1+b(x_2-(a/b)y_2)$

对比等号左右两边可以马上得到下面的递推公式：

$$\{\begin{array}{l}{x}_1=y_2\\ y_1=x_2-(a/b)y_2\end{array}$$

由此我们便可以通过$x_2$和$y_2$来进行反推出$x_1$和$y_1$了。代码如下：

int exGcd(int a, int b, int &x, int &y)   //x和y使用引用
{
    if(b == 0)
    {
        x = 1;
        y = 0;
        return a;
    }
    int g = exGcd(b, a%b, x, y);    //递归计算exGcd（b,a%b）
    int temp = x;             //存放x的值
    x = y;
    y = temp - (a/b)*y;       //更新y = x(old) - a/b*y(old)
    return g;                 //g是gcd
}

由于这里我们使用了引用，因此当exGcd函数结束时x和y就是所求的解。

显然，在得到这样一组解之后，就可以通过下面的式子得到全部解（其中K为任意整数）：

$$\{\begin{array}{l}{x}^{\prime }=x+\frac{b}{gcd}\times K\\ y^{\prime }=y-\frac{a}{gcd}\times K\end{array}$$

简单证明一下：

假设新的解为$x+s_1$、$y-s_2$，即有$a\ast (x+s_1)+b\ast (y-s_2)=gcd$成立，

通过代入$ax+by=gcd$可以得到$a{s}_1=b{s}_2$，

于是$\frac{s_1}{s_2}=\frac{b}{a}$成立。

为了让和$s_1$和$s_2$尽可能小，可以让分子和分母同时除以一个尽可能大的数，同时保证它们仍然是整数。显然，由于$\frac{b}{gcd}$与$\frac{a}{gcd}$互质，因此gcd是允许作为除数的最大数，于是$\frac{s_1}{s_2}=\frac{b}{a}=\frac{b/gcd}{a/gcd}$

得$s_1$和$s_2$的最小取值就是$\frac{b}{gcd}$与$\frac{a}{gcd}$。

证毕

也就是说，x和y的所有解分别以$\frac{b}{gcd}$与$\frac{a}{gcd}$为周期。

那么其中x的最小非负整数解是什么呢？直观来说就是$x\%\frac{b}{gcd}$。但是由于通过exGcd函数计算出来的x、y可以为负数，因此实际上$x\%\frac{b}{gcd}$会得到一个负数，例如$(-15)$。考虑到即便x是负数，$x\%\frac{b}{gcd}$的范围也是在($-\frac{b}{gcd}$,0)，因此对于任意整数来说，$(x\%\frac{b}{gcd}+\frac{b}{gcd})\%\frac{b}{gcd}$才是对应的最小非负整数解。

如果gcd=1，即ax+by=1时，全部解的公式简化为下式，且x的最小非负整数解也可以简化为$(x\%b+b)\%b$

$$\{\begin{array}{l}{x}^{\prime }=x+b\times K\\ y^{\prime }=y-a\times K\end{array}$$

2.方程ax+by=c的求解

当我们知道了ax+by=gcd的解后，最主要的应用就是：求解ax+by=c，其中c为任意整数。

首先，假设ax+by=gcd有一组解$(x_0,y_0)$，现在在其等号的左右两边同时乘以$\frac{c}{gcd}$，

即有$a\frac{c{x}_0}{gcd}+b\frac{c{y}_0}{gcd}=c$成立，因此$(x,y)=(\frac{c{x}_0}{gcd},\frac{c{y}_0}{gcd})$是ax+by=c的一组解。

但是显然这样做的充要条件是c%gcd==0（因为：$a\frac{x}{gcd}+b\frac{y}{gcd}=\frac{c}{gcd}$，其中左边为整数因此右边也为整数）

否则第一步将无法做到。

为了获取全部解的公式，我们可以模仿前面的做法，假设新的解为$x+s_1$，$y-s_2$，然后将$a\ast (x+s_1)+b\ast (y-s_2)=c$与$ax+by=c$联立，发现同样可得$\frac{s_1}{s_2}=\frac{b}{a}$成立。于是因为同样的原因，$s_1$和$s_2$的最小取值仍然是$\frac{b}{gcd}$与$\frac{a}{gcd}$。因此ax+by=c的全部解的公式为：

$$\{\begin{array}{l}{x}^{\prime }=x+\frac{b}{gcd}\ast K=\frac{c{x}_0}{gcd}+\frac{b}{gcd}\ast K\\ y^{\prime }=y-\frac{a}{gcd}\ast K=\frac{c{y}_0}{gcd}-\frac{a}{gcd}\ast K\end{array}$$

我们发现，这与ax+by=gcd全部解的公式是一样的，唯一不同的是初始解（x，y）不同。因此，对于ax+by=c来说，其解同样分别以$\frac{b}{gcd}$与$\frac{a}{gcd}$为周期。

但是我们要注意的是：在ax+by=gcd的全部解的基础上都乘以$\frac{c}{gcd}$只能获取一部分解。

这是因为：

由于$c\ge gcd$（根据c%gcd==0得到），因此ax+by=gcd的全部解乘以$\frac{c}{gcd}$会导致周期放大为原先的$\frac{c}{gcd}$倍（x和y的周期会边为$\frac{bc}{gc{d}^2}$和$\frac{ac}{gc{d}^2}$），所以导致漏解。

3.同余式$ax\equiv c(modm)$的解

既然已经解决了ax+by=c，不得不提的就是同余式$ax\equiv c(modm)$的求解。

什么式同余式：对于a、b、m来说，如果m整除a-b（即（a-b）%m=0），那么就说a与b模m同余，对应的同余式为$ax\equiv c(modm)$，m称为同余式的模。

那么我们解决同余式$ax\equiv c(modm)$的解呢？

根据同余式的定义，有$(ax-c)\%m=0$成立

因此存在整数y，使$ax-c=my$成立

移项并令y=-y后即得$ax+my=c$

由前面的结论我们知道，当$c\%gcd(a,m)==0$时方程才有解，且解的形式如下，

其中（x，y）是ax+my=c的其中一组解，可以先通过求解ax+my=gcd(a,m)得到$(x_0,y_0)$，然后由公式$(x,y)=(\frac{c{x}_0}{gcd(a,m)},\frac{c{y}_0}{gcd(a,m)})$直接得到。

$$\{\begin{array}{l}{x}^{\prime }=x+\frac{m}{gcd(a,m)}\times K\\ y^{\prime }=y-\frac{a}{gcd(a,m)}\times K\end{array}（K为任意整数）$$

虽然对方程ax+my=c来说，K可以取任意整数，但是对同余式来说会有很多解在模m意义下是相同的（由于只关心x，因此下面只考虑x）。对同余式来说没只需要找出那些在模m意义下不同的解。因此考虑$x^{\prime }=x+\frac{m}{gcd(a,m)}\ast K$，会发现当K分别取0、1、2、…、gcd(a,m)-1时，所得到的解在模m意义下是不同的，而其他解都可以对应到K取这gcd(a,m)个数值之一。由此可得到结论：

设a，c，m是整数，其中$m\ge 1$，则

①若$c\%gcd(a,m)\ne 0$，则同余式方程$ax\equiv c(modm)$无解。

②若$c\%gcd(a,m)=0$，则同余式方程$ax\equiv c(modm)$恰好由gcd(a,m)个模m意义下不同的解，且解的形式为：

$x^{\prime }=x+\frac{m}{gcd(a,m)}\ast K$

4.逆元的求解以及(b/a)%m的计算

接着我们来解决最后一个问题，假设a、m是整数，求a模m的逆元。

什么是逆元（此处特指乘法的逆元）：

假设a、b、m是整数，m>1，且有$ab\equiv 1(modm)$成立，

那么就说a和b互为模m的逆元，一般记作$a\equiv \frac{1}{b}$或$b\equiv \frac{1}{a}$

通俗的讲就是，如果两个整数的乘积模m后等于1，那么就称它们互为m的逆元

那么逆元有什么用处呢？

对乘法来说有$(b\ast a)\%m=((b\%m)\ast (a\%m))\%m$成立，但是对除法来说，$(b/a)\%m=((b\%m)/(a\%m))\%m$却不成立，$(b/a)\%m=((b\%m)/a)\%m$也不成立，

例如：如果要对12/4对2取模，采用((12%2)/4)%2会得到错误的结果0，而实际上结果是1。

这时就需要逆元来计算(b/a)%m:

通过找到a模m的逆元x，就有(b/a)%m=(b*x)%m成立（只考虑整数取模，也即假设b%a=0，即b是a的整数倍），于是就把除法取模转化为乘法取模，这对于解决被除数b非常大（使得b已经取过模，不是原始值）的问题来说是非常实用的。

由定义知，求a模m的逆元，就是求解同余式$ax\equiv 1(modm)$，并且在实际使用中，一般把x的最小正整数解称为a模m的逆元，因此下文中提到的逆元都默认为x的最小正整数解。显然，同余式$ax\equiv 1(modm)$是否有解取决于1%gcd(a,m)是否为0，而这等价于gcd(a,m)是否为1：

①如果$gcd(a,m)\ne 1$，那么同余式$ax\equiv 1(modm)$无解，a不存在模m的逆元。

②如果$gcd(a,m)=1$，那么同余式$ax\equiv 1(modm)$在(0,m)上有唯一解，可以通过求解ax+my=1得到。

注意：由于gcd(a,m)=1，因此ax+my=1=gcd(a,m)，直接使用扩展欧几里得算法解出x之后就可以用$(x\%m+m)$得到(0,m)范围内的解，也就是所需要的逆元。

下面代码使用了扩展欧几里得算法来求解a模m的逆元，使用条件是gcd(a,m)=1,当然如果m是素数，就肯定成立了，可以放心使用。

int inverse(int a,int m)
{
	int x, y;
	int g = exGcd(a, m, x, y);   //求解ax+my=1
    return (x%m+m)%m;            //a模m的逆元为(x%m+m)%m
}

补充说明：

如果m是素数，且a不是m的倍数，则还可以直接使用费马小定理来得到逆元，这种做法不需要使用扩展欧几里得算法。

费马小定理：设m是素数，a是任意整数且a不能整除m，则$a^{m-1}\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}m)$。

使用费马小定理来推到逆元的方法非常简单：由$a^{m-1}\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}m)$可知$a\ast a^{m-2}\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}m)$，直接由逆元的定义便可以知道，$a^{m-2}\%m$就是a模m的逆元，而这可以通过快速幂的方法来求解。