2023.4.7【模板】快速沃尔什变换FWT
题目概述
给定长度为 \(2^n\) 两个序列 \(A,B\)，设
\(C_i=\sum_{j\oplus k = i}A_j \times B_k\)
分别当 \(\oplus\) 是 or,and,xor 时求出 \(C\)
我们通常将这个操作，叫做“位运算卷积”，因为它的卷积是按照位运算法则“卷”起来的。
算法流程
或
当运算符是或时，我们类似于FFT设计一个函数，使得原来的\(O(n^2)\)卷积可以使用\(O(n)\)对应点乘来实现。
而在或运算中，这个函数可以简单地设计为
证明过程如下：
题目中的\(C_i\)要求\(C_i = \Sigma_{j|k = i}a_jb_k\)
(\(F_c[i]\)指以c为基的函数，构造与\(F_i\)，\(G_i\)一致)
如何由a快速转换成\(F_a\)呢？
考虑分治，考虑0~\(2^{n - 2} - 1\)和\(2^{n - 2}\)~\(2^{n - 1} - 1\)两个区间，考虑与左边的数k(小于\(2^{n - 2}\))或起来等于k的数字，其与k + \(2^{n - 2}\)或起来一定是k + \(2^{n - 2}\)，所以每层中将左半边的相应位置累加到右半边即可，然而左半边的数或上右半边的数，结果肯定是右半边的数，不会对左半边产生贡献，所以左半边的值仍然等于原来的值。
如何将\(F_a\)快速转换成a呢？
刚刚我们发现，右半边的答案被累加了一个左半边，所以减掉就好了。
这两个过程可以合并：
inline void OR(int *f,int op)
{
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
f[j + mid] = (f[j + mid] + f[j] * op) % MOD;
}
不仅用于变换，我们发现这个过程相当于以“二进制下哪些位为1”为条件对i进行了一次子集求和，即f[i]代表了i及其子集的和，这在很多方面都有极大应用。
与
仿照或运算，我们定义（不一样）：
所以
我们发现，分治的时候，如果右半边的一个数k与上另一个数等于k，那么另一个数与上k - \(2^{n -2}\)（最高位清空，其余不变）也一定等于k - \(2^{n - 2}\)。所以左半边的答案加上右半边，而右半边不变即可。
inline void AND(int *f,int op)
{
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
f[j] = (f[j] + f[j + mid] * op) % MOD;
}
异或
要求做到\(C_i=\sum_{j\oplus k = i}A_j \times B_k\)
我们设
需要做到\(\Sigma_{i = 0}^{n}g(x,i)C_i = \Sigma_{j = 0}^{n}g(x,j)a[j] \ * \ \Sigma_{i = 0}^ng(x,i)b[i]\)
即\(\Sigma_{j = 0}^n \Sigma_{k = 0}^ng(x,j\oplus k)a[j]b[k] = \Sigma_{j = 0}^n \Sigma_{k = 0}^ng(x,j)g(x,k)a[j]b[k]\)
所以参数g需要满足\(g(x,j \oplus k) = g(x,j)g(x,k)\)
有一个结论，即异或前后1的个数的奇偶性不会发生改变（分类讨论推一下）
设\(popcount(x)\)为x二进制表示下1的个数
则\(g(x,i) = (-1)^{popcount(i \cap x)}\)
此处\(i \cap x\)为将i中x相应的位数是1的位，可以说是\(i \& x\)
这个地方选择-1这个数字因为它最好体现“奇偶性”这个概念
所以我们得到
对于新加入的一位，我们进行分类讨论：
1.若k为左半边的点，左半边的答案仍然为原来的答案，右半边的点在最高位上0&1仍然是0，没有改变popcount，所以加上右半边答案。
2.若k为右半边的点，左半边的点在最高位上&后等于0，直接加上，但是右半边的点在最高位1&1=1,所有数popcount都加了1，全部反号，所以减去右半边答案。
综上，在每次向上一层时：
左 = 左 + 右
右 = 左 - 右
在逆过程中，为了得到原来的左、右：
左 = (左 + 右) / 2
右 = (左 - 右) / 2
这两个过程也可以合并
inline void XOR(int *f,int op)
{
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
{
int X = f[j],Y = f[j + mid];
f[j] = (X + Y) % MOD;
f[j + mid] = (X - Y + MOD) % MOD;
f[j] = 1ll * f[j] * op % MOD;
f[j + mid] = 1ll * f[j + mid] * op % MOD;
}
}
如果是逆过程，op就是2关于MOD的逆元，否则就是1
Code
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N = (1 << 17) + 5,MOD = 998244353,inv = 499122177;
int a[N],b[N],c[N],n;
inline void OR(int *f,int op)
{
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
f[j + mid] = (f[j + mid] + f[j] * op) % MOD;
}
inline void AND(int *f,int op)
{
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
f[j] = (f[j] + f[j + mid] * op) % MOD;
}
inline void XOR(int *f,int op)
{
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
{
int X = f[j],Y = f[j + mid];
f[j] = (X + Y) % MOD;
f[j + mid] = (X - Y + MOD) % MOD;
f[j] = 1ll * f[j] * op % MOD;
f[j + mid] = 1ll * f[j + mid] * op % MOD;
}
}
signed main()
{
cin>>n;
n = (1 << n) - 1;
for(int i = 0;i <= n;i++) cin>>a[i];
for(int j = 0;j <= n;j++) cin>>b[j];
memset(c,0,sizeof(c));
OR(a,1);OR(b,1);
for(int i = 0;i <= n;i++) c[i] = 1ll * a[i] * b[i] % MOD;
OR(a,MOD - 1);OR(b,MOD - 1);OR(c,MOD - 1);
for(int i = 0;i <= n;i++) cout<<c[i]<<" ";
cout<<endl;
memset(c,0,sizeof(c));
AND(a,1);AND(b,1);
for(int i = 0;i <= n;i++) c[i] = 1ll * a[i] * b[i] % MOD;
AND(a,MOD - 1);AND(b,MOD - 1);AND(c,MOD - 1);
for(int i = 0;i <= n;i++) cout<<c[i]<<" ";
cout<<endl;
memset(c,0,sizeof(c));
XOR(a,1);XOR(b,1);
for(int i = 0;i <= n;i++) c[i] = 1ll * a[i] * b[i] % MOD;
XOR(a,inv);XOR(b,inv);XOR(c,inv);
for(int i = 0;i <= n;i++) cout<<c[i]<<" ";
cout<<endl;
return 0;
}
一道好的练习题
洛谷P3175 HAOI2015 按位或 P3175 [HAOI2015]按位或 - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)
题目描述
刚开始你有一个数字 \(0\)，每一秒钟你会随机选择一个 \([0,2^n-1]\) 的数字，与你手上的数字进行或（C++,C 的 |，pascal 的 or）操作。选择数字 \(i\) 的概率是 \(p_i\)。保证 \(0\leq p_i \leq 1\)，\(\sum p_i=1\) 。问期望多少秒后，你手上的数字变成 \(2^n-1\)。
算法流程
这道题需要用到一个结论：min-max容斥原理
对于一个集合S，有：
"下面我们尝试着给出证明，这里只证明第一个等式好了，后边的可以自行推出。
其实只需要证明一件事，就是除了min(T)=max(S)的那个值，其他的min值都被消掉了就可以了(这里说明一下，我们假定集合中的元素两两相异，如果存在相同的值的话，我们给其中几个加上一些eps扰动一下即可，反正不影响最值就是了)。
先来说明max(S)的系数为什么是1，假设中S最大的元素是a，那么我们会发现只有min(a)=max(S)所以max(S)的系数必须是1。
然后再说明为什么别的min都被消掉了，假设某个元素b的排名是k，那么min(T)=b当且仅当我们选出的集合是后n-k个的元素构成的集合的子集然后并上{b}得到的，我们会发现显然这样的集合有2n−k种，而显然这其中恰有2n−k−1中是有奇数个元素的，恰有2n−k−1种是有偶数个元素的，两两相消自然就成0了，当然上述等式在k=n的时候不成立，但是此时剩下的刚好是最大值。"
(选自洛谷上shadowice1984的题解)
更好的一点是，这个等式在期望意义下仍然成立。
设\(E(max(S))\)为最晚的一位（二进制）出现时间的期望值。
套用公式得到：\(max(S) = \sum_{T \subseteq S}(-1)^{|T| + 1}min(T)\)，考虑min(T)为只要T中任何一位出现就好了，所以
那么和T有交集的集合如何算呢？
有交集的不好算，我们就算没有交集的，没有交集的集合一定属于T关于全集的补集，即\(T \oplus (2^n - 1)\)，所以
至于\(T \oplus (2^n - 1)\)的子集和，用FWT的or类型计算即可
注意每次加答案是要判断\({1 - \sum_{G \cap T = \emptyset}p[G]}\)是否为0
#include<bits/stdc++.h>
using namespace std;
const int N = 1 << 20;
double p[N],ans = 0.0000,eps = 1e-8;
int n,sta,num[N];
int main()
{
ios::sync_with_stdio(0);cin.tie(0);cout.tie(0);
cin>>n;
n = (1 << n) - 1;
for(int i = 0;i <= n;i++) cin>>p[i],num[i] = num[i >> 1] + (i & 1);
for(int mid = 1;mid <= n;mid <<= 1)
for(int i = 0,r = mid << 1;i + mid <= n;i += r)
for(int j = i;j < i + mid;j++)
p[mid + j] += p[j];
for(int i = 1;i <= n;i++)
if(1 - p[n ^ i] > eps)
ans += ((num[i] & 1) ? 1 : -1) * (1.00 / (1 - p[n ^ i]));
if(ans < eps) cout<<"INF";
else cout<<fixed<<setprecision(9)<<ans;
return 0;
}