多项式学习笔记

1. 前置知识

1.1 基础

\(f(x) = \sum_{i=0}^na_ix^i\) 被称为一个 \(n\) 次多项式。

\(\deg f(x)\) 表示多项式的次数。

\(f(x)g(x) = h(x)\) 称为多项式乘法，也叫多项式卷积，满足 \(h_n = \sum_{i + j = n}f_ig_j\)。

1.2 点值

给定一个多项式 \(f(x)\)，再给 \(m\) 个点 \(x_1, \dots, x_m\)，求 \(f(x_1), \dots, f(x_m)\)。

1.3 插值

1.3.1 定义

对于一个多项式 \(f(x)\)，给定 \(n\) 横坐标不相同的个点 \((x_0, y_0), \dots, (x_{n - 1}, y_{n - 1})\)， 求一个多项式 \(f(x)\) 使得对于所有 \(0 \le i < n\)，\(f(x_i) = y_i\)。

不难证明如果 \(\deg f(x) < n\) 则 \(f(x)\) 是唯一的。

考虑反证法，存在 \(f(x)\) 和 \(g(x)\) 满足要求。

则设 \(r(x) = f(x) - g(x)\)，显然 \(r(x)\) 有 \(n\) 个根。但是 \(\deg r(x) < n\)，所以 \(r(x) = 0\)，矛盾！

所以唯一性得证。

1.3.2 拉格朗日插值

Lagrange 插值法，是一种可以计算 \(f(x)\) 的各项系数的方法。

我们考虑构造 \(n\) 个多项式满足 \(P_i(x_i) = y_i, P_i(x_j) = 0\)，显然答案就是 \(\sum P_i(x)\)。

不难得出：

\[P_i(x) = y_i\frac{\prod_{i \not = j}(x - x_j)}{\prod_{i \not =j}(x_i - x_j)} \]

则最终的答案：

\[f(x) = \sum y_i \prod_{i \not = j}\frac{x - x_j}{x_i - x_j} \]

这个用代码实现复杂度是 \(O(n^2)\) 的。

1.4 复数

1.4.1 定义

定义 \(i^2 = -1\)，称为虚数。

对于所有 \(a, b \in \mathbb{R}\) 的 \(z = a + bi\) 称为复数。

1.4.2 基础知识

复数可以用笛卡尔坐标系来表示，横坐标表示实数，纵坐标表示虚树，\(a + bi\) 对应 \((a, b)\)。

一个复数可以对应到一个向量上。

不妨设 \(r = \sqrt{a^2 + b^2}\)，\(\theta\) 为 \(\overrightarrow{z} = a + bi\) 的夹角，则我们也可以表示为：\(\overrightarrow{z} = r(\cos \theta + i \sin \theta)\)。

1.4.3 欧拉公式

\[e^{ix} = \cos x + i \sin x \]

证明可以用泰勒展开来证明。

1.4.4 单位根

在复数域上，方程 \(x^n = 1\) 会有 \(n\) 个不同的解，这 \(n\) 个不同的解记作 \(\varepsilon_n^0, \dots, \varepsilon_n^{n-1}\)。

根据欧拉公式，我们有：

\[\varepsilon_n^k = e^{i \frac{2\pi k}{n}} = \cos \frac{2\pi k}{n} + i \sin \frac{2\pi k}{n} \]

这样可以快速计算单位根。

画在单位圆上，单位根会将整个单位圆 \(n\) 等分。

单位根有一些性质和定理：

折半定理： \(\varepsilon^2_{2n} = \varepsilon^1_{n}\)

2. 多项式运算基础

2.1 快速傅里叶变换 （FFT）

2.1.1 大致思路

考虑到多项式卷积的朴素算法是 \(O(n^2)\) 的。我们需要优化到 \(O(n \log n)\)。

根本原理在于只需知道 \(n\) 个点就能确定 \(h(x)\)，我们通过选取特殊的点来优化计算。

整个过程分为三步：

1. 选取一些点 \(x_0, \dots, x_{n -1}\)， 计算 \(f(x)\) 和 \(g(x)\) 的点值。

2. 根据 \(h(x) = f(x)g(x)\) 计算出 \(h(x)\) 的点值。

3. 插值得到 \(h(x)\)。

2.1.2 DFT 点值

我们考虑选取 \(n\) 次单位根 \(\varepsilon_n^0, \varepsilon_n^1, \dots, \varepsilon_n^{n-1}\) 作为点值。

考虑利用分治来求解。我们要求 \(n = 2^k\)，如果不是就补齐即可。

不妨设 \(f(x) = \sum_{i = 0}^{n - 1}a_ix^i\)，我们根据系数下标的奇偶性分开成两个多项式：

\[f^{(0)} (x) = a_0x^0 + a_2x^2 + \dots + a^{n-2}x^{\frac{n - 2}{2}} \]

\[f^{(1)} (x) = a_1x^0 + a_3x^2 + \dots + a^{n-1}x^{\frac{n - 2}{2}} \]

不难发现 \(f(x) = f^{(0)}(x^2) + xf^{(1)}(x^2)\)，我们将问题转化成计算这两个多项式在 \(\varepsilon_{\frac{n}{2}}^0, \dots, \varepsilon_{\frac{n}{2}}^{\frac{n}{2}-1}\) 的点值。

根据折半定理，\(\varepsilon_n^k = \varepsilon_\frac{n}{2}^{k \bmod \frac{n}{2}}\)，所以我们可以通过这两个多项式推出原来的多项式。

2.1.3 IDFT 插值

不妨考虑写成矩阵的形式：

\[\begin{bmatrix} y_0\\ y_1\\ \vdots\\ y_{n-1} \end{bmatrix} V = \begin{bmatrix} a_0\\ a_1\\ \vdots\\ a_{n-1} \end{bmatrix} \]

其中 \(V(i,j) = \varepsilon_n^{ij \bmod n}\)。

则我们现在知道 \(V, Y\)，需要求出 \(A\)。

我们只需要求出 \(V^{-1} . Y\) 即可。

根据观察可以发现 \(V^{-1} = \frac{1}{n}\varepsilon_n^{-ij}\)。

现在相当于我们对于函数 \(F(x) = \sum y_ix^i\) 求得其在 \(\frac{1}{n}\varepsilon_{n}^{-i}\) 这 \(n\) 个点的系数了，做法和之前的 DFT 类似。

所以我们就可以计算出系数了。

2.1.4 代码实现

考虑到 DFT 和 IDFT 除了再代入树枝上有细小区别其他都一样，所以我们用一个函数来实现，传递参数来表示是 DFT 还是 IDFT。

