[数据结构与算法]欧几里德算法及其扩展

欧几里德算法其实就是辗转相除法，求最大公因数。

具体做法是：用较大数除以较小数，再用所得的余数（第一余数）去除除数，再用所得的余数（第二余数）去除第一余数，如此反复，直到最后余数是0为止。最后的除数就是这两个数的最大公约数。

记gcd（a，b）为整数a和b的最大公因数。

public static int gcd(int m,int n){
    return n == 0 ? m :gcd(n,m%n);
}

欧几里德算法的扩展——裴蜀（贝祖）等式

（最下方有相关知识点的代码实现模板，请自行选择是否越过文字解释）

对任何整数a,b和它们的最大公约数d，关于未知数x和y的线性方程（裴蜀等式）：

ax+by=m（m为d的倍数）当且仅当gcd（a,b）=d时有整数解

裴蜀等式有解时必然有无穷多个整数解，每组解x、y都称裴蜀数

特别地，方程ax+by=1有整数解当且仅当整数a,b互素。

x，y可以通过扩展欧几里德求解：

欧几里德算法停止的状态是：a=gcd ,b=0 ,那么，a'x+b'y=gcd 此时x=1,y为任意数。

因为，这时候只要a=gcd的系数是1，那么只要b的系数是0或者其他值（无所谓是多少，反正任意数乘以0都等于0但是a 的系数一定是1），这时，我们可以得到：a*1+b*0=gcd。

当然这是最终状态，但是我们可以以此反推到最初状态。

假设当前我们要处理的是求出a和b的最大公约数，并求出x和y使得ax+by=gcd，而我们已经求出了下一个状态：b和a%b的最大公约数，并求出了一组X1和Y1使得:b*X1+(a%b)*Y1=gcd,那么这两个相邻的状态之间是否存在一种关系呢？

令a%b=k，==> a = b*(a/b)+k 即 a%b?= a - (a/b)*b

那么:?gcd=b*X1 + (a-(a/b)*b)*Y1

? ?? ?? ? ? ? =b*X1 + a*Y1 - (a/b)*b*Y1

? ?? ?? ? ? ? =a*Y1 + b*(X1 - a/b*Y1)

对比 a*Y1 + b*(X1 - a/b*Y1) 和 ax+by=gcd ，求一组x和y使得ax+by=gcd

由此得到递推公式：

? ? ? ? x=Y1

? ? ? ? y=X1-a/b*Y1

但是需要注意的是，根据这个递推所求的x是ax+by=d的解，而扩展欧几里德算法就是在求a,b的最大公约数d=gcd(a,b)的同时，求出贝祖等式ax+by=m的一个解（X0,Y0）

所以实际上，X0=x*m*（d^（-1））

通项：（t是一个倍数）

? ? ? ? x=X0 + (b/gcd) * t? ? ? ? （所有的x对b同模）? ? ? ? 求得的x1、x2、x3、x4之间是k倍的（b/d）

? ? ? ? y=Y0 - (a/gcd）* t? ? ? ? ?(所有的y对a同模）????????求得的y1、y2、y3、y4之间是k倍的（a/d)

相关知识点代码实现模板如下：

/**
*扩展欧几里德算法
*x=y1
*y=x1-a/b*y1
*/
public class Main{
    static long x;
    static long y;
    public static long gcd(long m,long n){
        return n == 0 ? m : gcd(n,m%n);
    }
    //最小公倍数
    public staic long min(long a,long b){
        return a*b/gcd(a,b);
    }
    /**
    *扩展欧几里德
    *调用完成后x y是ax+by=gcd（a，b）的解
    */
    public static long ext_gcd(long a,long b){
        if(b==0){
            x=1;
            y=0;
            return a;
        }
        long res=ext_gcd(b,a%b);
        long x1=x;
        x=y;
        y=x1-a/b*y;
        return res;
    }
    /**
    *线性方程
    *ax+by=m 当m是gcd（a，b）倍数时有解
    *等价于ax=m mod b
    */
    public static long linearEquation(long a,long b,long m)throws Exception{
        long d=ext_gcd(a,b);
        //m不是gcd（a，b）的倍数，这个方程无解
        if（m%d！=0）throw new Exception（"无解"）；
        long n=m/d;
        x*=n;
        y*=n;
        return d;
    }
    /**
    *求逆元
    *ax = 1 (%mo),gcd(a,mo)=1
    *ax+mo*y=1
    */
    public static long inverseElement(long a,long mo) throws Exception{
        long d=linearEquantion(a,mo,1);//ax+mo*y=1
        x=(x%mo+mo)%mo;//保证x>0
        return d;
    }
}

?如果想要得到x大于0的第一个解：

b/=d;? ? x = (x0%b+b)%b