[SDOI2015]序列统计 题解
题目地址:洛谷:【P3321】[SDOI2015]序列统计 – 洛谷、BZOJ:Problem 3992. — [SDOI2015]序列统计
题目描述
小C有一个集合S,里面的元素都是小于M的非负整数。他用程序编写了一个数列生成器,可以生成一个长度为N的数列,数列中的每个数都属于集合S。小C用这个生成器生成了许多这样的数列。但是小C有一个问题需要你的帮助:给定整数x,求所有可以生成出的,且满足数列中所有数的乘积mod M的值等于x的不同的数列的有多少个。小C认为,两个数列{Ai}和{Bi}不同,当且仅当至少存在一个整数i,满足Ai≠Bi。另外,小C认为这个问题的答案可能很大,因此他只需要你帮助他求出答案mod 1004535809的值就可以了。
输入输出格式
输入格式:
一行,四个整数,N、M、x、|S|,其中|S|为集合S中元素个数。第二行,|S|个整数,表示集合S中的所有元素。
输出格式:
一行,一个整数,表示你求出的种类数mod 1004535809的值。
输入输出样例
输入样例#1:
4 3 1 2 1 2
输出样例#1:
8
题解
参考资料:bzoj 3992: [SDOI2015]序列统计 (NTT+快速幂+DP) – CSDN博客
本题需要的数学姿势:数学笔记:数论(欧拉函数、阶、原根) | KSkun’s Blog
首先,计数自然是一个动态规划,我们考虑设计状态dp[i][j]为选了i个数且结果mod M的值是j,转移如下
dp[i][j] = \sum_{k \in S} dp[i-1][j'] \ (j'k \bmod M = j)
直观上感觉直接做是一个O(n|S|)的,好像只能得10分。
接下来我们要考虑一个问题,模数很特殊,并不是常见的1000000007,而是1004535809,它的原根为3,是常用的NTT模数,难道出题人在指引我们使用NTT?
想用NTT做,首先得搞出卷积,卷积下标是加减的事情没办法像上面那样乘吧。那么考虑利用原根把乘法转换为加法,dp[i][j]表示选到第i个数结果mod M的值为g^j,改变以后的转移方程就成这样了
dp[i][j] = \sum_{g^k \bmod M \in S} dp[i-1][j'] \ ((j' + k) \bmod (M - 1) = j)
或者这样?用f[i]表示g^i是否在集合中。
dp[i][j] = \sum_{k=1}^{M-1} dp[i-1][j-k] \cdot f[k]
咦,后面那部分现在可以卷积了,而且你会发现每一层的转移中f数组都是不变的,也就是说这是个线性变换,参照矩阵快速幂的操作,我们可以做多项式的快速幂,即dp[n] = dp[0] \times f^n。答案便是 dp[n][i] \ (g^i=x) 。
于是我们的思路就很明确了:求原根、算出f数组、多项式快速幂。
没完呢!这里有一个坑点,集合内的元素可以有0,0显然对答案没贡献,所以你要跳过它处理。
代码
// Code by KSkun, 2018/4
#include <cstdio>
#include <cmath>
#include <cstring>
#include <algorithm>
typedef long long LL;
inline char fgc() {
static char buf[100000], *p1 = buf, *p2 = buf;
return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2) ? EOF
: *p1++;
}
inline LL readint() {
register LL res = 0, neg = 1;
char c = fgc();
while(c < '0' || c > '9') {
if(c == '-') neg = -1;
c = fgc();
}
while(c >= '0' && c <= '9') {
res = res * 10 + c - '0';
c = fgc();
}
return res * neg;
}
inline LL fpow(LL n, LL k, LL p) {
n %= p;
LL res = 1;
while(k) {
if(k & 1) res = res * n % p;
n = n * n % p;
k >>= 1;
}
return res;
}
const int MAXN = 1 << 15, G = 3, MO = 1004535809;
int n, m, x, s, t;
LL f[MAXN], g[MAXN];
bool vis[MAXN];
LL rev[MAXN], dft1[MAXN], dft2[MAXN], dftr[MAXN];
inline void calrev(int n, int len) {
for(int i = 0; i < n; i++) {
rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (len - 1));
}
}
inline void dft(int n, LL *arr, int f) {
for(int i = 0; i < n; i++) {
if(i < rev[i]) std::swap(arr[i], arr[rev[i]]);
}
for(int i = 1; i < n; i <<= 1) {
LL gn = fpow(G, (MO - 1) / (i << 1), MO);
if(f == -1) gn = fpow(gn, MO - 2, MO);
for(int j = 0; j < n; j += i << 1) {
LL g = 1;
for(int k = 0; k < i; k++) {
LL x = arr[j + k], y = g * arr[j + k + i] % MO;
arr[j + k] = (x + y) % MO;
arr[j + k + i] = ((x - y) % MO + MO) % MO;
g = g * gn % MO;
}
}
}
if(f == -1) {
LL invn = fpow(n, MO - 2, MO);
for(int i = 0; i < n; i++) {
arr[i] = arr[i] * invn % MO;
}
}
}
inline void ntt(LL *a, LL *b) {
memset(dft1, 0, sizeof(dft1));
memset(dft2, 0, sizeof(dft2));
memset(dftr, 0, sizeof(dftr));
int n, len = 0;
for(n = 1; n <= (m - 1) << 1; n <<= 1) len++;
for(int i = 0; i < n; i++) {
dft1[i] = a[i];
dft2[i] = b[i];
}
calrev(n, len);
dft(n, dft1, 1); dft(n, dft2, 1);
for(int i = 0; i < n; i++) {
dftr[i] = dft1[i] * dft2[i] % MO;
}
dft(n, dftr, -1);
for(int i = 0; i < m - 1; i++) {
a[i] = (dftr[i] + dftr[i + m - 1]) % MO;
}
}
inline LL calg(LL n) {
int fact[8005], tot = 0, x = n - 1;
for(int i = 2; i * i <= x; i++) {
if(x % i == 0) {
fact[tot++] = (n - 1) / i;
while(x % i == 0) x /= i;
}
}
if(x > 1) fact[tot++] = (n - 1) / x;
bool success;
for(int i = 2;; i++) {
success = true;
for(int j = 0; j < tot; j++) {
if(fpow(i, fact[j], n) == 1) {
success = false;
break;
}
}
if(success) return i;
}
}
int main() {
n = readint(); m = readint(); x = readint(); s = readint();
for(int i = 1; i <= s; i++) {
t = readint();
vis[t] = 1;
}
int xpow, G = calg(m), gn = 1;
for(int i = 0; i < m - 1; i++) {
if(vis[gn]) f[i] = 1;
if(gn == x) xpow = i;
gn = gn * G % m;
}
g[0] = 1;
for(int i = n; i; i >>= 1) {
if(i & 1) ntt(g, f);
ntt(f, f);
}
printf("%lld", g[xpow]);
return 0;
}