递归一般会比较慢，我们发现最后一个数 \(x\) 的底层位置刚好为其二进制位反转后的位置，我们通过位运算计算出这个位置，直接自下往上合并即可。

这里给出 P3803 【模板】多项式乘法（FFT） 模板：

点击查看代码

#include <iostream>
#include <cstdio>
#include <cmath>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;
const int N = 3e6 + 5;
const double PI = acos(-1);

struct cnum {//复数 z = a + bi 
	double a, b;
	cnum (double _a = 0, double _b = 0) :
		a(_a), b(_b) {} 
};
cnum operator+(cnum x, cnum y) {
	return cnum(x.a + y.a, x.b + y.b);
}
void operator+=(cnum &x, cnum y) {
	x = x + y;
}
cnum operator-(cnum x, cnum y) {
	return cnum(x.a - y.a, x.b - y.b);
}
void operator-=(cnum &x, cnum y) {
	x = x - y;
}
cnum operator*(cnum x, cnum y) {
	return cnum(x.a * y.a - x.b * y.b, x.a * y.b + x.b * y.a);
}
void operator*=(cnum &x, cnum y) {
	x = x * y;
} 

struct poly {//多项式 
	int n;//次数 
	vector<int> a;//系数
	poly (int _n = 0) {
		n = _n;
		a = vector<int>(n, 0);
	} 
};
poly operator+(poly x, poly y) {
	poly ans = poly(max(x.n, y.n));
	for (int i = 0; i < x.n; i++)
		ans.a[i] += x.a[i];
	for (int j = 0; j < y.n; j++)
		ans.a[j] += y.a[j];
	return ans;
} 
poly operator-(poly x, poly y) {
	poly ans = poly(max(x.n, y.n));
	for (int i = 0; i < x.n; i++)
		ans.a[i] += x.a[i];
	for (int j = 0; j < y.n; j++)
		ans.a[j] -= y.a[j];
	return ans;
}
int rev[N] = {0};

void FFT(cnum *a, int lim, int op) {//op 表示是 DFT 还是 IDFT 
	for (int i = 0; i < lim; i++)
		if (i < rev[i])
			swap(a[i], a[rev[i]]);
	for (int mid = 1; mid < lim; mid <<= 1) {
		cnum wn = cnum(cos(PI / mid), op * sin(PI / mid));
		for (int j = 0; j < lim; j += (mid << 1)) {
			cnum w = cnum(1, 0);
			for (int k = 0; k < mid; k++, w = w * wn) {
				cnum x = a[j + k], y = w * a[j + k + mid];
				a[j + k] = x + y;
				a[j + k + mid] = x - y;
			}
		}
	}
}
cnum a[N], b[N];
poly operator*(poly f, poly g) {
	int lim = 1, pw = 0, len = f.n + g.n - 1;
	while (lim < len)
		lim <<= 1, pw++;
	for (int i = 0; i < lim; i++) {
		a[i] = cnum((i < f.n ? f.a[i] : 0), 0);
		b[i] = cnum((i < g.n ? g.a[i] : 0), 0);
		rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (pw - 1)));
	}
	FFT(a, lim, 1), FFT(b, lim, 1);
	for (int i = 0; i < lim; i++)
		a[i] *= b[i];
	FFT(a, lim, -1);
	poly ans = poly(lim);
	for (int i = lim - 1; i >= 0; i--)
		ans.a[i] = (int)(a[i].a / lim + 0.5);
	while (ans.n > 0 && ans.a.back() == 0)
		ans.n--, ans.a.pop_back();
	return ans;
}

int main() {
	int n, m;
	cin >> n >> m;
	++n, ++m;
	poly f = poly(n), g = poly(m);
	for (int i = 0; i < n; i++)
		cin >> f.a[i];
	for (int i = 0; i < m; i++)
		cin >> g.a[i];
	poly h = f * g;
	for (int i = 0; i < h.n; i++)
		cout << h.a[i] << " ";
	for (int i = 1; i <= n + m - 1 - h.n; i++)
	    cout << 0 << " ";
	cout << endl; 
	return 0;
}

2.2. 快速数论变换 （NTT）

考虑到 FFT 需要用到复数和三角函数，精度很低，而且不持支取模。

回顾整个过程，关键在于需找合适的数字进行 DFT，要满足其若干次放构成一个完整的系，我们不难想到原根。

\(g\) 是模 \(m\) 意义下的原根当且仅当 \(1 \le i \le \varphi(m), g^{i} \bmod m\) 各不相同。

不难发现如果有了这个性质，我们就可以进行快速的点值和插值。

然而还有问题，对于质数 \(p\) 和原根 \(g\)，模 \(p\) 意义下的 \(n\) 次单位根是 \(g^{\frac{p - 1}{n}}\)，这就要求 \(n | p - 1\)，而分治过程中所有的 \(n\) 都是 \(2\) 的幂，所以 \(p-1\) 要含有尽可能多的 \(2\)。

而 \(998244353 = 2^{23} \times 7 \times 17 + 1\) 刚好满足这个性质，其最小原根是 \(3\)，所以我们就在这个模数的意义下进行 \(NTT\)。

这里给出代码：

点击查看代码

#include <iostream>
#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;
const int mod = 998244353;

int fpow(int a, int b, int p) {//快速幂取模 
	if (b == 0)
		return 1;
	int ans = fpow(a, b / 2, p);
	ans = 1ll * ans * ans % p;
	if (b % 2 == 1)
		ans = 1ll * a * ans % p;
	return ans;
}

int mmi(int a, int p) {//求逆元 
	return fpow(a, p - 2, p);
}

struct poly {//多项式类 
	int n;//表示项数，次数为 n - 1 
	vector<int> a;//系数 
	poly (int _n = 0) {//初始化一个 _n 项多项式 
		n = _n;
		a = vector<int>(n, 0);
	}
};
poly operator+(poly x, poly y) {//重载多项式加法 
	poly ans = poly(max(x.n, y.n));
	for (int i = 0; i < x.n; i++)
		ans.a[i] = (ans.a[i] + x.a[i]) % mod;
	for (int i = 0; i < y.n; i++)
		ans.a[i] = (ans.a[i] + y.a[i]) % mod;
	return ans;
}

poly operator-(poly x, poly y) {//重载多项式减法
 	poly ans = poly(max(x.n, y.n));
	for (int i = 0; i < x.n; i++)
		ans.a[i] = (ans.a[i] + x.a[i]) % mod;
	for (int i = 0; i < y.n; i++)
		ans.a[i] = (ans.a[i] - y.a[i] + mod) % mod;
	return ans;
}

