20200601数论总结

最新推荐文章于 2022-08-24 21:17:50 发布

cqbzcsq

最新推荐文章于 2022-08-24 21:17:50 发布

阅读量747

点赞数

分类专栏：总结数学文章标签： c++ 算法

本文链接：https://blog.csdn.net/C20180602_csq/article/details/106480070

版权

数学同时被 2 个专栏收录

43 篇文章 0 订阅

订阅专栏

总结

40 篇文章 1 订阅

订阅专栏

exgcd

解不定方程ax+by=gcd(a,b)

bx+(a%b)y=gcd(b,a%b)=gcd(a,b)

bx+(a-(a/b)*b)y=gcd(a,b)

ay+bx-(a/b)*by=gcd(a,b)

ay+b(x-(a/b)*y)=gcd(a,b)

递归即可

excrt

$x\equiv x_1(mod X)$

$x\equiv x_2(modY)$

$x=x_1+t_1X$

$x=x_2+t_2Y$

$x_1+t_1X=x_2+t_2Y$

$t_1X=t_2Y+x_2-x_1$

有贝祖定理可知，gcd(X,Y)|（x2-x1）

两边同时除一个g=gcd(X,Y)

$t_1\frac{X}{g}=t_2\frac{Y}{g}+\frac{x_2-x_1}{g}$

写成mod Y/g的形式

$t_1\frac{X}{g}=\frac{x_2-x_1}{g}(mod\frac{Y}{g})$

此时X/g与Y/g互质，存在X/g的逆元（用exgcd求逆元），则可以求出t1

$t_1=\frac{x_2-x_1}{g}*(\frac{X}{g})^{-1}(mod\frac{Y}{g})$

$t_1=\frac{x_2-x_1}{g}*(\frac{X}{g})^{-1}+k*\frac{Y}{g}$

那么带入x=x1+t1*X

$x=x_1+X*\frac{x_2-x_1}{g}*(\frac{X}{g})^{-1}+k*\frac{XY}{g}$

再写成mod (XY)/g的形式

$x=x_1+X*\frac{x_2-x_1}{g}*(\frac{X}{g})^{-1}(mod\frac{XY}{g})$

这样我们就合并了两个同余方程

注意这里的逆元是在modY/g意义下的

例题：[NOI2018]屠龙勇士

先用exgcd手动解出kx=a(mod p)的解，再用excrt合并即可

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<set>
using namespace std;
#define N 100005
#define LL long long
inline LL gi()
{
	char c;LL num=0,flg=1;
	while((c=getchar())<'0'||c>'9')if(c=='-')flg=-1;
	while(c>='0'&&c<='9'){num=num*10+1ll*c-48;c=getchar();}
	return num*flg;
}
LL a[N],nw[N],p[N];
multiset<LL> S;
multiset<LL>::iterator it;
void exgcd(LL a,LL b,LL &gcd,LL &x,LL &y)
{
	if(!b){x=1;y=0;gcd=a;return;}
	exgcd(b,a%b,gcd,y,x);
	y-=x*(a/b);
}
LL mul(LL x,LL y,LL mod)
{
	if(x<0)x=(x%mod+mod)%mod;
	if(y<0)y=(y%mod+mod)%mod;
	if(x<y)swap(x,y);
	LL ret=0;
	while(y){
		if(y&1){ret=ret+x;if(ret>=mod)ret-=mod;}
		y>>=1;x<<=1;if(x>=mod)x-=mod;
	}
	return ret;
}
LL cans,cmod;
int main()
{
	int T,n,m,i;
	LL mx,b;
	T=gi();
	while(T--){
		n=gi();m=gi();S.clear();
		for(i=1;i<=n;i++)a[i]=gi();
		for(i=1;i<=n;i++)p[i]=gi();
		for(i=1;i<=n;i++)nw[i]=gi();
		for(i=1;i<=m;i++)S.insert(gi());
		cans=0;cmod=1;mx=0;
		bool flg=0;
		for(i=1;i<=n;i++){
			it=S.upper_bound(a[i]);
			if(it!=S.begin())it--;
			b=*it;S.erase(it);S.insert(nw[i]);
			mx=max(mx,(a[i]-1)/b+1);
			b%=p[i];a[i]%=p[i];//b*x=a[i](mod p[i])
			if(b==0&&a[i]==0)continue;
			if(b==0&&a[i]!=0){flg=1;break;}
			LL g,x,y,t1,t2;
			exgcd(b,p[i],g,x,y);//b*x+p[i]*y=g
			if(a[i]%g!=0){flg=1;break;}
			p[i]/=g;
			x=mul(a[i]/g,x,p[i]);
			
			exgcd(cmod,p[i],g,t1,t2);//cmod*t1+p[i]*t2=g
			if((x-cans)%g!=0){flg=1;break;}
			LL A=cmod/g*p[i];
			cans=(cans+mul(t1,mul(cmod,(x-cans)/g,A),A))%A;
			cmod=A;
		}
		if(flg==1) printf("-1\n");
		else{
			if(cans>=mx)printf("%lld\n",cans);
			else printf("%lld\n",cans+cmod*((mx-cans-1)/cmod+1));
		}
	}
}

