【LG3321】[SDOI2015]序列统计
题面
题解
前置芝士：原根
我们先看一下对于一个数(p)，它的原根(g)有什么性质(好像就是定义)：
(g^0\%p,g^1\%p,g^2\%p...g^{p-2}\%p) 恰好等于 ([0,p - 1])中所有数。
那么怎么求呢？
对(varphi(p))分解质因数，得到(varphi(p)=p_1^{a_1}p_2^{a_2}p_3^{a_3}...p_n^{a_n})
从(2sim p-1) 枚举 (g)，如果满足 (g) 对于 (forall p_i) ，有(g^{frac {varphi(p)}{p_i}} eq1;mod;p)
则该数是个原根,( ext{break})，否则( ext{continue})
关于此题：
有了上面的铺垫，我们想一想这题怎么做。
设(f[i][j])表示选了(i)个数，乘积(\%m)为(j)的方案数，
则有转移：
[f[2*i][c]=sum_{a*b\%m=c}f[i][a]*f[i][b]
]
这时候我们复杂度是(O(m^2log n))的，跑不过去，而转移次数已经无法优化了，想办法优化转移。
观察这个转移，如果它的判断条件为((a+b)\%m=c)，我们不就可以卷起来了吗？
想想什么能把乘法换成加法？对数！！！
但是因为是模意义下的对数，所以我们取个原根就行了。
( herefore) 设(C=log_gc\%m,A=log_ga\%m,B=log_gb\%m)。
则有转移：
[f[2*i][C]=sum_{(A+B)\%varphi (m)=c}f[i][A]*f[i][B]
]
那么就可以用(NTT)搞了，
注意最后要(f[i][j]+=f[i][j+varphi (m)])，因为你每次卷起来后(varphi (m))~((2varphi (m)-2))项都是要算贡献的。
注意：集合中(\%m=0)的数也要判一下。
代码
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <map>
using namespace std;
inline int gi() {
register int data = 0, w = 1;
register char ch = 0;
while (!isdigit(ch) && ch != '-') ch = getchar();
if (ch == '-') w = -1, ch = getchar();
while (isdigit(ch)) data = 10 * data + ch - '0', ch = getchar();
return w * data;
}
int fpow(int x, int y, int mod) {
int res = 1;
while (y) {
if (y & 1) res = 1ll * res * x % mod;
x = 1ll * x * x % mod;
y >>= 1;
}
return res;
}
const int Mod = 1004535809, G = 3, iG = fpow(G, Mod - 2, Mod);
int GetRoot(int x) {
int fact[10000], tot = 0;
int phi = x - 1;
for (int i = 2; i * i <= phi; i++) {
if (phi % i == 0) {
fact[++tot] = i;
while (phi % i == 0) phi /= i;
}
}
if (phi > 1) fact[++tot] = phi;
phi = x - 1;
for (int i = 2; i <= phi; i++) {
bool flg = 1;
for (int j = 1; j <= tot && flg; j++)
if (fpow(i, phi / fact[j], x) == 1) flg = 0;
if (flg) return i;
}
return -1;
}
const int MAX_M = 2.4e4 + 5;
int Limit, rev[MAX_M];
void NTT(int *p, int op) {
for (int i = 0; i < Limit; i++) if (i < rev[i]) swap(p[i], p[rev[i]]);
for (int i = 1; i < Limit; i <<= 1) {
int rot = fpow(op == 1 ? G : iG, (Mod - 1) / (i << 1), Mod);
for (int j = 0; j < Limit; j += (i << 1)) {
int w = 1;
for (int k = 0; k < i; k++, w = 1ll * w * rot % Mod) {
int x = p[j + k], y = 1ll * w * p[i + k + j] % Mod;
p[j + k] = (x + y) % Mod, p[i + j + k] = (x - y + Mod) % Mod;
}
}
}
if (op == -1) {
int inv = fpow(Limit, Mod - 2, Mod);
for (int i = 0; i < Limit; i++) p[i] = 1ll * p[i] * inv % Mod;
}
}
map<int, int> mp;
int N, M, X, S, F[MAX_M], H[MAX_M];
void mul(int *A, int *B, int *C) {
static int res[MAX_M], a[MAX_M], b[MAX_M];
for (int i = 0; i < Limit; i++) a[i] = A[i], b[i] = B[i];
NTT(a, 1), NTT(b, 1);
for (int i = 0; i < Limit; i++) a[i] = 1ll * a[i] * b[i] % Mod;
NTT(a, -1);
for (int i = 0; i < M - 1; i++) res[i] = (a[i] + a[i + M - 1]) % Mod;
for (int i = 0; i < M - 1; i++) C[i] = res[i];
}
int main () {
N = gi(), M = gi(), X = gi(), S = gi();
int g = GetRoot(M); for (int i = 0; i < M - 1; i++) mp[fpow(g, i, M)] = i;
for (int i = 1, x; i <= S; i++) {
x = gi() % M;
if (x) F[mp[x % M]]++;
}
H[mp[1]] = 1;
int p = 0;
for (Limit = 1; Limit <= 2 * M; Limit <<= 1, ++p) ;
for (int i = 0; i < Limit; i++) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (p - 1));
while (N) {
if (N & 1) mul(H, F, H);
mul(F, F, F);
N >>= 1;
}
printf("%d
", H[mp[X]]);
return 0;
}