我是靠谱客的博主 小巧帆布鞋,最近开发中收集的这篇文章主要介绍多项式各种操作,觉得挺不错的,现在分享给大家,希望可以做个参考。

概述

PS:洛谷上有大部分以下操作的模板题。

感谢 cyk 巨神 的 多项式算法合集。

多项式求逆

problem:

给定一个 ( n − 1 ) (n-1) (n1) 次多项式 A ( x ) A(x) A(x) ,求多项式 B ( x ) B(x) B(x) 满足:

A ( x ) B ( x ) ≡ 1 (   m o d    x n ) A(x)B(x)≡1(,mathrm{mod}; x^n) A(x)B(x)1(modxn)

solution:

为了方便,我们设 n = 2 k ( k > 0 ) n=2^k(k>0) n=2k(k>0)(若 n n n 不是 2 2 2 的整数次幂,可以像 f f t fft fft 那样强行扩大到 2 k 2^k 2k)。

用类似于倍增的方法,即假设我们已知

A ( x ) C ( x ) ≡ 1 (   m o d    x n 2 ) A(x)C(x)equiv1(,mathrm{mod};x^{frac n 2}) A(x)C(x)1(modx2n)

因为 A ( x ) B ( x ) ≡ 1 (   m o d    x n ) A(x)B(x)≡1(,mathrm{mod}; x^n) A(x)B(x)1(modxn),那么有 A ( x ) B ( x ) ≡ 1 (   m o d    x n 2 ) A(x)B(x)≡1(,mathrm{mod}; x^{frac n 2}) A(x)B(x)1(modx2n),两式相减得:

A ( x ) ( B ( x ) − C ( x ) ) ≡ 0 (   m o d    x n 2 ) A(x)(B(x)-C(x))equiv 0(,mathrm{mod};x^{frac n 2}) A(x)(B(x)C(x))0(modx2n)

那么消掉 A ( x ) A(x) A(x),两边平方得:

B 2 ( x ) + C 2 ( x ) − 2 B ( x ) C ( x ) ≡ 0 (   m o d    x n ) B^2(x)+C^2(x)-2B(x)C(x)equiv0(,mathrm{mod};x^n) B2(x)+C2(x)2B(x)C(x)0(modxn)

这个时候再把 A ( x ) A(x) A(x) 乘回去,得到:

A ( x ) ( B 2 ( x ) + C 2 ( x ) − 2 B ( x ) C ( x ) ) ≡ 0 (   m o d    x n ) A(x)(B^2(x)+C^2(x)-2B(x)C(x))equiv0(,mathrm{mod};x^n) A(x)(B2(x)+C2(x)2B(x)C(x))0(modxn)

然后再综合 A ( x ) B ( x ) ≡ 1 (   m o d    x n ) A(x)B(x)equiv1(,mathrm{mod};x^n) A(x)B(x)1(modxn) 化简,最终得到:

B ( x ) = 2 C ( x ) − A ( x ) C 2 ( x ) (   m o d    x n ) B(x)=2C(x)-A(x)C^2(x)(,mathrm{mod};x^n) B(x)=2C(x)A(x)C2(x)(modxn)

时间复杂度 O ( n log ⁡ n ) O(nlog n) O(nlogn)

code:

注意代码中用的变量和推导中用的不一样。

typedef vector<int> poly;
poly Inv(poly A,int len){
poly C,B(1,Power(A[0],P-2));
for(int limit=4;limit<(len<<2);limit<<=1){
init(limit);
C=A,C.resize(limit>>1);
C.resize(limit),NTT(C,limit,1);
B.resize(limit),NTT(B,limit,1);
for(int i=0;i<limit;++i)
B[i]=mul(B[i],dec(2,mul(B[i],C[i])));
NTT(B,limit,-1),B.resize(limit>>1);
}
B.resize(len);return B;
}

多项式开根

problem:

已知一个 ( n − 1 ) (n-1) (n1) 次多项式 A ( x ) A(x) A(x),求多项式 B ( x ) B(x) B(x) 满足

B 2 ( x ) ≡ A ( x ) (   m o d    x n ) B^2(x)≡A(x)(,mathrm{mod} ;x^n) B2(x)A(x)(modxn)

保证 A ( 0 ) = 1 A(0)=1 A(0)=1

solution:

和求逆一样,为了方便,记 n = 2 k ( k > 0 ) n=2^k(k>0) n=2k(k>0)

假设我们知道

C 2 ( x ) ≡ A ( x ) (   m o d    x n 2 ) C^2(x)equiv A(x)(,mathrm{mod} ;x^{frac n 2}) C2(x)A(x)(modx2n)

因为 B 2 ( x ) ≡ A ( x ) (   m o d    x n ) B^2(x)≡A(x)(,mathrm{mod}; x^n) B2(x)A(x)(modxn),那么有 B 2 ( x ) ≡ A ( x ) (   m o d    x n 2 ) B^2(x)≡A(x)(,mathrm{mod}; x^{frac n 2}) B2(x)A(x)(modx2n),那么

