再探扩展欧几里得算法的应用

今天又复习了一下扩展欧几里得算法，并且了解了一些扩展欧几里得算法可以解决的新东西。。。

对欧几里得算法一点都不会的请移驾本空间另一篇文章：传送门。

下面以三个题为例总结相关知识：

POJ2142 The Balance

这道题要求用扩展欧几里得算法求解a*x+b*y=c的一组解，要求此解满足|x|+|y|最小，|x|+|y|相同时a*|x|+b*|y|最小。

不妨设a>b，若a<b，交换a、b即可。

首先令d=gcd(a,b)，用扩展欧几里得求解a*x+b*y=d，然后令x=x*(c/d)，y=y*(c/d)，就求出了方程的任意一组解(x,y)。

根据相关定理，方程的全部解可以表示为：x+i*(b/d)，y-i*(a/d)，其中i为常数。观察函数：|x+i*(b/d)|+|y-i*(a/d)|，当y-i*(a/d)在0附近时，函数取最小值，证明略。因此我们只需要枚举y-i*(a/d)在0附近的两个整数解，二者中选更优的即可。

还有一种讲起来更通俗的方法：|x|+|y|取最小值，要么x为方程的最小正整数解，要么y为方程的最小正整数解。不过此法我没有写程序提交。

#include<iostream>

using namespace std;

int a,b,c,d,x,y,t1,t2,x1,y1,x2,y2,flag;

 

int abso(int x)

{

return x>0?x:-x;

}

 

int exgcd(int a,int b,int &x,int &y)

{

if(!b) {x=1; y=0; return a;}

int c=exgcd(b,a%b,x,y),z=x;

x=y; y=z-y*(a/b);

return c;

}

 

int main()

{

while(cin>>a>>b>>c&&a)

{

if(a<b) {swap(a,b); flag=1;} else flag=0;

d=exgcd(a,b,x,y);

x*=c/d; y*=c/d;

t1=y*d/a;

if(y-a/d*t1<0) t1--;

t2=t1+1;

x1=abso(x+b/d*t1); y1=abso(y-a/d*t1);

x2=abso(x+b/d*t2); y2=abso(y-a/d*t2);

if(x1+y1<x2+y2||(x1+y1==x2+y2&&a*x1+b*y1<a*x2+b*y2)) {x=x1; y=y1;}

else {x=x2; y=y2;}

if(flag) cout<<y<<' '<<x<<endl; else cout<<x<<' '<<y<<endl;

}

return 0;

}

 

POJ2115 C Looooops 此题就是求一次同余方程c*a mod p=b的最小正整数解a。令d=gcd(c,p)，若d|b，方程有解，并且所有的解关于p/d同余；否则无解。做法是用扩展欧几里得算法求方程c*x+p*y=d的一组解，则原方程的一组解为x*(b/d)。根据解关于p/d同余的性质可以取模求出最小正整数解。

#include<iostream>

using namespace std;

typedef long long LL;

LL f[35],a,b,c,k,p,d,i,x,y;

LL exgcd(LL a,LL b,LL &x,LL &y)

{

if(!b) {x=1; y=0; return a;}

LL c=exgcd(b,a%b,x,y),z=x;

x=y; y=z-y*(a/b);

return c;

}

int main()

{

for(i=1,f[0]=1;i<=32;i++) f[i]=f[i-1]*2;

while(cin>>a>>b>>c>>k&&k)

{

p=f[k];

b=(b-a+p)%p;

d=exgcd(c,p,x,y);

if(b%d==0) cout<<(x*b/d%p+p)%(p/d)<<endl; else cout<<"FOREVER\n";

}

return 0;

}

POJ2891 Strange Way to Express Integers 求最小正整数T，满足T mod ai=bi。这道题看第一眼：裸的中国剩余定理。再一看，没保证所有ai互质。。。先说说如果ai互质，用中国剩余定理的做法：

记所有ai的最小公倍数为mul。对于每个ai，求出一个ci，使得ci能被数列a中除了ai以外的所有数整除，并且ci mod ai=1。

ci可以这样求：首先用扩展欧几里得定理求出u*x+v*y=1的任意一组解(x,y)，其中u=mul div ai，v=ai，则u*x就是ci。可以这样理解：ci一定是u的倍数，因为它能被数列中除了ai的所有数整除；不妨设ci是u的x倍。由于ci mod ai=1,也就是说ci-1能被ai整除，所以ci-1是v的倍数，不妨设ci-1是v的y倍，所以ci-1=v*y，也就是u*x+v*y=1。

求出ci后，令di=ci*bi。把所有di相加就得到了一个满足要求的T值。由于让t最小，只需要mod mul即可，但是要注意mod出来负数时要处理一下。

但是这题并没有保证所有ai互质，因此我们使用扩展欧几里得算法求解。

对于第i个方程T mod ai=bi，假设之前求得的解为T，之前所有ai的最小公倍数为lcm，那么需要让T加上一个值 k*lcm，其中k为常数。并且使得：(T+k*lcm)mod ai=bi。这个方程就是上面刚讲过的一次线性同余方程了。若该方程无解，那么T不存在。

#include<iostream>

#include<cstdio>

using namespace std;

typedef long long LL;

int n,i,fail;

LL a,b,lcm,now,k,d,x,y;

 

LL exgcd(LL a,LL b,LL &x,LL &y)

{

if(!b) {x=1; y=0; return a;}

LL c=exgcd(b,a%b,x,y),z=x;

x=y; y=z-y*(a/b);

return c;

}

 

int main()

{

while(cin>>n)

{

scanf("%lld%lld",&a,&b);

lcm=a; now=b; fail=0;

for(i=1;i<n;i++)

{

scanf("%lld%lld",&a,&b);

b=(b-now%a+a)%a;

d=exgcd(lcm,a,x,y);

if(b%d==0) k=x*(b/d)%a; else fail=1;

now+=k*lcm;

lcm=lcm/d*a;

now=(now%lcm+lcm)%lcm;

}

if(fail) cout<<"-1\n"; else cout<<now<<endl;

}

return 0;

}