快速沃尔什变换

Published: by Creative Commons Licence

  • Tags:

引言

快速沃尔什变换(Fast Walsh Transform)适用于计算一类多项式卷积的高效算法。定义

\[C=A\otimes B=\sum_{k=0}^nX^k\sum_{i\oplus j=k}A_iB_j\]

其中$A$,$B$,$C$均为$n$阶多项式。而$A_i$,$B_j$分别表示多项式$A$中$X^i$项的系数和多项式$B$中$X^j$项的系数。

看着眼熟吧,这不是和多项式的乘法非常像吗,当$i\oplus j=i+j$时,上面就是我们熟知的多项式乘法公式。我们知道解决多项式乘法不仅有$O(n^2)$的朴素算法,还有$O(n\log_2n)$的快速傅里叶变换。快速傅里叶变换实际上是将两个多项式$A$,$B$转换为点值表示法DFT(A)和DFT(B),之后利用向量乘法得到DFT(C)=DFT(A)DFT(B),最后从DFT(C)中还原多项式$C$。我们不免想到,是否能利用同样的技术进行加速呢?

或运算

如果$i\oplus j=i\mid j$,那么我们可以定义:

\[FWT(A)=\sum_{i=0}^{n}X^i\sum_{j\in i}A_j\]

其中$j\in i$表示j中$1$出现的二进制位置是$i$的子集,换言之$i=i\mid j$。我们注意到

\[FWT(C)_k=\sum_{t\in k}C_t=\sum_{t\in k}\sum_{i\mid j=t}A_iB_j=\sum_{i\mid j\in k}A_iB_j\\ =\sum_{i\in k}A_i\sum_{j\in k}B_j=(\sum_{i\in k}A_i)(\sum_{i\in k}B_i)=FWT(A)_k\cdot FWT(B)_k\]

因此我们发现了$FWT(C)=FWT(A)\cdot FWT(B)$。下面我们讨论如何快速计算$FWT(A)$。记$m=\lceil \frac{n}{2} \rceil$,$L(A)=\sum_{i=0}^{m-1}A_iX^i$,$R(A)=\sum_{i=m}^nA_iX^i$。很显然$A=L(A)+R(A)\cdot X^m$。

\[FWT(A)=\sum_{i=0}^{n}X^i\sum_{j\in i}A_j=\sum_{i=0}^{n}X^i(\sum_{j\in i,j<m}A_j+\sum_{j\in i,j\geq m}A_j)\\ =\sum_{i=0}^{m-1}X^i(\sum_{j\in i,j<m}A_j+\sum_{j\in i,j\geq m}A_j)+\sum_{i=m}^{n}X^i(\sum_{j\in i,j<m}A_j+\sum_{j\in i,j\geq m}A_j)\\ =\sum_{i=0}^{m-1}X^i\sum_{j\in i}L(A)_j+X^m\sum_{i=0}^{m-1}X^i(\sum_{j\in i}L(A)_j+\sum_{j+m\in i+m}R(A)_j)\\ =FWT(L(A))+X^m(FWT(L(A))+FWT(R(A)))\]

我们可以简单地记作

\[FWT(A)=(FWT(L(A)), FWT(L(A))+FWT(R(A)))\]

这样我们可以每次都将$A$拆分为前部和后部进行分治计算。时间复杂度为$O(n\log_2n)$。而$FWT$的逆也很容易得到:

\[FWT^{-1}(A)=FWT^{-1}(L(A))+X^mFWT^{-1}(R(A)-L(A))\\ =(FWT^{-1}(L(A)),FWT^{-1}(R(A)-L(A)))\]

同样可以采用分治算法解决。

且运算

如果$\oplus$代表且运算,那么$FWT$函数定义如下:

\[FWT(A)=\sum_{i=0}^nX^i\sum_{i\in j}A_j\]

其余推理类似,不再赘述,直接给出公式:

