多项式全家桶
一、前置芝士
1. 基本概念
- 多项式： 有限项相加的求和式 \(\sum a_nx^n\)，记作 \(f(x) = \sum a_nx^n\)
- 多项式的度(次数)：对于一个多项式 \(f(x)\)，称其最高次项的次数为该多项式的度（\(degree\)），也称次数，记作 \(\deg f\)
- 级数：将数列的项依次用加号连接起来的函数
- 幂级数：每项均为非负整数次幂函数乘以常数系数的级数称为幂级数
更多&拓展
对于多项式，我们也可以将其直接定义为一个系数序列，表示为：
此处我们认为，\(x\) 只是一个形式符号，即一个对于系数位置的标识符
若支持无穷项，则称 \(f(x)\) 为形式幂级数：
Q：多项式和级数的区别？
A：多项式是有限的,级数是无限的
\(e.g.~\sum\limits_{i = 0}^{\infty} a_nx^n\) 是幂级数但不是多项式
2. 多项式的表示方法
① 系数表示法
存储方式
存储时只存储每一位的系数即可
② 点值表示法
即把多项式放入平面直角坐标系中，得到一个函数图像
存储方式
通常地，带入 \(n\) 个 \(x\)，得到 \(n\) 个 \(y\)
计算多项式乘法时：点值表示法可以 \(O(n)\) 求，即每一位对应乘起来即可
而这两种表示法可以相互转化，但是暴力转化是 \(O(n^2)\) 的，那么我们就要使用一个神奇的 FFT 优化多项式乘法，并且它是 \(O(nlog_n)\) 的
3. 复指数函数 & 复三角函数
- 复指数函数：对于复数 \(z = x + i \times y\)，定义 \(e^z\) 为其复指数函数。特别地，我们通常将 \(e^z\) 记作 \(\exp z\)
在代数上 \(\exp z = e^x(\cos y + i\sin y)\)
证明
前置知识：欧拉定理
- 复三角函数：
4. 复数的三种形式
详细见 复数 - ricky_lin
① 代数形式
适用范围：计算复数的加减乘除
② 三角形式
③ 指数形式
适用范围：这两种形式用于计算复数的乘除两个运算以及后面的运算较为方便
三角形式 知识层面要求更低，但 指数形式 会更加地方便
5. 单位根
本小节需要掌握的前置知识点
- 复数
- 欧拉定理
① 定义和约定
我们定义 \(x^n = 1\) 在复数意义下的解是 \(n\) 次复根。
这 \(n\) 个解都为 \(n\) 次单位根
根据学过的知识 \(n\) 次单位根将单位圆分成 \(n\) 等份
那么 \(x^n = 1\) 解集表示为 \(\{\omega_n^k | k = 0,1,\cdots,n-1\}\).
一般地，我们默认 \(n\) 次单位根 \(\omega_n\) 表示：从 \(1\) 开始逆时针方向的第一个解，即上方的 \(\omega_n\)
② 性质和证明
- 性质1：\(\omega_n^n = 1\)
证明
由定义 $x^n = 1$ 可证- 性质2：\(\omega_n^k = \omega_{2n}^{2k}\)
证明
$$\omega_{2n}^{2k} = \exp \frac {2\times 2k\times \pi i} {2n} = \exp \frac {2k\pi i} n = \omega_{n}^{k}$$- 性质3：\(\omega_{2n}^{k+n} = -\omega_{2n}^k\)
证明
我们先需要知道一件事：
然后：
二、多项式乘法
上面我们讲到多项式若用点值表示，那么乘法是 \(O(n)\) 的，但是显然地，我们常用系数进行表示
然而,对于系数转点值，需要代入 \(n\) 个点，显然复杂度是 \(O(n^2)\) 的
然后我们考虑，找一些 \(x\) 使得 \(x^m = 1\)，这样的话就会计算 \(m\) 次乘方之后进入循环，不需要做全部的 \(m\) 次运算了
但是，显然地，这样的 \(x\) 在实数域内只有可能是 \(\pm 1\)，所以我们要找到一些特殊性质的数。
1. FFT——快速傅里叶变换
① 递归
第一种方法是将 \(x\) 拓展到复数域
令 \(f_1(x) = a_0 + a_2x + \cdots + a_{n-2}x^{\frac n 2 -1}\)，\(f_2(x) = a_1 + a_3x + \cdots + a_{n-1}x^{\frac n 2 -1}\)（保证 \(x\) 为偶数）
可以得到：\(f(x) = f_1(x^2) + xf_2(x^2)\)
代入 \(x = \omega_n^k\) 得：
代入 \(\omega_{n}^{k + \frac n 2} = -\omega_{n}^k\) 得：
然后便可以递归求解，复杂度 \(O(n\log n)\)
重要提示
考虑上面的式子的适用条件是每一次多项式的长度都是 \(2\) 倍数，那么原先的多项式的长度就必须是 \(2\) 的次幂，不满足的在后面补 \(0\) 即可
② 倍增
但是考虑递归常数过大，我们考虑倍增怎么做
递归的时候需要把系数按照 \(\{a_0,a_2,a_4,a_6,\cdots\}, \{a_1,a_3,a_5,a_7,\cdots\}\)来分组，并交换位置，我们将其称为 位逆序变换：
例子
以 \(8\) 项多项式为例：
| 状态 | 序列 |
|---|---|
| 初始序列 | \(\{a_0,a_1,a_2,a_3,a_4,a_5,a_6,a_7\}\) |
| 第一次二分 | \(\{a_0,a_2,a_4,a_6\},\{a_1,a_3,a_5,a_7\}\) |
| 第二次二分 | \(\{a_0,a_4\},\{a_2,a_6\},\{a_1,a_5\},\{a_3,a_7\}\) |
| 第三次二分 | \(\{a_0\},\{a_4\},\{a_2\},\{a_6\},\{a_1\},\{a_5\},\{a_3\},\{a_7\}\) |
根据上面的例子，我们可以找到一个规律：
初始序列的下标在二进制下 翻转便可以得到其在最终序列中的下标
例子
以 \(x_1\) 为例，\(1\) 在二进制下是 \(001(1)\) 最终序列中是 \(100(4)\)，同样的 \(x_4\) 也跑到了下标为 \(1\) 的位置
位逆序变换实现
位逆序置换可以 \(O(n)\) 从小到大递推实现：
我们设 \(R(x)\) 表示长度为 \(k\) 的二进制数 \(x\) 翻转后的数（高位补 \(0\)）。我们要求的是 \(R(0), R(1), \cdots , R(n − 1)\)。
首先 \(R(0) = 0\)。
我们从小到大求 \(R(x)\)。因此在求 \(R(x)\) 时，\(R(\lfloor\frac x 2\rfloor)\) 的值是已知的。因此我们把 \(x\) 右移一位（除以 2），然后翻转，再右移一位，就得到了 \(x\) 除了（二进制）个位之外
其它位的翻转结果，把最后一位 \(0/1\) 拼到最高位，就得到了 \(R(x)\)。
数学语言：
③ 蝶形运算优化
我们从上面知道，我们需要通过下面两个式子进行推导：
使用 位逆序置换 后，对于给定的 \(n,k\)：
- \(f_1(\omega_{\frac n 2}^k)\) 的值存储在下标为 \(k\) 的位置上，\(f_2(\omega_{\frac n 2}^k)\) 的值存储在下标为 \(\frac n 2 + k\) 的位置上
- \(f(\omega_n^k)\) 的值存储在数组下标为 \(k\)，\(f(\omega_n^{k + \frac n 2})\) 的值将存在数组下标为 \(k + \frac n 2\) 的位置
所以说我们每次变换的时候不需要再新开一个数组，只需要在原数组上进行覆写，即可。
这样我们就完成了 DFT——离散傅里叶变换
④ IDFT——快速傅里叶逆变换
在前面，我们将两个多项式用 \(O(nlog_n)\) 的复杂度，变成了点值表示法，并用 \(O(n)\) 的时间，求出了乘积。
现在，我们需要将点值表示法在 \(O(nlog_n)\) 的时间内，变为系数表示法。
方法即为：将 \(\omega_n^k\) 的共轭复数 \(\omega_n^0\) 带回到 DFT 中，再对其所有系数都除以 len(多项式长度)
证明
我们知道 \(\omega_n^k\) 的共轭复数为 \(\omega_{n}^{-k}\)
同时我们也知道了一个函数 \(f(x)\) 在 \(\omega_n^0,\omega_n^1,\cdots,\omega_n^{n-1}\) 的点值 \(g(x)\)，并设为 \(y_0,y_1,\cdots,y_{n-1}\),我们需要求出原函数的系数 \(a_0,a_1,\cdots,a_{n-1}\)
代入 \(x = \omega_n^{-k}\) 得到：
记 \(S(x) = \sum\limits_{i = 0}^{n-1}(\omega_n^x)^i\)，由等比数列求和得：
由上式得到，当 \(\omega_n^x = 1\) 即 \(x = 0\) 时，\(S(0) = \sum\limits_{i=0}^{n-1} 1^i = n\)，否则 \(S(x) = 0\)
代回原式：
code
#include<cstdio>
#include<vector>
#include<cmath>
using namespace std;
typedef long long ll;
const double PI = acos(-1);
struct Complex{
double x,y;
Complex(double x = 0,double y = 0): x(x),y(y){};
Complex operator + (const Complex &B){return Complex(x+B.x,y+B.y);}
Complex operator - (const Complex &B){return Complex(x-B.x,y-B.y);}
Complex operator * (const Complex &B){return Complex(x*B.x-y*B.y,x*B.y+y*B.x);}
};
typedef vector<Complex> Poly;
Poly A,B;
int len = 1,llog;
void FFT(vector<int> &rev,Poly &v,int inv){
for(int i = 0; i < len; ++i)
if(i < rev[i]) swap(v[i],v[rev[i]]);
for(int k = 1; k < len; k <<= 1){
Complex omega(cos(PI/k), inv * sin(PI/k));
for(int i = 0; i < len; i += k * 2){
Complex w(1,0);
for(int j = 0; j < k; ++j, w = w * omega){
Complex s = v[i+j], t = v[i+j+k] * w;
v[i+j] = s + t; v[i+j+k] = s - t;
}
}
}
// for(int i = 0; i < len; ++i) printf("%lf ",v[i].x);puts("");
if(inv == -1) for(int i = 0; i < len; ++i) v[i].x /= len;
}
int n,m;
int main() {
scanf("%d%d",&n,&m);
while(len <= n + m) len <<= 1,++llog;
A.resize(len);B.resize(len);
for(int i = 0; i <= n; ++i) scanf("%lf",&A[i].x);
for(int i = 0; i <= m; ++i) scanf("%lf",&B[i].x);
vector<int> rev(len);
for(int i = 0; i < len; ++i) rev[i] = (rev[i>>1] >> 1) | ((i&1) << (llog-1));
FFT(rev, A, 1), FFT(rev, B, 1);
for(int i = 0; i < len; ++i) A[i] = A[i] * B[i];
FFT(rev, A, -1);
for(int i = 0; i <= n + m; ++i) printf("%lld ",(ll)(A[i].x + 1e-3));
return 0;
}
⑤ 从线性代数角度理解FFT
DFT 本身是个线性变换，可以理解为将目标多项式当作向量，左乘一个矩阵得到变
换后的向量，以模拟把单位复根代入多项式的过程：
现在我们已经得到最左边的矩阵 \(A\) 了，中间的 \(x\) 值在目标多项式的点值表示中也是一一对应的
所以，根据矩阵的基础知识，我们只要在式子两边左乘中间那个大矩阵 \(B\) 的逆矩阵就行了。
由于这个矩阵的元素非常特殊，它的逆矩阵也有特殊的性质，就是每一项取倒数，再除以变换的长度 n，就能得到它的逆矩阵。
注：这里的长度 \(n\) 也指一个 \(2\) 的幂，不足则补 \(0\)。
2. NTT——快速数论变换
由于 FFT 算 \(\omega\) 的时候精度极有可能爆炸，导致一些神奇问题，所以我们需要一些精度更高的算法，而且我们更多时候是考虑的取模意义下的运算。
于是我们便有了 FFT
① 补充亿点芝士
阶
- 若 \(a,p\) 互素，且 \(p > 1\)，对于 \(a^n \equiv 1 (\bmod p)\) 最小的 \(n\)，我们称之为 \(a\) 模 \(p\) 的阶，记作 \(\delta_p(a)\)，例如 \(\delta_7(2) = 3\)
性质：\(a,a^2,\cdots,a^{\delta_p(a)}\) 模 \(p\) 两两不同余
证明
反证法：
若 \(\exists~i\neq j~且~i,j\leq \delta_p(a)~且~a^i \equiv a^j (\bmod p)\)，则有 \(a^{|i-j|} \equiv 1 (\bmod p)\)
显然地，\(0 < |i-j| < \delta_p(a)\)，这和阶的定义矛盾，故原命题成立。
原根
设 \(p\in\mathbb{Z^+},a\in\mathbb Z\)，若 \(\delta_p(a) = \varphi(p)\)，则称 \(a\) 为模 \(p\) 的一个原根，记为 \(g\)
\(\delta_7(3) = 6 = \varphi(7)\)，因此 \(3\) 是模 \(7\) 的一个原根。
换一个角度理解
就是说 \(a^{\varphi(p)} \equiv 1 \bmod p\) 且 \(\forall i \in \mathbb{Z},1 \leq i < \varphi(p)~\text{满足}~a^i\not\equiv 1\bmod p\)，则 \(a\) 为 \(p\) 的原根
可以得到一个结论
原根判定定理 & 原根计算
设 \(m \geq 3\)，\(\gcd(g, m) = 1\)，则 \(g\) 是模 \(m\) 的原根的充要条件是：
对于 \(\varphi(m)\) 的每个素因数 \(p\)，都有原根。
没什么快速求法，只能暴力枚举判断。若素数 \(p\) 有原根，其最小原根 \(g_p = O (p^{0.25+\epsilon})\)，
其中 \(\epsilon > 0\)，这保证了我们暴力找一个数的最小原根，复杂度是可以接受的。
通常模数常见的有 \(998244353, 1004535809, 469762049\)。这几个的原根都是 \(3\)。
我们用原根可以代替 FFT 中的 \(\omega\)，因为它满足虚数 \(\omega\) 的所有性质（模意义下）
证明
- 性质1：\(\omega_n^0 = \omega_n^n = 1\)
显然地 \(\omega_n^0 = g^0 = 1\)，\(\omega_n^n = 1\)
根据费马小定理 \(g^{p-1} \equiv 1(\bmod p)\)，得证 \(\omega_n^n \equiv g^{p-1} (\bmod p)\)，得证。
- 性质2：\(\omega_n^k = \omega_{2n}^{2k}\)
由定义得：\(\omega_{2n}^{2k} = (g^{\frac {p-1} {2n}})^{2k} = g^{\frac {k(p-1)} {n}} = (g^{\frac {p-1} n})^k = \omega_n^k\)
- 性质3：\(\omega_n^{k + \frac n 2} = -\omega_n^k\)
\[(\omega_n^\frac n 2) ^ 2 = g^{p-1} \equiv 1 \bmod p \]又因为 \(p\in prime\)，所以 \(g^{0\sim p-1} \bmod p\) 两两不同
所以 \(g^{\frac {p-1} 2} \equiv -1 \bmod p\)
即证：\(g^{k+\frac {p-1} 2} = -g^k\)
在实现 NTT 的时候，把 DFT 中的 \(\omega_n^k = \cos (\frac k n 2\pi) + \sin (\frac k n 2\pi)\) 替换为 \(g^{\frac {k(p−1)} n} \bmod p\)，IDFT 就是将 \(g \rightarrow g^{-1}\) 即可。
例题：P3723 [AH2017/HNOI2017] 礼物
给定两个长为 n 的数列 a, b，a, b 均可循环移动，你需要选择一个整数 c，最小化
\[\sum_{i=1}^n (a_i-b_i+c)^2 \]\(n\leq 5\times 10^4,a_i,b_i\leq 100\)
solution
把柿子拆一拆： $$ \begin{aligned} & \sum_{i=1}^n (a_i-b_i+c)^2\\ = & \sum_{i=1}^n a_i^2+b_i^2+c^2-2a_ib_i + 2c(a_i-b_i)\\ = & \sum_{i=1}^n (a_i^2+b_i^2) + 2c\sum_{i=1}^n(a_i-b_i) + nc^2 - 2\sum_{i=1}^na_ib_i \end{aligned} $$ 前三项都是关于一个 $c$ 的二次函数，那么前三项通过 $c$ 最小化即可，现在我们需要做的，便是最大化$\sum\limits_{i=1}^n a_ib_i$我们先倍长 \(b\)，断环为链：
\[\sum\limits_{i=1}^n a_ib_{i+k} \]因为此处 \(a,b\) 之差不变且为 \(k\)，所以我们使用减法卷积
将 \(a\) 数组逆序，得到：
\[A_{n-k+1} = \sum\limits_{i=1}^n a_{n-i+1}b_{i+k} \]然后直接
NTT即可
三、多项式基本运算
1. 前置知识
本节只列出定义和常见公式，详细请自行查询相关资料学习
① 求导
② 泰勒展开
对于一个函数 \(f(x)\) 以及一个点 \(x_0\)，我们在 \(x_0\) 处对函数 \(f\) 进行一个拟合，设拟合函数为 \(T\)，那么泰勒展开的一般形式如下：
2. 多项式牛顿迭代
牛顿迭代解决了以下 problem:
给定多项式 \(g(x)\)，已知有 \(f(x) 满足：\)
\[g(f(x)) \equiv 0~~(\bmod x^n) \]求出模 \(x^n\) 意义下的 \(f(x)\)
solution
考虑倍增首先当 \(n = 1\) 是时 \([x^0]g(f(x)) = 0\) 的解需要单独求出。
假设现在已经得到了模 \(x^{\lceil\frac n 2\rceil}\) 意义下的解 \(f_0(x)\)，要求模 \(x^n\) 意义下的解 \(f(x)\)
将 \(g(f(x))\) 在 \(f_0(x)\) 处进行泰勒展开，有：
\[\sum\limits_{i=0}^{+\infty}\frac {g^{(i)}(f_0(x))} {i!}(f(x)-f_0(x))^i \equiv 0~~(\bmod x^n) \]因为 \(f(x) - f_0(x)\) 的最低非零项次数最低为 \(\lceil\frac n 2 \rceil\)，所以：
\[\forall i\geq2,(f(x)-f_0(x))^i\equiv 0~~~(\bmod x^n) \]那么：
\[\begin{aligned} & \sum\limits_{i=0}^{+\infty}\frac {g^{(i)}(f_0(x))} {i!}(f(x)-f_0(x))^i \\ \equiv & g'(f_0(x))\cdot(f(x)-f_0(x))\\ \equiv & g(f_0(x)) + g'(f_0(x))(f(x)-f_0(x))\\ \equiv & 0~~(\bmod x^n) \end{aligned} \]\[f(x) \equiv f_0(x) - \frac {g(f_0(x))} {g'(f_0(x))}~~(\bmod x^n) \]
3. 多项式求逆
设给定函数为 \(h(x)\)，求 \(f(x)\) 使得 \(f(x) \cdot h(x) \equiv 1~~(\bmod x^n)\)
solution
由上式可得： $$ g(f(x)) = \frac 1 {f(x)} - h(x) \equiv 0~~(\bmod x^n) $$ 使用牛顿迭代得到： $$ \begin{aligned} f(x) & \equiv f_0(x) - \frac {\frac 1 {f_0(x)} - h(x)} { - \frac 1 {f_0^2(x)}} & (\bmod x^n) \\ & \equiv 2f_0(x) - f_0^2(x)h(x) & (\bmod x^n) \end{aligned} $$ 时间复杂度： $$ T(n) = T(\frac n 2) + O(n\log n) = O(n\log n) $$
4. 多项式求ln、exp
① 定义
我们先了解 \(e^{f(x)}/exp(f(x))\) 在多项式中的含义是什么
我们先泰勒展开，可以得到：
\(exp\) 代表任意多个相应组合对象的组合。（\(f(x)^i\) 即为将 \(i\) 个对象组合；\(i!\) 消除了各组 之间的差异，避免重复计数）
\(ln\) 便是 \(exp\) 的逆运算。
② 多项式求ln
设给定函数为 \(h(x)\)，求 \(f(x)\) 使得 \(f(x) = ln~h(x)~~(\bmod x^n)\)
两边同时求导，得：
最后使用多项式求逆，并且积分积回去即可：
③ 多项式exp
设给定函数为 \(h(x)\)，求 \(f(x)\) 使得 \(f(x) \equiv e^{h(x)}~~(\bmod x^n)\)
由上式可得：
使用牛顿迭代得：
时间复杂度：
5. 多项式开方
设给定函数为 \(h(x)\)，求 \(f(x)\) 使得 \(f^2(x) \equiv h(x)~~(\bmod x^n)\)
由上式可得：
使用牛顿迭代得：
时间复杂度：
边界条件需要求一个数在模意义下的二次剩余，需使用 \(Cipolla\) 算法
6. 多项式除法&取模
给定多项式 \(f(x),g (x)\)，求 \(g(x)\) 除 \(f(x)\) 的商 \(Q(x)\) 和余数 \(R(x)\)，即
\(f(x) = g(x)\cdot Q(x) + R(x)\)。
发现若能消除 \(R(x)\) 的影响则可直接多项式求逆解决。
考虑构造变换：
即，\(f^R(x)\) 表示反转 \(f(x)\) 的系数，\(\deg f\) 表示多项式 \(f(x)\) 的度
设 \(n = \deg f\)
本文来自博客园，作者：ricky_lin，转载请注明原文链接：https://www.cnblogs.com/rickylin/p/blog_summary_polynomial.html
