BZOJ3992 [SDOI2015]序列统计

本文版权归ljh2000和博客园共有，欢迎转载，但须保留此声明，并给出原文链接，谢谢合作。

 

 

 

本文作者：ljh2000
作者博客：http://www.cnblogs.com/ljh2000-jump/
转载请注明出处，侵权必究，保留最终解释权！

 

题目链接：BZOJ3992

正解：$NTT$+原根

解题报告：　　

　　考虑朴素$DP$：令$f[i][j]$表示序列的前$i$位的乘积取模为$j$的方案数，那么转移就很$simple$了，就是每次枚举集合内的一个数然后直接转出就好了，这个的复杂度是$O(nm^2)$的。

　　考虑这个乘法的过程很不资瓷啊，不能优化呀，如果变成加法了就可以变成卷积的形式了呀。

　　这是一类巧妙转化：把乘法变成加法卷积->利用原根的性质。考虑$m$为质数，那么必然存在原根$g$，而$m$以内的每个数都可以表示成$g$的次幂，然后就变成加法辣！因为要取模，直接上$NTT$，但还是$O(nmlogm)$的呀，$n$巨大无比...

　　每次转移的多项式都是一样的，直接快速幂就好了。

 

　 　补一发原根的求法：

　　暴力的话就是枚举原根然后按照定义$check$就好了，

　　还有一种做法，就是把$p-1$质因数分解之后，对于$g^{(p-1)/质因子}$$check$一下，如果存在一个为$1$就不是原根，这样的话因为原根一般比较小而且复杂度很优秀，所以能有理有据的跑过去...

 

//It is made by ljh2000
//有志者，事竟成，破釜沉舟，百二秦关终属楚；苦心人，天不负，卧薪尝胆，三千越甲可吞吴。
#include <algorithm>
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <complex>
#include <vector>
#include <cstdio>
#include <string>
#include <bitset>
#include <queue>
#include <cmath>
#include <ctime>
#include <map>
#include <set>
#define lc root<<1
#define rc root<<1|1
#define pr pair<int,int>
#define MP make_pair
#define fr first
#define sc second
#define rep(i,j,k) for(int i=j;i<=k;++i)
#define per(i,j,k) for(int i=j;i>=k;--i)
#define reg(i,x) for(int i=first[x];i;i=next[i])
using namespace std;
typedef long long LL;
typedef long double LB;
typedef complex<double> C;
const double pi = acos(-1);
const double eps = 1e-9;
const int MAXN = 20011;
const int mod = 1004535809;
int middle[MAXN];
inline int fast_pow(int x,int y,int mm){ int r=1; while(y>0) { if(y&1) r=1LL*r*x%mm; x=1LL*x*x%mm; y>>=1; } return r; }
inline int getint(){
    int w=0,q=0; char c=getchar(); while((c<'0'||c>'9') && c!='-') c=getchar();
    if(c=='-') q=1,c=getchar(); while (c>='0'&&c<='9') w=w*10+c-'0',c=getchar(); return q?-w:w;
}

namespace NTT{
	const int G = 3;
	int n,m,a[MAXN],b[MAXN],L,R[MAXN];
	inline void DFT(int *a,int n,int f){
		rep(i,0,n-1) if(i<R[i]) swap(a[i],a[R[i]]);
		int wn,w,x,t;
		for(int i=1;i<n;i<<=1) {
			wn=fast_pow(G,(mod-1)/(i<<1),mod);
			for(int j=0;j<n;j+=(i<<1)) {
				w=1;
				for(int k=0;k<i;k++,w=1LL*w*wn%mod) {
					x=a[j+k]; t=1LL*a[j+i+k]*w%mod;
					a[j+k]=(x+t)%mod;
					a[j+i+k]=(x-t+mod)%mod;
				}
			}
		}
		if(f==1) return ;
		reverse(a+1,a+n);
	}

	inline void work(){
		int savn=n,savm=m;
		m+=n; for(L=0,n=1;n<=m;n<<=1) L++;
		rep(i,0,n-1) R[i]=(R[i>>1]>>1) | ((i&1) << (L-1));
		rep(i,savn+1,n) a[i]=0;
		rep(i,savm+1,n) b[i]=0;
		rep(i,savm+1,n) middle[i]=0;
		DFT(a,n,1); DFT(b,n,1);
		rep(i,0,n) a[i]=1LL*a[i]*b[i]%mod;
		DFT(a,n,-1); int ni=fast_pow(n,mod-2,mod);
		rep(i,0,m) middle[i]=1LL*a[i]*ni%mod;
		rep(i,savm+1,m) middle[i-savm]+=middle[i],middle[i-savm]%=mod;
	}
}

int n,m,final,S,mi[MAXN],rt,match[MAXN];
int pri[MAXN],cc,a[MAXN],ans[MAXN];
bool hav0;
inline void getrt(){
	if(m==2) { rt=1; return ; }
	int x=m-1;
	for(int i=2;i*i<m;i++) {
		if(x%i==0) {
			pri[++cc]=i;
			while(x%i==0) x/=i;
			if(x==1) break;
		}
	}
	if(x>1) pri[++cc]=x;

	bool flag;
	for(rt=2;;rt++) {
		flag=true;
		for(int i=1;i<=cc;i++)
			if(fast_pow(rt,(m-1)/pri[i],m)==1) { flag=false; break; }
		if(flag) return ;
	}
}

inline void work(){
	n=getint(); m=getint(); final=getint(); S=getint(); int x;
	getrt(); mi[0]=1; rep(i,1,m-1) mi[i]=1LL*mi[i-1]*rt%m,match[mi[i]]=i;
	rep(i,1,S) { x=getint(); if(x!=0) a[match[x]]=1; else hav0=true; }
	if(final==0) { if(!hav0) puts("0"); else printf("%d",(fast_pow(S,n,mod)-fast_pow(S-1,n,mod)+mod)%mod); return ; }

	ans[0]=1; int nn=0;
	while(n>0) {
		if(n&1) {
			rep(i,0,nn) NTT::a[i]=ans[i];
			rep(i,0,m-1) NTT::b[i]=a[i];
			NTT::n=nn;
			NTT::m=m-1;
			NTT::work();
			nn=m-1; nn%=m;
			rep(i,0,m-1) ans[i]=middle[i];
		}
		NTT::n=m-1; NTT::m=m-1;
		rep(i,0,m-1) NTT::a[i]=a[i];
		rep(i,0,m-1) NTT::b[i]=a[i];
		NTT::work();
		rep(i,0,m-1) a[i]=middle[i];
		n>>=1;
	}

	printf("%d",ans[match[final]]);
}

int main()
{
#ifndef ONLINE_JUDGE
	freopen("3992.in","r",stdin);
	freopen("3992.out","w",stdout);
#endif
    work();
    return 0;
}
//有志者，事竟成，破釜沉舟，百二秦关终属楚；苦心人，天不负，卧薪尝胆，三千越甲可吞吴。