求一个形如(xequiv a_i (mod p_i)) 的线性方程组的解
我们考虑合并两个方程:
[egin{cases}
x=a_1 (mod p_1)\
\
x=a_2 (mod p_2)\
end{cases}
]
先考虑将(x)设为(a_1),再不断的加(p_1)使得(x)在模(p_2)意义下为(a_2)
我们可以用(exgcd)求得一组(mp_1+np_2= gcd(p_1,p_2))的解,即(mp_1equiv gcd(p1,p2) (mod p2))
也就是说我们每次增加(mp_1)即可让(x)在模(p_1)意义下增加(gcd(p_1,p_2))
令(d=gcd(p_1,p_2)),那么我们应该增加((a_2-a_1)/d) 次,即(x)应该加上((a_2-a_1)/d*mp_1),然后我们知道方程合并后的模数应该是(P=p_1*p_2/gcd(p_1,p_2)),所以再将结果对(P)取模即可。
代码:
LL exgcd(LL x,LL y,LL &a,LL &b){
LL z;
return y?(z=exgcd(y,x%y,b,a),b-=x/y*a,z):(a=1,b=0,x);
}
LL mult(LL a,LL b,LL p){
LL ans=0;
while(b){
if(b&1)ans=(ans+a)%p;
b>>=1; a=(a<<1)%p;
}
return ans;
}
struct CRT{
LL a[N],p[N],P,x;
bool excrt(LL a1,LL p1,LL a2,LL p2){
LL m,n,d=exgcd(p1,p2,m,n);
LL a=(a2-a1%p2+p2)%p2;
if(a%d!=0)return false;
P=p1/d*p2; m=(m%p2+p2)%p2;
x=(x+mult(a/d*p1,m,P))%P; //注意这里一般先计算a/d*p1,然后再乘m,因为题目一般保证了所有p_i的最小公倍数小于某个值,这样前面的运算就不会溢出,而后面的乘m可能会溢出,所以还要用快速乘。
return true;
}
bool run(){
x=a[1],P=p[1];
for(int i=2;i<=n;i++)
if(!excrt(x,P,a[i],p[i]))
return false;
return true;
}
}crt;
题目要把一个形如(Axequiv a (mod p))的方程变成(xequiv b(mod p))的形式。
肯定有无解的情况,我们求(d=gcd(A,p)),然后将方程每一项除以(d),如果(d)不是(a)的约数,那么方程无解。
然后我们得到了等价方程(A'xequiv a' (mod p')),其中(A'=A/d, a'=a/d, p'=p/d)
我们用(exgcd)求得(mA+np=gcd(A,p))的一组解,则这组解满足(mA'+np'=1),即(mA'equiv1 (mod p')) ,(m)即为(A')在模(p')意义下的逆元。
这样上面的方程就可以转化为(xequiv a'm (mod p'))
这样就是一个(xequiv b(mod p))的形式了。
代码:
bool flag=0;
for(int i=1;i<=n;i++){
LL A1=A[i],a1=a[i],p1=p[i];
LL m,n,d=exgcd(A1,p1,m,n);
if(a1%d!=0){flag=1; break;}
A1/=d,a1/=d,p1/=d; m=(m%p1+p1)%p1;
crt.a[i]=mult(a1,m,p1),crt.p[i]=p1;
}