[總結]多項式求逆代替分治 $\text{FFT}$

因爲我懶得不想學蠢得學不會分治 \(\text{FFT}\) ,發現能夠用多項式求逆來完整地代替...ui

文章節選自分治 FFT 與多項式求逆,轉載方便本身查看。更多多項式求逆和分治 \(\text{FFT}\) 的內容與聯繫,可參見原博客。spa

問題提出

給定 \(\forall i\in[1,n),g[i]\),求遞推式code

\[f[i]=\begin{cases}1 & \text{ if } i=0 \\ \sum_{j=1}^if[i-j]g[j] & \text{ otherwise }\end{cases}\]get

的前 \(n\) 項,即 \(\forall i\in[0,n),f[i]\) 。答案模費馬質數。博客

求逆代替分治

這是分治 \(\text{FFT}\) 的板子,咱們考慮用多項式求逆來解決這個問題。it

考慮翻轉 \(f\)\(g\) 的地位ast

\[\begin{aligned}f[i]&=\sum_{j=1}^if[i-j]g[j]\\\Rightarrow g[i]&=f[i]-\sum_{j=1}^{i-1}f[j]g[i-j]\end{aligned}\]class

im

\[h[i]=\begin{cases}1 & \text{ if } i=0\\ -g[i] & \text{ otherwise }\end{cases}\]

\[\begin{aligned}g[i]&=f[i]-\sum_{j=1}^{i-1}f[j]g[i-j]\\g[i]&=f[i]+\sum_{j=1}^{i-1}f[j](-g[i-j])\\g[i]&=f[i]h[0]+\sum_{j=1}^{i-1}f[j]h[i-j]\\g[i]&=\sum_{j=1}^if[j]h[i-j]\end{aligned}\]

咱們須要讓 \(j\) 的下界變成 \(0\),才能真正變成卷積。所以咱們強行將 \(f[0]\) 置爲 \(0\) 便可。爲何在這裏不會出現問題?由於在這裏,\(f[0]\) 根本不會和 \(h\) 作卷積(咱們初始的下界就是 \(j=1\))。

\[\begin{aligned}g[i]&=\sum_{j=0}^if[j]h[i-j]\\G&=F\otimes H\end{aligned}\]

咱們用多項式求逆算出 \(H^{-1}\),而後和 \(G\) 相乘便可獲得 \(F\)

代碼實現

#include <bits/stdc++.h>
using namespace std;
const int N = (100000+5)<<2, yzh = 998244353;

int n, g[N], h[N], f[N], tmp[N], R[N], len, L;

int quick_pow(int a, int b) {
    int ans = 1;
    while (b) {
        if (b&1) ans = 1ll*ans*a%yzh;
        b >>= 1, a = 1ll*a*a%yzh;
    }
    return ans;
}
void NTT(int *A, int o) {
    for (int i = 0; i < len; i++) if (i < R[i]) swap(A[i], A[R[i]]);
    for (int i = 1; i < len; i <<= 1) {
        int gn = quick_pow(3, (yzh-1)/(i<<1));
        if (o == -1) gn = quick_pow(gn, yzh-2);
        for (int j = 0; j < len; j += (i<<1))
            for (int k = 0, x, y, g = 1; k < i; k++, g = 1ll*g*gn%yzh) {
                x = A[j+k], y = 1ll*g*A[i+j+k]%yzh;
                A[j+k] = (x+y)%yzh, A[i+j+k] = (x-y)%yzh;
            }
    }
    if (o == 1) return;
    for (int i = 0, inv = quick_pow(len, yzh-2); i < len; i++)
        A[i] = 1ll*A[i]*inv%yzh;
}
void poly_inv(int *A, int *B, int deg) {
    if (deg == 1) {B[0] = quick_pow(A[0], yzh-2); return; }
    poly_inv(A, B, (deg+1)>>1);
    for (L = 0, len = 1; len <= (deg<<1); len <<= 1) ++L;
    for (int i = 0; i < len; i++) R[i] = (R[i>>1]>>1)|((i&1)<<L-1);
    for (int i = 0; i < deg; i++) tmp[i] = A[i];
    for (int i = deg; i < len; i++) tmp[i] = 0;
    for (int i = (deg+1)>>1; i < len; i++) B[i] = 0;
    NTT(tmp, 1), NTT(B, 1);
    for (int i = 0; i < len; i++) B[i] = 1ll*B[i]*(2ll-1ll*B[i]*tmp[i]%yzh)%yzh;
    NTT(B, -1);
}
void work() {
    scanf("%d", &n);
    for (int i = 1; i < n; i++) scanf("%d", &g[i]);
    h[0] = 1;
    for (int i = 1; i < n; i++) h[i] = -g[i];
    poly_inv(h, f, n);
    for (int i = n; i < len; i++) f[i] = 0;
    NTT(f, 1), NTT(g, 1);
    for (int i = 0; i < len; i++) f[i] = 1ll*f[i]*g[i]%yzh;
    NTT(f, -1); f[0] = 1;
    for (int i = 0; i < n; i++) printf("%d ", (f[i]+yzh)%yzh);
}
int main() {work(); return 0; }
相關文章
相關標籤/搜索