void DFT(vector<int> &a, vector<int> &rev, int lim, int op) {//op = 1 DFT 否则是 IDFT 
	for (int i = 0; i < lim; i++)
		if (i < rev[i])
			swap(a[i], a[rev[i]]);
	for (int i = 1; i < lim; i <<= 1) {
		int wn = fpow((op == 1 ? 3 : (mod + 1) / 3), (mod - 1) / (i << 1), mod);
		for (int j = 0; j < lim; j += (i << 1)) {
			int w = 1;
			for (int k = 0; k < i; k++, w = 1ll * w * wn % mod) {
				int x = a[j + k], y = 1ll * w * a[i + j + k] % mod;
				a[j + k] = (x + y) % mod, a[i + j + k] = (x - y + mod) % mod;
			}
		}
	}
	if (op == -1) {
		int inv = mmi(lim, mod);
		for (int i = 0; i < lim; i++)
			a[i] = 1ll * a[i] * inv % mod;
	}
}

poly NTT(poly f, poly g) {//多项式NTT，返回两个多项式模意义下的卷积 
	int len = f.n + g.n - 1;
	int lim = 1, pw = 0;
	while (lim < len)
		lim <<= 1, pw++;
	f.n = g.n = lim;
	f.a.resize(lim);
	g.a.resize(lim);
	vector<int> rev = vector<int>(lim, 0);
	for (int i = 0; i < lim; i++)
		rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (pw - 1)));
	DFT(f.a, rev, lim, 1), DFT(g.a, rev, lim, 1);
	for (int i = 0; i < lim; i++)
		f.a[i] = 1ll * f.a[i] * g.a[i] % mod;
	DFT(f.a, rev, lim, -1);
	while (f.n > 0 && f.a.back() == 0)
		f.n--, f.a.pop_back();
	return f;
}

int main() {
	int n, m;
	cin >> n >> m;
	++n, ++m;
	poly f = poly(n), g = poly(m);
	for (int i = 0; i < n; i++)
		cin >> f.a[i];
	for (int i = 0; i < m; i++)
		cin >> g.a[i];
	poly h = NTT(f, g);
	for (int i = 0; i < h.n; i++)
		cout << h.a[i] << " ";
	for (int i = 1; i <= n + m - 1 - h.n; i++)
	    cout << 0 << " ";
	cout << endl; 
	return 0;
} 

2.3 MTT

MTT 其实很简单，解决任意模数 NTT 问题，我们只需要选取三个满足 NTT 的模数，然后用 CRT 合并得出解即可。考虑到值域不大，我们可以找到唯一的解。

2.4 bluestein (Chirp Z 变换)

2.4.1 原理

不妨考虑任意模长 DFT：求 \(A(x) = \sum_{i = 0}^{n-1}a_ix^i\) 在 \(\omega_{n}^0 \sim \omega_{n}^{n-1}\) 的取值。不保证 \(n\) 是 \(2\) 的幂。

我们要求的就是 \(b_m = \sum_{k = 0}^{n-1}a_k\omega_{n}^{mk}\)，我们需要进行一些变化。

技巧： 将 \(mk\) 转化成下阶乘幂，有：

\[mk = \frac{1}{2}[m^{\underline{2}} + (k+1)^{\underline{2}} - (m - k)^{\underline{2}}] \]

注意我们最好用二次下阶乘幂而不是平方，这样就不会设计开方等比较复杂且不一定存在的操作。

我们代入得到：

\[b_m = \omega_n^{\frac{1}{2}m^{\underline{2}}}\sum_{k=0}^{n-1}a_k\omega_n^{\frac{1}{2}(k+1)^{\underline{2}} - \frac{1}{2}(m - k)^{\underline{2}}} \]

我们考虑平移从而变成卷积的形式：令 \(B_{m+n} = \omega_n^{-\frac{1}{2}m^{\underline{2}}}b_m\)。

我们尝试构造 \(B = A * C\) 从而求出 \(b_m\)。

不妨设 \(A_k = a_k\omega_n^{\frac{1}{2}{(k+1)}^{\underline{2}}}\)，则：

\[B_{m + n} = \sum_{k=0}^{n-1}A_k \omega_n^{-\frac{1}{2}(m-k)^{\underline{2}}} \]

为了构造卷积，我们需要 \(C_{m + n - k} = \omega_n^{-\frac{1}{2}(m-k)^{\underline{2}}}\)， \(C_j = \omega_{n}^{-\frac{1}{2}(j-n)^{\underline{2}}}\)。

所以我们得到：

\[B_{m + n} = \sum_{k = 0}^{n - 1}A_kC_{m + n - k} \]

但是这还不是标准的卷积，为了达到效果，我们令 \(A_{n} = A_{n + 1} = \dots = A_{n + m} = 0\)，这样就可以写成：

\[B_{m + n} = \sum_{k=0}^{m + n}A_{k}C_{m + n - k} \]

我们运用一次 FFT 就可以求得 \(B\)，从而推出 \(b\)。

回顾整个过程不难发现我们不一定需要单位根，对于任意得等比数列我们都可以在 \(O(n \log n)\) 的时间内算出其 DFT 的结果。

2.4.2 变种

不妨考虑 DFT 中我们不求 \(n\)，而是求某个小于 \(n\) 的 \(k\)，我们求 \(\omega_{k}^0 \sim \omega_{k}^{k-1}\) 的所有取值需要怎么办。

我们借助降幂多项式，具体的，我们设：

\[B(x) = \sum_{i = 0}^{n-1}a_ix^{i \bmod k} \]

则 \(B\) 是一个 \(k\) 项多项式，我们有结论：对于 \(\omega_{k}^0 \sim \omega_{k}^{k-1}\)，\(B(x)\) 在这些点的取值与 \(A(x)\) 相等。

这个很好证明，由于单位根的周期性，我们会发现指数对 \(k\) 取模不受影响，所以得出这个结论。

所以我们借助这个 \(B(x)\) 对其进行 bluestein 即可。

2.4.3 应用

HDU4656 Evaluation

由于模数很小，我们转化为求出 \(F(0) \sim F(p-1)\) 即可。

可以直接用多项式快速点值，但是太慢了，我们尝试用 bluestein 来做。

考虑到模意义下原根的良好性质，及 \(g^{0} \sim g^{p-1}\) 构成模 \(p\) 意义下的简化剩余系，则我们转化为求一个等比数列的值。

这样就可以用 bluestein 来解决了。

如果用 FFT 精度可能不够，而如果用 \(NTT\)，考虑到 \(p - 1 = 2 \times 2 \times 16667\)，没办法分治，所以这道题只能用 MTT 来做。

