[SDOI 2015] 序列统计

目录

• $ ext{Description}$
• $ ext{Solution}$
• $ ext{Code}$

( ext{Description})

传送门

( ext{Solution})

首先能想到朴素 (mathtt{DP})：枚举数列长度，长度 (-1) 数列乘积，新加入的数，这是 (mathcal O(n imes m^2)) 的。

有一个优化类似快速幂，就是把数与数列合并变成数列和数列合并，可以做到 (mathcal O(m^2log n))。

观察题目性质，所有数值都是属于 ([0,m)) 的，我们 自然 地想到了用 (m) 的原根的幂来表示数值（原根的 ([0,m)) 次幂正好包含 ([0,m))（注意 (m) 是质数））。

好处就是将数值乘积转化成指数相加，令 (g) 为当前状态，(f) 为转移方式，就有一个卷积：

[g(x)=sum_{i=0}^{m-1}g(i)f(x-i) ]

这个用 (mathtt{NTT}) 处理即可。

详见代码。

( ext{Code})

#include <cstdio>

#define rep(i,_l,_r) for(register signed i=(_l),_end=(_r);i<=_end;++i)
#define fep(i,_l,_r) for(register signed i=(_l),_end=(_r);i>=_end;--i)
#define erep(i,u) for(signed i=head[u],v=to[i];i;i=nxt[i],v=to[i])
#define efep(i,u) for(signed i=Head[u],v=to[i];i;i=nxt[i],v=to[i])
#define print(x,y) write(x),putchar(y)

template <class T> inline T read(const T sample) {
    T x=0; int f=1; char s;
    while((s=getchar())>'9'||s<'0') if(s=='-') f=-1;
    while(s>='0'&&s<='9') x=(x<<1)+(x<<3)+(s^48),s=getchar();
    return x*f;
}
template <class T> inline void write(const T x) {
    if(x<0) return (void) (putchar('-'),write(-x));
    if(x>9) write(x/10);
    putchar(x%10^48);
}
template <class T> inline T Max(const T x,const T y) {if(x>y) return x; return y;}
template <class T> inline T Min(const T x,const T y) {if(x<y) return x; return y;}
template <class T> inline T fab(const T x) {return x>0?x:-x;}
template <class T> inline T gcd(const T x,const T y) {return y?gcd(y,x%y):x;}
template <class T> inline T lcm(const T x,const T y) {return x/gcd(x,y)*y;}
template <class T> inline T Swap(T &x,T &y) {x^=y^=x^=y;}

#include <iostream>
using namespace std;

const int mod=1004535809,maxn=8005,ig=334845270;

int len,m,goal,S,rt,lim,bit,rev[maxn<<2],f[maxn<<2],g[maxn<<2],co[maxn],a[maxn],ilim;

int qkpow(int x,int y,int Mod) {
	int r=1;
	while(y) {
		if(y&1) r=1ll*r*x%Mod;
		x=1ll*x*x%Mod; y>>=1;
	}
	return r;
}

void GetRoot() { // 求原根
	int x; bool flag;
	for(rt=2;rt<m;++rt) {
		x=m-1; flag=0;
		for(int i=2;i*i<=x;++i) {
			if(x%i) continue;
			while(!(x%i)) x/=i;
			if(qkpow(rt,(m-1)/i,m)==1) {flag=1; break;}
		}
		if(x>1&&qkpow(rt,(m-1)/x,m)==1) continue;
		if(!flag) return;
	}
}

void NTT(int *t,bool op) {
    int wn,w,tmp;
    rep(i,0,lim-1) if(i<rev[i]) swap(t[i],t[rev[i]]);
    for(int mid=1;mid<lim;mid<<=1) {
        wn=qkpow(op?3:ig,(mod-1)/(mid<<1),mod);
        for(int i=0;i<lim;i+=(mid<<1)) {
            w=1;
            for(int j=0;j<mid;++j,w=1ll*w*wn%mod) {
                tmp=1ll*t[i+j+mid]*w%mod;
                t[i+j+mid]=(t[i+j]-tmp+mod)%mod,t[i+j]=(t[i+j]+tmp)%mod;
            }
        }
    }
    if(!op) rep(i,0,lim-1) t[i]=1ll*t[i]*ilim%mod;
}

void Goose() {
    g[0]=1;
    while(len) {
        if(len&1) {
            NTT(f,1),NTT(g,1);
            rep(i,0,lim-1) g[i]=1ll*g[i]*f[i]%mod;
            NTT(f,0),NTT(g,0);
            rep(i,m-1,lim-1) g[i%(m-1)]=(g[i%(m-1)]+g[i])%mod,g[i]=0;
            // 由于 f,g 的范围是 [0,m)，做卷积时会超出范围，将超出范围的加回去再归零
        }
        len>>=1;
        NTT(f,1);
        rep(i,0,lim-1) f[i]=1ll*f[i]*f[i]%mod;
        NTT(f,0);
        rep(i,m-1,lim-1) f[i%(m-1)]=(f[i%(m-1)]+f[i])%mod,f[i]=0;
    }
}

int main() {
	int x;
	len=read(9),m=read(9),goal=read(9),S=read(9);
	GetRoot();
	rep(i,0,m-2) co[qkpow(rt,i,m)]=i;
	rep(i,1,S) {
		x=read(9);
		if(!x) continue;
		a[i]=co[x],++f[a[i]];
        // 初始化转移方式
	}
	lim=1;
	while(lim<=2*m-2) lim<<=1,++bit;
    // 由于 f,g 的范围是 [0,m)，做卷积时会有 m-1+m-1 项
    ilim=qkpow(lim,mod-2,mod);
	rep(i,0,lim-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<bit-1);
    Goose();
    print(g[co[goal]],'
');
	return 0;
}