• 多项式的各类计算(多项式的逆/开根/对数/exp/带余除法/多点求值)


    • 预备知识:FFT/NTT

    • 多项式的逆

      给定一个多项式 F(x)F(x),请求出一个多项式 G(x)G(x),满足 F(x)G(x)1(mod xn)F(x)*G(x) equiv 1(mod x^n)
      系数对 998244353998244353 取模,1n1051≤n≤10^5

      首先将多项式的长度拓展至22的次幂,然后我们要求的是
      G(x)F(x)1 (mod xn)G(x)*F(x) equiv 1 (mod x^n)假设已经求出了
      H(x)F(x)1 (mod xn2)H(x)*F(x) equiv 1 (mod x^{⌈frac n2⌉})又因为有
      G(x)F(x)1 (mod xn2)G(x)*F(x) equiv 1 (mod x^{⌈frac n2⌉})两式相减有
      (G(x)H(x))F(x)0 (mod xn2)(G(x)-H(x))*F(x) equiv 0 (mod x^{⌈frac n2⌉})G(x)H(x)0 (mod xn2)G(x)-H(x) equiv 0 (mod x^{⌈frac n2⌉})因为左边得到的多项式前n2{⌈frac n2⌉}项都是00,平方前nn项都是00,所以有
      (G(x)H(x))20 (mod xn)(G(x)-H(x))^2 equiv 0 (mod x^{n})G2(x)2G(x)H(x)+H2(x)0 (mod xn)G^2(x)-2G(x)H(x)+H^2(x)equiv 0 (mod x^n)两边同时乘以F(x)F(x)
      G(x)2H(x)+H2(x)F(x)0 (mod xn)G(x)-2H(x)+H^2(x)F(x)equiv 0 (mod x^n)G(x)2H(x)H2(x)F(x) (mod xn)G(x)equiv2H(x)-H^2(x)F(x) (mod x^n)
      我们就得到了递推式,时间复杂度如下T(n)=T(n/2)+O(nlog2n)=O(nlog2n)T(n)=T(n/2)+O(nlog_2n)=O(nlog_2n)

    • 多项式开根

      给定一个 n1n−1 次多项式 F(x)F(x),求一个在 mod xnmod x^n 意义下的多项式 G(x)G(x),使得 G2(x)F(x) (mod xn)G^2(x) equiv F(x) (mod x^n)(多项式的系数在 mod 998244353mod 998244353 的意义下进行运算,1n1051≤n≤10^5
      同样假设求出了H2(x)F(x) (mod xn2)H^2(x)equiv F(x) (mod x^{⌈frac n2⌉})且有
      G2(x)F(x) (mod xn2)G^2(x)equiv F(x) (mod x^{⌈frac n2⌉})所以
      G(x)H(x)0 (mod xn2)G(x)-H(x) equiv 0 (mod x^{⌈frac n2⌉})(G(x)H(x))20 (mod xn)(G(x)-H(x))^2 equiv 0 (mod x^{n})G2(x)2G(x)H(x)+H2(x)0 (mod xn)G^2(x)-2G(x)H(x)+H^2(x)equiv 0 (mod x^n)F(x)2G(x)H(x)+H2(x)0 (mod xn)F(x)-2G(x)H(x)+H^2(x)equiv 0 (mod x^n)
      所以2G(x)F(x)H(x)+H(x) (mod xn)2G(x)equiv frac{F(x)}{H(x)}+H(x) (mod x^n)
      此处除法转化为乘上多项式的逆,就能递推了,时间复杂度也是O(nlog2n)O(nlog_2n)

    • 多项式求自然对数

      给出 n1n−1 次多项式 F(x)F(x),求一个 mod xnmod x^n 下的多项式 G(x)G(x),满足 G(x)lnF(x) (mod xn)G(x) equiv ln F(x) (mod x^n).
      多项式的系数在 mod 998244353mod 998244353 的意义下进行运算,1n1051≤n≤10^5
      G(x)=ln F(x)=F(x)F(x)dx=F(x)F1(x)dxegin{aligned} G(x)&=ln F(x)\ &=int frac{F'(x)}{F(x)}dx\ &=int F'(x)F^{-1}(x)dx end{aligned}所以G=FF1G=int F'F^{-1},直接多项式求逆+多项式求导积分就行了
      求逆复杂度O(nlog2n)O(nlog_2n),求导积分复杂度O(n)O(n),总时间复杂度仍是O(nlog2n)O(nlog_2n)

    • 多项式exp

      多项式exp是用牛顿迭代做的,我也不会证明牛顿迭代,只会背公式.
      本来牛顿迭代公式是这个样子的xn+1=xnf(xn)f(xn)x_{n+1}=x_n-frac{f(x_n)}{f'(x_n)}
      对于多项式的牛顿迭代长这样un+1=unf(u,v)fu(u,v)u_{n+1}=u_n-frac{f(u,v)}{f_u(u,v)}

      • 其中f(u,v)f(u,v)表示已知多项式为vv,未知多项式为uu的多项式方程,就拿exp为例,原方程为uev (mod xn)uequiv e^v (mod x^n)两边取对数为ln uv (mod xn)ln uv0 (mod xn)egin{aligned}ln uequiv v (mod x^n)\ln u-vequiv 0 (mod x^n)end{aligned}那么f(u,v)=ln uvf(u,v)=ln u-v.
      • fu(u,v)f_u(u,v)表示的是f(u,v)f(u,v)uu取导.注意是对uu,不是对多项式里的xx.
        举个栗子,ln vln vvv求导是1vfrac 1v,但是对xx求导是vvfrac{v'}{v}(此处vv'是对xx求导)

      那么我们就可以推推柿子un+1un(ln unv(ln unv)) (mod xn)un(ln unv1un) (mod xn)un(ln unv)un (mod xn)un(1ln un+v) (mod xn)egin{aligned}u_{n+1}&equiv u_n-left(frac{ln u_n-v}{(ln u_n-v)'} ight) &(mod x^n)\&equiv u_n-left(frac{ln u_n-v}{frac1{u_n}} ight) &(mod x^n)\&equiv u_n-(ln u_n-v)u_n &(mod x^n)\&equiv u_n(1-ln u_n+v) &(mod x^n)end{aligned}然后就套用前面的模板迭代做了

      TIP:TIP:其实多项式开根也可以用这个方法推出来,有兴趣可以去试试看.

    代码
    #include <cmath>
    #include <cstdio>
    #include <cstring>
    #include <algorithm>
    using namespace std;
    namespace READ {
    	inline void read(int &num) {
    		char ch; int flg=1;
    		while((ch=getchar())<'0'||ch>'9')if(ch=='-')flg=-flg;
    		for(num=0;ch>='0'&&ch<='9';num=num*10+ch-'0',ch=getchar());
    		num*=flg;
    	}
    }
    using namespace READ;
    const int mod = 998244353, G = 3, inv2 = 499122177;
    const int MAXN = 1<<18; // 注意MAXN要大于2*n
    int rev[MAXN], Wn[MAXN+1][2], inv[MAXN+1];
    inline void reset(int len) { for(int i = 0; i < len; ++i) rev[i] = (rev[i>>1]>>1)|((i&1)?(len>>1):0); }
    inline int qmul(int a, int b) {
    	int res = 1;
    	while(b) {
    		if(b&1) res = 1ll * res * a % mod;
    		a = 1ll * a * a % mod; b>>=1;
    	}
    	return res;
    }
    inline int add(int x, int y) { return x+y>mod ? x+y-mod : x+y; }
    namespace Poly { 
    	int A[MAXN],B[MAXN],C[MAXN],D[MAXN],E[MAXN],F[MAXN],GG[MAXN],H[MAXN];//每个不同函数开不同的临时数组,否则可能重复使用
    	inline void clear(int *b, int n) { for(int i = 0; i < n; ++i) b[i] = 0; } //多清零
    	inline void NTT(int *arr, int n, int flg) {
    		int wn, w, A0, wA1, i, j, k;
    		for(i = 0; i < n; ++i) if(i < rev[i]) swap(arr[i], arr[rev[i]]);
    		for(i = 2; i <= n; i<<=1) {
    			wn = Wn[i][(~flg)?0:1];
    			for(j = 0, w = 1; j < n; j += i, w = 1)
    				for(k = j; k < j + i/2; ++k, w = 1ll * w * wn % mod) {
    					A0 = arr[k], wA1 = 1ll * w * arr[k + i/2] % mod;
    					arr[k] = add(A0, wA1);
    					arr[k + i/2] = add(A0, -wA1+mod);
    				}
    		}
    		if(!(~flg)) for(i = 0; i < n; ++i) arr[i] = 1ll * arr[i] * inv[n] % mod;
    	}
    	
    	inline void INV(int *a, int *b, int n) {
    		clear(b, n<<1), b[0] = qmul(a[0], mod-2);
    		clear(A, n<<1), clear(B, n<<1);
    		int len, lim;
    		for(len = 2; len < (n<<1); len<<=1) {
    			lim = len<<1;
    			for(int i = 0; i < len; ++i) A[i] = a[i], B[i] = b[i];
    			reset(lim); NTT(A, lim, 1); NTT(B, lim, 1);
    			for(int i = 0; i < lim; ++i) b[i] = ((2ll - 1ll * A[i] * B[i] % mod) * B[i] % mod + mod) % mod;
    			NTT(b, lim, -1);
    			for(int i = len; i < lim; ++i) b[i] = 0;
    			}
    		for(int i = 0; i < len; ++i) A[i] = B[i] = 0;
    		for(int i = n; i < len; ++i) b[i] = 0;
    	}
    	
    	inline void SQRT(int *a, int *b, int n) {
            clear(b, n<<1); b[0] = int(sqrt(a[0])+0.5);
            int len, lim, *A = GG, *B = H;
            clear(A, n<<1), clear(B, n<<1);
            for(len = 2; len < (n<<1); len<<=1) {
                lim = len<<1;
                for(int i = 0; i < len; ++i) A[i] = a[i];
                INV(b, B, len);
                reset(lim); NTT(A, lim, 1); NTT(B, lim, 1);
                for(int i = 0; i < lim; ++i) A[i] = 1ll * A[i] * B[i] % mod;
                NTT(A, lim, -1);
                for(int i = 0; i < len; ++i) b[i] = 1ll * (b[i] + A[i]) % mod * inv2 % mod;
                for(int i = len; i < lim; ++i) b[i] = 0;
            }
            for(int i = 0; i < len; ++i) A[i] = B[i] = 0;
            for(int i = n; i < len; ++i) b[i] = 0;
        }
    
    	inline void INT(int *a, int *b, int n) {
    		clear(b, n<<1);
    		for(int i = n-1; i; --i)
    			b[i] = 1ll * a[i-1] * qmul(i, mod-2) % mod;
    		b[0] = 0;
    	}
    	
    	inline void DIF(int *a, int *b, int n) {
    		clear(b, n<<1);
    		for(int i = 1; i < n; ++i)
    			b[i-1] = 1ll * a[i] * i % mod;
    		b[n-1] = 0;
    	}
    	
    	inline void LN(int *a, int *b, int n) {
    		int *A = E, *B = F, len = 1;
    		clear(A, n<<1), clear(B, n<<1);
    		INV(a, A, n); DIF(a, B, n);
    		while(len < (n<<1)) len<<=1;
    		reset(len); NTT(A, len, 1); NTT(B, len, 1);
    		for(int i = 0; i < len; ++i) A[i] = 1ll * A[i] * B[i] % mod;
    		NTT(A, len, -1); INT(A, b, n);
    	}
    	
    	inline void EXP(int *a, int *b, int n) {
    		clear(b, n<<1), b[0] = 1;
    		int len, lim, *A = C, *B = D;
    		clear(A, n<<1), clear(B, n<<1);
    		for(len = 2; len < (n<<1); len<<=1) {
    			lim = len<<1;
    			LN(b, A, len);
    			for(int i = 0; i < len; ++i) A[i] = ((i==0)-A[i]+a[i]+mod) % mod, B[i] = b[i];
    			for(int i = len; i < lim; ++i) A[i] = B[i] = 0;
    			reset(lim); NTT(A, lim, 1); NTT(B, lim, 1);
    			for(int i = 0; i < lim; ++i) A[i] = 1ll * A[i] * B[i] % mod;
    			NTT(A, lim, -1);
    			for(int i = 0; i < len; ++i) b[i] = A[i];
    		}
    		for(int i = n; i < len; ++i) b[i] = 0;
    	}
    }
    using namespace Poly;
    inline void Pre() { //预处理原根及逆元
    	for(int i = 1; i <= MAXN; i<<=1)
    		Wn[i][0] = qmul(G, (mod-1)/i),
    		Wn[i][1] = qmul(G, (mod-1)-(mod-1)/i),
    		inv[i] = qmul(i, mod-2);
    }
    int main () {
    	Pre();
    }
    

    Upd:Upd:现在看,发现我以前sb了…不用开不同的数组来取地址,直接static定义就行了…

    • 带余除法(待更…)

    • 多项式多点求值(待更…)

  • 相关阅读:
    [bzoj4653] [NOI2016]区间
    [bzoj5285] [HNOI2018]寻宝游戏
    [bzoj4071] [Apio2015]巴邻旁之桥
    [bzoj1146] [CTSC2008]网络管理Network
    [bzoj3004] [SDOi2012]吊灯
    [bzoj5321] [Jxoi2017]加法
    [bzoj5010] [FJOI2017]矩阵填数
    [bzoj3504] [CQOI2014]危桥
    ASP.NET
    ASP.NET
  • 原文地址:https://www.cnblogs.com/Orz-IE/p/12039439.html
Copyright © 2020-2023  润新知