多項式不全家桶


多項式不全家桶(暫不更新)

前置知識需要FFT和NTT。

多項式求逆

板子題鏈接
求逆指的是給定一個多項式 \(F(x)\),你需要求出一個多項式 \(G(x)\),使其滿足 \(F(x) * G(x)\equiv1\pmod {x^n}\)

考慮為什么要在 \(\pmod {x^n}\) 情況下來求解,其實多項式可以近似理解為級數,顯然級數的項數是無窮的,所以我們用 \(\pmod {x^n}\) 來使得這個多項式的項數有窮。(大致理解就行)

首先大致的思路是分治解決問題。

顯然的如果多項式只有 \(1\) 項,那么解就是其逆元。

我們考慮從原式在 mod \(x^{\lceil\frac n 2\rceil}\) 的推到 mod \(x^ n\) 的情 況:

\[F(x)* G_1(x)\equiv1\pmod {x^{\lceil\frac n 2\rceil}} \]

其中 \(G_1(x)\) 指在 \(\pmod {x^{\lceil\frac n 2\rceil}}\) 的情況下的解,\(G(x)\) 指在 \(\pmod {x^n}\) 的情況下的解。

那么也有:

\[F(x)* G(x)\equiv1\pmod {x^{\lceil\frac n 2\rceil}} \]

兩個式子相減可以得到:

\[G_1(x)-G(x)\equiv0\pmod{x^{\lceil\frac n 2\rceil}} \]

平方得 :

\[G_1^{2}(x)-2G_1(x)G(x)+G^2(x)\equiv0\pmod{x^n} \]

兩邊同乘多項式 \(F(x)\) 得:

\[F(x) * G_1^{2}(x) - 2G_1(x) + G(x)\equiv0\pmod{x^n} \]

移項顯然得:

\[G(x)\equiv(2-F(x)G_1(x))G_1(x)\pmod{x^n} \]

然后你發現右邊的東西已經可以計算了。

復雜度非常顯然的 \(O(nlogn)\)

但是需要注意幾個問題,我會在代碼里寫出來。

代碼
#include <bits/stdc++.h>
using namespace std;
const int N = 4e5 + 5, mod = 998244353, g1 = 3, g2 = (mod + 1) / 3;
int n, a[N], b[N], d[N], A[N], B[N];
inline int qpow(int a, int b, int res = 1) {
    for(; b; b >>= 1, a = 1ll * a * a % mod) 
        if(b & 1) res = 1ll * res * a % mod;
    return res;
}
void NTT(int *a, int flag, int len) {
    for(register int i = 0; i < len; ++i) if(i < d[i]) swap(a[i], a[d[i]]);
    for(register int l = 2, m = 1; l <= len; l <<= 1, m <<= 1) {
        int w1 = qpow(flag == 1 ? g1 : g2, (mod - 1) / l);
        for(register int *p = a; p != a + len; p += l) {
            int w = 1;
            for(register int i = 0; i < m; ++i, w = 1ll*w*w1%mod) {
                int t = 1ll * w * p[i + m] % mod;
                p[i + m] = (p[i] - t + mod) % mod, p[i] = (p[i] + t) % mod;
            }
        }
    }
}
void niyuan(int *a, int *b, int n) {
    memset(A, 0, n * 16), memset(B, 0, n * 16);//如果需要多次計算的話,別忘了清零
    b[0] = qpow(a[0], mod - 2);
    for(register int l = 2, m = 1; l <= n << 2; l <<= 1, m <<= 1) {//注意我們需要跑到 l < n * 4,即4倍
        for(register int i = 0; i < m; ++i) A[i] = a[i], B[i] = b[i];
        for(register int i = 0; i < l; ++i) d[i] = (d[i>>1] >> 1) | ((i & 1) ? m : 0);
        NTT(A, 1, l), NTT(B, 1, l);
        for(register int i = 0; i < l; ++i) b[i] = ((2ll - 1ll*A[i]*B[i]%mod) * B[i] % mod + mod) % mod;
        NTT(b, -1, l);
        int ni = qpow(l, mod - 2);
        for(register int i = 0; i < l; ++i) b[i] = 1ll * b[i] * ni % mod;
        for(register int i = m; i < l; ++i) b[i] = 0;//為什么每次需要清零,主要是由推導的式子決定的。
    }
}
int main() {
    scanf("%d", &n);
    for(register int i = 0; i < n; ++i) scanf("%d", &a[i]);
    niyuan(a, b, n);
    for(register int i = 0; i < n; ++i) printf("%d ", b[i]);
    return 0;
}

