多项式学习笔记

多项式学习笔记…

多项式求逆

给定一个多项式 $g(x)$ ，求一个多项式 $f(x)$ 满足：

$$f(x)\ast g(x)\equiv1 \pmod {x^n}$$

也记作 $f(x)\equiv g^{-1}(x)\pmod{x^n}$。

令 $f_\ast(x)\equiv g^{-1}(x)\pmod {x^{n/2}}$：

$$\begin{aligned} f_\ast(x)&\equiv f(x) &\pmod{x^{n/2}}\\ f_\ast(x)- f(x)&\equiv0 &\pmod{x^{n/2}}\\ \big(f_\ast(x)- f(x)\big)^2&\equiv0 &\pmod{x^{n}}\\ f_\ast^2(x)-2f_\ast(x) f(x)+f^2(x)&\equiv0 &\pmod{x^{n}}\\ f_\ast^2(x)g(x)-2f_\ast(x)+f(x)&\equiv0 &\pmod{x^{n}}\\ \end{aligned}$$

$$f(x)\equiv f_\ast(x)\big(2-f_\ast(x)g(x)\big)\pmod{x^n}$$

倍增即可。
时间复杂度 $O(n\log n)。$

多项式开方

给定一个多项式 $g(x)$ ，求一个多项式 $f(x)$ 满足：

$$f^2(x)\equiv g(x) \pmod {x^n}$$

令 $f_\ast^2(x)\equiv g(x)\pmod {x^{n/2}}$：

$$\begin{aligned} f^2_{\ast}(x)&\equiv g(x) &\pmod{x^{n/2}}\\ f^2_{\ast}(x)-g(x)&\equiv 0 &\pmod{x^{n/2}}\\ \big(f^2_{\ast}(x)-g(x)\big)^2&\equiv 0 &\pmod{x^{n}}\\ \big(f^2_{\ast}(x)+g(x)\big)^2&\equiv 4f_\ast^2(x)g(x) &\pmod{x^{n}}\\ \left(\frac{f^2_{\ast}(x)+g(x)}{2f_\ast(x)}\right)^2&\equiv g(x) &\pmod{x^{n}}\\ \frac{f^2_{\ast}(x)+g(x)}{2f_\ast(x)}&\equiv f(x) &\pmod{x^{n}}\\ \end{aligned}$$

$$f(x)\equiv 2^{-1}f_\ast(x)+2^{-1}g(x)f_\ast^{-1}(x)\pmod {x^n}$$

倍增即可。
时间复杂度 $O(n\log n)$。

多项式带余除法

给定一个 $n$ 次多项式 $f(x)$，一个次数为 $m$ 的多项式 $g(x)$ ，求出多项式 $Q(x),R(x)$ 满足：

$$\deg Q=n-m,\deg R <m$$

$$f(x)=Q(x)\ast g(x)+R(x)$$

即 $Q(x)$ 为商，$R(x)$ 为余数。

定义一种变换（这里的 $R$ 更 $R(x)$ 没关系）：

$$f^R(x)=f\left(\frac{1}{x}\right)x^{\deg f}$$

其实就是把系数 reverse 一下。
考虑把几个多项式都搞成变换过的形式：

$$\begin{aligned} f(x)&=Q(x)\ast g(x)+R(x)\\[1ex] f\left(\frac{1}{x}\right)&=Q\left(\frac{1}{x}\right)\ast g\left(\frac{1}{x}\right)+R\left(\frac{1}{x}\right)\\ x^n f\left(\frac{1}{x}\right)&=x^{n-m}Q\left(\frac{1}{x}\right)\ast x^mg\left(\frac{1}{x}\right)+x^{n-m+1}x^{m-1}\left(\frac{1}{x}\right)\\[1ex] f^R(x)&=Q^R(x)\ast g^R(x)+x^{n-m+1}R^R(x)\\ \end{aligned}$$

放在模 $x^{n-m+1}$ 在 $R(x)$ 就没了，然后 $\deg Q = n-m < n-m+1$ 不会受到影响。

$$f^R(x)=Q^R(x)\ast g^R(x) \pmod{x^{n-m+1}}$$

多项式求逆可求 $Q(x)$ ，再求 $R(x)=f(x)-Q(x)\ast g(x)$ 即可。

