【学习笔记】多项式
前面的内容先鸽掉了，之后再补。
快速傅里叶变换 FFT
Code。
快速数论变换 NTT
NTT 模数原根表
Code。
任意模数多项式乘法
9NTT 做法
NTT+CRT。
Code。
5FFT 做法
Code。
快速莫比乌斯变换 FMT
快速沃尔什变换 FWT
多项式全家桶
前置：泰勒展开
如果函数足够光滑的话，在已知函数在某一点的各阶导数值的情况之下，泰勒公式可以用这些导数值做系数构建一个多项式来近似函数在这一点的邻域中的值。
在某一点使用泰勒公式即在该点（在下方公式中即为 \(x_0\)）处泰勒展开（\(f^{(k)}(x)\) 表示多项式 \(f(x)\) 的 \(n\) 阶导）：
多项式牛顿迭代
多项式牛顿迭代用于求多项式的零点，即对于函数 \(G(x)\)，求满足 \(G(F(z))\equiv 0 (\bmod\ x^n)\) 的多项式 \(F(z)\)。
考虑倍增。假设已经求出了
扩展到 \(\bmod\ x^n\) 意义下，将 \(G(F(z))\) 在 \(F_0(z)\) 处泰勒展开，有：
由于 \(F(z)\) 的前 \(\frac{n}{2}\) 项与 \(F_0(z)\) 相同，所以 \(F(z)-F_0(z)\) 的最低次项次数为 \(\frac{n}{2}\)，在 \(k\ge 2\) 时，\((F(z)-F_0(z))^k\) 的最低次项次数均大于等于 \(2\cdot \frac{n}{2}=n\)，所以该式在 \(\bmod\ x^n\) 意义下表示为：
又因为 \(G(F(z))\equiv 0 (\bmod\ x^n)\)，则
多项式求逆
已知多项式 \(A(x)\)，且 \(A(x)\times B(x) \equiv 1(\bmod\ x^n)\)（定义 \(\times\) 为多项式乘法），求多项式 \(B(x)\)。
手动模拟样例 \(A(x)=-x^2-x+1\)，发现 \(B^{-1}(x)=1+x+2x^2+3x^3+5x^4 \cdots\)，其系数为斐波那契数列且有无穷多项。对 \(x^n\) 取模相当于只保留前 \(n\) 项。
模拟求解时间复杂度 \(O(n^2)\) 较劣，考虑倍增。
假设已得出模 \(x^m\) 意义下的解 \(B_m(x)\)，有（\(B_m(x)\) 为 \(B_{2m}(x)\) 的前 \(m\) 项）：
两边平方得（模数可同时平方）：
两边同乘 \(A(x)\) 得：
此时就可以使用多项式乘法快速求解了。
每次倍增的时间复杂度是 \(O(m\log m)\)，总时间复杂度 \(O(n\log n)\)。
（也可用牛顿迭代推导。）
多项式带余除法
已知多项式 \(A(x),C(x)\) ，且 \(A(x)=B(x)\times C(x)+R(x)\)，求 \(C(x),R(x)\)。
已知 \(A(x)\) 最高次为 \(n\)，\(B(x)\) 最高次为 \(m\)，不难得出 \(C(x)\) 最高次为 \(n-m\)，\(R(x)\) 最高次为 \(m-1\)。
定义 \(f^T(x)=x^nf(x^{-1})\)，有：
两边同乘 \(x^n\) 得：
在模 \(x^{n-m+1}\) 意义下：
不难发现 \(f^T(x)\) 的本质是把多项式 \(f(x)\) 的系数翻转，所以 \(A^T(x),B^T(x)\) 易得，且：
\(C^T(x)\) 易求，则 \(C(x),R(x)\) 易求。
时间复杂度 \(O(n\log n)\)。
code。
多项式开根
已知多项式 \(A(x)\)，且 \(A(x)=B(x)^2\)，求 \(B(x)\)。
考虑倍增。
假设已得出模 \(x^m\) 意义下的解 \(B_m(x)\)，有（\(B_m(x)\) 为 \(B_{2m}(x)\) 的前 \(m\) 项）：
两边平方得（模数可同时平方）：
此时就可以使用多项式求逆快速求解了。
每次倍增的时间复杂度是 \(O(m\log m)\)，总时间复杂度 \(O(n\log n)\)。
（也可用牛顿迭代推导。）
code。
多项式对数函数 \(\ln\)
已知 \(F(x)\)，求 \(G(x)\equiv \ln F(x)(\bmod\ x^n)\)。
令函数 \(f(x)=\ln x\)，则 \(G(x)\equiv f(F(x))(\bmod\ x^n)\)。
已知复合函数求导公式：\(f(g(x))'=f'(g(x))g'(x)\)，则
又因为 \(\ln'(x)=\frac{1}{x}\)，则
先对 \(F(x)\) 的求导、求逆，二者相乘得到 \(G'\)，再对 \(G'\) 积分得到 \(G\) 即可。
求导公式：\({x^y}'=y\cdot x^{y-1}\)
积分公式：\(\displaystyle \int x^ydx =\frac{1}{y+1}x^{y+1}\)
多项式指数函数 \(\exp\)
已知 \(A(x)\)，求 \(F(x)\equiv e^{A(x)}(\bmod\ x^n)\)。
对原式变形得
设
即要求多项式 \(G(F(x))\) 的零点。
代入牛顿迭代可得
用多项式 \(\ln\)、乘法、求导等解决即可。
时间复杂度 \(O(n\log n)\)。
多项式指数快速幂
已知 \(F(x)\)，求 \(G(x)\equiv F(x)^k(\bmod\ x^n)\)。
取对数
取指数
用多项式 \(\ln\) 和 \(exp\) 快速求解即可。
时间复杂度 \(O(n\log n)\)。
由于 \(\exp\) 模板里保证 \(A_0=1\)，可以直接求解，但在更一般化的情况下 \(A_0\) 可能不等于 \(1\) 甚至等于 \(0\)。在加强版中需要找到序列中第一个非 \(0\) 数，对多项式降次并把所有系数乘上该非 \(0\) 数的逆元，以此把首项化为 \(1\)。最后将次数升回并将所有系数乘上原 \(A_0\) 的 \(k\) 次方。需要注意的是，这里的 \(k\) 需要对 \(\varphi(p)\) 取模。
代码
终于卡出了一份常数还可以的代码。内容包含求逆、\(\ln\)、\(\exp\)。
之后找时间把开根和除法也卡好常加进去。
#include<bits/stdc++.h>
using namespace std;
#define add(a,b) ((a+b>=mod)?(a+b-mod):(a+b))
#define dec(a,b) ((a-b<0)?(a-b+mod):(a-b))
#define mul(a,b) (1ll*a*b%mod)
#define swap(a,b) (a^=b,b^=a,a^=b)
namespace FastIO{
const int SIZE=1<<16;
char buf[SIZE],obuf[SIZE],str[60];
int bi=SIZE,bn=SIZE,opt;
inline int read(register char *s){
while(bn){
for(;bi<bn&&buf[bi]<=' ';bi=-~bi);
if(bi<bn) break;
bn=fread(buf,1,SIZE,stdin);
bi&=0;
}
register int sn=0;
while(bn){
for(;bi<bn&&buf[bi]>' ';bi=-~bi) s[sn++]=buf[bi];
if(bi<bn) break;
bn=fread(buf,1,SIZE,stdin);
bi&=0;
}
s[sn]&=0;
return sn;
}
inline bool read(register int &x){
int n=read(str),bf=0;
if(!n) return 0;
register int i=0;
(str[i]=='-')&&(bf=1,i=-~i);
(str[i]=='+')&&(i=-~i);
for(x=0;i<n;i=-~i) x=(x<<3)+(x<<1)+(str[i]^48);
bf&&(x=~x+1);
return 1;
}
inline bool read(register long long &x){
int n=read(str),bf=1;
if(!n) return 0;
register int i=0;
(str[i]=='-')&&(bf=-1,i=-~i);
for(x=0;i<n;i=-~i) x=(x<<3)+(x<<1)+(str[i]^48);
(bf<0)&&(x=~x+1);
return 1;
}
inline void write(register int x){
if(!x) obuf[opt++]='0';
else{
(x<0)&&(obuf[opt++]='-',x=~x+1);
register int sn=0;
while(x) str[sn++]=x%10+'0',x/=10;
for(register int i=sn-1;i>=0;i=~-i) obuf[opt++]=str[i];
}
(opt>=(SIZE>>1))&&(fwrite(obuf,1,opt,stdout),opt&=0);
}
inline void write(register long long x){
if(!x) obuf[opt++]='0';
else{
(x<0)&&(obuf[opt++]='-',x=~x+1);
register int sn=0;
while(x) str[sn++]=x%10+'0',x/=10;
for(register int i=sn-1;i>=0;i=~-i) obuf[opt++]=str[i];
}
(opt>=(SIZE>>1))&&(fwrite(obuf,1,opt,stdout),opt&=0);
}
inline void write(register unsigned long long x){
if(!x) obuf[opt++]='0';
else{
register int sn=0;
while(x) str[sn++]=x%10+'0',x/=10;
for(register int i=sn-1;i>=0;i=~-i) obuf[opt++]=str[i];
}
(opt>=(SIZE>>1))&&(fwrite(obuf,1,opt,stdout),opt&=0);
}
inline void write(register char x){
obuf[opt++]=x;
(opt>=(SIZE>>1))&&(fwrite(obuf,1,opt,stdout),opt&=0);
}
inline void Fflush(){
opt&&fwrite(obuf,1,opt,stdout);
opt&=0;
}
}
#define FI(n) FastIO::read(n)
#define FO(n) FastIO::write(n)
#define Flush FastIO::Fflush()
const int N=4e5+10,mod=998244353;
inline int _pow(int a,int b){
int ans=1;
while(b){
if(b&1) ans=mul(ans,a);
a=mul(a,a),b>>=1;
}
return ans;
}
namespace Poly{
const int N=4e5+10,G3=3,invG=332748118,MAXN=4e5;
int inv[N],gw[N],f[N],g[N],F[N],G[N],A[N],B[N];
inline void init(){
inv[0]=inv[1]=1;
for(int i=2;i<=MAXN;i++) inv[i]=mul((mod-mod/i),inv[mod%i]);
return ;
}
inline void NTT(int *a,int n,int type){
for(register int i=0,j=0;i<n;i++){
if(i<j) swap(a[i],a[j]);
for(register int k=n>>1;(j^=k)<k;k>>=1) ;
}
for(register int w=2,k=1,r=1;w<=n;w<<=1,k<<=1,r++){
int gn=_pow(type?G3:invG,(mod-1)/w);gw[0]=1;
for(int i=1;i<k;i++) gw[i]=mul(gw[i-1],gn);
for(register int i=0;i<n;i+=w)
for(register int j=0;j<k;j++){
register int x=a[i+j],y=mul(a[i+j+k],gw[j]);
a[i+j]=add(x,y),a[i+j+k]=dec(x,y);
}
}
if(!type)
for(int i=0;i<n;i++) a[i]=mul(a[i],inv[n]);
return ;
}
inline void MUL(int *a,int *b,int *c,int n){
for(int i=0;i<n;i++) f[i]=a[i],g[i]=b[i];
int len=1;while(len<(n<<1)) len<<=1;
NTT(f,len,1),NTT(g,len,1);
for(int i=0;i<len;i++) f[i]=mul(f[i],g[i]);
NTT(f,len,0);
for(int i=0;i<n;i++) c[i]=f[i];
for(int i=0;i<len;i++) f[i]=g[i]=0;
return ;
}
inline void inversion(int *a,int *b,int n){
if(n==1){b[0]=_pow(a[0],mod-2);return ;}
inversion(a,b,(n+1)>>1);
for(int i=0;i<n;i++) f[i]=a[i],g[i]=b[i];
int len=1;while(len<(n<<1)) len<<=1;
NTT(f,len,1),NTT(g,len,1);
for(int i=0;i<len;i++) f[i]=mul(f[i],mul(g[i],g[i]));
NTT(f,len,0);
for(int i=0;i<n;i++) b[i]=dec(add(b[i],b[i]),f[i]);
for(int i=0;i<len;i++) f[i]=g[i]=0;
return ;
}
inline void derivation(int *a,int *b,int n){
for(int i=1;i<n;i++) b[i-1]=mul(i,a[i]);
b[n-1]=0;
return ;
}
inline void integration(int *a,int *b,int n){
for(int i=1;i<n;i++) b[i]=mul(inv[i],a[i-1]);
b[0]=0;return ;
}
inline void ln(int *a,int *b,int n){
derivation(a,F,n),inversion(a,G,n);
int len=1;while(len<(n<<1)) len<<=1;
NTT(F,len,1),NTT(G,len,1);
for(int i=0;i<len;i++) F[i]=mul(F[i],G[i]);
NTT(F,len,0),integration(F,b,n);
for(int i=0;i<len;i++) F[i]=G[i]=0;
return ;
}
inline void exp(int *a,int *b,int n){
if(n==1){b[0]=1;return ;}
exp(a,b,(n+1)>>1),ln(b,A,n),A[0]=dec(1+a[0],A[0]);
for(int i=1;i<n;i++) A[i]=dec(a[i],A[i]);
int len=1;while(len<(n<<1)) len<<=1;
NTT(A,len,1),NTT(b,len,1);
for(int i=0;i<len;i++) b[i]=mul(b[i],A[i]);
NTT(b,len,0);
for(int i=n;i<len;i++) A[i]=b[i]=0;
return ;
}
}
int n;
int f[N],g[N];
int main(){
FI(n),Poly::init();
for(int i=0;i<n;i++) FI(f[i]);
Poly::exp(f,g,n);
for(int i=0;i<n;i++) FO(g[i]),FO(' ');
Flush;
return 0;
}
应用
常系数齐次线性递推
把线性递推式写作 \(\sum_{i=0}^{k-1} f_ix^i\)，需要求解的答案即为 \(x^n\)。
由 \(x^k=\sum_{i=0}^{k-1} f_ix^i\)，我们可以把 \(x^n\) 降次，最后可以得到一个形如 \(\sum_{i=0}^{k-1} f'_ix^i\) 的多项式，就可以直接求解。
不难发现每次会从原多项式中减去 \(a\) 倍的 \(x^k\) 并替换为 \(a\) 倍的 \(\sum_{i=0}^{k-1} f_ix^i\)，本质上等于多项式取模的过程。即对 \(x^k-\sum_{i=0}^{k-1} f_ix^i\) 取模，每减去一个 \(x^k\) 就会加上一个 \(\sum_{i=0}^{k-1} f_ix^i\)。
快速幂同时多项式取模计算即可。
时间复杂度 \(O(n\log^2n)\)。
这一篇题解讲得详细易懂。为节约时间这里就不再详解了。
任意无向连通图数量
设 \(f(n)\) 表示 \(n\) 个点无向连通图的数量，\(g(n)\) 表示 \(n\) 个点连通图的数量。
显然，有
（钦定 \(1\) 在 \(f\) 所代表的的连通块中）
将第一个式子代入第二个得
这是一个卷积形式。写出生成函数：
则
套一个多项式求逆板子即可。
时间复杂度 \(O(n\log n)\)。
浙公网安备 33010602011771号