2.5 循环卷积

循环卷积的定义：

\[c_r = \sum_{p + q \equiv r \pmod n}a_p \times b_q \]

解决这个问题，不难想到构造 \(2n\) 长度的 \(d = a \times b\)，然后 \(c_r = d_r + d_{r + n}\)，但是这样太笨了。

不难发现，\(D(x)\) 对 \(n\) 降幂后就是 \(C(x)\)，我们考虑将其联系在一起。

我们希望求出 \(C(x)\) 的系数，这就相当于计算 \(D(x)\) 在 \(\omega_n^0 \sim \omega_n^{n-1}\) 的取值，考虑到 \(D = A \times B\)，于是我们就有 \(C = IDFT(DFT(A) \cdot DFT(B))\)，这里是指点乘。

于是我们就可以转化成计算 IDFT 与 DFT 了。

这个可以推广到 \(k\) 次循环卷能得到：\(C = IDFT(DFT^k(A))\)，也就是一次 DFT 后对每个数取 \(k\) 次方。

3.多项式运算综合

注意一些代码实现的时候关于对 \(x^m\) 取模需要考虑，这个很关键。

3.1 多项式乘法

问题

给定多项式 \(A(x), B(x)\)，求 \(C(x) = A(x)B(x)\)。

方法

直接利用 FTT 或 NTT 即可，时间复杂度 \(O(n \log n)\)。

3.2 多项式取模

问题

给定多项式 \(A(x)\) 和整数 \(n\)，求 \(B(x) \bmod x^n\)。

方法

首先我们需要定义单项式的取模：

\[cx^k \bmod x^n = \begin{cases} 0& k \ge n\\ cx^k & k < n \end{cases} \]

而多项式取模就是对于每一项分别取模后相加。

其实就是保留所有次数小于 \(n\) 的项。可以看成去掉高阶无穷小量。

3.3 多项式求逆

问题

不妨设 \(n = \deg A(x) + 1\)。

一般形式： 给定多项式 \(A(x)\) 和 \(C(x)\)，求 \(B(x)\) 使得 \(A(x)B(x) \equiv C(x) \pmod {x^n}\)。

特殊形式： 当 \(C(x) = 1\) 时，求 \(B(x)\) 使得 \(A(x)B(x) \equiv 1 \pmod {x^n}\)。

方法

首先，如果推出了特殊形式，我们可以通过两边同时乘上 \(C(x)\) 得到一般形式的解，直接使用多项式乘法时间复杂度 \(O(n \log n)\)，所以我们考虑如何求解特殊形式。

考虑常数项，如果 \(A_0 = 0\)，显然无解，因为 \(0\) 不可逆。

否则，我们根据乘法的定义不难推出 \(B_0 = 1\)，然后根据 \(A_0B_1 + A_1B_0 = 0\) 推出 \(B_1\)，以此类推。

但是这样时间复杂度是 \(O(n^2)\) 考虑更快。

考虑倍增，假设我们知道了 \(B_1(x) = B(x) \bmod x^m\)，我们希望通过一些变化推出 \(B_2(x) = B(x) \bmod x^{2m}\) 从而最终推出 \(B(x)\)。

然后就是推式子了：

\[\begin{aligned} A(x)B_1(x) \equiv 1 \pmod {x^m}\\ A(x)B_1(x) - 1\equiv 0 \pmod {x^m}\\ \end{aligned} \]

考虑平方后会导致 \(\bmod x^{2m}\) 依然是 \(0\)，所以：

\[(A(x)B_1(x) - 1)^2\equiv 0 \pmod {x^{2m}}\\ \]

为了方便，去掉 \(x\)：

\[\begin{aligned} A^2B_1^2 - 2AB_1 + 1 \equiv 0 \pmod{x^{2m}}\\ 2AB_1 - A^2B_1^2 \equiv 1 \pmod{x^{2m}}\\ A(2B_1 - AB_1^2)\equiv 1 \pmod{x^{2m}}\\ \end{aligned} \]

所以 \(B_{2} = B_1(2 - AB_1)\)，可以通过多项式乘法在 \(O(n \log n)\) 内计算，则总时间复杂度为：

\[\sum_{k}\frac{n}{2^k} \log n = 2 n \log n = O(n \log n) \]

所以我们就解决了多项式求逆，注意要判断常数项是否是 \(0\)。

3.4 多项式除法

问题

给定 \(A(x), B(x)\)，求 \(C(x), R(x)\) 满足 \(A(x) = B(x)C(x) + R(x)\) 且 \(\deg R(x) < \deg B(x)\)。

方法

如果 \(\deg A(x) < \deg B(x)\)，则 \(R(x) = A(x), C(x) = 0\)。

否则，我们定义运算 \(A^R(x)\) 表示：

\[A^R(x) = x^nA(\frac{1}{x}) \]

其中 \(n = \deg A(x)\)，其实本质就是反转系数。

不妨设 \(n = \deg A(x), m = \deg B(x)\)， 显然有 \(\deg C(x) = n - m\)，则：

\[\begin{aligned} A(x) &= B(x)C(x) + R(x)\\ A(\frac{1}{x}) &= B(\frac{1}{x})C(\frac{1}{x}) + R(\frac{1}{x})\\ x^nA(\frac{1}{x}) &= x^mB(\frac{1}{x})x^{n-m}C(\frac{1}{x}) + x^nR(\frac{1}{x})\\ A^R(x) &= B^R(x)C^R(x) + R(\frac{1}{x})x^{n}\\ \end{aligned} \]

显然 \(\deg R(x) \le n - 1\)，我们把它视作 \(m-1\) 次多项式，然后进行反转得到：

\[A^R(x) = B^R(x)C^R(x) + R^R(x)x^{n - m + 1} \]

由于 \(\deg C(x) = n - m\)，所以我们只要求出 \(C(x) \bmod x^{n - m + 1}\) 即可。所以两边同时取模消去 \(R\) 得到：

\[A^R(x) \equiv B^R(x)C^R(x) \pmod {x^{n - m + 1}} \]

从而：

\[C^R(x) \equiv A^R(x){B^R}(x)^{-1} \pmod {x^{n - m + 1}} \]

利用多项式求逆在 \(O(n \log n)\) 解决，然后再反转系数即可。

\(R(x)\) 可以直接通过 \(A(x) - B(x)C(x)\) 配合多项式乘法计算，也是 \(O(n \log n)\) 的。

3.5 多项式对数函数

问题

设 \(n = \deg A(x) + 1\)，求 \(B(x)\) 使得 \(B(x) \equiv \ln(A(x)) \pmod {x^n}\)。

方法

