crt,excrt学习总结

(crt,Chinese Remainder Theorem)

概述

• 前置技能:同余基础性质,(exgcd).
• (crt),中国剩余定理.用于解决模数互质的线性同余方程组.大概长这样:

[egin{equation} left{ egin{array}{lr} xequiv a_1(mod m_1),\ xequiv a_2(mod m_2),\ xequiv a_3(mod m_3),\ ......\ xequiv a_n(mod m_n).\ end{array} ight. end{equation} ]

• 其中,所有的 (m) 两两互质,否则需要使用 (excrt),此处不讨论.
• 记(N=prod m_i),对于每个方程(xequiv a_i(mod m_i)),尝试找一个 (y_i) 满足((N/m_i)*y_iequiv 1(mod m_i).)
• 这里的(N/m_i)就是其他所有模数的积,因为所有 (m) 互质,所以(N/m_i)与(m_i)互质.那么上面的线性同余方程一定可以通过 (exgcd) 找到这样的 (y_i).
• 再记(x_i=(N/m_i)*y_i),那么(x_i*a_iequiv a_i(mod m_i),x_i*a_iequiv 0(mod m_j,forall j ot= i)).
• 将各个方程合并起来即可得出原方程组的解(x=sum x_ia_i),容易验证,这个构造出的 (x) 满足原方程组,且对于(forall x'equiv x(mod N),x')也是合法的解.
• 时间复杂度为(O(ncdot(logprod m_i))).

应用

• 解线性同余方程组,拆分处理,合并答案.

题目

曹冲养猪

• 题意:求解一个线性同余方程组的最小正整数解,其中模数互质.
• 直接使用 (crt) 求解.因为要求最小正整数解,所以在合并的时候对 (N) 取模,且最后确保为正.

#include"bits/stdc++.h"
#define int long long
using namespace std;
typedef long long LoveLive;
inline int read()
{
	int out=0,fh=1;
	char jp=getchar();
	while ((jp>'9'||jp<'0')&&jp!='-')
		jp=getchar();
	if (jp=='-')
		{
			fh=-1;
			jp=getchar();
		}
	while (jp>='0'&&jp<='9')
		{
			out=out*10+jp-'0';
			jp=getchar();
		}
	return out*fh;
}
const int MAXN=11;
int a[MAXN],m[MAXN];
int n;
void exgcd(int a,int b,int &x,int &y)
{
	if(!b)
		{
			x=1;
			y=0;
			return;
		}
	exgcd(b,a%b,x,y);
	int t=x;
	x=y;
	y=t-(a/b)*y;
}
int crt(int a[],int m[],int n)
{
	int res=0,N=1,x,y;
	for(int i=1;i<=n;++i)
		N*=m[i];
	for(int i=1;i<=n;++i)
		{
			int tmp=N/m[i];
			exgcd(m[i],tmp,x,y);
			res=(res+y*tmp*a[i])%N;
		}
	return (res+N)%N;
}
signed main()
{
	n=read();
	for(int i=1;i<=n;++i)
		m[i]=read(),a[i]=read();
	int ans=crt(a,m,n);
	printf("%lld
",ans);
	return 0;
}

(Excrt,Extended crt)

概述

• 前置技能:(crt)(其实不是很必要),(exgcd).
• 在 (crt) 中,要求所有模数互质.如果模数不互质, (crt) 就解决不了问题.因为构造过程中会出现(N/m_i)与(m_i)不互质的情况,同余方程((N/m_i)*y_iequiv 1(mod m_i))会出现无解的情况.
• 这时候就需要使用 (excrt) 进行处理.
• (excrt) 的基本思路与 (crt) 不同,(crt) 是直接通过构造得出解, (excrt) 是通过将方程两两合并,直到最后只剩下一个方程,也就得到了通解.
• 对于如下方程组:

[egin{equation}left{ egin{array}{lr} xequiv a_1(mod m_1),\ xequiv a_2(mod m_2),\ xequiv a_3(mod m_3),\ ......\ xequiv a_n(mod m_n).\ end{array} ight.end{equation} ]

• 先考虑合并前两个方程,将其写成一个等价的线性同余方程.
• 由同余的基本性质可知,一个合法解 (x_0) 满足 (x_0=a_1+k_1*m_1=a_2+k_2*m_2.)
• 将后两个式子移项整理可得(k_2*m_2-k_1*m_1=a_1-a_2.)
• 这里就是线性不定方程的形式了,尝试使用 (exgcd) 求解.
• 可记(gcd(m_1,m_2)=g,a_1-a_2=c.)那么当(g ot| c)时,显然无解,这两个方程无法同时满足,那么整个方程组也就无解.
• 否则,先通过 (exgcd) 求出满足方程 (k_2*m_2-k_1*m_1=g)的一组解(k_2',-k_1').
• 方程两边乘上(c/g),即(-k_1')乘上(c/g)就得到了(-k_1).这里可以让(-k_1)对(m_2)取模,防止其过大.
• 然后代入第一个方程求出(x_0=a_1+k_1*m_1),注意上一步求出的是(-k_1),这里要取负号.
• 这两个方程的通解即为(x=x_0+lcm(m_1,m_2)),写成同余式为(xequiv x_0(mod lcm(m_1,m_2)).)
• 新得到的同余式和原来的两个同余式并在一起是等价的,那么将新得到的同余式继续按照上述步骤与其他同余式合并,直到最后只剩下一个同余式,就是要求的通解了.
• 时间复杂度为(O(ncdot log(prod m_i)).)

应用

• 解线性同余方程组,合并答案.

题目

Strange Way to Express Integers

• 题意:求解一个线性同余方程组的最小正整数解,不保证模数互质.
• 使用 (excrt) 处理即可.注意调整解.

#include"bits/stdc++.h"
#define int long long
using namespace std;
typedef long long LoveLive;
inline int read()
{
	int out=0,fh=1;
	char jp=getchar();
	while ((jp>'9'||jp<'0')&&jp!='-')
		jp=getchar();
	if (jp=='-')
		{
			fh=-1;
			jp=getchar();
		}
	while (jp>='0'&&jp<='9')
		{
			out=out*10+jp-'0';
			jp=getchar();
		}
	return out*fh;
}
const int MAXN=1e5+10;
int a[MAXN],m[MAXN],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,x,y);
	int t=x;
	x=y;
	y=t-(a/b)*y;
	return d;
}
int excrt(int a[],int m[],int n)
{
	int k1,k2;
	for(int i=2;i<=n;++i)
		{
			int g=exgcd(m[i],m[1],k2,k1);
			int c=a[1]-a[i];
			if(c%g)
				return 0;
			k1*=c/g;
			k1%=m[i];
			a[1]=a[1]+(-k1)*m[1];//注意修改的先后顺序 
			m[1]=m[1]*m[i]/g;	
		}
	return 1;
} 
signed main()
{
	while(~scanf("%d",&n))
		{
			for(int i=1;i<=n;++i)
				m[i]=read(),a[i]=read();
			int flag=excrt(a,m,n);
			if(!flag)
				puts("-1");
			else
				printf("%lld
",(a[1]%m[1]+m[1])%m[1]);
		}
	return 0;
}

屠龙勇士

• 题面比较长...直接看链接内容即可(注意数据表格).
• 每轮剑的攻击力可以用一个 (multiset) 维护,可以直接预处理出来,记每轮攻击的攻击力为(atk_i).
• 考虑记要求的攻击次数为 (x) ,则每次攻击为使恢复若干次后生命值恰好为 (0) ,只需满足(atk_i*xgeq a_i,atk_i*x equiv a_i(mod p_i).)
• 注意到数据中要么满足特性 (1) , 即(a_ileq p_i),要么满足所有的 (p_i=1).
• 那么可以判断后分情况处理,对于满足特性 (1) 的部分,显然只需处理同余式,同余式成立的情况下一定有(atk_i*xgeq a_i).对于所有 (p_i=1) 的部分,只需求出将每只龙的生命值打到非正数的次数,取最大值即可,即(max lceil a_i/atk_i ceil),容易处理.
• 对于 (a_ileq p_i) 的部分,我们得到若干个同余式,考虑使用 (excrt) 进行处理.
• 同余式形式是 (atk_i*x equiv a_i(mod p_i)) ,需要先转化为 (excrt) 的标准形式,即

[egin{equation}left{ egin{array}{lr} xequiv a_1(mod m_1),\ xequiv a_2(mod m_2),\ xequiv a_3(mod m_3),\ ......\ xequiv a_n(mod m_n).\ end{array} ight.end{equation} ]

• 考虑将同余式两边乘上逆元 (atk_i^{-1}),可得到(x equiv atk_i^{-1}*a_i(mod p_i)),就化为了标准形式.
• 但 (atk_i) 与 (p_i) 可能不互质,此时 (atk_i) 是没有逆元的.可以根据同余的消去律,将(atk_i,a_i,p_i)都约去三者的 (gcd) .若此时 (atk_i,p_i) 仍不互质,则整个方程组也就无解了.否则求出逆元后,转化为标准形式.最后对所有转化后的同余式使用 (excrt) 求解.
• 注意中间运算会出现$long long * long long % long long $的情况,使用快速乘,即将快速幂中的乘法运算改为加法运算,中途取模.

#include"bits/stdc++.h"
#define int long long
using namespace std;
typedef long long LoveLive;
inline int read()
{
	int out=0,fh=1;
	char jp=getchar();
	while ((jp>'9'||jp<'0')&&jp!='-')
		jp=getchar();
	if (jp=='-')
		{
			fh=-1;
			jp=getchar();
		}
	while (jp>='0'&&jp<='9')
		{
			out=out*10+jp-'0';
			jp=getchar();
		}
	return out*fh;
}
multiset<int> sword;
multiset<int>::iterator it;
int n,m;
const int MAXN=1e5+10;
int fmul(int a,int b,int mod)
{
	assert(mod>0);
	a=(a%mod+mod)%mod;
	b=(b%mod+mod)%mod;
	int res=0;
	while(b)
		{
			if(b&1)
				res=(res+a)%mod;
			a=(a+a)%mod;
			b>>=1;
		}
	return res;
}
int exgcd(int a,int b,int &x,int &y)
{
	if(!b)
		{
			x=1;
			y=0;
			return a;
		}
	int d=exgcd(b,a%b,x,y);
	int t=x;
	x=y;
	y=t-(a/b)*y;
	return d;
}
int a[MAXN],p[MAXN],atk[MAXN];
int award[MAXN];
bool check()
{
	for(int i=1;i<=n;++i)
		if(p[i]!=1)
			return false;
	return true;
}
void pre()
{
	for(int i=1;i<=n;++i)
		{
			it=sword.upper_bound(a[i]);
			if(it!=sword.begin())
				--it;
			atk[i]=*it;
			sword.erase(it);
			sword.insert(award[i]);
		}
}
void solve_sp()
{
	int ans=-1;
	for(int i=1;i<=n;++i)
		{
			ans=max(ans,(a[i]+atk[i]-1)/atk[i]);
		}
	printf("%lld
",ans);
}
int A[MAXN],M[MAXN];
int inv(int k,int mod)
{
	int x,y;
	exgcd(k,mod,x,y);
	return (x%mod+mod)%mod;
}
bool build_equations()
{
	int x,y;
	for(int i=1;i<=n;++i)
		{
			int g=exgcd(a[i],p[i],x,y);
			g=exgcd(g,atk[i],x,y);
			atk[i]/=g;
			a[i]/=g;
			p[i]/=g;
			int flag=exgcd(atk[i],p[i],x,y);
			if(flag!=1)
				return false;
			A[i]=fmul(a[i],inv(atk[i],p[i]),p[i]);
			M[i]=p[i];
		}
	return true;
}
void excrt()
{
	int k1,k2;
	for(int i=2;i<=n;++i)
		{
			int g=exgcd(M[i],M[1],k2,k1);
			int c=A[1]-A[i];
			if(c%g)
				{
					puts("-1");
					return;
				}
			k1=fmul(k1,c/g,M[i]);
			/*k1*=c/g;
			k1%=M[i];*/
			A[1]=A[1]-fmul(k1,M[1],M[1]*(M[i]/g));//打括号,防止爆long long 
			M[1]=M[1]*(M[i]/g);
			A[1]%=M[1];
		}
	printf("%lld
",(A[1]%M[1]+M[1])%M[1]);
}
signed main()
{
//	freopen("dragon8.in","r",stdin);
	int T=read();
	while(T--)
		{
			sword.clear();
			n=read();
			m=read();
			for(int i=1;i<=n;++i)
				a[i]=read();
			for(int i=1;i<=n;++i)
				p[i]=read();
			for(int i=1;i<=n;++i)
				award[i]=read();
			for(int i=1;i<=m;++i)
				sword.insert(read());
			pre();
			if(check())
				{
					solve_sp();
					continue;
				}
			if(build_equations())
				excrt();
			else
				puts("-1");
		}
	return 0;
}