多项式算法

Published: by Creative Commons Licence

  • Tags:

快速数论变换

了解快速傅里叶变换的人应该知道,快速傅里叶变换之所以能成立源于对任意正整数n,存在这样的一个单位根$w(n)$,满足四个条件:

  1. $w(n)^0,w(n)^1, \ldots, w(n)^{n-1}$各不相同
  2. $w(n)^n = 1$
  3. $w(dn)^{dk}=w(n)^k$。
  4. $w(2n)^n=-1$

在快速傅里叶变换中,我们选取$w(n)=e^{(2\pi /n)i}$作为单位根。

快速傅里叶变换会存在一定的精度损失问题,我们可以使用快速数论变换来替代。但是快速数论变换实际上还是沿用了快速傅里叶变换的整体框架,只是替换了单位根。我们需要找一个素数p,并找到它的一个原根g。可以发现当我们取$w(n)=g^{\frac{p-1}{n}}$时,可以满足上面所有条件。不是所有p都适合用于快速数论变换,我们需要找一个比较特殊的p,使得$n|p-1$。由于n一般是2的幂次,因此我们希望$p-1$包含尽可能多的因子2,下面介绍两个常用的选择:

\[998244353=2^{23}\cdot17 \cdot 7 +1\\ 1004535809=2^{21}\cdot 479+1\]

剩下的和快速傅里叶变换没有什么区别,直接套用就好了。

多项式求逆

给出多项式A,找到A的逆多项式B,满足$AB=1\pmod{x^n}$。

当n为1时,问题非常简单,此时A的逆多项式一定是$1/A[0]$。假设我们已经找到了某个多项式C,使得$AC=1\pmod{x^{n/2}}$。那么,可以推导出:

\[AC=AB=1\pmod{x^{n/2}}\\ \Rightarrow A(B-C)=0\pmod{x^{n/2}}\\ \Rightarrow (B-C)=0\pmod{x^{n/2}}\\ \Rightarrow (B-C)^2=0 \pmod{x^n}\]

最后一步需要说明下,对于$B-C$的后n项都是0,因此对于$B-C$的$x^i$的系数,若$i<n/2$,那么系数一定是$B-C$的后$n/2$项中两项的乘积的累加,一定是0,若$i\geq n/2$,那么一定是后$n/2$项和前$n/2$项的系数乘积的累和,也是0。

将二次方程展开:

\[(B-C)^2=B^2+C^2-2BC=0\pmod{x^n}\\ \Rightarrow B^2=2BC-C^2\pmod{x^n}\]

之后两端都乘上A,借助$AB=1$可知:

\[AB^2=2ABC-AC^2\pmod{x^n}\\ \Rightarrow B=2C-AC^2\pmod{x^n}\]

这样我们就得到了多项式B,并且可以发现,整个流程我们仅假设$A[0]$可逆,即一个多项式可逆当且仅当其$x^0$项的系数可逆。

多项式求逆的时间复杂度为$T(n)=T(\frac{n}{2})+\alpha n\log_2n=O(2\alpha n\log_2n)=O(n\log_2n)$,其中$\alpha$是固定常数。

多项式除法取模

对于阶数为n的多项式$f(x)$,记$f^R(x)=x^nf(\frac{1}{x})$,换句话说$f^R(x)$是$f(x)$的系数前后交换后得到的。

接下来考虑阶数为$n$的多项式$A$和阶数为$m$的多项式$B$,以及阶数小于$m$的多项式$R$,我们希望得到如下形式等式$A(x)=B(x)C(x)+R(x)$。推一下公式:

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

最后的结果可以通过对多项式求逆得到。

而取模可以通过$R(x)=A(x)-B(x)C(x)$得到。

多项式除法和取模都仅涉及到一次的多项式求逆和乘法,因此时间复杂度为$O(n\log_2n)$。

多项式多点求值

给出$n$阶多项式$f$,以及$m$个不同的数$x_1,\ldots, x_m$,要求计算$f(x_1),\ldots,f(x_m)$。

很显然可以设计一套暴力算法,$O(nm)$时间复杂度解决。但是这个问题存在丧心病狂的$O(n\log_2n+m(\log_2m)^2)$的解法。

我们可以$m$个数值分成相同大小的两部分。第一部分是$X_1={x_1,\ldots,x_{\frac{m}{2}}}$,第二部分是$X_2={x_{\frac{m}{2}+1},\ldots, x_m}$。

