[SDOI2015]序列统计

题目链接

一句话题意：给定一个集合$S$，求用其中的元素生成一个长度为$N$的序列，求满足$prod_{i=1}^{n}a_i=X a_iin S$的序列的方案数,其中$X$给定，且保证集合$S$中所有数都小于$m$且$m$为质数。

题解：

　　1.首先考虑DP

　　　　设dp[i][j]表示前i个数已经选完他们的乘积为j的方案数，转移方程很好想，显然是$dp[i][j]=sum _{l=0}^{m} sum_{r*l=j(mod m),rin S}dp[i-1][l]$

　　　　复杂度：$O(nm^2)$显然空间时间都吃不消。

　　2.然后考虑优化DP

　　　　首先看一眼DP转移式第一眼发现dp[i][j]的转移只与i-1有关，因此滚动，这样空间问题就会解决（时间？？？）

　　　　然后再看一眼DP转移方程，发现转移方程很灵性，它存在$dp[i*j (mod m)]=sum dp[i]*1 jin S$的形式，接下来是这道题最神的地方，因为m为质数，可以证明m必定存在一个原根，而质数原根的0到m-2次方在模m意义下分别对应着1到m-1这些数，也就是说我们可以把1到m-1这些数分别用原根的0到m-2次幂表示出来（有什么用呢？？？）

　　　　

　　　　这样表示完之后，根据$g^{n}*g^{m}=g^{n+m}$的公式我们就可以把转移式中的乘法转化为加法啦。

　　　　注意这里有一个问题，就是当i+j大于m-2时怎么办，解决办法其实很简单，直接把i+j模上m-1就可以啦。证明 根据扩展欧拉定理$g^a=g^{a(mod phi(m))}(mod m) phi(m)=m-1$

　　　　这样转移方程变为$dp[i+j(mod m-1)]=sum _{G^jin S}dp[i]$,设计一个数组g，当且仅当$G^iin S$时g[i]=1。这样方程又变成了$dp[i+j (mod m-1)]=sum dp[i]*g[j]$显然是卷积的形式，可以将转移复杂度优化到$mlogm$但总时间复杂度为$nmlongm$还是过不了（难受啊）

　　　　最后我们可以发现，对于n次转移每次乘的g数组都是相同的因此我们可以直接求出g数组n次卷积后的结果再乘上dp数组，对于出g数组n次卷积的求法我们可以直接使用快速幂进行处理这样复杂度就变成了$mlogmlogn$可过

具体细节还是看代码吧（觉得这道题特别神奇，fft优化dp，处理i*j形式的卷积）（第一次没照题解打的NTT题.......）

#include<bits/stdc++.h>
#define ll long long 
using namespace std;
ll mod=1004535809,n=1,lim,m,mm,x,slen,s[10000],G=3,inv;int rev[1000000],mp[1000000],j,l;
ll a[1000000],b[1000000],dp[1000000],g[1000000],ttt[1000000],invv[1000000];
inline ll read(){
    ll x=0,f=1;char s=0;
    while(s<'0'||s>'9'){if(s=='-')f=-1;s=getchar();}
    while(s>='0'&&s<='9'){x=x*10+s-'0';s=getchar();}
    return x*f;
}
inline ll ksm(ll a,ll b,ll p){a%=p;ll tmp=1;for(;b;b>>=1,a=a*a%p)if(b&1)tmp=tmp*a%p;return tmp;}
inline void init(){
    for(register int i=0;i<n;++i){
        int t=0;for(j=0;j<lim;++j)
        if((i>>j)&1)t|=(1<<(lim-1-j));
        rev[i]=t;
    }
    for(l=2;l<=n;l*=2){
        ttt[l]=ksm(G,(mod-1)/l,mod);
        invv[l]=ksm(ttt[l],mod-2,mod);
    }
}
inline bool check(int x){
    for(register int i=2;i*i<=m-1;++i)
    if(ksm(x,(m-1)/i,m)==1)return false;
    return true;
}
inline ll findroot(){
    if(m==2)return 1;
    for(register int i=2;i<=m;++i)
    if(check(i))return i;   
}
inline void NTT(ll *a,int flag){
    for(register int i=0;i<n;i++)if(rev[i]>i)swap(a[i],a[rev[i]]);
    for(l=2;l<=n;l*=2){
        int m=l/2;ll tmp=ttt[l],omg;
        if(flag==-1)tmp=invv[l];
        for(ll *p=a;p!=a+n;p+=l){
            omg=1;for(j=0;j<m;++j,omg=omg*tmp%mod){
                ll t=p[j+m]*omg%mod;
                p[j+m]=p[j]-t;p[j+m]=(p[j+m]+mod)%mod;
                p[j]=p[j]+t;p[j]%=mod;
            }
        }
    }
}
inline void ksm2(ll tt){
    for(;tt;tt>>=1){
        if(tt&1){
            for(register int i=0;i<n;++i)b[i]=dp[i],a[i]=g[i];
            NTT(a,1),NTT(b,1);for(register int i=0;i<n;++i)a[i]=a[i]*b[i]%mod;
            NTT(a,-1);for(register int i=0;i<m-1;++i)dp[i]=(a[i]*inv%mod+a[i+m-1]*inv%mod)%mod;
        }
        for(register int i=0;i<n;i++)b[i]=g[i],a[i]=g[i];
        NTT(a,1),NTT(b,1);for(register int i=0;i<n;++i)a[i]=a[i]*b[i]%mod;
        NTT(a,-1);for(register int i=0;i<m-1;++i)g[i]=(a[i]*inv%mod+a[i+m-1]*inv%mod)%mod;
    }
}
int main(){ll tt;
    mm=read(),m=read(),x=read(),slen=read();
    ll tmp=findroot();for(register int i=0;i<m-1;++i)mp[ksm(tmp,i,m)]=i;
    while(n<(m-1)*2)n<<=1,lim++;for(register int i=1;i<=slen;++i){
        tt=read();if(tt==0)continue;
        g[mp[tt]]++;
    }dp[0]=1;
    init();inv=ksm(n,mod-2,mod);ksm2(mm);printf("%lld
",dp[mp[x]]);
}