最大公约数

暴力求解法

int gcd(int x,int y)
{
    if(!x)
        return y;
    if(!y)
        return x;
    int result = x > y?y:x;
    while(x%result != 0 || y%result != 0)
    {
        --result;
    }
    return result;
}

辗转相除法

简介

又称欧几里得算法(Euclid's Algorithm)

公式

( gcd(x,y)= egin{cases} x& ext{y=0}\ gcd(y,x \% y)& ext{y>0} end{cases} )

c++代码

int Euclid(int x,int y)
{
    if(y)
        Euclid(y,x % y);
    else
        return x;
}

证明

符号说明

gcd(x,y)：x和y的最大公约数
cd(x,y)：x和y所有公约数的集合
d(x)：x所有公约数的集合

即证(gcd(x,y) = gcd(y,x\%y))

1. 当(x=y)时，左边(= y = gcd(y,0) =)右边
2. 当(x<y)时，左边(= gcd(x,y) = gcd(y,x) =)右边
3. 当(x>y)时，不妨设(x = ky + r(k,rin N^+))，则(r=x\%y)
   • (forall tin cd(x,y))，对(r = x - ky)，两边÷t,得(frac{r}{t}=frac{x}{t}-kfrac{y}{t} = m)，易得(min N^+)，所以(tin c(r))，又(tin c(y))，所以(tin cd(y,r) = cd(y,x\%y)),所以(cd(x,y) subseteq cd(y,x\%y))
   • (forall t in cd(y,x\%y))，对(x = ky + r)，两边÷t,得(frac{x}{t}=kfrac{y}{t}+kfrac{r}{t} = m)，易得(min N^+)，所以(t in c(x)))，又(t in c(y))，所以(t in cd(x,y)),所以(cd(y,x\%y) subseteq cd(x,y))
   • 综上，(cd(x,y) = cd(y,x\%y))，故(gcd(x,y) = gcd(y,x\%y))

第三步的另一种证法
设(t=gcd(x,y),x=mt,y=nt)，则(r=x-ky=mt-knt=(m-kn)t)，故(t in c(r))
下证(m-kn ext{与} n ext{互质})
反证法：若(m-kn ext{与} n ext{不互质})，则存在因子w，使得(m-kn=aw,n=bw)，由(egin{cases}x=mt=(aw+kn)t=(aw+kbw)t=(a+kb)tw\y=nt=bwtend{cases})得(gcd(x,y)=tw)，矛盾
又(egin{cases}r=(m-kn)t\y=ntend{cases})，故(t=gcd(y,r))
综上(gcd(x,y) = gcd(y,r) = gcd(y,x\%y))

非递归的c++代码

通过上面的式子和证明可以知道，(gcd(x,y))不断保证(x>=y)

int gcd(int x,int y)
{
    int r;
    if(x < y)
    {
        r = x;
        x = y;
        y = r;
    }
    do
    {
        r = x % y;
        x = y;
        y = r;
    }while(r);
    return x;
}

缺陷

欧几里德算法是计算两个数最大公约数的传统算法，无论从理论还是从实际效率上都是很好的。但是却有一个致命的缺陷，这个缺陷在素数比较小的时候一般是感觉不到的，只有在大素数时才会显现出来。
一般实际应用中的整数很少会超过64位（当然现在已经允许128位了），对于这样的整数，计算两个数之间的模是很简单的。对于字长为32位的平台，计算两个不超过32位的整数的模，只需要一个指令周期，而计算64位以下的整数模，也不过几个周期而已。但是对于更大的素数，这样的计算过程就不得不由用户来设计，为了计算两个超过64位的整数的模，用户也许不得不采用类似于多位数除法手算过程中的试商法，这个过程不但复杂，而且消耗了很多CPU时间。对于现代密码算法，要求计算128位以上的素数的情况比比皆是，设计这样的程序迫切希望能够抛弃除法和取模。

更相减损术

简介

来源于《九章算术》

可半者半之，不可半者，副置分母、子之数，以少减多，更相减损，求其等也。以等数约之。

白话文
第一步：任意给定两个正整数；判断它们是否都是偶数。若是，则用2约简；若不是则执行第二步。
第二步：以较大的数减较小的数，接着把所得的差与较小的数比较，并以大数减小数。继续这个操作，直到所得的减数和差相等为止。
则第一步中约掉的若干个2与第二步中等数的乘积就是所求的最大公约数。
其中所说的“等数”，就是最大公约数。求“等数”的办法是“更相减损”法。

公式

