多项式全家桶

Galx_Trail · 2026-08-13 15:44:17 · 算法·理论

前言：lyy还是太强了/bx

前置芝士：

代数基本定理：

1. 一个 d 次多项式可以被 d+1 个点唯一确定。

2. 还有要用的后面补吧

虚数，复数杂烩：

基础定义：i^2=-1。

复平面：想象一下，就是在一维实数轴上扩展一个虚数轴。

复数的表示：

1. 代数法：a+bi，其中 a 为实部，bi 为虚部，这个应该好理解。

2. 极角表示，想象一下，复平面上一个点到 O 为一个向量，其模长为 r，辐角为 \theta，则可以表示为：r\times \cos\theta + r\times i\sin\theta。

在复数计算上，加法就相当于向量相加，乘法相当于模长相乘，辐角相加，这在后面理解 FFT 有帮助。

欧拉公式

还是从复平面的视角理解，物理意义就是当 \theta \in [0,2\pi]，e^{i\theta} 刚好画出一个单位圆。

单位复数根

形如，求解一个方程：

在实数域上只有 1（偶次方还有 -1），根据代数基本定理，在复数域上有 n 个根，我们怎么找这 n 个根？

既然 1 在复平面上是 (1,0)，对应着 2\pi 的整数倍，记为 2k\pi，对上面展开：

这 n 个根均匀分布在单位圆上，令逆时针第一个点为主 n 次单位根。 记为 \omega_n。 剩下的所有根都是 \omega_n 多次旋转相同角度得到的，即：\omega_n^0,\omega_n^1,\omega_n^2\dots\omega_n^{n-1}。

一些引理

代数证明读者自证，应该挺简单的（？

1. 消去引理：

复平面上本质角度 \theta 是不变的，所以它们也是相等的。

2. 折半引理

如果 n 是偶数，那么这 n 个单位根的平方，恰好对应 \frac{n}{2} 个 \frac{n}{2} 次单位根，每个出现 2 次。

为什么？ 想一想平方意味着什么？ 是不是角度相乘？

那么以 n=8 为例子，上下部分各 4 个点，平方后，上面 4 个点旋转到下面 4 个点的位置完全重合，下面同理，所以得证。 这是最核心的，是 FFT 能够做到 O(n\log n) 的原理。

3. 求和引理

当然可以直接等比数列做，只是说需要证明其能作用在复数域上。 这里还是以复平面的视角看：

根据 n 次单位根的定义，n 个根构成了一个正 n 边形，把它们想象成 n 个与 O 构成的向量，因为它是一个完美对称的旋转体，所以向量和一定为 0。

拉格朗日插值

一般公式

已知 m 个点 (x_1,y_1),(x_2,y_2)\dots(x_m,y_m)（其中 x_i 互不相同），寻找一个最高次数为 m-1 次的多项式 f(x) 使得：\forall i, f(x_i)=y_i。

不考虑证明 （其实我不会），直接给出构造：

感性的理解，就是只有当 x=x_i 才有 y_i 的值，这么构造能完全满足所有条件。

有了这个，就可以做模板题了：

