一.生成函数与多项式乘法.
众所周知,生成函数是OI中的计数利器,而很多生成函数的题目需要快速计算式子中的多项式运算,所以本篇主要是介绍多项式运算的.
值得注意的是,这里很多奇怪的多项式运算本来是不存在模 x n x^n xn意义下的定义的,但是由于生成函数本质是形式幂级数,所以直接取泰勒展开后的结果就是了.
多项式乘法作为最基础的多项式运算,这里不再介绍,具体参见快速傅里叶变换FFT与快速数论变换NTT入门.
二.多项式求逆.
多项式的逆:对于一个多项式
F
(
x
)
F(x)
F(x),其逆
G
(
x
)
G(x)
G(x)定义为:
G
(
x
)
F
(
x
)
=
1
G(x)F(x)=1
G(x)F(x)=1
但很多情况下
G
(
x
)
G(x)
G(x)是一个无穷级数的形式,所以一般我们只需要取其最后
n
n
n项(次数低于
n
n
n次),记为:
G
(
x
)
F
(
x
)
≡
1
(
m
o
d
x
n
)
G(x)F(x)\equiv 1\,\,(mod\,\,x^{n})
G(x)F(x)≡1(modxn)
当然可以暴力求解,但是有一种倍增的思路:我们假设已经求出了模 x ⌈ n 2 ⌉ x^{\left\lceil\frac{n}{2}\right\rceil} x⌈2n⌉意义下 F ( x ) F(x) F(x)的逆 H ( x ) H(x) H(x),推出 G ( x ) G(x) G(x)与 H ( x ) H(x) H(x)之间的关系.
现在有两个条件:
G
(
x
)
F
(
x
)
≡
1
(
m
o
d
x
⌈
n
2
⌉
)
H
(
x
)
F
(
x
)
≡
1
(
m
o
d
x
⌈
n
2
⌉
)
G(x)F(x)\equiv 1 \,\,(mod\,\,x^{\left\lceil\frac{n}{2}\right\rceil})\\ H(x)F(x)\equiv 1\,\,(mod\,\,x^{\left\lceil\frac{n}{2}\right\rceil})
G(x)F(x)≡1(modx⌈2n⌉)H(x)F(x)≡1(modx⌈2n⌉)
将两式相减可得:
(
G
(
x
)
−
H
(
x
)
)
F
(
x
)
≡
0
(
m
o
d
x
⌈
n
2
⌉
)
G
(
x
)
−
H
(
x
)
≡
0
(
m
o
d
x
⌈
n
2
⌉
)
(G(x)-H(x))F(x)\equiv 0\,\,(mod\,\,x^{\left\lceil\frac{n}{2}\right\rceil})\\ G(x)-H(x)\equiv 0\,\,(mod\,\,x^{\left\lceil\frac{n}{2}\right\rceil})
(G(x)−H(x))F(x)≡0(modx⌈2n⌉)G(x)−H(x)≡0(modx⌈2n⌉)
再平方一下可以得到:
(
G
(
x
)
−
H
(
x
)
)
2
≡
0
(
m
o
d
x
n
)
G
2
(
x
)
+
H
2
(
x
)
−
2
G
(
x
)
H
(
x
)
≡
0
(
m
o
d
x
n
)
(G(x)-H(x))^{2}\equiv 0\,\,(mod\,\,x^{n})\\ G^{2}(x)+H^{2}(x)-2G(x)H(x)\equiv 0\,\,(mod\,\,x^{n})
(G(x)−H(x))2≡0(modxn)G2(x)+H2(x)−2G(x)H(x)≡0(modxn)
两边同时乘上
F
(
x
)
F(x)
F(x)得到:
F
(
x
)
G
2
(
x
)
+
F
(
x
)
H
2
(
x
)
−
2
F
(
x
)
G
(
x
)
H
(
x
)
≡
0
(
m
o
d
x
n
)
G
(
x
)
+
F
(
x
)
H
2
(
x
)
−
2
H
(
x
)
≡
0
(
m
o
d
x
n
)
G
(
x
)
≡
2
H
(
x
)
−
F
(
x
)
H
2
(
x
)
(
m
o
d
x
n
)
G
(
x
)
≡
H
(
x
)
(
2
−
F
(
x
)
H
(
x
)
)
(
m
o
d
x
n
)
F(x)G^{2}(x)+F(x)H^{2}(x)-2F(x)G(x)H(x)\equiv 0\,\,(mod\,\,x^{n})\\ G(x)+F(x)H^{2}(x)-2H(x)\equiv 0\,\,(mod\,\,x^{n})\\ G(x)\equiv 2H(x)-F(x)H^{2}(x)\,\,(mod\,\,x^{n})\\ G(x)\equiv H(x)(2-F(x)H(x))\,\,(mod\,\,x^{n})
F(x)G2(x)+F(x)H2(x)−2F(x)G(x)H(x)≡0(modxn)G(x)+F(x)H2(x)−2H(x)≡0(modxn)G(x)≡2H(x)−F(x)H2(x)(modxn)G(x)≡H(x)(2−F(x)H(x))(modxn)
然后就可以快乐倍增了.
时间复杂度 O ( n log n ) O(n\log n) O(nlogn).
代码如下:
int ta[N+9];
void Poly_inv(int *a,int n,int *res){
if (!n) {res[0]=Get_inv(a[0]);return;}
Poly_inv(a,n>>1,res);
for (int i=0;i<=n;++i) ta[i]=a[i];
Get_len(n<<1);
NTT(res,len,0);
NTT(ta,len,0);
for (int i=0;i<len;++i) smul(res[i],sub(2,mul(res[i],ta[i]))),ta[i]=0;
NTT(res,len,1);
for (int i=n+1;i<=n<<1;++i) res[i]=0;
}
三.多项式对数.
多项式对数主要是指多项式
ln
\ln
ln,即给定多项式
F
(
x
)
F(x)
F(x),求
G
(
x
)
G(x)
G(x)满足:
e
G
(
x
)
≡
F
(
x
)
(
m
o
d
x
n
)
e^{G(x)}\equiv F(x)\,\,(mod\,\,x^{n})
eG(x)≡F(x)(modxn)
对数函数在求导后就变成除法了,所以尝试先求导再积分:
ln
F
(
x
)
=
∫
(
ln
F
(
x
)
)
′
d
x
=
∫
ln
′
F
(
x
)
F
′
(
x
)
d
x
=
∫
F
′
(
x
)
F
(
x
)
d
x
\ln F(x)=\int(\ln F(x))'\mathrm{d}x=\int\ln'F(x)F'(x)\mathrm{d}x=\int\frac{F'(x)}{F(x)}\mathrm{d}x
lnF(x)=∫(lnF(x))′dx=∫ln′F(x)F′(x)dx=∫F(x)F′(x)dx
现在的式子就可以直接做了.
不过有一个问题是,这样子多项式 ln \ln ln必须得满足 F ( x ) F(x) F(x)的常数项为 1 1 1,最终求出的 G ( x ) G(x) G(x)的常数项为 0 0 0.
至于原因,主要是这样按照泰勒展开定义的多项式 ln \ln ln必须满足这样的条件.
时间复杂度 O ( n log n ) O(n\log n) O(nlogn).
代码如下:
void Poly_ln(int *a,int n,int *res){
if (!n) return;
Poly_inv(a,n-1,res);
for (int i=1;i<=n;++i) ta[i-1]=mul(i,a[i]);
Get_len(n-1<<1);
NTT(res,len,0);
NTT(ta,len,0);
for (int i=0;i<len;++i) smul(res[i],ta[i]),ta[i]=0;
NTT(res,len,1);
for (int i=n;i>=1;--i) res[i]=mul(inv[i],res[i-1]);
res[0]=0;
for (int i=n+1;i<n<<1;++i) res[i]=0;
}
四.牛顿迭代.
牛顿迭代是一种估算方程解的方法,不过这里我们直接应用到多项式上了.
对于给定的规则函数
T
(
x
)
T(x)
T(x),我们有方程:
T
(
F
(
x
)
)
≡
0
(
m
o
d
x
n
)
T(F(x))\equiv 0\,\,(mod\,\,x^{n})
T(F(x))≡0(modxn)
现在我们要求 F ( x ) F(x) F(x).
考虑倍增求解,如果我们已经求出了在模
x
⌈
n
2
⌉
x^{\left\lceil\frac{n}{2}\right\rceil}
x⌈2n⌉意义下的解
G
(
x
)
G(x)
G(x),即:
T
(
G
(
x
)
)
≡
0
(
m
o
d
x
⌈
n
2
⌉
)
T(G(x))\equiv 0\,\,(mod\,\,x^{\left\lceil\frac{n}{2}\right\rceil})
T(G(x))≡0(modx⌈2n⌉)
那么我们将
T
(
F
(
x
)
)
T(F(x))
T(F(x))在
G
(
x
)
G(x)
G(x)处泰勒展开:
T
(
F
(
x
)
)
=
∑
i
=
0
+
∞
T
(
i
)
(
G
(
x
)
)
i
!
(
F
(
x
)
−
G
(
x
)
)
i
=
T
(
G
(
x
)
)
+
T
′
(
G
(
x
)
)
(
F
(
x
)
−
G
(
x
)
)
+
1
2
T
′
′
(
G
(
x
)
)
(
F
(
x
)
−
G
(
x
)
)
2
+
⋯
T(F(x))=\sum_{i=0}^{+\infty}\frac{T^{(i)}(G(x))}{i!}(F(x)-G(x))^{i}\\ =T(G(x))+T'(G(x))(F(x)-G(x))+\frac{1}{2}T''(G(x))(F(x)-G(x))^{2}+\cdots
T(F(x))=i=0∑+∞i!T(i)(G(x))(F(x)−G(x))i=T(G(x))+T′(G(x))(F(x)−G(x))+21T′′(G(x))(F(x)−G(x))2+⋯
由于我们求的是模
x
n
x^{n}
xn意义下的结果,所以对于
(
F
(
x
)
−
G
(
x
)
)
2
(F(x)-G(x))^{2}
(F(x)−G(x))2有:
(
F
(
x
)
−
G
(
x
)
)
2
≡
0
(
m
o
d
n
)
(F(x)-G(x))^{2}\equiv 0\,\,(mod\,\,n)
(F(x)−G(x))2≡0(modn)
也就是说:
T
(
F
(
x
)
)
≡
0
(
m
o
d
x
n
)
T
(
G
(
x
)
)
+
T
′
(
G
(
x
)
)
(
F
(
x
)
−
G
(
x
)
)
≡
0
(
m
o
d
x
n
)
F
(
x
)
≡
G
(
x
)
−
T
(
G
(
x
)
)
T
′
(
G
(
x
)
)
(
m
o
d
x
n
)
T(F(x))\equiv 0\,\,(mod\,\,x^{n})\\ T(G(x))+T'(G(x))(F(x)-G(x))\equiv 0\,\,(mod\,\,x^{n})\\ F(x)\equiv G(x)-\frac{T(G(x))}{T'(G(x))}\,\,(mod\,\,x^{n})
T(F(x))≡0(modxn)T(G(x))+T′(G(x))(F(x)−G(x))≡0(modxn)F(x)≡G(x)−T′(G(x))T(G(x))(modxn)
推导出来的结果在做多项式操作的时候非常有效.
五.多项式指数.
多项式指数主要是指多项式
exp
\exp
exp,即给定多项式
F
(
x
)
F(x)
F(x),求
G
(
x
)
G(x)
G(x)满足:
G
(
x
)
≡
e
F
(
x
)
(
m
o
d
x
n
)
G(x)\equiv e^{F(x)}\,\,(mod\,\,x^{n})
G(x)≡eF(x)(modxn)
我们先对两边求
ln
\ln
ln,得到:
ln
G
(
x
)
≡
F
(
x
)
(
m
o
d
x
n
)
ln
G
(
x
)
−
F
(
x
)
≡
0
(
m
o
d
x
n
)
\ln G(x)\equiv F(x)\,\,(mod\,\,x^{n})\\ \ln G(x)-F(x)\equiv 0\,\,(mod\,\,x^{n})
lnG(x)≡F(x)(modxn)lnG(x)−F(x)≡0(modxn)
设在模
x
⌈
n
2
⌉
x^{\left\lceil\frac{n}{2}\right\rceil}
x⌈2n⌉意义下的答案为
H
(
x
)
H(x)
H(x),那么根据牛顿迭代公式有:
G
(
x
)
≡
H
(
x
)
−
ln
H
(
x
)
−
F
(
x
)
1
H
(
x
)
(
m
o
d
x
n
)
G
(
x
)
≡
H
(
x
)
(
1
−
ln
H
(
x
)
+
F
(
x
)
)
(
m
o
d
x
n
)
G(x)\equiv H(x)-\frac{\ln H(x)-F(x)}{\frac{1}{H(x)}}\,\,(mod\,\,x^{n})\\ G(x)\equiv H(x)\left(1-\ln H(x)+F(x)\right)\,\,(mod\,\,x^{n})
G(x)≡H(x)−H(x)1lnH(x)−F(x)(modxn)G(x)≡H(x)(1−lnH(x)+F(x))(modxn)
然后就可以快乐倍增了.
时间复杂度 O ( n log n ) O(n\log n) O(nlogn).
与 ln \ln ln相对的, exp \exp exp必须满足 F ( x ) F(x) F(x)常数项为 0 0 0, G ( x ) G(x) G(x)常数项为 1 1 1的条件.
代码如下:
int tb[N+9];
void Poly_exp(int *a,int n,int *res){
if (!n) {res[0]=1;return;}
Poly_exp(a,n>>1,res);
Poly_ln(res,n,tb);
tb[0]=sub(add(a[0],1),tb[0]);
for (int i=1;i<=n;++i) tb[i]=sub(a[i],tb[i]);
Get_len(n<<1);
NTT(res,len,0);
NTT(tb,len,0);
for (int i=0;i<len;++i) smul(res[i],tb[i]),tb[i]=0;
NTT(res,len,1);
for (int i=n+1;i<=n<<1;++i) res[i]=0;
}
六.多项式快速幂.
给定多项式
F
(
x
)
F(x)
F(x),求
G
(
x
)
G(x)
G(x)满足:
G
(
x
)
≡
F
k
(
x
)
(
m
o
d
x
n
)
G(x)\equiv F^{k}(x)\,\,(mod\,\,x^{n})
G(x)≡Fk(x)(modxn)
当然可以直接拿快速幂写,但是这样时间复杂度为 O ( n log n log k ) O(n\log n\log k) O(nlognlogk),然而还有更加优秀的做法.
我们考虑对其
ln
\ln
ln后
exp
\exp
exp,即:
G
(
x
)
≡
e
k
ln
F
(
x
)
(
m
o
d
x
n
)
G(x)\equiv e^{k\ln F(x)}\,\,(mod\,\,x^{n})
G(x)≡eklnF(x)(modxn)
问题是这个做法只能用在常数项为 1 1 1的情况.
对于常数项不为 1 1 1的情况,我们取其最低的不为 0 0 0的项,设为 a m x m a_mx^{m} amxm,那么直接将这个单项式作为因式提出来后再做快速幂,最后将 ( a m x m ) k (a_mx^{m})^{k} (amxm)k乘回去就好了.
时间复杂度降到了 O ( n log n ) O(n\log n) O(nlogn).
代码如下:
int tc[N+9];
void Poly_power(int *a,int n,int k,int *res){
if (!k) {res[0]=1;return;}
int m=0;
for (;!a[m]&&m<=n;++m);
if (m>n/k) return;
int mk=m*k,t0=Power(a[m],k),t1=Get_inv(a[m]);
for (int i=0;i<=n-mk;++i) tb[i]=mul(a[i+m],t1);
Poly_ln(tb,n-mk,tc);
for (int i=0;i<=n-mk;++i) smul(tc[i],k),tb[i]=0;
Poly_exp(tc,n-mk,res);
for (int i=n-mk;i>=0;--i) res[i+mk]=mul(res[i],t0),tc[i]=0;
for (int i=0;i<mk;++i) res[i]=0;
}
七.多项式开根.
给定多项式
F
(
x
)
F(x)
F(x),求
G
(
x
)
G(x)
G(x)满足:
G
2
(
x
)
≡
F
(
x
)
(
m
o
d
x
n
)
G^{2}(x)\equiv F(x)\,\,(mod\,\,x^{n})
G2(x)≡F(x)(modxn)
当然可以直接 ln \ln ln后 exp \exp exp,但是还有一种常数更加优秀的牛顿迭代做法.
直接进行移项,得到:
G
2
(
x
)
−
F
(
x
)
≡
0
(
m
o
d
x
n
)
G^{2}(x)-F(x)\equiv 0\,\,(mod\,\,x^{n})
G2(x)−F(x)≡0(modxn)
设在模
x
⌈
n
2
⌉
x^{\left\lceil \frac{n}{2}\right\rceil}
x⌈2n⌉意义下的答案为
H
(
x
)
H(x)
H(x),那么根据牛顿迭代公式有:
G
(
x
)
≡
H
(
x
)
−
H
2
(
x
)
−
F
(
x
)
2
H
(
x
)
(
m
o
d
x
n
)
G
(
x
)
≡
1
2
(
H
(
x
)
+
F
(
x
)
H
(
x
)
)
(
m
o
d
x
n
)
G(x)\equiv H(x)-\frac{H^{2}(x)-F(x)}{2H(x)}\,\,(mod\,\,x^{n})\\ G(x)\equiv \frac{1}{2}\left(H(x)+\frac{F(x)}{H(x)}\right)\,\,(mod\,\,x^{n})
G(x)≡H(x)−2H(x)H2(x)−F(x)(modxn)G(x)≡21(H(x)+H(x)F(x))(modxn)
至于常数项写个二次剩余就行了.
时间复杂度 O ( n log n ) O(n\log n) O(nlogn).
代码如下:
void Poly_sqrt(int *a,int n,int *res){
if (!n) {int t=Cipolla(a[0]);res[0]=min(t,mod-t);return;}
Poly_sqrt(a,n>>1,res);
Poly_inv(res,n,tb);
for (int i=0;i<=n;++i) ta[i]=a[i];
Get_len(n<<1);
NTT(ta,len,0);
NTT(tb,len,0);
for (int i=0;i<len;++i) smul(ta[i],tb[i]),tb[i]=0;
NTT(ta,len,1);
for (int i=0;i<=n;++i) res[i]=mul(inv2,add(res[i],ta[i])),ta[i]=0;
for (int i=n+1;i<=n<<1;++i) ta[i]=0;
}
八.多项式三角函数.
给定多项式
F
(
x
)
F(x)
F(x),求
G
(
x
)
G(x)
G(x)和
H
(
x
)
H(x)
H(x)满足:
G
(
x
)
≡
sin
F
(
x
)
(
m
o
d
x
n
)
H
(
x
)
≡
cos
F
(
x
)
(
m
o
d
x
n
)
G(x)\equiv \sin F(x)\,\,(mod\,\,x^{n})\\ H(x)\equiv \cos F(x)\,\,(mod\,\,x^{n})
G(x)≡sinF(x)(modxn)H(x)≡cosF(x)(modxn)
生成函数本质上是形式幂级数,所以将
e
x
,
sin
(
x
)
,
cos
(
x
)
e^{x},\sin(x),\cos(x)
ex,sin(x),cos(x)三者泰勒展开:
e
x
=
∑
i
=
0
+
∞
x
i
i
!
=
1
+
x
+
x
2
2
!
+
x
3
3
!
+
⋯
sin
(
x
)
=
∑
i
=
0
+
∞
(
−
1
)
i
(
2
i
+
1
)
!
x
2
i
+
1
=
x
−
x
3
3
!
+
x
5
5
!
−
x
7
7
!
+
⋯
cos
(
x
)
=
∑
i
=
0
+
∞
(
−
1
)
i
(
2
i
)
!
x
2
i
=
1
−
x
2
2
!
+
x
4
4
!
−
x
6
6
!
+
⋯
e^{x}=\sum_{i=0}^{+\infty}\frac{x^{i}}{i!}=1+x+\frac{x^{2}}{2!}+\frac{x^{3}}{3!}+\cdots\\ \sin(x)=\sum_{i=0}^{+\infty}\frac{(-1)^{i}}{(2i+1)!}x^{2i+1}=x-\frac{x^{3}}{3!}+\frac{x^{5}}{5!}-\frac{x^{7}}{7!}+\cdots\\ \cos(x)=\sum_{i=0}^{+\infty}\frac{(-1)^{i}}{(2i)!}x^{2i}=1-\frac{x^{2}}{2!}+\frac{x^{4}}{4!}-\frac{x^{6}}{6!}+\cdots
ex=i=0∑+∞i!xi=1+x+2!x2+3!x3+⋯sin(x)=i=0∑+∞(2i+1)!(−1)ix2i+1=x−3!x3+5!x5−7!x7+⋯cos(x)=i=0∑+∞(2i)!(−1)ix2i=1−2!x2+4!x4−6!x6+⋯
用这些东西可以推导出用三角函数表示
exp
\exp
exp的欧拉公式:
e
ϕ
i
=
cos
ϕ
+
i
sin
ϕ
e^{\phi i}=\cos \phi+i\sin \phi
eϕi=cosϕ+isinϕ
那么同样可以用这些东西推导出用
exp
\exp
exp表示三角函数的公式:
sin
ϕ
=
i
e
−
ϕ
i
−
e
ϕ
i
2
cos
ϕ
=
e
ϕ
i
+
e
−
ϕ
i
2
\sin \phi=i\frac{e^{-\phi i}-e^{\phi i}}{2}\\ \cos \phi=\frac{e^{\phi i}+e^{-\phi i}}{2}
sinϕ=i2e−ϕi−eϕicosϕ=2eϕi+e−ϕi
然后我们知道 i i i在模一个NTT模数 p p p意义下为 g p − 1 4 g^{\frac{p-1}{4}} g4p−1,其中 g g g为一个原根.
剩下就是套板子的事情了,时间复杂度 O ( n log n ) O(n\log n) O(nlogn).
注意在实现的时候,可以先求出 e ϕ i e^{\phi i} eϕi,再通过求逆得到 e − ϕ i e^{-\phi i} e−ϕi,比直接用调用两次 exp \exp exp常数小一些.
代码如下:
void Poly_sin(int *a,int n,int *res){
int im=wn[0][N>>2];
for (int i=0;i<=n;++i) res[i]=mul(im,a[i]);
Poly_exp(res,n,tc);
for (int i=0;i<=n;++i) res[i]=0;
Poly_inv(tc,n,res);
for (int i=0;i<=n;++i) res[i]=mul(mul(inv2,im),sub(res[i],tc[i])),tc[i]=0;
}
void Poly_cos(int *a,int n,int *res){
int im=wn[0][N>>2];
for (int i=0;i<=n;++i) res[i]=mul(im,a[i]);
Poly_exp(res,n,tc);
for (int i=0;i<=n;++i) res[i]=0;
Poly_inv(tc,n,res);
for (int i=0;i<=n;++i) res[i]=mul(inv2,add(res[i],tc[i])),tc[i]=0;
}