B ( x ) − C ( x ) ≡ 0 (   m o d    x n 2 ) B(x)-C(x)equiv0(,mathrm{mod} ;x^{frac n 2}) B(x)C(x)0(modx2n)

平方,得:

B 2 ( x ) + C 2 ( x ) − 2 B ( x ) C ( x ) ≡ 0 (   m o d    x n ) B^2(x)+C^2(x)-2B(x)C(x)equiv0(,mathrm{mod};x^n) B2(x)+C2(x)2B(x)C(x)0(modxn)

又由于 B 2 ( x ) ≡ A ( x ) (   m o d    x n ) B^2(x)equiv A(x)(,mathrm{mod};x^n) B2(x)A(x)(modxn),所以:

B ( x ) = A ( x ) 2 C ( x ) + C ( x ) 2 (   m o d    x n ) B(x)=frac{A(x)}{2C(x)}+frac{C(x)}{2}(,mathrm{mod};x^n) B(x)=2C(x)A(x)+2C(x)(modxn)

不难发现需要用到上面所说的多项式求逆

时间复杂度 O ( n log ⁡ n ) O(nlog n) O(nlogn)好像常数有点大

code:

typedef vector<int> poly;
poly Sqrt(poly A,int len){
poly B(1,1),C,D;
for(int limit=4;limit<(len<<2);limit<<=1){
C=A,C.resize(limit>>1);
init(limit),D=Inv(B,limit>>1);
C.resize(limit),NTT(C,limit,1);
D.resize(limit),NTT(D,limit,1);
for(int i=0;i<limit;++i)
C[i]=mul(C[i],D[i]);
NTT(C,limit,-1),B.resize(limit>>1);
for(int i=0;i<(limit>>1);++i)
B[i]=mul(inv2,add(B[i],C[i]));
}
B.resize(len);return B;
}

多项式除法 & & & 取模

problem:

给定一个 n n n 次多项式 A ( x ) A(x) A(x) 和一个 m m m 次多项式 B ( x ) B(x) B(x),求一个 ( n − m ) (n-m) (nm) 次多项式 C ( x ) C(x) C(x) ( m − 1 ) (m-1) (m1) 次多项式 D ( x ) D(x) D(x) 满足:

A ( x ) = B ( x ) C ( x ) + D ( x ) A(x)=B(x)C(x)+D(x) A(x)=B(x)C(x)+D(x)

solution:

对于一个最高项次数为 n n n 的多项式 A ( x ) A(x) A(x),我们考虑对它进行一个操作,即令 A R ( x ) = x n A ( 1 x ) A^R(x)=x^nA(frac 1 x) AR(x)=xnA(x1)

发现 A R ( x ) A^R(x) AR(x) 就是将 A ( x ) A(x) A(x) 的系数反转的式子(可以代几个 n n n 比较小的多项式来帮助自己理解)。

现在的式子是:

A ( x ) = B ( x ) C ( x ) + D ( x ) A(x)=B(x)C(x)+D(x) A(x)=B(x)C(x)+D(x)

其中, A ( x ) A(x) A(x) 最高项为 n n n B ( x ) B(x) B(x) 最高项为 m m m C ( x ) C(x) C(x) 最高项为 n − m n-m nm D ( x ) D(x) D(x) 的最高项为 m − 1 m-1 m1

那么两边同时乘 x n x^n xn,并将 1 x frac 1 x x1 代入,得到:

x n A ( 1 x ) = x m B ( 1 x ) x n − m C ( 1 x ) + x n − m + 1 × x m − 1 D ( 1 x ) x^nA(frac 1 x)=x^mB(frac 1 x)x^{n−m}C(frac 1 x)+x^{n−m+1}times x^{m−1}D(frac 1 x) xnA(x1)=xmB(x1)xnmC(x1)+xnm+1×xm1D(x1)

通过上面的公式,即 A R ( x ) = x n A ( 1 x ) A^R(x)=x^nA(frac 1 x) AR(x)=xnA(x1),得到:

A R ( x ) = B R ( x ) C R ( x ) + x n − m + 1 D R ( x ) A^R(x)=B^R(x)C^R(x)+x^{n-m+1}D^R(x) AR(x)=BR(x)CR(x)+xnm+1DR(x)

考虑在 m o d    x n − m + 1 mathrm{mod}; x^{n−m+1} modxnm+1 意义下, A , B A,B A,B 已知, C C C 最高为 n − m n-m nm 不影响,而 D D D 被消去。

所以我们有:

A R ( x ) ≡ B R ( x ) C R ( x ) (   m o d    x n − m + 1 ) A^R(x)equiv B^R(x)C^R(x)(,mathrm{mod};x^{n-m+1}) AR(x)BR(x)CR(x)(modxnm+1)

用多项式求逆即可解出 C R ( x ) C^R(x) CR(x),也就解出了 C ( x ) C(x) C(x)

最后把 C ( x ) C(x) C(x) 代回去就可以解出 D ( x ) D(x) D(x)

复杂度 O ( n log ⁡ n ) O(nlog n) O(nlogn)

code:

PS:代码中多项式的 ∗ * − - 是重载过的。