首先我们需要利用一个 \(\ln\) 的等式：

\[\ln(1 - x) = -\sum_{i \ge 1}\frac{x^i}{i} \]

证明这个等式，我们考虑通过证明其导函数相同并且在 \(0\) 点取值相同来证。

首先在 \(0\) 点取值都为 \(0\)，考虑对两边求导（注意需要用链式法则）：

\[\begin{aligned} \ln(1-x)' &= -\frac{1}{1-x}\\ (-\sum_{i \ge 1}\frac{x^i}{i})' &= -\sum_{i \ge 1}\frac{ix^{i-1}}{i}\\ &= -\sum_{i\ge 0}x^i\\ &= -\frac{1}{1-x} \end{aligned} \]

从而等式得证。

然后开始正文。

我们依然从导数的角度考虑，我们通过对两边计算导数得到：

\[B'(x) \equiv \frac{A'(x)}{A(x)} \pmod {x^n} \]

由于我们可以轻松地从 \(B'(x)\) 积分推出 \(B(x)\)，所以我们需要多项式积分，求导和求逆，前两个都很简单，但是第三个求逆要求 \(A_0 \not = 0\)。

假设 \(A_0 = 0\)，则我们将最小的系数不为 \(0\) 的项提出来，得到：

\[B(x) \equiv \ln(A_kx^k) + \ln(C(x)) \]

其中 \(C(x)\) 的常数项显然是 \(1\)。

然后处理一下前一项的影响就可以正常做了。

时间复杂度是多项式求逆的 \(O(n \log n)\)。

这里给出保证 \(a_0 = 1\) 的代码。

3.6 多项式指数函数

问题

设 \(n = \deg A(x) + 1\)，求 \(B(x) \equiv \exp(A(x)) \pmod {x^n}\)

方法

依然考虑倍增。不妨设 \(F(x) \equiv B(x) \pmod {x^m}\)，求 \(B(x) \bmod x^{2m}\)。

不妨考虑 \(\ln(B(x))\) 这个函数在 \(f_0(x)\) 这一点的 taylor 展开式对 \(x^{2m}\) 取模：

\[\ln(B(x)) \equiv \ln(f_0(x)) + \ln'(f_0(x))(B(x) - f_0(x)) \pmod {x^{2m}} \]

后面的无穷项由于 \((B(x) - f_0(x))\) 平方后最小的非零项至少是 \(x^m\)，对 \(x^{2m}\) 取模均会变成 \(0\)，故全部被忽略掉了。

由于 \(A(x) = \ln(B(x))\)，我们整理一下可以得到：

\[A(x) \equiv \ln(f_0(x)) + f_0^{-1}(x)(B(x) - f_0(x)) \pmod {x^{2m}} \]

最后整理得到：

\[B(x) \equiv f_0(x)[A(x) - \ln(f_0(x)) + 1] \pmod {x^{2m}} \]

多项式乘法和多项式求逆均是 \(O(n \log n)\)，所以总复杂度还是 \(O(n \log n)\)。

这里给出保证 \(a_0 = 0\) 的代码：

3.7 多项式快速幂

问题

给定 \(A(x)\) 和 \(k\)，求 \(A^k(x) \bmod {x^n}\)。

方法

我们考虑先求 \(\ln\) 再求 \(\exp\)：

\[\begin{aligned} A^k(x) &\equiv \exp(\ln(A^k(x))) \pmod {x^n}\\ A^k(x) &\equiv \exp(k\ln(A(x))) \pmod {x^n}\\ \end{aligned} \]

做一次 \(\ln\) 再做一次 \(\exp\) 即可。

3.8 多项式开根

问题

给定 \(A(x)\) 和 \(k\)，求 \(A^{\frac{1}{k}}(x) \bmod {x^{n}}\)。

方法

首先不一定总是存在，我们有一个充要条件：设最小非零项是 \(d\)，则 \(k | d\) 且 \(a^{\frac{1}{k}}\) 存在。

然后我们沿用上一个快速幂的套路可以得到：

\[A^{\frac{1}{k}}(x)\equiv \exp(\frac{1}{k}\ln(A(x))) \pmod {x^n} \]

和快速幂本质上是相同的。

3.9 多项式开方

问题

上一个问题中 \(k = 2\) 的特例。

方法

依然考虑倍增，不妨设 \(B(x)^2 \equiv A(x) \pmod {x^n}\)，假设我们已经知道了 \(f(x) \equiv B(x) \pmod {x^{m}}\)，则：

\[\begin{aligned} f^2 &\equiv A \pmod{x^m}\\ f^2 - A &\equiv 0 \pmod{x^m}\\ (f^2 - A)^2 &\equiv 0 \pmod{x^{2m}}\\ (f^2 + A)^2 &\equiv 4f^2A \pmod {x^{2m}}\\ (\frac{f^2 + A}{2f})^2 &\equiv A \pmod {x^{2m}}\\ \frac{f^2 + A}{2f} &\equiv B \pmod {x^{2m}}\\ \frac{1}{2}f + {\frac{1}{2}}Af^{-1} &\equiv B \pmod {x^{2m}} \end{aligned} \]

然后就可以倍增计算，时间复杂度 \(O(n \log n)\)。

4. 多项式全家桶代码

目前的代码：

点击查看代码

#include <iostream>
#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;
const int N = 1e5 + 5;
const int mod = 998244353;


int fpow(int a, int b, int p) {//快速幂取模 
	if (b == 0)
		return 1;
	int ans = fpow(a, b / 2, p);
	ans = 1ll * ans * ans % p;
	if (b % 2 == 1)
		ans = 1ll * a * ans % p;
	return ans;
}

int mmi(int a, int p) {//求逆元 
	return fpow(a, p - 2, p);
}

namespace Poly { 
	typedef vector<int> poly;
	int rev[N * 4] = {0};

	poly operator+(poly x, int y) {
		if (x.size() == 0)
			x.resize(1);
		x[0] += y;
		return x;
	}
	
	poly operator+(poly x, poly y) {//重载多项式加法 
		poly ans = poly(max(x.size(), y.size()));
		for (int i = 0; i < (int)x.size(); i++)
			ans[i] = (ans[i] + x[i]) % mod;
		for (int i = 0; i < (int)y.size(); i++)
			ans[i] = (ans[i] + y[i]) % mod;
		while (ans.size() > 0 && ans.back() == 0)
			ans.pop_back();
		return ans;
	}
	
	poly operator-(poly x, poly y) {//重载多项式加法 
		poly ans = poly(max(x.size(), y.size()));
		for (int i = 0; i < (int)x.size(); i++)
			ans[i] = (ans[i] + x[i]) % mod;
		for (int i = 0; i < (int)y.size(); i++)
			ans[i] = (ans[i] - y[i] + mod) % mod;
		while (ans.size() > 0 && ans.back() == 0)
			ans.pop_back();
		return ans;
	}
	
	void operator+=(poly &x, poly y) {
		x.resize(max(x.size(), y.size()));
		for (int i = 0; i < (int)y.size(); i++)
			x[i] = (x[i] + y[i]) % mod;
	}
	
	void operator-=(poly &x, poly y) {
		x.resize(max(x.size(), y.size()));
		for (int i = 0; i < (int)y.size(); i++)
			x[i] = (x[i] - y[i] + mod) % mod;
	}
	
	void DFT(poly &a, int lim, int op) {//op = 1 DFT 否则是 IDFT 
		for (int i = 0; i < lim; i++)
			if (i < rev[i])
				swap(a[i], a[rev[i]]);
		for (int i = 1; i < lim; i <<= 1) {
			int wn = fpow((op == 1 ? 3 : (mod + 1) / 3), (mod - 1) / (i << 1), mod);
			for (int j = 0; j < lim; j += (i << 1)) {
				int w = 1;
				for (int k = 0; k < i; k++, w = 1ll * w * wn % mod) {
					int x = a[j + k], y = 1ll * w * a[i + j + k] % mod;
					a[j + k] = (x + y) % mod, a[i + j + k] = (x - y + mod) % mod;
				}
			}
		}
		if (op == -1) {
			int Inv = mmi(lim, mod);
			for (int i = 0; i < lim; i++)
				a[i] = 1ll * a[i] * Inv % mod;
		}
	}
	
	poly NTT(poly f, poly g) {//多项式NTT，返回两个多项式模意义下的卷积 
		int len = (int)f.size() + (int)g.size() - 1;
		int lim = 1, pw = 0;
		while (lim < len)
			lim <<= 1, pw++;
		f.resize(lim);
		g.resize(lim);
		for (int i = 0; i < lim; i++)
			rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (pw - 1)));
		DFT(f, lim, 1), DFT(g, lim, 1);
		for (int i = 0; i < lim; i++)
			f[i] = 1ll * f[i] * g[i] % mod;
		DFT(f, lim, -1);
		while (f.size() > 0 && f.back() == 0)
			f.pop_back();
		return f;
	}
	
	poly poly_inv(poly f, int M) {//求多项式的逆，需要处理一下2的幂 
		if (M == 1) {
			poly ans = poly(1);
			ans[0] = mmi(f[0], mod);
			return ans;
		}
		poly ans = poly_inv(f, (M + 1) / 2);
		int lim = 1, pw = 0;
		while (lim < (M << 1))
			lim <<= 1, pw++;
		for (int i = 0; i < lim; i++)
			rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (pw - 1)));
		poly g = f;
		g.resize(M), g.resize(lim), ans.resize(lim);
		DFT(g, lim, 1), DFT(ans, lim, 1);
		for (int i = 0; i < lim; i++)
			ans[i] = 1ll * ans[i] * (2 - 1ll * g[i] * ans[i] % mod + mod) % mod;
		DFT(ans, lim, -1);
		ans.resize(M);
		return ans;
	}
	
	poly poly_sqr(poly f, int M) {//多项式开根号 
		if (M == 1) {
			poly ans = poly(1);
			ans[0] = 1;
			return ans;
		}
		poly ans = poly_sqr(f, (M + 1) / 2);
		int lim = 1, pw = 0;
		while (lim < (M << 1))
			lim <<= 1, pw++;
		for (int i = 0; i < lim; i++)
			rev[i] = ((rev[i >> 1] >> 1) | ((i & 1) << (pw - 1)));
		poly inv = poly_inv(ans, M);
		poly g = f;
		g.resize(M), g.resize(lim);
		inv.resize(lim), ans.resize(lim);
		DFT(g, lim, 1), DFT(inv, lim, 1);
		for (int i = 0; i < lim; i++)
			inv[i] = 1ll * inv[i] * g[i] % mod;
		DFT(inv, lim, -1);
		for (int i = 0; i < lim; i++)
			ans[i] = 1ll * (mod + 1) / 2 * (ans[i] + inv[i]) % mod;
		ans.resize(M);
		return ans;
	}
	
	poly poly_rev(poly f) {//反转系数
		reverse(f.begin(), f.end());
		return f;
	}
	
	pair<poly, poly> poly_div(poly f, poly g) {//多项式带余除法 
		int n = (int)f.size(), m = (int)g.size();
		if (n < m) 
			return make_pair(poly(0), f);
		poly q, r;
		q = NTT(poly_rev(f), poly_inv(poly_rev(g), n - m + 1));
		q.resize(n - m + 1);
		q = poly_rev(q);
		r = f - NTT(g, q);
		return make_pair(q, r);
	}
	int inv[N * 4] = {0};
	
	poly deriv(poly f) {//求导 
		poly ans = poly(f.size() - 1);
		for (int i = 0; i < (int)ans.size(); i++)
			ans[i] = 1ll * f[i + 1] * (i + 1) % mod;
		return ans;
	}
	
	poly intrg(poly f) {//积分 
		poly ans = poly(f.size() + 1);
		inv[1] = 1;
		for (int i = 2; i < (int)ans.size(); i++) 
			inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod; 
		for (int i = 1; i < (int)ans.size(); i++) 
			ans[i] = 1ll * f[i - 1] * inv[i] % mod;
		return ans;
	} 
	
	poly poly_ln(poly f, int M) {//多项式ln 
		poly b = intrg(NTT(deriv(f), poly_inv(f, M)));
		b.resize(M);
		return b;
	}
	
	poly poly_exp(poly f, int M) {//求多项式exp 
		if (M == 1) {
			poly ans = poly(1);
			ans[0] = 1;
			return ans;
		}
		poly ans = poly_exp(f, (M + 1) / 2);
		poly g = f;
		g.resize(M);
		poly Ln = poly_ln(ans, M);
		ans = NTT(ans, g - Ln + 1);
		ans.resize(M);
		return ans;
	}
	
	poly poly_fpw(poly f, int k) {
		f = poly_ln(f, (int)f.size());
		for (int i = 0; i < (int)f.size(); i++)
			f[i] = 1ll * f[i] * k % mod;
		f = poly_exp(f, (int)f.size());
		return f;
	}
} 