exlucas

求C(n,m)mod X，时间复杂度为O(Σp^k)

问题可以转化为求

n!%(p^k)

由于p的倍数在mod p^k下无逆元

我们可以维护二元组(x,y)表示x*p^y，其中x中不含有p因子，可以直接模

先把p的倍数单独拿出来计算到y中，发现可以化为子问题

考虑怎么求x，预处理出p^k以内的所有数除去p倍数的阶乘

由于循环节大小为p^k，循环的部分可以直接快速幂，多出来的部分直接用预处理的阶乘值

最后用excrt合并一下就好了，合并的时候注意要把每一个答案的所有pi因子去掉，再乘回来

例题：Kapitał

（最后使用的手动crt合并，细节略多。。。调了我好久。。。）

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int mod1=512;
const int mod2=1953125;
const int MOD=1000000000;
int fac1[mod1+5],fac2[mod2+5];
#define LL long long
int ksm(int x,LL y,int mod)
{
	int ret=1;
	while(y){
		if(y&1)ret=1ll*ret*x%mod;
		y>>=1;x=1ll*x*x%mod;
	}
	return ret;
}
struct node{int x;LL y;node(int _x,LL _y){x=_x;y=_y;}};
node mul(node x,node y,int mod){return node(1ll*x.x*y.x%mod,x.y+y.y);}
node dvd(node x,node y,int mod,int phi){return node(1ll*x.x*ksm(y.x,phi-1,mod)%mod,x.y-y.y);}
node getfac(LL x,int p,int mod,int fac[])
{
	if(!x)return node(1,0);
	LL tmp=1ll*fac[x%mod]*ksm(fac[mod-1],x/mod,mod)%mod;
	return mul(node(tmp,x/p),getfac(x/p,p,mod,fac),mod);
}
int main()
{
	LL n,m;int k,i;
	scanf("%lld%lld%d",&n,&m,&k);
	fac1[0]=fac2[0]=1;
	for(i=1;i<mod1;i++)fac1[i]=1ll*fac1[i-1]*((i%2)?i:1)%mod1;
	for(i=1;i<mod2;i++)fac2[i]=1ll*fac2[i-1]*((i%5)?i:1)%mod2;
	node ans1=getfac(n+m,2,mod1,fac1);
	ans1=dvd(ans1,getfac(n,2,mod1,fac1),mod1,mod1/2);
	ans1=dvd(ans1,getfac(m,2,mod1,fac1),mod1,mod1/2);
	node ans2=getfac(n+m,5,mod2,fac2);
	ans2=dvd(ans2,getfac(n,5,mod2,fac2),mod2,mod2/5*4);
	ans2=dvd(ans2,getfac(m,5,mod2,fac2),mod2,mod2/5*4);
	LL tmp=min(ans1.y,ans2.y);
	LL x=ksm(mod2,mod1/2-1,mod1),y=ksm(mod1,mod2/5*4-1,mod2);
	
	ans1.x=1ll*ans1.x*ksm(ksm(5,mod1/2-1,mod1),tmp,mod1)%mod1;
	ans2.x=1ll*ans2.x*ksm(ksm(2,mod2/5*4-1,mod2),tmp,mod2)%mod2;
	ans1.y-=tmp;ans2.y-=tmp;
	ans1.x=1ll*ans1.x*ksm(2,ans1.y,mod1)%mod1;
	ans2.x=1ll*ans2.x*ksm(5,ans2.y,mod2)%mod2;
	int ans=(1ll*ans1.x*x%MOD*mod2%MOD+1ll*ans2.x*y%MOD*mod1%MOD)%MOD;
	for(i=ksm(10,k-1,1000000007);i>=1;i/=10)
		printf("%d",(ans/i)%10);
}

