• bzoj 3527: [Zjoi2014]力 快速傅里叶变换


    题意:

    给出n个数qi,给出Fj的定义如下: 
     
    令Ei=Fi/qi,求Ei

      fft的那一堆东西还是背不到啊。。。这次写虽说完全自己写的,但是还是在参见了以前fft程序的情况下调了很久,主要在如下几点写错:1、非递归中内层数组调用中下表忘掉加k 2、每次转换乘的那个数是cos(...)+isin(...),不要记混了,且里面是(a/b*2*PI) 3、pp[]没有每次清零这一些逗B错误。

    #include<iostream>
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<cmath>
    using namespace std;
    #define MAXN 1100000
    #define PI 3.1415926535897932384
    int m,n;
    struct Complex
    {
            double x,y;
            Complex(){};
            Complex(double x,double y=0):x(x),y(y){};
            Complex operator +(Complex a)
            {
                    return Complex(x+a.x,y+a.y);
            }
            Complex operator -(Complex a)
            {
                    return Complex(x-a.x,y-a.y);
            }
            Complex operator *(Complex a)
            {
                    return Complex(x*a.x-y*a.y,x*a.y+y*a.x);
            }
    };
    Complex ww[MAXN][2];
    void dft(Complex g[],int len,bool d)
    {
            Complex t;
            for (int i=1;i<len;i<<=1)
            {
                    for (int j=0;j<len;j+=(i<<1))
                    {
                            for (int k=0;k<i;k++)
                            {
                                    t=g[k+j];
                                    g[k+j]=g[k+j]+g[k+j+i]*ww[k * (n/(i<<1))][d];
                                    g[k+j+i]=t-g[k+j+i]*ww[k * (n/(i<<1))][d];
                            }
                    }
            }
    }
    int pp[MAXN];
    Complex g1[MAXN],g2[MAXN];
    void fft(double s1[],double s2[],int m,double res[])
    {
            int i,j,k,x;
            n=m;
            while (n != (n&(-n)))n-=n&(-n);
            n<<=2;
            memset(pp,0,sizeof(pp));
            for (i=0;i<n;i++)
            {
                    for (x=1,j=n>>1;j;j>>=1,x<<=1)
                    {
                            pp[i]+=((i&j)!=0)*x;
                    }
            }
            for (i=0;i<n;i++)g1[pp[i]]=s1[i];
            for (i=0;i<n;i++)g2[pp[i]]=s2[i];
            for (i=0;i<=n;i++)
            {
                    ww[i][0]=Complex(cos(2*PI*i/n),-sin(2*PI*i/n));
                    ww[i][1]=Complex(ww[i][0].x,-ww[i][0].y);
            }
            dft(g1,n,0);
            dft(g2,n,0);
            for (i=0;i<n;i++)g2[i]=g1[i]*g2[i];
            for (i=0;i<n;i++)g1[pp[i]]=g2[i];
            dft(g1,n,1);
            for (i=n;i>=0;i--)
                    res[i]=g1[i].x/n;
    }
    double q1[MAXN],q2[MAXN],a[MAXN];
    double r1[MAXN],r2[MAXN];
    double f[MAXN];
    double s1[MAXN],s2[MAXN];
    int main()
    {
            //freopen("input.txt","r",stdin);
            //freopen("output.txt","w",stdout);
            int n;
            int i,j,k;
            scanf("%d",&n);
            int x,y,z;
            for (i=0;i<n;i++)
                    scanf("%lf",q1+i),q2[n-i-1]=q1[i];
            for (i=1;i<n;i++)
                    a[i]=1.0/i/i;
            fft(q1,a,n,r1);
            fft(q2,a,n,r2);
            for (i=0;i<n;i++)
            {
                    f[i]+=r1[i];
                    f[i]-=r2[n-i-1];
            }
            for (i=0;i<n;i++)
            {
                    printf("%.5lf ",f[i]);
            }
    
    }
    by mhy12345(http://www.cnblogs.com/mhy12345/) 未经允许请勿转载

    本博客已停用,新博客地址:http://mhy12345.xyz

  • 相关阅读:
    nginx防止盗链
    Nginx防盗链详细设置
    [bzoj2127]happiness
    [bzoj2400]Optimal Marks
    [bzoj1738]发抖的牛
    [bzoj1741]穿越小行星群
    [bzoj3123]森林
    [bzoj2588]Count on a tree
    [bzoj3144]切糕
    [bzoj1787]紧急集合
  • 原文地址:https://www.cnblogs.com/mhy12345/p/4174996.html
Copyright © 2020-2023  润新知