typedef vector<int> poly;
poly operator/(poly A,poly B){
int limit=1,len=A.size()-B.size()+1;
reverse(A.begin(),A.end());
reverse(B.begin(),B.end());
while(limit<=len)
limit<<=1;
B=Inv(B,limit),B.resize(len);
A=A*B,A.resize(len);
reverse(A.begin(),A.end());
return A;
}
poly operator%(poly A,poly B){
poly C=A-(A/B)*B;
C.resize(B.size()-1);
return C;
}

多项式微积分(求导 & & & 积分)

solution:

f ( x ) = x n f(x)=x^n f(x)=xn,那么它的导数 f ′ ( x ) = n x n − 1 f'(x)=nx^{n-1} f(x)=nxn1,根据这个可以快速求多项式的导数。

证明如下:
首先根据导数定义,有:

f ′ ( x ) = lim ⁡ Δ x → 0 ( x + Δ x ) n − x n Δ x f'(x)=lim_{Delta xrightarrow 0}frac{(x+Delta x)^n-x^n}{Delta x} f(x)=Δx0limΔx(x+Δx)nxn

二项式定理展开,得到:

f ′ ( x ) = lim ⁡ Δ x → 0 ∑ i = 0 n ( n i ) x i Δ x n − i − x n Δ x = lim ⁡ Δ x → 0 ∑ i = 0 n − 1 ( n i ) x i Δ x n − i + x n Δ x 0 − x n Δ x = lim ⁡ Δ x → 0 ∑ i = 0 n − 1 ( n i ) x i Δ x n − i − 1 = lim ⁡ Δ x → 0 ∑ i = 0 n − 2 ( n i ) x i Δ x n − i − 1 + ( n n − 1 ) x n − 1 begin{aligned} f'(x)&=lim_{Delta xrightarrow 0}frac{sum_{i=0}^nbinom n ix^iDelta x^{n-i}-x^n}{Delta x}\ &=lim_{Delta xrightarrow 0}frac{sum_{i=0}^{n-1}binom n ix^iDelta x^{n-i}+x^nDelta x^0-x^n}{Delta x}\ &=lim_{Delta xrightarrow 0}sum_{i=0}^{n-1}binom n ix^iDelta x^{n-i-1}\ &=lim_{Delta xrightarrow 0}sum_{i=0}^{n-2}binom n ix^iDelta x^{n-i-1}+binom n {n-1}x^{n-1} end{aligned} f(x)=Δx0limΔxi=0n(in)xiΔxnixn=Δx0limΔxi=0n1(in)xiΔxni+xnΔx0xn=Δx0limi=0n1(in)xiΔxni1=Δx0limi=0n2(in)xiΔxni1+(n1n)xn1

又由于有限个无穷小的乘积依旧是无穷小,所以得到:

f ′ ( x ) = ( n n − 1 ) x n − 1 = n x n − 1 f'(x)=binom n {n-1}x^{n-1}=nx^{n-1} f(x)=(n1n)xn1=nxn1

证毕。

积分好像可以简单理解为导数的逆运算,还不是很会。。。

code:

PSinv 数组是预处理的逆元。

typedef vector<int> poly;
poly Deriv(poly A){
for(int i=0;i<A.size()-1;++i)
A[i]=mul(A[i+1],i+1);
A.pop_back();return A;
}
poly Integ(poly A){
A.push_back(0);
for(int i=A.size()-1;i;--i)
A[i]=mul(A[i-1],inv[i]);
A[0]=0;return A;
}

多项式 Ln

problem:

给出 ( n − 1 ) (n-1) (n1) 次多项式 A ( x ) A(x) A(x),求一个   m o d     x n bmod{:x^n} modxn 下的多项式 B ( x ) B(x) B(x),满足:

B ( x ) ≡ ln ⁡ A ( x ) (   m o d     x n ) B(x) equiv ln A(x)(bmod ;x^n) B(x)lnA(x)(modxn)

保证 A ( 0 ) = 1 A(0)=1 A(0)=1

solution:

有一个结论,即若 B ( x ) = ln ⁡ A ( x ) B(x) = ln A(x) B(x)=lnA(x),则 A ′ ( x ) = A ( x ) B ′ ( x ) A'(x)=A(x)B'(x) A(x)=A(x)B(x)