\[FWT(C)=FWT(A)\cdot FWT(B)\\ FWT(A)=(FWT(L(A))+FWT(R(A)),FWT(R(A)))\\ FWT^{-1}(A)=(FWT^{-1}(L(A)-R(A)),FWT^{-1}(R(A)))\]

异或运算

如果$\oplus$代表亦或运算。记$E$为所有二进制形式中1的个数为偶数个的正整数集合,而$O$表示所有二进制形式中1的个数为奇数个的正整数集合。那么$FWT$函数定义如下:

\[FWT(A)=\sum_{i=0}^nX^i(\sum_{i\&j\in E}A_j-\sum_{i\&j\in O}A_j)\]

由于(可以通过枚举证明)

\[(i\&j)\oplus k=(i\oplus k)\&(j\oplus k)\]

因此:

\[FWT(C)_k=\sum_{k\&t\in E}C_t-\sum_{k\&t\in O}C_t\\ =\sum_{k\&t\in E}\sum_{i\oplus j=t}A_iB_j-\sum_{k\&t\in O}\sum_{i\oplus j=t}A_iB_j\\ =\sum_{k\&(i\oplus j)\in E}A_iB_j-\sum_{k\&(i\oplus j)\in O}A_iB_j\\ =\sum_{(k\&i)\oplus (k\& j)\in E}A_iB_j-\sum_{(k\&i)\oplus (k\& j)\in O}A_iB_j\\\]

记$b(i)$表示$i$的二进制形式中$1$的个数。记$r$表示$i\oplus j$中可能的$1$的数量,初始值为$b(i)+b(j)$。遍历每一位,若浏览某一位时,$i$与$j$均为$0$,则$r$不变,若一者为$1$,$r$依旧不变,否则二者都为$1$,则$r$减少2。很显然这个算法最终得到的就是$b(i\&j)$。因此我们证明了$b(i\&j)$与$b(i)+b(j)$的二进制形式中1出现的次数拥有相同的奇偶性。

\[=\sum_{k\&i \in E}A_i\sum_{k\&j\in E}B_j+\sum_{k\&i \in O}A_i\sum_{k\&j\in O}B_j-\sum_{k\&i \in E}A_i\sum_{k\&j\in O}B_j-\sum_{k\&i \in O}A_i\sum_{k\&j\in E}B_j\\ =(\sum_{k\&i \in E}A_i-\sum_{k\&i \in O}A_i)(\sum_{k\&i \in E}B_i-\sum_{k\&i \in O}B_i)\\ =FWT(A)_kFWT(B)_k\]

因此我们可以推出:

\[FWT(C)=FWT(A)\cdot FWT(B)\]

接着我们考虑如何拆分多项式。

\[FWT(A)=\sum_{i=0}^nX^i(\sum_{i\&j\in E}A_j-\sum_{i\&j\in O}A_j)\\ =\sum_{i=0}^nX^i(\sum_{i\&j\in E,j<m}A_j+\sum_{i\&j\in E,j\geq m}A_j-\sum_{i\&j\in O,j<m}A_j-\sum_{i\&j\in O,j\geq m}A_j)\\ =\sum_{i=0}^{m-1}X^i(\sum_{i\&j\in E,j<m}A_j+\sum_{i\&j\in E,j\geq m}A_j-\sum_{i\&j\in O,j<m}A_j-\sum_{i\&j\in O,j\geq m}A_j)\\ +\sum_{i=m}^nX^i(\sum_{i\&j\in E,j<m}A_j+\sum_{i\&j\in E,j\geq m}A_j-\sum_{i\&j\in O,j<m}A_j-\sum_{i\&j\in O,j\geq m}A_j)\\ =\sum_{i=0}^{m-1}X^i(\sum_{i\&j\in E}L(A)_j+\sum_{i\&j\in E}R(A)_j-\sum_{i\&j\in O}L(A)_j-\sum_{i\&j\in O}R(A)_j)\\ +X^m\sum_{i=0}^{m-1}X^i(\sum_{i\&j\in E}L(A)_j+\sum_{i\&j\in O}R(A)_j-\sum_{i\&j\in O}L(A)_j-\sum_{i\&j\in E}R(A)_j)\\ =FWT(L(A))+FWT(R(A))+X^m(FWT(L(A))-FWT(R(A)))\\ =(FWT(L(A))+FWT(R(A)),FWT(L(A))-FWT(R(A)))\]

