浅谈多项式 & 生成函数科技（一）

快速数论变换 · 2025-10-11 15:16:52 · 算法·理论

前言

感觉没有什么能说的。

也许是新式颓废方式？

前排广告：

• 从斐波那契通项公式到卡特兰数通项公式
• 莫比乌斯反演学习笔记
• 基于 A 思想的 SPFA 优化-SPFA（已被 hack）（怎么被当成题解误杀了）

引子

P10780 BZOJ3028 食物

P4841 [集训队作业2013] 城市规划

很容易看出来这是两道计数题，但是怎么计数呢？不知道。

:::align{center} :::

这个时候就需要多项式与生成函数闪亮登场了。

大部分人应该知道多项式，但对于生成函数却很陌生。没关系，接下来我会讲解这些科技。

章一·初识生成函数

生成函数，可以简单理解为一个有着无穷项的多项式，一般可以记作 F(x)=\sum\limits_{i\ge0}f_ix^i。

（接下来将省略 F(x) 的 (x)）

其中，f_i 是第 i 项的系数，但是，第 i 项系数更多情况下写为 [x^i]F。

然而，我们有一些数并不会用到，例如 [x^{1145141919810}]F，而且由于空间限制，我们一般也存不下。

于是，便引入了“截断”操作。我们通过舍弃 F(x) 中 x^{n} 及更高次的所有项，得到一个项数不多于 n 项的多项式。

将生成函数 F 从 x^n 处截断记为 F\pmod {x^n}。一定要注意，这是截断，不要当成取模 （虽然与取模差不了多少）！

截断后的生成函数，就是我们日常常见的多项式了，写作 \sum\limits_{i=0}^{n-1}f_ix^i。

而多项式可以像整数一样进行加减乘除与取模，还可以求 \exp 与 \ln，甚至可以求其（反）三角函数多项式，这些我们后面再讲。

生成函数又分为两类：

• 普通型生成函数（OGF）
  • 通常写作 F=\sum\limits_{i\ge0}f_ix^i。
• 指数型生成函数（EGF）
  • 通常写作 F=\sum\limits_{i\ge0}f_i\dfrac{x^i}{i!}。

OGF 通常用来解决无标号计数问题，EGF 通常用来解决有标号计数问题。

而这两种计数问题又指的是什么呢？

章二·两类计数问题

一、无标号计数问题

引子中的第一道题食物就是一道经典的无编号计数。

但何为无标号计数？

无标号计数，利用群论知识，可简单理解为有多少个等价类。

但是这样表达太过严谨，我用一个例子来告诉你是什么：

例如，你的同一份代码，在洛谷上提交了两遍，一次由于测评机波动被卡常导致 TLE 了，而另一次刚好在一个优秀的时间提交，极限 AC；那么如果从代码的角度看，这两次提交测评时相同的；但如果从结果的角度来看，就是不同的。

没错，这个例子就是我。

在这个例子中，从代码的角度来看，就是无标号计数；而从测评结果、提交顺序来看，就是有标号计数。

好了，你已经理解了这两者的内涵，接下来讲解遇到无标号计数时，如何求解。

我们以食物为例。

一共带 n 件物品，要求如下：

• 承德汉堡：偶数个，则它的生成函数为 A=1+x^2+x^4+x^6+\cdots；
• 可乐：0 或 1 个，则它的生成函数为 B=1+x；
• 鸡腿：0、1 或 2 个，则它的生成函数为 C=1+x+x^2；
• 蜜桃多：奇数个，则为 D=x+x^3+x^5+x^7+\cdots；
• 鸡块：4 的倍数个，则为 E=1+x^4+x^8+x^{12}+\cdots；
• 包子：0、1、2 或 3 个，则为 F=1+x+x^2+x^3；
• 土豆片炒肉：不超过一个，则为 G=1+x；
• 面包：3 的倍数个，则为 H=1+x^3+x^6+x^9+\cdots。

然后，只要将这 8 个生成函数相乘得到 K，然后 [x^n]K 就能得到答案了。

但是，我并没有讲如何求多项式的乘积，一些大蛇口中的答案可能马上就要呼之欲出。但是，这题根本用不到多项式乘法。

让我们对 A,D,E,H 进行简单的数学推导，得到它们的封闭形式：

然后再将这八个相乘，得到 K=\dfrac{x}{(x-1)^4}。

但这个时候就不知道如何取 [x^n]K 了，怎么办呢？

此时就要请出广义二项式定理了！

广义二项式定理

公式如下：

其中 k\in\N,\alpha\in\R。

然后呢？这跟我们的式子有什么关系吗？

有的，兄弟有的。K=\dfrac{x}{(x-1)^4}=x(x-1)^{-4}，接下来利用广义二项式展开就可以了。

最后得到 K=\sum\limits_{k=1}^\infty\dbinom{k+2}{k-1}x^k=\sum\limits_{k=1}^\infty\dfrac{k(k+1)(k+2)}{6}x^k。

则 [x^n]K=\dfrac{n(n+1)(n+2)}{6}。

从这道题可以看出，解决无标号计数的方法如下：