int main() {
	
	return 0;
} 

5. 多项式进阶操作

5.1 多项式平移

5.1.2 普通多项式

有一个多项式 \(f(x) = \sum_{i = 0}^na_ix^i\)，求 \(f(x + c)\) 的系数。

我们直接暴力拆项得到：

\[\begin{aligned} f(x + c) &= \sum_{i=0}^na_i(x + c)^i\\ &= \sum_{i=0}^na_i\sum_{j=0}^i\binom{i}{j}x^jc^{i-j}\\ &= \sum_{j=0}^n\frac{x^j}{j!}\sum_{i=j}^ni!a_i\frac{c^{i-j}}{(i - j)!}\\ &= \sum_{i=0}^n\frac{x^i}{i!}\sum_{j=i}^nj!a_j\frac{c^{j-i}}{(j - i)!}\\ \end{aligned} \]

这是一个变形的卷积，设后面的和式为 \(ans_i\)，我们要求 \(ans_0 \sim ans_n\)，考虑化成卷积的形式，不妨令 \(k = j - i\)，则我们得到：

\[\begin{aligned} ans_i &= \sum_{k=0}^{n-i}(k+i)!a_{k+i}\frac{c^k}{k!}\\ &= \sum_{j + k = n- i}(n - j)!a_{n - j}\frac{c^k}{k!} \end{aligned} \]

