多项式学习笔记
- 基础操作
- FFT & NTT
- 分治 FFT & NTT
- 拉格朗日插值
- 多项式求导与积分
- 多项式求逆
- 多项式开根
- 多项式带余除法
- 多项式对数函数
- 多项式指数函数
- 参考代码
基础操作
FFT & NTT
基础知识。
考虑如何快速计算多项式乘法,发现常规的做法是没有前途的,考虑换种方法计算。
有一种计算方法是这样的(设要计算 \(A\times B\) ):
- 将 \(A\) 转化为点值表示。
- 将 \(B\) 转化为点值表示。
- \(O(n)\) 地计算点值表示的乘积 \(C\)。
- 将 \(C\) 转化为数值表示。
复杂度瓶颈为点值表示与数值表示的转化,朴素做是 \(O(n^2)\)。
为了加速这个过程,FFT 引入了单位复根 \(\omega\),NTT 引入了原根 \(g\)。
单位复根 \(\omega\) 是 \(x^n=1\) 的解集,其满足如下性质:
- \(\omega^n_n=1\);
- \(\omega^k_n=\omega^{2k}_{2n}\);
- \(\omega^{k+n}_{2n}=-\omega^k_{2n}\)。
考虑分治,对一个多项式 \(f(x)\) 奇偶性分组处理,比如一个多项式:
\[f(x)=a_0+a_1x+a_2x^2+a_3x^3+a_4x^4+a_5x^5+a_6x^6+a_7x^7 \]考虑奇偶性分组后的状态:
\[f(x)=a_0+a_2x^2+a_4x^4+a_6x^6+x(a_1+a_3x^2+a_5x^4+a_7x^6) \]建立新函数:
\[G(x)=a_0+a_2x+a_4x^2+a_6x^3\\H(x)=a_1+a_3+a_5x^2+a_7x^3 \]那么待求解的 \(f(x)\) 可改为:
\[f(x)=G(x^2)+xH(x^2) \]带入 \(\omega\) 可得:
\[\begin{aligned} \operatorname{DFT}(f(\omega_n^k)) &=\operatorname{DFT}(G((\omega_n^k)^2)) + \omega_n^k \times \operatorname{DFT}(H((\omega_n^k)^2))\\ &=\operatorname{DFT}(G(\omega_n^{2k})) + \omega_n^k \times \operatorname{DFT}(H(\omega_n^{2k}))\\ &=\operatorname{DFT}(G(\omega_{n/2}^k)) + \omega_n^k \times \operatorname{DFT}(H(\omega_{n/2}^k)) \end{aligned} \\ \begin{aligned} \operatorname{DFT}(f(\omega_n^{k+n/2})) &=\operatorname{DFT}(G(\omega_n^{2k+n})) + \omega_n^{k+n/2} \times \operatorname{DFT}(H(\omega_n^{2k+n}))\\ &=\operatorname{DFT}(G(\omega_n^{2k})) - \omega_n^k \times \operatorname{DFT}(H(\omega_n^{2k}))\\ &=\operatorname{DFT}(G(\omega_{n/2}^k)) - \omega_n^k \times \operatorname{DFT}(H(\omega_{n/2}^k)) \end{aligned} \]然后就可以分治做了。
诶等等,递归常数特别大呢!
为了不递归,我们考虑一开始就预处理出来每一项最后到了什么地方,这个做法是位逆序置换(又名蝴蝶置换)。
考虑找规律,发现将下标二进制表示后翻转即可,设 \(R_i\) 为 \(i\) 变到的地方,那么可以考虑这样求:
- 右移一位。
- 翻转。
- 再右移一位。
- 处理最高位。
那么 \(R_i\) 的递推式可以轻易的写成:
\[R_i=\lfloor\frac{R_{\lfloor\frac{i}{2}\rfloor}}{2}\rfloor+(i\bmod2)\times\frac{len}{2} \]以上就是如何将数值表示转化为点值表示了,所以怎样转回来呢?
直接反转后除以 \(len\),或是取单位根的倒数跑即可,具体推导过程见 OI-Wiki。
模板:
void FFT(int n,cp* x,int op) {
for(int i=0;i
NTT 呢?
首先解释 \(g\),\(g\) 是原根,在 \(\bmod p\) 的情况下 \(\forall i\not =j,g^i\not\equiv g^j\pmod p\)。
我们用 \(g^{\frac{\varphi(p)}{n}}\) 代替 \(\omega\) 进行计算,由 \(\varphi(p)=p-1\) 且 \(n=2^q\),所以强制要求 \(p=k\times 2^q+1\),其中 \(2^q\ge n\)。
较 FFT 精度高,但却对模数有着致命的要求。
时间复杂度 \(O(n\log n)\)。
分治 FFT & NTT
问题:求 \(f_n\),其中:
\[f_i=\sum^{i=1}_{j=0} f_jg_{i-j} \]考虑 CDQ 分治,对于一个区间 \([l,r]\),将其划分成 \([l,mid],[mid+1,r]\) 两段。
先处理 \([l,mid]\),再计算 \([l,mid]\) 对 \([mid+1,r]\) 的贡献,再处理 \([mid+1,r]\)。
中间的计算使用 NTT/FFT 解决,时间复杂度 \(O(n\log^2n)\)
拉格朗日插值
问题,给定 \(n\) 个点 \((x_i,y_i)\),确定一个多项式 \(y=f(x)\),求出 \(f(k)\) 的值。
显然的事实:
\[f(x)\equiv f(y)\pmod{x-y} \]证明:\(f(x)-f(y)=(a_0-a_0)+\color{red}a_1(x-y)\color{black}+a_2(x^2-y^2)+a_3(x^3-y^3)+\cdots\)。
那么我们就可以列一个线性同余方程组:
\[\begin{cases}f(x)\equiv y_1\pmod {x-y_1}\\f(x)\equiv y_2\pmod{x-y_2}\\\ldots\\f(x)\equiv y_n\pmod{x-y_n}\end{cases} \]直接 CRT,所有模数的积为:
\[\prod_{i=1}^n (x-y_i) \]那么 \(m_i=\prod _{j\not=i} (x-y_i),m_i^{-1}=\prod _{j\not=i} \frac {1}{x_i-x_j}\)。
所以:
\[f(x)=\sum^n_{i=1} y_i \prod_{j\not=i} \frac {x-x_j}{x_i-x_j} \]朴素实现复杂度是 \(O(n^2)\) 的。
在所取点值连续的时候,即 \(x_i=i\) 的时候,可以直接用前缀积优化到 \(O(n)\)。
在需要动态加点的时候,可以使用重心拉格朗日插值法:
回顾之前的式子:
\[f(x)=\sum^n_{i=1} y_i \prod_{j\not=i} \frac {x-x_j}{x_i-x_j} \]设 \(g=\prod^n_{i=1} k-x_i\),那么:
\[f(x)=g\sum^n_{i=0}\prod\frac {y_i}{(k-x_i)(x_i-x_j)} \]设 \(t_i=\frac {y_i}{\prod_{j\not=i} x_i-x_j}\),那么:
\[f(x)=g\sum^n_{i=0}\frac {t_i}{k-x_i} \]加入一个点只需计算 \(t_i\) 即可。
求自然数幂和可以直接拉格朗日插值,时间复杂度 \(O(n)\)。
多项式求导与积分
由于求导和积分满足可加性,所以拆成单项式分别求导或积分。
如果你学过微积分相关理论,会知道:
\[A'(x)=\sum^n_{i=1}ia_ix^{i-1}\\\int A(x)dx=(\sum^n_{i=0}\frac {a_i}{i+1}x^{i+1})+c \]直接做即可,为数不多的 \(O(n)\) 时间复杂度。
多项式求逆
给定 \(G(x)\),求一个 \(F(x)\) 使得 \(F(x)G(x)\equiv0\pmod {x^n}\)。
首先假设我们已经求出了一个 \(f(x)\) 使得 \(f(x)G(x)\equiv 0\pmod{x^{n/2}}\),考虑如何求出 \(F(x)\)。
因为 \(F(x)G(x)\equiv 0\pmod{x^n}\),所以必然 \(F(x)G(x)\equiv 0 \pmod{x^{n/2}}\),两式相减,得 \(F(x)-f(x)\equiv0\pmod{x^{n/2}}\),此处 \(G(x)\) 被约掉了。
两边同时平方,得:
\[F^2(x)-2F(x)f(x)+f^2(x)\equiv0\color{red}{\pmod{x^{n}}} \]注意标红的部分,想一想,为什么?
考虑卷积的定义,在 \([0,n/2]\) 的部分每一项都是 \(0\),而右边 \([n/2,n]\) 的部分由 \(a_i=\sum^i_{j=1} a_ja_{i-j}\) 得来,容易发现,\(j\) 与 \(i-j\) 中必有一个小于零。
两边同时乘以 \(G(x)\) 并移项,得:
\[F(x)\equiv 2f(x)-G(x)f^2(x)\pmod{x^n} \]直接使用 FFT 计算即可,时间复杂度 \(O(n\log n)\)。
多项式开根
给定 \(g(x)\),求 \(f(x)\),使得:\(f^2(x)\equiv g(x)\pmod{x^n}\)。
类似求逆,假设已经求出了 \(f_0^2(x)\equiv g(x)\pmod{x^{n/2}}\)。
移项得:
\[f_0^2(x)-g(x)\equiv 0\pmod{x^{n/2}} \]两边平方得:
\[(f_0^2(x)-g(x))^2\equiv 0\pmod {x^n} \]使用简单的数学公式 \((x-y)^2=(x+y)^2-4xy\) 得:
\[(f_0^2(x)+g(x))^2\equiv 4f_0^2(x)g(x)\pmod {x^n} \]将 \(4f_0^2(x)\) 移到左边,得:
\[(\frac{f_0^2(x)+g(x)}{2f_0(x)})^2\equiv g(x) \pmod{x^n} \]左右两边同时开根,得:
\[\frac{f_0^2(x)+g(x)}{2f_0(x)} \equiv f(x)\pmod{x^n} \]\[2^{-1}f_0(x)+2^{-1}f_0^{-1}(x)g(x)\equiv f(x)\pmod{x^n} \]直接计算即可,时间复杂度 \(O(n\log n)\)。
多项式带余除法
问题:给定 \(n\) 次多项式 \(F(x)\) 和 \(m\) 次多项式 \(G(x)\),求 \(Q(x)\) 与 \(R(x)\),满足:
- \(\deg Q(x)=n-m\),\(\deg R(x)
; - \(F(x)=Q(x)\times G(x)+R(x)\)。
对于一个多项式 \(f(x)\),定义 \(f_R(x)=x^nf(\frac 1 x)\),等价于将多项式 \(f\) 的系数翻转一次。
那么:
\[F(x)=Q(x)\times G(x)+R(x)\\F(\frac 1 x)=Q(\frac 1 x)\times G(\frac 1 x)+R(\frac 1 x)\\x^nF(\frac 1 x)=x^{n-m} Q(\frac 1 x)\times x^mG(\frac 1 x)+x^{n-m+1}x^{m-1}R(\frac 1 x)\\F_R(x)=Q_R(x)\times G_R(x)+x^{n-m+1}R_R(x)\\F_R(x)\equiv Q_R(x)\times G_R(x)\pmod {x^{n-m+1}}\\Q_R(x)\equiv F_R(x)\times G_R^{-1}(x) \pmod{x^{n-m+1}} \]直接求出 \(Q(x)\),那么回代可得 \(R(x)\),时间复杂度 \(O(n\log n)\)。
多项式对数函数
求 \(B(x)\equiv \ln A(x)\pmod{x^n}\)。
设 \(F(x)=\ln x\),那么:
\[B(x)=F(A(x))\\B'(x)=F'(A(x))A'(x)\\B'(x)=\frac {A'(x)}{A(x)}\\B(x)=\int \frac {A'(x)dx}{A(x)} \]直接模拟即可,时间复杂度 \(O(n\log n)\)。
多项式指数函数
求 \(F(x)\equiv e^{A(x) \pmod{x^n}}\)。
考虑多项式牛顿迭代。
多项式牛顿迭代解决的问题是,给定函数 \(G(x)\),求 \(G(F(x))\equiv 0\pmod {x^n}\)。
既然是迭代,那么当 \(n=1\) 的时候是要单独求出的,这里不讲。
类似求逆,假设已经求出 \(G(F_0(x))\equiv 0 \pmod{x^{\lceil \frac n x\rceil}}\),考虑如何扩展。
把 \(G(F(x))\) 在 \(F_0(x)\) 处泰勒展开:
\[G(F(x))=G(F_0(x))+\frac {G'(F_0(x))}{1!}(F(x)-F_0(x))+\frac {G''(F_0(x))}{2!}(F(x)-F_0(x))+\ldots \]因为 \(F(x)\) 和 \(F_0(x)\) 的最后几项相同,所以 \((F(x)-F_0(x))^2\) 的最低非 \(0\) 次数大于 \(n\),所以可以得到:
\[G(F(x))\equiv G(F_0(x))+G'(F_0(x))(F(x)-F_0(x))\pmod{x^n}\\G(F_0(x))+G'(F_0(x))(F(x)-F_0(x))\equiv 0\pmod{x^n}\\G(F_0(x))+G'(F_0(x))F(x)-G'(F_0(x))F_0(x)\equiv 0\pmod{x^n}\\G'(F_0(x))F(x)\equiv G'(F_0(x))F_0(x)-G(F_0(x))\pmod {x^n}\\F(x)\equiv F_0(x)-\frac {G(F_0(x))}{G'(F_0(x))} \pmod {x^n} \]那么,将 \(F(x)=e^{A(x)}\) 变形得 \(\ln F(x)-A(x)=0\)。
设 \(G(F(x))=\ln F(x)-A(x)\),那么 \(G'(F(x))=\frac 1 {F(x)}\),带入牛顿迭代公式得:
\[F(x)\equiv F_0(x)-\frac{\ln F_0(x)-A(x)}{F_0(x)}\pmod{x^n}\\F(x)\equiv F_0(x)(1-\ln F_0(x)+A(x))\pmod{x^n} \]直接做,时间复杂度 \(O(n\log n)\)。
参考代码
不想太长,所以放一份在剪贴板。