1. 构造，使得对于每件物品，其生成函数系数不为 0 的项指数符合物品数目要求。
2. 推导，将无穷形式化为收敛形式。
3. 相乘，由多项式的意义可得，相乘表示选择。
4. 取系数，取得 [x^n]F，作为答案。

二、有标号计数问题

你说得对，但是城市规划并不是一道经典的有标号计数问题，主要是找不到能直接推式子的原题了 （更主要的原因是因为不想手写题面）。

但是，在此之前，我们要先解释：为什么 EGF 每一项会多一个 \dfrac{1}{i!}。只有理解了这个，遇到有标号计数时才能更得心应手一点。

我们在处理有标号计数问题时，常常会携带一些组合意义，如果使用 OGF，会出现组合数，乘法时的转移常常会变成这样：F(i)\cdot G(n-i)\cdot \dbinom{n}{i}\to H(n)，处理时会变得麻烦。此时使用 EGF，分母的阶乘就会自动的帮你处理好这个组合数的问题。

然后，我们来解决这个题。

我们不妨设 f_i 为 i 个点时，有标号无向连通图的个数，而 g_i 为 i 个点时无向图的个数。显然，g_i=2^{\small\dbinom{i}{2}}。

此时，设 F 为 \langle f_i\rangle 的生成函数，G 为 \langle g_i\rangle 的生成函数。

接下来，试图找到 F 与 G 的关系。

我们枚举 i，表示有 i 个连通块，将 F 连乘 i 次，以表示此时方案数。但是要给 F^i 乘上 \dfrac{1}{i!}，来表示连通块的排列是任意的。

最终得到 G=\sum\dfrac{F^i}{i!}，发现这就是 e^F 的泰勒展开，于是 G=e^F，两边求 \ln 就好了。

由于 G 是好求的，所以直接上多项式 \ln 模板。

但是大部分人应该不会多项式全家桶，所以知道做法即可，马上会讲到这些东西。

章三·多项式并不全家桶

提示：

1. 这部分篇幅可能稍长，如果你会多项式全家桶了，可以快进到后面。
2. 这里仅给出模数为 998244353 的 NTT 全家桶。

::::warning[警告]{open} 本人 NTT 常数略大！详见记录：https://www.luogu.com.cn/record/240720129 。 ::::

我们定义多项式最高次幂次数（即多项式的阶）+1 定义为多项式的长度。

令 n 为 F 的长度，m 为 G 的长度。

一、多项式加减

最简单的东西了，根据小学二年级所学过的知识，能得到 H=F\pm G\Longleftrightarrow [x^n]H=[x^n]F\pm[x^n]G。

时间复杂度：O(\max(n,m))。

二、多项式乘法

求 H=F\times G 中 x^n 的系数，由小学二年级知识可得 [x^n]H=\sum\limits_{i=0}^n[x^i]F\cdot[x^{n-i}]G。

但是，这样子求复杂度显然是 O(nm) 的。

如果我们局限于这样普通的乘法显然是没有办法进行优化的，于是考虑一些其他做法。

我们可以从拉格朗日插值中得到启发：l+1 个点可唯一确定一个最高次项为 l 次的多项式。

所以，我们可以考虑对于每个多项式，求出 n+m 个点的值。但是，求点值是 O\!\left((n+m)^2\right) 的，求点值积是 O(n+m) 的，朴素插值是 O\!\left((n+m)^2\right) 的，总的复杂度就是 O\!\left((n+m)^2\right) 的。

真的是太棒了！恭喜我们成功发现了比暴力还慢的做法！

那接下来该如何优化呢？

如果我们是随便选点的话，那么复杂度肯定降不下来。但是如果我们选择的点有一些优美的性质，那就说不准了。

想一下，我们求点值的瓶颈在哪里？是不是就是因为我们已知的信息难以重复利用导致更多的开销？！我们看看对称轴在 y 轴上的二次函数，当你求出 f(x) 时，是不是同时也能知道 f(-x) 的值？

容易发现，这类任意一项次数为偶数的函数，都是偶函数，都有这样优美的性质；同样发现任意一项次数为奇数的函数是奇函数，也拥有这样优美的性质。

但是，当这两者混合起来时，这优美的性质就荡然无存，于是，我们要先对函数做一件事：奇偶分离。即奇函数部分放到一边，偶函数部分放到另一边。

这样子，奇函数部分对称中心永远是 (0,0)，偶函数部分对称轴永远是 x=0。且只要将奇函数部分提取一个 x 出来，就又变成偶函数了。

这样我们就只需要求 \dfrac{(n+m)}{2} 个点了！

但是，复杂度好像并没有任何优化……

瓶颈仍然是在于求值与插值。

注意到我们奇偶分离后，若令 t=x^2，就会发现偶函数部分变成了一个普通关于 t 的函数，奇函数部分先提出一个 x 后，也变成了一个关于 t 的普通函数。

那样的话，我们只需要对奇偶两部分再不断进行换元与奇偶分离，就变成了可爱的分治，然后递归解决，复杂度就会降到 O\!\left((m+n)\log(m+n)\right)。

由于需要一直分，除递归中最后一层外每一层必须是偶数项，所以要将两个多项式补齐到 2^{\lceil\log_2(n+m)\rceil} 项。

