逆元:
对于$ax equiv 1 (mod ~m)$这样的同余方程。
将x称为a的逆元,满足乘以a即可为乘以a的逆元x。
逆元的求解:
对单个逆元有三种求法:
1:费马小定理:快速幂求解。
2:欧拉定理:x = a^(phi(m)-1)
可以发现,就是费马小定理的扩展。
3:扩展欧几里得。
对于$ax equiv 1 (mod ~m)$。
我们将商设为y,即可转化为:$ax = my+1$
$ax = my+1 ightarrow ax+(-y)*m = 1 $
那么就是扩展欧几里得,求解满足条件的x即可,然后把x取模到区间(0~m-1)。
LL exgcd(LL a,LL b,LL &x,LL &y) { if(b == 0) { x = 1,y = 0; return a; } int t = exgcd(b,a%b,y,x); y -= a/b*x; return t; } int main() { int n,p;n = read(),p = read(); for(rg int i = 1;i <= n;++i) { LL x,y; LL gcd = exgcd(i,p,x,y); x = (x+p)%p; printf("%lld ",x); } system("pause"); }
对于多个逆元的递推:
这里给出线性递推的式子:inv[i] = (mod-mod/i)*inv[mod%i]%mod;
中国剩余定理:
不证明还是很简单的~
中国剩余定理用于求解,对于多组的a[i]和b[i],在满足$x equiv ai~ mod ~(bi)$时,求解最小的非负整数x。
条件:必须满足所有b[i]互质
解法:定义$M = prod_{i = 1}^{n} b[i]$
$Mi = frac{M}{b[i]]}$
可以构造出解$x = sum_{i = 1}^{n} ai * Mi * inv[Mi]$
x要对M取模
对于Mi的逆元,这里使用扩展欧几里得求解。
// Author: levil #include<bits/stdc++.h> using namespace std; typedef long long LL; typedef pair<int,int> pii; const int N = 2e3+5; const int M = 2e6+5; const LL Mod = 1e9+7; #define rg register #define pi acos(-1) #define INF 1e9 #define CT0 cin.tie(0),cout.tie(0) #define IO ios::sync_with_stdio(false) #define dbg(ax) cout << "now this num is " << ax << endl; namespace FASTIO{ inline LL read(){ LL x = 0,f = 1;char c = getchar(); while(c < '0' || c > '9'){if(c == '-') f = -1;c = getchar();} while(c >= '0' && c <= '9'){x = (x<<1)+(x<<3)+(c^48);c = getchar();} return x*f; } void print(int x){ if(x < 0){x = -x;putchar('-');} if(x > 9) print(x/10); putchar(x%10+'0'); } } using namespace FASTIO; void FRE(){/*freopen("data1.in","r",stdin); freopen("data1.out","w",stdout);*/} LL a[15],b[15]; LL exgcd(LL a,LL b,LL &x,LL &y) { if(b == 0) { x = 1,y = 0; return a; } LL t = exgcd(b,a%b,y,x); y -= a/b*x; return t; } int main() { int n;n = read(); LL M = 1; for(rg int i = 1;i <= n;++i) a[i] = read(),b[i] = read(),M *= a[i]; LL ans = 0; for(rg int i = 1;i <= n;++i) { LL Mi = M/a[i],x,y; exgcd(Mi,a[i],x,y); x = (x+a[i])%a[i]; ans += Mi*x*b[i]; } printf("%lld ",ans%M); system("pause"); }
扩展中国剩余定理:
我们知道,中国剩余定理在求解同余方程组$x equiv ai~ mod ~(bi)$时,必须要满足模数b[i]全部互质。
当b[i]不满足全部互质时,就有了扩展中国剩余定理。
我们假设当前已经求出了前k-1个同余方程组的解为X,$M = prod_{i = 1}^{k-1} b[i]$
那么他们的通解为X+t*M。
那么对于第k个式子的解,xk = X+t*M.
那么因为$xk equiv a[i]~ mod ~(b[i])$
进一步转化$xk equiv a[i]~ mod ~(b[i]) ightarrow X+t*M equiv a[i]~mod~(b[i]) ightarrow X + t*m = b[i]*y+a[i]$
那么我们就是要求解一个满足条件的t使得$X + t*m = b[i]*y+a[i]$成立。
这里显然是用扩展欧几里得
对上式移项得$t * M + b[i]*(-y) = a[i]-X$
令d = a[i]-X
我们知道扩展欧几里得可以求解出$t * M + b[i]*(-y) = gcd(M,b[i])$
当gcd(M,b[i]) | d时,这里也就是判断同余方程是否有解的关键。
当不满足这个条件时,说明t * M + b[i]*(-y) = d无解。
当满足条件时,对于满足$t * M + b[i]*(-y) = d$的x就是上面的x的d/gcd倍。
求出这个x然后xk就有了,然后X = xk。
注意的是这里的X,最后还需要对M取模。
且M = M*b[i]/gcd。本质上我们应该是对M*b[i]。但是这里求出的gcd对前k个同余方程都满足,所以我们再对M/gcd,防止M爆longlong。
细节:X*d/gcd会爆longlong,需要快速乘。
M = M*(b[i]/gcd) 而不是 M = M*b[i]/gcd,后者会爆longlong
Code:洛谷P4777 【模板】扩展中国剩余定理(EXCRT)
// Author: levil #include<bits/stdc++.h> using namespace std; typedef long long LL; typedef pair<int,int> pii; const int N = 1e5+5; const int M = 2e6+5; const LL Mod = 1e9+7; #define rg register #define pi acos(-1) #define INF 1e9 #define CT0 cin.tie(0),cout.tie(0) #define IO ios::sync_with_stdio(false) #define dbg(ax) cout << "now this num is " << ax << endl; namespace FASTIO{ inline LL read(){ LL x = 0,f = 1;char c = getchar(); while(c < '0' || c > '9'){if(c == '-') f = -1;c = getchar();} while(c >= '0' && c <= '9'){x = (x<<1)+(x<<3)+(c^48);c = getchar();} return x*f; } void print(int x){ if(x < 0){x = -x;putchar('-');} if(x > 9) print(x/10); putchar(x%10+'0'); } } using namespace FASTIO; void FRE(){/*freopen("data1.in","r",stdin); freopen("data1.out","w",stdout);*/} LL a[N],b[N]; LL exgcd(LL a,LL b,LL &x,LL &y) { if(b == 0) { x = 1,y = 0; return a; } LL t = exgcd(b,a%b,y,x); y -= a/b*x; return t; } LL Mul(LL a,LL b,LL p) { LL re = 0; while(b) { if(b&1) re = (re+a)%p; a = (a+a)%p; b >>= 1; } return re; } int main() { int n;n = read(); for(rg int i = 1;i <= n;++i) b[i] = read(),a[i] = read(); LL X = a[1],M = b[1]; for(rg int i = 2;i <= n;++i) { LL d = ((a[i]-X)%b[i]+b[i])%b[i],x,y; LL gcd = exgcd(M,b[i],x,y); //if(d%gcd != 0){break;}//无解 LL t = Mul(x,d/gcd,b[i]); X = X+t*M; M = M*(b[i]/gcd); X = (X+M)%M; } printf("%lld ",X); system("pause"); }