不考虑第一步，及x,y不都为偶数
( gcd(x,y)= egin{cases} x& x = y\ gcd(y,x - y)& x > y\ gcd(y,x)& x < y end{cases} )

公式证明(x>y)

方法1：
(forall tin cd(x,y))，设(h = x - y)，两边÷t,得(frac{h}{t}=frac{x}{t}-kfrac{y}{t} = m)，易得(min N^+)，所以(tin c(x-y))，又(tin c(y))，所以(tin cd(y,x-y)),所以(cd(x,y) subseteq cd(y,x-y))
(forall t in cd(y,x-y))，对(x = y + h)，两边÷t,得(frac{x}{t}=kfrac{y}{t}+kfrac{h}{t} = m)，易得(min N^+)，所以(t in c(x))，又(t in c(y))，所以(t in cd(x,y)),所以(cd(y,x-y) subseteq cd(x,y))
综上，(cd(x,y) = cd(y,x-y))，故(gcd(x,y) = gcd(y,x\%y))
方法2:
设(t=gcd(x,y),x=mt,y=nt)，则(r=x-y=mt-nt=(m-n)t)，故(t in c(x-y))
下证(m-n ext{与} n ext{互质})
反证法：若(m-n ext{与} n ext{不互质})，则存在因子w，使得(m-n=aw,n=bw)，由(egin{cases}x=mt=(aw+n)t=(aw+bw)t=(a+b)tw\y=nt=bwtend{cases})得(gcd(x,y)=tw)，矛盾
又(egin{cases}x-y=(m-n)t\y=ntend{cases})，故(t=gcd(y,r))
综上(gcd(x,y) = gcd(y,x-y))

c++代码

int gcd(int x, int b)
{
    int times = 0;
    while(!(x&1) && !(y&1))
    {
        ++times;
        x >>= 1;
        y >>= 1;
    }
    while(x != y)
    {
        if(x > y)
            x -= y;
        else
            y -= x;
    }
    return x<<times;
}

缺陷

更相减损术和辗转相除法的主要区别在于前者所使用的运算是“减”，后者是“除”。从算法思想上看，两者并没有本质上的区别，但是在计算过程中，如果遇到一个数很大，另一个数比较小的情况，可能要进行很多次减法才能达到一次除法的效果，从而使得算法的时间复杂度退化为O(N)，其中N是原先的两个数中较大的一个。相比之下，辗转相除法的时间复杂度稳定于O(logN)

Stein算法

简介

与更相减损术相似，差别在于:
在更相减损法中，若两个是偶数则同除以2，结果乘以2。如果增加一个判断，若为一奇一偶则偶数除以2，结果不变，若为两个奇数才相减，这样就变成了目前计算大整数最大公约数的非常好的一个算法。

证明

若能证明当k与b互为质数，gcd(kx,y)=gcd(x,y)，便可推出一般情况，即当k=2时，计算一个偶数和一个奇数的最大公约数时，可以先将偶数除以2。
设(t=gcd(x,y),kx=mt,y=nt,x=frac{mt}{k})
下证(frac{m}{k})与(n)互质
反证法：设(frac{m}{k}=pw,n=qw)，由(egin{cases}kx=kpwt\y=qwtend{cases})，得(gcd(x,y)=wt)，矛盾
由(egin{cases}x=frac{m}{k}t\y=ntend{cases})，得(t=gcd(x,y))
所以(gcd(kx,y)=gcd(x,y))

代码实现

int gcd(int x, int y)
{
    if(!x)
        return y;
    if(!y)
        return x;
    int times = 0;
    while(!(x&1) && !(y&1))
    {
        ++times;
        x >>= 1;
        y >>= 1;
    }
    while(!(x&1))
        x >>= 1;
    while(!(y&1))
        y >>= 1;
    while(x != y)
    {
        if(x > y)
            x -= y;
        else
            y -= x;
    }
    return x<<times;
}

优化版：优化的根据，两个奇数的差又是偶数，可以再除以2

void divide2(int &n)
{
    while(!(n&1))
        n >>= 1;
}
int gcd(int x, int y)
{
    if(!x)
        return y;
    if(!y)
        return x;
    int times = 0;
    while(!(x&1) && !(y&1))
    {
        ++times;
        x >>= 1;
        y >>= 1;
    }
    divide2(x);
    divide2(y);
    while(x != y)
        if(x > y)
        {
            x -= y;
            divide2(x);
        }
        else
        {
            y -= x;
            divide2(y);
        }
    return x<<times;
}