但是插值的复杂度仍没有变，而且求点值时会出现一点问题。

因为我们会用到 x=\pm\sqrt{t}，已经默认了 t\ge0，但是这样子就没有相反数供我们在递归过程中优化求点值了！

那该怎么办呢？

根据《数学·人教 B 版·必修四》的知识，我们发现复数就可以很好的取代我们现在所用的实数，因为它可以保证平方后能存在负数，而且两个共轭复数的平方也互为共轭。

什么？你没学复数？该罚！

复数

定义与形式

复数，常用 z 表示，一般记为 z=a+b\text{i}，其中，a,b\in\R，\text{i} 为虚数单位，定义为 \text{i}^2=-1。

复数并不只有这一种形式，它还有三角形式，写法为：z=\rho(\cos\theta+\text{i}\sin\theta)。

由欧拉公式可得：z=\rho\exp(\text{i}\theta)。

这两者是可以转化的，有如下公式：

主播主播，e 指数上带了个 \text{i} 是什么东西啊？

提出这样问题的人，推荐去观看 3b1b 的视频，也推荐大家关注这位 UP。

复数部分到此为止，回归主线。

我们选择了复数来求点值，但是要是随便选的话，显然也无法进行简单的求解与后续的优化。

为了能够使得求解更简单，我们不妨令递归过程中最后一层所求的点值为 x=1 的值。

这样子，我们取的点只要满足 z^n=1 就可以了。

这个方程显然是有 n 个解的。

这 n 个解，称为 n 次单位根，记为 w。

借一张大佬的图。

:::align{center}

上图是 16 次单位根在复平面的表示 :::

单位根具有对称性与周期性，在这里我们不多做证明，读者自证不难。

另外，单位根有一个引理：折半引理，w^{2k}_{2n}=w^k_n。同样，读者自证不难。

那该怎么求呢？

显然，这 n 个点将单位圆平均的分成了 n 份，如图中逆时针排布的话，则第 k 个点与 x 轴的夹角为 \dfrac{2k\pi}{n}。

则由复数的三角形式可得，第 k 个 n 次单位根 w_n^k=\cos\frac{2k\pi}{n}+\text{i}\sin\frac{2k\pi}{n}。

讲了这么多，我们会发现我们仅仅优化了求点值，而插值仍然是朴素的 O\!\left((n+m)^2\right)。

我们希望插值也可以利用上这些性质优美的数，将复杂度降下来。但该怎么将他们统一为一种相似的形式呢？

我们再求完点值之后，能得到这样一组等式：

不妨将其写为矩阵形式，如下：

由于我们使用了单位根，所以矩阵最终形式如下

然后我们要求矩阵 \begin{bmatrix}h_0\\h_1\\h_2\\\vdots\\h_s\end{bmatrix}。

则有：

既然有着相似的形式，那么插值就也可以跟求点值一样方便的求出来。

这种多项式乘法，有一名，曰：快速傅里叶变换。

但实际上，将系数表示转化为点值表示才是傅里叶变换，全名离散傅里叶变换，简写为 DFT，将其变换回去是其逆变换，全名逆离散傅里叶变换，简写为 IDFT。

模板

#include<bits/stdc++.h>
#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/tree_policy.hpp>
#include<ext/pb_ds/hash_policy.hpp>
#include<ext/rope>
using namespace std;
using namespace __gnu_pbds;
using namespace __gnu_cxx;
//#define int long long
#define FileIn(file_name) freopen(#file_name".in","r",stdin)
#define FileOut(file_name) freopen(#file_name".out","w",stdout)
#define File(file_name) FileIn(file_name),FileOut(fail_name)
//tree< ,null_type,less< >,rb_tree_tag,tree_order_statistics_node_update> tr;
//gp_hash_table< , > ht;
constexpr int N=6e6+6;
constexpr double PI=acos(-1.0);
int rev[N];
inline void init(int len){
    int b=log2(len);
    for(int i=1;i<len;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline void FFT(complex<double>* a,int n,int type=1){
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=1;len<=(n>>1);len<<=1){
        complex<double> wn={cos(PI/len),type*sin(PI/len)};
        for(int i=0;i<n;i+=(len<<1)){
            complex<double> w={1,0};
            for(int j=0;j<len;++j,w*=wn){
                complex<double> a0=a[i+j],a1=w*a[i+j+len];
                a[i+j]=a0+a1,a[i+j+len]=a0-a1;
            }
        }
    }
    return;
}
inline void IFFT(complex<double>*a,int n){
    FFT(a,n,-1);
    for(int i=0;i<n;++i) a[i]={a[i].real(),a[i].imag()/2/n+0.5};
    return;
}
complex<double> a[N];
signed main(){
    cin.tie(nullptr)->sync_with_stdio(false);
    int n,m,len=1;cin>>n>>m;
    for(int i=0,k;i<=n;++i) cin>>k,a[i]={k,0};
    for(int i=0,k;i<=m;++i) cin>>k,a[i]={a[i].real(),k};
    while(len<=n+m) len<<=1;
    init(len);
    FFT(a,len);
    for(int i=0;i<len;++i) a[i]*=a[i];
    IFFT(a,len);
    for(int i=0;i<=n+m;++i) cout<<(int)(a[i].imag())<<' ';
    return 0;
}

注意到我只做了 1 次 FFT 与 1 次 IFFT 还用到了 rev 等并没有提到的东西。

首先讲这个 rev 是怎么回事。

这是 FFT 的位逆序置换，随着向下递归，系数会这样划分，如下所示：

我们比较第一层与最后一层的 a 下标的二进制：

那我们就可以与预处理出其变换，然后直接调换其顺序，能够再缩小常数的同时保证正确性。

那程序中的递推为什么是正确的呢？

分类讨论，先考虑偶数的情况。此时有 rev[i<<1]=rev[i]>>1;，则 rev[i]=rev[i>>1]>>1;；再考虑奇数的情况。此时有 rev[i<<1|1]=(rev[i]>>1)+1;，则 rev[i]=(rev[i>>1]>>1)+(1<<(b-1));。

总结一下，浓缩成 rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));