线段树+线性基

维护一个序列，支持区间异或上某个值k，求一个区间的数的线性基大小

线段树+线性基

注意，一个线性基中所有的数异或上k再重构线性基，不等于把所有的数单独异或上k，再一个一个插入线性基

我们可以维护一个区间数两两异或构成的线性基，这样一次修改对整个线性基是没有影响的

考虑把所有线段树节点的最右边的点拿出来，维护他们的值，合并的时候把左右儿子的线性基合并，再把这两个特征值加入父亲的线性基即可，这样就维护出了所有点两两异或的线性基

查询的时候任意将一个单点值插入答案线性基，再把所有的线性基两两合并，相邻节点的特征值合并，就可以得到奇数个数与偶数个数的异或线性基（因为只要存在一个单点，它异或进一个由偶数个数异或的方案，就可以得到所有奇数个数的方案）

代码：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
inline int gi()
{
	char c;int num=0,flg=1;
	while((c=getchar())<'0'||c>'9')if(c=='-')flg=-1;
	while(c>='0'&&c<='9'){num=num*10+c-48;c=getchar();}
	return num*flg;
}
#define N 200005
#define LOG 29
#define lc i<<1
#define rc i<<1|1
int bastmp[LOG+2],cnttmp;
struct node{
	int l,r,bas[LOG+1],cnt,la;int val;//bool flg;
	void add(int x){
		//if(!x){flg=1;return;}
		for(int i=LOG;i>=0;i--)if(x&(1<<i)){
			if(!bas[i]){bas[i]=x,cnt++;return;}
			else{x^=bas[i];if(!x)return;}
		}
		//flg=1;
	}
	/*void modify(int k){// force O(30^2)
		cnttmp=0;cnt=0;
		for(int i=LOG;i>=0;i--)if(bas[i])
			bastmp[++cnttmp]=bas[i]^k,bas[i]=0;
		for(int i=1;i<=cnttmp;i++)add(bastmp[i]);
		if(flg)add(k);
	}*/
}a[N<<2];
int val[N];
void pushup(int i)
{
	int LC=lc,RC=rc;
	if(a[LC].cnt<a[RC].cnt)swap(LC,RC);
	a[i].cnt=a[LC].cnt;//a[i].flg=a[LC].flg|a[RC].flg;
	for(int j=LOG;j>=0;j--)a[i].bas[j]=a[LC].bas[j];
	for(int j=LOG;j>=0;j--)if(a[RC].bas[j])
		a[i].add(a[RC].bas[j]);
	a[i].add(a[LC].val^a[RC].val);
	a[i].val=a[rc].val;
}
void pushdown(int i)// O(900*2)
{
	if(a[i].la){
		//a[lc].modify(a[i].la);
		a[lc].val^=a[i].la;a[lc].la^=a[i].la;
		//a[rc].modify(a[i].la);
		a[rc].val^=a[i].la;a[rc].la^=a[i].la;
		a[i].la=0;
	}
}
void build(int i,int l,int r)
{
	a[i].l=l;a[i].r=r;
	if(l==r){a[i].val=val[l];return;}
	int mid=(l+r)>>1;
	build(lc,l,mid);build(rc,mid+1,r);
	pushup(i);
}
void insert(int i,int l,int r,int k)
{
	if(a[i].l>r||a[i].r<l)return;
	if(l<=a[i].l&&a[i].r<=r){
		//a[i].modify(k);
		a[i].val^=k;a[i].la^=k;
		return;
	}
	pushdown(i);//O(logn*1800)
	insert(lc,l,r,k);insert(rc,l,r,k);
	pushup(i);
}
void query(int i,int l,int r)
{
	if(a[i].l>r||a[i].r<l)return;
	if(l<=a[i].l&&a[i].r<=r){
		for(int j=LOG;j>=0;j--)if(a[i].bas[j])
			a[0].add(a[i].bas[j]);
		a[0].add(a[0].val^a[i].val);
		a[0].val=a[i].val;
		return;
	}
	pushdown(i);
	query(lc,l,r);query(rc,l,r);
}
int main()
{
	int n,Q,i,op,k,l,r;
	n=gi();Q=gi();
	for(i=1;i<=n;i++)val[i]=gi();
	build(1,1,n);
	while(Q--){
		op=gi();l=gi();r=gi();
		if(op==1){k=gi();insert(1,l,r,k);}
		else{
			memset(a[0].bas,0,sizeof(a[0].bas));a[0].cnt=a[0].val=0;
			query(1,l,r);
			printf("%d\n",1<<a[0].cnt);
		}
	}
}

