hdu4609 3-idiots FFT

     跟n条边，求任选三条边可以组成三角形的概率。求出n条边组成三角形的方案数，再除以C（n,3）就可以了，所以这题转化成n条边求可组成的三角形数。

首先用一个数组num[i]记录长度为i的边出现了几次，之后求num的卷积，即可得到n条边中任取两条边的和的长度各出现了多少次，拿第一组样例来说

a[]={1 3 3 4}转化成num={0,1,0,2,1}，{0,1,0,2,1}*{0,1,0,2,1}得到{0,0,1,0,4,2,4,4,1},乘积即两边之和出现的次数，但是我们在选取的时候一条边不能用两次，而a+b和b+a只记一次，所以要对这个乘积处理一下，即num[a[i]+a[i]]--, num[i]/=2, 处理后的结果{0,0,0,0,2,1,1,2,0}就是C（n,2）得到的长度和出现的次数.之后枚举a[i]，假设a[i]为三边的最大值，cnt=sum[max]-sum[a[i]]就是任选两边长度符合的种数，其中sum[i]为num[]的前缀和，而这个cnt中还有非法的情况:

cnt-=(ll)(n-i-1)*i;//b,c一个大于a[i],一个小于
cnt-=(ll)(n-1);//b,c中一个是a[i]
cnt-=(ll)(n-i-1)*(n-i-2)/2;//b,c都大于a[i]

最后累加所有的cnt就是最后的答案。num[]的卷积可以用快速傅里叶变换（FFT）以nlogn的复杂度完成，没做过fft的同学可以参考
http://blog.csdn.net/night_raven/article/details/20546435

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <cstring>
using namespace std;
typedef long long ll;
const double PI = acos(-1.0);
struct comp
{
double r,i;
comp(double rt=0,double it=0)
{
r=rt;
i=it;
}
comp operator +(const comp& b)
{
return comp(r+b.r,i+b.i);
}
comp operator -(const comp &b)
{
return comp(r-b.r,i-b.i);
}
comp operator *(const comp &b)
{
return comp(r*b.r-i*b.i,r*b.i+i*b.r);
}
};

void change(comp y[],int len)//二进制转置--雷德算法
{
int i,j,k;
for(i = 1, j = len/2;i < len-1;i++)
{
if(i < j)swap(y[i],y[j]);
k = len/2;
while( j >= k)
{
j -= k;
k /= 2;
}
if(j < k)j += k;
}
}

void fft(comp y[],int len,int on)
/* on=1 DFT 把一个多项式的系数向量转化为点集表示；
on=-1,IDFT 把一个点集转化成多项式的系数向量*/
{
change(y,len);
for(int h = 2;h <= len;h <<= 1)
{
comp wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
for(int j = 0;j < len;j += h)
{
comp w(1,0);
for(int k = j;k < j+h/2;k++)
{
comp u = y[k];
comp t = w*y[k+h/2];
y[k] = u+t;
y[k+h/2] = u-t;
w = w*wn;
}
}
}
if(on == -1)
for(int i = 0;i < len;i++)
y[i].r /= len;
}

void conv(comp f[],int len)//求f的卷积
{
fft(f,len,1);
for (int i=0; i<len; i++)
f[i]=f[i]*f[i];
fft(f,len,-1);
}
int n,m;
const int maxn=404000;
comp x1[maxn];
int a[maxn];
ll num[maxn],sum[maxn];
int p,q,k,t;
int len1,len;
int main()
{
// freopen("in.txt","r",stdin);
int tt;
scanf("%d",&tt);
while(tt--)
{
memset(num,0,sizeof num);
scanf("%d",&n);
for (int i=0; i<n; i++)
scanf("%d",&a[i]),num[a[i]]++;
sort(a,a+n);
len1=a[n-1]+1;
len=1;
while (len<len1*2) len<<=1;
for (int i=0; i<len1; i++)
x1[i]=comp(num[i],0);
for (int i=len1; i<len; i++)
x1[i]=comp(0,0);

conv(x1,len);
for (int i=0; i<len; i++)
num[i]=(ll)(x1[i].r+0.5);//四舍五入

len=2*a[n-1];

//去重
for (int i=0; i<n; i++)
num[a[i]+a[i]]--;//相同的组合减点一个
for (int i=1; i<=len; i++)
num[i]>>=1;//选择有序，除以2

memset(sum,0,sizeof sum);
for (int i=1; i<=len; i++)
sum[i]=sum[i-1]+num[i];

ll cnt=0;
ll tot=(ll)n*(n-1)*(n-2)/6; //C(n,3);
ll ans=0;
for (int i=0; i<n; i++)
{
//假设a[i]为三角形的最长边，枚举

cnt=sum[len]-sum[a[i]];//b+c>a[i]的取法

cnt-=(ll)(n-i-1)*i;//b,c一个大于a[i],一个小于
cnt-=(ll)(n-1);//b,c中一个是a[i]
cnt-=(ll)(n-i-1)*(n-i-2)/2;//b,c都大于a[i]
ans+=cnt;
}
double out=(double)ans/(double)tot;
printf("%.7lf\n",out);

}
return 0;
}