模板:Chirp-Z 变换 / ARIS0_0 - 12

· · 题解

Chirp-Z 变换

给长度为 n 的多项式 f(x) 和正整数 p, q, m,在 O(n \log n) 的时间内求 f(p), f(p \cdot q), f(p \cdot q^2), \dots, f(p \cdot q^{m - 1})

有重要等式:xy = \dbinom{x + y}{2} - \dbinom{x}{2} - \dbinom{y}{2}

考虑代入点值 x_i = p \cdot q^i。则有

\begin{align*} y_i &= \sum_{j = 0}^{n - 1} f_jx_i^j \\ &= \sum_{j = 0}^{n - 1} f_j \cdot (p \cdot q^i)^j \\ &= \sum_{j = 0}^{n - 1} f_j \cdot p^j \cdot q^{ij} \\ &= \sum_{j = 0}^{n - 1} f_j \cdot p^j \cdot q^{\binom{i + j}{2} - \binom{i}{2} - \binom{j}{2}} \\ &= q^{-\binom{i}{2}}\sum_{j = 0}^{n - 1} \left(f_j \cdot p^j \cdot q^{-\binom{j}{2}}\right) \cdot q^{\binom{i + j}{2}} \end{align*}

注意到后面那个求和符号是个差卷积的形式,所以直接做一个 O(n \log n) 的卷积就行了。

#include <cstdio>
#include <vector>
#include <algorithm>
#define bit_ceil(x) ((x) <= 1 ? 1 : (1 << (32 - __builtin_clz((x) - 1))))
typedef std::vector<int> poly;
namespace FastIn
{
    const int S = 4e6 + 7;
    char *p1, *p2, buf[S];
    #define gc() (p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, S, stdin), p1 == p2) ? EOF : *p1++)
    inline void read(int &x)
    {
        x = 0;
        char c = gc();
        while (c < '0' || c > '9')
            c = gc();
        for (; c >= '0' && c <= '9'; c = gc())
            x = (x << 3) + (x << 1) + (c & 15);
    }
    #undef gc
}
using FastIn::read;
const int MAXN = 1e6 + 7, mod = 998244353;
inline int M(const int &x) { return x >= mod ? x - mod : x; }
inline int qpow(int a, int b)
{
    int res = 1;
    for (; b; a = (long long)(a) * a % mod, b >>= 1)
        (b & 1) && (res = (long long)(res) * a % mod);
    return res;
}
inline void NTT(poly &a, const int &rvs)
{
    static std::vector<int> rev;
    int n = a.size(); (int)(rev.size()) < n && (rev.resize(n), true);
    for (int i = 0; i < n; ++i)
        if (i < (rev[i] = (rev[i >> 1] >> 1) | ((i & 1) * (n >> 1))))
            std::swap(a[i], a[rev[i]]);
    for (int len = 1, W; len < n; len <<= 1)
    {
        W = qpow(rvs == 1 ? 3 : 332748118, (mod - 1) / (len << 1));
        for (int i = 0; i < n; i += (len << 1))
            for (int j = i, f0, f1, w = 1; j < i + len; ++j, w = (long long)(w) * W % mod)
                f0 = a[j], f1 = (long long)(w) * a[j + len] % mod, a[j] = M(f0 + f1), a[j + len] = M(f0 + mod - f1);
    }
    if (rvs == -1)
        for (int i = 0, invn = qpow(n, mod - 2); i < n; ++i)
            a[i] = (long long)(a[i]) * invn % mod;
}
inline poly operator*(const poly &A, const poly &B)
{
    int n = A.size() + B.size() - 1, len = bit_ceil(n);
    poly a = A, b = B;
    a.resize(len), b.resize(len), NTT(a, 1), NTT(b, 1);
    for (int i = 0; i < len; ++i)
        a[i] = (long long)(a[i]) * b[i] % mod;
    return NTT(a, -1), a.resize(n), a;
}
inline poly operator*=(poly &a, const poly &b) { return a = a * b; }
int n, m, p, q, iq, qi[MAXN << 1], qi2[MAXN << 1], iqi[MAXN], iqi2[MAXN];
poly a, f, g;
int main()
{
    read(n), p = 1, read(q), read(m), a.resize(n), iq = qpow(q, mod - 2);
    for (int i = 0; i < n; ++i)
        read(a[i]);
    f.resize(n), g.resize(n + m - 1);
    qi[0] = qi2[0] = iqi[0] = iqi2[0] = 1; // 为了卡常,记 qi[i] = q^i, qi2[i] = q^{\binom{i}{2}}, iqi[i] = q^{-i}, iqi2[i] = q^{-\binom{i}{2}},这些都可以线性递推求出
    for (int i = 1; i < n + m - 1; ++i)
        qi[i] = (long long)(qi[i - 1]) * q % mod, qi2[i] = (long long)(qi2[i - 1]) * qi[i - 1] % mod;
    for (int i = 1; i < std::max(n, m); ++i)
        iqi[i] = (long long)(iqi[i - 1]) * iq % mod, iqi2[i] = (long long)(iqi2[i - 1]) * iqi[i - 1] % mod;
    for (int i = 0, pi = 1; i < n; ++i, pi = (long long)(pi) * p % mod)
        f[i] = (long long)(a[i]) * pi % mod * iqi2[i] % mod; // 题解中的 f_j p^j q^{-\binom{j}{2}}
    for (int i = 0; i < n + m - 1; ++i)
        g[i] = qi2[i]; // 题解中的 q^{\binom{i + j}{2}}
    std::reverse(f.begin(), f.end()), f *= g;
    for (int i = 0; i < m; ++i)
        printf("%lld ", (long long)(iqi2[i]) * f[n + i - 1] % mod);
    return 0;
}
// 尘埃覆写往日的档案,因果重构奇迹的日常
// 12