洛谷p3803 FFT入门

最新推荐文章于 2022-12-22 15:00:56 发布

weixin_30615767

最新推荐文章于 2022-12-22 15:00:56 发布

阅读量90

点赞数

CC 4.0 BY-SA版权

原文链接：http://www.cnblogs.com/kekxy/p/7526182.html

本文深入浅出地讲解了快速傅立叶变换（FFT）的基本原理，并通过洛谷p3803题目实例展示了如何利用FFT进行多项式卷积运算，实现从O(n^2)到O(nlogn)的时间复杂度优化。

洛谷p3803 FFT入门

ps：花了我一天的时间弄懂fft的原理，感觉fft的折半很神奇！
大致谈一谈FFT的基本原理：
　对于两个多项式的卷积，可以O(n^2)求出来（妥妥的暴力）
　显然一个多项式可以用a0+a1X+a2X^2+a3X^3+a4X^4……表示。
　也可以用(x1,y1),(x2,y2),(x3,y3),(x4,y4)的点集来表示。
　用点值表示有一个好处：两个多项式的卷积可以直接取相同的x值，y值相乘得到。
　那么，怎么转化为点值表示呢？
　直接代进去？显然也是O(n^2)，没用……
　设\(A(X)=a0+a1X+a2X^2+a3X^3+a4X^4……\)
　\(A[0](X)=a0+a2X+a4X^2……\)
　\(A[1](X)=a1+a3X+a5X^2……\)
　那么A(X)=A0(X^2)+X*A1(X^2);
　然而，X^2仍会有N种不同的取值。
　引入一个复数根可以很好的解决这个问题：设i=\(\sqrt{-1}\)
　\(cosθ+i∗sinθ\)有很多很好的性质：
　例如\((cosθ+i∗sinθ)^k=coskθ+i*sinkθ\)
　记\(Wn=cos\frac{2π}{n}+i∗sin\frac{2π}{n}\)
　\(A(W{n}^{k})=A[0](((W{n}^{k})^2)+W{n}^{k}*A[1]((W{n}^{k})^2)\)
　因为\[(W{n}^{k})^2=W{n}^{2k}=cos\frac{2π*2k}{n}+i∗sin\frac{2π*2k}{n}=cos\frac{2πk}{\frac{n}{2}}+i∗sin\frac{2πk}{\frac{n}{2}}=(W{\frac{n}{2}^{k}})\]
　所以\[A(W{n}^{k})=A[0]( (W{n}^{k})^2)+W{n}^{k}*A[1]((W{n}^{k})^2)=A[0]W{\frac{n}{2}^{k}}+W{n}^{k}*A[1]W{\frac{n}{2}}^{k}\]
　根据三角函数的周期性：2π为一个周期。
　\(W{n}^{k}=W{n}^{k　mod　n}\)
　那么\(A(W{n}^{k})=A[0]W{\frac{n}{2}}^{k}+W{n}^{k}*A[1]W{\frac{n}{2}^{k}}\)可以拆成两部分,设\(t<\frac{n}{2}\)
　\(A(W{n}^{t})=A[0]W{\frac{n}{2}}^{t}+W{n}^{t}*A[1]W{\frac{n}{2}^{t}}\)
　\(A(W{n}^{t+\frac{n}{2}})=A[0]W{\frac{n}{2}}^{t}-W{n}^{t}*A[1]W{\frac{n}{2}^{t}}\)
　因为\(W{n}^{t+\frac{n}{2}}=W{n}^{t}*W{n}^{\frac{n}{2}}\)
　\(W{n}^{\frac{n}{2}}=W{2}=cos\frac{2π}{2}+i*sin\frac{2π}{2}=cosπ+i*sinπ=-1\)
　证毕，啦啦啦！
　一直折半下去算就可以将时间优化到O(nlogn)
　 SouthEast
　上图摘自算法导论。算法导论下载，提取码qcfs
　一题洛谷的裸题p3803
【代码】

#include<cstdio>
#include<algorithm>
#include<cmath>

using namespace std;

typedef long long ll;

const int N=1e6+10;
const double pi=acos(-1.0);
int n,m,nn,k;
struct Complex {
    double real,image;   //real+i*image
    Complex(){}
    Complex(double _real,double _image)
    {
        real=_real; image=_image;
    }
    friend Complex operator + (Complex A,Complex B) { return (Complex(A.real+B.real,A.image+B.image)); }
    friend Complex operator - (Complex A,Complex B) { return (Complex(A.real-B.real,A.image-B.image)); }
    friend Complex operator * (Complex A,Complex B) { return (Complex(A.real*B.real-A.image*B.image,A.real*B.image+A.image*B.real)); }
} A[N<<2],B[N<<2],C[N<<2];
int rev[N<<2];

void FFT(Complex *A,int n,int DFT)
{
    for(int i=0;i<n;i++) if(i<rev[i]) swap(A[i],A[rev[i]]);
    for(int s=1;(1<<s)<=n;s++)
    {
        int mi=(1<<s);
        Complex wn=Complex(cos(2*pi/mi),sin(2*pi/mi)*DFT);
        for(int t=0;t<n;t+=mi)
        {
            Complex w=Complex(1,0);
            for(int j=0;j<(mi>>1);j++)
            {
                Complex u=A[t+j],v=w*A[t+j+(mi>>1)];
                A[t+j]=u+v; A[t+j+(mi>>1)]=u-v;
                w=w*wn;
            }
        }
    }
    if(DFT==-1) for(int i=0;i<n;i++) A[i].real/=n,A[i].image/=n;
}

int main()
{
    scanf("%d%d",&n,&m);
    for(int i=0;i<=n;i++) { int xx; scanf("%d",&xx); A[i]=Complex(xx,0); }
    for(int i=0;i<=m;i++) { int xx; scanf("%d",&xx); B[i]=Complex(xx,0); }
    m+=n; k=0;
    for(nn=1;nn<=m;nn<<=1) k++;
    for(int i=0;i<nn;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(k-1));
    FFT(A,nn,1); FFT(B,nn,1);
    for(int i=0;i<=nn;i++) C[i]=A[i]*B[i];
    FFT(C,nn,-1);
    for(int i=0;i<=m;i++) printf("%d ",(int)(C[i].real+0.5));
    printf("\n");
    return 0;
}

ps：这篇总结写得真心累……

转载于:https://www.cnblogs.com/kekxy/p/7526182.html