存在整数对 x,y ,使得 gcd(a,b)=ax+by。
求解 x,y的方法的理解
设 a>b。
1,显然当 b=0,gcd(a,b)=a。此时 x=1,y=0;
2,a>b>0 时
设 ax
1+ by
1= gcd(a,b);
bx
2+ (a mod b)y
2= gcd(b,a mod b);
根据朴素的
欧几里德原理有 gcd(a,b) = gcd(b,a mod b);
则:ax
1+ by
1= bx
2+ (a mod b)y
2;
即:ax
1+ by
1= bx
2+ (a - [a / b] * b)y
2=ay
2+ bx
2- [a / b] * by
2;
也就是ax
1+ by1 == ay
2+ b(x
2- [a / b] *y
2);
根据恒等定理得:x
1=y
2; y
1=x
2- [a / b] *y
2;
这样我们就得到了求解 x
1,y
1 的方法:x
1,y
1 的值基于 x
2,y
2.
上面的思想是以递归定义的,因为 gcd 不断的递归求解一定会有个时候 b=0,所以递归可以结束。
c++实现:
- int exgcd(int a,int b,int &x,int &y){
- if (b==0){
- x=1,y=0;
- return a;
- }
- int q=gcd(b,a%b,y,x);
- y-=a/b*x;
- return q;
- }
或者
- 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);
- int t=x;x=y;y=t-a/b*y;
- return r;
- }
这样求出来的x0和y0是gcd(a,b)=ax+by的一个特解,通解可以用下面的方法表示:
x = x0 + (b/gcd)*t
y = y0 – (a/gcd)*t 证明如下:
首先把 a/gcd(a,b)化为a',b/gcd(a,b)化为b',这样ax+by=1,gcd(a,b)=1(注意这里a,b实际上写作a'和b',为简单起见省掉撇号)
由于
ax+by = ax0+by0
所以
a(x-x0) = -b(y-y0) = P
所以
a|P,b|P(“|”表示整除关系)
因为a,b互素,所以
ab|P
即
P=kab
所以
x-x0=kb,y-y0=-ka
整理得到
x=x0+kb,y=y0-ka
再由gcd(a,b)=ax+by的通解可以推导出更为一般的ax+by=c的通解,此处若c mod gcd(a,b)==0则有整数解,否则无解,证明如下:
由于我们知道,存在一组x与y使得a*x+b*y=gcd(a,b)。
将等式两边同时乘以整数k,即a*x*k+b*y*k=gcd(a,b)*k。如果c mod gcd(a,b)=f,则0<=f<gcd(a,b)。
那么可以令c=gcd(a,b)*k+f。这样一来,就有a*x*k+b*y*k+f=c。
若f 0,由于f<gcd(a,b)<=a<=b(假设a<=b),所以不存在f=a*m(m为整数),也就不存在a*(x*k+m)+b*y*k=c。也就是说,不存在a*x+b*y=c的整数解x与y。
所以f=0,即只有当c mod gcd(a,b)=0时,a*x+b*y=c有正整数解。得证。
而ax+by=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的所有整数解。
至于这个通解为什么是这样的证明,我也没搞太懂,个人认为可以类比上面求gcd(a,b)=ax+by的通解证明,还望各路数论大神指点。
扩展欧几里得定理常用于求解中国剩余定理相关问题和求乘法逆元的问题:
求解方程组c≡x(mod m),c≡b(mod n)的c,即
求解同余方程a≡c(mod m),b
≡c(mod n)
将其化为a=c+mx,b=c+ny,作差得
mx-ny=a-b,从而化为线性不定方程利用扩展欧几里得定理解决。
模m乘法逆元:
定义:对于整数a,m,如果存在整数b,满足ab ≡ 1(mod m),则说,b是a的模m乘法逆元。
定理:a存在模m的乘法逆元的充要条件是gcd(a,m) = 1
证明:充分性:
因为
gcd(a,m) = 1
根据欧拉定理,有
a^φ(m) ≡ 1(mod m)
因此
a * a^(φ(m)-1) mod m = 1
所以存在a的模m乘法逆元,即a^(φ(m)-1)
必要性:
假设存在a模m的乘法逆元为b,则
ab ≡ 1 (mod m)
所以
ab = km +1
所以
1 = ab - km
由欧几里得定理,有
gcd(a,m) = 1
由定理知:
对于ax + my = 1,可以看出x是a模m的乘法逆元,y是m模a的乘法逆元。
反过来,要计算a模m的乘法逆元,就相当于求ax + my = 1的x的最小正整数解,从而化为线性不定方程解决。
下面是如何求解最小正整数解(转):一般,我们能够找到无数组解满足条件,但是一般是让你求解出最小的那组解,怎么做?我们求解出来了一个特殊的解 x’ 那么,我们用 x' % m其实就得到了最小的解了(更广泛情况应为x0%(m/gcd))。为什么?
可以这样思考:
x 的通解不是 x' + m*t (t为任意整数)吗?(更广泛情况应为x0+(m/gcd)*t)(个人理解:x的通解可以由任何特解表示出来,因此也可以用最小正整数特解表示出来,假设最小正整数解为x0,则x=x0+m*t,那么x’=x0+m*t)
那么,也就是说, a 关于 m 的所有逆元是关于 m 同余的,那么根据最小整数原理,一定存在一个最小的正整数,它是 a 关于m 的逆元,而最小的逆元肯定是在(0 , m)之间的,而且只有一个,这就好解释了(更广泛情况应为a 关于 m 的所有逆元是关于 m/gcd 同余的,因为x%(m/gcd)==x0,x为a的任意逆元,x0为最小正整数解,x0在(0,m/gcd)之间,括号外面最后一句话解释了x0的存在性且唯一性)
。
可能有人注意到了,这里,我写通解的时候并不是 x0 + (m/gcd)*t ,但是想想一下就明白了,gcd = 1,所以写了跟没写是一样的,但是,由于问题的特殊性,有时候我们得到的特解 x0 是一个负数,还有的时候我们的 m 也是一个负数这怎么办?
当 m 是负数的时候,我们取 m 的绝对值就行了,当 x0 是负数的时候,他模上 m 的结果仍然是负数(在计算机计算的结果上是这样的,虽然定义的时候不是这样的),这时候,我们仍然让 x0 对abs(m) 取模,然后结果再加上abs(m) 就行了,于是,我们不难写出下面的代码求解一个数 a 对于另一个数 m 的乘法逆元:
- LL cal(LL a,LL m)
- {
- LL x,y;
- LL gcd=exgcd(a,m,x,y);
- if(1%gcd!=0) return -1;
- x*=1/gcd;
- if(m<0) m=-m;
- LL ans=x%m;
- if(ans<=0) ans+=m;
- return ans;
- }
这个地方可能会有人和我刚看到这段代码时一样的疑问,既然明知道要判断是否gcd==1了,那干嘛要写1%gcd!=0啊,直接写gcd!=1不就好了,其实这是因为这个子函数是一个模板,一般我们要求的等式为
a*x+b*y=c(此处b就相当于上面代码中m),c==1时我们可以像上面这么写(当然可以更简洁),如果c为任意数,代码就变成了如下的样子:
- LL equation(LL a, LL b, LL c, LL &x, LL &y)
- {
- LL g = extend_Euclid(a, b, x, y);
-
- if(c % g) return -1;
- LL ran = b / g;
- x *= c/g;
- if(ran < 0) ran = -ran;
- x = (x%ran + ran) % ran;
- return 0; }
这才是模板最基本的样子。
下面附上本题代码:
- #include <iostream>
- #include <cstdio>
- #include <cstring>
- #include <math.h>
- #include <algorithm>
- #define ll long long
- using namespace std;
- int gcd(int a,int b)
- {
- return b?gcd(b,a%b):a;
- }
- ll n,m;
- ll M;
- ll extend_Euclid(ll a,ll b,ll &x,ll &y)
- {
- if(b==0)
- {
- x=1;y=0;
- return a;
- }
- ll ans=extend_Euclid(b,a%b,y,x);
- y-=(a/b)*x;
- return ans;
- }
- ll equation(ll a,ll b,ll &x,ll &y,ll c)
- {
- ll gcd=extend_Euclid(a,b,x,y);
- if(c%gcd) return -1;
- ll k=abs(b/gcd);
- x*=c/gcd;
- x=(x%k+k)%k;
- return 0;
- }
- ll judge(ll x,ll y)
- {
- ll p,q;
- if(equation(2*n,-2*m,p,q,-x+y)==-1)
- return M+1;
- ll t=2*p*n+x;
- if(t<0||t>M)
- return M+1;
- return t;
- }
- ll solve(ll x,ll y)
- {
- ll g=gcd(n,m);
- M=1ll*m/g*n;
- ll ans=M+1;
- ans=min(ans,judge(x,y));
- ans=min(ans,judge(x,-y));
- ans=min(ans,judge(-x,y));
- ans=min(ans,judge(-x,-y));
- if(ans==M+1)
- return -1;
- return ans;
- }
- int main()
- {
- int k;
- scanf("%d%d%d",&n,&m,&k);
- ll x,y;
- for(int i=0;i<k;i++)
- {
- scanf("%lld%lld",&x,&y);
- printf("%lld\n",solve(x,y));
- }
- return 0;
- }
以上内容部分转载于:
http://blog.csdn.net/synapse7/article/details/9901195;
百度百科;
http://blog.csdn.net/zhjchengfeng5/article/details/7786595