最大公约数
暴力求解法
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)= \begin{cases} x& \text{y=0}\\ gcd(y,x \% y)& \text{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)\)
- 当\(x=y\)时,左边\(= y = gcd(y,0) =\)右边
- 当\(x<y\)时,左边\(= gcd(x,y) = gcd(y,x) =\)右边
- 当\(x>y\)时,不妨设\(x = ky + r(k,r\in N^+)\),则\(r=x\%y\)
- \(\forall t\in cd(x,y)\),对\(r = x - ky\),两边÷t,得\(\frac{r}{t}=\frac{x}{t}-k\frac{y}{t} = m\),易得\(m\in N^+\),所以\(t\in c(r)\),又\(t\in c(y)\),所以\(t\in 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}=k\frac{y}{t}+k\frac{r}{t} = m\),易得\(m\in 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 \text{与} n \text{互质}\)
反证法:若\(m-kn \text{与} n \text{不互质}\),则存在因子w,使得\(m-kn=aw,n=bw\),由\(\begin{cases}x=mt=(aw+kn)t=(aw+kbw)t=(a+kb)tw\\y=nt=bwt\end{cases}\)得\(gcd(x,y)=tw\),矛盾
又\(\begin{cases}r=(m-kn)t\\y=nt\end{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)=
\begin{cases}
x& x = y\\
gcd(y,x - y)& x > y\\
gcd(y,x)& x < y
\end{cases}
\)
公式证明(x>y)
方法1:
\(\forall t\in cd(x,y)\),设\(h = x - y\),两边÷t,得\(\frac{h}{t}=\frac{x}{t}-k\frac{y}{t} = m\),易得\(m\in N^+\),所以\(t\in c(x-y)\),又\(t\in c(y)\),所以\(t\in 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}=k\frac{y}{t}+k\frac{h}{t} = m\),易得\(m\in 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 \text{与} n \text{互质}\)
反证法:若\(m-n \text{与} n \text{不互质}\),则存在因子w,使得\(m-n=aw,n=bw\),由\(\begin{cases}x=mt=(aw+n)t=(aw+bw)t=(a+b)tw\\y=nt=bwt\end{cases}\)得\(gcd(x,y)=tw\),矛盾
又\(\begin{cases}x-y=(m-n)t\\y=nt\end{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\),由\(\begin{cases}kx=kpwt\\y=qwt\end{cases}\),得\(gcd(x,y)=wt\),矛盾
由\(\begin{cases}x=\frac{m}{k}t\\y=nt\end{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;
}