令 \(a_j' = a_{n-j}(n-j)!\)，并且 \(ans_i' = ans_{n-i}\)，则：

\[ans_i' = \sum_{j + k = i}a_j' \frac{c^k}{k!} \]

这样就是标准的卷积形式，可以再 \(O(n \log n)\) 内求解了。

5.1.2 下降幂多项式

有一个下降幂多项式 \(f(x) = \sum_{i=0}^nb_ix^{\underline{i}}\)，求 \(f(x + c)\) 的下降幂系数。

我们需要用到下降幂的二项式定理：

\[(x + y)^{\underline{n}} = \sum_{i=0}^n\binom{n}{i}x^{\underline{i}}y^{\underline{n-i}} \]

所以我们还是类似上面的暴力拆项得到：

\[\sum_{i=0}^n\frac{x^{\underline{i}}}{i!}\sum_{j = i}^nj!b_j\frac{c^{\underline{j - i}}}{(j - i)!} \]

所以我们设后面的和式是 \(ans_i\)，并且 \(b_j' = (n - j)!b_{n - j}, ans_i' = ans_{n-i}\)，则：

\[ans_i' = \sum_{j + k = i}b_j' \frac{c^{\underline{k}}}{k!} \]

还是卷积的形式。

5.2 下降幂多项式的操作

5.2.1 点值

考虑如何快速求下降幂多项式再 \(c \sim n + c\) 的点值。

首先，只要我们可以求出再 \(0 \sim n\) 的点值，就可以通过平移转化到这个问题，于是我们考虑求 \(0 \sim n\) 的值。

推式子：

\[\begin{aligned} y_k &= f(k)\\ &= \sum_{i=0}^nb_ik^{\underline{i}}\\ &= \sum_{i=0}^kb_i\frac{k!}{(k-i)!} \end{aligned} \]

两边同时除以 \(k!\) 得到：

\[\frac{y_k}{k!} = \sum_{i=0}^kb_i\frac{1}{(k-i)!} \]

刚好是卷积的形式，\(O(n \log n)\) 计算。

5.2.2 插值

假设知道了 \(c \sim n + c\)，求 \(f(x + c)\)。

同样还是转化成 \(0 \sim n\) 的插值，考虑上面推出的式子，令 \(Y(x) = \sum_{i}\frac{y_i}{i!}x^i,B(x) = \sum_{i}b_ix^i\)，可以得到：

\[Y(x) \equiv B(x)e^{x} \pmod {x^{n+1}} \]

则：

\[B(x) \equiv Y(x)e^{-x}\pmod {x^{n+1}} \]

所以还是卷积就可以了。

5.2.3 乘法

两个下降幂多项式的乘法可以先点值再插值即可。时间复杂度 \(O(n \log n)\)

6. 题目

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

显然用 EGF 来求解。不难发现对于不连通的情况方案数为 \(2^{\binom{n}{2}}\)，则设 \(f_n = 2^{\binom{n}{2}}\)，\(g_n\) 为答案。

则对于它们的 EGF 满足 \(F(x) = \exp(G(x))\)，运用一次多项式 exp 即可求出答案。

当然这道题也可以用别的方法来做，考虑它们之间满足的关系式：

\[g_n = \sum_{i=1}^{n}\binom{n-1}{i-1}f_{i}g_{n-i} \]

则我们拆开组合数：

\[\frac{g_n}{(n-1)!} = \sum_{i=1}^n\frac{f_{i}}{(i - 1)!}\frac{g_{n - i}}{(n - i)!} \]

所以我们不妨设 \(G(x) = \sum\frac{g_n}{n!}, H(x) = \sum\frac{g_n}{(n-1)!}, F(x) = \sum\frac{f_n}{(n-1)!}\)，则我们希望求出 \(F\)，可以用 \(F = H \times G^{-1}\) 求出。

这样就只用写一个多项式求逆即可。

P4389 付公主的背包

显然答案就是 \(\prod_{i=1}^n\frac{1}{1-x^{v_i}}\) 的系数，但是如果用多项式求逆和卷积直接算会爆炸，我们要考虑优化。

首先一个很重要的技巧就是看到乘积的时候想到用 \(\ln, \exp\) 转化成求和。

所以：

\[\begin{aligned} \prod_{i=1}^n\frac{1}{1-x^{v_i}} &= \exp(\ln(\prod_{i=1}^n\frac{1}{1-x^{v_i}}))\\ &= \exp(\sum_{i=1}^n\ln(\frac{1}{1-x^{v_i}}))\\ &= \exp(\sum_{i=1}^n\ln(\frac{1}{1-x^{v_i}}))\\ &= \exp(-\sum_{i=1}^n\ln(1-x^{v_i}))\\ &= \exp(\sum_{i=1}^n\sum_{j \ge 1}\frac{x^{v_ij}}{j})\\ \end{aligned} \]

不妨记录每种体积的物品个数，然后直接用 \(O(m \ln m)\) 的时间枚举倍数计算内部这个多项式，然后再求一次 \(\exp\) 即可。

CF438E The Child and Binary Tree

不妨设 \(f_n\) 表示权值和为 \(n\) 的答案，\(g_n\) 表示 \(n\) 权值可不可行，可以得到：

\[f_n = \sum_{i=1}^ng_i\sum_{j=0}^{n - i}f_if_{n-i-j} \]

考虑到 \(g_0 = 0, f_0 = 1\)，对于其生成函数我们可以得到：

\[F = F^2 \times G + 1 \]

考虑这个方程的解：

