FFT_深入浅出解释fft

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<ctime>
#include<cmath>
#include<iostream>
#include<climits>
#include<string>
#include<deque>
#include<queue>
#include<vector>
#include<algorithm>
#include<map>
#include<stack>
#include<set>
using namespace std;
#define LL long long
#define pi (acos(-1.0))
 
 
struct complex
{
	double r,i;
	complex(double _r=0,double _i=0){ r=_r ;  i=_i;  }
	complex operator + (const complex &b) { return complex(r+b.r,i+b.i); }
	complex operator - (const complex &b) { return complex(r-b.r,i-b.i); }
	complex operator * (const complex &b) { return complex(r*b.r-i*b.i,r*b.i+i*b.r);}
};
 
 
void change(complex 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(complex y[],int len,int on)
{
	change(y,len);
	for(int h=2; h<=len; h<<=1)
	{
		complex wn(cos(-on*2*pi/h),sin(-on*2*pi/h));
		for(int j=0; j<len;j+=h)
		{
			complex w(1,0);
			for(int k=j;k<j+h/2;k++)
			{
				complex u=y[k];
				complex 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;
}
 
 
const int maxn=810010;
complex x1[maxn];
LL num[maxn];
int a[maxn];
 
 
void dofft(LL len1)
{
	int len=1;
	while(len < 2*len1) len<<=1;
	for(int i=0;i<len1;i++) x1[i]=complex(num[i],0);
	for(int i=len1;i<len;i++) x1[i]=complex(0,0);
	fft(x1,len,1);
	for(int i=0;i<len;i++)  x1[i]=x1[i]*x1[i];
	fft(x1,len,-1);
	for(int i=0;i<len;i++) num[i]=(LL)(x1[i].r+0.5);
}
 
 
int main()
{
#ifndef ONLINE_JUDGE
	freopen("in.txt","r",stdin);
	freopen("out.txt","w",stdout);
#endif
	int t; scanf("%d",&t);
	int n;
	while(t--)
	{
		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);
		dofft(a[n-1]+1);
		LL 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]/=2;
		for(int i=1;i<=len;i++) num[i]+=num[i-1];
		double ans=0.0;
		for(int i=0;i<n;i++)
		{
			ans+=num[len]-num[a[i]];
			ans-=(n-i-1)*(double)i;
			ans-=(double)(n-1);
			ans-=(n-i-1)*(double)(n-i-2)/2.0;
		}
		double all=(double)n*(n-1)*(n-2)/6.0;
		printf("%.7f\n",ans/all);
	}
	return 0;
}