【学习笔记】快速乘，CRT，exCRT，BSGS

CRT，exCRT，BSGS，玄学算法：MillerRabin判大质数，PollardRho大数因子寻找

这篇很大程度都是学习这位博主qaq：https://www.cnblogs.com/812-xiao-wen/

讲的真好

后面两个玄学算法的链接：

MillerRabin：https://www.cnblogs.com/812-xiao-wen/p/10543928.html

PollardRho：https://www.cnblogs.com/812-xiao-wen/p/10544546.html

第一遍看的懵懵懂懂，我有时间一定会去学的

之所以把这些放到一起，是因为都用到了快速乘这个东西，所以先介绍下快速乘，之后主要总结exCRT和BSGS

快速乘

快速乘主要用于用longlong乘法也会溢出的时候（模数也会是longlong），此时就需要不爆longlong算乘法取模的算法，所以就有了快速乘

原理参考：https://www.cnblogs.com/812-xiao-wen/p/10543023.html

这里就放个板子吧：快速乘和用快速乘的快速幂

typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
//O(1)快速乘 
inline ll ksc(ll x,ll y,ll p)
{
	ll z=(ld)x/p*y;
	ll res=(ull)x*y-(ull)z*p;
	return (res+p)%p;
}
//O(log)快速幂
inline ll ksm(ll x,ll y,ll p)
{
	ll res=1;
	while(y){
		if(y&1) res=ksc(res,x,p);
		x=ksc(x,x,p);y>>=1;
	}
	return res;
} 

中国剩余定理：

求解形如：

x=a1(mod m1)

x=a1(mod m2)

…

x=an(mod mn)

其中mi两两互质

这样的线性同余方程组的一个特解

具体算法：

m=∏mi，Mi=m/mi，ti为Mi·ti=1(mod mi)的一个解

对于方程组，有解：x=Σai·Mi·ti，通解为：x+km k为整数

板子：（求的是最小非负整数解）

typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
//扩展欧几里得
//输入a,b,任意x,y 
//返回ax+by=gcd(a,b)的一组特解x0,y0 
ll exgcd(ll a,ll b,ll &x,ll &y)
{
	if(b==0){
		x=1,y=0;
		return a;
	}
	ll d=exgcd(b,a%b,x,y); 
	ll z=x;x=y;y=z-y*(a/b);
	return d;
}
//O(1)快速乘 
inline ll ksc(ll x,ll y,ll p)
{
	ll z=(ld)x/p*y;
	ll res=(ull)x*y-(ull)z*p;
	return (res+p)%p;
}
const int maxn=1e6;
ing n;//n个方程 
ll a[maxn];//n个方程对应的a
ll m[maxn];//n个方程对应的m 
ll crt()
{
	ll ans=0;//解 
	ll M=1; //对应原理中所说的m,即mi连乘 
	ll Mi;
	ll X,Y;
	for(int i=1;i<=n;i++) M*=m[i];
	for(int i=1;i<=n;i++){
		Mi=M/m[i];
		exgcd(Mi,m[i],X,Y);
		ans+=(ans+ksc(ksc(a[i],Mi,M),X,M))%M; 
	}
	return ans;
}

扩展中国剩余定理：

exCRT用于解决mi不两两互质的情况

使用数学归纳法：首先假设求出前k-1个方程构成的方程组的一个解x，记m=lcm(m1,m2,…,mk-1)

考虑第k个方程，求出一个整数t，使得x+t·m=ak(mod mk)

该方程等价于m·t=ak-x (mod mk)，相当于再解这个线性同余方程，用exgcd

判断有没有解，若无解，则方程组无解

若有解，则对于前k个方程，x’=x+t·m为其一个解

循环使用exgcd，最终就能求出整个方程组的解

板子：（求的是最小非负整数解）

typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
//扩展欧几里得
//输入a,b,任意x,y 
//返回ax+by=gcd(a,b)的一组特解x0,y0 
ll exgcd(ll a,ll b,ll &x,ll &y)
{
	if(b==0){
		x=1,y=0;
		return a;
	}
	ll d=exgcd(b,a%b,x,y); 
	ll z=x;x=y;y=z-y*(a/b);
	return d;
}
//O(1)快速乘 
inline ll ksc(ll x,ll y,ll p)
{
	ll z=(ld)x/p*y;
	ll res=(ull)x*y-(ull)z*p;
	return (res+p)%p;
}
const int maxn=1e6;
int n;//n个方程 
ll a[maxn];//n个方程对应的a
ll m[maxn];//n个方程对应的m
ll exCRT()
{
	ll ans=a[1]%m[1];//前1个方程的解ans 
	ll M=m[1];//前1个方程mi的lcm,对应原理里的m
	for(int i=2;i<=n;i++){
		ll d;//gcd
		ll X,Y;
		ll c=(a[i]%m[i]-ans%m[i]+m[i])%m[i];//对应原理的ak-x
		d=exgcd(M,m[i],X,Y);
		if(c%d!=0) return -1;
		X=ksc(X,c/d,m[i]/d);//对应求解第k个该求的线性同余方程的一个特解
		ll tempM=M; 
		M=M/d*m[i];//前k个mi的lcm 
		ans=((ans+ksc(X,tempM,M))%M+M)%M; 
	}
	return ans; 
}

BSGS：O(sqrt§)

小大步算法，用于解决高次同余方程

高次同余方程有两种形式a^x=b(mod p) 和 x^a=b(mod p)，BSGS用于解决前者，其中a,p互质，求这个同余方程的最小整数解

算法原理：看算法竞赛进阶指南

板子：

typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
//O(1)快速乘 
inline ll ksc(ll x,ll y,ll p)
{
	ll z=(ld)x/p*y;
	ll res=(ull)x*y-(ull)z*p;
	return (res+p)%p;
}
//O(log)快速幂
inline ll ksm(ll x,ll y,ll p)
{
	ll res=1;
	while(y){
		if(y&1) res=ksc(res,x,p);
		x=ksc(x,x,p);y>>=1;
	}
	return res;
}
//计算同余方程的最小非负整数解，无解时返回-1
ll bsgs(ll a,ll b,ll p)
{
	map<ll,ll> mp;mp.clear();
	ll t=(ll)sqrt((double)p)+1;
	b%=p;
	for(ll j=0;j<t;j++){
		ll val=ksc(b,ksm(a,j,p),p);
		mp[val]=j;
	}
	a=ksm(a,t,p);
	if(a==0) return b==0?1:-1;
	for(ll i=0;i<=t;i++){
		ll val=ksm(a,i,p);
		ll j=mp.find(val)==mp.end()?-1:mp[val];
		if(j>=0&&i*t-j>=0) return i*t-j;//不会爆longlong 
	}
	return -1;
}