Note: 多项式卷积
参考资料：[学习笔记&教程] 信号, 集合, 多项式, 以及各种卷积性变换 (FFT,NTT,FWT,FMT) , 快速傅里叶变换(FFT)详解
前言
不会数学。菜完了
FWT
FWT 通过线性变换将问题转化为点乘形式，从而将 \(O(n^2)\) 的复杂度降到 \(O(n\log n)\) ——DeepSeek
1.或卷积
构造 \(T(A)_i=\sum_{j\subseteq i}A_j\) ，那么有
那么太好了，只要求出 \(T(A)\) 和 \(T(B)\) 就可以求出 \(T(A\ *_{\oplus} B)\) 了。\(T(A)\) 可以分治求出：将 \(A\) 分为最高位为 0、最高位为 1 的两部分 \(A_0,A_1\)（这也相当于对半砍）。前后两部分按位比较可以发现他们除了最高位不同剩下啥都相同 是这样吗，是这样吧，太好了现在只需要看最高位的贡献了。显然，在或运算下最高位为 0 的可以给最高位为 1 的贡献，反之不能，那么就能有递推式
（人话：前半拉不变，后半拉每项加上对应的前半拉）
现在需要构造搭配的逆变换由 \(T(A\ * B)\) 求出原卷积。很显然啊，加上的现在全给它减掉就行了，也就是
（懒得翻译成人话了）
那么这样经过两次变换、一次计算和逆变换后就可以求出卷积数组了
void fwt_or(int *a,int op)
{
for (int i=2;i<=n+1;i<<=1)
for (int p=i>>1,j=0;j<n;j+=i)
for (int k=j;k<j+p;k++) a[k+p]=(a[k+p]+op*a[k]+MOD)%MOD;
//第一维枚举分治的区间长度 i，第二维枚举区间首项 j，第三维枚举区间前半拉的下标 k，后半拉的下标就是 k+p
}
2.与卷积
和或卷积同工异曲
构造 \(T(A)_i=\sum_{i\subseteq j}A_j\) ，那么肯定也有\(T(A)_i\times T(B)_i=T(A\ * B)_i\)。
现在考虑怎么计算这个变换。和或运算相反，与运算最高位为 1 的数可以给最高位 0 的数贡献，反之不能，所以有转移式
逆变换也就直接出来了
void fwt_and(int *a,int op)
{
for (int i=2;i<=n+1;i<<=1)
for (int p=i>>1,j=0;j<n;j+=i)
for (int k=j;k<j+p;k++) a[k]=(a[k]+op*a[k+p]+MOD)%MOD;
}
3.异或卷积
异或有个小性质：\((i \&k)\oplus(j\&k)=(i\&j)and\ k\)。设 \(d(i)\)表示 \(popcount(i)\) 的奇偶性，那么好，\(j\) 每多个 1 的出现 \(d(i\oplus j)\) 就会翻转一次
然后我也不道啊，这里的变换为 \(T(A)_i=\sum_{d(i \&j)=0}A_j-\sum_{d(i \&j)=1}A_j\) ，然后也有 \(T(A)_i\times T(B)_i=T(A\ * B)_i\)
这里变换的计算不像前两个那样好算，因为左右半拉会互相贡献。虽然但是，最高位为 0 时与运算不会对其造成影响，所以可以直接加；反之，最高位为 1 时与运算会翻转 \(d(i)\) 的奇偶性，所以要取反
逆变换直接逆逆完了
void fwt_xor(int *a,int op)
{
for (int i=2;i<=n+1;i<<=1)
for (int p=i>>1,j=0;j<n;j+=i)
for (int k=j;k<j+p;k++)
{
int p1=a[k],p2=a[k+p];
a[k]=((p1+p2)%MOD*op)%MOD,a[k+p]=((p1-p2+MOD)%MOD*op)%MOD;
}
}
点击查看完整代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N=(1<<17)+3;
const int MOD=998244353;
int n,A[N],B[N];
int a[N],b[N],inv2;
void ea() { for (int i=0;i<=n;i++) a[i]=A[i],b[i]=B[i]; }
void eaa() { for (int i=0;i<=n;i++) a[i]=(a[i]*b[i])%MOD; }
void out()
{
for (int i=0;i<=n;i++) cout<<a[i]<<" ";
cout<<"\n";
}
void fwt_or(int *a,int op)
{
for (int i=2;i<=n+1;i<<=1)
for (int p=i>>1,j=0;j<n;j+=i)
for (int k=j;k<j+p;k++) a[k+p]=(a[k+p]+op*a[k]+MOD)%MOD;
}
void fwt_and(int *a,int op)
{
for (int i=2;i<=n+1;i<<=1)
for (int p=i>>1,j=0;j<n;j+=i)
for (int k=j;k<j+p;k++) a[k]=(a[k]+op*a[k+p]+MOD)%MOD;
}
void fwt_xor(int *a,int op)
{
for (int i=2;i<=n+1;i<<=1)
for (int p=i>>1,j=0;j<n;j+=i)
for (int k=j;k<j+p;k++)
{
int p1=a[k],p2=a[k+p];
a[k]=((p1+p2)%MOD*op)%MOD,a[k+p]=((p1-p2+MOD)%MOD*op)%MOD;
}
}
int qsm(int a,int b)
{
int res=1;
while (b)
{
if (b&1) res=(res*a)%MOD;
a=(a*a)%MOD,b>>=1;
}
return res;
}
signed 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>>A[i];
for (int i=0;i<=n;i++) cin>>B[i];
ea(),fwt_or(a,1),fwt_or(b,1),eaa(),fwt_or(a,-1),out();
ea(),fwt_and(a,1),fwt_and(b,1),eaa(),fwt_and(a,-1),out();
ea(),fwt_xor(a,1),fwt_xor(b,1),eaa(),fwt_xor(a,inv2=qsm(2,MOD-2)),out();
return 0;
}
FFT
前置知识
复数
设实数 \(a,b\)，有 \(i^2=-1\)，形如 \(a+bi\) 的数成为复数。其中 \(i\) 被称为虚数单位。
在复平面中，从原点到 \((a,b)\) 的向量表示复数 \(a+bi\)。复数的乘法与向量的乘法相同，几何意义下模长相乘、幅角相加。
单位圆与单位根
在平面直角坐标系上，圆心为原点、半径为单位长度的圆为单位圆。
在复数域下，方程 \(x^n=1\) 的 \(x\) 被称为 \(n\) 次单位根，在几何意义下，这 \(n\) 个解位于单位圆上，构成正 \(n\) 边形的顶点，其中一个顶点为 \(1\)。
下文设幅角最小且为正的 \(n\) 次单位根为 \(\omega_n\)，那么其余单位根则为 \(\omega_n^2,\omega_n^3...\omega_n^n\)。根据欧拉公式有
显然有 \(\omega_{2n}^{2k}=\omega_{n}^k,\omega_{2n}^{k}=-\omega_{2n}^{k+n},\omega_n^{k+n}=\omega_{n}^k\)
FFT
由于单位根有着很好的性质，所以可以分治求出 \(F(x)\) 所有单位根的点值。设多项式 \(A(x)=\sum_{i=0}^n a_i x^i\)，根据次数的奇偶性分类，令
那么显然有
将 \(\omega _n^k\) 代入 \(x\)，有
单有这个式子无法缩小问题的规模。考虑令 \(k\leqslant\frac{n}{2}\)，有
有了上面两个式子，我们发现在下一层递归求出 \(\frac{n}{2}\) 个点值后就得到这层的 \(n\) 个点值了。在整个递归过程中，需要求出 \(n\log n\) 个点值，每个值都可 \(O(1)\) 转移，所以求值的复杂度为 \(O(n\log n)\)
IDFT
求出 \(H(x)\) 的点值后还需插值求出原函数。改天再写 \(O\omega O\)
浙公网安备 33010602011771号