【SDOI 2015】序列统计

在遭受了巨大的精神打击后，我决定了，要做一个压行的好孩子！！！

博客如果有不对之处，请一定要指出啊啊啊啊啊谢谢谢谢谢大佬！！！

我们首先可以写出一个递推式：

$f[a+c][b\ast d取模mod]=f[a][b]\ast f[c][d]$

可以发现，这个格式并不是卷积。那我们怎么办呢？

我们知道，有$a^i\ast a^j=a^{i+j}$，那我们可以将式子转化为原根i次方。然后我们就可以用卷积来做了。

然后我再来说一下打法：solve函数类似快速幂，每一次乘法相当于添加一个数（就像dp外层的1至n的循环）。而g数组则表示相应的余数对应的方案数。

最后，因为运算时不带取模，我们需要再加上$a[i+m-1]$。为什么加m-1而不是加m？因为此处保证m为质数，也就是说m的欧拉函数是m-1。（如有不懂欧拉函数的童鞋们可以去康康 ）那么，我们可以保证的是，0到m-2构成的数，在m-1及之后相同的长度还有一毛一样的。而这取模m-1就是相同的余数。为什么只加一次呢？这要看我们变换的范围。

再最后，我们知道m-1并不一定是阶。这就可能有0到m-2之间有重复余数但不相加的情况。为什么可以直接输出呢？仔细看看我们的初始化，我们将输入的值为余数处理。实际上，这些不同的编号（i）代表的是同一个数。

awsl感觉自己在口胡。不知所言。

所以又到了激动人心的上代码环节。

#include<cstdio>
#include<iostream>
using namespace std;

const int mod = 1004535809, N = 50002;
int n, m, X, len, inv, lim = 1, rev[N], ori, L, pos = -1, a[N], b[N], g[N], f[N];
bool vis[N];

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

int cal(const int p) {//找原根
	if(p == 2)
		return 1;
	bool flag;
	for(int i = 2; ; ++ i) {
		flag = 0;
		for(int j = 2; j * j < p; ++ j)
			if(qkpow(i, (p - 1) / j, p) == 1) {
				flag = 1;
				break;
			}
		if(! flag)
			return i;
	}
}

void Swap(int &x, int &y) {
	x ^= y ^= x ^= y;
}

void NTT(int *t, const int op) {
	for(int i = 0; i < lim; ++ i) if(i < rev[i]) Swap(t[i], t[rev[i]]);
	for(int mid = 2; mid <= lim; mid <<= 1) {
		int per = qkpow(3, (mod - 1) / mid, mod);
		if(op == -1) per = qkpow(per, mod - 2, mod);
		for(int i = 0; i < lim; i += mid) {
			int w = 1;
			for(int k = 0; k < (mid >> 1); ++ k, w = 1ll * w * per % mod) {
				int X = t[i + k], Y = 1ll * w * t[i + k + (mid >> 1)] % mod;
				t[i + k] = (X + Y) % mod;
				t[i + k + (mid >> 1)] = (X - Y + mod) % mod; 
			}
		}
	}
}

void mul(int *g, int *f) {
	for(int i = 0; i < lim; ++ i) a[i] = g[i] % mod, b[i] = f[i] % mod;
	NTT(a, 1);
	NTT(b, 1);
	for(int i = 0; i < lim; ++ i) a[i] = 1ll * a[i] * b[i] % mod;
	NTT(a, -1);
	for(int i = 0; i < lim; ++ i) a[i] = 1ll * a[i] * inv % mod;//因为在模意义之下，不能直接除以lim
	for(int i = 0; i < m - 1; ++ i)	g[i] = (a[i] + a[i + m - 1]) % mod;
}

void solve(int y) {
	inv = qkpow(lim, mod - 2, mod);
	g[0] = 1;
	while(y) {
		if(y & 1) mul(g, f);
		mul(f, f);
		y >>= 1;
	}
}

void init(const int Lim) {
	while(lim < Lim) {
		lim <<= 1;
		++ L;
	}
	for(int i = 0; i < lim; ++ i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (L - 1));
}

int main() {
	ios::sync_with_stdio(false);
	cin.tie(NULL);
	int num;
	cin >> n >> m >> X >> len;
	ori = cal(m);
	for(int i = 1; i <= len; ++ i) cin >> num, vis[num] = 1;
	for(int i = 0, j = 1; i < m - 1; ++ i, j = 1ll * j * ori % m) {
		if(vis[j]) f[i] = 1;//有这种形式的g^i的数，其中g为原根
		if(X == j) pos = i;
	}
	init(m - 1 << 1);
	solve(n);
	if(pos != -1) cout << g[pos] % mod;
	else cout << "0";
	putchar('\n');
	return 0;
}

谢谢！