【JSOI 2012】分零食（倍增 + FFT）

95 阅读 0 评论 63 点赞

我是靠谱客的博主粗心自行车，这篇文章主要介绍【JSOI 2012】分零食（倍增 + FFT），现在分享给大家，希望可以做个参考。

题目链接：【JSOI 2012】分零食

题目大意：将 $m$ 颗糖依次分给 $n$ 个小朋友，每个小朋友的快乐度为 $f(x=得到的糖果数)=ox^2+sx+u$ 。特别地， $f(0)=1$ 。糖果分完后剩下的小朋友将得不到糖果。求每种情况下小朋友快乐度乘积之和 $mod 998244353$ 的结果。

令 $g[i][j]$ 表示前 $i$ 个小朋友得到 $j$ 颗糖的答案。则： $g[i][j]=sum_{k=1}^{j}{f(j-k)g[i-1][k]}$ 。发现这是个卷积的形式。想到 $FFT$ ，发现 $g[n]=f^n$ 。

令 $A[n]=sum_{j=1}^{n}{g[i]}$ 。我们使用倍增的方法求出 $A[n]$ 。

$A[n+1]=sum_{i=1}^{n+1}{g[i]}$
$A[n+1]=A[n]+g[n+1]$

$A[2n]=sum_{i=1}^{2n}{g[i]}$
$A[2n]=A[n]+sum_{i=n+1}^{2n}{g[i]}$
$A[2n]=A[n]+f^nsum_{i=1}^n{f^i}$
$A[2n]=A[n]+g[n]*A[n]$

然后用 $FFT$ 做卷积即可。

复制代码

#include <cmath>
#include <cstdio>
#include <complex>
#include <iostream>
using namespace std;
typedef complex<double> cd;
const double pi = acos(-1.);
const int maxn = 1 << 15 | 5;
cd f[maxn], g[maxn], dp[maxn], tmp[maxn];
int n, m, o, s, u, mod, lim, bit, r[maxn];
void fft(cd *a, int dft) {
    for (int i = 0; i < lim; i++) {
        if (i < r[i]) {
            swap(a[i], a[r[i]]);
        }
    }
    for (int k = 1; k < lim; k <<= 1) {
        cd wn0(cos(pi / k), dft * sin(pi / k));
        for (int i = 0; i < lim; i += k << 1) {
            cd wnk(1, 0);
            for (int j = i; j < i + k; j++, wnk *= wn0) {
                cd x = a[j], y = wnk * a[j + k];
                a[j] = x + y, a[j + k] = x - y;
            }
        }
    }
    if (dft == -1) {
        for (int i = 0; i < lim; i++) {
            a[i] /= lim;
            a[i] = int(a[i].real() + 0.5) % mod;
        }
    }
}
void solve(cd *dp, cd *g, int n) {
    if (!n) {
        g[0] = 1;
        return;
    }
    if (n & 1) {
        solve(dp, g, n - 1);
        fft(g, 1);
        for (int i = 0; i < lim; i++) {
            g[i] *= f[i];
        }
        fft(g, -1);
        for (int i = 0; i <= m; i++) {
            dp[i] += g[i];
            dp[i] = int(dp[i].real() + 0.5) % mod; 
        }
        for (int i = m + 1; i < lim; i++) {
            g[i] = 0;
        }
    } else {
        solve(dp, g, n >> 1);
        for (int i = 0; i < lim; i++) {
            tmp[i] = dp[i];
        }
        fft(tmp, 1);
        fft(g, 1);
        for (int i = 0; i < lim; i++) {
            tmp[i] *= g[i];
        }
        fft(tmp, -1);
        for (int i = 0; i < lim; i++) {
            g[i] *= g[i];
        }
        fft(g, -1);
        for (int i = 0; i <= m; i++) {
            dp[i] += tmp[i];
            dp[i] = int(dp[i].real() + 0.5) % mod; 
        }
        for (int i = m + 1; i < lim; i++) {
            g[i] = 0;
        }
    }
}
int main() {
    scanf("%d %d", &m, &mod);
    scanf("%d", &n);
    scanf("%d", &o);
    scanf("%d", &s);
    scanf("%d", &u);
    for (lim = 1; lim <= 2 * m; lim <<= 1)  bit++;
    for (int i = 0; i < lim; i++) {
        r[i] = (r[i >> 1] >> 1) | ((i & 1) << (bit - 1));
    }
    for (int i = 1; i <= m; i++) {
        f[i] = (o * i * i + s * i + u) % mod;
    }
    fft(f, 1);
    solve(dp, g, n);
    printf("%dn", int(dp[m].real() + 0.5));
    return 0;
}

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
#include <cmath>
#include <cstdio>
#include <complex>
#include <iostream>
using namespace std;
typedef complex<double> cd;
const double pi = acos(-1.);
const int maxn = 1 << 15 | 5;
cd f[maxn], g[maxn], dp[maxn], tmp[maxn];
int n, m, o, s, u, mod, lim, bit, r[maxn];
void fft(cd *a, int dft) {
    for (int i = 0; i < lim; i++) {
        if (i < r[i]) {
            swap(a[i], a[r[i]]);
        }
    }
    for (int k = 1; k < lim; k <<= 1) {
        cd wn0(cos(pi / k), dft * sin(pi / k));
        for (int i = 0; i < lim; i += k << 1) {
            cd wnk(1, 0);
            for (int j = i; j < i + k; j++, wnk *= wn0) {
                cd x = a[j], y = wnk * a[j + k];
                a[j] = x + y, a[j + k] = x - y;
            }
        }
    }
    if (dft == -1) {
        for (int i = 0; i < lim; i++) {
            a[i] /= lim;
            a[i] = int(a[i].real() + 0.5) % mod;
        }
    }
}
void solve(cd *dp, cd *g, int n) {
    if (!n) {
        g[0] = 1;
        return;
    }
    if (n & 1) {
        solve(dp, g, n - 1);
        fft(g, 1);
        for (int i = 0; i < lim; i++) {
            g[i] *= f[i];
        }
        fft(g, -1);
        for (int i = 0; i <= m; i++) {
            dp[i] += g[i];
            dp[i] = int(dp[i].real() + 0.5) % mod; 
        }
        for (int i = m + 1; i < lim; i++) {
            g[i] = 0;
        }
    } else {
        solve(dp, g, n >> 1);
        for (int i = 0; i < lim; i++) {
            tmp[i] = dp[i];
        }
        fft(tmp, 1);
        fft(g, 1);
        for (int i = 0; i < lim; i++) {
            tmp[i] *= g[i];
        }
        fft(tmp, -1);
        for (int i = 0; i < lim; i++) {
            g[i] *= g[i];
        }
        fft(g, -1);
        for (int i = 0; i <= m; i++) {
            dp[i] += tmp[i];
            dp[i] = int(dp[i].real() + 0.5) % mod; 
        }
        for (int i = m + 1; i < lim; i++) {
            g[i] = 0;
        }
    }
}
int main() {
    scanf("%d %d", &m, &mod);
    scanf("%d", &n);
    scanf("%d", &o);
    scanf("%d", &s);
    scanf("%d", &u);
    for (lim = 1; lim <= 2 * m; lim <<= 1)  bit++;
    for (int i = 0; i < lim; i++) {
        r[i] = (r[i >> 1] >> 1) | ((i & 1) << (bit - 1));
    }
    for (int i = 1; i <= m; i++) {
        f[i] = (o * i * i + s * i + u) % mod;
    }
    fft(f, 1);
    solve(dp, g, n);
    printf("%dn", int(dp[m].real() + 0.5));
    return 0;
}