之后我们定义$g(x)=\prod_{i=1}^{\frac{m}{2}} (x-x_i)$,很显然对于任意$x'\in X_1$,满足$g(x')=0$。我们利用多项式除法可以得到:$f(x)=g(x)A(x)+R(x)$,其中$R$的阶数小于$\frac{m}{2}$。利用$x'\in X_1 \Rightarrow f(x')=R(x')$,因此我们现在只需要计算$X_1$在某个不超过$m$阶的多项式下的取值。重复这个过程,直到点集不可再分。

分治的时间复杂度为:

\[T(m)=2T(\frac{m}{2})+\alpha m\log_2m = \alpha m(\log_2m)^2=O(m(\log_2m)^2)。\]

总的时间复杂度为$O(n\log_2n+m(\log_2m)^2)$。

提供一道题目

等差卷积

考虑计算下面多项式:

\[H(x)=\sum_{i=0}^nx^i\sum_{j=i}^nf(j)g(j-i)\]

它非常类似于一般卷积,只是内部是等差,而一般卷积是等和。我们定义$f'(x)=f(n-x)$。那么代入得到

\[H'(x)=\sum_{i=0}^nx^i\sum_{j=0}^if'(j)g(i-j)\]

很显然$H'$可以利用一般的快速卷积在$O(n\log_2n)$时间复杂度内得到,下面我们来看$H'$和$H$的关系。

\(H'_i=\sum_{j=0}^if'(j)g(i-j)=\sum_{j=0}^if(n-j)g(i-j)\) 这里包含了所有差为$n-i$的项的和。因此可以得到:

\[H_i'=H_{n-i}\]

即最终我们可以通过翻转$H'$得到$H$。

注意这里我们要在$f$前面填充0之前翻转,且$H'$也仅翻转后面的$n$项目。

求序列等幂和

假设提供一个序列${a_i}$,记$f(k)=\sum_{i=1}^na_i^k$,要求我们计算出$f(1),f(2),\ldots, f(m)$。

这个问题可以用生成函数$G(x)$来求解。

\[G(x)=\sum_{i=0}\sum_{j=1}^na_{j}^ix^i\\ =\sum_{j=1}^n\sum_{i=0}(a_{j}x)^i\\ =\sum_{j=1}^n\frac{1}{1-a_jx}\\ =\frac{\sum_{j=1}^n\prod_{j\neq i}(1-a_jx)}{\prod_{j}(1-a_jx)}\]

假如我们记$D(x)=\prod_{j}(1-a_jx)$,记$U(x)=\sum_{j=1}^n\prod_{j\neq i}(1-a_jx)$,那么我们可以发现$D_i$和$U_i$存在某种联系。事实上,考虑所有对$x^j$有贡献的系数一定是$j$个${-a_i}$序列中的数的乘积,而这个乘积会在$U_i$中出现贡献$n-j$次。因此我们可以直接断定$U_i=(n-i)D_i$。

于是乎我们得到了$U(x)$和$D(x)$后,就可以用多项式求逆的技术得出$G(x)=U(x)D^{-1}(x)$。

这里我们计算$D(x)$的时候,可以利用分治+快速卷积的技术将时间复杂度优化到$O(n(\log_2n)^2 )$。

多项式差与多项式前缀和

定义1:对于一个任意多项式$f(x)$,我们记其前缀和函数为$g(x)=\sum_{i=0}^xf(i)$。

命题1:对于任意$n$阶多项式$A(x)$,我们始终可以找到另外一个$n+1$阶多项式$B(x)$,满足$B(x)-B(x-1)=A(x)$

证明:

我们可以利用归纳法进行证明,当$A(x)$是0阶多项式的时候,即$A(x)=c$,那么我们可以得到$B(x)=cx$,很显然$B(x)-B(x-1)=cx-c(x-1)=c=A(x)$。

假设对于任意不大于$k-1$阶的多项式$A(x)$,命题都成立,那么当$A(x)$为$k$阶多项式的时候,记$A(x)=a_kx^k+\ldots+a_0$,注意到$x^k-(x-1)^k$是$k-1$阶多项式,因此我们可以通过为$B(x)$选择合适的最高项系数$b_{k+1}=\frac{a_k}{k}$,将问题化简为$(B(x)-b_{k+1}x^{k+1})-(B(x-1)-b_{k+1}(x-1)^{k+1})=A(x)-(b_{k+1}x^{k+1}-b_{k+1}(x-1)^{k+1})$,我们重新标记后得到公式:$B'(x)-B'(x-1)=A'(x)$,其中我们要找的$B'(x)$是$k$阶多项式,而$A'(x)$是$k-1$阶多项式,而之前的归纳法已经证明了$B'(x)$的存在性,因此$B(x)=B'(x)+b_{k+1}x^{k+1}$也是存在的。

命题2:如果多项式$f(x)$是$n$阶多项式,那么其前缀和函数$g(x)$一定是$n+1$阶多项式

证明:

首先注意到$g(x)-g(x-1)=f(x)$,按照命题1我们一定可以刻画出这样一个$n+1$阶多项式$g(x)$,现在只需要保证$g(0)=f(0)$即可。注意到修改$g(x)$的常数项不会影响$g(x)-g(x-1)=f(x)$等式的成立,因此我们可以通过修改$g(x)$的常数项为$f(0)$即可找到与前缀和函数完全相同的一个$n+1$阶多项式。

Chirp-Z Transform

考虑这样一个问题:给定一个$n$阶多项式$f(x)=\sum_{i=0}^na_ix^i$,以及某个数$c$和$m$,要求计算$f(c^0),f(c^1),\ldots,f(c^m)$。其中$1\leq n, m\leq 3\times 10^5$。

暴力做法$n^2$显然不可行。

可以发现

\[f(c^j)=\sum_{i=0}^na_ic^{ij}\]

由于$ij={i+j\choose 2}-{i\choose 2}-{j\choose 2}$,因此有

\[\begin{aligned} f(c^j)&=\sum_{i=0}^na_ic^{ij}\\ &=\sum_{i=0}^na_ic^{ {i+j\choose 2}-{i\choose 2}-{j\choose 2} }\\ \end{aligned}\]

稍加转换可以得到:

\[\begin{aligned} &f(c^j)c^{j\choose 2}=\sum_{i=0}^na_ic^{ {i+j\choose 2}-{i\choose 2} S}\\ \Rightarrow &f(c^j)c^{j\choose 2}=\sum_{i=0}^n[c^{i+j\choose 2}][a_ic^{-{i\choose 2}}] \end{aligned}\]

可以发现右边是两个多项式的等差卷积,因此可以直接用快速卷积算法$O((n+m)\log_2(n+m))$快速求解。

提供一道题目,卡常题,慎做。

多元多项式卷积

非常神奇的技术,从这个博客学到的。

考虑给定两个$k$元多项式$f(x_1,\ldots,x_k)$和$g(x_1,\ldots,x_k)$,求它们的卷积对$(x_1^{n_1},\ldots,x_k^{n_k})$取模的结果(注意这里取模并不是不进位的意思,而是舍弃发生进位的项的意思)。记$N=\prod_{i=1}^kn_i$

我们可以用一个$0$到$N-1$的下标来表示多项式中的一项。

考虑如何移除进位项带来的影响,我们可以引入一个新的哈希函数$\chi(i)=\sum_{j=1}^{k-1} \lfloor \frac{i}{\prod_{t\leq j} n_t} \rfloor$。考虑到$0\leq \lfloor \frac{a+b}{\prod_{t\leq j} n_t}\rfloor -\lfloor \frac{a}{\prod_{t\leq j} n_t} \rfloor - \lfloor \frac{b}{\prod_{t\leq j} n_t} \rfloor\leq 1$,因此有$0\leq \chi(a+b)-\chi(a)-\chi(b)<k$,其中如果$\chi(a+b)-\chi(a)-\chi(b)=0$,则表示$a+b$过程中不发生进位。

之后我们将$f$重新写作如下的形式$f(x,z)=\sum_{i=0}^{N-1}f_ix^iz^{\chi(i)}$,同理$g(x,z)=\sum_{i=0}^{N-1}g_ix^iz^{\chi(i)}$。我们可以从它们中还原出原先的$k$元多项式。

考虑$f$和$g$中两项的乘积:

\[f_ix^iz^{\chi(i)}\cdot g_jx^jz^{\chi(j)}=f_ig_jx^{i+j}z^{\chi(i)+\chi(j)}[\chi(i+j)=\chi(i)+\chi(j)]\]

其中它们的乘积会被统计到我们结果中,当且仅当$\chi(i)+\chi(j)=\chi(i+j)$。由于$0\leq \chi(a+b)-\chi(a)-\chi(b)<k$以及$\chi(i)$可能比较大,因此我们可以让结果对$z^{k}-1$取模,可以发现信息是不会丢失的。

接下来我们只需要计算$f*g$,并舍弃那些$x$项指数$i$与$t$项指数$j$不匹配的项,即$\chi(i)\neq j$,这些项是发生了进位的。

最后我们要考虑如何快速计算两个二元多项式的卷积。事实上我们可以按$t$的指数对两边进行展开,令$F_i=\sum_{i=0}^{N-1}f_ix^i[\chi(i)=t]$,$G_i=\sum_{i=0}^{N-1}g_ix^i[\chi(i)=t]$:

\[\begin{aligned} f*g=&(\sum_{i=0}^{N-1}f_ix^iz^{\chi(i)}) * (\sum_{i=0}^{N-1}g_ix^iz^{\chi(i)})\\ =&(\sum_{t=0}^{k-1}z^t\sum_{i=0}^{N-1}f_ix^i[\chi(i)=t]) * (\sum_{t=0}^{k-1}z^t\sum_{i=0}^{N-1}g_ix^i[\chi(i)=t])\\ =&(\sum_{t=0}^{k-1}z^tF_t)*(\sum_{t=0}^{k-1}z^tG_t)\\ =&\sum_{t=0}^{k-1}\sum_{p=0}^{k-1}z^{t+p}F_t*G_p \end{aligned}\]

现在问题就变成了$k^2$个多项式的乘积了,一种方案就是直接暴力FFT,时间复杂度为$O(k^2N\log_2N)$,当然我们可以先对每个$F_i$和$G_i$预处理它们的DFT形式,那么公式中卷积就可以简单的用点积来替代,之后把同类项合并,在一同取逆即可,时间复杂度可以减少到$O(kN\log_2N+k^2N)$,考虑到$k\leq \log_2 N$,因此时间复杂度可以简写为$O(kN\log_2N)$。

提供几道题目:

话说子集卷积其实也是多元多项式的卷积,对$(x_1^2,\ldots,x_k^2)$进行取模而已。所以同样可以用这个算法解子集卷积,不过由于FFT操作会慢于大部分FWT,因此最好还是用FWT做子集卷积。

参考资料