rev 的意义搞定，那我为什么只需要进行一次正变换，一次逆变换就可以了呢？

这是 FFT 中的一个能大幅减少常数的操作：三次便两次优化。

那它的原理是怎么回事呢？

注意到我们是将多项式以复数存储，实部为 F，虚部为 G。

而对于一个复数 z=a+b\text{i}，z^2=(a^2-b^2)+2ab\text{i}。

那么虚部就是我们要求的东西了，不过要再除以二，而且原来的 \dfrac{1}{n} 不能丢。

而 NTT 与 FFT 相差不大，仅仅是将单位根换为了原根。

附上 NTT 代码（已经经过了卡常，可放心食用）：

#include<bits/stdc++.h>
using namespace std;
namespace IO{
    static constexpr int BUFSIZE=1<<21;
    char ibuf[BUFSIZE],*is=ibuf,*it=ibuf,obuf[BUFSIZE];
    static int cnt=0,top,wr[50];
    inline void flush(){
        fwrite(obuf,1,cnt,stdout),cnt=0;
        return;
    }
    inline char get(){
        if(is==it) it=(is=ibuf)+fread(ibuf,1,BUFSIZE,stdin);
        return is==it? EOF:*is++;
    }
    inline void put(char c){
        obuf[cnt++]=c;
        if(cnt==BUFSIZE) flush();
        return;
    }
    struct AutoFlush{~AutoFlush(){flush();}}flusher;
    inline int read(){
        char c=get();int x=0,neg=0;
        while(!isdigit(c)){
            if(c=='-') neg^=1;
            c=get();
        }
        do x=(x<<1)+(x<<3)+(c&15); while(isdigit(c=get()));
        return neg? -x:x;
    }
    template<typename Tp>
    inline void write(Tp x,const char& c='\0'){
        if(x<0) put('-'),x=-x;
        do wr[++top]=x%10; while(x/=10);
        while(top) put(wr[top--]|48);
        if(c) put(c);
        return;
    }
}
using IO::read;
using IO::write;
constexpr int N=6e6+6,p=998244353;
int rev[N];
inline void init(int len){
    int b=__lg(len);
    for(int i=1;i<len;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline int qpw(int a,int T){
    int res=1;
    for(;T;T>>=1,a=1ll*a*a%p) if(T&1) res=1ll*res*a%p;
    return res;
}
inline void NTT(int* a,int n,int g){
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=2,wn,l,a0,a1;len<=n;len<<=1){
        l=len>>1,wn=qpw(g,(p-1)/len);
        for(int i=0;i<n;i+=len){
            for(int j=0,w=1;j<l;++j,w=1ll*w*wn%p){
                a0=a[i+j],a1=1ll*w*a[i+j+l]%p;
                a[i+j]=(0ll+a0+a1)%p,a[i+j+l]=((a0-a1)%p+p)%p;
            }
        }
    }
    return;
}
inline void INTT(int* a,int n,int l1,int l2){
    NTT(a,n,332748118);
    for(int i=0,invn=qpw(n,p-2);i<=l1+l2;++i) a[i]=1ll*a[i]*invn%p;
    return;
}
int a[N],b[N],n,m,len=1;
signed main(){
    n=read(),m=read();
    while(len<=n+m) len<<=1;
    for(int i=0;i<=n;++i) a[i]=read();
    for(int i=0;i<=m;++i) b[i]=read();
    init(len);
    NTT(a,len,3),NTT(b,len,3);
    for(int i=0;i<len;++i) a[i]=1ll*a[i]*b[i]%p;
    INTT(a,len,n,m);
    for(int i=0;i<=n+m;++i) write(a[i],' ');
    return 0;
}

三、多项式乘法逆

注意：

1. 本段已经默认存在逆，若 [x^0]F=0 则不存在逆；
2. 本段多项式乘法逆都在 \!\!\!\pmod{x^n} 意义下，否则许多的多项式乘法逆都是无穷项。

这里讲解一个较常用的倍增法。

很显然的一个东西就是 [x^0]F^{-1}=([x^0]F)^{-1}。

假设我们现在已经求出了 F_0^{-1}\pmod{x^{\lceil\frac{n}{2}\rceil}}，则有：

将 [3] 式平方得：

再同时乘以 F，并移项，可得：

发现可以递归求解，时间复杂度是 O(n\log n) 的。

模板

/*我们所可以自慰的，想来想去，也还是所谓对于将来的希望。

希望是附丽于存在的，有存在，便有希望，有希望，便是光明。*/
#include<bits/stdc++.h>
using namespace std;
constexpr int N=5e5+5,p=998244353;
int rev[N];
inline void init(int len){
    int b=__lg(len);
    for(int i=1;i<len;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline int qpw(int a,int T){
    int res=1;
    for(;T;T>>=1,a=1ll*a*a%p) if(T&1) res=1ll*res*a%p;
    return res;
}
inline void NTT(int* a,int n,int g){
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=2,wn,l,a0,a1;len<=n;len<<=1){
        l=len>>1,wn=qpw(g,(p-1)/len);
        for(int i=0;i<n;i+=len){
            for(int j=0,w=1;j<l;++j,w=1ll*w*wn%p){
                a0=a[i+j],a1=1ll*w*a[i+j+l]%p;
                a[i+j]=(0ll+a0+a1)%p,a[i+j+l]=((a0-a1)%p+p)%p;
            }
        }
    }
    return;
}
inline void INTT(int* a,int n,int l){
    NTT(a,n,332748118);
    for(int i=0,invn=qpw(n,p-2);i<=l;++i) a[i]=1ll*a[i]*invn%p;
    return;
}
int a[N],inv[N],f[N],n,len=1;
signed main(){
    cin.tie(nullptr)->sync_with_stdio(false);
    cin>>n;while(len<=n) len<<=1;
    for(int i=0;i<n;++i) cin>>a[i];
    f[0]=qpw(a[0],p-2);
    for(int t=2,k;t<=len;t<<=1){
        k=t<<1;init(k);
        copy(a,a+t,inv),fill(inv+t,inv+k,0);
        NTT(f,k,3),NTT(inv,k,3);
        for(int i=0;i<k;++i) f[i]=1ll*f[i]*((2-1ll*inv[i]*f[i]%p)%p+p)%p;
        INTT(f,k,t);
        fill(f+t,f+k,0);
    }
    for(int i=0;i<n;++i) cout<<f[i]<<' ';
    return 0;
}

四、多项式 \ln & \exp

多项式 \ln

首先，多项式 \ln 要保证 [x^0]F=1，否则常数项不收敛，但如果你是 double 玩家的话不用管它。

首先对 \ln F 求导然后再积分回去，显然结果不会改变。

在这过程中，我们可以得到：

求导与积分显然是 O(n) 的，求逆是 O(n\log n)，总复杂度就是 O(n\log n)。

模板

/*我们所可以自慰的，想来想去，也还是所谓对于将来的希望。

希望是附丽于存在的，有存在，便有希望，有希望，便是光明。*/
#include<bits/stdc++.h>
using namespace std;
constexpr int N=5e5+5,p=998244353;
int rev[N],len;
inline void init(int n){
    int b=__lg(n);
    for(int i=1;i<n;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline int qpw(int a,int T){
    int res=1;
    for(;T;T>>=1,a=1ll*a*a%p) if(T&1) res=1ll*res*a%p;
    return res;
}
inline void NTT(int* a,int n,int g){
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=2,wn,l,a0,a1;len<=n;len<<=1){
        l=len>>1,wn=qpw(g,(p-1)/len);
        for(int i=0;i<n;i+=len){
            for(int j=0,w=1;j<l;++j,w=1ll*w*wn%p){
                a0=a[i+j],a1=1ll*w*a[i+j+l]%p;
                a[i+j]=(0ll+a0+a1)%p,a[i+j+l]=((a0-a1)%p+p)%p;
            }
        }
    }
    return;
}
inline void INTT(int* a,int n){
    NTT(a,n,332748118);
    for(int i=0,invn=qpw(n,p-2);i<n;++i) a[i]=1ll*a[i]*invn%p;
    return;
}
int a[N],n,f[N],g[N],h[N],inv[N],fac[N];
inline void derivative(int* a,int* res,int n){
    for(int i=1;i<n;++i) res[i-1]=1ll*a[i]*i%p;
    return res[n-1]=0,void();
}
inline void integrate(int* a,int* res,int n){
    for(int i=1;i<n;++i) res[i]=1ll*a[i-1]*inv[i]%p;
    return res[0]=0,void();
}
inline void polyinv(int* a,int* f,int n){
    static int temp[N];
    fill(f,f+(n<<1),0);
    f[0]=qpw(a[0],p-2);
    for(int t=2;t<=(n<<1);t<<=1){
        int k=t<<1;
        init(k);
        copy(a,a+t,temp);
        fill(temp+t,temp+k,0);
        NTT(f,k,3);NTT(temp,k,3);
        for(int i=0;i<k;++i) 
            f[i]=1ll*f[i]*((2-1ll*temp[i]*f[i]%p)%p+p)%p;
        INTT(f,k);
        fill(f+t,f+k,0);
    }
    return;
}
signed main(){
    cin.tie(nullptr)->sync_with_stdio(false);
    cin>>n;
    for(int i=0;i<n;++i) cin>>a[i];
    int t=inv[1]=1;
    for(int i=2;i<=n;++i) inv[i]=1ll*(p-p/i)*inv[p%i]%p;
    derivative(a,g,n);
    polyinv(a,f,n);
    while(t<(n<<1)) t<<=1;
    init(t);
    NTT(f,t,3);NTT(g,t,3);
    for(int i=0;i<t;++i) f[i]=1ll*f[i]*g[i]%p;
    INTT(f,t);
    integrate(f,g,n);
    for(int i=0;i<n;++i) cout<<g[i]<<' ';
    return 0;
}

多项式 \exp

与 \ln 一样，\exp 也要求你的常数项收敛，则必须有 [x^0]F=0。

首先对 \exp F 求导，得：

然后，我们计算系数。

再把两边积分回去，得：

注意到这个式子的形式非常猎奇与众不同，其他式子一般都是直接卷积就能得到，但这个经过类似卷积的形式后只能得到一项系数。

暴力显然是 O(n^2) 的，如果为每个系数都暴力 FFT，一方面会常数飞起，一方面复杂度是 O(n^2\log n) 的，更劣了。接下来就是考虑如何让时间复杂度更优。

而这，属于另一项科技。接下来就以多项式 \exp 与它的模板题为例，讲解这项科技。

五、分治 FFT

当我们计算完 [x^n]\exp F 时，注意到会对后面的每一项产生贡献，而产生了贡献，就意味着我们可以在处理前面系数的同时，计算对后面的贡献。

但是逐个向前递推求解显然是不可能的。

注意到 CDQ 分治刚好能够应用于我们现在的贡献处理问题：对于一个区间，我们计算其前半段，然后计算前半段对于后半段的贡献，然后再进行后半段的计算。

什么，你告诉我你不会 CDQ 分治！孩子赶紧重开吧。

CDQ 分治

上面的话已经将 CDQ 分治的思想说的比较明白了，就是在分治求解的同时计算对后面的贡献。接下来以【模板】三维偏序来详解 CDQ 分治。

三个维度显然不好直接做，尝试消去一维。

因为 a_j\le a_i 可以看作第一个要满足的条件（当然你也可以选别的，只要你不怕记错），则可以先将其按照 a 这一维排序，这样保证了 a 这一维单调，同时后面的贡献也变得易于计算。

然后按照一般的、正常的、格式化的处理方法，再用分治来解决 b 这一维。

最后思考如何解决 c 这一维。注意到处理完前两维本质上就是查询前半段有多少个同时满足前两维约束的元素，其 c 也能满足限制。很容易转化为前缀和来解决，由于会动态修改，所以要用树状数组来解决。

简单而言，CDQ 分治解决问题的方法就是：排序解决一维，分治解决许多维，数据结构解决一维（当然你也可以用好多层树套树来解决（此方法没试过，有能力者可尝试））。

然后用主定理分析此题中复杂度：

得到时间复杂度为 O(n\log^2 n) 的。

而对于更广泛的情况，则有：

可得 T(n)=O(n\log^{k-1}n)。

但这个算法是离线算法，而且维度越多，表现越差。在解决 k 维偏序问题时，有这么一句话：

十维以上，不如 n^2。

暴力进队，乱搞 AC。

解决方法详见 FHR 的课件，这里不再过多赘述。

好了，经过了短暂的跑题前置知识讲解，你应该就会这项科技的原理与思想了，然后考虑如何计算对后面产生的贡献。

假设我们已经求出了 \{[x^k]F\mid k\in[l,m)\}，要求解其对 \{[x^k]F\mid k\in[m,r)\} 的贡献。

容易发现，贡献就是 F_{[l,m)}*G_{[0,r-l)}。

把模板题代码扔一下吧。

模板

/*我们所可以自慰的，想来想去，也还是所谓对于将来的希望。

希望是附丽于存在的，有存在，便有希望，有希望，便是光明。*/
#include<bits/stdc++.h>
using namespace std;
constexpr int N=5e5+5,p=998244353;
int rev[N],f[N],g[N],d[N],q[N],n;
inline void init(int len){
    int b=__lg(len);
    for(int i=1;i<len;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline int qpw(int a,int T){
    int res=1;
    for(;T;T>>=1,a=1ll*a*a%p) if(T&1) res=1ll*res*a%p;
    return res;
}
inline void NTT(int* a,int n,int g){
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=2,wn,l,a0,a1;len<=n;len<<=1){
        l=len>>1,wn=qpw(g,(p-1)/len);
        for(int i=0;i<n;i+=len){
            for(int j=0,w=1;j<l;++j,w=1ll*w*wn%p){
                a0=a[i+j],a1=1ll*w*a[i+j+l]%p;
                a[i+j]=(0ll+a0+a1)%p,a[i+j+l]=((a0-a1)%p+p)%p;
            }
        }
    }
    return;
}
inline void INTT(int* a,int n,int l1,int l2){
    NTT(a,n,332748118);
    for(int i=0,invn=qpw(n,p-2);i<=l1+l2;++i) a[i]=1ll*a[i]*invn%p;
    return;
}
inline void CDQ(int l,int r){
    if(l==r) return;
    int mid=(l+r)>>1,len=1;
    CDQ(l,mid);
    for(;len<=r-(l<<1)+mid+1;len<<=1);
    memset(d,0,len<<2),memset(q,0,len<<2);
    memcpy(d,f+l,(mid-l+1)<<2),memcpy(q,g,(r-l+1)<<2);
    init(len);
    NTT(d,len,3),NTT(q,len,3);
    for(int i=0;i<len;++i) d[i]=1ll*d[i]*q[i]%p;
    INTT(d,len,r-l+1,mid-l+1);
    for(int i=mid+1;i<=r;++i) f[i]=(0ll+f[i]+d[i-l])%p;
    return CDQ(mid+1,r);
}
signed main(){
    cin.tie(nullptr)->sync_with_stdio(false);
    cin>>n;f[0]=1;
    for(int i=1;i<n;++i) cin>>g[i];
    CDQ(0,n-1);
    for(int i=0;i<n;++i) cout<<f[i]<<' ';
    return 0;
}

让我们再回到多项式 \exp，发现只需要略微做一些修改，就可以直接套用分治 NTT 的板子了！

模板

/*我们所可以自慰的，想来想去，也还是所谓对于将来的希望。

希望是附丽于存在的，有存在，便有希望，有希望，便是光明。*/
#include<bits/stdc++.h>
using namespace std;
constexpr int N=5e5+5,p=998244353;
int rev[N],f[N],g[N],d[N],q[N],n,inv[N],fac[N];
inline void init(int len){
    int b=__lg(len);
    for(int i=1;i<len;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline int qpw(int a,int T){
    int res=1;
    for(;T;T>>=1,a=1ll*a*a%p) if(T&1) res=1ll*res*a%p;
    return res;
}
inline void NTT(int* a,int n,int g){
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=2,wn,l,a0,a1;len<=n;len<<=1){
        l=len>>1,wn=qpw(g,(p-1)/len);
        for(int i=0;i<n;i+=len){
            for(int j=0,w=1;j<l;++j,w=1ll*w*wn%p){
                a0=a[i+j],a1=1ll*w*a[i+j+l]%p;
                a[i+j]=(0ll+a0+a1)%p,a[i+j+l]=((a0-a1)%p+p)%p;
            }
        }
    }
    return;
}
inline void INTT(int* a,int n,int l1,int l2){
    NTT(a,n,332748118);
    for(int i=0,invn=qpw(n,p-2);i<=l1+l2;++i) a[i]=1ll*a[i]*invn%p;
    return;
}
inline void CDQ(int l,int r){
    if(l==r) return f[l]=1ll*f[l]*inv[l]%p,void();
    int mid=(l+r)>>1,len=1;
    CDQ(l,mid);
    for(;len<=r-(l<<1)+mid+1;len<<=1);
    memset(d,0,len<<2),memset(q,0,len<<2);
    memcpy(d,f+l,(mid-l+1)<<2),memcpy(q,g,(r-l+1)<<2);
    init(len);
    NTT(d,len,3),NTT(q,len,3);
    for(int i=0;i<len;++i) d[i]=1ll*d[i]*q[i]%p;
    INTT(d,len,r-l+1,mid-l+1);
    for(int i=mid+1;i<=r;++i) f[i]=(0ll+d[i-l]+f[i])%p;
    return CDQ(mid+1,r);
}
signed main(){
    cin.tie(nullptr)->sync_with_stdio(false);
    cin>>n;f[0]=fac[0]=inv[0]=1;
    for(int i=0;i<n;++i) cin>>g[i],g[i]=1ll*i*g[i]%p;
    for(int i=1;i<=n;++i) fac[i]=1ll*fac[i-1]*i%p;
    inv[n]=qpw(fac[n],p-2);
    for(int i=n-1;i;--i) inv[i]=1ll*inv[i+1]*(i+1)%p;
    for(int i=1;i<=n;++i) inv[i]=1ll*inv[i]*fac[i-1]%p;
    CDQ(0,n-1);
    for(int i=0;i<n;++i) cout<<f[i]<<' ';
    return 0;
}

六、多项式除法 & 取模

式子长成这个模样：F=G*Q+R。

若 R=0 的话，那直接做一个多项式乘法逆就完成了。但这个情况并不普遍。更多情况下还是会出现余多项式的。

我们构造一个变换 F^R，使得对于多项式 H，若它的阶为 h，则 H^R=x^hH(\frac{1}{x})。

显然，这个变换就是将系数进行翻转。

然后，我们将原来式子中的 x 全部替换为 \frac{1}{x}，然后两边同乘以 x^n，可得：

然后两边都丢到 \pmod{x^{n-m+1}} 意义下，R^R 的影响就被消除了，然后直接多项式求逆就好了。余多项式就是在做一个减法的事。

模板

/*我们所可以自慰的，想来想去，也还是所谓对于将来的希望。

希望是附丽于存在的，有存在，便有希望，有希望，便是光明。*/
#include<bits/stdc++.h>
using namespace std;
constexpr int N=5e5+5,p=998244353;
int rev[N],len=1,m,n,f[N],g[N],ff[N],gg[N],q[N],r[N],tmp[N],inv[N];
inline int qpw(int a,int T){
    int res=1;
    for(;T;T>>=1,a=1ll*a*a%p) if(T&1) res=1ll*res*a%p;
    return res;
}
inline void init(int n){
    int b=__lg(n);
    for(int i=1;i<n;++i) rev[i]=(rev[i>>1]>>1)+((i&1)<<(b-1));
    return;
}
inline void NTT(int*const& a,int n,int type){
    const int g=(~type)? 3:332748118;
    for(int i=0;i<n;++i) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int len=2,wn,l,a0,a1;len<=n;len<<=1){
        l=len>>1,wn=qpw(g,(p-1)/len);
        for(int i=0;i<n;i+=len){
            for(int j=0,w=1;j<l;++j,w=1ll*w*wn%p){
                a0=a[i+j],a1=1ll*w*a[i+j+l]%p;
                a[i+j]=(0ll+a0+a1)%p,a[i+j+l]=((a0-a1)%p+p)%p;
            }
        }
    }
    if(!~type) for(int i=0,invn=qpw(n,p-2);i<n;++i) a[i]=1ll*a[i]*invn%p;
    return;
}
signed main(){
    cin.tie(nullptr)->sync_with_stdio(false);
    cin>>n>>m;
    for(int i=0;i<=n;++i) cin>>f[i],ff[i]=f[i];
    for(int i=0;i<=m;++i) cin>>g[i],gg[i]=g[i];
    reverse(f,f+n+1),reverse(g,g+m+1);
    for(len=1;len<=n-m;len<<=1);
    for(int i=n-m+1;i<=n;++i) f[i]=0;
    for(int i=n-m+1;i<=m;++i) g[i]=0;
    inv[0]=qpw(g[0],p-2);
    for(int t=2,k;t<=len;t<<=1){
        k=t<<1;init(k);
        copy(g,g+t,tmp),fill(tmp+t,tmp+k,0);
        NTT(inv,k,1),NTT(tmp,k,1);
        for(int i=0;i<k;++i) inv[i]=1ll*inv[i]*((2-1ll*inv[i]*tmp[i]%p)%p+p)%p;
        NTT(inv,k,-1);
        fill(inv+t,inv+k,0);
    }
    for(len=1;len<=(n-m+1)<<1;len<<=1);
    init(len);
    NTT(f,len,1),NTT(inv,len,1);
    for(int i=0;i<len;++i) q[i]=1ll*f[i]*inv[i]%p;
    NTT(q,len,-1),reverse(q,q+n-m+1);
    for(int i=0;i<=n-m;++i) cout<<q[i]<<' ';
    cout<<'\n';
    for(len=1;len<=n+1;len<<=1);
    for(int i=n-m+1;i<len;++i) q[i]=0;
    init(len),NTT(gg,len,1),NTT(q,len,1);
    for(int i=0;i<len;++i) gg[i]=1ll*q[i]*gg[i]%p;
    NTT(gg,len,-1);
    for(int i=0;i<m;++i) r[i]=((ff[i]-gg[i])%p+p)%p;
    for(int i=0;i<m;++i) cout<<r[i]<<' ';
    return 0;
}

七、多项式开根 & 幂函数

多项式幂函数

显然，我们有 O(n\log n\log k) 的优秀做法。

然而，这道题可以做到更优。

两边取 \ln 得

再 \exp 回去，可得

但是我们忽略了一件事：[x^0]F 可能为 0。此时我们只要将前面为 0 的项先“删掉”，最后在补回来就好了。

多项式开根

直接跳过二次根，讨论开 k 次根。

可以和上面一样，直接两边求 \ln，然后再 \exp 回去，但是我们没保证 [x^0]F=1。

那还是和上面一样，两边都乘上 [x^0]F 的逆。

最难的就是要求一个 f_0 的 k 次剩余。

八、多项式三角函数

我恨三角函数！我恨三角恒等变换！我恨解三角形！

话虽这么说，但是多项式三角函数还是比较优美的。

首先我们可以将 \sin 与 \cos 的写另一种形式出来。

然后将 F 带入，可得。

而 \tan F，可以由 \dfrac{\sin F}{\cos F} 求出。

九、多项式反三角函数

一开始看，觉得没有什么办法。

此时上万能的微积分就好了。先对其求导，在积分回去。

然后将 F 代进去得。

章四·后话

就先写到这里吧，已经带不动专栏了，而且也累了。

在 NOIP 之后，要是没退役，会更新《浅谈多项式 & 生成函数科技（二）》，会将部分没给出的模板题代码、类包装的全家桶、应用等坑给填完，敬清期待。

参考资料

• 指数型生成函数 & 下降幂 学习笔记——by @C_liar；
• 快速傅里叶变换（FFT）学习笔记——by @KobeBeanBryantCox；
• 多项式初等函数——by OI Wiki；
• 《算法导论·第三版》——by 机械工业出版社；
• 题解 P3803 【【模板】多项式乘法（FFT）】——by @NaCly_Fish
• [可能有用科技]多项式任意次幂&任意次根——by @望月Asta