多项式求导/积分

对于多项式 $f(x)$ 。

求导：

$$f'(x)=\sum_{i=0} [x^i]f(x) \cdot ix^{i-1}$$

$$[x^n]g(x)=(n+1) [x^{n+1}]f(x)$$

积分：

$$\int f(x)\,\mathrm{d}x = C+\sum_{i=0}[x^i]f(x)\cdot \frac{1}{i+1}x^{i+1}$$

$$[x^n]g(x)=\frac{1}{n}[x^{n-1}]f(x)$$

复杂度都是 $O(n)$ 的。

多项式牛顿迭代

给定一个多项式 $g(x)$ ，求模 $x^n$ 意义下的 $f(x)$ 满足：

$$g\big(f(x)\big)\equiv 0\pmod{x^n}$$

结论：

$$f(x)\equiv f_\ast(x)-\frac{g\big(f_\ast(x)\big)}{g'\big(f_\ast(x)\big)}\pmod{x^n}$$

$f_\ast (x)$ 是已经找到的在模 $x^{n/2}$ 意义下的解，上式是一个倍增。

干啥的，方便推导一些东西。

求逆

$$f(x)^{-1}\equiv h(x)\pmod{x^n}$$

设 $g\big(f(x)\big)=f^{-1}(x)-h(x)$ ，有方程 $g\big(f(x)\big)\equiv 0 \pmod{x^n}$ 。牛顿迭代一下：

$$\begin{aligned} f(x)&\equiv f_\ast (x) -\frac{f_\ast^{-1}(x)-h(x)}{-f_\ast^{-2}(x)}&\pmod{x^n}\\ f(x)&\equiv 2f_\ast(x)-f_\ast^2(x)h(x) &\pmod{x^n} \end{aligned}$$

开方

$$f^2(x)\equiv h(x)\pmod{x^n}$$

设 $g\big(f(x)\big)=f^{2}(x)-h(x)$ ，有方程 $g\big(f(x)\big)\equiv 0 \pmod{x^n}$ 。牛顿迭代一下：

$$\begin{aligned} f(x)&\equiv f_\ast (x)-\frac{f_\ast^2(x)-h(x)}{2 f_\ast(x)} &\pmod{x^n}\\ f(x)&\equiv 2^{-1}f_\ast(x)+2^{-1}h(x)f_\ast^{-1}(x)&\pmod {x^n}\\ \end{aligned}$$

EXP

$$\ln f(x)\equiv h(x) \pmod{x^n}$$

设 $g\big(f(x)\big)=\ln f(x)-h(x)$ ，有方程 $g\big(f(x)\big)\equiv 0 \pmod{x^n}$ 。牛顿迭代一下：

$$\begin{aligned} f(x)&\equiv f_\ast(x)-\frac{\ln f_\ast(x)-h(x)}{f^{-1}_\ast(x)} \pmod{x^n}\\ f(x)&\equiv f_\ast(x)\big(1-\ln f_\ast (x)+h(x)\big) \pmod{x^n}\\ \end{aligned}$$

---

以上推倒的式子复杂度都为：

$$T(n)=T\left(\frac{n}{2}\right)+O(n\log n)=O(n\log n)$$

话说这里推导的求逆和开方和之前的是一样的，但显然用牛顿迭代推导更方便。

多项式对数函数

给定一个多项式 $g(x)$ ，求一个多项式 $f(x)$ 满足：

$$f(x)\equiv \ln g(x) \pmod{x^n}$$

具体的由麦克劳林级数定义：

$$\ln f(x)=-\sum_{i=1}^{+\infty} \frac{\big(1-f(x)\big)^i}{i}$$

若 $\ln g(x)$ 存在需满足 $[x^0]g(x)=1$
有 $\ln' x=\frac{1}{x}$，两边同时求导，注意链式法则。再积分就完事了。

$$\begin{aligned} f(x)&\equiv \ln g(x) &\pmod{x^n}\\ f'(x)&\equiv \frac{g'(x)}{g(x)} &\pmod{x^n}\\ f(x)&\equiv \int \frac{g'(x)}{g(x)}\;\mathrm{d} x &\pmod{x^n}\\ \end{aligned}$$

一次求导，一次求逆，一次乘法，一次积分完事了。
复杂度 $O(n \log n)$ 。

递推法：

$$\begin{aligned} f(x)&\equiv \ln g(x) &\pmod{x^n}\\ f'(x)g(x)&\equiv g'(x) &\pmod{x^n}\\ \end{aligned}$$

$$\begin{aligned}{} \sum_{i=0}^n [x^i]f'(x)[x^{n-i}]g(x) &= [x^n]g'(x)\\ \sum_{i=0}^n (i+1)[x^{i+1}]f'(x)[x^{n-i}]g(x) &= n[x^{n+1}]g'(x) \end{aligned}$$

$$[x^n]f(x)=[x^n]g(x)-\frac{1}{n}\sum_{i=1}^{n-1} i[x^{n-i}]g(x)[x^i]f(x)$$

多项式指数函数

给定一个多项式 $g(x)$ ，求一个多项式 $f(x)$ 满足：

$$f(x)\equiv \exp g(x) \pmod{x^n}$$

具体的由麦克劳林级数定义：

$$\exp f(x)=\sum_{i=0}^{+\infty} \frac{f^i(x)}{i!}$$

若 $\exp g(x)$ （即 $e^{g(x)}$）存在需满足 $[x^0]g(x)=0$。

注：$\exp' x= \exp x$

通过牛顿迭代给出了一种 $O(n\log n)$ 的倍增做法: EXP 。再来考虑一下普通做法：

$$\begin{aligned} f(x)&\equiv \exp g(x) &\pmod{x^n}\\ f'(x)&\equiv f(x)g'(x)&\pmod{x^n} \end{aligned}$$

提取系数：

$$\begin{aligned}{} [x^n]f'(x)&=\sum_{i=0}^n[x^{n-i}]f(x)\cdot [x^{i}]g'(x)\\ (n+1) [x^{n+1}]f(x)&=\sum_{i=0}^n [x^{n-i}]f(x)\cdot (i+1) [x^{i+1}]g(x)\\ [x^n]f(x)&=\frac{1}{n}\sum_{i=1}^{n} [x^{n-i}]f(x)\cdot i[x^{i}]g(x) \end{aligned}$$

注：$i=1$ 枚举因为 $[x^0]g(x)=0$

分治 FFT ，复杂度 $O(n\log^2 n)$ ，实际上常数小，比倍增的快。

多项式快速幂

给定一个多项式 $g(x)$ ，求一个多项式 $f(x)$ 满足：

$$f(x)\equiv g^k(x) \pmod{x^n}$$

如果只是 niave 多项式乘法+快速幂，复杂度 $O(n\log n \log k)$。

显然有更好的做法：

$$f(x)=g^k(x)=e^{(\ln g(x)) k}$$

会 EXP,Ln 便完事了，复杂度 $O(n\log n)$ ，但需要满足 $[x^0]g(x)=1$ 。

注意复杂度与 $k$ 无关，$k$ 只要 $\bmod 998244353$ 即可。（$k$ 跟 $\ln g(x)$ 是数乘）

代码实现

/*
    NTT 多项式全家桶
    p_mul 乘法; p_inv 求逆; p_div 带余数除法; p_sqrt 开方;
    p_ln Ln； p_exp EXP; p_int 积分; p_der 求导; p_pow 快速幂;
    DCFFT 分治 FFT 板子;

    to be continue ...
    多项式三角函数，多项式反三角函数，多项式多点求值，多项式快速差值.....
*/

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#define clr(f,n) memset(f,0,sizeof(long long)*(n))
#define cpy(f,g,n) memcpy(f,g,sizeof(long long)*(n))
#define Outarr(x,n) cerr<<#x<<" : "; for(int i=0;i<n;++i) cerr<<x[i]<<" ";cout<<endl;
#define outarr(x,n) for(int i=0;i<n;++i) printf("%lld%c",x[i],(i==n-1)?'\n':' ');
#define MOD(x) ((x)<mod?(x):((x)%mod))
using namespace std;
template<typename _Tp>
inline void red(_Tp &x) {
    x = 0; bool fg = 0; char ch = getchar();
    while (ch < '0' || ch > '9') {if (ch == '-') fg ^= 1; ch = getchar();}
    while (ch >= '0' && ch <= '9') x = (x << 1) + (x << 3) + (_Tp)(ch ^ 48), ch = getchar();
    if (fg) x = -x;
}
typedef long long ll;

namespace poly {
    const int mod = 998244353;
    const int N = (1 << 19);
    const int _G = 3;
    const int _iG = 332748118;
    const int inv2 = 499122177;

    ll fpow(ll a, ll b, ll p) {
        ll r = 1;
        for (; b; a = a * a % p, b >>= 1) if (b & 1) r = r * a % p;
        return r;
    }

    int rev[N], rev_n;
    void prerev(int n) {
        if (n == rev_n) return;
        rev_n = n;
        for (int i = 0; i < n; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (n >> 1) : 0);
    }
    // NTT : fg=1 DFT fg=-1 IDFT
    void NTT(ll f[], int n, int fg) {
        prerev(n);
        for (int i = 0; i < n; ++i) if (i < rev[i]) swap(f[i], f[rev[i]]);
        for (int h = 2; h <= n; h <<= 1) {
            ll Dt = fpow((fg == 1) ? _G : _iG, (mod - 1) / h, mod), w;
            int len = h >> 1;
            for (int j = 0; j < n; j += h) {
                w = 1;
                for (int k = j; k < j + len; ++k) {
                    ll tmp = MOD(f[k + len] * w);
                    f[k + len] = f[k] - tmp; (f[k + len] < 0) &&(f[k + len] += mod);
                    f[k] = f[k] + tmp; (f[k] >= mod) &&(f[k] -= mod);
                    w = MOD(w * Dt);
                }
            }
        }
        if (fg == -1) {
            ll invn = fpow(n, mod - 2, mod);
            for (int i = 0; i < n; ++i) f[i] = MOD(f[i] * invn);
        }
    }
    // f(x) = f*g(x) n = def f ; m = def g ; len = 最终长度(保留几位)
    void p_mul(ll f[], ll g[], int n, int m, int len) {
        static ll a[N], b[N];
        int nn = 1 << (int)ceil(log2(n + m - 1));
        clr(a, nn); clr(b, nn); cpy(a, f, n); cpy(b, g, m);
        NTT(a, nn, 1); NTT(b, nn, 1);
        for (int i = 0; i < nn; ++i) a[i] = MOD(a[i] * b[i]);
        NTT(a, nn, -1);
        for (int i = 0; i < len; ++i) f[i] = a[i];
    }
    // f(x) = g^-1(x)  f(x) 为 g(x) 模 x^n 意义下的逆
    void p_inv(ll g[], int n, ll f[]) {
        static ll sav[N];
        int nn = 1 << (int)ceil(log2(n));
        clr(f, n * 2);
        f[0] = fpow(g[0], mod - 2, mod);
        for (int h = 2; h <= nn; h <<= 1) {
            cpy(sav, g, h); clr(sav + h, h);
            NTT(sav, h << 1, 1); NTT(f, h << 1, 1);
            for (int i = 0; i < (h << 1); ++i) {
                f[i] = f[i] * (2ll - f[i] * sav[i] % mod + mod) % mod;
            }
            NTT(f, h << 1, -1); clr(f + h, h);
        }
        clr(f + n, nn * 2 - n);
    }
    // f^2(x) = g(x)  f(x) 为 g(x) 模 x^n 意义下的开方
    void p_sqrt(ll g[], int n, ll f[]) {
        static ll sav[N], r[N];
        int nn = 1 << (int)ceil(log2(n));
        clr(f, n * 2); f[0] = 1; // g[0] should be 1 otherwise
        for (int h = 2; h <= nn; h <<= 1) {
            cpy(sav, g, h); clr(sav + h, h); p_inv(f, h, r);
            NTT(sav, h << 1, 1); NTT(r, h << 1, 1);
            for (int i = 0; i < (h << 1); ++i) sav[i] = MOD(sav[i] * r[i]);
            NTT(sav, h << 1, -1);
            for (int i = 0; i < h; ++i) f[i] = MOD((f[i] + sav[i]) * inv2);
            clr(f + h, h);
        }
        clr(f + n, nn * 2 - n);
    }
    // f(x) = g(x) * q(x) + r(x) : q(x) 为商 r(x) 为余数
    void p_div(ll f[], ll g[], int n, int m, ll q[], ll r[]) {
        static ll sav1[N], sav2[N];
        int nn;
        for (nn = 1; nn < n - m + 1; nn <<= 1);
        clr(sav1, nn); clr(sav2, nn); cpy(sav1, f, n); cpy(sav2, g, m);
        reverse(sav1, sav1 + n); reverse(sav2, sav2 + m);
        p_inv(sav2, n - m + 1, q); p_mul(q, sav1, n - m + 1, n, n - m + 1);
        reverse(q, q + n - m + 1); cpy(r, g, m);
        p_mul(r, q, m, n - m + 1, m - 1);
        for (int i = 0; i < m - 1; ++i) r[i] = MOD(f[i] - r[i] + mod);
    }
    // 预处理乘法逆元
    ll inv[N];
    void Initinv(int n) {
        inv[0] = inv[1] = 1;
        for (int i = 2; i <= n; ++i) inv[i] = (mod - mod / i) * inv[mod % i] % mod;
    }
    // 对 f(x) 进行积分  Initinv() first
    void p_int(ll f[], int n) {
        for (int i = n - 1; i; --i) f[i] = MOD(f[i - 1] * inv[i]);
        f[0] = 0;
    }
    // 对 f(x) 进行求导
    void p_der(ll f[], int n) {
        for (int i = 1; i < n; ++i) f[i - 1] = MOD(f[i] * i);
        f[n - 1] = 0;
    }
    //  f(x) <- ln f(x)  f[0] should be 1
    void p_ln(ll f[], int n) {
        static ll g[N];
        p_inv(f, n, g); p_der(f, n);
        p_mul(f, g, n, n, n + n);
        p_int(f, n);
    }

    // f(x) <- exp f(x) (倍增版) f[0] should be 0
    void p_exp(ll f[], int n) {
        static ll g[N], sav[N];
        clr(g, n * 2); clr(sav, n * 2); g[0] = 1;
        for (int h = 2; h <= n; h <<= 1) {
            cpy(sav, g, h); p_ln(sav, h);
            for (int i = 0; i < h; ++i) sav[i] = MOD(f[i] - sav[i] + mod);
            sav[0] = MOD(sav[0] + 1);
            p_mul(g, sav, h, h, h);
        }
        cpy(f, g, n);
    }

    /*
    void _p_exp(ll f[],ll g[],int l,int r) {
        static ll A[N],B[N];
        if(r-l==1) {if(l>0) f[l]=MOD(f[l]*fpow(l,mod-2,mod));return ;}
        int mid=(l+r)>>1,len=mid-l;
        _p_exp(f,g,l,mid);
        cpy(A,f+l,len); clr(A+len,len); cpy(B,g,len<<1);
        p_mul(A,B,len<<1,len<<1,len<<1);
        for(int i=mid;i<r;++i) f[i]=MOD(f[i]+A[i-l]);
        _p_exp(f,g,mid,r);
    }
    // f(x) <- exp f(x) (分治 FFT 版) f[0] should be 0
    void p_exp(ll f[],int n) {
        static ll g[N];
        cpy(g,f,n); clr(f,n); f[0]=1;
        for(int i=0;i<n;++i) g[i]=MOD(g[i]*i);
        _p_exp(f,g,0,n);
    }
    */

    // f(x) <- f^k(x)  f(x) 模 x^n 意义下的 k 次
    void p_pow(ll f[], int n, ll k) {
        p_ln(f, n);
        for (int i = 0; i < n; ++i) f[i] = MOD(f[i] * k);
        p_exp(f, n);
    }

    // 分治FFT [l,r)   F[n] = sum(0<i<=n) F[n-i]G[i]
    void DCFFT(ll f[], ll g[], int l, int r) {
        static ll A[N], B[N];
        if (r - l == 1) return ;
        int mid = (l + r) >> 1, len = mid - l;
        DCFFT(f, g, l, mid);
        cpy(A, f + l, len); clr(A + len, len); cpy(B, g, len << 1);
        p_mul(A, B, len << 1, len << 1, len << 1);
        for (int i = mid; i < r; ++i) f[i] = MOD(f[i] + A[i - l]);
        DCFFT(f, g, mid, r);
    }
}