给你一个模 \(m\) 意义下的数集，需要用这个数集生成一个数列，使得这个数列在的乘积为 \(x\)。
问方案数模 \(1004535809\)。 保证 \(m\) 是质数。
首先要知道 \(1004535809\) 是个 \(NTT\) 常用模数，原根为 \(3\)。
然后推一推 \(DP\) 式子，设 \(f[i][j]\) 表示填好了前 \(i\) 个数，乘积为 \(j\) 的方案数。
转移显然 \(f[i][j]=\sum\limits_{p\times k=j}f[i-1][p]\times sum[k]\) 。其中 \(sum[i]\) 表示数集中 \(i\) 是否出现。
然鹅这样做是 \(O(m^2n)\) 的
观察到如果把转移写成矩阵，那每一次的转移矩阵都是相等的，这启发我们进行快速幂，复杂度变为 \(O(m^2\log n)\)。
发现这个 \(p\times k=j\) 不好处理，如果把它变成 \(p+k=j\) 的形式就可以 \(NTT\) 优化了。
\(m\) 是质数，它的原根数量为 \(\phi(m-1)\) 个。
\(m\ge 3\)。
相当于保证 \(m\) 一定有原根。
根据原根的性质，\(g^0,g^1,g^2\dots g^{m-2}\) 两两互不相同。
那 \(DP\) 式子可以改写为 \[f[i][g^a]=\sum_{g^b\times g^c=g^a} f[i-1][g^b]\times sum[g^c]\]
\[f[i][g^a]=\sum\limits_{g^{b+c}=g^a} f[i-1][g^b]\times sum[g^c]\]
\[f[i][a]=\sum\limits_{b+c=a} f[i-1][b]\times sum[c]\]
快速幂的时候 \(NTT\) 优化即可。
哦对了，因为这个式子只有 \(g^0\sim g^{m-2}\) 总共 \(m-1\) 项，所以下标实际上都应该对 \(m-1\) 取模，因为 \(0\) 是没有对应的项的。
时间复杂度 \(O(m\log n\log 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\n",a[pos[X]]); return 0; }