扩展欧几里得
扩展欧几里得就是要解决ax+by=c的解得算法,先将原式转化成ax+by=gcd(a,b);因为辗转相除法gcd(a,b)=gcd(b,a%b),所以最后一定能转化成bx+0*y=gcd(a,b)=b;所以有一组解x=1,y=0(y为什么都行,反正系数为0,乘什么都为0)。因为转化的这个式子一定有解,所以只要c%gcd(a,b)==0,原式就一定有解,但一般都要求x最小的解,怎么求呢?我们先假设ax+by=gcd(a,b)有一组解是x1,y1.任取另一组解是(x2,y2),则ax1+by1=ax2+by2=gcd(a,b);变形得
a1(x1-x2)=b1(y1-y2),其中a1=a/gcd,b1=b/gcd;注意,此时a1和b1互素,因此x1-x2一定是b1的整数倍,设它为kb1,同理,y1-y2=ka1;最后发现这推导并没用到ax+by右边是什么。
最后可以得出结论对于ax+by=c的一组解(x0,y0),则它的任意整数解都可以写成(x0+kb1,y0+ka1)其中b1=b/gcd(a,b),a1=a/gcd(a,b).
这样我们就可以用(x0%b1+b1)%b1得到它的最小正整数解了。注意此时的x0是求出的x乘上c/gcd(a,b),因为你求出的数缩小了几倍。
现在再说说乘法逆元,乘法逆元是形如ax=1(modb),最小的正整数x称为amodb的逆元。先说一下,ax=c(modb)可以转化成ax+by=c的形式,为什么呢?转化一下形式,ax%b=c%b => ax+k1b=c+k2b => ax+(k1-k2)b=c=>ax+by=c;
所以ax=1(mod)b可以转化成ax+by=1;这样就可以用扩欧解了。
这有一个讲的非常详细的博客详细的博客。
下面贴几个写的题。
HDU 1576 A/B
要求(A/B)%9973,但由于A很大,我们只给出n(n=A%9973)(我们给定的A必能被B整除,且gcd(B,9973) = 1)。
Input
数据的第一行是一个T,表示有T组数据。
每组数据有两个数n(0 <= n < 9973)和B(1 <= B <= 10^9)。
Output
对应每组数据输出(A/B)%9973。
Sample Input
Sample Output
每组数据有两个数n(0 <= n < 9973)和B(1 <= B <= 10^9)。
2 1000 53 87 123456789
7922 6060
裸的扩欧。A+9973y=n;把A写成Bx就行了
include<iostream> #include<cstdio> #include<cstring> #include<algorithm> using namespace std; #define ll long long const int M=9973; ll n,B; ll ex_gcd(ll a,ll b,ll &x,ll &y){ if(b==0){ x=1;y=0;return a; } ll k=ex_gcd(b,a%b,x,y); ll tmp=x; x=y; y=tmp-a/b*y; return k; } ll work(){ ll x,y; ll d=ex_gcd(B,M,x,y); ll c=n/d; ll mod=M/d; x=(x*c%mod+mod)%mod; return x; } int main(){ int t; scanf("%d",&t); while(t--){ scanf("%lld%lld",&n,&B); printf("%lld\n",work()); } return 0; }
POJ 1061 青蛙的约会
直接贴链接吧博客
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define ll long long
ll x,y,m,n,L;
ll ex_gcd(ll a,ll b,ll &x1,ll &y1){
if(b==0){
x1=1;y1=0;return a;
}
ll k=ex_gcd(b,a%b,x1,y1);
ll tmp=x1;
x1=y1;
y1=tmp-a/b*y1;
return k;
}
ll work(){
ll x1,y1;
ll c=y-x;
ll d=ex_gcd(m-n,L,x1,y1);
if(c%d!=0) return -1;
return (x1*c/d%(L/d)+(L/d))%(L/d);
}
int main(){
while(scanf("%lld%lld%lld%lld%lld",&x,&y,&m,&n,&L)==5){
if(m==n) {
puts("Impossible");
continue;
}
if(m<n) swap(m,n),swap(x,y);
if(work()==-1) puts("Impossible");
else printf("%lld\n",work());
}
return 0;
}