而逆函数容易得到:

\[FWT^{-1}(A)=(FWT^{-1}(\frac{L(A)+R(A)}{2}),FWT^{-1}(\frac{L(A)-R(A)}{2}))\]

整体的时间复杂度为$O(n\log_2n)$。

题目1:给定长度为$2^m$的向量$x$,要求计算$n$个连续的$x$的异或卷积:$x\oplus x\oplus \ldots \oplus x$。其中$1\leq m\leq 20$,$n\leq 10^18$,结果对素数$p$取模。

考虑到$FWT(x^n)=FWT(x)^n$,因此可以先对$x$求一次FWT,之后将每个元素$x_i$替换为$x_i^n$,最后IFWT恢复。

题目2:给定长度为$2^m$的向量$x$,要求计算$n$个连续的$x$的异或卷积:$x\oplus x\oplus \ldots \oplus x$。其中$1\leq m\leq 20$,$n\leq 10^18$,结果对任意数$p$取模,其中$1\leq p\leq 10^9$。

由于IFWT的时候需要做除2操作,因此这里不能直接使用。有一种叫做K进制FWT的技术,其只需要在做IFWT的最后阶段,将整个向量除上$2^m$即可。但是这里依旧不能保证$gcd(2^m,p)=1$,因此我们需要将$p$先乘上$2^m$做取模,之后每个向量元素能整除$2^m$的性质就不会被影响。最后的阶段才将整个向量对$p$取模。

题目3:给定两个序列$A_0,b$,其中$A_0$的长度为$2^m$,$b$的长度为$m+1$。定义$A_{i,j}=\sum_{k}A_{i-1,k}b[bitcount(k\oplus j)]$。要求输出$A_n$,结果对$p$取模。其中$1\leq m\leq 20$,$1\leq n\leq 10^{18}$

提供一道题目

我们将$b$展开,定义$c[i]=b[bitcount(i)]$,那么我们要求的就是$A_{i,j}=\sum_{k}A_{i-1,k}c[k\oplus j]$。可以发现$A_i=A_{i-1}\oplus c$,因此结果为$A_n=A_{0}\oplus c^n$。

我们可以用题目$2$提到的方法计算出最终结果。

快速子集卷积

考虑这样一个问题,给定序列$a$和$b$,要求找到另外一个序列$c$,满足:

\[c_i=\sum_{j\mid k=i\land j\&k=0}a_jb_k\]

我们可以定义$A_{i,j}=a_j[bitcount(j)=i]$,同理定义$B$和$C$。那么考虑$T_k=\sum_{i=0}^kA_{i}\oplus B_{k-i}$,其与$C_k$的区别非常小,实际上$C_{k,i}=T_{k,i}[bitcount(i)=k]$。

因此如果我们完成上面共$n^2$次卷积,就可以求出$C_1,C_2,\ldots$,而我们要求的结果为$c=\sum_iC_i$。

上面的过程时间复杂度为$O(n^32^n)$。

我们考虑到可以预先处理出$A_1,A_2,\ldots$以及$B_1,B_2,\ldots$的FWT形式,然后两两点积即可,最后再整体做一次IFWT即可(FWT和IFWT都是线性变换)。即$T_{k}=IFWT(\sum_{i=0}^kFWT(A_{i})\oplus FWT(B_{k-i}))$。这样就只需要做$n^2$次的点积以及$O(n)$次的FWT和IFWT操作即可,时间复杂度为$O(n^22^n)$。