\[\frac{1 \pm \sqrt{1 - 4G}}{2G} \]

如何确定应该是哪个？我们考虑 \(x \to 0\) 时，如果取正号则 \(F \to \infty\)，显然不行；取负号则 \(F \to 1\)，满足 \(f_0 = 1\)，所以取负号。再整理一下得到：

\[\frac{2}{1 + \sqrt{1 - 4G}} \]

多项式求逆和开根即可。

AVL Trees Gym - 100341C

考虑 dp，设 \(f_{n,m}\) 表示总共 \(n\) 个节点，高度为 \(m\) 的 AVL 个数，我们可以得到：

\[f_{n,m} = \sum_{i + j = n - 1}(f_{i, m - 1}f_{j, m - 1} + f_{i, m - 1}f_{j, m - 2} + f_{i, m - 2}f_{j, m - 1}) \]

设生成函数 \(F_m(x) = \sum_nf_{n,m}x^n\)，则我们可以得到：

\[F_m(x) = x(F_{m - 1}^2(x) + 2F_{m-1}(x)F_{m-2}(x)) \]

考虑到 \(m \le 16\)，所以最大的情况下 \(F_m(x)\) 是一个 \(65536\) 项的多项式，我们求 \(n\) 的系数可以求 \(\omega^0 \sim \omega^{65535}\) 的值，然后用 IDFT 公式 \(a_t = \frac{1}{n}\sum_{0\le i < n}\omega^{-tj}y_j\) 求出 \(n\) 的系数即可。

每一个求值可以从 \(F_0(x)\) 开始推，所以时间复杂度是 \(O(h2^h)\)。

注意当 \(n \ge 2^h\) 是需要特判掉这种情况。

提交记录

HDU7057 Buying Snacks

同样，设 \(f_{n,m}\) 表示前 \(n\) 个花 \(m\) 块钱，转移是：

\[f_{n,m} = f_{n-1,m} + f_{n-1,m-1} + f_{n-1,m-2} + f_{n-2, m - 1} + 2f_{n-2, m - 2} + f_{n-2, m - 3} \]

由于这道题 \(m\) 比较小 \(n\) 比较大，考虑设 \(F_n(x) = \sum_mf_{n,m}x^m\)，则：

\[F_n(x) = (1 + x + x^2)F_{n-1}(x) + (x + 2x^2 + x^3)F_{n-2}(x) \]

由于最终答案可能是个次数很大的多项式，我们没办法只通过 \(m\) 个值来还原，所以我们考虑矩阵快速幂。

将矩阵存的数改成多项式，然后用 NTT 优化可以做到 \(O(m \log m \log n)\)，每次只保留后 \(m\) 位。

但是一次矩阵乘法不加优化会做 24 次 DFT 和 IDFT，所以我们先用 8 次将这八项的 DFT 算出来，剩下的直接用点值来算，最后再用 4 次 IDFT 还原就可以减小一半的常数。

提交记录

P6620 [省选联考 2020 A 卷] 组合数问题

考虑先将 \(f(x)\) 转化成下阶乘幂多项式 \(f(x) = \sum_{i=0}^mb_ik^{\underline{i}}\)。

直接推式子得到：

\[\begin{aligned} \sum_{k=0}^n\sum_{i=0}^mb_ik^{\underline{i}}x^k\binom{n}{k} \end{aligned} \]

对于这种既有下阶乘幂的也有组合数的，一般把两者全部拆成阶乘来处理，有一个好用的公式：

\[k^{\underline{i}}\binom{n}{k} = \binom{n - i}{k - i}n^{\underline{i}} \]

所以得到：

\[\begin{aligned} \sum_{k=0}^n\sum_{i=0}^mb_ik^{\underline{i}}x^k\binom{n}{k} &= \sum_{k=0}^n\sum_{i=0}^mb_ix^k\binom{n - i}{k - i}n^{\underline{i}}\\ &= \sum_{i=0}^mb_in^{\underline{i}}\sum_{k=i}^nx^k\binom{n - i}{k - i}\\ &= \sum_{i=0}^mb_in^{\underline{i}}x^i\sum_{k=0}^{n - i}x^k\binom{n - i}{k}\\ &= \sum_{i=0}^mb_in^{\underline{i}}x^i(x + 1)^{n - i}\\ \end{aligned} \]

用斯特林数 \(O(m^2)\) 求出下降幂多项式然后求值即可。

P6667 [清华集训2016] 如何优雅地求和

考虑到题目给的是点值，意味着我们可以求出多项式的下降幂表示，现在考虑如果 \(f(x) = x^{\underline{m}}\) 的情况，其他都是这个的线性组合。

\[\begin{aligned} Q &= \sum_{k=0}^n\binom{n}{k}x^k(1-x)^{n-k}k^{\underline{m}}\\ &= \sum_{k=0}^n\binom{n - m}{k - m}n^{\underline{m}}x^k(1-x)^{n-k}\\ &= n^{\underline{m}}\sum_{k=m}^n\binom{n - m}{k - m}x^k(1-x)^{n-k}\\ &= n^{\underline{m}}x^m\sum_{k=0}^{n - m}\binom{n - m}{k}x^k(1-x)^{(n - m)-k}\\ &= n^{\underline{m}}x^m(x + 1 - x)^{n-m}\\ &= n^{\underline{m}}x^m \end{aligned} \]

这就做完了。

但是这道题也可以用 PGF 概率生成函数来做，考虑随机变量 \(X\) 表示 \(n\) 次试验成功的次数，\(p\) 是一次实验成功的概率，则：

\[P(X = k) = p^k(1-p)^{n-k}\binom{n}{k} \]

所以：

\[Q = \sum_{k=0}^nf(k)P(X = k) = E(f(X)) \]

考虑 \(f(x) = x^{\underline{m}}\)，根据 PGF 相关知识就能得到答案是 \(n^{\underline{m}}x^{m}\)。

CF622F The Sum of the k-th Powers

太典了。连续点值求另一个点值。

考虑函数 \(f(x) = \sum_{i=1}^x x^k\)，由于其一阶差分是 \(f(x) = x^k\)，所以这个函数应该是一个 \(k + 1\) 次多项式。

我们显然希望通过插值来求出这个多项式，但是由于直接拉格朗日插值太慢了，我们要一些优化。

注意到点的选取不受限制，所以我们选取 \(1 \sim k + 2\)，作为点值。

而带入原来的插值公式就会发现后面的乘积的分母可以化成阶乘，而由于我们只求 \(n\) 点的值，所以分子可以用前后缀积预处理计算，这样时间复杂度就是 \(O(k \log k)\) 了。