首页 > 技术文章 > 快速傅里叶变换应用之二 hdu 4609 3-idiots

holyprince 2014-03-12 17:50 原文

快速傅里叶变化有不同的应用场景,hdu4609就比较有意思。题目要求是给n个线段,随机从中选取三个,组成三角形的概率。

初始实在没发现这个怎么和FFT联系起来,后来看了下别人的题解才突然想起来:组合计数问题可以用多项式的卷积来解决。于是将给的数据进行卷积相乘,利用FFT即可求出三角形任意两条线段组合的可能数目。

然后遍历初始数据,将其作为最长边(这里一开始也没想明白,其实就是只要最长边大于短边之和,其他两个不等式也自然可以满足)。那么理论上说比它长的所有两边组合可能都可以。当然在这里要考虑三种特殊情况:(即在两边组合数目中减去这些情况)

1.这两个边有可能一个边比最长边长,一个边小于最长边

2.其中一个边就是要选的这个边

3.两个边其实都比最长边长,这种情况要除以二

 

PS:G++使用的是longlong类型,C++是_int64,好久没写忘记了。

longlong在代码中间乘的运算也要加上,否则还是会出错。

#include <iostream>
#include <cmath>
#include <algorithm> //spell!
#include <string.h>
#define MAXN 400040 
#define PI acos(-1.0)
using namespace std;

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[MAXN];
 


void bitrev(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 *in,int n,int flag)
{
    int i,j,k;
    complex u,t;
    bitrev(in,n);
    for(int i=2;i<=n;i=i*2)
    {
        complex wn(cos((2*PI*flag)/i),sin((2*PI*flag)/i));//初始化单位复根
        for(j=0;j<n;j=j+i)
        {
            complex w(1,0);
            for(k=j;k<j+i/2;k++)
            {
                u=in[k];
                t=w*in[k+i/2];
                in[k]=u+t;
                in[k+i/2]=u-t;
                w=w*wn;
            }
        }
    }
    if(flag==-1)
        for(int i=0;i<n;i++)
            in[i].r=in[i].r/n;
}


int a[100003];
long long res[MAXN]; 
long long sum[MAXN];
long long num[MAXN];
int main() {
    int T;
    scanf("%d",&T);
    while(T--)
    {
        int n,i;
        scanf("%d",&n);
        memset(res,0,sizeof(res));
        memset(sum,0,sizeof(sum));
        memset(num,0,sizeof(num));
        for(int j=0;j<n;j++)
        {
            scanf("%d",&a[j]);
            num[a[j]]++;
        }
        sort(a,a+n);
        for(i = 0;i <=a[n-1];i++)
        {
            x1[i].r=num[i];
            x1[i].i=0;
        }
        int expandn=1;
        while(expandn<2*(a[n-1]+1))expandn=expandn*2;

        for(i = a[n-1]+1;i<expandn;i++)
        {
            x1[i].r=0;
            x1[i].i=0;
        }
        fft(x1,expandn,1);
        for(i=0;i<expandn;i++)
            x1[i]=x1[i]*x1[i];
        fft(x1,expandn,-1);
        for(i=0;i<expandn;i++)
        {
            res[i]=(long long)(x1[i].r+0.5);
        }
        //去除本身
        for(i=0;i<n;i++)
            res[a[i]+a[i]]--;
        //变为组合
        for(i=0;i<expandn;i++)
            res[i]=res[i]/2;
        //求出两边之和为i的所有可能
        //expandn=(a[n-1]+1)*2;
        sum[0]=res[0];
        for(i=1;i<expandn;i++)
            sum[i]=res[i]+sum[i-1];
        long long ans=0;
        for(i=0;i<n;i++)
        {
            ans+=sum[expandn-1]-sum[a[i]];//比长度为a[i]大的所有可能
            //去除一个大于a[i],一个小于a[i]
            ans=ans-(long long)(n-1-i)*i;
            //去除一个取自己
            ans=ans-(n-1);
            //去除取两个都大
            ans=ans-(long long)(n-1-i)*(n-2-i)/2;
        }
        long long  all = (long long)n*(n-1)*(n-2)/6;
        printf("%.7lf\n",(double)ans/all);
    }
}
hdu 4609

 

推荐阅读