[SDOI2015] 序列统计

Description

给你一个模 (m) 意义下的数集，需要用这个数集生成一个数列，使得这个数列在的乘积为 (x)。
问方案数模 (1004535809)。 保证 (m) 是质数。

Solution

首先要知道 (1004535809) 是个 (NTT) 常用模数，原根为 (3)。

然后推一推 (DP) 式子，设 (f[i][j]) 表示填好了前 (i) 个数，乘积为 (j) 的方案数。

转移显然 (f[i][j]=sumlimits_{p imes k=j}f[i-1][p] imes sum[k]) 。其中 (sum[i]) 表示数集中 (i) 是否出现。

然鹅这样做是 (O(m^2n)) 的

观察到如果把转移写成矩阵，那每一次的转移矩阵都是相等的，这启发我们进行快速幂，复杂度变为 (O(m^2log n))。

发现这个 (p imes k=j) 不好处理，如果把它变成 (p+k=j) 的形式就可以 (NTT) 优化了。

(m) 是质数，它的原根数量为 (phi(m-1)) 个。

(mge 3)。

相当于保证 (m) 一定有原根。

根据原根的性质，(g^0,g^1,g^2dots g^{m-2}) 两两互不相同。

那 (DP) 式子可以改写为 $$f[i][g^a]=sum_{g^b imes g^c=g^a} f[i-1][g^b] imes sum[g^c]$$

[f[i][g^a]=sumlimits_{g^{b+c}=g^a} f[i-1][g^b] imes sum[g^c] ]

[f[i][a]=sumlimits_{b+c=a} f[i-1][b] imes sum[c] ]

快速幂的时候 (NTT) 优化即可。

哦对了，因为这个式子只有 (g^0sim g^{m-2}) 总共 (m-1) 项，所以下标实际上都应该对 (m-1) 取模，因为 (0) 是没有对应的项的。

时间复杂度 (O(mlog nlog m))

求原根的时候对 (m-1) 因数分解，然后暴力判是否可行即可。

Code

#include<bits/stdc++.h>
using std::min;
using std::max;
using std::swap;
using std::vector;
typedef double db;
typedef long long ll;
#define pb(A) push_back(A)
#define pii std::pair<int,int>
#define all(A) A.begin(),A.end()
#define mp(A,B) std::make_pair(A,B)
#define int long long
const int G=3;
const int N=8005*8;
const int mod=1004535809;

int lim,rev[N];
int n,m,X,g;
int p[N],pc;
int pos[N],in;
int a[N],b[N],c[N];

int ksm(int a,int b,int mod,int ans=1){
    while(b){
        if(b&1) ans=ans*a%mod;
        a=a*a%mod;b>>=1;
    } return ans;
}

int getg(){
    int Phi=m-1,sq=sqrt(Phi);
    for(int i=2;i<=sq;i++){
        if(Phi%i==0) p[++pc]=i,i*i!=Phi?p[++pc]=Phi/i:0;
    } for(int i=2;i<=Phi;i++){
        int flag=0;
        for(int j=1;j<=pc;j++){
            if(ksm(i,p[j],m)==1) {
                flag=1;
                break;
            }
        } if(!flag) return i;
    }
}

int inv(int x){
    return ksm(x,mod-2,mod);
}

int getint(){
    int X=0,w=0;char ch=getchar();
    while(!isdigit(ch))w|=ch=='-',ch=getchar();
    while( isdigit(ch))X=X*10+ch-48,ch=getchar();
    if(w) return -X;return X;
}

void ntt(int *f,int opt){
    for(int i=0;i<lim;i++) if(i<rev[i]) swap(f[i],f[rev[i]]);
    for(int mid=1;mid<lim;mid<<=1){
        int tmp=ksm(G,(mod-1)/(mid<<1),mod);
        if(opt==-1) tmp=inv(tmp);
        for(int R=mid<<1,j=0;j<lim;j+=R){
            int w=1;
            for(int k=0;k<mid;k++,w=w*tmp%mod){
                int x=f[j+k],y=f[j+k+mid]*w%mod;
                f[j+k]=(x+y)%mod;f[j+k+mid]=(mod+x-y)%mod;
            }
        }
    }
}

void mul(int *a,int *b){
    for(int i=0;i<m;i++) c[i]=b[i];
    ntt(a,1),ntt(c,1);
    for(int i=0;i<lim;i++) a[i]=a[i]*c[i]%mod;
    ntt(a,-1);
    for(int i=0;i<lim;i++) a[i]=a[i]*in%mod;
    for(int i=0;i<m;i++) a[i]=(a[i]+a[i+m])%mod;
    for(int i=m;i<lim;i++) a[i]=c[i]=0;
}

void sqr(int *a){
    ntt(a,1);
    for(int i=0;i<lim;i++) a[i]=a[i]*a[i]%mod;
    ntt(a,-1);
    for(int i=0;i<lim;i++) a[i]=a[i]*in%mod;
    for(int i=0;i<m;i++) a[i]=(a[i]+a[i+m])%mod;
    for(int i=m;i<lim;i++) a[i]=0;
}

void ksm(int *ans,int *a,int b){
    while(b){
        if(b&1) mul(ans,a);
        sqr(a);b>>=1;
    }
}

signed main(){
    n=getint(),m=getint(),X=getint();
    g=getg();int now=1;
    pos[1]=0;for(int i=1;i<m-1;i++) pos[(now*=g)%=m]=i;
    lim=1;while(lim<=m) lim<<=1;lim<<=1;in=inv(lim);
    for(int i=0;i<lim;i++) rev[i]=(rev[i>>1]>>1)|(i&1?lim>>1:0);
    int len=getint();
    for(int i=1;i<=len;i++){
        int x=getint();
        if(x) b[pos[x]]++;
    } a[0]=1;m--;
    ksm(a,b,n);
    printf("%lld
",a[pos[X]]);
    return 0;
}