给定正整数n1,n2,...,nk(未必两两互质),要求找到x,满
足x≡ai(mod ni) (i=1,2...k)
x ≡ a1( mod n1)
足x≡ai(mod ni) (i=1,2...k)
x ≡ a1( mod n1)
x ≡ a2( mod n2)
...........
即n1,....nk之间不一定互质
那便是解k个线性方程,如下
x + u*n1 = a1
x - v*n2 = a2
.............
先取前两个方程
x + u*n1 = a1
x - v*n2 = a2
关于此方程,在gcd(n1,n2)|(a1-a2)时有解联立方程,得n1*u + n2*v = (a1-a2) ,利用扩展欧几里得算法,可得一解(u0,v0)
那么此方程组解集为:
u = u0 + t*(n2/gcd(n1,n2)) t为任意整数
带入原方程,得x=a1 - u0*n1 + t*n1*n2/gcd(n1,n2) t为任意整数
因为n1*n2/gcd(n1,n2)=lcm(n1,n2)即最小公倍数
所以x=a1 - u0*n1 + t*lcm(n1,n2)
所以可得x ≡ a1 - u0*n1 (mod lcm(n1,n2))
即
x ≡ a1( mod n1)
x ≡ a2( mod n2)
==>合并为
x ≡ a1 - u0*n1 (mod lcm(n1,n2))
所以a1==> a1 - u0*n1
n1==>lcm(n1,n2)
合并到最后只剩下一个时,便可求出x的集合
代码如下:
#include<iostream>
using namespace std;
int gcdEx(int a,int b,int& x,int& y)
{ //求ax+by=gcd(a,b) 的整数解 返回gcd(a,b);
if(b==0) {
x=1; y=0; return a;
}
int x1,y1;
int gcd=gcdEx(b,a%b,x1,y1);
x=y1;
y=x1-a/b*y1;
return gcd;
}
int CRTEx(int *a,int *n,int k)
{
int a1=a[0],n1=n[0],lcm,u0,x,y,a2,n2,d;
for(int i=1;i<k;i++)
{
a2=a[i]; n2=n[i];
d=gcdEx(n1,n2,x,y);
if((a2-a1)%d) return -1; //无解
lcm=n1/d*n2;
u0=(a1-a2)/d*x%n2; //n1*x + n2*y = (a-a1) 的一个解u0 ,别忘记%n1
a1=(a1-u0*n1)%lcm; //合并a
n1=lcm; //合并n
}
return (a1%n1+n1)%n1; //只剩下方程x≡a1(mod n1)
}