快速沃尔什变换
其实是不太想接触这类东西的。
但是万一他就考你个子集卷积呢?万一呢万一呢?万一 Θ(3n) 比 Θ(2n) 少很多昏呢?
可能大概或许是把 command_block’s blog 抄了一遍?
不过也加上了一些自己的理解和规范化的符号吧。
目录:
- 快速变换的本质
- FWT 的推导
- 位矩阵理论
- OR 运算
- AND 运算
- XOR 运算
- 子集卷积
- FMT
- 集合幂级数
- 性质
- 代码模板
- 对进制的拓展
- OR - 每维取 max
- AND - 每维取 min
- 数论领域:将每个质数视作一维
- 例题
快速变换的本质
首先要搞明白快速____变换都是一类什么东西。
简单来说就是我们要求卷积,卷积就是如下操作:
设所有数组大小是 Θ(n) 的,你有两个数组 f,g,要求一个新数组 h。
定义 h 是 f 与 g 的卷积 h=f∘g,当且仅当:
h_k=\sum_i\sum_j[k=i\oplus j]f_ig_j$$($\oplus$ 表示一种对下标的运算,例如 $\text{FFT}$ 中 $\oplus$ 表示 $+$。)
我一般喜欢将上式简记成 $f_i\times g_j\xrightarrow+h_{i\oplus j}$。
朴素做显然是 $\Theta(n^2)$ 的,但**假如存在一种快速的变换**,我们称 $f$ 的变换是 $f'$,$f'$ 再变换一次是 $f$。
> 注:这里为了符号的简洁性我们假装变换两次之后回到原序列,而事实上这种变换并不一定有这种性质,我们还需要一种**逆变换**。
>
> 这里能这么简记的主要原因是:通常情况下,对一个序列**变换两次并无意义**,且变换后序列总是和变换后序列运算,变换前序列总是和变换前序列运算,且除了变换和逆变换之外的运算均**不会改变一个序列是否处于变换态**这个性质。
>
> 简单来说就是一个正常人只需要数 $'$ 的个数就能确定这个数到底是变换态还是非变换态,就像数负号个数判断数的符号一样。
这种变换应存在如下性质:
- 其变换 $f\to f'$ 和逆变换 $f'\to f$ 都是 $\Theta(ns)$ 的($s$ 表示低于 $\Theta(n)$ 的复杂度)。
- 定义 $f*g$ 表示 $f$ 与 $g$ 的**点积**,$h=f*g\iff h_i=f_ig_i$,则 $(f\circ g)'=f'*g'$。
则我们求 $f\circ g$ 只需要求 $(f'*g')'$ 即可,复杂度 $\Theta(ns)$。
快速傅里叶变换正是基于如上思想使得我们能方便地对多项式进行乘法。
## FWT 的推导
### 位矩阵理论
> 注意:位矩阵的部分我的定义与 command_block 的定义并不相同,我定义的矩阵相当于其 blog 的矩阵的转置,因为我将其视作原序列到变换序列的转移系数,而 command_block 将其视为变换序列接受原序列的转移系数。
>
> 我认为我的定义要稍微好一点,至少解释 $\text{AND/OR-FWT}$ 可以做高维前后缀和的地方符号更通顺一点(
>
> 所以在 $\text{AND},\text{OR}$ 的部分请不要记错位矩阵的具体部分。
回到主题,现在 $f_i\times g_j\xrightarrow+h_{i\oplus j}$ 中的 $\oplus$ 是一种位运算了,那我们如何构造一种变换来完成如上操作呢?
虽然我们不知道这种变换是什么,但我们首先可以假装它是一种**线性变换**,即存在一个 $n\times n$ 的变换矩阵 $c$ 满足 $f'_j=\sum_{i}f_ic_{i,j}$。
(题外话:你发现逆变换就形如一个类似反演的东西,不考虑复杂度的情况下,你只需要对矩阵求个逆即可实现逆变换,当然这也就必须要求我们的变换矩阵存在逆。)
我们对定义式 $(f\circ g)'=f'*g'$进行恒等变换:
$$\forall l,\left(\newcommand{\e}{\iff}\begin{aligned}&\sum_k(f\circ g)_kc_{k,l}=\left(\sum_if_ic_{i,l}\right)\times\left(\sum_ig_ic_{i,l}\right)\\\e&\sum_k\left(\sum_i\sum_j[k=i\oplus j]f_ig_j\right)c_{k,l}=\sum_i\sum_jf_ig_jc_{i,l}c_{j,l}\\\e&\sum_k\sum_{i,j,i\oplus j=k}f_ig_jc_{i\oplus j,l}=\sum_i\sum_jf_ig_jc_{i,l}c_{j,l}\\\e&\sum_i\sum_jf_ig_jc_{i\oplus j,l}=\sum_{i}\sum_jf_{i}g_jc_{i,l}c_{j,l}\end{aligned}\right)
观察等式两边,发现你只需要满足 ∀i,j,l,ci⊕j,l=ci,lcj,l 就可以了。
首先由于你得保证矩阵有逆,所以收起你全 0 矩阵的想法。
以下在不引起歧义的情况下,将 xi 定义为 x 在二进制表示下的第 i 位。
我们考虑二进制数的特殊性,你发现如果你知道了 c0/1,0/1,那么你可以这样构造:
ci,j=k∏cik,jk
通过验证可以发现其满足我们的要求:
∀i,j,l,ci⊕j,l=k∏c(i⊕j)k,lk=k∏cik⊕jk,lk=k∏cik,lkcjk,lk=(k∏cik,lk)(k∏cjk,lk)=ci,lcj,l
所以现在就是要找一个初始矩阵,我们就可以找到这样一个很大很大的矩阵了。
但我们不能显式地使用矩阵来实现 f→f′,因为还是 Θ(n2) 的。
我们尝试采用分治算法,设 i∗=imod2m:
∀j,fj′=i=0∑2m−1fici,j=i=0∑2m−1−1fici,j+i=2m−1∑2m−1fici,j=i=0∑2m−1−1ficim,jmci∗,j∗+i=2m−1∑2m−1ficim,jmci∗,j∗=c0,jmi=0∑2m−1−1fici∗,j∗+c1,jmi=2m−1∑2m−1fici∗,j∗
关键的一步在第二行到第三行的变化,这里把第 m 位拆开,拆成两个子问题。
我们发现 ∑i=02m−1−1fici∗,j∗,∑i=2m−12m−1fici∗,j∗ 是两个规模减半的子问题,于是我们的复杂度分析形如普通分治,Θ(nlogn)。
不难发现,位运算的不同影响的仅仅是 c0/1,0/1,我们将这样一个 2×2 的矩阵称作该位运算的**「位矩阵」**。
OR 运算
由于代码里的 |,&,^,popcount(S) 在 LATEX 里面看着不是很好看,我们用 ∪,∩,⊕,∣S∣ 代替。
构造是一种人类智慧,由于构造的推导并不适合公式呈现,这里仅给出构造的结果:
c∪=[1011],c∪−1=[10−11]
注:仅仅将一个序列进行 OR-FWT 而不逆变换回去可以求高维前缀和,仅仅逆变换是高维前缀差分(仅限这一种矩阵)。
理解:观察得 ci,j=[i∪j=j],所以 fj′=∑ifici,j=∑i,i∪j=jfi=∑i⊆jfi,也就是子集和。
理解性记忆技巧:右上三角,逆矩阵将顶点变为 −1。
AND 运算
c∩=[1101],c∩−1=[1−101]
同理,仅 AND-FWT 可以求高维后缀和,仅逆变换是高维后缀差分。
观察得 ci,j=[i∩j=j],所以 fj′=∑ifici,j=∑i,i∩j=jfi=∑j⊆ifi,也就是超集和。
理解性记忆技巧:左下三角,逆矩阵将顶点变为 −1。
XOR 运算
c⊕=[111−1],c⊕−1=[0.50.50.5−0.5]
理解性记忆技巧:左上三角,右下角是负数。
注意 0.5 要换成模意义下 21。
c⊕,i,j=(−1)∣i∩j∣,化简时可能会用到。
题外话,在 i,j∈Z2n 意义下,∣i∩j∣mod2 其实就是 i 与 j 的点积,本文记为 i⊙j,∣i∩j∣mod2=0 其实就是 i 与 j 正交。
我们敏锐地注意到 i 与 j 的向量意义加法就是 i⊕j,那么我们就可以把加法乘法的结合律拿过来:
i⊙(j⊕k)=(i⊙j)⊕(i⊙k)=i⊙j⊕i⊙k
那 c⊕,i,j=(−1)i⊙j,很有用的一点是异或运算的奇偶性变化与加法一致,也就是 (−1)a⊕b=(−1)a+b=(−1)a(−1)b。
所以我们可以得到二级结论(其实也可以从位矩阵的定义推来):
c⊕,i,⨁jbj=(−1)i⊙⨁jbj=(−1)⨁ji⊙bj=j∏(−1)i⊙bj=j∏c⊕,i,bj
后文会用到此处的定义。
子集卷积
等价于求:
hk=i,j∑[i∪j=k][i∩j=∅]figj
考虑 i∪j=k,这就是朴素的 OR-FWT,那如何处理 i∩j=∅ 呢?
变变形, i∩j=∅⟺∣i∣+∣j∣=∣i∪j∣=∣k∣。
所以我们多开一维来记录元素个数即可。
FMT
实现难度低于 FWT,代码速度可能会快一点?
但是只能做 AND,OR。
一言以蔽之:将上述的“高维(前缀/后缀)(和/差分)”写一个真的那种东西上去,看起来实现难度就低很多。
似乎也有自己的一套推导过程,不过可以规约到 FWT 上嘛,那就略过吧。
集合幂级数
可能是对于一些记法的定义?
我对幂级数并没有很深刻的理解,但似乎许多多项式能做的集合幂级数都能做。
如果没有特殊指明,下文均用 x 表示自变量,∑cxn 表示幂级数,幂级数相乘表示其做卷积,通常卷积的类型是比较明确的。
若用数组表示幂级数,则使用 fi 下标形式来表示某一项,若使用 ∑cxn 来表示幂级数,则使用 [xi]∑cxn 来表示某一项。
性质
由于这是一个线性变换,所以集合幂级数跟线性代数是一家的。
- a′+b′=(a+b)′。
- ca′=(ca)′。
这意味着如果你把两个序列相加,那可以变换后直接加而不必变换回去。
所以我们通常把所有序列都 FWT,然后一顿操作变换回去。
代码模板
这是代码。
贴一下关键部分:
1 2 3 4 5 6 7 8 9 10 11 12 13
| #define ci const int ci iv2=p-p/2, OR[2][2]={1,1,0,1}, IOR[2][2]={1,p-1,0,1}, AND[2][2]={1,0,1,1}, IAND[2][2]={1,0,p-1,1}, XOR[2][2]={1,1,1,p-1}, IXOR[2][2]={iv2,iv2,iv2,p-iv2}; void fwt(ll a[],int n,ci c[2][2]) { for(int m=1;m<n;m<<=1) for(int o=0;o<n;o+=m<<1) F(i,o,o+m-1) tie(a[i],a[i+m])=make_pair((a[i]*c[0][0]+a[i+m]*c[1][0])%p, (a[i]*c[0][1]+a[i+m]*c[1][1])%p); }
|
注意 m 的上限不要写成 <=n,要不然 ⊕ 卷积会炸掉。
代码。
关键部分:
1 2 3 4 5 6
| F(i,0,(n=1<<(l=read()))-1) a[ppct(i)][i]=read(); F(i,0,n-1) b[ppct(i)][i]=read(); F(i,0,l) fwt(a[i],n,OR),fwt(b[i],n,OR); F(i,0,l) F(j,0,l) if(i+j<=l) F(k,0,n-1) c[i+j][k]=(c[i+j][k]+a[i][k]*b[j][k])%p; F(i,0,l) fwt(c[i],n,IOR); F(i,0,n-1) printf("%lld ",c[ppct(i)][i]);
|
需要注意的是,不要会错意,ci+j,k+ai,kbj,k 这个柿子是 FWT 意义下才成立的,所以不要 IFWT 之后再加上,要等到最后的时候再 IFWT。
对进制的拓展
将快速沃尔什变换视作对 Z2n 的变换,那我们能不能将其拓展到 Zkn 意义下呢?
事实证明是可以的,且这些拓展具有某种共通的规律。
由于精力有限,此处仅给出位矩阵的形式,FWT 的过程推导不再展开。
OR - 每维取 max
我们以 k=4 为例:
c∪=1000110011101111,c∪−1=1000−11000−11000−11
相信你可以发现某种规律。
原矩阵 c 是标准的右上三角,逆矩阵仅保留了主对角线和紧贴着右上的一侧。
AND - 每维取 min
k=4:
c∪=1111011100110001,c∪−1=1−10001−10001−10001
原矩阵 c 是标准的右上三角,逆矩阵仅保留了主对角线和紧贴着左下的一侧。
尝试构造 ⊕ 卷积?
但我们遗憾地发现,可爱的 ⊕ 的拓展是不进位加法,在任意模数的意义下需要分圆多项式。
但我们发现 k 维 FMT 还是那么的恰到好处,你可以直接高维前缀和/差分实现 ∩/∪ 卷积,十分舒适。
数论领域:将每个质数视作一维
考虑把每个数视作无穷维的向量,将每个质数视作一维,于是我们就有了一些向量,向量的加法是数的乘法,向量的 ∩ 是 gcd,向量的 ∪ 是 lcm。
同理我们可以定义此意义下的 FMT,其实就是狄利克雷(前缀/后缀)(和/差分)。
可以快速(Θ(n))求出 ck=∑(i,j)=kaibj。
例题
给定一个序列,问有多少种选择子序列的方式使得子序列 and 和为 0。
n≤106,ai≤106。
先来看传统容斥做法:
在集合意义下的容斥通常是把状态由 S 的答案变成包含 S 的答案和,即令 gT=∑S⊇TfS,求完 g 再高维差分。
设 fS 为包含 S 的方案数,sS=∑i[ai⊇S],则 fS=2sS−1,较为明显。
FWT 做法:
最最朴素的想法是背包,设 fi,j 表示前 i 个物品考虑完了,and 和是 j 的方案数,转移就是 fi,j←fi−1,j+∑ai∩k=jfi−1,k。
你考虑把他写成 and 卷积的形式,那就是每次 gS=[ai=S],fi←fi−1+fi−1∘g。
然后你发现 gS′ 就是 [S⊆ai],所有 gS=1 的位置都会让 fS′ 翻倍。
所以你直接算出那个容斥做法的 s,也就是每个位置都被翻过多少次倍,最后逆变换即可。
殊途同归,换了一个新角度理解。
给定 n 个可重集和 x,y,z,ai,bi,ci,每个可重集可以用 (a,b,c) 描述,表示这个可重集里有 x 个 a,y 个 b,z 个 c,求在每个数组里选一个数使得它们 ⊕ 和为 t 的方案数,对每个 t<2k 求出以上值。
n≤105,k≤17,x,y,z≤109。
还是考虑暴力设状态 fi,S 表示考虑了前 i 个数组异或出了 S 的方案数,转移就是 fi,S←∑T⊕ai=Sfi−1,Sgi,T,其中 gi,[ai,bi,ci]=[x,y,z]。
考虑写成异或卷积的形式,则最终结果就是 ∏(∘)igi。(∏(∘) 表示以 ∘ 为乘法的连乘。)
也就是我们要求 (∏(⋅)gi′)′。
回顾我们的位矩阵理论(我猜你们一定分得清位矩阵数组 ci,j 和原题的数组 ci),gi 最大的特点就是有效值少,我们把它拆开(第一维不参与 FWT):
gi,j′=k∑gi,kck,j=xcai,j+ycbi,j+zcci,j=x(−1)ai⊙j+y(−1)bi⊙j+z(−1)ci⊙j
(最后一步运用了 ⊕ 位矩阵中的理论。)
然后你就发现 gi,j′ 一共只有 23=8 种可能的取值,那你如果能分别求出这些取值的系数,不就可以直接快速幂做了吗。
那我们现在就是要求类似于 pj,u,v,w=∑i[ai⊙j=u][bi⊙j=v][ci⊙j=z]。
之后由于 3 这个数很小,存在讨论解方程的做法,但我想介绍更本质的可拓展的一种做法。
我们考虑拓展原题的限制,得到更普遍的解法,将题面改成这样:
给定 n 个可重集和 wj,ai,j(i∈[1,n],j∈[0,m)),表示第 i 个可重集里有 wj 个 ai,j,求在每个数组里选一个数使得它们 ⊕ 和为 t 的方案数,对每个 t<2k 求出以上值。
看到上面的最后一步,那我们现在要求的就是 pT,S=∑i∏j[ai,j⊙T=Sj](Si=0/1 则 wi 系数为 1/(−1)),其意义为“以 S 为 1/(−1) 系数合成的数对 (∏(⋅)gi′)T 贡献了几次方”。
大体可以看出 T 是一个最外层的变量,不能动什么手脚,那我们就考虑先固定这个 T,看怎么快速求解 2m 个 S 对应的 pT,S。
如果固定 T 的话,那 pT,S 就相当于具有 S 这种特质的可重集的个数,可重集就是 ai,j 中的 i 这一维,这对于理解后文有帮助。
求解这 2m 个未知数自然需要 2m 个线性无关的方程,这引导我们把 A⊆[0,m)∩N 的 2m 个取值变到一个幂级数 hA=∑ix⨁j∈Aai,j,求解 hA 对应的 FWT 数组 hA′,考虑:
hA′[xT]hA′=(i∑x⨁j∈Aai,j)′=i∑(x⨁j∈Aai,j)′=i∑c(j∈A⊕ai,j,T)=i∑j∈A∏(−1)ai,j⊙T
最后一步的推导过程在 ⊕ 位矩阵部分。
你考虑能不能把 [xT]hA′=∑i∏j∈A(−1)ai,j⊙T 和 pT,S=∑i∏j[ai,j⊙T=Sj] 建立什么对应关系,这样就可以有方程了。
(其实他俩长得挺像的。)
前文说过,pT,S 指的是具有 S 这种特质的可重集的个数(i 的个数),你考虑如果 S⊙A=0(偶数个 j 在 A 中起了作用),即有偶数个 1 在 −1 的指数上产生了贡献,那这 pT,S 个数对 [xT]hA′ 的贡献就是 (−1)2n=1,同理如果 S⊙A=1,则 pT,S 个数对 [xT]hA′ 的贡献系数是 −1。
即:
[xT]hA′=i∑j∈A∏(−1)ai,j⊙T=S∑(−1)S⊙ApT,S
这个关系就很明确了,你发现右边就是 [xA](∑SpT,SxS)′,所以你令 sA=[xT]hA′,再给 s 求个 IFWT,你就得到了 pT,S,这个过程写成代码就是 hA,T→hA,T′→pT,A′→pT,A 即 p=((h′)T)′,形如一个 FWT − 转置 − IFWT,就挺神奇的。
(本来的解方程是高消三方,但是你发现关系就是 IFWT,所以直接 2mm 搞定。)
梳理:
- ∀A,求出 hA=∑ix⨁j∈Aai,j 和它的 FWT 幂级数 hA′。
- 令 sA=[xT]hA′,再给 s 求个 IFWT,则 pT,S=sS′。
- 利用 pT,S 将以 S 为法则用 w 合成的数快速幂求出 (∏(⋅)gi′)T,再次 IFWT 即得到答案数组。
复杂度 Θ(n2m+(m+k)2m+k)。
口胡题
给定一个集合,求对于 k∈[1,n],至少选多少个数才能使它们的 gcd=k,并求出有多少种情况满足这个至少的情况。
ai,n≤3×105。
首先枚举这个 k,把所有倍数拿出来,相当于套了一个调和级数的 ln。
考虑到答案其实不会超过 max(ω)+1,所以我们直接枚举答案并计数。
然后你发现如果答案是 k,那你就看 (i=1∏k−1(∘gcd)s)1 有无值即可,其中 si 表示数 i 的出现次数。
给定长为 n 的数组 s,求对于每个满足以下条件的五元组 (a,b,c,d,e) 的权值和。
条件:
- sa∩sb=∅。
- ∣(sa∪sb)∩sc∩(sd⊕se)∣=1。
五元组 (a,b,c,d,e) 的权值定义为 f(sa∪sb)⋅f(sc)⋅f(sd⊕se),其中 f(0)=0,f(1)=1.f(i)=f(i−1)+f(i−2)。
n≤106,si<217。
我们考虑把
a,b,c,d,e,sa∩sb=∅,∣(sa∪sb)∩sc∩(sd⊕se)∣=1∑f(sa∪sb)⋅f(sc)⋅f(sd⊕se)
转化为
i,j,k,∣i∩j∩k∣=1∑f(i)a,b,sa∩sb=∅,sa∪sb=i∑1f(j)(c,sc=j∑1)f(k)d,e,sd⊕se=k∑1
接下来就都好处理了,三个括号分别用子集卷积,直接处理,异或卷积,最后全都与卷积起来。
题意冗长难以简化建议自己去看题。
设 eS 为 S 是否合法,fS 为选了 S 中的点的划分权值和,oS=∑i∈Swi,则转移:
fS←oSp1T⊊S∑f∁SToTpeT
我疑惑了许久究竟初值怎么设,没想到是 f0=1,返璞归真了属于是。
直接做是 Θ(3n) 的,考虑这个东西长得很像子集卷积,但是是子集卷自己,也就是传说中的半在线卷积。
其实半在线子集卷积跟原来没什么区别,只不过第一维枚举的是答案数组,这样保证了算到每个数的时候没有后效性。
对于前面的系数,我们需要每次算完了 IFWT 回来乘上系数再 FWT 回去。
请注意不要觉得乘完了 0 就万事无忧了,你之后有可能对这个数进行 0 次幂,所以请务必最后置为 0。