生成树求和加强版

先拆位，然后用变元矩阵树定理

把边权为0的看做1，为1的看做x，为2的看做x^2

如果我们用矩阵树定理求出最后行列式（是一个多项式）

再求它在mod x^3意义下的多项式，那么权值和就是(2*a2+a1)*3^t（ai表示最后x^i的系数）

我们可以带入三次单位根，用二元组(x,y)表示x+y*w31（w30=1,w32=-w31-1）

这就是一个自然的3次循环卷积

我们通过矩阵树定理可以求出三个二元组，分别表示带入w30,w31,w32时的答案多项式

现在的任务就是插值求出原多项式

我们有三个方程

$a+bw_3^0+c(w_3^0)^2=a_0+b_0w_3^1$

$a+bw_3^1+c(w_3^1)^2=a_1+b_1w_3^1$

$a+bw_3^2+c(w_3^2)^2=a_2+b_2w_3^1$

化简一下

$a+b+c=a_0+b_0w_3^1$

$a+bw_3^1+c(-w_3^1-1)=a_1+b_1w_3^1$

根据待定系数法，b0必定等于0，其实剩下就有了3个方程，所以是可以解出a,b,c的

$a+b+c=a_0$

$a-c=a_1$

$b-c=b_1$

稍微解一下就好了

这样就可以不用写IDFT（当然只在循环卷积次数很小的时候适用，因为要手解方程）

只用做两次矩阵树，目前是LOJ最快(20200601)

