zoukankan      html  css  js  c++  java
  • FWT(快速沃尔什变换)小结

    在多项式卷积的处理中,我们实际上实现的是下面的一个式子

    [C_k=sum_{i+j=k}A_iB_j ]

    然而事实上有些和(sang)蔼(xin)可(bing)亲(kuang)的出题人,并不会让你直接求这样的式子,比如把+换成-什么的

    但是更多时候,你拿到手上的却是这样一个式子

    [C_k=sum_{iigoplus j=k}A_iB_j ]

    其中(igoplus)可能是(or,and,xor)中的一种或多种

    那么现在我们就会使用FWT(快速沃尔什变化)来优化时间复杂度

    前置芝士

    FFT

    理解思想即可,不在赘述

    正文

    (p.s. 注意下文的乘号应结合具体语境理解,是表示按位相乘还是位运算卷积)

    考虑我们在实现FFT的时候干了什么?

    我们搞出了一个点值表示法,然后直接对点值进行求乘积,再把它转回去

    那么这对位运算是否成立呢?

    即我们希望有(FWT(C)=FWT(A) imes FWT(B)),此时可以直接按位相乘

    直觉告诉我们是有的,但是是什么样的呢?

    or

    你要问我是怎么找的我也无从说起,十分诡异

    想要了解更多的可以参考:https://www.cnblogs.com/ACMLCZH/p/8022502.html

    这里直接给出结论及其正确性的证明

    对于一个有(2^n)项的(A_i)

    [FWT(A)=egin{cases}A & n=0 \ (FWT(A_0),FWT(A_0+A_1)) & n>0 end{cases} ]

    什么是(A_0)(A_1)

    我们把(A)想象成一个(2^n)维的向量,那么前(2^{n-1})项是(A_0),后(2^{n-1})项是(A_1)

    这之后的括号表示的意义同样(在之后会阐述(FWT(A))的本质)

    如果只是感性理解这玩意还是比较容易的:对于(A_0)中的每一项,它们的二进制表示最高位均是(0),当它们与(A_1)的某一项进行了(or)卷积的话,那么最高位就变成了1,对后面的(2^{n-1})项产生了贡献

    那么问题就变成了如何证明这时的(FWT(C)=FWT(A) imes FWT(B))

    我们考虑使用数学归纳法证明,设此时(A,B,C)均为(2^n)

    (n=0),结论显然成立

    (n>0),假设结论在(n=k)时成立

    那么此时就应该有(注意这里用(*)表示(or)卷积)

    [FWT(A_0*B_0)=FWT(A_0) imes FWT(B_0)\ FWT(A_1*B_0)=FWT(A_1) imes FWT(B_0)\ FWT(A_0*B_1)=FWT(A_0) imes FWT(B_1)\ FWT(A_1*B_1)=FWT(A_1) imes FWT(B_1) ]

    (理解一下,在括号里面的表示位运算卷积,外面的则为按位相乘)

    同时显然(FWT(A+B)=FWT(A)+FWT(B))(+表示按位相加)

    首先是等式左边

    [C=(A_0*B_0,A_1*B_0+A_0*B_1+A_1*B_1)\ ]

    这也很好理解,和上面(FWT)的公式一样理解

    [egin{aligned} FWT(C)=&(FWT(A_0*B_0),FWT(A_0*B_0)+FWT(A_0*B_1+A_1*B_0+A_1*B_1))\ =&(FWT(A_0*B_0),FWT(A_0*B_0)+FWT(A_0*B_1)+FWT(A_1*B_0)+FWT(A_1*B_1))\ =&(FWT(A_0) imes FWT(B_0),FWT(A_0) imes FWT(B_0)+FWT(A_1) imes FWT(B_1)+FWT(A_1) imes FWT(B_0)+FWT(A_1) imes FWT(B_1)\ end{aligned} ]

    再来看等式左边

    [egin{aligned} FWT(A) imes FWT(B)=&(FWT(A_0),FWT(A_0+A_1)) imes(FWT(B_0),FWT(B_0+B_1))\ =&(FWT(A_0) imes FWT(B_0),(FWT(A_0)+FWT(A_1)) imes(FWT(B_0)+FWT(B_1)))\ =&(FWT(A_0) imes FWT(B_0),FWT(A_0) imes FWT(B_0)+FWT(A_0) imes FWT(B_1)+FWT(A_1) imes FWT(B_0)+FWT(A_1) imes FWT(B_1)) end{aligned} ]

    证毕

    好像结束了?让我们来看一下(FWT(A))的本质

    回到最开始的式子(c_i=sum_{j|k=i}a_jb_k),假设我们记(C_i=sum _{jsubseteq i}c_j),则

    [egin{aligned} C_i=&sum_{j|ksubseteq i}a_ib_j\ =&sum_{jsubseteq i,ksubseteq i}a_ib_j\ =&A_jB_k end{aligned} ]

    我们发现(A,B,C)就起到了上文中(FWT)的作用,于是当我们要求某个集合的子集中的元素和时,除了(O(3^n))的枚举子集以外,还有(O(n2^n))的求出数组(A)(FWT(A))的方法,(A_i)(i)的子集和就是(FWT(A)_i)

    and

    证明和上面是类似的,在这里直接给出结论

    [FWT(A)=egin{cases}A & n=0 \ (FWT(A_0+A_1),FWT(A_1)) & n>0 end{cases} ]

    以及我们注意到,当求出(FWT(A))时,我们有(FWT(A)_i=sum_{isubseteq j}A_j),它和上面的(or)卷积在其他除了FWT的题目中也有很好的应用

    xor

    依然是只给出结论,证明留给读者

    [FWT(A)=egin{cases}A & n=0 \ (FWT(A_0)+FWT(A_1),FWT(A_0)-FWT(A_1)) & n>0 end{cases} ]

    IDFT相关

    我们实现了FFT中类似于DFT的一步,那么IDFT呢?

    直接按照上面的逆运算做回去就可以了,直接上结论:
    对于(or)卷积

    [IFWT(A)=egin{cases}A& n=0 \ (IFWT(A_0),IFWT(A_1)-IFWT(A_0)) & n>0 end{cases} ]

    对于(and)卷积

    [IFWT(A)=egin{cases}A& n=0 \ (IFWT(A_0)-IFWT(A_1),IFWT(A_1))& n>0end{cases} ]

    对于(xor)卷积

    [IFWT(A)=egin{cases} A& n=0\ (frac{IFWT(A_0)+IFWT(A_1)}{2},frac{IFWT(A_0)-IFWT(A_1)}{2}end{cases} ]

    code

    模板题:luogu 4717

    #include<iostream>
    #include<string.h>
    #include<string>
    #include<stdio.h>
    #include<algorithm>
    #include<math.h>
    #include<vector>
    #include<queue>
    #include<map>
    using namespace std;
    typedef long long ll;
    const ll maxd=998244353,inv2=499122177;
    const double pi=acos(-1.0); 
    int n,N;
    ll a[1001000],b[1001000],c[1001000];
    
    int read()
    {
        int x=0,f=1;char ch=getchar();
        while ((ch<'0') || (ch>'9')) {if (ch=='-') f=-1;ch=getchar();}
        while ((ch>='0') && (ch<='9')) {x=x*10+(ch-'0');ch=getchar();}
        return x*f;
    }
    
    void fwt_or(int lim,ll *a,int typ)
    {
        int mid;
        for (mid=1;mid<lim;mid<<=1)
        {
            int len=(mid<<1),sta,j;
            for (sta=0;sta<lim;sta+=len)
            {
                for (j=0;j<mid;j++)
                {
                    if (typ==1) a[sta+mid+j]=(a[sta+j+mid]+a[sta+j])%maxd;
                    else a[sta+mid+j]=(a[sta+j+mid]+maxd-a[sta+j])%maxd;
                }
            }
        }
    }
    
    void fwt_and(int lim,ll *a,int typ)
    {
        int mid;
        for (mid=1;mid<lim;mid<<=1)
        {
            int len=(mid<<1),sta,j;
            for (sta=0;sta<lim;sta+=len)
            {
                for (j=0;j<mid;j++)
                {
                    if (typ==1) a[sta+j]=(a[sta+j]+a[sta+mid+j])%maxd;
                    else a[sta+j]=(a[sta+j]+maxd-a[sta+j+mid])%maxd;
                }
            }
        }
    }
    
    void fwt_xor(int lim,ll *a,int typ)
    {
        int mid;
        for (mid=1;mid<lim;mid<<=1)
        {
            int len=(mid<<1),sta,j;
            for (sta=0;sta<lim;sta+=len)
            {
                for (j=0;j<mid;j++)
                {
                    int x=a[sta+j],y=a[sta+j+mid];
                    a[sta+j]=(x+y)%maxd;
                    a[sta+j+mid]=(x+maxd-y)%maxd;
                    if (typ==-1) {a[sta+j]=(a[sta+j]*inv2)%maxd;a[sta+j+mid]=(a[sta+j+mid]*inv2)%maxd;}
                }
            }
        }
    }
    
    int main()
    {
        n=read();N=1<<n;int i;
        for (i=0;i<N;i++) a[i]=read();
        for (i=0;i<N;i++) b[i]=read();
        fwt_or(N,a,1);fwt_or(N,b,1);
        for (i=0;i<N;i++) c[i]=(a[i]*b[i])%maxd;
        fwt_or(N,a,-1);fwt_or(N,b,-1);fwt_or(N,c,-1);
        for (i=0;i<N;i++) printf("%lld ",c[i]);printf("
    ");
        fwt_and(N,a,1);fwt_and(N,b,1);
        for (i=0;i<N;i++) c[i]=(a[i]*b[i])%maxd;
        fwt_and(N,a,-1);fwt_and(N,b,-1);fwt_and(N,c,-1);
        for (i=0;i<N;i++) printf("%lld ",c[i]);printf("
    ");
        fwt_xor(N,a,1);fwt_xor(N,b,1);
        for (i=0;i<N;i++) c[i]=(a[i]*b[i])%maxd;
        fwt_xor(N,a,-1);fwt_xor(N,b,-1);fwt_xor(N,c,-1);
        for (i=0;i<N;i++) printf("%lld ",c[i]);printf("
    ");
        return 0;
    }
    
  • 相关阅读:
    hdu 1269 迷宫城堡(强联通分量,基础)
    hdu 2102 A计划(BFS,基础)
    python 变量命名规范
    rpm常用选项
    memcached
    session共享
    Nginx高级使用
    nginx 反向代理
    Nginx基本使用
    github 建立博客
  • 原文地址:https://www.cnblogs.com/encodetalker/p/10783732.html
Copyright © 2011-2022 走看看