洛咕 P3321 [SDOI2015]序列统计

显然dp就是设\(f[i][j]\)表示dp了i轮,对m取膜是j的方案数

\(f[i][xy\mod m]=f[i-1][x]\times f[i-1][y]\)

这是\(O(nm^2)\)

像我这样的蒟蒻都能想到用类似快速幂一样的东西来转移是吧,那么就\(O(log_2 nm^2)\)

非常难受,还是过不去

如果可以优化一下dp转移就好了,比如把乘改成加,就能用NTT了

然后就要用到一个叫做原根的东西,学NTT的时候只是记了一下不知道这货有啥用

质数\(m\)原根\(g\)的性质:对\(m-1\)进行唯一分解,\(m-1=p_1^{a_1}p_2^{a_2}\cdots p_k^{a_k}\)\(g\)对任何一个\(i\)都满足\(g^{\frac{m-1}{p_i}}\neq 1\)

然后,这个题要用到的就是,\(g^0,g^1,\cdots,g^{m-2}\)互不相同

那么,把原来f[1][i]变成f[1][g^{i-1}],就可以一一对应了,又把乘变成了加

#include<bits/stdc++.h>
#define il inline
#define vd void
#define int ll
#define mod 1004535809
typedef long long ll;
il int gi(){
    int x=0,f=1;
    char ch=getchar();
    while(!isdigit(ch)){
        if(ch=='-')f=-1;
        ch=getchar();
    }
    while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
    return x*f;
}
il int Pow(int x,int y,int mo){
    int ret=1;
    while(y){
        if(y&1)ret=ret*x%mo;
        x=x*x%mo;y>>=1;
    }
    return ret;
}
const int G=3,iG=Pow(G,mod-2,mod);
int n,m,qx,k,F[8000];
struct yyb{
    ll f[8000];
    yyb(){memset(f,0,sizeof f);}
};
ll N,lg,rev[1<<14],P[1<<14],iP[1<<14],invN;
il vd ntt(ll*A,int n,int t){
    for(int i=0;i<n;++i)if(rev[i]>i)std::swap(A[i],A[rev[i]]);
    for(int o=1;o<n;o<<=1){
        int W=t?P[o]:iP[o];
        for(int*p=A;p!=A+n;p+=o<<1)
            for(int i=0,w=1;i<o;++i,w=w*W%mod){
                int t=w*p[i+o]%mod;
                p[i+o]=(p[i]-t+mod)%mod;p[i]=(p[i]+t)%mod;
            }
    }
    if(!t)for(int i=0;i<n;++i)A[i]=A[i]*invN%mod;
}
il yyb mul(const yyb&a,const yyb&b){
    static int A[1<<14],B[1<<14];
    memset(A,0,N*8),memset(B,0,N*8);
    memcpy(A,a.f,sizeof a.f),memcpy(B,b.f,sizeof b.f);
    ntt(A,N,1),ntt(B,N,1);
    for(int i=0;i<N;++i)A[i]=A[i]*B[i]%mod;
    ntt(A,N,0);
    yyb ret;
    for(int i=0;i<=m*2;++i)(ret.f[i%(m-1)]+=A[i])%=mod;
    return ret;
}
signed main(){
    n=gi(),m=gi(),qx=gi(),k=gi();	
    int gM;
    {//get gM
        int p[10],tot=0;
        for(int x=m-1,i=2;i<=x;++i)
            if(x%i==0){
                p[++tot]=i;
                while(x%i==0)x/=i;
            }
        for(int i=2;i<=m;++i){
            for(int j=1;j<=tot;++j)if(Pow(i,(m-1)/p[j],m)==1)goto gg;
            gM=i;break;gg:;
        }
        for(int i=0,j=1;i<m-1;++i,j=j*gM%m)F[j]=i;
    }
    N=1,lg=0;while(N<m<<1)N<<=1,++lg;
    for(int i=0;i<N;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<lg-1);
    for(int i=1;i<N;i<<=1)P[i]=Pow(G,(mod-1)/(i<<1),mod),iP[i]=Pow(iG,(mod-1)/(i<<1),mod);
    invN=Pow(N,mod-2,mod);
    yyb X,ans;
    int t;
    while(k--){
        t=gi();
        if(t)X.f[F[t]]=1;
    }
    ans.f[F[1]]=1;
    while(n){
        if(n&1)ans=mul(ans,X);
        X=mul(X,X);n>>=1;
    }
    printf("%lld\n",ans.f[F[qx]]);
    return 0;
}
posted @ 2018-12-03 20:42  菜狗xzz  阅读(171)  评论(0编辑  收藏  举报