P3321 [SDOI2015]序列统计

题意

首先考虑暴力DP怎么做:

(f_{i,j})表示选了(i)个数,乘积为(j)的方案数。
(f_{i,j}=sumlimits_{a*b\%m= j}f_{i-1,a}*f_{1,b})
这是(O(nm^2))的,我们考虑优化:

1.首先注意到(n)很大,于是考虑倍增地转移:
(f_{2*i,j}=sumlimits_{a*b\%m= j}f_{i,a}*f_{i,b})
于是就变为(O(m^2logn))的了。

2.发现下标很有意思,似乎是卷积的形式,但是这是乘积的形式,我们考虑将它化为加法的形式,即取(log)。由于是模意义下的,我们考虑指标。

我们先求出(m)的一个原根(g)。根据原根的性质,我们对于一个(xin [0,m-1]),只有唯一一个(k),满足(g^kequiv xpmod{p}),因此此我们能将上式变为:
(f_{2*i,j}=sumlimits_{g^a*g^b\%m= g^j}f_{i,a}*f_{i,b})
(f_{2*i,j}=sumlimits_{a+b\%(m-1)= j}f_{i,a}*f_{i,b})

于是就可以(NTT)加速了。

注意我们对下标取了模,但实际上(NTT)时是不会对下标取模的,也就是说(f_{i,j})的答案会存在第(i)位和第(i+m-1)位,我们要将它们加起来。

code:

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxm=32010;
const ll mod=1004535809;
const ll G=3;
const ll invG=334845270;
int n,m,K,S,lim=1,len;
int pos[maxm],pri[maxm];
ll f[maxm],g[maxm],a[maxm],b[maxm],c[maxm];
unordered_map<int,int>mp;
inline ll power(ll x,ll k,ll mod)
{
	ll res=1;
	while(k)
	{
		if(k&1)res=res*x%mod;
		x=x*x%mod;k>>=1;
	}
	return res;
}
inline int getroot(int mod)
{
	int tmp=mod-1,tot=0;
	for(int i=2;i*i<=tmp;i++)
	{
		if(tmp%i)continue;
		pri[++tot]=i;
		while(tmp%i==0)tmp/=i;
	}
	if(tmp>1)pri[++tot]=tmp;
	for(int i=2;i<=mod-1;i++)
	{
		bool flag=1;
		for(int j=1;j<=tot;j++)
			if(power(i,(mod-1)/pri[j],mod)==1){flag=0;break;}
		if(flag)return i;
	}
	return -1;
}
inline void NTT(ll* a,int op)
{
    for(int i=0;i<lim;i++)if(i<pos[i])swap(a[i],a[pos[i]]);
    for(int l=1;l<lim;l<<=1)
    {
        ll wn=power(op==1?G:invG,(mod-1)/(l<<1),mod);
        for(int i=0;i<lim;i+=l<<1)
        {
           	ll w=1;
            for(int j=0;j<l;j++,w=w*wn%mod)
            {
                int x=a[i+j],y=w*a[i+l+j]%mod;
                a[i+j]=(x+y)%mod;a[i+l+j]=(x-y+mod)%mod;
            }
        }
    }
    if(op==1)return;
    ll inv=power(lim,mod-2,mod);
    for(int i=0;i<lim;i++)a[i]=a[i]*inv%mod;
}
inline void mul(ll* f,ll* g,ll* h)
{
	for(int i=0;i<lim;i++)a[i]=f[i],b[i]=g[i];
	NTT(a,1);NTT(b,1);
	for(int i=0;i<lim;i++)a[i]=a[i]*b[i]%mod;
	NTT(a,-1);
	for(int i=0;i<m-1;i++)c[i]=(a[i]+a[i+m-1])%mod;
	for(int i=0;i<m-1;i++)h[i]=c[i];
}
int main()
{
	scanf("%d%d%d%d",&n,&m,&K,&S);
	int gi=getroot(m);
	for(int i=0;i<m-1;i++)mp[power(gi,i,m)]=i;
	for(int i=1;i<=S;i++)
	{
		int x;scanf("%d",&x);x%=m;
		if(x)f[mp[x]]++;
	}
	while(lim<=2*m)lim<<=1,len++;
	for(int i=0;i<lim;i++)pos[i]=(pos[i>>1]>>1)|((i&1)<<(len-1));
	g[mp[1]]=1;
	while(n)
	{
		if(n&1)mul(g,f,g);
		mul(f,f,f);n>>=1;
	}
	printf("%lld
",g[mp[K]]);
	return 0;
}
原文地址:https://www.cnblogs.com/nofind/p/12142957.html