欧几里得
感谢博主,文章参考该博文写的,里面还有一些详细的证明,此文已省https://www.cnblogs.com/frog112111/archive/2012/08/19/2646012.html
https://blog.csdn.net/azx59285/article/details/101093214
欧几里得算法O(logn)
欧几里德算法又称辗转相除法,用于计算两个整数a,b的最大公约数。
基本算法:设a=qb+r,其中a,b,q,r都是整数,则gcd(a,b)=gcd(b,r),即gcd(a,b)=gcd(b,a%b)。
实现方法:
int gcd(int a,int b)
{
if(b==0)
return a;
return
gcd(b,a%b);
}
扩展欧几里得算法
基本算法(裴蜀定理):对于不完全为 0 的非负整数 a,b,gcd(a,b)表示 a,b 的最大公约数,必然存在整数对 x,y ,使得 gcd(a,b)=ax+by。
求解x、y:
- 显然当 b=0,gcd(a,b)=a。此时 x=1,y=0;
- ab!=0 时
设 ax1+by1=gcd(a,b);
bx2+(a mod b)y2=gcd(b,a mod b);
根据朴素的欧几里德原理有 gcd(a,b)=gcd(b,a mod b);
则:ax1+by1=bx2+(a mod b)y2;
即:ax1+by1=bx2+(a-(a/b)*b)y2=ay2+bx2-(a/b)*by2;
根据恒等定理得:x1=y2; y1=x2-(a/b)*y2;
这样我们就得到了求解 x1,y1 的方法:x1,y1 的值基于 x2,y2.
上面的思想是以递归定义的,因为 gcd 不断的递归求解一定会有个时候 b=0,所以递归可以结束。
递归实现代码:(基于欧几里得代码写)
只求x和y:
void exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1;
y=0;
return;
}
int x1,y1;
exgcd(b,a%b,x1,y1);
x=y1;
y=x1-(a/b)*y1;
}
求x和y,并且求gcd(a,b):
int exgcd(int a,int b,int &x,int &y)
{
if(b==0) //递归结束跟欧几里得算法一样
{
x=1;
y=0;
return a;
}
int r=exgcd(b,a%b,x,y); //计算上一层的结果,得到上一层的x和y
int t=x;
x=y;
y=t-a/b*y;
return r;
}
题目:扩展欧几里得
刷题链接:https://www.acwing.com/problem/content/description/879/
AC代码:
#include<iostream>
using namespace std;
int n,a,b;
void exgcd(int a,int b,int &x,int &y)
{
if(b==0)
{
x=1; y=0; return ;
}
int x1,y1;
exgcd(b,a%b,x1,y1);
x=y1;
y=x1-(a/b)*y1;
return ;
}
int main()
{
cin>>n;
while(n--)
{
cin>>a>>b;
int x,y;
exgcd(a,b,x,y);
cout<<x<<" "<<y<<endl;
}
return 0;
}
扩展欧几里得算法的应用
(1)求解不定方程;
(2)求解模线性方程(线性同余方程);
(3)求解模的逆元;
求解不定方程
对于不定整数方程a * x+b * y = c,若c mod Gcd(a,b) =0,则该方程存在整数解,否则不存在整数解。
-
先求x * a+y * b = Gcd(a, b)的解
上面已经列出找一个整数解的方法,在找到x * a+y * b = Gcd(a, b)的一组解x0,y0后,x * a+y * b = Gcd(a, b)的其他整数解满足:
x = x0 + b/Gcd(a, b) * t
y = y0 - a/Gcd(a, b) * t (其中t为任意整数) -
再求a * x+b * y = c的解
只需将x * a+y * b = Gcd(a, b)的每个解乘上 c/Gcd(a, b) 即可。
在找到x * a+y * b = Gcd(a, b)的一组解x0,y0后,应该是得到x * a+y * b = c的一组解x1 = x0*(c/Gcd(a,b)) , y1 = y0*(c/Gcd(a,b))
x * a+y * b = c的其他整数解满足:
x = x1 + b/Gcd(a, b) * t
y = y1 - a/Gcd(a, b) * t(其中t为任意整数)
x 、y就是x * a+y * b = c的所有整数解.
相关证明可参考:http://www.cnblogs.com/void/archive/2011/04/18/2020357.html
用扩展欧几里得算法解不定方程ax+by=c;
代码如下:
bool linear_equation(int a,int b,int c,int &x,int &y)
{
int gcd=exgcd(a,b,x,y);
if(c%gcd) //c不能被gcd(a,b)整除,方程无解
return false;
int k=c/gcd;
x*=k;y*=k; //求得的只是其中一组解,其他的解通过列举t得到
return true; //表示方程有解
}
套用exgcd模板求得的是一组特殊解,但其实这一个方程式是有一个解系,在很多问题中是要你求得最小整数解
方法一、exgcd求得x和y后
x = x%b;(求出的是最小整数解,有可能是正数或者负数)
x=(x%b+b)%b;(求出的就是最小正整数解)
方法二、可以说这是求最小正整数的模板
LL e_gcd(LL a,LL b,LL&x,LL&y)
{
if(b==0)
{
x=1;
y=0;
return a;
}
int r=e_gcd(b,a%b,x,y);
LL t=x;
x=y;
y=t-a/b*y;
return r;
}
LL cal(LL a,LL b,LL c)
{
LL x,y;
LL gcd=e_gcd(a,b,x,y);
if(c%gcd) return -1;
x*=c/gcd; //转化为x*a+y*b=c的解
y*=c/gcd;
b/=gcd; //约去c后原来的b就变成了b/gcd;
if(b<0) b=-b; //如果b为负数就去绝对值
LL ans=x%b;
if(ans<0) ans+=b; //求最小正整数解
return ans;
}
这两种本质上没啥区别,只是在一些问题中a,b等系数可能为负,第一种需要预处理,而第二种则可以直接用
求解模线性方程(线性同余方程)
同余方程 ax≡b (mod n) 对于未知数x有解,当且仅当gcd(a,n) | b时方程有解(b是gcd(a,n)的倍数),方程有gcd(a,n)个解
求解方程ax≡b (mod n) 相当于求解方程 ax+ny=b (x,y为整数),ax=ny+b
- 由上述解方程ax+ny=b的方法,我们可以得到
ax≡b (mod n) 的最小正整数解为:x0 = (ans%n+n)%n - 接下来求其所有的解集
首先看一个简单的例子:5x=4(mod3)
解得x = 2,5,8,11,14…
由此可以发现一个规律,就是解的间隔是3.
那么这个解的间隔是怎么决定的呢?
如果可以设法找到第一个解,并且求出解之间的间隔,那么就可以求出模的线性方程的解集了.
我们设解之间的间隔为dx.
那么有
ax ≡ b(mod n);
a(x+dx) ≡ b(mod n);
两式子相减,得到:a* dx (mod n) = 0;
也就是说a * dx是a的倍数(显而易见嘛),同时也是n的倍数,即a * dx是a和n的公倍数,为了求出dx,我们应该求出a和n的最小公倍数,此时对应的dx是最小的(解相互之间的最小间隔)
设a和n的最大公约数为gcd,那么a和n的最小公倍数为(a*n)/gcd
即a *dx =a *n/gcd
所以dx = n/gcd;
因此解之间的间隔就求出来了,也就得到ax≡b (mod n)的解集了
这个证明其实就响应了上面对于不定方程通解的定义,二者是等价的,间隔其实就是n/gcd(a,n),即x = x0 + n/Gcd(a, n) * t ,t为0-Gcd(a, n) 的一个值
代码实现:
bool modular_linear_equation(int a,int b,int n)
{
int x,y,x0,i;
int d=exgcd(a,n,x,y); //求得ax+ny=gcd(a,n)的一个解
if(b%d)
return false; //不存在解x
x0=x*(b/d)%n; //得ax+ny=b的最小整数解,x*(b/d)是其中一个整数解,x*(b/d)%n是最小整数解
for(i=1;i<d;i++) //按照间隔输出所有的解,均保证在n以内
printf("%d\n",(x0+i*(n/d))%n);
return true;
}
求模的逆元
-
首先引入逆元的概念:
逆元是模运算的一个概念,我们通常说A是B模C的逆元,实际上是指A * B ≡ 1 (mod C) ,也就是说A与B的乘积模C的余数为1。可表示为A= B^(-1) mod C.
通常用来求A或者B,将其中一个变量看做a,进而另一个变量就能转化为上述的x,求解x
打个比方,7 模 11 的逆元,即:7^(-1) mod 11 = 8,这是因为 7 × 8 = 5 × 11 + 1,所以说 7 模 11 的逆元是 8 -
同余方程ax≡ b(mod n) ,如果gcd(a,n) = 1,即a和n互质,则方程只有唯一解。
在这种情况下,如果b==1,同余方程就是ax ≡ 1(mod n),gcd(a,n)=1;这时求出的x为a在模n情况下的逆元
对于同余方程 ax= 1(mod n ), gcd(a,n)= 1 的求解就是求解方程
ax+ ny= 1,x, y 为整数。这个可用扩展欧几里德算法求出,原同余方程的唯一解就是用扩展欧几里德算法得出的 x 。
经典例题:RSA解密技术
该题核心代码:
求pq两个质数直接循环即可
#include<iostream>
using namespace std;
#include<algorithm>
typedef long long LL;
LL n=1001733993063167141,d=212353,C=20190324;
//LL n=55,d=3,C=19; //小的测试数据
void exgcd(LL a,LL b,LL &x,LL &y)
{
if(b==0)
{
x=1;y=0;return ;
}
LL x1,y1;
exgcd(b,a%b,x1,y1);
x=y1;
y=x1-(a/b)*y1;
return ;
}
LL quickMul(LL a,LL b)
{
LL di=a,ans=0;
while(b)
{
if(b%2)
{
ans=(ans+di)%n;
}
di=(di+di)%n;
b>>=1;
}
return ans;
}
LL quickPower(LL a,LL b)
{
LL di=a,ans=1;
while(b)
{
if(b%2)
{
ans=quickMul(ans,di)%n; //注意这个地方二者乘积可能会爆出范围,因此需要用快速乘来边乘边取模
}
di=quickMul(di,di)%n;
b>>=1;
}
return ans;
}
int main()
{
for(LL p=1;p<=n/p;p++) //只枚举p然后直接得到q,可以减少循环层数
{
if(n%p) continue;
LL q=n/p;
LL b=(p-1)*(q-1);
if(__gcd(d,b)==1)
{
cout<<p<<" "<<q<<endl;
LL x,y;
exgcd(d,b,x,y);
LL e=(x%b+b)%b;
cout<<"e : "<<e<<endl;
cout<<quickPower(C,e)<<endl;
}
}
return 0;
}
竞赛解题方式:巧用python大数解题
https://blog.csdn.net/weixin_43914593/article/details/112979612
完结,在知识的海洋中遨游~加油吧