HDU 4609 快速傅里叶变换

题意：给出10^5个边长，问任取三个组成三角形的概率。

这题没遇到还不知道什么是傅里叶变换FFT，当时怎么想时间复杂度都不能降到n*log(n)。实际是这样的，把每个记录长度个数的数组存在num数组中。

直接摘大神的吧 网址http://www.cnblogs.com/kuangbin/archive/2013/07/24/3210565.html

学会了FFT。

这题就很容易了。

其实题目是给了n条线段。问随机取三个，可以组成三角形的概率。

其实就是要求n条线段，选3条组成三角形的选法有多少种。

首先题目给了a数组，

如样例一：

4

1 3 3 4

把这个数组转化成num数组，num[i]表示长度为i的有num[i]条。

样例一就是

num = {0   1   0    2    1}

代表长度0的有0根，长度为1的有1根，长度为2的有0根，长度为3的有两根，长度为4的有1根。

使用FFT解决的问题就是num数组和num数组卷积。

num数组和num数组卷积的解决，其实就是从{1 3 3 4}取一个数，从{1 3 3 4}再取一个数，他们的和每个值各有多少个

例如{0 1 0 2 1}*{0 1 0 2 1} 卷积的结果应该是{0 0  1  0  4  2  4  4  1 }

长度为n的数组和长度为m的数组卷积，结果是长度为n+m-1的数组。

{0 1 0 2 1}*{0 1 0 2 1} 卷积的结果应该是{0 0  1  0  4  2  4  4  1 }。

这个结果的意义如下：

从{1 3 3 4}取一个数，从{1 3 3 4}再取一个数

取两个数和为 2 的取法是一种：1+1

           和为 4 的取法有四种：1+3， 1+3  ，3+1 ，3+1

           和为 5 的取法有两种：1+4 ，4+1；

           和为 6的取法有四种：3+3,3+3,3+3,3+3,3+3

           和为 7 的取法有四种： 3+4,3+4,4+3,4+3

           和为 8 的取法有 一种：4+4

利用FFT可以快速求取循环卷积，具体求解过程不解释了，就是DFT和FFT的基本理论了。

总之FFT就是快速求到了num和num卷积的结果。只要长度满足>=n+m+1.那么就可以用循环卷积得到线性卷积了。

弄完FFT得到一个num数组，这个数组的含义在上面解释过了。

#include <iostream>#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>using namespace std;#define N 400005#define pi acos(-1.0) // PI值struct complex{    double r,i;    complex(double real=0.0,double image=0.0)    {        r=real;        i=image;    }    // 以下为三种虚数运算的定义    complex operator + (const complex o)    {        return complex(r+o.r,i+o.i);    }    complex operator - (const complex o)    {        return complex(r-o.r,i-o.i);    }    complex operator * (const complex o)    {        return complex(r*o.r-i*o.i,r*o.i+i*o.r);    }} x1[N],x2[N];void brc(complex *y,int l) // 二进制平摊反转置换 O(logn){    register int i,j,k;    for(i=1,j=l/2; i<l-1; i++)    {        if(i<j)    swap(y[i],y[j]); // 交换互为下标反转的元素        // i<j保证只交换一次        k=l/2;        while(j>=k) // 由最高位检索，遇1变0，遇0变1，跳出        {            j-=k;            k/=2;        }        if(j<k)    j+=k;    }}void fft(complex *y,int l,double on) // FFT O(nlogn)// 其中on==1时为DFT，on==-1为IDFT{    register int h,i,j,k;    complex u,t;    brc(y,l); // 调用反转置换    for(h=2; h<=l; h<<=1) // 控制层数    {        // 初始化单位复根        complex wn(cos(on*2*pi/h),sin(on*2*pi/h));        for(j=0; j<l; j+=h) // 控制起始下标        {            complex w(1,0); // 初始化螺旋因子            for(k=j; k<j+h/2; k++) // 配对            {                u=y[k];                t=w*y[k+h/2];                y[k]=u+t;                y[k+h/2]=u-t;                w=w*wn; // 更新螺旋因子            } // 据说上面的操作叫蝴蝶操作…        }    }    if(on==-1)    for(i=0; i<l; i++)    y[i].r/=l; // IDFT}long long num[N];long long sum[N];int a[N>>1];int main(){    int t,n,l,maxn;    scanf("%d",&t);    while(t--)    {        memset(num,0,sizeof(num));        maxn=0;        scanf("%d",&n);        for(int i=0; i<n; i++)        {            scanf("%d",&a[i]);            num[a[i]]++;        }        sort(a,a+n);        maxn=a[n-1];        l=1;        while(l<2*(maxn+1)) l<<=1;        for(int i=0; i<=maxn; i++)            x1[i].r=num[i],x1[i].i=0;        for(int i=maxn+1; i<=l; i++)            x1[i].r=x1[i].i=0;        fft(x1,l,1);        for(int i=0; i<=l; i++) x1[i]=x1[i]*x1[i];        fft(x1,l,-1);        for(int i=0; i<=2*maxn; i++) num[i]=(long long)(x1[i].r+0.5);        for(int i=0; i<n; i++) num[a[i]+a[i]]--;        for(int i=1; i<=2*maxn; i++) num[i]/=2;        sum[0]=0;        for(int i=1; i<=2*maxn; i++)            sum[i]=sum[i-1]+num[i];        long long cnt=0;        for(int i=0; i<n; i++) //枚举每条边并把这条边当成最大边        {            cnt+=sum[2*maxn]-sum[a[i]];            cnt-=n-1;    //减去包括这条边的            cnt-=(long long)i*(n-i-1);//减去有一条边比这条边大的情况            cnt-=(long long)(n-i-1)*(n-i-2)/2;//减去有两条都比该条边大的情况        }        long long tol=(long long)n*(n-1)*(n-2)/6;        printf("%.7f\n",double((double)cnt/(double)tol));    }    return 0;}