中国剩余定理CRT

更新日志 2024/12/09:开工。

2024/12/10:完成证明、添加模板和例题


概念

就不写引入部分了,直接进入正题。

中国剩余定理用来求解形如下式的方程的一组特解:

\[\begin{cases} x\equiv r_1\pmod {m_1}\\ x\equiv r_2\pmod {m_2}\\ x\equiv r_3\pmod {m_3}\\ \dots\\ x\equiv r_n\pmod {m_n} \end{cases} \]

其中 \(m\) 两两互质。

如果不互质,见拓展中国剩余定理

思路

考虑构造法,构造出一组可行解即可。

事实上,先人已经为我们构造出了一种简单的 \(O(n)\) 方案,所以我们直接记即可。

为了简洁,这里记 \(M\gets\prod\limits_{j=1}^nm_i\)

\[x=\sum_{i=1}^n[r_i\cdot\frac{M}{m_i}\cdot(\frac{M}{m_i})^{-1}] \]

其中 \((\frac{M}{m_i})^{-1}\) 表示 \(\frac{M}{m_i}\) 在模 \(m_i\) 意义下的逆元。

证明

我们证明 \(\forall i\in[1,n],x\equiv r_i\pmod {m_i}\) 即可。

\[\because\forall j\ne i,\frac{M}{m_j}\equiv 0\pmod {m_i} \]

\[\begin{align*} \forall i\in[1,n],x&\equiv \sum_{j=1}^n[r_j\cdot\frac{M}{m_j}\cdot(\frac{M}{m_j})^{-1}]&\pmod {m_i}\\ x&\equiv r_i\cdot\frac{M}{m_i}\cdot(\frac{M}{m_i})^{-1}&\pmod {m_i}\\ x&\equiv r_i&\pmod{m_i} \end{align*} \]

证毕。

常见要求

求通解

不难发现 \(\forall i\in[1,n],M\equiv0\pmod{m_i}\)

故而通解为 \(x_0+kM,k\in Z\)

求最小正整数解

我们根据通解知道 \(x_0\) 加减任意个 \(M\) 仍然为解。

直接 \(x_0\bmod M\) 即可。若 \(\le0\)\(+M\) 就行了。

模板

int n;
int r[N],m[N];
int exgcd(int a,int b,int &x,int &y){
	if(!b){
		x=1,y=0;
		return a;
	}
	int d=exgcd(b,a%b,y,x);
	y-=a/b*x;
	return d;
}
int CRT(){
	int M=1;rep(i,1,n)M*=m[i];
	int res=0,inv,tmp;
	rep(i,1,n){
		exgcd(M/m[i],m[i],inv,tmp);
		res+=r[i]*(M/m[i])*inv;
	}
	return res;
}

请注意数据范围,如数据范围很大,开启long long且计算答案时实时%M,保险起见乘法过程中强制转换__int128

下面是一个修改示例:

int n;
ll r[N],m[N];
ll exgcd(ll a,ll b,ll &x,ll &y){
	if(!b){
		x=1,y=0;
		return a;
	}
	ll d=exgcd(b,a%b,y,x);
	y-=a/b*x;
	return d;
}
ll CRT(){
	ll M=1;rep(i,1,n)M*=m[i];
	ll res=0,inv,tmp;
	rep(i,1,n){
		exgcd(M/m[i],m[i],inv,tmp);
		res=((res+(__int128)r[i]*(M/m[i])*inv%M)%M+M)%M;
	}
	return res;
}

例题

LG1495

代码

前注:非题解,不做详细讲解

求中间有一段开了__int128

#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
typedef unsigned long long ull;
typedef __int128 i128;
typedef double db;
typedef long double ld;
typedef pair<int,int> pii;
typedef pair<ll,ll> pll;
typedef pair<int,ll> pil;
typedef pair<ll,int> pli;
template <typename Type>
using vec=vector<Type>;
template <typename Type>
using grheap=priority_queue<Type>;
template <typename Type>
using lrheap=priority_queue<Type,vector<Type>,greater<Type> >;
#define fir first
#define sec second
#define pub push_back
#define pob pop_back
#define puf push_front
#define pof pop_front
#define chmax(a,b) a=max(a,b)
#define chmin(a,b) a=min(a,b)
#define rep(i,x,y) for(int i=x;i<=y;i++)
#define per(i,x,y) for(int i=x;i>=y;i--)

const int inf=0x3f3f3f3f;
const ll INF=0x3f3f3f3f3f3f3f3f;
const int mod=1e9+7/*998244353*/;

const int N=15;

int n;
ll r[N],m[N];
ll exgcd(ll a,ll b,ll &x,ll &y){
	if(!b){
		x=1,y=0;
		return a;
	}
	ll d=exgcd(b,a%b,y,x);
	y-=a/b*x;
	return d;
}
ll CRT(){
	ll M=1;rep(i,1,n)M*=m[i];
	ll res=0,inv,tmp;
	rep(i,1,n){
		exgcd(M/m[i],m[i],inv,tmp);
		res=((res+(__int128)r[i]*(M/m[i])*inv%M)%M+M)%M;
	}
	return res;
}

int main(){
    ios::sync_with_stdio(false);
    cin.tie(0);cout.tie(0);
    cin>>n;
    rep(i,1,n)cin>>m[i]>>r[i];
    cout<<CRT();
    return 0;
}
posted @ 2024-12-09 14:00  HarlemBlog  阅读(5)  评论(0编辑  收藏  举报