代码：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define N 105
#define M 5005
const int mod=1000000007;
const int inv3=333333336;
int Pow(int a,int b){int s=1;for(;b;b>>=1,a=1ll*a*a%mod) b&1&&(s=1ll*s*a%mod);return s;}
struct cp{
	int r,i;
	cp(){}
	cp(int x,int y){r=x;i=y;}
	cp operator + (const cp &t)const{return cp((r+t.r)%mod,(i+t.i)%mod);}
	cp operator - (const cp &t)const{return cp((r+mod-t.r)%mod,(i+mod-t.i)%mod);}
	cp operator * (const cp &t)const{
		return cp((1ll*r*t.r+mod-1ll*i*t.i%mod)%mod,(1ll*r*t.i+1ll*i*t.r+mod-1ll*i*t.i%mod)%mod);
	}
	cp getinv()const{
		int inv=Pow((1ll*r*r+1ll*i*i+mod-1ll*r*i%mod)%mod,mod-2);
		return cp(1ll*(r+mod-i)*inv%mod,1ll*(mod-i)*inv%mod);
	}
	bool empty(){return !r&&!i;}
}a[N][N],w[4],F[3],G[3];
int n,m,ans,e[M][3],mx;
cp det()
{
	int i,j,k;
	cp ret=cp(1,0);
	for(i=1;i<n;i++){
		if(a[i][i].empty()){
			for(j=i+1;j<n;j++){
				if(!a[j][i].empty()){
					for(k=1;k<=n;k++)
						swap(a[i][k],a[j][k]);
					ret=cp(0,0)-ret;
					break;
				}
			}
		}
		if(a[i][i].empty()) return cp(0,0);
		ret=ret*a[i][i];
		cp inv=a[i][i].getinv();
		for(j=i+1;j<n;j++) if(!a[j][i].empty()){
			cp t=a[j][i]*inv;
			for(k=i;k<n;k++) a[j][k]=a[j][k]-a[i][k]*t;
		}
	}
	return ret;
}
int main()
{
	freopen("sum.in","r",stdin);
	freopen("sum.out","w",stdout);
	int i,j,k,pw;
	w[0]=w[3]=cp(1,0);//1
	w[1]=cp(0,1);//w3
	w[2]=cp(mod-1,mod-1);//w3^2=-1-w3
	scanf("%d%d",&n,&m);
	for(i=1;i<=m;i++){
		scanf("%d%d%d",&e[i][0],&e[i][1],&e[i][2]);
		mx=max(mx,e[i][2]);
	}
	for(pw=1;pw<=mx;pw*=3){
		for(k=0;k<2;k++){
			memset(a,0,sizeof a);
			for(i=1;i<=m;i++){
				cp c=w[e[i][2]/pw%3*k%3];
				int u=e[i][0],v=e[i][1];
				a[u][v]=a[u][v]-c;a[v][u]=a[v][u]-c;
				a[u][u]=a[u][u]+c;a[v][v]=a[v][v]+c;
			}
			F[k]=det();
		}
		//printf("%d %d\n",F[0].r,F[0].i);
		//printf("%d %d\n",F[1].r,F[1].i);
		//if(F[0].i!=0){printf("problem\n");}
		int A=1ll*(F[0].r+2ll*mod-F[1].r-F[1].i)%mod*inv3%mod;
		int B=(A+F[1].i)%mod;
		ans=(ans+(B+A*2ll)%mod*pw)%mod;
	}
	printf("%d\n",ans);
}

线性基与贪心结合

1、离线+线性基，在线性基中维护每个值的删除时间，贪心地更新线性基，这样可以做到线性基中删除，例题:ZROI或许

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>
using namespace std;
inline int gi()
{
	char c;int num=0,flg=1;
	while((c=getchar())<'0'||c>'9')if(c=='-')flg=-1;
	while(c>='0'&&c<='9'){num=num*10+c-48;c=getchar();}
	return num*flg;
}
#define N 2000005
int n,Q,x[N],op[N],ed[N];
map<int,int> mp;
int bas[35],id[35],cnt;
void add(int x,int pos)
{
	for(int i=n-1;i>=0;i--){
		if(x&(1<<i)){
			if(!bas[i]){bas[i]=x;id[i]=pos;cnt++;return;}
			else{
				if(id[i]<pos)swap(x,bas[i]),swap(pos,id[i]);
				x^=bas[i];
			}
		}
	}
}
void del(int pos)
{
	for(int i=n-1;i>=0;i--){
		if(id[i]==pos){
			bas[i]=id[i]=0;cnt--;
			return;
		}
	}
}
int ans;
int main()
{
	freopen("A.in","r",stdin);
	freopen("A.out","w",stdout);
	int i;
	n=gi();Q=gi();
	for(i=1;i<=Q;i++){
		op[i]=gi();x[i]=gi();
		if(op[i]==1)
			mp[x[i]]=i,ed[i]=Q+1;
		else
			ed[mp[x[i]]]=i;
	}
	for(i=1;i<=Q;i++){
		if(op[i]==1)add(x[i],ed[i]);
		else del(i);
		ans^=(1<<(n-cnt));
	}
	printf("%d\n",ans);
}

