FFT & NTT & FWT
只是学习笔记，真心推荐 cmd ，他讲的真的细到把所有的前置知识都讲了一遍。\
FFT & NTT & FWT 大杂烩
首先我们引入卷积的概念，对于两个多项式进行卷积，形如：
\
有
\
其中 $ F_n , f_i , g_j $ 均表示多项式 $ F , f, g $ 的某一项系数。
一般来说两个多项式进行朴素实现是 $ O(mn) $ ， $ m ,n $ 分别是两个多项式的次数。
那么就下来介绍几种对于特殊的卷积可以加快速度的算法。
<hr>一、FFT
FFT 是针对于加法卷积（也就是多项式相乘）的快速计算方法。
大体思路
首先我们知道对于一个 $ n $ 次多项式可以由 $ n+1 $ 个不同的点来表示，而对于两个多项式 $ f(x) , g(x) $ 相乘，我们可以通过 $ n+1 $ 个点在 $ O(n) $ 时间复杂度内完成，然后在通过插值将得到的 $ n+1 $ 个点还原成一个多项式。这就是 FFT 的大体思路。
但是朴素求值和插值都是 $ O(n^2) $ 的，所以我们现在要优化求值和插值的过程。
加速
前置：单位根
首先大家应该学过复数，接下来介绍单位根 。
有一个单位圆，我们通过 $ n $ 条线将它均分成了 $ n $ 份，我们依次编号 $ 0,1,2,\cdots ,n-1 $ ，他们对应的复数就是单位根 $ w_n^{i} $
如图是 $ n = 4 $ 的情况，他恰好对应着坐标轴。
有一天FFT的发明者 傅里叶 突然想把 $ w_n^i $ 代入到多项式中，然后发现这玩意有很好的性质以至于可以快速求值和插值。
首先由几个性质在几何角度看来非常容易证明：
1、 $ w_n^0 = 1 $
2、 $ w_n^k = w_n^{ k \mod n } $
3、 $ (w_n^k) ^j = w_n^{jk} $
4、 $ w_n^k \times w_n^j = w_n^{k+j} $
5、 $ w_{2n}^{2k} = w_n^k $
6、 $ n $ 是偶数， $ w_n^{k + \frac{n}{2}} = - w_n^k $
<hr>主体
那么接下来讲解 FFT 的主要过程：
求值：
对一个长度为 $ n $ 的序列，它对应着一个 $ n-1 $ 次多项式 $ F(x) $ 。
对于 $ F(x) $ ，我们按奇偶把他一分为二，下标为偶数的在前面，称为 $ Fl(x) $ ，下标为奇数的在后面，称为 $ Fr(x) $ 。
那么我们尝试写出 $ F(x) $ 与 $ Fl(x) , Fr(x) $ 的关系式：
\
接下来就是单位根出场的时候了，我们将 $ w_n^k $ 代入 $ F(x) $ 和上述式子：
对于 $ k \lt \frac{n}{2} $ 有
\[\begin{aligned}F( w_n^k ) &= Fl( w_n^{2k} ) + w_n^k Fr( w_n^{2k} ) \\&= Fl(w_{\frac{n}{2}}^k) + w_n^k Fr( w_{\frac{n}{2}}^k )\end{aligned}\]
对于 $ k \lt \frac{n}{2} $ ，我们代入 $ w_n^{k + \frac{n}{2}} $ 有
\[\begin{aligned}F(w_n^{k+\frac{n}{2}}) &= Fl(w_n^{2k+n}) + w_n^{k+\frac{n}{2}} Fr(w_n^{2k+n}) \\&= Fl(w_{\frac{n}{2}}^k) - w_n^k Fr(w_{\frac{n}{2}}^k)\end{aligned}\]
上面的式子意味着什么，它意味着我们可以通过分治去 $ O(n \log n) $ 计算 $ F(w_n^k) $ 的值。
好，我们现在已经可以 $ O(n \log n) $ 的将一个多项式转化为 $ ( w_n^k , F( w_n^k ) ) $ 的点值表达式了，然后可以 $ O(n) $ 点乘出最终的多项式 ， 然后我们现在需要快速将这个多项式插回去。
插值：
设我们刚才把 $ F(x) $ 转成了 $ G(x) $ 。那么我们直接不加证明的给出：
\
算了还是证明一下：
我们将 $ G(i) $ 代入
\[\begin{aligned}右边 &= \sum_{i=0}^{n-1} w_n^{-ik} \sum_{j=0}^{n-1} w_n^{ij} F(j) \\&= \sum_{j=0}^{n-1} F(j) w_n^{ij} \sum_{i=0}^{n-1} w_n^{-ik} \\&= \sum_{j=0}^{n-1} F(j) \sum_{i=0}^{n-1} w_n^{i(j-k)}\end{aligned}\]
分类讨论，对于 $ j=k $ 的情况：
贡献为 $ F(k) \sum_{i=0}^{n-1} w_n^0 = n F(k) $
对于 $ j \not = k $ 的情况：
贡献为 $ F(j) \sum_{i=0}^{n-1} w_n^{i p} $ ， 等比数列求和， $ \sum_{i=0}^{n-1} w_n^{i p} = \frac{ w_n^{ n p } - 1 }{ w_n^p - 1 } = 0 $ ，也就是没有贡献！
所以最终就有
\[右边 = n F(k) = 左边\]
所以我们将插值也转成了求值，同样可以 $ O(n \log n) $ 求。
实现
朴素实现
现在可以再来看一开始的流程图
整个过程就比较清晰了吧。
给出模板题： P3803 【模板】多项式乘法（FFT）
码(未卡常版本)#include<bits/stdc++.h>using namespace std;typedef long long ll;typedef pair<int,int> pii;typedef unsigned long long ull;#define mk make_pair#define ps push_back#define fi first#define se secondconst int N=5e6+10,inf=0x3f3f3f3f;const ll mod=1e9+7,linf=0x3f3f3f3f3f3f3f3f;inline ll read(){ char c=getchar();ll x=0;bool f=0; while(!isdigit(c))f=c=='-'?1:0,c=getchar(); while(isdigit(c))x=(x<<1)+(x<<3)+(c^48),c=getchar(); return f?-x:x;}const double PI=acos(-1);int n,m;struct jj{ // 复数结构体 double x,y; inline jj operator +(const jj&k){return {x+k.x,y+k.y};} inline jj operator -(const jj&k){return {x-k.x,y-k.y};} inline jj operator *(const jj&k){return {x*k.x-y*k.y,x*k.y+y*k.x};} inline jj operator /(const jj&k){return {(x*k.x-y*k.y)/(k.x*k.x+k.y*k.y),(y*k.x+x*k.y)/(k.x*k.x+k.y*k.y)};}}f,g,tp;inline void fft(jj *f,int l,int r,bool fl){ if(l==r)return; int mid=l+r>>1,len=r-l+1; for(int i=l;i<=r;++i) tp=f; for(int i=l;i<=mid;++i){ f=tp[(i-l)*2+l],f=tp[(i-l)*2+l+1];// 按照奇偶分裂 } fft(f,l,mid,fl),fft(f,mid+1,r,fl);// 分治计算 jj op={cos(2*PI/len),sin(2*PI/len)},now={1,0};// op=w_n^1 if(fl)op.y*=-1;// 如果 fl 为 1 表示是插值，此时 op=w_n^-1 for(int i=l;i<=mid;++i){ tp=f+now*f; tp=f-now*f; now=now*op; } for(int i=l;i<=r;++i) f=tp;}signed main(){ #ifndef ONLINE_JUDGE freopen("in.in","r",stdin); freopen("out.out","w",stdout); #endif ios::sync_with_stdio(0),cin.tie(0),cout.tie(0); n=read(),m=read(); for(int i=0;i<=n;++i) f={read(),0}; for(int i=0;i<=m;++i) g={read(),0}; m+=n;n=1; while(n<=m)n<<=1; --n; fft(f,0,n,0),fft(g,0,n,0); for(int i=0;i<=n;++i) f=f*g; fft(f,0,n,1); for(int i=0;i<=m;++i){ cout<<(int)(f.x/(n+1)+0.499)<<' ';//有精度问题，四舍五入 }}常数非常大，luogu最大点跑了 800 ms ，考虑一些优化。
优化
首先可以把递归改为非递归版的迭代。其次一个重要优化是 蝴蝶变换 。他说的是我们递归的过程会把 $ F(x) $ 按奇偶分裂，我们最终分裂好的整个序列原来在 $ i $ （下标从 0 开始）位置的数，现在在 $ i $ 按照二进制翻转后的数为下标的地方，比如说 $ 1010111_2 (87_{10}) $ ，在翻转后是 $ 1110101_2 (117_{10}) $ 。
首先证明这个事：
我们脑动模拟一下递归过程，每次按照偶在左半边，奇在右半边去分，相当于看二进制下最后一位，按照大小确定了最终二进制下第一位的大小。
有点抽象，但是有图：
以下将经过 $ i $ 次分裂的序列称为 $ F_i $
这是一开始的 $ F_0(x) $
然后我们按照最后一位的大小，分开
这时可以看到 $ F_1 $ 的下标的最后一位和 $ F_0 $ 下标的第一位相同，因为一开始我们是按第一位大小排的，后来我们按最后一位大小排的。同时前一半和后一半没有关系了，他们各自内部都是按大小排好的，而且此时所有数的最后一位是啥不重要了，相当于所有数最后一位没了，倒数第二位顺延为最后一位，继续递归处理。
那么下一层也就根据 $ F_1 $ 的最后一位，也就是倒数第二位确定了第二位的是啥，然后 $ \cdots $ ，也就相当于二进制翻转。
而事实也确实如此。
好，设 $ pos_i $ 表示 $ i $ 翻转后的答案，那么我们怎么求 $ pos_i $ ? 可以考虑先翻转最后一位之前的，由 $ pos_{i/2} $ 得到，同时翻转最后最后一位，于是有 pos=(pos>>1)|(i&1?n>>1:0) ， $ n $ 是序列长度。
所以我们就获得了一个常数较小的 FFT：
码(卡常版)#include<bits/stdc++.h>using namespace std;typedef long long ll;typedef pair<int,int> pii;typedef unsigned long long ull;#define mk make_pair#define ps push_back#define fi first#define se secondconst int N=1e6+10,inf=0x3f3f3f3f;const ll mod=1e9+7,linf=0x3f3f3f3f3f3f3f3f;inline ll read(){ char c=getchar();ll x=0;bool f=0; while(!isdigit(c))f=c=='-'?1:0,c=getchar(); while(isdigit(c))x=(x<<1)+(x<<3)+(c^48),c=getchar(); return f?-x:x;}const double PI=acos(-1);struct jj{ double x,y; inline jj operator +(const jj&k){return {x+k.x,y+k.y};} inline jj operator -(const jj&k){return {x-k.x,y-k.y};} inline jj operator *(const jj&k){return {x*k.x-y*k.y,x*k.y+y*k.x};}}f,g,ji;int pos;inline void fft(jj *f,int n,bool fl){ for(int i=0;i<n;++i) if(i<pos)swap(f,f]);//pos 表示预处理二进制翻转后的位置 for(int len=2,mid;len<=n;len<<=1){ mid=len>>1; jj op=ji,lp;//ji表示预处理 w_n^1 ，因为 sin,cos 太慢了，所以预处理会少调用几次 if(fl)op.y*=-1; for(int i=0;i<n;i+=len){ jj now={1,0}; for(int j=i;j<i+mid;++j){ lp=now*f; f=f-lp;f=f+lp; now=now*op; } } }}int n,m;signed main(){ #ifndef ONLINE_JUDGE freopen("in.in","r",stdin); freopen("out.out","w",stdout); #endif ios::sync_with_stdio(0),cin.tie(0),cout.tie(0); n=read(),m=read(); for(int i=0;i<=n;++i) f={read(),0}; for(int i=0;i<=m;++i) g={read(),0}; m+=n;n=1;ji={1,0}; while(n<=m)n<<=1,ji={cos(2*PI/n),sin(2*PI/n)}; for(int i=0;i<n;++i) pos=(pos>>1)|(i&1?n>>1:0); fft(f,n,0),fft(g,n,0); for(int i=0;i<=n;++i) f=f*g; fft(f,n,1); for(int i=0;i<=m;++i){ cout<<(int)(f.x/n+0.499)<<' '; }}luogu最大点 469 ms ，而且出奇的好写。
例题
P1919 【模板】高精度乘法 | A*B Problem 升级版
另一个板子题，FFT 加速即可。
P3338 力
给出序列 $ q $ ，有
\
求所有的 $ E_i $
我们考虑如何把 $ E_i $ 转成卷积的形式。首先把烦人的除法去掉，令 $ G(x) = \frac{1}{x^2} $ ，特殊的， $ G(0) = 0 $ ，同时令 $ q_0 = 0 $ ，那么有
\
前面已经是卷积的形式了，无需转化，后面的不是，我们考虑给他转成卷积的形式。
画个图看看：
他们是这样相乘的，但是卷积是这样相乘的：
所以我们考虑给 $ q(x) $ 转过来，变成 $ Q(x) $ ，所以后面的式子就变成了
\[\sum_{i=0}^{k=n-j} Q(i) G(k-i)\]
对两者 FFT 然后加起来即可。
二、 NTT
NTT ，中文 快速数论变换 是用来解决卷积过程中需要取模的快速加法卷积，说白了就是可以取模的 FFT。
这个东西就是直接把 $ w_n^1 $ 替换成了 $ g^{ \lfloor \frac{mod-1}{n} \rfloor } $ ，其中 $ mod $ 是模数，$ g $ 是模数的一个原根，$ n $ 是序列长度 ，但是原根这个东西我也不会，目前只知道 998244353 的原根是 3 ，所以 NTT 到此为止。
码#include<bits/stdc++.h>using namespace std;typedef long long ll;#define int lltypedef pair<int,int> pii;typedef unsigned long long ull;#define mk make_pair#define ps push_back#define fi first#define se secondconst int N=1e6+10,inf=0x3f3f3f3f;const ll mod=998244353,linf=0x3f3f3f3f3f3f3f3f,g=3;inline ll read(){ char c=getchar();ll x=0;bool f=0; while(!isdigit(c))f=c=='-'?1:0,c=getchar(); while(isdigit(c))x=(x<<1)+(x<<3)+(c^48),c=getchar(); return f?-x:x;}inline ll qpow(ll x,ll y){ ll ans=1; while(y){ if(y&1)ans=ans*x%mod; x=x*x%mod;y>>=1; } return ans;}ll ng=qpow(g,mod-2);// g 的逆元int pos;inline void ntt(ll *f,int n,bool fl){ for(int i=0;i<n;++i) if(i<pos)swap(f,f]); for(int len=2,mid;len<=n;len<<=1){ int op=qpow(fl?g:ng,(mod-1)/len),tp;mid=len>>1; for(int i=0;i<n;i+=len){ ll now=1; for(int j=i;j<i+mid;++j){ tp=f*now%mod; f=f-tp; if(f<0)f+=mod; f+=tp; if(f>=mod)f-=mod; now=now*op%mod; } } }}int n,m;ll f,ff;signed main(){ #ifndef ONLINE_JUDGE freopen("in.in","r",stdin); freopen("out.out","w",stdout); #endif ios::sync_with_stdio(0),cin.tie(0),cout.tie(0); n=read(),m=read(); for(int i=0;i<=n;++i) f=read(); for(int i=0;i<=m;++i) ff=read(); m+=n;n=1; while(n<=m)n<<=1; for(int i=0;i<n;++i) pos=(pos>>1)|(i&1?n>>1:0); ntt(f,n,0),ntt(ff,n,0); for(int i=0;i<n;++i) f=f*ff%mod; ntt(f,n,1); ll ny=qpow(n,mod-2); for(int i=0;i<=m;++i) cout<<f*ny%mod<<' ';}但是我的 NTT 好像没有 FFT 快？？？
三、 FWT
在 OI 中，FWT是用来解决有关位运算卷积的快速卷积算法。
即对于 $ c_i = \sum_{i = j \oplus k } a_j b_k $ ，在知道 $ a,b $ 的情况下，快速计算 $ c $ 。
或运算
此时 $ \oplus \to | $ 。
我们设 $ A,B,C $ 分别对应 $ a,b,c $ 进行了 FWT 后的序列，我们定义 $ A_i = \sum_{i=i|j} a_j $ ，那么自然有 $ C_i = A_i \cdot B_i $ 。
考虑如何快速求 $ A $ ，我们考虑初始序列 $ a $ ，它的长度 $ n $ 是 2 的整数次幂，我们将 $ a $ 一分为二，左边为 $ a_0 $ ，右边为 $ a_1 $ ，两者没有任何关系 ，直接递归求出 $ A_0 , A_1 $ ，此时需要需要合并 $ A_0,A_1 $ 为 $ A $ ，即，把 $ a_0 $ 向 $ A_1 $ 的贡献考虑进来， $ a_1 $ 向 $ A_0 $ 中的贡献考虑进来。
注意到将 $ a_1 $ 向 $ A_0 $ 考虑时，其中 $ a_1 $ 的目前最高位上一定是 1 ，而 $ a_0 $ 中的数一定是 0 ，所以 $ a_1 $ 不可能向 $ A_0 $ 做贡献。
再看 $ a_0 $ 向 $ A_1 $ 做贡献，注意到 $ \large a_{0_j} | a_{1_i} = a_{1_i} $ 当且仅当 $ \large a_{0_j} | a_{0_i} = a_{0_i} $ ，即能向 $ a_{0_i} $ 做贡献的都能向 $ a_{1_i} $ 做贡献，所以有
\
其中 $ + $ 表示多项式相加， $ he(x,y) $ 表示将 $ x,y $ 两个序列顺序连接。所以我们可以分治去求 $ A $ 。
类似的可以退出逆运算为
\
也就可以在得到 $ C $ 之后反推出 $ c $ 。
<hr>注
感觉就和 FFT ， NTT 一样， FWT(|) 找了一种特殊运算，使得可以分治 $ O(n \log n) $ 的去求值和插值。同时分治的过程可以看成逐渐考虑二进制上每一位的影响。
<hr>与运算
类比或运算，这里直接给出式子
求值：
\
插值：
\
<hr>异或运算
异或运算是 FWT 中最难理解的一个，是一种非常神奇的构造。
以下 $ \oplus \to xor $
首先我们定义 $ x \circ y = popcount( x \oplus y ) \mod 2 $ ，然后我们先给出 $ ( x \circ y ) \oplus ( x \circ z ) = x \circ ( y \oplus z ) $ 。
证明：
我们一位一位的去考虑，在一开始，两边都是 0，然后考虑新的一位，如果要做出改变，首先 $ x $ 在这位上得是 1，否则不用考虑。然后看 $ y \oplus z $ 的值，如果为 1 ，说明两个数一个是 1，一个是 0，所以此时左边右边的值都会改变，不影响等式成立。如果 $ y \oplus z = 0 $ ，说明两者一样，此时两边的值都不会改变，也不影响等式成立。
综上， $ ( x \circ y ) \oplus ( x \circ z ) = x \circ ( y \oplus z ) $ 。
接下来我们定义 $ A_i = \sum_{i \circ j = 0} a_j - \sum_{i \circ j = 1} a_j $ ，那么有
\[\begin{aligned}A_i \cdot B_i &= ( \sum_{i \circ j = 0} a_j - \sum_{i \circ j = 1} a_j ) \cdot ( \sum_{i \circ j = 0} b_j - \sum_{i \circ j = 1} b_j ) \\&= (\sum_{i \circ j = 0} a_j \sum_{i \circ k =0} b_k + \sum_{i \circ j = 1} a_j \sum_{i \circ k =1} b_k) - ( \sum_{i \circ j = 0} a_j \sum_{i \circ k =1} b_k + \sum_{i \circ j = 1} a_j \sum_{i \circ k =0} b_k ) \\&= \sum_{ i \circ (j \oplus k) =0 } a_j b_k - \sum_{ i \circ (j \oplus k) =1 } a_j b_k \\&= C_i\end{aligned}\]
接下来考虑如何快速计算 $ A $ ，我们仍然尝试将它分为 $ A_0,A_1 $ 递归求解，接下来考虑合并。
首先举一个例子：
假设 $ A $ 长度为 8，我们写成二进制形式
然后我们一分为二，递归求解，但是注意递归到下一层后，他们在本层的第一位将不再考虑，也就是
这也就是为什么之前注里面说分治的过程是逐渐考虑二进制位的影响。
理解了这里，下面就很好理解了。
对于 $ A $ ，分为 $ A_0 ,A_1 $ 之后，我们新考虑一位，定义 $ x' $ 为 考虑新一位的 $ x $ ，那么有 $ a_{0_i} ' \circ a_{1_j} ' = a_{1_i} \circ a_{1_j} , a_{0_i} ' \circ a_{0_j} ' = a_{0_i} \circ a_{0_j} , a_{1_i} ' \circ a_{1_j} ' = (a_{1_i} \circ a_{1_j}) \oplus 1 , a_{1_i} ' \circ a_{0_j} ' = a_{0_i} \circ a_{0_j} $ ，于是我们有
\
同时有逆推
\
模板题
P4717 【模板】快速莫比乌斯/沃尔什变换 (FMT/FWT)
码#include<bits/stdc++.h>using namespace std;typedef long long ll;#define int lltypedef pair<int,int> pii;typedef unsigned long long ull;#define mk make_pair#define ps push_back#define fi first#define se secondconst int N=1e6+10,inf=0x3f3f3f3f;const ll linf=0x3f3f3f3f3f3f3f3f,mod=998244353;inline ll read(){ char c=getchar();ll x=0;bool f=0; while(!isdigit(c))f=c=='-'?1:0,c=getchar(); while(isdigit(c))x=(x<<1)+(x<<3)+(c^48),c=getchar(); return f?-x:x;}inline void OR(int *f,int n,int op){ for(int len=2;len<=n;len<<=1){ int mid=len>>1; for(int j=0;j<n;j+=len){ for(int k=j;k<j+mid;++k){ f=(f+f*op)%mod; } } }}inline void AND(int *f,int n,int op){ for(int len=2;len<=n;len<<=1){ int mid=len>>1; for(int j=0;j<n;j+=len){ for(int k=j;k<j+mid;++k) f=(f+f*op)%mod; } }}inline void XOR(int *f,int n,int op){ for(int len=2;len<=n;len<<=1){ int mid=len>>1,tp; for(int j=0;j<n;j+=len){ for(int k=j;k<j+mid;++k){ tp=f; f=(f+f)*op%mod; f=(tp-f)*op%mod; } } }}inline ll qpow(ll x,ll y){ ll ans=1; while(y){ if(y&1)ans=ans*x%mod; x=x*x%mod;y>>=1; } return ans;}int a,b,c,n;signed main(){ #ifndef ONLINE_JUDGE freopen("in.in","r",stdin); freopen("out.out","w",stdout); #endif ios::sync_with_stdio(0),cin.tie(0),cout.tie(0); n=read();int m=1<<n; for(int i=0;i<m;++i) a=read(); for(int i=0;i<m;++i) b=read(); OR(a,m,1),OR(b,m,1); for(int i=0;i<m;++i) c=a*b%mod; OR(c,m,-1);OR(a,m,-1),OR(b,m,-1); for(int i=0;i<m;++i) cout<<(c+mod)%mod<<' '; cout<<'\n'; AND(a,m,1);AND(b,m,1); for(int i=0;i<m;++i) c=a*b%mod; AND(a,m,-1),AND(b,m,-1),AND(c,m,-1); for(int i=0;i<m;++i) cout<<(c+mod)%mod<<' '; cout<<'\n'; XOR(a,m,1),XOR(b,m,1); for(int i=0;i<m;++i) c=a*b%mod; XOR(c,m,qpow(2,mod-2)); for(int i=0;i<m;++i) cout<<(c+mod)%mod<<' ';}小结
感觉三者都是通过构造找了一些特殊的点或特殊的运算，使得我们可以通过这写特殊的性质达到快速计算的效果，当然这一切的前提都是因为我们用到的 $ O(n) $ 点乘对点是没有要求的，这就支撑着我们可以通过构造特殊点来完成快速卷积。
页:
[1]