#include <iostream>
#include <cstring>
#include <vector>
#include <algorithm>
#include <iomanip>
#include <unordered_map>
#include <set>
#include <functional>
#include <numeric>
using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
    using ll = long long;
    using ull = unsigned long long;
    using i128 = __int128;
    void init() {
//      freopen("1.in","r",stdin);
//      freopen("my.out","w",stdout);
        ios::sync_with_stdio(false);
        cin.tie(nullptr);
        cout.tie(nullptr);
    }
    constexpr ll mod = 998244353;
    ll fpow(ll a, ll b){
        ll res = 1;
        a %= mod;
        while (b){
            if (b & 1) res = res * a % mod;
            a = a * a % mod;
            b >>= 1;
        }
        return res;
    }
    ll inv(ll x){
        return fpow(x, mod - 2);
    }
    void Main() {
        ll k;
        int n; cin >> n >> k;
        vec<pair<ll,ll>> a(n + 1);
        for (int i = 1; i <= n; ++i) cin >> a[i].first >> a[i].second;
        ll ans = 0;
        for (int i = 1; i <= n; ++i) {
            ll s1 = a[i].second % mod, s2 = 1;
            for (int j = 1; j <= n; ++j){
                if (i == j) continue;
                s1 = (s1 * ((k - a[j].first) % mod + mod)) % mod;
                s2 = (s2 * ((a[i].first - a[j].first) % mod + mod)) % mod % mod;
            }
            ans = (ans + s1 * inv(s2) % mod) % mod;
        }
        cout << ans;
    }
}
int main() {
    CZW::init();
    int Test = 1;
//  cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

拓展问题

例如：

定义函数 S_k(n)=\sum_{i=1}^{n} i^k，求 S_k(n) \pmod p 的值，（其中 n\le 10^{18},k\le 10^6）。

设 f(n)=S_k(n)，则前向差分 \Delta f(n)=f(n + 1) - f(n) = (n + 1) ^ k。

根据有限差分定理，既然 \Delta f(n) 为 k 次多项式，则 S_k(n) 为一个 k+1 次多项式。

根据前置芝士，我们知道需要 k + 2 个点来唯一确认这个多项式。

再看数据，显然 O(k^2) 做不了，需要优化。

回头看过程，我们发现点是可以自己取的，那我们可以取特殊点值来优化，例如：令 x_i=i。

代入原式：

考虑分母，分子分开优化：

• 分子：

为了避免出现 0 的情况，我们可以预处理前缀积与后缀积即可。

• 分母：

显然这个可以用阶乘逆元处理。

但要注意，我们求插值是对一个前缀和求，所以我们构造的纵坐标是对 1 到 i 的前缀 j^k 和。

#include <iostream>
#include <cstring>
#include <vector>
#include <algorithm>
#include <iomanip>
#include <unordered_map>
#include <set>
#include <functional>
#include <numeric>
using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
    using ll = long long;
    using ull = unsigned long long;
    using i128 = __int128;
    constexpr ll mod = 1e9 + 7, N = 1e6 + 5;
    ll fpow(ll a, ll b){
        ll res = 1;
        a %= mod;
        while(b){
            if (b & 1) res = res * a % mod; 
            a = a * a % mod;
            b >>= 1;
        }
        return res;
    }
    int primes[N], cnt;
    ll pk[N];
    bool is[N];
    void init() {
//      freopen("1.in","r",stdin);
//      freopen("my.out","w",stdout);
        ios::sync_with_stdio(false);
        cin.tie(nullptr);
        cout.tie(nullptr);
        pk[1] = 1;
    }
    void Main() {
        int n, k;
        cin >> n >> k;
        for (int i = 2; i < N; ++i){
            if (!is[i]){
                primes[++cnt] = i;
                pk[i] = fpow(i, k);
            }
            for (int j = 1; j <= cnt && i * primes[j] < N; ++j){
                is[i * primes[j]] = true;
                pk[i * primes[j]] = pk[i] * pk[primes[j]] % mod;
                if (i % primes[j] == 0){
                    break;
                }
            }
        }
        vec<ll> y(k + 5, 0);
        y[0] = 0;
        for (int i = 1; i < k + 5; ++i) {
            y[i] = (y[i - 1] + pk[i]) % mod;
        }
        if (n <= k + 2) {
            cout << y[n];
            return ;
        }
        int m = k + 2;
        vec<ll> pre(m + 3, 0), suf(m + 3, 0);
        pre[0] = 1, suf[m + 1] = 1;
        for (int i = 1; i <= m; ++i) {
            pre[i] = pre[i - 1] * (n - i) % mod;
        }
        for (int i = m; i >= 1; --i) {
            suf[i] = suf[i + 1] * (n - i) % mod;
        }
        vec<ll> inv(m + 3, 0), fact(m + 3, 0);
        fact[0] = 1;
        for (int i = 1; i <= m; ++i) fact[i] = fact[i - 1] * i % mod;
        inv[m] = fpow(fact[m], mod - 2);
        for (int i = m - 1; i >= 1; --i) inv[i] = inv[i + 1] * (i + 1) % mod;
        inv[0] = 1;
        ll ans = 0;
        for (int i = 1; i <= m; ++i){
            ll sum1 = pre[i - 1] * suf[i + 1] % mod;
            ll sum2 = inv[i - 1] * inv[m - i] % mod;
            ll sum3 = y[i] * sum1 % mod * sum2 % mod;
            if ((m - i) % 2 == 0) {
                ans = (ans + sum3) % mod;
            }else {
                ans = ((ans - sum3) % mod + mod) % mod;
            }
        }
        cout << ans;
    }
}
int main() {
    CZW::init();
    int Test = 1;
//  cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

P5364 [SNOI2017] 礼物

由题意得：

求：s_n=2\times s_{n-1}+n^k。

注意到 n 很大，k 很小，而且这是一个递推式，考虑用矩阵加速优化。

考虑对后面的常数项进行二项式展开：

所以，我们可以将 n^0,n^1,n^2\dots n^k 加入矩阵转移：

如何构造 base 辅助转移？ 考虑到转移后第 0 行 s_{i+1}=2\times s_i+(i+1)^k。

令 base_{0,0}=2 且按照上面转化 \forall i \in [0,k], base_{0,i+1}=\binom{k}{i}，其余为 0。

对于 i\in [0,k]，有 (n+1)^i=\sum_{j=0}^{i}{\binom{i}{j}\times n^j}，所以：

令 base_{i+1,j+1}=\binom{i}{j}，其余为 0。

最后就是矩阵快速幂和组合数的板子了。

#include <iostream>
#include <cstring>
#include <vector>
#include <algorithm>
#include <iomanip>
#include <unordered_map>
#include <set>
#include <functional>
#include <numeric>
using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
    using ll = long long;
    using ull = unsigned long long;
    using i128 = __int128;
    constexpr ll mod = 1e9 + 7;
    void init() {
//      freopen("1.in","r",stdin);
//      freopen("my.out","w",stdout);
        ios::sync_with_stdio(false);
        cin.tie(nullptr);
        cout.tie(nullptr);
    }
    int sz;
    ll fpow(ll a, ll b) {
        ll res = 1;
        a %= mod;
        while (b) {
            if (b & 1) res = res * a % mod;
            a = a * a % mod;
            b >>= 1;
        }
        return res;
    }
    struct Mat {
        ll mat[15][15];
        Mat(){
            for (int i = 0; i <= sz; ++i) for (int j = 0; j <= sz; ++j) mat[i][j] = 0;
        }
    };
    Mat mul(const Mat& a, const Mat& b){
        Mat c;
        for (int i = 0; i < sz; ++i){
            for (int kk = 0; kk < sz; ++kk){
                if (!a.mat[i][kk]) continue;
                for (int j = 0; j < sz; ++j){
                    c.mat[i][j] = (c.mat[i][j] + a.mat[i][kk] * b.mat[kk][j] % mod) % mod;
                }
            }
        }
        return c;
    }
    Mat Mat_fpow(Mat a, ll b) {
        Mat res;
        for (int i = 0; i < sz; ++i) res.mat[i][i] = 1;
        while(b){
            if (b & 1) res = mul(res, a);
            a = mul(a, a);
            b >>= 1;
        }
        return res;
    }
    void Main() {
        ll n; int k;
        cin >> n >> k;
        sz = k + 2;
        vec<vec<ll>> C(15, vec<ll> (15, 0));
        for (int i = 0; i < 15; ++i){
            C[i][0] = 1;
            for (int j = 1; j <= i; ++j){
                C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % mod;
            }
        }
        Mat base;
        base.mat[0][0] = 2;
        for (int j = 0; j <= k; ++j) base.mat[0][j + 1] = C[k][j];
        for (int i = 0; i <= k; ++i){
            for (int j = 0; j <= i ;++j){
                base.mat[i + 1][j + 1] = C[i][j];
            }
        }
        Mat res = Mat_fpow(base, n - 1);
        ll sum = res.mat[0][1];
        cout << (sum + fpow(n, k)) % mod;
    }
}
int main() {
    CZW::init();
    int Test = 1;
//  cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

FFT

这个玩意我觉得是学 NTT 和 FWT 的基础，它们本质几乎一样的，当然理解方法很多，这给出一种我觉得比较好理解的方法（？

*先考虑一个问题：给定两个多项式 A 和 B 让你求其卷积，即求 $C=AB$。**

显然如果我们暴力按照系数去每位相乘是 O(n^2) 的。 这时 FFT 的思想就是将其转化为点值表达式对应相乘再还原。

形式化的讲，找一个 n\times n 的可逆变换矩阵 M，使得：

这样的话，C=M^{-1}(M_A \circ M_B)。

那 FFT 怎么处理的呢？

FFT 将变换空间引入到复数域，利用 n 次单位根的性质，分治求出点值，以达到快速求解。

具体怎么求点值？ 举一个例子：我们有一个次数界为 n 的多项式 A，表示为 A(x)=a_0+a_1x+a_2x^2+\dots a_{n-1}x^{n-1}，求出 n 个单位复数根。

这时绝妙的点来了！ 将其拆成奇数次项与偶数次项组成的多项式：

于是我们可以重组原多项式：A(x)=A^{[0]}(x^2)+x\times A^{[1]}(x^2)。

现在我们需要求解原多项式前半段 \omega_n^k 和后半段 \omega_n^{k+n/2}。

先求前半段（y_k）：

根据折半引理：

再求后半段（y_{k+n/2}）：

因为 \omega_n^n=1 且 \omega_n^{n/2}=-1，代入得：

这时观察两个式子，本质差别就在正负号上！ 这说明只要我们把子问题求出即可求出原问题。

合并就是著名的蝴蝶操作：

有了：

即可：

可以结合图来理解。 应该是长得像蝴蝶所以叫蝴蝶操作？

现在我们有一些点值，怎么做逆变换呢？

观察发现，FFT 的可逆变换矩阵 M 就是有名的范德蒙德矩阵 (V_n)_{j,k}=\omega_n^{jk}，其逆矩阵 (V_n)^{-1}_{j,k}=\frac{1}{n}\omega_{n}^{-jk}（显然还需要证明范德蒙德矩阵有逆矩阵，只需证明 \det(V_n) 不为 0 即可）。

考虑证明，令 U=(V_n)^{-1}：

所以当 j=k 时为 1，其余根据求和引理始终为 0。

这时我们就有了一个递归的写法了，求逆只需按照上面改一下即可。

当然实际上我们实现中并不会使用递归，常数大且容易爆栈。

这就要说到迭代实现 FFT 了

给一个例子：

观察最后一层我们所处理出来的最小子问题，若它们成为一个序列，与原序列二进制下标作对比：

观察到二进制下标是反转的！ 所以我们 O(n) 按照底层样式排序，进行从低到高做 FFT。

大概长这样：

这样我们就学完 FFT 了！ 可以先做一道模板题。

#include <algorithm>
#include <cmath>
#include <cstring>
#include <functional>
#include <iomanip>
#include <iostream>
#include <numeric>
#include <set>
#include <unordered_map>
#include <vector>

using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
    using ll = long long;
    using ull = unsigned long long;
    using i128 = __int128;
    using lb = long double;
    const lb pi = acos(-1.0);
    void init() {
        //      freopen("1.in","r",stdin);
        //      freopen("my.out","w",stdout);
        ios::sync_with_stdio(false);
        cin.tie(nullptr);
        cout.tie(nullptr);
    }
    constexpr int N = 4e6 + 5;
    struct Complex {
        lb r, i;
        Complex(lb rr = 0, lb ii = 0) : r(rr), i(ii) {}
        Complex operator+(const Complex &other) {
            return {r + other.r, i + other.i};
        }
        Complex operator-(const Complex &other) {
            return {r - other.r, i - other.i};
        }
        Complex operator*(const Complex &other) {
            return {r * other.r - i * other.i, r * other.i + i * other.r};
        }
    };
    Complex a[N], b[N];
    int rev[N];
    void FFT(Complex *A, int lim, int op) {
        for (int i = 0; i < lim; ++i)
            if (i < rev[i])
                swap(A[i], A[rev[i]]);
        for (int m = 1; m < lim; m <<= 1) {
            Complex wn(cos(pi / m), op * sin(pi / m));
            for (int j = 0; j < lim; j += (m << 1)) {
                Complex w = 1;
                for (int k = 0; k < m; ++k) {
                    Complex x = A[k + j], y = A[j + k + m] * w;
                    A[j + k] = x + y, A[j + k + m] = x - y;
                    w = w * wn;
                }
            }
        }
        if (op == -1) {
            for (int i = 0; i < lim; ++i)
                A[i].r /= lim;
        }
    }
    void Main() {
        int n, m;
        cin >> n >> m;
        for (int i = 0; i <= n; ++i)
            cin >> a[i].r;
        for (int j = 0; j <= m; ++j)
            cin >> b[j].r;
        int lim = 1, len = 0;
        int sz = n + m + 1;
        while (lim < sz) {
            lim <<= 1;
            ++len;
        }
        // cout << lim << endl;
        for (int i = 0; i < lim; ++i)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (len - 1));
        FFT(a, lim, 1), FFT(b, lim, 1);
        for (int i = 0; i < lim; ++i)
            a[i] = a[i] * b[i];
        FFT(a, lim, -1);
        for (int i = 0; i < sz; ++i)
            cout << (int)round(a[i].r) << ' ';
    }
} // namespace CZW
int main() {
    CZW::init();
    int Test = 1;
    // cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

还有一个高精度乘法的问题，其实就是这个的延伸，只用最后处理一下进位和前导 0 即可。

NTT

因为 FFT 用到浮点数，可能在取模难处理或者求答案时会被挂精度，我们要找一个东西去完全替代 \omega_n。

这直接给结论：用原根代替。

考虑证明：

对于一个质数 p 和任意非零整数 a，有 g^{p-1} \equiv 1 \pmod p 且 p-1 是满足这个条件的最小正整数。

这样的话 \forall i \in [1,p-2],g^i \not \equiv 1 \pmod p。

我们构造 NTT 单位根：

考虑是否满足三大引理：

1. 满足周期性 / 消去引理：(\omega_n)^n=(g^{\frac{p-1}{n}})^n=g^{p-1} \equiv 1 \pmod p。

2. 满足折半引理：(\omega_n^{k+\frac{n}{2}})^2 = \omega_n^{2k+n} = \omega_n^{2k} \cdot \omega_n^n \equiv \omega_n^{2k} \cdot 1 \pmod p。

3. 满足求和引理：\sum_{j=0}^{n-1} (\omega_n^k)^j \equiv \frac{1 - (\omega_n^k)^n}{1 - \omega_n^k} \pmod p。 其中分子 1 - (\omega_n^{n})^k \equiv 1 - 1 = 0 \pmod p，又因为 g 为原根，所以 1 - \omega_n^{k} \not\equiv 1 - 1 = 0 \pmod p，所以得证。

这样我们只用在 FFT 基础上套用原根即可（注意 NTT 对模数要求很苛刻，可以去搜一下）。

神级（秘）引用

分治 NTT 应用

套用 cdq 分治的思想，将前一半 NTT 先求出来，加入到右边去一起计算即可。

#include <algorithm>
#include <cmath>
#include <cstddef>
#include <cstring>
#include <functional>
#include <iomanip>
#include <iostream>
#include <numeric>
#include <set>
#include <unordered_map>
#include <vector>

using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
using ll = long long;
using ull = unsigned long long;
using i128 = __int128;
using lb = long double;
void init() {
    //      freopen("1.in","r",stdin);
    //      freopen("my.out","w",stdout);
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cout.tie(nullptr);
}
constexpr ll mod = 998244353, G = 3, GI = (mod + 1) / 3, N = 4e5 + 5;
ll fpow(ll a, ll b) {
    ll res = 1;
    a %= mod;
    while (b) {
        if (b & 1)
            res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}
int rev[N];
ll f[N], g[N];
ll a[N], b[N];
void NTT(ll *A, int n, int op) {
    for (int i = 0; i < n; ++i)
        if (i < rev[i])
            swap(A[i], A[rev[i]]);
    for (int m = 1; m < n; m <<= 1) {
        ll wn = fpow(op == 1 ? G : GI, (mod - 1) / (m << 1));
        for (int j = 0; j < n; j += m << 1) {
            ll w = 1;
            for (int k = 0; k < m; ++k) {
                ll x = A[j + k], y = w * A[j + k + m] % mod;
                A[j + k] = (x + y) % mod;
                A[j + k + m] = ((x - y) % mod + mod) % mod;
                w = w * wn % mod;
            }
        }
    }
    if (op == -1) {
        ll inv = fpow(n, mod - 2);
        for (int i = 0; i < n; ++i)
            A[i] = A[i] * inv % mod;
        ;
    }
}
void cdq(int l, int r) {
    if (l == r)
        return;
    int mid = (ll)(l + r) >> 1;
    cdq(l, mid);
    int lim = 1, len = 0;
    while (lim <= (mid - l + 1) + (r - l)) {
        lim <<= 1;
        ++len;
    }
    for (int i = 0; i < lim; ++i) {
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (len - 1));
    }
    for (int i = 0; i < lim; ++i)
        a[i] = b[i] = 0;
    for (int i = l; i <= mid; ++i)
        a[i - l] = f[i];
    for (int i = 1; i <= r - l; ++i)
        b[i - 1] = g[i];
    NTT(a, lim, 1), NTT(b, lim, 1);
    for (int i = 0; i < lim; ++i)
        a[i] = a[i] * b[i] % mod;
    NTT(a, lim, -1);
    for (int i = mid + 1; i <= r; ++i)
        f[i] = (f[i] + a[i - l - 1]) % mod;
    cdq(mid + 1, r);
}
void Main() {
    int n;
    cin >> n;
    for (int i = 1; i < n; ++i)
        cin >> g[i];
    f[0] = 1;
    cdq(0, n - 1);
    for (int i = 0; i < n; ++i)
        cout << f[i] << ' ';
}
} // namespace CZW
int main() {
    CZW::init();
    int Test = 1;
    // cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

期望转化为系数

看了好几遍才看懂 QwQ。。

设 X 为我们抽取袜子的总次数，要找期望 E[X]。

常规定义是：E[X]=\sum_k{k\times P(X=k)}，这非常难算。

但是我们可以做一个转化，对于非负整数随机变量，有一个等价公式：

意思是抽了 k 次还没结束，换句话讲就是前 k 次抽出来的袜子全是不同颜色。

那如何算 P(X>k) 的值？

假设有 S=\sum A_i，无顺序取出 k 只袜子，共有 \binom{S}{k} 种。

所以：P(X>k)=\frac{\text{选出 } k \text{ 只袜子颜色互不相同的方案数}}{\binom{S}{k}}。

那怎么求上面的东西？ 想象一下，对于第 i 种袜子，我们要不拿一个（A_i 种可能），要不就不拿（1 种可能），所以单个 i 的状态为：(1+A_ix)。

所以 P(x)=\prod_{i=1}^n(1+A_ix)。 暴力展开，对于第 k 次方的系数就代表着挑了 k 个 A_ix 和 n-k 个 1，换句话讲就是从 n 种颜色里挑出 k 种不同的颜色的方案数！

所以记 x^k 的系数为 [x^k]P(x)，则：

所以我们只需用分治 NTT 求出 P(x) 每项系数，即可解决这道题！

#include <algorithm>
#include <bit>
#include <cmath>
#include <cstring>
#include <functional>
#include <iomanip>
#include <iostream>
#include <numeric>
#include <set>
#include <unordered_map>
#include <vector>

using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
using ll = long long;
using ull = unsigned long long;
using i128 = __int128;
using lb = long double;
constexpr ll mod = 998244353, G = 3, GI = (mod + 1) / 3;
ll fpow(ll a, ll b) {
    ll res = 1;
    a %= mod;
    while (b) {
        if (b & 1)
            res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}
void init() {
    //      freopen("1.in","r",stdin);
    //      freopen("my.out","w",stdout);
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cout.tie(nullptr);
}
namespace Poly {
vec<int> rev;
void NTT(vec<ll> &A, int lim, int op) {
    if ((int)rev.size() != lim) {
        rev.resize(lim);
        int len = __builtin_ctz(lim);
        for (int i = 0; i < lim; ++i)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (len - 1));
    }
    for (int i = 0; i < lim; ++i)
        if (i < rev[i])
            swap(A[i], A[rev[i]]);
    for (int m = 1; m < lim; m <<= 1) {
        ll wn = fpow(op == -1 ? GI : G, (mod - 1) / (m << 1));
        for (int j = 0; j < lim; j += (m << 1)) {
            ll w = 1;
            for (int k = 0; k < m; ++k) {
                ll x = A[j + k] % mod, y = A[j + k + m] * w % mod;
                A[j + k] = (x + y) % mod;
                A[j + k + m] = ((x - y) % mod + mod) % mod;
                w = w * wn % mod;
            }
        }
    }
    if (op == -1) {
        ll inv = fpow(lim, mod - 2);
        for (int i = 0; i < lim; ++i)
            A[i] = A[i] * inv % mod;
    }
}
vec<ll> multiply(vec<ll> A, vec<ll> B) {
    int n = A.size(), m = B.size(), sz = n + m - 1;
    int lim = 1, len = 0;
    while (lim <= sz) {
        lim <<= 1;
        ++len;
    }
    if (lim == 1)
        return A;
    A.resize(lim, 0), B.resize(lim, 0);
    NTT(A, lim, 1), NTT(B, lim, 1);
    for (int i = 0; i < lim; ++i)
        A[i] = A[i] * B[i] % mod;
    NTT(A, lim, -1);
    A.resize(sz);
    return A;
}
} // namespace Poly
vec<ll> a;
vec<ll> cdq(int l, int r) {
    if (l == r) {
        return {1, a[l]};
    }
    int mid = (ll)(l + r) >> 1;
    return Poly::multiply(cdq(l, mid), cdq(mid + 1, r));
}
void Main() {
    int n;
    cin >> n;
    a.resize(n + 1);
    ll sum = 0;
    for (int i = 0; i < n; ++i) {
        cin >> a[i];
        sum += a[i];
    }
    vec<ll> d = cdq(0, n - 1);
    vec<ll> fact(n + 1, 1);
    for (int i = 1; i <= n; ++i)
        fact[i] = fact[i - 1] * i % mod;
    vec<ll> down(n + 1, 1), inv(n + 1, 1);
    for (int i = 1; i <= n; ++i) {
        ll val = (sum - i + 1) % mod;
        if (val < 0)
            val += mod;
        down[i] = down[i - 1] * val % mod;
    }
    inv[n] = fpow(down[n], mod - 2);
    for (int i = n - 1; i >= 0; --i) {
        ll val = (sum - i) % mod;
        if (val < 0)
            val += mod;
        inv[i] = inv[i + 1] * val % mod;
    }
    ll ans = 0;
    for (int i = 0; i <= n; ++i) {
        if (i < (int)d.size()) {
            ll num = d[i] * fact[i] % mod;
            num = num * inv[i] % mod;
            ans = (ans + num) % mod;
        }
    }
    cout << ans;
}
} // namespace CZW
int main() {
    CZW::init();
    int Test = 1;
    // cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

FWT

本质上其实这玩意跟高维前缀和很像，因为位运算的性质以至于每一位互不影响，就像二维前缀和，我们可以对 i 维先求再将求完的答案更新对 j 维再求一次即可（当然你也可以用克罗内克积和转移矩阵来理解，当然我是不会的）。

所以我们可以从一维开始考虑，以异或为例子：

给定 A=[A_0,A_1] 和 B=[B_0,B_1]，求异或卷积 C=A \oplus B。

根据异或规则：

跟上面类似，FWT 还是想处理出一个 A' 和 B' 使得 C' 等于对应位置乘积。

那么就开始凑数（可以自己尝试一下相加，相减），这直接给结论：

逆变换也非常简单，已知 C'_0=C_0+C_1 和 C'_1=C_0-C_1，求 C_0,C_1，显然是一个二元一次方程组，随便解。

既然一维会了，有了上面的性质，就可以对整体进行分治，因为高 / 低位互相影响不到，我们可以从高位到低位进行分治：

假设数组长度为 4（下标 00, 01, 10, 11）：

1. 按最高位分治：最高位为 0 的一半是：A_{00}, A_{01}（简称左半边）最高位为 1 的一半是：A_{10}, A_{11}（简称右半边）。

因为最高位的异或不受低位影响，我们直接对这两半整体使用“魔法”：

新左半边 = 老左半边 + 老右半边。 新右半边 = 老左半边 - 老右半边。

2. 按次高位（最低位）分治：现在左右半边各自的内部，再根据次高位继续做同样的加减法。

AND 和 OR 的凑数比 XOR 简单，可以自己推一下。 （反正给代码）

这样就可以写模板了。 毕竟像我这么蒻的都能一次写对

#include <algorithm>
#include <cmath>
#include <cstring>
#include <functional>
#include <iomanip>
#include <iostream>
#include <numeric>
#include <set>
#include <unordered_map>
#include <vector>

using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
    using ll = long long;
    using ull = unsigned long long;
    using i128 = __int128;
    using lb = long double;
    void init() {
        //      freopen("1.in","r",stdin);
        //      freopen("my.out","w",stdout);
        ios::sync_with_stdio(false);
        cin.tie(nullptr);
        cout.tie(nullptr);
    }
    constexpr ll mod = 998244353, inv2 = (mod + 1) / 2;
    void FWT(vec<ll> &A, int lim, int op, int type) {
        for (int m = 1; m < lim; m <<= 1) {           // 枚举处理到二进制第几位
            for (int j = 0; j < lim; j += (m << 1)) { // 以2*m为一个块跳着处理
                for (int k = 0; k < m; ++k) {
                    if (type == 1) {
                        if (op == 1)
                            (A[j + k] += A[j + k + m]) %= mod;
                        else
                            (A[j + k] += mod - A[j + k + m]) %= mod;
                    } else if (type == 2) {
                        if (op == 1)
                            (A[j + k + m] += A[j + k]) %= mod;
                        else
                            (A[j + k + m] += mod - A[j + k]) %= mod;
                    } else {
                        ll x = A[j + k], y = A[j + k + m];
                        if (op == 1) {
                            A[j + k] = (x + y) % mod;
                            A[j + k + m] = ((x - y) % mod + mod) % mod;
                        } else {
                            A[j + k] = (x + y) % mod * inv2 % mod;
                            A[j + k + m] =
                                ((x - y) % mod + mod) % mod * inv2 % mod;
                        }
                    }
                }
            }
        }
    }
    vec<ll> solve(vec<ll> a, vec<ll> b, int op, int lim) {
        FWT(a, lim, 1, op), FWT(b, lim, 1, op);
        for (int i = 0; i < lim; ++i)
            a[i] = a[i] * b[i] % mod;
        FWT(a, lim, -1, op);
        return a;
    }
    void Main() {
        int n;
        cin >> n;
        int lim = 1 << n;
        vec<ll> a(lim + 5, 0), b(lim + 5, 0);
        for (int i = 0; i < lim; ++i)
            cin >> a[i];
        for (int i = 0; i < lim; ++i)
            cin >> b[i];
        vec<ll> ans_and = solve(a, b, 1, lim), ans_or = solve(a, b, 2, lim);
        vec<ll> ans_xor = solve(a, b, 3, lim);
        for (int i = 0; i < lim; ++i)
            cout << ans_or[i] << ' ';
        cout << endl;
        for (int i = 0; i < lim; ++i)
            cout << ans_and[i] << ' ';
        cout << endl;
        for (int i = 0; i < lim; ++i)
            cout << ans_xor[i] << ' ';
    }
} // namespace CZW
int main() {
    CZW::init();
    int Test = 1;
    // cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

FWT好题

看到 n 很小，考虑枚举行是否翻转的状态，通过神级观察（为什么都说很版啊，我蒻成区了），翻转任意行所对列产生的影响：翻转行的状态 i \oplus 某个列状态 j。

这么看挺抽象的，举个例子就懂了，假设是行状态 i=101 即翻转第 0 行和第 1 行，某个列状态 j=110 那我们朴素翻转后 j \rightarrow j'=011，这时候再看上面结论，你会发现，i \oplus j=011=j'。

这下一切都简单了，令 A_x 为原始列向量的值，B_x 表示某个行状态所产生的最小贡献。

所以：

转化一下式子：

显然是一个FWT板子，这样就做完了。

#include <algorithm>
#include <cmath>
#include <cstring>
#include <functional>
#include <iomanip>
#include <iostream>
#include <numeric>
#include <set>
#include <unordered_map>
#include <vector>

using namespace std;
namespace CZW {
#define endl "\n"
#define vec std::vector
#define pb push_back
#define eb emplace_back
using ll = long long;
using ull = unsigned long long;
using i128 = __int128;
using lb = long double;
void init() {
    //      freopen("1.in","r",stdin);
    //      freopen("my.out","w",stdout);
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cout.tie(nullptr);
}
constexpr int N = 1 << 20, M = 1e5 + 5;
ll a[N << 1 | 1], b[N << 1 | 1];
int val[M];
void FWT(ll *A, int lim, int op) {
    for (int m = 1; m < lim; m <<= 1) {
        for (int j = 0; j < lim; j += (m << 1)) {
            for (int k = 0; k < m; ++k) {
                ll x = A[j + k], y = A[j + k + m];
                if (op == 1) {
                    A[j + k] = x + y;
                    A[j + k + m] = x - y;
                } else {
                    A[j + k] = (x + y) / 2;
                    A[j + k + m] = (x - y) / 2;
                }
            }
        }
    }
}
void Main() {
    int n, m;
    cin >> n >> m;
    string s;
    for (int i = 0; i < n; ++i) {
        cin >> s;
        for (int j = 0; j < m; ++j) {
            if (s[j] == '1')
                val[j] |= 1 << i;
        }
    }
    int lim = 1 << n;
    for (int i = 0; i < m; ++i)
        ++a[val[i]];
    for (int i = 0; i < lim; ++i) {
        int x = __builtin_popcount(i);
        b[i] = min(x, n - x);
    }
    FWT(a, lim, 1), FWT(b, lim, 1);
    for (int i = 0; i < lim; ++i)
        a[i] *= b[i];
    FWT(a, lim, -1);
    ll ans = 4e18;
    for (int i = 0; i < lim; ++i)
        ans = min(ans, a[i]);
    cout << ans;
}
} // namespace CZW
int main() {
    CZW::init();
    int Test = 1;
    // cin >> Test;
    while (Test--)
        CZW::Main();
    return 0;
}

给一道比较经典用FWT的树上异或Trick，祭奠我逝去的一整个下午。。。。

提示一下：题目可以转化为求满足条件路径异或和两个点的方案数。而图上异或路径有一个 Trick 就是，对于任意两个点 (u,v) 其异或路径为生成树上简单路径 + 任意环张成的线性空间。