P4774 [NOI2018] 屠龙勇士 - 数论、中国剩余定理

题意

给出 \(n\) 个线性同余方程构成的方程组 \(\begin{cases} a_1x\equiv b_1 \pmod{p_1} \\ a_2x\equiv b_2 \pmod{p_2} \\ \dots \\ a_nx\equiv b_n \pmod{p_n} \end{cases}\)

问其大于等于某个数的最小解，如果无解输出 \(-1\)。

题解

首先考虑某个线性同余方程 \(a_ix\equiv b_i \pmod{p_i}\) 怎么解。它其实等价于不定方程 \(a_ix+p_iy=b_i\)，设其特解为 \(\begin{cases} x=x_0 \\ y=y_0 \end{cases}\)，那么其通解就是 \(\begin{cases} x=x_0+k\frac{p_i}{\gcd(p_i,a_i)} \\ y=y_0+k\frac{a_i}{\gcd(p_i,a_i)} \end{cases},k\in \mathbb{Z}\)。转化成同余式就是 \(x\equiv x_0 \pmod{\frac{p_i}{\gcd(p_i,a_i)}}\)。

然后套 exCRT 即可。注意 exCRT 过程中一定要多取模！！！

代码

#include <cstdio>
#include <cstring>
#include <cctype>
#include <set>
#include <cassert>
#include <algorithm>
using namespace std;
#define For(Ti,Ta,Tb) for(int Ti=(Ta);Ti<=(Tb);++Ti)
#define Dec(Ti,Ta,Tb) for(int Ti=(Ta);Ti>=(Tb);--Ti)
template<typename T> void Read(T &x){
	x=0;int _f=1;
	char ch=getchar();
	while(!isdigit(ch)) _f=(ch=='-'?-1:_f),ch=getchar();
	while(isdigit(ch)) x=x*10+(ch^48),ch=getchar();
	x=x*_f;
}
template<typename T,typename... Args> void Read(T &x,Args& ...others){
	Read(x);Read(others...);
}
typedef long long ll;
ll Mul(ll a,ll b,ll p){
	a=(a%p+p)%p,b=(b%p+p)%p;
	ll res=0;
	while(b){
		if(b&1) res=(res+a)%p;
		b>>=1,a=a*2%p;
	}return res;
}
ll Lcm(ll a,ll b){
	if(!a||!b) return a^b;
	return a/__gcd(a,b)*b;
}
ll ExEuclid(ll a,ll b,ll &x,ll &y){
	if(!b){x=1,y=0;return a;}
	ll t,g=ExEuclid(b,a%b,t,y);
	x=y,y=t-a/b*y;return g;
}
bool SolveEq(ll a,ll b,ll c,ll &x){//ax=b(mod c)
	ll y,g=ExEuclid(a,c,x,y);
	if(b%g) return 0;
	x=Mul(x,b/g,c/g);return 1;
}
const int N=1e5+5;
int T,n,m;ll a[N],p[N],rew[N],atk[N],sword[N];
int main(){
	Read(T);
	while(T--){
		Read(n,m);
		For(i,1,n) Read(a[i]);
		For(i,1,n) Read(p[i]);
		For(i,1,n) Read(rew[i]);
		For(i,1,m) Read(atk[i]);
		multiset<ll> sw(atk+1,atk+m+1);
		For(i,1,n){
			auto it=sw.upper_bound(a[i]);
			if(it!=sw.begin()) --it;
			sword[i]=*it;sw.erase(it);
			sw.insert(rew[i]);
		}
		ll lc=1,x=0,mx=0,succ=1;
		For(i,1,n){
			ll sol1;
			if(!SolveEq(sword[i],a[i],p[i],sol1)){
				puts("-1");succ=0;break;
			}
			ll cur=sol1;p[i]/=__gcd(p[i],sword[i]);
			if(!SolveEq(lc,((cur-x)%p[i]+p[i])%p[i],p[i],sol1)){
				puts("-1");succ=0;break;
			}
//			assert(lc*sol1+x==p[i]*sol2+cur);
			ll mod=Lcm(lc,p[i]);
			x=(Mul(lc,sol1,mod)+x+mod)%mod;lc=mod;
			mx=max(mx,(a[i]+sword[i]-1)/sword[i]);
		}
		if(succ){
			ll k=(mx-x+lc-1)/lc;
			printf("%lld\n",x+k*lc);
		}
	}
	return 0;
}