多項式開根

板子題鏈接

多項式開根是已知多項式 \(A(x)\) ,在 \(\pmod {x^n}\) 的條件下,求得一個多項式 \(B(x)\),使得 \(B^2(x) = A(x) \equiv1\pmod {x^n}\)

我們依舊采用分治的思想進行解決。
一般我們習慣設 \(B(x)\) 的常數項為 1,這樣主要是為了方便計算。

我們考慮從原式在 mod \(x^{\lceil\frac n 2\rceil}\) 的推到 mod \(x^ n\) 的情 況:

\[B^2(x)\equiv A(x) (mod\ x^{n}) \]

設多項式 \(B'^2(x)\) 為在 \((mod\ x^{\frac{n}{2}})\) 的情況下的解,
那么:

\[B'^2(x)\equiv A(x) (mod\ x^{\frac{n}{2}}) \]

然后兩個式子相減得到 :

\[B(x)-B'(x)\equiv 0(mod\ x^{\frac{n}{2}}) \]

然后平方得 :

\[A(x)-2B(x)B'(x)+B'^2(x)\equiv 0(mod\ x^n) \]

最后移項得 :

\[B(x)\equiv \dfrac{A(x)+B'^2(x)}{2B'(x)} (mod\ x^n) \]

非常顯然的右邊的式子可以用求逆和乘法解決掉。

復雜度 \(O(nlogn)\)

代碼
#include <bits/stdc++.h>
using namespace std;
const int N = 8e5 + 5, mod = 998244353, g1 = 3, g2 = (mod + 1) / 3, ni2 = (mod + 1) / 2;
int n, a[N], b[N], d[N], A[N], B[N], inv[N], c[N];
inline int qpow(int a, int b, int res = 1) {
    for(; b; b >>= 1, a = 1ll * a * a % mod) 
        if(b & 1) res = 1ll * res * a % mod;
    return res;
}
void NTT(int *a, int flag, int len) {
    for(register int i = 0; i < len; ++i) if(i < d[i]) swap(a[i], a[d[i]]);
    for(register int l = 2, m = 1; l <= len; l <<= 1, m <<= 1) {
        int w1 = qpow(flag == 1 ? g1 : g2, (mod - 1) / l);
        for(register int *p = a; p != a + len; p += l) {
            int w = 1;
            for(register int i = 0; i < m; ++i, w = 1ll*w*w1%mod) {
                int t = 1ll * w * p[i + m] % mod;
                p[i + m] = (p[i] - t + mod) % mod, p[i] = (p[i] + t) % mod;
            }
        }
    }
}
void niyuan(int *a, int *b, int n) {
    memset(A, 0, n * 16), memset(B, 0, n * 16);
    b[0] = qpow(a[0], mod - 2);
    for(register int l = 2, m = 1; l < n << 2; l <<= 1, m <<= 1) {
        for(register int i = 0; i < m; ++i) A[i] = a[i], B[i] = b[i];
        for(register int i = 0; i < l; ++i) d[i] = (d[i>>1] >> 1) | ((i & 1) ? m : 0);
        NTT(A, 1, l), NTT(B, 1, l);
        for(register int i = 0; i < l; ++i) b[i] = ((2ll - 1ll*A[i]*B[i]%mod) * B[i] % mod + mod) % mod;
        NTT(b, -1, l);
        int ni = qpow(l, mod - 2);
        for(register int i = 0; i < l; ++i) b[i] = 1ll * b[i] * ni % mod;
        for(register int i = m; i < l; ++i) b[i] = 0;
    }
}
void kaigen(int *a, int *b, int n) {
    memset(c, 0, n * 16), memset(inv, 0, n * 16);
    b[0] = 1;//默認
    for(register int l = 2, m = 1; l < n << 2; l <<= 1, m <<= 1) {
        for(register int i = 0; i < m; ++i) c[i] = a[i];
        niyuan(b, inv, m);
        for(register int i = 0; i < l; ++i) d[i] = (d[i>>1] >> 1) | ((i&1) ? m : 0);
        NTT(c, 1, l), NTT(inv, 1, l);
        for(register int i = 0; i < l; ++i) c[i] = 1ll * c[i] * inv[i] % mod;
        NTT(c, -1, l);
        int ni = qpow(l, mod - 2);
        for(register int i = 0; i < l; ++i) c[i] = 1ll * c[i] * ni % mod;
        for(register int i = 0; i < m; ++i) b[i] = 1ll * ni2 * ((b[i] + c[i])%mod) % mod;
    }
}
int main() {
    scanf("%d", &n);
    for(register int i = 0; i < n; ++i) scanf("%d", &a[i]);
    kaigen(a, b, n);
    for(register int i = 0; i < n; ++i) printf("%d ", b[i]);
    return 0;
}

多項式求導與積分

求導 : 多項式求導的話一定會使它的次數降低。
多項式求導可以近似理解為對於多項式的每一項分別求導,冪函數的求導公式就是 \(f'(x^n) = n* x^{n-1}\)
所以多項式求導的代碼也就大致出來了:

void qiudao(int *a, int *b, int n) {
     for(register int i = 1; i < n; ++i) b[i-1] = 1ll * i * a[i] % mod; b[n-1] = 0;
}

積分
反正我只知道積分可以大概理解為是求導的逆運算,而且我也不會算,所以直接告訴你代碼吧(其實可以從求逆的角度理解):

void jifen(int *a, int *b, int n) {
     for(register int i = 1; i < n; ++i) b[i] = 1ll * a[i-1] * qpow(i, mod - 2) % mod; b[0] = 0;
}

注意積分不能把常數項求出來,所以我們默認為0。

多項式對數函數(多項式 ln)

已知多項式 \(A(x)\),我們需要求出一個多項式使其滿足 \(B(x) = lnA(x)\)
下面開始推式子,首先我們對兩邊分別求導:

\[B(x) = lnA(x) \]

\[B'(x) = \frac{A'(x)}{A(x)} \]

右邊的式子我們可以理解為復合函數的求導展開。
然后我們發現右邊的式子很可做。
然后我們就可以得到 \(B'(x)\),最后在積分一下就得到了 \(B(x)\)

代碼
#include <bits/stdc++.h>
using namespace std;
const int N = 4e5 + 5, mod = 998244353, g1 = 3, g2 = (mod + 1) / 3;
int n, a[N], b[N], A[N], B[N], c[N], inv[N], d[N];
int qpow(int a, int b, int res = 1) {
    for(; b; b >>= 1, a = 1ll * a * a % mod) 
        if(b & 1) res = 1ll * res * a % mod;
    return res;
}
void NTT(int *a, int flag, int len) {
    for(register int i = 0; i < len; ++i) if(i < d[i]) swap(a[i], a[d[i]]);
    for(register int l = 2, m = 1; l <= len; l <<= 1, m <<= 1) {
        int w1 = qpow(flag == 1 ? g1 : g2, (mod-1) / l);
        for(register int *p = a; p != a + len; p += l) {
            int w = 1;
            for(register int i = 0; i < m; ++i, w = 1ll*w*w1%mod) {
                int t = 1ll * w * p[i + m] % mod;
                p[i + m] = (p[i] - t + mod) % mod;
                p[i] = (p[i] + t) % mod;
            }
        }
    }
}
void qiuni(int *a, int *b, int n) {
    memset(A, 0, n * 16), memset(B, 0, n * 16);
    b[0] = qpow(a[0], mod - 2);
    for(register int l = 2, m = 1; l < n << 2; l <<= 1, m <<= 1) {
        for(register int i = 0; i < m; ++i) A[i] = a[i], B[i] = b[i];
        for(register int i = 0; i < l; ++i) d[i] = (d[i>>1] >> 1) | ((i&1) ? m : 0);
        NTT(A, 1, l), NTT(B, 1, l);
        for(register int i = 0; i < l; ++i) b[i] = ((2ll - 1ll*A[i]*B[i]%mod) * B[i] % mod + mod) % mod;
        NTT(b, -1, l);
        int ni = qpow(l, mod - 2);
        for(register int i = 0; i < l; ++i) b[i] = 1ll * b[i] * ni % mod;
        for(register int i = m; i < l; ++i) b[i] = 0;
    }
}
void qiudao(int *a, int *b, int n) {
    for(register int i = 1; i < n; ++i) b[i-1] = 1ll * i * a[i] % mod; b[n-1] = 0;
}
void jifen(int *a, int *b, int n) {
    for(register int i = 1; i < n; ++i) b[i] = 1ll * a[i-1] * qpow(i, mod - 2) % mod; b[0] = 0;
}
void loyin(int *a, int *b, int n) {
    memset(c, 0, n * 16), memset(inv, 0, n * 16);
    qiudao(a, c, n), qiuni(a, inv, n);
    int len = 1;
    while(len <= n*2) len <<= 1;
    for(register int i = 0; i < len; ++i) d[i] = (d[i>>1] >> 1) | ((i&1) ? len >> 1 : 0);
    NTT(c, 1, len), NTT(inv, 1, len);
    for(register int i = 0; i < len; ++i) c[i] = 1ll * c[i] * inv[i] % mod;
    NTT(c, -1, len);
    int ni = qpow(len, mod - 2);
    for(register int i = 0; i < len; ++i) c[i] = 1ll * c[i] * ni % mod;
    jifen(c, b, n);
}
int main() {
    scanf("%d", &n);
    for(register int i = 0; i < n; ++i) scanf("%d", &a[i]);
    loyin(a, b, n);
    for(register int i = 0; i < n; ++i) printf("%d ", b[i]); 
    return 0;
}

分治FFT

板子題鏈接

分治FFT主要是解決這樣一類自己會轉移到自己的問題,比如給定多項式 \(A(x)\),讓你求出 \(B_i = \sum^{i}_{j = 1} B_{i-j} * A_j\),邊界 \(B_0 = 1\)

這就是分治FFT主要解決的問題,主要的思想是利用CDQ分治的思想,每次計算左半變對於右半邊的共享,這樣復雜度是 \(O(nlog^2n)\)的。

其實其主要思想就是cdq分治,每次做一個卷積(乘法),好像不知道怎么解釋啊,自己試一下,當然主要的細節問題是邊界問題,就是每次做多項式乘法的邊界問題,這個在代碼里標注。

代碼
#include <bits/stdc++.h>
using namespace std;
const int N = 4e5 + 5, mod = 998244353, g1 = 3, g2 = (mod + 1) / 3;
int n, a[N], b[N], d[N], cdq1[N], cdq2[N];
int qpow(int a, int b, int res = 1) {
    for(; b; b >>= 1, a = 1ll * a * a % mod) 
        if(b & 1) res = 1ll * res * a % mod;
    return res;
}
void NTT(int *a, int flag, int len) {
    for(register int i = 0; i < len; ++i) if(i < d[i]) swap(a[i], a[d[i]]);
    for(register int l = 2, m = 1; l <= len; l <<= 1, m <<= 1) {
        int w1 = qpow(flag == 1 ? g1 : g2, (mod-1) / l);
        for(register int *p = a; p != a + len; p += l) {
            int w = 1;
            for(register int i = 0; i < m; ++i, w = 1ll*w*w1%mod) {
                int t = 1ll * w * p[i + m] % mod;
                p[i + m] = (p[i] - t + mod) % mod;
                p[i] = (p[i] + t) % mod;
            }
        }
    }
}
void fenzhi(int l, int r) {
    if(l == r) return;
    int mid = (l + r) / 2, len = r - l + 1;
    fenzhi(l, mid);
    int *A = cdq1, *B = cdq2;
    memset(A, 0, len * 8), memset(B, 0, len * 8);//注意每次清空數組
    for(register int i = l; i <= mid; ++i) A[i-l] = b[i];
    for(register int i = 0; i < r-l+1; ++i) B[i] = a[i+1];
    int n = 1;
    while(n <= len) n <<= 1;//之所以不用到 n <= len << 1,是因為這是一個循環卷積。
    for(register int i = 0; i < n; ++i) d[i] = (d[i>>1] >> 1) | ((i&1) ? n >> 1 : 0);
    NTT(A, 1, n), NTT(B, 1, n);
    for(register int i = 0; i < n; ++i) A[i] = 1ll * A[i] * B[i] % mod;
    NTT(A, -1, n);
    int ni = qpow(n, mod - 2);
    for(register int i = 0; i < n; ++i) A[i] = 1ll * A[i] * ni % mod;
    for(register int i = mid + 1; i <= r; ++i) b[i] = (b[i] + A[i-l-1]) % mod;//注意取值與邊界
    fenzhi(mid + 1, r);
}
int main() {
    scanf("%d", &n), b[0] = 1;
    for(register int i = 1; i < n; ++i) scanf("%d", &a[i]);
    fenzhi(0, n-1);
    for(register int i = 0; i < n; ++i) printf("%d ", b[i]); 
    return 0;
}

多項式指數函數(多項式 exp)

我為什么要先講一個分治FFT,主要是因為我的exp是寫的分治FFT,復雜度是 \(O(nlog^2n)\),雖然復雜度比正解的 \(O(nlogn)\) 要慢一點,但是代碼好寫,比正解的代碼碼量少一半,而且因為正解的常數過於巨大,如果卡常不好的話,我這個更快,洛谷上的提交記錄顯示,我的代碼是比大多數人的正解都快的,所以我選擇寫分治FFT。

顯然 exp 是給你一個多項式 \(A(x)\),讓你求得一個多項式 \(B(x)\),使其滿足 \(e^{A(x)} = B(x)\)

顯然我們接着推式子:
我們依舊對於左右兩邊求導:

\[B'(x) = A'(x) * e^{A(x)} \]

替換得到 :

\[B'(x) = A'(x) * B(x) \]

我們直接積分求解得:

\[B(x) = \int A'(x) * B(x) \]

然后我們可能覺得右邊的式子沒有什么思路,但是我們發現這個式子與分治FFT非常像,只不過是帶了一個積分,但是多項式的積分公式也很好求,所以我們可以選擇采用分治FFT的思路,注釋在代碼里。

代碼
#include <bits/stdc++.h>
using namespace std;
const int N = 4e5 + 5, mod = 998244353, g1 = 3, g2 = (mod + 1) / 3;
int n, a[N], b[N], d[N], cdq1[N], cdq2[N];
int qpow(int a, int b, int res = 1) {
    for(; b; b >>= 1, a = 1ll * a * a % mod) 
        if(b & 1) res = 1ll * res * a % mod;
    return res;
}
void NTT(int *a, int flag, int len) {
    for(register int i = 0; i < len; ++i) if(i < d[i]) swap(a[i], a[d[i]]);
    for(register int l = 2, m = 1; l <= len; l <<= 1, m <<= 1) {
        int w1 = qpow(flag == 1 ? g1 : g2, (mod-1) / l);
        for(register int *p = a; p != a + len; p += l) {
            int w = 1;
            for(register int i = 0; i < m; ++i, w = 1ll*w*w1%mod) {
                int t = 1ll * w * p[i + m] % mod;
                p[i + m] = (p[i] - t + mod) % mod;
                p[i] = (p[i] + t) % mod;
            }
        }
    }
}
void fenzhi(int l, int r) {
    if(l == r) {
        if(l) b[l] = 1ll * b[l] * qpow(l, mod - 2) % mod;//積分公式,乘一個逆元
        else b[l] = 1;
        return;
    }
    int mid = (l + r) / 2, len = r - l + 1;
    fenzhi(l, mid);
    int *A = cdq1, *B = cdq2;
    int n = 1;
    while(n < len) n <<= 1;
    memset(A, 0, n * 4), memset(B, 0, n * 4);//其他的與分治FFT無任何區別
    for(register int i = 0; i < n; ++i) d[i] = (d[i>>1] >> 1) | ((i&1) ? n >> 1 : 0);
    for(register int i = l; i <= mid; ++i) A[i-l] = b[i];
    for(register int i = 0; i < len; ++i) B[i] = a[i];
    NTT(A, 1, n), NTT(B, 1, n);
    for(register int i = 0; i < n; ++i) A[i] = 1ll * A[i] * B[i] % mod;
    NTT(A, -1, n);
    int ni = qpow(n, mod - 2);
    for(register int i = 0; i < n; ++i) A[i] = 1ll * A[i] * ni % mod;
    for(register int i = mid+1; i <= r; ++i) b[i] = (b[i] + A[i-l-1]) % mod;
    fenzhi(mid + 1, r);
}
int main() {
    scanf("%d", &n);
    for(register int i = 0; i < n; ++i) scanf("%d", &a[i]);
    for(register int i = 1; i < n; ++i) a[i-1] = 1ll * a[i] * i % mod; a[n-1] = 0;//求導
    fenzhi(0, n-1);
    for(register int i = 0; i < n; ++i) printf("%d ", b[i]);
    return 0;
}


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM