BZOJ3992 [SDOI2015]序列统计

4 篇文章 0 订阅
1 篇文章 0 订阅

额,跟xuruifan学习了一下午相关知识终于搞过去了

首先我们得知道原根的相关性质,对于一个素数p,他的一个原根g的0~p-2次幂在模p意义下取遍1~p-1的所有值

求原根可以先找出p-1的所有质因数,然后从2开始枚举,如果对任意的p-1的质因数p[i],x的(p-1)/p[i]次幂模p都不等于1,则x是p的一个原根

原根一般比较小,可以暴力找

这样在1~p-1内的每一个值i就唯一对应一个g的幂g^i',这样我们就可以把i*j=k转换成i'+j'=k',然后就变成多项式乘法的形式了

当然这是要求集合里没有0的,不过这题数据错了里边有0,直接无视掉就好了

于是相当于我们要求一个多项式乘法,结果的系数取模,运算过程中下标也要取模

先说结果如何取模,我们可以用NTT,不过对模数有一些要求,具体见网上资料

于是管xuruifan要了一份NTT模板,NTT里的3就是模数的原根

然后下标嘛,如果原来下标都在限制之内,那么乘之后下标都不会超过下标模数的二倍,所以每次乘法的时候再O(n)把超过的地方加到他模之后的位置就行了

于是就NTT+快速幂就过了

#include<iostream>
#include<cstdlib>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<ctime>
#include<algorithm>
#include<iomanip>
#include<vector>
#include<stack>
#include<queue>
#include<map>
#include<set>
#include<bitset>
using namespace std;
#define MAXN 200010
#define MAXM 1010
#define ll long long
#define INF 1000000000
#define MOD 1004535809
#define eps 1e-8
ll G;
ll p[MAXN];
ll a[MAXN],b[MAXN];
ll t1[MAXN],t2[MAXN];
int tot;
int n,X,s;
int L;
int R[MAXN];
int N;
ll mi(ll x,ll y,ll mod){
	if(y<0){
		return mi(mi(x,-y,mod),mod-2,mod);
	}
	ll re=1;
	while(y){
		if(y&1){
			(re*=x)%=mod;
		}
		(x*=x)%=mod;
		y>>=1;
	}
	return re;
}
ll getG(){
	int i,j;
	ll x=n-1;
	for(i=2;x!=1;i++){
		if(!(x%i)){
			p[++tot]=i;
			while(!(x%i)){
				x/=i;
			}
		}
	}
	for(i=2;;i++){
		bool flag=1;
		for(j=1;j<=tot;j++){
			if(mi(i,(n-1)/p[j],n)==1){
				flag=0;
				break;
			}
		}
		if(flag){
			return i;
		}
	}
}
void ntt(ll *a,int f){
	int i,j,k;
	for(i=0;i<n;i++){
		if(i<R[i]){
			swap(a[i],a[R[i]]);
		}
	}
	for(i=1;i<n;i<<=1){
		ll wn=mi(3,f*(MOD-1)/(i<<1),MOD);
		for(j=0;j<n;j+=(i<<1)){
			ll w=1;
			for(k=0;k<i;k++,(w*=wn)%=MOD){
				ll x=a[j+k],y=w*a[j+k+i]%MOD;
				a[j+k]=(x+y)%MOD;
				a[j+k+i]=(x-y+MOD)%MOD;
			}
		}
	}
	if(f==-1){
		ll ni=mi(n,-1,MOD);
		for(i=0;i<n;i++){
			(a[i]*=ni)%=MOD;
		}
	}
}
void mul(ll *c,ll *x,ll *y){
	int i;
	for(i=0;i<n;i++){
		t1[i]=x[i];
		t2[i]=y[i];
	}
	ntt(t1,1);
	ntt(t2,1);
	for(i=0;i<n;i++){
		c[i]=t1[i]*t2[i]%MOD;
	}
	ntt(c,-1);
	for(i=n-1;i>=N;i--){
		(c[i-N]+=c[i])%=MOD;
		c[i]=0;
	}
}
int main(){
	int i,x,y;
	scanf("%d%d%d%d",&y,&n,&X,&s);
	ll G=getG();
	ll now=1;
	for(i=0;i<n-1;i++){
		p[now]=i;
		(now*=G)%=n;
	}
	for(i=1;i<=s;i++){
		scanf("%d",&x);
		if(x){
			a[p[x]]++;
		}
	}
	N=(n-1)<<1;
	for(n=1;n<=N;n<<=1){
		L++;
	}
	N>>=1;
	for(i=0;i<n;i++){
		R[i]=R[i>>1]>>1|((i&1)<<(L-1));
	}
	b[0]=1;
	while(y){
		if(y&1){
			mul(b,a,b);
		}
		mul(a,a,a);
		y>>=1;
	}
	printf("%lld\n",b[p[X]]);
	return 0;
}

/*
2 3 0 3
0 1 2
*/


  • 1
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值