Math

发布时间 2023-07-29 18:28:47作者: EnderWave

数论

Part1 线性筛

int p[N],n,k,q,cnt;
bool vis[N];
void init(){
	for(int i=2;i<=n;++i){
		if(!vis[i]){
			p[++cnt]=i;
			vis[i]=1;
		}
		for(int j=1;j<=cnt;j++){
			if(i*p[j]>n)break;
			if(i%p[j]==0){
				vis[i*p[j]]=1;
				break;
			}
			vis[i*p[j]]=1;
		}
	}
}

Part2 gcd

GCD

int gcd(int a,int b){
	return b==0?a:gcd(b,a%b);
}

EXGCD

void exgcd(ll a,ll b,ll&x,ll&y){
	if(b==0){
		x=1;y=0;
		return;
	}
	exgcd(b,a%b,x,y);
	ll z=x;x=y;y=z-(a/b)*y;
}

Part3 中国剩余定理

CRT

for(int i=1;i<=n;i++){
		M[i]=mul/m[i];ll x=0,y=0;
		exgcd(M[i],m[i],x,y);if(x<0)x=x+m[i];
		ans=(ans+(a[i]*M[i]*x)%mul)%mul;
	}

EXCRT

ll mul(ll a,ll b,ll p){
	ll res=0;
	while(b){
		if(b&1){
			res=(res+a)%p; 
		}
		b>>=1;
		a=(a+a)%p;
	}
	return res;
} 
ll excrt(){
	ll M,ans,x,y;
	M=a[1];ans=b[1];
	for(int i=2;i<=n;i++){
		ll c=(b[i]-ans%a[i]+a[i])%a[i];
		exgcd(M,a[i],x,y);
		ll ag=a[i]/d;
		x=mul(x,c/d,ag);
		if(c%d!=0)return -1;
		ans+=x*M;M*=ag;
		ans=(ans%M+M)%M;
	}
	return ans;
}

Part4 Lucas定理

ll pow(ll a,ll b){
	res=1;
	while(b){
		if(b&1)res=(res*a)%p;
		a=(a*a)%p;
		b>>=1;
	}
	return res;
}
ll C(ll n,ll m){
	if(m>n)return 0;
	return ((a[n]*pow(a[m],p-2))%p*pow(a[n-m],p-2))%p;
}
ll Lucas(ll n,ll m){
	if(!m)return 1;
	return C(n%p,m%p)*Lucas(n/p,m/p)%p;
}

多项式

FFT

struct Complex{
	double x,y;
	Complex(double x=0,double y=0):x(x),y(y){
	}
}a[N],b[N];
Complex operator *(Complex a,Complex b){
	return Complex{a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};
}
Complex operator +(Complex a,Complex b){
	return Complex{a.x+b.x,a.y+b.y};
}
Complex operator -(Complex a,Complex b){
	return Complex{a.x-b.x,a.y-b.y};
}
void FFT(Complex*a,int type){
	for(int i=0;i<lim;i++)
		if(i<r[i])
			std::swap(a[i],a[r[i]]);
	for(int mid=1;mid<lim;mid<<=1){
		Complex wn(cos(PI/mid),(double)type*sin(PI/mid));
		for(int pos=0,len=mid<<1;pos<lim;pos+=len){
			Complex w(1,0);
			for(int k=0;k<mid;k++,w=w*wn){
				Complex x=a[pos+k];
				Complex y=w*a[pos+k+mid];
				a[pos+k]=x+y;
				a[pos+k+mid]=x-y;
			}
		}
	}
	if(type==1)return;
	for(int i=0;i<=lim;i++)
		a[i].x/=lim;
}