传送门
中国剩余定理(CRT)
我的第一反应是小学奥数题——韩信点兵。
转化成数学语言,就是给你 (n) 个关于 (x) 的同余方程(保证 (a_i) 互质),求最小整数解。
[egin{cases}
xequiv b_1pmod {a_1}\
xequiv b_2pmod {a_2}\
xequiv b_3pmod {a_3}\
cdotscdotscdotscdots\
end{cases}]
怎么做呢?
可以仿照拉格朗日插值法,构造一组通解:
[ans=sum_{i=1}^{n}b_iprod_{j=1,j
e i}^{n}(a[j]*a[j]^{-1})
]
其中 (a[j]^{-1}) 表示的是 (a[j]) 在模 (a[i]) 下的逆元。
因为 (a[i]) 之间两两互质,所以我们可以放心的用欧拉定理:
[a^{psi(m)}equiv 1pmod m
]
所以
[a^{psi(m)-1}equiv frac{1}{a}pmod m
]
可见 (a[j]^{psi(a[i])-1}) 即是逆元。
总结一下思路
预处理出 (prod a_i) 和 (psi(1) ext{ 到 }psi(max\_a)),带入公式中,利用快速幂求出逆元,把答案累加起来。
注意事项
- 数据有锅,第一个点的 (a) 高达 (1093258),注意数组开的范围
- 模数最大为 (10^{18}) ,所以需要使用快
(gui)速乘
AC代码
#include<iostream>
using namespace std;
bool vis[2000005];
int n,cnt,prime[2000005];
long long a[15],b[15],prod_a=1,phi[2000005],maxa,ans;
inline long long ksc(long long x,long long y,long long mod)
{
return (x*y-(long long)((long double)x/mod*y)*mod+mod)%mod;
}
long long q_p(long long a,long long b){
if(b==0) return 1;
if(b==1) return a%prod_a;
long long res=q_p(a,b/2);
if(b&1) return ksc(ksc(res,res,prod_a),a,prod_a);
return ksc(res,res,prod_a);
}
int main()
{
cin>>n;
for(int i=1;i<=n;i++) cin>>a[i]>>b[i];
for(int i=1;i<=n;i++) prod_a*=a[i],maxa=max(maxa,a[i]);
phi[1]=1;
for(int i=2;i<=maxa;i++){
if(!vis[i]) prime[++cnt]=i,phi[i]=i-1;
for(int j=1;j<=cnt&&i*prime[j]<=maxa;j++){
vis[i*prime[j]]=1;
if(i%prime[j]){
phi[i*prime[j]]=phi[i]*phi[prime[j]];
}else{
phi[i*prime[j]]=phi[i]*prime[j];
break;
}
}
}
for(int i=1;i<=n;i++){
ans=(ans+ksc(ksc(prod_a/a[i],q_p(prod_a/a[i],phi[a[i]]-1),prod_a),b[i],prod_a))%prod_a;
}
cout<<ans;
return 0;
}