BZOJ2194. 快速傅立叶之二

清疚 2021-10-24 02:22 257阅读 0赞

[传送门][Link 1]

求 $c\[k\]=\\sum\_\{i=k\}^\{n\}(a\[i\]\*b\[i-k\])$

为了搞 $FFT$，考虑把 $a,b$ 的下标变成相加等于 $k$ 这样的形式

设 $d$ 为 $b$ 翻转后的数组，即 $d\[i\]=b\[n-1-i\]$，或者说 $b\[i\]=d\[n-1-i\]$

原式 $c\[k\]=\\sum\_\{i=k\}^\{n\}(a\[i\]\*b\[i-k\])$ 可变成 $c\[k\]=\\sum\_\{i=k\}^\{n\}(a\[i\]\*b\[n-1-(n-1-i+k)\])$

把 $d\[i\]=b\[n-1-i\]$ 代入得 $c\[k\]=\\sum\_\{i=k\}^\{n\}(a\[i\]\*d\[n-1-i+k\])$

发现 $a,d$ 的下标相加恒为 $n-1+k$ 所以可以 $FFT$

为了更好理解，不妨设 $e\[n-1+i\]=c\[i\]$，则 $e\[n-1+k\]=\\sum\_\{i=k\}^\{n\}(a\[i\]\*d\[n-1-i+k\])$

此时 $a,d$ 下标相加恰好等于 $e$ 的下标，就是裸的 $FFT$

最后要求输出 $c\[i\]$ ，其实就是输出 $e\[n-1+i\]$

#include<iostream>
    #include<cstdio>
    #include<algorithm>
    #include<cstring>
    #include<cmath>
    using namespace std;
    typedef long long ll;
    inline 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<<1)+(x<<3)+(ch^48); ch=getchar(); }
        return x*f;
    }
    const int N=4e5+7;
    const double pi=acos(-1.0);
    struct CP {
        double x,y;
        CP (double xx=0,double yy=0) { x=xx,y=yy; }
        inline CP operator + (const CP &tmp) const {
            return CP(x+tmp.x,y+tmp.y);
        }
        inline CP operator - (const CP &tmp) const {
            return CP(x-tmp.x,y-tmp.y);
        }
        inline CP operator * (const CP &tmp) const {
            return CP(x*tmp.x-y*tmp.y,x*tmp.y+y*tmp.x);
        }
    }A[N],D[N];
    int n,p[N];
    void FFT(CP *A,int len,int type)
    {
        for(int i=0;i<len;i++) if(i<p[i]) swap(A[i],A[p[i]]);
        for(int mid=1;mid<len;mid<<=1)
        {
            CP wn(cos(pi/mid),type*sin(pi/mid));
            for(int R=mid<<1,j=0;j<len;j+=R)
            {
                CP w(1,0);
                for(int k=0;k<mid;k++,w=w*wn)
                {
                    CP x=A[j+k],y=w*A[j+mid+k];
                    A[j+k]=x+y;
                    A[j+mid+k]=x-y;
                }
            }
        }
    }
    int main()
    {
        n=read();
        for(int i=0;i<n;i++) A[i].x=read(),D[n-1-i].x=read();
        int len=1,tot=0;
        while(len<=2*(n-1)) len<<=1,tot++;
        for(int i=0;i<len;i++) p[i]=(p[i>>1]>>1) | ((i&1)<<(tot-1));
        FFT(A,len,1); FFT(D,len,1);
        for(int i=0;i<=len;i++) A[i]=A[i]*D[i];
        FFT(A,len,-1);
        for(int i=0;i<n;i++) printf("%d\n",int(A[n-1+i].x/len+0.5));
        return 0;
    }

转载于:https://www.cnblogs.com/LLTYYC/p/11254677.html

[Link 1]: https://www.lydsy.com/JudgeOnline/problem.php?id=2194