本文章由 WyOJ Shojo 从洛谷专栏拉取，原发布时间为 2022-06-08 20:10:58
多项式模板大杂烩
本文中 $f(x),F,F(x)$ 均指同一多项式。 $f(x)= \{a_0,a_1,...,a_n\}$ 表示以系数表示法表示的 $n$ 次多项式 $f(x)$，$[x^k]F=a_k$。
$1)$ FFT(NTT)
1. FFT：
FFT依靠分治实现。
具体来说，每次我们将 $f(x)$ 的系数按照奇偶次项分开，假设我们有一个最高次为 $n-1(n=2^k)$ 次的多项式：
$$f(x)=\{a_0,a_1,...,a_{n-1}\}$$
我们的目标是 $\forall i \in [0,n-1]$，求出 $f(\omega_{n}^{i})$
设
$$ \begin{aligned} f_{0}(x)&=\{a_0,a_2,...,a_{n-2}\} \ f_{1}(x)&=\{a_1,a_3,...,a_{n-1}\} \end{aligned} $$
那么我们有
$$f(x)=f_{0}(x^2)+x \cdot f_{1}(x^2)$$
根据单位根的性质有 $k \in [0,\frac{n}{2})$ 时：
$$ \begin{aligned} f(\omega_{n}^{k})&=f_0((\omega_{n}^{k})^2)+\omega_{n}^{k} \cdot f_1((\omega_{n}^{k})^2) \ &= f_0(\omega_{\frac{n}{2}}^{k})+\omega_{n}^{k} \cdot f_1(\omega_{\frac{n}{2}}^{k})\ f(\omega_{n}^{k+\frac{n}{2}})&=f_0((\omega_{n}^{k+\frac{n}{2}})^2)+\omega_{n}^{k+\frac{n}{2}} \cdot f_1((\omega_{n}^{k+\frac{n}{2}})^2) \ &= f_0(\omega_{\frac{n}{2}}^{k})-\omega_{n}^{k} \cdot f_1(\omega_{\frac{n}{2}}^{k}) \end{aligned} $$
特别的，当 $n=1$ 时， $f(\omega_{1}^{0})=\{(\omega_{1}^{0},a_0)\}$ .
也就是说，如果已知 $f_0(x),f_1(x)$ 的点值表示法，那么我们就可以枚举 $k$，在 $O(n)$ 时间内得到 $f(x)$ 的点值表示法，又因为每次 $f_0,f_1$ 长度减半，最多进行 $O(\log n)$ 次，所以总复杂度为 $O(n\log n)$ .
2. 蝴蝶变换：
在FFT的实际实现中，我们不采用递归，而是使用了叫做“蝴蝶变换”的东西。
观察分治过程 $(n=8)$：
第一层： $\{a_{(000)2},a{(001)2},a{(010)2},a{(011)2},a{(100)2},a{(101)2},a{(110)2},a{(111)2}\}$
第二层： $\{a{(000)2},a{(010)2},a{(100)2},a{(110)2}\},\{a{(001)2},a{(011)2},a{(101)2},a{(111)2}\}$
第三层： $\{a{(000)2},a{(100)2}\},\{a{(010)2},a{(110)2}\},\{a{(001)2},a{(101)2}\},\{a{(011)2},a{(111)2}\}$
第四层： $\{a{(000)2}\},\{a{(100)2}\},\{a{(010)2}\},\{a{(110)2}\},\{a{(001)2}\},\{a{(101)2}\},\{a{(011)2}\},\{a{(111)_2}\}$
我们发现个位置叶子层上第 $i$ 个位置对应的 $a$ 的下标恰好是 $i$ 二进制反转后的结果（其实这是必然），所以我们可以直接把 $\{a_{i}\}$ 变换为最后一层，然后一层一层地向上合并。
另外，FFT要求多项式必须是 $2^k-1$ 次（即长度是 $2^k$），在系数表示法后补 $0$ 即可。
3. IFFT：
IFFT即FFT的逆变换，将点值表示法转换回系数表示法。
现在我们已知 $f(x)=\{(\omega_{n}^{0},f(\omega_{n}^{0})),...,(\omega_{n}^{n-1},f(\omega_{n}^{n-1}))\}$，假设它对应的系数表示法是 $\{a_i\}$，那么有：
$$a_k=\frac{1}{n} \sum_{i=0}^{n-1}f(\omega_{n}^i)(\omega_{n}^{-k})^i$$
可以发现其实是构造了一个多项式 $B(x)=\{f(\omega_{n}^{0}),...,f(\omega_{n}^{n-1})\}$，使得 $a_k=\frac{1}{n}B(\omega_{n}^{-k})$，将刚才FFT式子中的 $\omega_{n}^k$ 变为 $\omega_{n}^{-k}$，再做一遍FFT，乘以 $\frac{1}{n}$ 即可。
证一下：
$$ \begin{aligned} \sum_{i=0}^{n-1}f(\omega_{n}^i)(\omega_{n}^{-k})^i &= \sum_{i=0}^{n-1}(\omega_{n}^{-k})^i\sum_{j=0}^{n-1}a_j(\omega_{n}^{i})^j \ &= \sum_{j=0}^{n-1}a_j\sum_{i=0}^{n-1}(\omega_{n}^{j-k})^i \ &= na_k+ \sum_{j=0 & j\not=k}^{n-1}a_j\frac{\omega_{n}^{n(j-k)}-1}{\omega_{n}^{j-k}-1}\ &=na_k \ \ (\omega_{n}^{n(j-k)}-1=\omega_{1}^{j-k}-1=1-1=0)\ \end{aligned} $$
FFT代码就不放了，待会直接看NTT的就行。
4. NTT&INTT：
原根 $(g_n^k)$，就是模 $p$ 意义下的单位根，具有与单位根相似的性质，得以用来计算模 $p$ 意义下的多项式乘法，直接套用FFT的式子即可。同时因为要求 $n=2^k$ 且 $n \mid p-1$，所以 $p$ 都需要可以表示为 $a\cdot2^k+1$ （当然通过中国剩余定理可以做到任意模数多项式乘法(MTT)，这个后面再说）。
最常见的NTT模数是 $998244353$，原根为 $3$，题中见到不常见的模数时先检查是不是多项式模数。
模数的原根（如果存在，其中素数一定有原根）也不会太大，遇到没见过的直接从 $3$ 开始暴力枚举然后检查即可。
代码：
#include <bits\/stdc++.h>
using namespace std;
const int Mod = 998244353;\/\/Mod = 2 ^ 23 * 119 + 1
const int g = 3, ginv = 332748118;\/\/998244353的一个原根及其逆元
const int MAXN = 2.1e6;
int n, m;
int a[MAXN], b[MAXN];
int rev[MAXN], bit, tot;
inline int qpow(int a, int b) {
int base = a, ans = 1;
while (b) {
if (b & 1) ans = 1ll * ans * base % Mod;
base = 1ll * base * base % Mod;
b >>= 1;
}
return ans;
}
void NTT(int a[], bool inv) {
for (int i = 0; i < tot; i++)
if (i < rev[i])
swap(a[i], a[rev[i]]);
for (int len = 1; len < tot; len <<= 1) {
int gn = qpow(inv ? ginv : g, (Mod - 1) \/ (len << 1));
for (int i = 0; i < tot; i += len * 2) {
int g0 = 1;
for (int j = 0; j < len; j++) {
int x = a[i + j], y = 1ll * g0 * a[i + j + len] % Mod;
a[i + j] = (x + y) % Mod;
a[i + j + len] = (x - y + Mod) % Mod;
g0 = 1ll * g0 * gn % Mod;
}
}
}
}
int main() {
scanf("%d%d", &n, &m);
for (int i = 0; i <= n; i++) scanf("%d", &a[i]);
for (int i = 0; i <= m; i++) scanf("%d", &b[i]);
while ((1 << bit) <= n + m) bit++;
tot = 1 << bit;
for (int i = 1; i < tot; i++) rev[i] = rev[i >> 1] >> 1 | (i & 1) << bit - 1;
NTT(a, 0), NTT(b, 0);
for (int i = 0; i < tot; i++) a[i] = 1ll * a[i] * b[i] % Mod;
NTT(a, 1);
int inv = qpow(tot, Mod - 2);
for (int i = 0; i <= n + m; i++) printf("%d ", 1ll * a[i] * inv % Mod);\/\/ \/ tot 变为 * inv
return 0;
}
有了乘法，往后就是多项式的各种操作了，尽量只写式子。
以下运算均在模 $998244353$ 意义下进行，如无特殊说明，都有 $n \le 10^5$ .
$2)$ 多项式求逆
给定 $n-1$ 次多项式 $f(x)$，求多项式 $f^{-1}(x)$，使得 $f(x)*f^{-1}(x) \equiv 1 \pmod {x^n}$ .
显然在 $n=1$ 时，$f^{-1}(x) \equiv a_0^{-1} \pmod x$ .
假设我们现在已经求出了 $G$，满足 $FG \equiv 1 \pmod {x^{\frac{n}{2}}}$，目标是求出 $H$，满足 $FH \equiv 1 \pmod {x^n}$ .
则
$$ \begin{aligned} F(G-H) &\equiv 0 \pmod {x^{\frac{n}{2}}}\ G-H &\equiv 0 \pmod {x^{\frac{n}{2}}}\ (G-H)^2 &\equiv 0 \pmod {x^{n}} \ H^2 &\equiv 2GH - G^2 \pmod {x^{n}}\ H &\equiv 2G-FG^2 \pmod {x^{n}}\ \end{aligned} $$
$n=1$ 时的答案已经知道，从模 $x$ 开始倍增递推即可。
在 $O(n \log n)$ 时间内解决。
inline void calc_inv(int a[], int n) {\/\/a ← a^(-1) (mod x^n) 多项式求逆
int bas = 1, len = 4, bit = 2;
bool opt = 0;
memset(inv[0], 0, sizeof(inv[0]));
inv[0][0] = qpow(a[0], Mod - 2);
while (bas < n) {
calc_rev(bit, len);
opt ^= 1;
memset(inv[opt], 0, sizeof(inv[opt]));
for (int i = 0; i < bas; i++) inv[opt][i] = 2ll * inv[opt ^ 1][i] % Mod;\/\/bas后面都是0
mul(inv[opt ^ 1], inv[opt ^ 1], len);
mul(inv[opt ^ 1], a, len);
bas <<= 1, len <<= 1, bit++;
for (int i = 0; i < bas; i++) inv[opt][i] = (inv[opt][i] - inv[opt ^ 1][i] + Mod) % Mod;
}
for (int i = 0; i < n; i++) a[i] = inv[opt][i];
}
$3)$ 多项式开根
给定 $n-1$ 次多项式 $f(x)$，保证常数项为 $1$，求次数 $< n$ 的多项式 $g(x)$，使得 $g^2(x) \equiv f(x) \pmod {x^n}$，多解情况下取常数项较小的作为答案。
显然 $n=1$ 时，$\sqrt F \equiv 1 \pmod x$ .
假设已经求出 $G$，满足 $G^2 \equiv F \pmod {x^\frac{n}{2}}$，目标求出 $H$，使得 $H^2 \equiv F \pmod {x^n}$ .
则
$$ \begin{aligned} G^2-H^2 &\equiv 0 \pmod {x^\frac{n}{2}}\ (G-H)(G+H) &\equiv 0 \pmod {x^\frac{n}{2}}\ \end{aligned} $$
因为要求常数项最小，即 $1$，所以 $H$ 必然与 $G$ 同余。
即
$$ \begin{aligned} H-G &\equiv 0 \pmod {x^\frac{n}{2}}\ H^2-2GH+G^2 &\equiv 0 \pmod {x^{n}}\ H &\equiv (2G)^{-1}(F+G^2) \pmod {x^{n}}\ H &\equiv 2^{-1}(G^{-1}F+G) \pmod {x^{n}} \end{aligned} $$
结合多项式求逆在 $O(n \log n)$ 时间内解决。
inline void calc_sqrt(int a[], int n) {\/\/a ← sqrt(a) (mod x^n)
int bas = 1, bit = 2, len = 4;
bool opt = 0;
memset(sqr[0], 0, sizeof(sqr[0]));
sqr[0][0] = 1;
while (bas < n) {
opt ^= 1;
memset(sqr[opt], 0, sizeof(sqr[opt]));
for (int i = 0; i < bas; i++) sqr[opt][i] = sqr[opt ^ 1][i];
calc_inv(sqr[opt ^ 1], bas << 1);
calc_rev(bit, len);
mul(sqr[opt ^ 1], a, len);
bas <<= 1, len <<= 1, bit++;
for (int i = 0; i < bas; i++) sqr[opt][i] = 1ll * inv2 * (sqr[opt ^ 1][i] + sqr[opt][i]) % Mod;
}
for (int i = 0; i < n; i++) a[i] = sqr[opt][i];
}
$4)$ 多项式取模
先介绍一个概念，假设我们有多项式 $f(x)=\{a_0,a_1,...,a_{n-1},a_n\}$，那么我们将 $\{a_n,a_{n-1},...,a_1,a_0\}$ 称为 $f(x)$ 的反射多项式，记作 $F^R$。
显然我们有：
$$f(x)=x^nf(\frac{1}{x})$$
进入正题。
给定 $n-1$ 次多项式 $f(x)$，与 $m-1$ 次多项式 $g(x)$，求多项式 $h(x),r(x)$，满足：
- $h(x)$ 次数为 $n-m$，$r(x)$ 次数 $<m-1$ .
- $F=GH+R$ .
不妨设 $R$ 的次数为 $m-2$， 那么有
$$ \begin{aligned} F&=GH+R\ x^{n-1}F&=x^{m-1}G * x^{n-m}H+x^{n-m+1}x^{m-2}R\ F^R&=G^RH^R+x^{n-m+1}R^R\ F^R &\equiv G^RH^R \pmod {x^{n-m+1}}\ H^R &\equiv F^R(G^R)^{-1} \pmod {x^{n-m+1}}\ \end{aligned} $$
结合多项式求逆在 $O(n \log n)$ 时间内解决，最后 $R=F-GH$ 求出 $R$ 即可。
inline void division(int a[], int n, int b[], int m) {
for (int i = 0; i < n; i++) ta[i] = a[n - i - 1];
for (int i = 0; i < m; i++) tb[i] = b[m - i - 1];
calc_inv(tb, n - m + 1);
int len = 1, bit = 0;
while (len < n - m + 1 << 1) len <<= 1, bit++;
calc_rev(bit, len);
mul(ta, tb, len);
reverse(ta, ta + n - m + 1);
for (int i = 0; i < n - m + 1; i++) printf("%d ", ta[i]); puts("");
len = 1, bit = 0;
while (len < m << 1) len <<= 1, bit++;
for (int i = n - m + 1; i < len; i++) ta[i] = 0;
calc_rev(bit, len);
mul(b, ta, len);
for (int i = 0; i < m - 1; i++) printf("%d ", (a[i] - b[i] + Mod) % Mod);
}
$5)$ 多项式 $\ln$
首先我们要知道几个式子：
- 若 $f(x)= \ln x$，则 $f'(x)= \dfrac{1}{x}$ . （自然对数函数求导）
- 若 $f(x)= k \cdot x^a$，则 $f'(x)= ka \cdot x^{a-1}$ . （幂函数求导）
- 若 $f(x)= k \cdot x^a$，则 $\int f(x)dx= \frac{k}{a+1}x^{a+1}+C$ . （幂函数求积）
- 若 $f(x)= g(h(x))$，则 $f'(x)=g'(h(x))h'(x)$ . （链式法则）
另外，对多项式求导（求积）等于对各个单项式求导（求积）后再加起来。
进入正题：
给定 $n-1$ 次多项式 $f(x)$，求多项式 $g(x)$，满足 $g(x) \equiv \ln f(x) \pmod {x^n}$，保证 $[x^0]F=1$ .
$\ln (1-f(x)) = - \sum_{i=0}^{+ \infty} \frac{f^i(x)}{i}$
根据上面的式子，我们有：
$$ \begin{aligned} G &\equiv \ln F \pmod {x^n}\ G' &\equiv F'F^{-1} \pmod {x^n} \ G &\equiv \int F'F^{-1} \pmod {x^n} \end{aligned} $$
结合多项式求逆在 $O(n \log n)$ 时间内解决。
void ln(int a[], int n) {
memset(ta, 0, sizeof(ta));
for (int i = 1; i < n; i++) ta[i - 1] = 1ll * a[i] * i % Mod;
calc_inv(a, n);
int len = 1, bit = 0;
while (len < n << 1) len <<= 1, bit++;
calc_rev(bit, len);
mul(a, ta, len);
for (int i = n - 1; i; i--) a[i] = 1ll * a[i - 1] * qpow(i, Mod - 2) % Mod;
}
$5 \frac{1}{2} )$ $\ \text Newton's \ \text Method$
在 $\exp$ 之前，先介绍一个比较好用的东西——多项式牛顿迭代法。
给定以多项式为参数的函数 $f(x)$，求模 $x^n$ 意义下的多项式 $g(x)$，使得 $f(g(x)) \equiv 0 \pmod {x^n}$ .
依然是倍增的思想，假设我们已经求出来了模 $x^{\frac{n}{2}}$ 意义下的根 $g_0(x)$，那么在模 $x^n$ 意义下的根 $f(x)$ 满足：
$$f(x) \equiv f_0(x)-\frac{g(f_0(x))}{g'(f_0(x))} \pmod {x^n}$$
OI wiki 上给出了证明，需要知道 泰勒展开。
另外，在求导的时候，我们的数域拓展到了多项式，例如当 $f(F)=GF^2+A$ 时，有 $f'(F)=2GF$ .
有了牛顿迭代，上面提到的多项式求逆和多项式开根就可以用此方法求出，举一个开根的例子：
因为有 $g^2(x) \equiv f(x) \pmod {x^n}$，所以我们就是要求出 $h(g(x))=g^2(x)-f(x) \equiv 0 \pmod {x^n}$ 的根，那么有：
$$ \begin{aligned} G &\equiv G_0-\frac{H(G_0)}{H'(G_0)}\ &\equiv G_0-\frac{G_0^2-F}{2G_0}\ &\equiv \frac{G_0^2+F}{2G_0}\ &\equiv 2^{-1}(G_0+FG_0^{-1}) \pmod {x^n} \end{aligned} $$
与另一种方法推出来的式子一致。
$6)$ 多项式 $\exp$
给定 $n-1$ 次多项式 $f(x)$，求多项式 $g(x)$，满足 $g(x) \equiv e^{f(x)}$，保证 $[x^0]F=0$ .
其中 $\exp f(x) = \sum_{i=0}^{+ \infty} \frac{f^i(x)}{i!}$ .
就是求 $h(g(x))= \ln g(x)-f(x) \equiv 0 \pmod {x^n}$ 的根，根据牛顿迭代有：
$$ \begin{aligned} G &\equiv G_0- \frac{ \ln G_0 - F}{ \frac{1}{G_0}}\ &\equiv G_0(1+F- \ln G_0) \pmod {x^n} \end{aligned} $$
结合多项式 $\ln$ 在 $O(n \log n)$ 时间内解决。
void calc_exp(int a[], int n) {
int bas = 1, bit = 2, len = 4;
bool opt = 0;
memset(exp[0], 0, sizeof(exp[0]));
exp[0][0] = 1;
while (bas < n) {
opt ^= 1;
memset(exp[opt], 0, sizeof(exp[opt]));
for (int i = 0; i < bas; i++) exp[opt][i] = exp[opt ^ 1][i];
calc_ln(exp[opt ^ 1], bas << 1);
for (int i = 0; i < bas << 1; i++) exp[opt ^ 1][i] = ((i == 0) + a[i] - exp[opt ^ 1][i] + Mod) % Mod;
calc_rev(bit, len);
mul(exp[opt], exp[opt ^ 1], len);
for (int i = bas << 1; i < len; i++) exp[opt][i] = 0;
bas <<= 1, bit++, len <<= 1;
}
for (int i = 0; i < n; i++) a[i] = exp[opt][i];
}
$7)$ 多项式快速幂
给定 $n-1$ 次多项式 $f(x)$，求多项式 $g(x)$，满足 $g(x) \equiv f^k(x) \pmod {x^n} \ (k \le 10^{10^5})$ 保证 $[x^0]F=1$.
可以看到，$k$ 足足有 $1e5$ 位，$O(n \log n \log k)$ 不用想了。
……但是我们有 $\ln$ 和 $\exp$ 啊，所以有：
$$F^k \equiv e^{k \ln F} \pmod {x^n}$$
所以这玩意当 $i \equiv j \pmod {998244353}$ 时，$F^i \equiv F^j \pmod {x^n}$ .
然后……就没了。
void calc_qpow(int a[], int n, int k) {
calc_ln(a, n);
for (int i = 0; i < n; i++) a[i] = 1ll * a[i] * k % Mod;
calc_exp(a, n);
}
$8)$ 多项式多点求值
给定 $n-1$ 次多项式 $f(x)$ 以及 $m$ 个非负整数 $\{a_i\}$，$\forall i \in [1,m]$，求出 $f(a_i)$ .
默认小写字母表示列向量，大写字母表示矩阵。
矩阵转置：对于矩阵 $A$ 中任意一个元素 $a_{i,j}$，将 $i,j$ 交换，得到的矩阵称为 $A$ 的转置矩阵，记作 $A^T$，其满足以下性质：
- $(A^T)^T=A$
- $(A+B)^T=A^T+B^T$
- $(AB)^T=B^TA^T$
转置原理：又称特勒根原理，其指出当我们要求解 $b=Aa$ 时，可以找到 $b'=A^Ta$ 的一种高效求解方法，将方法中每一步对 $a$ 的操作拆解出来，即把 $A^T$ 分解为 $E_1E_2...E_n$，然后倒序执行每一个操作的转置，即将 $E_n^TE_{n-1}^T...E_1^T$ 施加于 $a$，这样我们就得到了 $b$ .
关于操作的转置，比如计算卷积 $\{c_0,c_1,...,c_{n+m}\}=\{b_0,b_1,...,b_m\}*\{a_0,a_1,...a_n\}$ 时，可以看作：
$$ \begin{bmatrix} c_0\ c_1\ c_2\ \ dots\ c_{n+m} \end{bmatrix}
\begin{bmatrix} b_0 &0 &0 &\cdots &0\ b_1 &b_0 &0 &\cdots &0\ b_2 &b_1 &b_0 &\cdots &0\ \ dots &\ dots &\ dots &\ddots &\ dots\ b_{n+m} &b_{n+m-1} &b_{n+m-2} &\cdots &b_0 \end{bmatrix} * \begin{bmatrix} a_0\ a_1\ a_2\ \ dots\ a_{n+m} \end{bmatrix} $$
将方阵转置，我们得到：
$$ \begin{bmatrix} c'0\ c'_1\ c'_2\ \ dots\ c'{n+m} \end{bmatrix}
\begin{bmatrix} b_0 &b_1 &b_2 &\cdots &b_{n+m}\ 0 &b_0 &b_1 &\cdots &b_{n+m-1}\ 0 &0 &b_0 &\cdots &b_{n+m-2}\ \ dots &\ dots &\ dots &\ddots &\ dots\ 0 &0 &0 &\cdots &b_0 \end{bmatrix} * \begin{bmatrix} a_0\ a_1\ a_2\ \ dots\ a_{n+m} \end{bmatrix} $$
此时 $c'i=\sum{j=i}^{n}a_{j}b_{j-i}$，即卷积的转置操作，记作 $\{b_i\}*^T\{a_i\}$。
回到正题，假设答案序列为 $\{b_1,b_2,...,b_n \}$，我们有 $(n=m)$：
$$ \begin{bmatrix} b_1\ b_2\ \ dots\ b_n \end{bmatrix}
\begin{bmatrix} 1 &a_1 &a_1^2 &\cdots &a_1^{n-1}\ 1 &a_2 &a_2^2 &\cdots &a_2^{n-1}\ \ dots &\ dots &\ dots &\ddots &\ dots\ 1 &a_{n} &a_{n}^2 &\cdots &a_{n}^{n-1} \end{bmatrix} * \begin{bmatrix} f_0\ f_1\ \ dots\ f_{n-1} \end{bmatrix} $$
转置之后得到：
$$ \begin{bmatrix} b'_1\ b'_2\ \ dots\ b'_n \end{bmatrix}
\begin{bmatrix} 1 &1 &1 &\cdots &1\ a_1 &a_2 &a_3 &\cdots &a_n\ \ dots &\ dots &\ dots &\ddots &\ dots\ a_1^{n-1} &a_2^{n-1} &a_3^{n-1} &\cdots &a_n^{n-1} \end{bmatrix} * \begin{bmatrix} f_0\ f_1\ \ dots\ f_{n-1} \end{bmatrix} $$
设 $\{b'_i\}$ 的 $\operatorname{OGF}$ 为 $B$，有：
$$ \begin{aligned} B &= \sum_{i=1}^{+\infty} b_ix^i\ &= \sum_{i=1}^{+\infty} \left(\sum_{j=1}^{n} a_j^{i-1}f_{j-1} \right)x^i\ &= \sum_{j=1}^{n} f_{j-1} \left(x\sum_{i=0}^{+\infty}a_j^ix^i \right)\ &= x\sum_{i=1}^{n}\frac{f_{i-1}}{1-a_ix} \end{aligned} $$
这个式子可以用分治 $\operatorname{FFT}$ 解决，维护 $M_{l,r}, P_{l,r}$ 分别表示区间 $[l,r]$ 的分子和分母，有：
$$ \begin{aligned} M_{l,l}&=1-a_lx\ P_{l,l}&=f_{l-1}\ M_{l,r}&=M_{l,mid}M_{mkd+1,r}\ P_{l,r}&=P_{l,mid}M_{mid+1,r}+P_{mid+1,r}*M_{l,mid} \end{aligned} $$
最后答案是 $M_{1,m}^{-1}P_{1,m}$ .
整个分治过程都可以看作是对 $\{f_i\}$ 进行的一系列初等变换，因此直接把所有操作倒序并执行其转置操作。
具体来讲，先计算出常数 $M_{l,r}$ ，然后令 $P_{1,m} = \{f_i\}^TM_{1,m}^{-1}$，接着从根向叶节点分治，将上述分治过程中的转移全部操作全部换为其转置操作，即：
$$ \begin{aligned} P_{l,mid}&=P_{l,r}^TM_{mid+1,r}\ P_{mid+1,r}&=P_{l,r}^TM_{l,mid} \end{aligned} $$
按顺序输出 $P_{l,l}$ 的常数项即可。
\/\/这里只放了第二次分治
void solve2(int p, int l, int r) {
if (l == r) return printf("%d\n", F[p][0]), void();
int mid = (l + r) >> 1;
F[p << 1] = mulT(F[p], M[p << 1 | 1]); F[p << 1].resize(mid - l + 1);
F[p << 1 | 1] = mulT(F[p], M[p << 1]); F[p << 1 | 1].resize(r - mid);
solve2(p << 1, l, mid);
solve2(p << 1 | 1, mid + 1, r);
}
$9)$ 线性递推
已知数列 $\{f_i\}$ 满足 $k$ 阶线性递推关系：
$$f_i=\sum_{j=1}^kg_jf_{i-j}$$
给定 $g_{1...k},f_{0...k-1}$，求 $f_n$ $(n \le 10^9,k \le 32000)$ .
（当 $k$ 较小时，一般直接使用矩阵快速幂解决）
考虑 $\operatorname{OGF}$，我们有：
$$ \begin{aligned} F(x)&=G(x)F(x)+R(x)\\ F(x)&=\frac{R(x)}{1-G(x)} \end{aligned} $$
其中 $\operatorname{deg}R < k$，用来调整 $0$ 到 $k-1$ 位。
那么只要题目满足 $F(x)=\dfrac{R(x)}{Q(x)}$ 我们就可以称其具有线性递推关系，下面介绍一下这种关系的机械处理方法。
$$ \begin{aligned} [x^n]F(x)&=[x^n]\frac{R(x)}{Q(x)}\ &= [x^n]\frac{R(x)Q(-x)}{Q(x)Q(-x)}\ \end{aligned} $$ 令 $U(x^2)=Q(x)Q(-x)$，$A(x^2)+xB(x^2)=R(x)Q(-x)$，有： $$ \begin{aligned} [x^n]F(x) &= [x^n]\frac{A(x^2)}{U(x^2)}+[x^n]x\frac{B(x^2)}{U(x^2)}\ &= \begin{cases} [x^{\frac{n}{2}}]\dfrac{A(x)}{U(x)}& n \bmod 2=0\ [x^{\frac{(n-1)}{2}}] \dfrac{B(x)}{U(x)}& n\bmod 2=1 \end{cases} \end{aligned} $$
$n$ 每次会减小一半，时间复杂度 $O(k\log k\log n)$，还可以进一步将常数 $\times \dfrac{1}{3}$，具体参考杜爷的博客 .
当然对于分母的处理方式不一定是乘上 $Q(-x)$，只要能将分母化为 $U(x^b)$ 的形式，就能从 $O(n)$ 降为 $O(\log_{b}n)$ .
模板还有很多，有时间再补（
鲁ICP备2025150228号