BZOJ3992[SDOI2015]序列统计

题目链接

洛谷

BZOJ

解析

头一回知道原根还可以这么考……

不难想到递推的做法\(dp[i][j]\)表示长度为\(i\)，乘积为\(j\)的答案，那么\(dp[i][j \cdot a[i] \ mod \ M] += dp[i - 1][j]\)，

首先我们发现\(0\)可以直接丢掉，因为包含\(0\)的序列对答案不产生任何贡献，然后题目说\(M\)是个质数，如果设\(G\)是\(M\)的原根，那么\([1, M - 1]\)和\(G^1 \ mod \ M, G^2 \ mod \ M, ..., G^{M - 1} \ mod \ M\)一一对应，那么可以用\(G\)头上的指数来代替\([1, M - 1]\)中的一个数，设\(dp[i][j]\)表示长度为\(i\)，乘积为\(G^j \ mod \ M\)的答案，则\(dp[i][(j + k)] += dp[i - 1][j]\)

用\(f[j]\)表示\(dp[i][j]\)，\(g[j]\)表示\(dp[i - 1][j]\)，\(h[j]\)表示给定的集合中是否存在\(G^j \ mod \ M\)这个数，那么
\[ f[j] = \sum_{k = 1}^{M - 1} g[k] \cdot h[j - k] \]
这是一个卷积

发现初始只有\(f[0] = 1\)，所以只需对\(h\)快速幂卷积即可

当然上面都是不考虑指数对\(M - 1\)取模的时候，考虑取模就每次卷积后把\(f[i](i > M - 1)\)加到\(f[i \ mod \ (M - 1)]\)的位置就可以了

代码

不知道为什么好像我的\(FFT\)和\(NTT\)都自带大常数……耗时比第一页最慢的都多了一倍加\(1s\)……

#include <cstring>
#include <iostream>
#include <cstdio>
#include <vector>
#define MAXN 8010
 
typedef long long LL;
const LL mod = 1004535809ll;
 
int qpower(int, int, int);
void pre_prime();
void divide(int, std::vector<int> &);
int get_g(int);
void pre_rev(int);
void NTT(int *, int, int);
void mul(int *, int *, int);
 
int N, M, X, SZ, n, g, pre[MAXN], f[MAXN << 2], ans[MAXN << 2], rev[MAXN << 2];
std::vector<int> prime, dvd;
 
inline void inc(int &x, int y) { x += y; if (x >= mod) x -= mod; }
inline void dec(int &x, int y) { x -= y; if (x < 0) x += mod; }
inline int add(int x, int y) { x += y; return x >= mod ? x - mod : x; }
inline int sub(int x, int y) { x -= y; return x < 0 ? x + mod : x; }
 
int main() {
    pre_prime();
    scanf("%d%d%d%d", &N, &M, &X, &SZ);
    g = get_g(M);
    for (int i = 1, j = g; i < M; ++i, j = (LL)j * g % M) pre[j] = i;
    for (int i = 0; i < SZ; ++i) {
        int t; scanf("%d", &t);
        if(t) ++f[pre[t]];
    }
    //debug
    //printf("%d\n", g);
    //for (int i = 0; i < M; ++i) printf("%d ", pre[i]);
    //puts("");
    ans[0] = 1;
    while ((1 << n) < (M << 1)) ++n;
    pre_rev(n);
    while (N) {
        if (N & 1) mul(ans, f, n);
        mul(f, f, n);
        N >>= 1;
    }
    printf("%d\n", ans[pre[X]]);
 
    return 0;
}
void pre_prime() {
    static bool isn_prime[MAXN];
    for (int i = 2; i < MAXN; ++i) {
        if (!isn_prime[i]) prime.push_back(i);
        for (int j = 0; j < prime.size(), i * prime[j] < MAXN; ++j) {
            isn_prime[i * prime[j]] = 0;
            if (i % prime[j] == 0) break;
        }
    }
}
int get_g(int x) {
    divide(x - 1, dvd);
    for (int i = 2; i < x; ++i)
        if (qpower(i, x - 1, x) == 1) {
            bool flag = 1;
            for (int j = 0; j < dvd.size(); ++j)
                if (qpower(i, (x - 1) / dvd[j], x) == 1) { flag = 0; break; }
            if (flag) return i;
        }
}
void divide(int x, std::vector<int> &res) {
    for (int i = 0; i < prime.size(); ++i)
        if (x % prime[i] == 0) {
            res.push_back(prime[i]);
            while (x % prime[i] == 0) x /= prime[i];
        }
}
int qpower(int x, int y, int p) {
    int res = 1;
    while (y) {
        if (y & 1) res = (LL)res * x % p;
        x = (LL)x * x % p; y >>= 1;
    }
    return res;
}
void NTT(int *arr, int sz, int tp) {
    for (int i = 0; i < (1 << sz); ++i)
        if (rev[i] > i) std::swap(arr[i], arr[rev[i]]);
    for (int len = 2, half = 1; len <= (1 << sz); len <<= 1, half <<= 1) {
        int wn = qpower(3, (mod - 1) / len, mod);
        if (tp == -1) wn = qpower(wn, mod - 2, mod);
        for (int i = 0; i < (1 << sz); i += len)
            for (int j = 0, w = 1; j < half; ++j, w = (LL)w * wn % mod) {
                int x = arr[i + j], y = (LL)arr[i + j + half] * w % mod;
                inc(arr[i + j], y); dec(arr[i + j + half] = x, y);
            }
    }
    if (tp == -1) {
        int inv = qpower(1 << sz, mod - 2, mod);
        for (int i = 0; i < (1 << sz); ++i) arr[i] = (LL)arr[i] * inv % mod;
    }
}
void mul(int *a, int *b, int sz) {
    static int tmp[MAXN << 2];
    for (int i = 0; i < (1 << sz); ++i) tmp[i] = b[i];
    NTT(a, sz, 1); NTT(tmp, sz, 1);
    for (int i = 0; i < (1 << sz); ++i) a[i] = (LL)a[i] * tmp[i] % mod;
    NTT(a, sz, -1);
    for (int i = 1; i < M; ++i) inc(a[i], a[i + M - 1]), a[i + M - 1] = 0; 
}
void pre_rev(int sz) {
    for (int i = 0; i < (1 << sz); ++i)
        rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << sz - 1));
}
//Rhein_E

转载于:https://www.cnblogs.com/Rhein-E/p/10528412.html