「快速沃尔什变换」学习笔记

发表于 | 分类于 |

想起来这玩意看懂了之后一直没写……
要不是 WJJ 今天问了我我还真就遗忘在任务列表里了（逃

思想

考虑类比 FFT 处理多项式乘法的方式，可以尝试设计一种变换，使得其能够将位运算卷积转化为对应位置乘积。

处理位运算卷积

按位或

对于序列 \(a,b\)，考虑序列 \(c\)，其中 \[ \newcommand\fwt[1]{\mathrm{FWT}_{\mathrm{ #1 }}} \newcommand\merge{\mathrm{merge}} \newcommand\xor{\operatorname{xor}} \newcommand\popcnt{\operatorname{popcount}} c_i = \sum\limits_{j|k=i} a_j b_k \]

考虑一个对于序列的变换 \(\fwt{or}(a)\)： \[ \fwt{or}(a)_i = \sum\limits_{j|i} a_j \]

则注意到 \[ \begin{aligned} \fwt{or}(c)_i &= \sum\limits_{(j|k)|i=i} a_j b_k \\ &= \sum\limits_{j|i}\sum\limits_{k|i}a_jb_k \\ &= \fwt{or}(a)_i \fwt{or}(b)_i \end{aligned} \]

这显然正是我们设计这个变换而所期望得到的。

接下来考虑如何快速实现这个变换（\(O(n^2)\) 肯定是接受不了的）。
考虑分治，每次将 \(a\) 数组按下标在二进制下的最高位分，即对半分。
设在 \(a\) 中下标最高位为 \(0\) 的部分为 \(a_0\)，为 \(1\) 的部分为 \(a_1\)。
则由于只有 \(0\) 对 \(1\) 产生贡献，而 \(0,1\) 均对 \(1\) 产生贡献，有 \[ \fwt{or}(a) = \merge(\fwt{or}(a_0),\fwt{or}(a_0)+\fwt{or}(a_1)) \]

其中两个序列的加法意义为对应位置相加，\(\merge\) 为将两个序列首尾相接形成一个新的序列。
如此即可分治。

按位与

对于序列 \(a,b\)，考虑序列 \(c\)，其中 \[ c_i = \sum\limits_{j\&k=i} a_j b_k \]

类似地，考虑 \[ \fwt{and}(a)_i = \sum\limits_{j\&i=i} a_j \]

显然 \(\fwt{and}(c)_i = \fwt{and}(a)_i \fwt{and}(b)_i\)。

类似地，\(0,1\) 都对 \(0\) 有贡献但只有 \(1\) 对 \(1\) 有贡献，故有 \[ \fwt{and}(a)_i = \merge(\fwt{and}(a_0)+\fwt{and}(a_1),\fwt{and}(a_1)) \]

按位异或

对于序列 \(a,b\)，考虑序列 \(c\)，其中 \[ c_i = \sum\limits_{j \xor k=i} a_j b_k \]

定义 \(a \oplus b = \popcnt(a\mathop\&b) \bmod 2\)，其中 \(\popcnt(x)\) 表示 \(x\) 在二进制表示下 \(1\) 的个数。
则显然有 \((i \oplus j) \xor (i \oplus k) = (j \xor k) \oplus i\)。

考虑 \[ \fwt{xor}(a)_i = \sum\limits_{i\oplus j=0}a_j - \sum\limits_{i\oplus j=1}a_j \]

则有 \[ \begin{aligned} \fwt{xor}(a)_i \fwt{xor}(b)_i &= \left(\sum\limits_{i\oplus j=0}a_j - \sum\limits_{i\oplus j=1} a_j\right)\left(\sum\limits_{i\oplus k=0}b_k - \sum\limits_{i\oplus k=1}b_k\right) \\ &= \left(\sum\limits_{i\oplus j=0}\sum\limits_{i\oplus k=0}a_jb_k + \sum\limits_{i\oplus j=1}\sum\limits_{i\oplus k=1}a_jb_k\right) - \left(\sum\limits_{i\oplus j=0}\sum\limits_{i\oplus k=1} a_jb_k + \sum\limits_{i\oplus j=1}\sum\limits_{i\oplus k=0}a_jb_k\right) \\ &= \sum\limits_{(i\oplus j)\xor(i\oplus k)=0}a_jb_k - \sum\limits_{(i\oplus j)\xor(i\oplus k)=1}a_jb_k \\ &= \sum\limits_{(j\xor k)\oplus i=0} a_jb_k - \sum\limits_{(j\xor k)\oplus i=1} a_jb_k \\ &= \fwt{xor}(c)_i \end{aligned} \]

分治合并时考虑最高位会贡献 \(1\) 还是 \(-1\)，有 \[ \fwt{xor}(a) = \merge(\fwt{xor}(a_0)+\fwt{xor}(a_1),\fwt{xor}(a_0)-\fwt{xor}(a_1)) \]

逆变换

逆变换什么的，把 FWT 倒过来做就好啦（
反正都是自底向上的，还是比较好做的。

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
#include <cstdio>
#define add(x,y) (x + y >= mod ? x + y - mod : x + y)
#define dec(x,y) (x < y ? x - y + mod : x - y)
using namespace std;
const int N = 1 << 17;
const int mod = 998244353;
const int inv = 499122177;
int n,a[N + 5],b[N + 5];
int f[N + 5],g[N + 5];
inline void fwt_or(int *a,int type)
{
    for(register int w = 2,m = 1;w <= n;w <<= 1,m <<= 1)
        for(register int i = 0;i < n;i += w)
            for(register int j = 0;j < m;++j)
                a[i | j | m] = type == 1 ? add(a[i | j | m],a[i | j]) : dec(a[i | j | m],a[i | j]);
}
inline void fwt_and(int *a,int type)
{
    for(register int w = 2,m = 1;w <= n;w <<= 1,m <<= 1)
        for(register int i = 0;i < n;i += w)
            for(register int j = 0;j < m;++j)
                a[i | j] = type == 1 ? add(a[i | j],a[i | j | m]) : dec(a[i | j],a[i | j | m]);
}
inline void fwt_xor(int *a,int type)
{
    for(register int w = 2,m = 1;w <= n;w <<= 1,m <<= 1)
        for(register int i = 0;i < n;i += w)
            for(register int j = 0;j < m;++j)
            {
                int t = a[i | j | m];
                a[i | j | m] = dec(a[i | j],t),
                a[i | j] = add(a[i | j],t),
                type == -1 && (a[i | j] = (long long)a[i | j] * inv % mod,a[i | j | m] = (long long)a[i | j | m] * inv % mod);
            }
}
int main()
{
    scanf("%d",&n),n = 1 << n;
    for(register int i = 0;i < n;++i)
        scanf("%d",a + i);
    for(register int i = 0;i < n;++i)
        scanf("%d",b + i);
    for(register int i = 0;i < n;++i)
        f[i] = a[i],g[i] = b[i];
    fwt_or(f,1),fwt_or(g,1);
    for(register int i = 0;i < n;++i)
        f[i] = (long long)f[i] * g[i] % mod;
    fwt_or(f,-1);
    for(register int i = 0;i < n;++i)
        printf("%d%c",f[i]," \n"[i == n - 1]);
    for(register int i = 0;i < n;++i)
        f[i] = a[i],g[i] = b[i];
    fwt_and(f,1),fwt_and(g,1);
    for(register int i = 0;i < n;++i)
        f[i] = (long long)f[i] * g[i] % mod;
    fwt_and(f,-1);
    for(register int i = 0;i < n;++i)
        printf("%d%c",f[i]," \n"[i == n - 1]);
    for(register int i = 0;i < n;++i)
        f[i] = a[i],g[i] = b[i];
    fwt_xor(f,1),fwt_xor(g,1);
    for(register int i = 0;i < n;++i)
        f[i] = (long long)f[i] * g[i] % mod;
    fwt_xor(f,-1);
    for(register int i = 0;i < n;++i)
        printf("%d%c",f[i]," \n"[i == n - 1]);
}