证明如下:
首先有链式法则,即对于复合函数 ( g ∘ f ) ( x ) = g ( f ( x ) ) (gcirc f)(x)=g(f(x)) (gf)(x)=g(f(x)),有 ( g ∘ f ) ′ ( x ) = g ′ ( f ( x ) ) f ′ ( x ) (gcirc f)'(x)=g'(f(x))f'(x) (gf)(x)=g(f(x))f(x)
那么对于 B ( x ) = ln ⁡ A ( x ) B(x)=ln A(x) B(x)=lnA(x),两边同时求导,得到 B ′ ( x ) = ( ln ⁡ A ( x ) ) ′ B'(x)=(ln A(x))' B(x)=(lnA(x))
等式右边可以看做一个复合函数,用链式法则化简后就是 ( ln ⁡ ′ A ( x ) ) A ′ ( x ) = A ′ ( x ) A ( x ) (ln' A(x))A'(x)=frac{A'(x)}{A(x)} (lnA(x))A(x)=A(x)A(x)
带回去,即可得到 A ′ ( x ) = A ( x ) B ′ ( x ) A'(x)=A(x)B'(x) A(x)=A(x)B(x)
证毕。

那么 B ′ ( x ) = A ′ ( x ) A ( x ) B'(x)=frac{A'(x)}{A(x)} B(x)=A(x)A(x),给 A ( x ) A(x) A(x) 求逆求导就可以算出 B ′ ( x ) B'(x) B(x)

然后再对 B ′ ( x ) B'(x) B(x) 积分就能解出 B ( x ) B(x) B(x)

code:

typedef vector<int> poly;
poly Ln(poly A,int len){
A=Integ(Deriv(A)*Inv(A,len)),A.resize(len);
return A;
}

多项式 Exp

前置知识 1 1 1:泰勒展开

泰勒展开是将一个在 x = x 0 x=x_0 x=x0 处具有 n n n 阶导数的函数 f ( x ) f(x) f(x) 利用关于 ( x − x 0 ) (x-x_0) (xx0) n n n 次多项式来逼近函数的方法。

若函数 f ( x ) f(x) f(x) 在包含 x 0 x_0 x0 的某个闭区间 [ a , b ] [a,b] [a,b] 上具有 n n n 阶导数,且在开区间 ( a , b ) (a,b) (a,b) 上具有 ( n + 1 ) (n+1) (n+1) 阶导数,则对闭区间 [ a , b ] [a,b] [a,b] 上任意一点 x x x,下式成立:

f ( x ) = f ( x 0 ) 0 ! + f ′ ( x 0 ) 1 ! ( x − x 0 ) + f ′ ′ ( x 0 ) 2 ! ( x − x 0 ) 2 + ⋯ + f ( n ) ( x 0 ) n ! ( x − x 0 ) n + ⋯ f(x)=frac{f(x_0)}{0!}+frac{f'(x_0)}{1!}(x-x_0)+frac{f''(x_0)}{2!}(x-x_0)^2+cdots+frac{f^{(n)}(x_0)}{n!}(x-x_0)^n+cdots f(x)=0!f(x0)+1!f(x0)(xx0)+2!f(x0)(xx0)2++n!f(n)(x0)(xx0)n+

我简单说一下我理解的泰勒展开。

等式右边的式子是一个多项式,所以我们的目标就是找一个多项式 g ( x ) g(x) g(x),使它与 f ( x ) f(x) f(x) 完全相同。

首先我们保证初始点 g ( x 0 ) = f ( x 0 ) g(x_0)=f(x_0) g(x0)=f(x0)

在此基础上,考虑曲线的变化趋势,即导数,保证在此处的导数相同,即 g ′ ( x 0 ) = f ′ ( x 0 ) g'(x_0)=f'(x_0) g(x0)=f(x0)

要进一步精确化,再考虑凹凸性,由于表征图像的凹凸性的参数为导数的导数,我们让二者导数的导数相等,即 g ′ ′ ( x 0 ) = f ′ ′ ( x 0 ) g''(x_0)=f''(x_0) g(x0)=f(x0)

那么以此类推,继续精确,让 x 0 x_0 x0 更高阶的导数相等,即 g ( n ) ( x 0 ) = f ( n ) ( x 0 ) g^{(n)}(x_0)=f^{(n)}(x_0) g(n)(x0)=f(n)(x0)

上面的式子归纳一下,就是 ∀ i , g ( i ) ( x 0 ) = f ( i ) ( x 0 ) forall i,g^{(i)}(x_0)=f^{(i)}(x_0) i,g(i)(x0)=f(i)(x0)

因为 x n x^n xn 的导数为 n x n − 1 nx^{n-1} nxn1,所以求 g ( x ) g(x) g(x) i i i 阶导数时为 a i × i ! = f ( i ) ( x 0 ) a_itimes i!=f^{(i)}(x_0) ai×i!=f(i)(x0),所以 a i = f ( i ) ( x 0 ) i ! a_i=frac{f^{(i)}(x_0)}{i!} ai=i!f(i)(x0)

可以看看这个,对理解泰勒展开会有帮助。

前置知识 2 2 2:牛顿迭代

摘自 Miskcoo’s Space 的牛顿迭代法在多项式运算的应用。

有这样一个问题,即已知一个函数 g ( x ) g(x) g(x),求一个多项式 f ( x )   m o d     x n f(x)bmod ;x^n f(x)modxn,满足方程

g ( f ( x ) ) ≡ 0 (   m o d     x n ) g(f(x))≡0(bmod;x^n) g(f(x))0(modxn)

首先在 n = 1 n=1 n=1 的时候, g ( f ( x ) ) ≡ 0 (   m o d     x ) g(f(x))≡0(bmod;x) g(f(x))0(modx),这是可以直接求出来的。

现在假设已经求出了 g ( f 0 ( x ) ) ≡ 0 (   m o d     x ⌈ n 2 ⌉ ) g(f_0(x))≡0(bmod;x^{⌈frac n 2⌉}) g(f0(x))0(modx2n),考虑如何扩展到   m o d     x n bmod;x^n modxn 下。

我们把 g ( f ( x ) ) g(f(x)) g(f(x)) f 0 ( x ) f_0(x) f0(x) 进行泰勒展开,得到:

g ( f ( x ) ) = g ( f 0 ( x ) ) + g ′ ( f 0 ( x ) ) 1 ! ( f ( x ) − f 0 ( x ) ) + g ′ ′ ( f 0 ( x ) ) 2 ! ( f ( x ) − f 0 ( x ) ) 2 + ⋯ g(f(x))=g(f_0(x))+frac{g'(f_0(x))}{1!}(f(x)−f_0(x))+frac{g''(f_0(x))}{2!}(f(x)−f_0(x))^2+⋯ g(f(x))=g(f0(x))+1!g(f0(x))(f(x)f0(x))+2!g(f0(x))(f(x)f0(x))2+

首先显然有 f ( x ) ≡ f 0 ( x ) (   m o d     x ⌈ n 2 ⌉ ) f(x)≡f_0(x)(bmod ;x^{⌈frac n 2⌉}) f(x)f0(x)(modx2n),所以 f ( x ) − f 0 ( x ) f(x)-f_0(x) f(x)f0(x) 最低项次数一定大于 ⌈ n 2 ⌉ lceil frac n 2rceil 2n。因此在   m o d     x n bmod ;x^n modxn 意义下,整个式子从 ( f ( x ) − f 0 ( x ) ) 2 (f(x)−f_0(x))^2 (f(x)f0(x))2 这一项开始都为 0 0 0 了,也即:

g ( f ( x ) ) ≡ g ( f 0 ( x ) ) + g ′ ( f 0 ( x ) ) ( f ( x ) − f 0 ( x ) ) (   m o d     x n ) g(f(x))≡g(f_0(x))+g'(f_0(x))(f(x)−f_0(x))(bmod ;x^n) g(f(x))g(f0(x))+g(f0(x))(f(x)f0(x))(modxn)

又因为 g ( f ( x ) ) ≡ 0 (   m o d     x n ) g(f(x))≡0(bmod;x^n) g(f(x))0(modxn),得到

f ( x ) ≡ f 0 ( x ) − g ( f 0 ( x ) ) g ′ ( f 0 ( x ) ) (   m o d     x n ) f(x)≡f_0(x)−frac{g(f_0(x))}{g'(f_0(x))}(bmod;x^n) f(x)f0(x)g(f0(x))g(f0(x))(modxn)

这样就大功告成了。


OK,下面就步入正题,多项式 Exp

problem:

给出 ( n − 1 ) (n-1) (n1) 次多项式 A ( x ) A(x) A(x),求一个   m o d      x n bmod{;x^n} modxn 下的多项式 B ( x ) B(x) B(x),满足:

B ( x ) ≡ e A ( x ) (   m o d     x n ) B(x) equiv e^{A(x)}(bmod;x^n) B(x)eA(x)(modxn)

保证 A ( 0 ) = 0 A(0)=0 A(0)=0

solution:

发现原问题可以转化一下,即求解:

ln ⁡ B ( x ) ≡ A ( x ) (   m o d     x n ) ln B(x)equiv A(x)(bmod ; x^n) lnB(x)A(x)(modxn)

还是令 n = 2 k n=2^k n=2k,假设我们已经知道了 ln ⁡ B 0 ( x ) − A ( x ) ≡ 0 (   m o d     x n 2 ) ln B_0(x)-A(x)equiv 0(bmod;x^{frac n 2}) lnB0(x)A(x)0(modx2n)

构造 P ( B ( x ) ) = ln ⁡ B ( x ) − A ( x ) P(B(x))=ln B(x)-A(x) P(B(x))=lnB(x)A(x),那么 P ′ ( B ( x ) ) = 1 B ( x ) P'(B(x))=frac{1}{B(x)} P(B(x))=B(x)1

那么按照牛顿迭代的推理可得:

B ( x ) ≡ B 0 ( x ) − P ( B 0 ( x ) ) P ′ ( B 0 ( x ) ) ≡ B 0 ( x ) ( 1 − ln ⁡ B 0 ( x ) + A ( x ) ) (   m o d     x n ) B(x)≡B_0(x)− frac{P(B_0(x))}{P'(B_0(x))}≡B_0(x)(1−ln B_0(x)+A(x))(bmod ;x^n) B(x)B0(x)P(B0(x))P(B0(x))B0(x)(1lnB0(x)+A(x))(modxn)

时间复杂度 O ( n log ⁡ n ) O(nlog n) O(nlogn)

code:

typedef vector<int> poly;
poly Exp(poly A,int len){
poly B(1,1),C;A.resize(len<<1);
for(int limit=2;limit<(len<<1);limit<<=1){
C=Ln(B,limit);
for(int i=0;i<limit;++i)
C[i]=dec(A[i],C[i]);
C[0]=add(C[0],1),B=B*C;
B.resize(limit);
}
B.resize(len);return B;
}

多项式快速幂

problem:

给定一个 ( n − 1 ) (n-1) (n1) 次多项式 A ( x ) A(x) A(x),求一个在   m o d     x n bmod; x^n modxn 意义下的多项式 B ( x ) B(x) B(x),满足:

B ( x ) ≡ A k ( x )   (   m o d     x n ) B(x) equiv A^k(x) (bmod; x^n) B(x)Ak(x) (modxn)

保证 A ( 0 ) = 1 A(0)=1 A(0)=1

solution:

如果直接用快速幂+ n t t ntt ntt的话,时间复杂度是 O ( n log ⁡ 2 n ) O(nlog^2n) O(nlog2n) 的,而且常数还比较大。

所以我们换一种方法,即发现:

A k ( x ) = ( e ln ⁡ A ( x ) ) k = e k × ln ⁡ A ( x ) A^k(x)=(e^{ln A(x)})^k=e^{ktimes ln A(x)} Ak(x)=(elnA(x))k=ek×lnA(x)

所以我们用多项式 Ln,Exp 就可以了。

时间复杂度 O ( n log ⁡ n ) O(nlog n) O(nlogn)

code:

typedef vector<int> poly;
poly Power(poly A,int len,int k){
A=Exp(Ln(A,len)*k,len),A.resize(len);
return A;
}

多项式多点求值

problem:

给定一个 n n n 次多项式 f ( x ) f(x) f(x),并给出 m m m 个整数 a i a_i ai,请你对于 ∀ i ∈ [ 1 , m ] forall i in [1,m] i[1,m],求出 f ( a i ) f(a_i) f(ai)

solution:

首先,根据余式定理,我们有:

f ( x ) ≡ f ( x 0 ) (   m o d     ( x − x 0 ) ) f(x)equiv f(x_0)(bmod;(x-x_0)) f(x)f(x0)(mod(xx0))

简单写一下它的证明:
f ( x ) = g ( x ) ( x − x 0 ) + A ( x ) f(x)=g(x)(x-x_0)+A(x) f(x)=g(x)(xx0)+A(x),那么 A ( x ) A(x) A(x) 为就是 f ( x ) f(x) f(x) ( x − x 0 ) (x-x_0) (xx0) 取模的结果,即 f ( x ) ≡ A ( x ) (   m o d   ( x − x 0 ) ) f(x)equiv A(x)(bmod (x-x_0)) f(x)A(x)(mod(xx0))
由于 ( x − x 0 ) (x-x_0) (xx0) 的最高项次数为 1 1 1,所以 A ( x ) A(x) A(x) 最高项次数为 0 0 0,即为常数项。
x 0 x_0 x0 代入 f ( x ) f(x) f(x),得到 f ( x 0 ) = g ( x 0 ) ( x 0 − x 0 ) + A ( x 0 ) = A ( x ) f(x_0)=g(x_0)(x_0-x_0)+A(x_0)=A(x) f(x0)=g(x0)(x0x0)+A(x0)=A(x),所以就有 f ( x ) ≡ f ( x 0 ) (   m o d     ( x − x 0 ) ) f(x)equiv f(x_0)(bmod;(x-x_0)) f(x)f(x0)(mod(xx0))

但是直接这样做是不行的,对于一个点,多项式取模是 O ( n log ⁡ n ) O(nlog n) O(nlogn) 的,甚至比 O ( n ) O(n) O(n) 的暴力还慢。

然后考虑用分治优化。

G l , r ( x ) = ∏ i = l r ( x − x i ) G_{l,r}(x)=prodlimits_{i=l}^{r}(x-x_i) Gl,r(x)=i=lr(xxi),我们现在已经算出了 f ( x ) ≡ A ( x ) ( m o d ∏ i = l r ( x − x i ) ) f(x)equiv A(x)pmod{prodlimits_{i=l}^r(x-x_i)} f(x)A(x)(modi=lr(xxi)),因此只需要用 A ( x ) A(x) A(x) 分别对 G l , ⌊ l + r 2 ⌋ ( x ) G_{l,lfloorfrac{l+r}{2}rfloor}(x) Gl,2l+r(x) G ⌊ l + r 2 ⌋ + 1 , r ( x ) G_{lfloorfrac{l+r}{2}rfloor+1,r}(x) G2l+r+1,r(x) 取模,再递归下去。分治到叶子的时候,就只剩下 f ( x ) ≡ A ( x ) ( m o d ( x − x i ) ) f(x)equiv A(x)pmod{(x-x_i)} f(x)A(x)(mod(xxi)),此时的 A ( x ) A(x) A(x) 就是 i i i 的答案。

每个 G ( x ) G(x) G(x) 可以提前分治 n t t ntt ntt预处理出来。

时间复杂度 T ( n ) = T ( n 2 ) + O ( n log ⁡ n ) = O ( n log ⁡ 2 n ) T(n)=T(frac{n}{2})+O(nlog n)=O(nlog^2 n) T(n)=T(2n)+O(nlogn)=O(nlog2n)

不知道为什么跑得贼慢,也许是我自带大常数吧。

code:

typedef vector<int> poly;
poly F,G[N];
void build(int root,int l,int r){
if(l==r){
G[root].clear();
G[root].push_back(P-a[l]),G[root].push_back(1);
return;
}
int mid=(l+r)>>1;
build(root<<1,l,mid),build(root<<1|1,mid+1,r);
G[root]=G[root<<1]*G[root<<1|1];
}
void calc(int root,int l,int r,poly now){
if(l==r)
{ans[l]=now[0];return;}
int mid=(l+r)>>1;
calc(root<<1,l,mid,now%G[root<<1]);
calc(root<<1|1,mid+1,r,now%G[root<<1|1]);
}

多项式快速插值

前置知识 1 1 1:拉格朗日插值

给定 ( n + 1 ) (n+1) (n+1) 个点 ( x i , y i ) (x_i,y_i) (xi,yi),其中每个 x x x 互不相同,那么存在唯一的次数不超过 n n n 的多项式 f ( x ) f(x) f(x),满足 ∀ i ∈ [ 1 , n ] forall iin [1,n] i[1,n] f ( x i ) = y i f(x_i)=y_i f(xi)=yi,即:

f ( x ) = ∑ i = 1 n y i ∏ j ≠ i x − x j x i − x j f(x)=sum_{i=1}^ny_iprod_{jne i}frac{x-x_j}{x_i-x_j} f(x)=i=1nyij=ixixjxxj

这个应该比较显然,因为对于 x k x_k xk,只有当上式中 i = k i=k i=k ∏ j ≠ i x − x j x i − x j = 1 prod_{jne i}frac{x-x_j}{x_i-x_j}=1 j=ixixjxxj=1,否则是 0 0 0

前置知识 2 2 2:洛必达法则

lim ⁡ x → a f ( x ) = 0 limlimits_{xrightarrow a}f(x)=0 xalimf(x)=0 lim ⁡ x → a g ( x ) = 0 limlimits_{xrightarrow a}g(x)=0 xalimg(x)=0,在点 a a a 的某去心邻域内两者都可导,且 g ′ ( x ) ≠ 0 g'(x)ne 0 g(x)=0,那么有:

lim ⁡ x → a f ( x ) g ( x ) = lim ⁡ x → a f ′ ( x ) g ′ ( x ) lim_{xrightarrow a}frac{f(x)}{g(x)}=lim_{xrightarrow a}frac{f'(x)}{g'(x)} xalimg(x)f(x)=xalimg(x)f(x)


接下来回到我们的问题。

problem:

给出 n n n 个点 ( x i , y i ) (x_i, y_i) (xi,yi),求一个 n − 1 n-1 n1 次的多项式 f ( x ) f(x) f(x),使得 f ( x i ) ≡ y i ( m o d 998244353 ) f(x_i)equiv y_ipmod{998244353} f(xi)yi(mod998244353)

solution:

考虑到直接用求拉格朗日插值是 O ( n 2 ) O(n^2) O(n2) 的,我们把它变一下形:

f ( x ) = ∑ i = 1 n y i ∏ j ≠ i ( x i − x j ) ∏ j ≠ i ( x − x j ) f(x)=sum_{i=1}^nfrac{y_i}{prod_{jne i}(x_i-x_j)}prod_{jne i} (x-x_j) f(x)=i=1nj=i(xixj)yij=i(xxj)

考虑先求出 G = ∏ j ≠ i ( x i − x j ) G=prod_{jne i}(x_i-x_j) G=j=i(xixj)

g ( x ) = ∏ j = 1 n ( x − x j ) g(x)=prod_{j=1}^n(x-x_j) g(x)=j=1n(xxj),那么 j ≠ i jne i j=i 的限制就相当于除了个 x − x i x-x_i xxi,所以 G = g ( x i ) x − x i G=frac{g(x_i)}{x-x_i} G=xxig(xi)

发现这时分子分母都是 0 0 0,根据洛必达法则,有 G = g ′ ( x i ) 1 = g ′ ( x i ) G=frac{g'(x_i)}{1}=g'(x_i) G=1g(xi)=g(xi)

接下来考虑分治优化,设 f l , r f_{l,r} fl,r 表示 [ l , r ] [l,r] [l,r] 的答案, m i d = ⌊ l + r 2 ⌋ mid=lfloorfrac{l+r}{2}rfloor mid=2l+r,那么有:

f l , r = ∑ i = l r y i g ′ ( x i ) ∏ j = l r [ j ≠ i ] ( x − x j ) = ∏ k = m i d + 1 r ( x − x k ) ∑ i = l m i d y i g ′ ( x i ) ∏ j = l m i d [ j ≠ i ] ( x − x j ) + ∏ k = l m i d ( x − x k ) ∑ i = m i d + 1 r y i g ′ ( x i ) ∏ j = m i d + 1 r [ j ≠ i ] ( x − x j ) = f l , m i d ∏ k = m i d + 1 r ( x − x k ) + f m i d + 1 , r ∏ k = l m i d ( x − x k ) begin{aligned} f_{l,r}&=sum_{i=l}^rfrac{y_i}{g'(x_i)}prod_{j=l}^{r}[jne i](x-x_j)\ &=prod _{k=mid+1}^r(x-x_k)sum_{i=l}^{mid}frac{y_i}{g'(x_i)}prod_{j=l}^{mid}[jne i](x-x_j)+prod _{k=l}^{mid}(x-x_k)sum_{i=mid+1}^rfrac{y_i}{g'(x_i)}prod_{j=mid+1}^{r}[jne i](x-x_j)\ &=f_{l,mid}prod _{k=mid+1}^r(x-x_k)+f_{mid+1,r}prod _{k=l}^{mid}(x-x_k) end{aligned} fl,r=i=lrg(xi)yij=lr[j=i](xxj)=k=mid+1r(xxk)i=lmidg(xi)yij=lmid[j=i](xxj)+k=lmid(xxk)i=mid+1rg(xi)yij=mid+1r[j=i](xxj)=fl,midk=mid+1r(xxk)+fmid+1,rk=lmid(xxk)

其实就是把 f l , r f_{l,r} fl,r 拆成 f l , m i d f_{l,mid} fl,mid f m i d + 1 , r f_{mid+1,r} fmid+1,r 的式子。

然后用分治 n t t ntt ntt预处理出 g g g,然后多点求值求 g ′ ( x i ) g'(x_i) g(xi),最后分治求 f 1 , n f_{1,n} f1,n 即可。

时间复杂度 O ( n log ⁡ n ) O(nlog n) O(nlogn)

code:

PS: h h h 数组就是事先处理出的 y i g ′ ( x i ) frac{y_i}{g'(x_i)} g(xi)yi

poly Get(int root,int l,int r){
if(l==r)
return poly(1,h[l]);
int mid=(l+r)>>1;
poly L=Get(root<<1,l,mid),R=Get(root<<1|1,mid+1,r);
return L*G[root<<1|1]+R*G[root<<1];
}

下降幂多项式乘法

problem:

给定一个 n n n 次下降幂多项式 A ( x ) A(x) A(x) m m m 次下降幂多项式 B ( x ) B(x) B(x),你要求出一个 ( n + m ) (n+m) (n+m) 次下降幂多项式 F ( x ) F(x) F(x) 满足 F ( x ) = A ( x ) B ( x ) F(x)=A(x)B(x) F(x)=A(x)B(x)

solution:

首先有一个重要的结论,即:

e x = ∑ i = 0 ∞ x i i ! e^x=sum_{i=0}^inftyfrac{x^i}{i!} ex=i=0i!xi

这个根据 e x e^x ex 0 0 0 处泰勒展开即可得到。

考虑对于 x n ‾ x^{underline n} xn 构建 EGF

x n ‾ = ∑ i = n ∞ i n ‾ i ! x i = ∑ i = n ∞ 1 ( i − n ) ! x i = ∑ i = n ∞ x i − n ( i − n ) ! x n = x n e x begin{aligned} x^{underline n}&=sum_{i=n}^{infty}frac{i^{underline n}}{i!}x^i\ &=sum_{i=n}^{infty}frac{1}{(i-n)!}x^i\ &=sum_{i=n}^{infty}frac{x^{i-n}}{(i-n)!}x^n\ &=x^ne^x end{aligned} xn=i=ni!inxi=i=n(in)!1xi=i=n(in)!xinxn=xnex

考虑对于 f ( x ) = ∑ i = 0 ∞ a i x i ‾ f(x)=sum_{i=0}^{infty}a_ix^{underline i} f(x)=i=0aixi 的点值构造 EGF:

g ( x ) = ∑ i = 0 ∞ f ( i ) i ! x i = ∑ i = 0 ∞ x i i ! ∑ j = 0 ∞ a j i j ‾ = ∑ j = 0 ∞ a j ∑ i = 0 ∞ i j ‾ i ! x i = ∑ j = 0 ∞ a j x j e x = e x ∑ i = 0 ∞ a i x i begin{aligned} g(x)&=sum_{i=0}^{infty}frac{f(i)}{i!}x^i\ &=sum_{i=0}^{infty}frac{x^i}{i!}sum_{j=0}^{infty}a_ji^{underline j}\ &=sum_{j=0}^{infty}a_jsum_{i=0}^{infty}frac{i^{underline j}}{i!}x^i\ &=sum_{j=0}^{infty}a_jx^je^x\ &=e^xsum_{i=0}^{infty}a_ix^i end{aligned} g(x)=i=0i!f(i)xi=i=0i!xij=0ajij=j=0aji=0i!ijxi=j=0ajxjex=exi=0aixi

因此只需要用普通多项式的系数卷个 e x e^x ex 就得到了点值的 EGF。

点值还原原多项式也只需要卷个 e − x e^{-x} ex

代码很简单,就是之前的操作,就不放了。

最后

以上就是小巧帆布鞋为你收集整理的多项式各种操作的全部内容,希望文章能够帮你解决多项式各种操作所遇到的程序开发问题。

如果觉得靠谱客网站的内容还不错,欢迎将靠谱客网站推荐给程序员好友。

本图文内容来源于网友提供,作为学习参考使用,或来自网络收集整理,版权属于原作者所有。
点赞(49)

评论列表共有 0 条评论

立即
投稿
返回
顶部