2、离线+线性基，对每个r维护线性基，在线性基中维护每个点的坐标，用靠后的点代替靠前的点，可以做到O(nlogn)区间线性基查询，例题：Ivan and Burgers

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
inline int gi()
{
	char c;int num=0,flg=1;
	while((c=getchar())<'0'||c>'9')if(c=='-')flg=-1;
	while(c>='0'&&c<='9'){num=num*10+c-48;c=getchar();}
	return num*flg;
}
#define N 500005
#define LOG 19
int a[N];
int bas[LOG+2],id[N];
struct node{
	int l,r,id;
	bool operator < (const node &t)const{return r<t.r||(r==t.r&&l<t.l);}
}q[N];
int ans[N];
void insert(int x,int pos)
{
	for(int i=LOG;i>=0;i--){
		if(x&(1<<i)){
			if(!bas[i]){bas[i]=x,id[i]=pos;break;}
			else{
				if(id[i]<pos)swap(x,bas[i]),swap(pos,id[i]);
				x^=bas[i];
			}
		}
	}
}
int query(int pos)
{
	int mx=0;
	for(int i=LOG;i>=0;i--)
		if(bas[i]&&id[i]>=pos)
			mx=max(mx,mx^bas[i]);
	return mx;
}
int main()
{
	int n,i,j,Q;
	n=gi();for(i=1;i<=n;i++)a[i]=gi();
	Q=gi();for(i=1;i<=Q;i++){q[i].l=gi();q[i].r=gi();q[i].id=i;}
	sort(q+1,q+Q+1);
	for(i=1,j=1;i<=n;i++){
		insert(a[i],i);
		while(j<=Q&&q[j].r==i)
			ans[q[j].id]=query(q[j].l),j++;
	}
	for(i=1;i<=Q;i++)printf("%d\n",ans[i]);
}

exBSGS

普通的BSGS可以解决模数p与底数a互质的情况

如果不互质的话可以设g=gcd(a,p)

原方程a^x=b(mod p)

写成不定方程的形式:a^x+kp=b

a*a^(x-1)+kp=b

发现a，p都是常数，则方程有解当且仅当g|b

所以我们两边除一个g

(a/g)*a^(x-1)+k(p/g)=b/g

再写回mod的形式

(a/g)*a^(x-1)=b/g (mod p/g)

这样我们再来判断gcd(a,p/g)是否整除b/g，递归下去直到a，p互质，就可以用普通BSGS来解决了

注意判断、输出答案的时候还要把除掉的东西乘回来

模板：

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
#include<map>
using namespace std;
int ksm(int x,int y,int mod)
{
    int ans=1;
    while(y){
        if(y&1) ans=1ll*ans*x%mod;
        y>>=1;x=1ll*x*x%mod;
    }
    return ans;
}
int gcd(int x,int y){return !y?x:gcd(y,x%y);}
int bsgs(int a,int b,int p)
{
    a%=p;b%=p;
    int cnt=0,t=1,g;
    if(b==t)return cnt;
    for(g=gcd(a,p);g!=1;g=gcd(a,p)){
        if(b%g) return -1;
        t=1ll*t*(a/g)%p;
        b/=g;p/=g;
        cnt++;
        if(b==t)return cnt;
    }
    map<int,int>mp;
    int m=int(sqrt(1.0*p)+0.999999),A,B,l,r,tmp;
    r=b;
    for(B=0;B<m;B++){
        mp[r]=B;
        r=1ll*r*a%p;
    }
    l=t;tmp=ksm(a,m,p);
    for(A=1;A<=m+1;A++){
        l=1ll*l*tmp%p;
        if(mp.count(l))
            return A*m-mp[l]+cnt;
    }
    return -1;
}
int main()
{
    int a,b,p,ans;
    while(~scanf("%d%d%d",&a,&p,&b)){
        if(!a&&!b&&!p) return 0;
        if(p==1)ans=0;
        else ans=bsgs(a,b,p);
        if(ans==-1) printf("No Solution\n");
        else printf("%d\n",ans);
    }
}

cqbzcsq

关注

0
点赞
踩
0

收藏

觉得还不错? 一键收藏
0
评论
20200601数论总结

exgcd解不定方程ax+by=gcd(a,b)bx+(a%b)y=gcd(b,a%b)=gcd(a,b)bx+(a-(a/b)*b)y=gcd(a,b)ay+bx-(a/b)*by=gcd(a,b)ay+b(x-(a/b)*y)=gcd(a,b)递归即可excrt有贝祖定理可知，gcd(X,Y)|（x2-x1）两边同时除一个g=gcd(X,Y)写成mod Y/g的形式此时X/g与Y/g互质，存在X/g的逆元（用exgc..
复制链接

扫一扫