HDU 4609 3-idiots(FFT计数)

#include<bits/stdc++.h>
using namespace std;
#define pie acos(-1.0)
const int maxn = 4e5 + 10;
typedef long long ll;
struct Cp{double x,y;}A[maxn],B[maxn],C[maxn],w,w0,x;
Cp operator +(Cp x,Cp y){Cp z;z.x=x.x+y.x;z.y=x.y+y.y;return z;}
Cp operator -(Cp x,Cp y){Cp z;z.x=x.x-y.x;z.y=x.y-y.y;return z;}
Cp operator *(Cp x,Cp y){Cp z;z.x=x.x*y.x-x.y*y.y;z.y=x.x*y.y+x.y*y.x;return z;}
int n;
ll sum[maxn],num[maxn],a[maxn],ans,b[maxn];
void Rader(Cp A[],int len){
    int j=len>>1;
    for(int i=1;i<len-1;i++){
        if(i<j) swap(A[i],A[j]);
        int k=len>>1;
        while(j>=k){j-=k,k>>=1;}
        if(j<k)j+=k;
    }
}
void FFT(Cp A[],int len,int f){
    Rader(A,len);
    for(int i=1;i<len;i<<=1){
        w.x=cos(pie/i);w.y=sin(pie/i)*f;
        for(int j=0;j<len;j+=i<<1){
            w0.x=1.0;w0.y=0.0;
            for(int k=0;k<i;k++){
                x=A[j+k];
                A[j+k]=x+w0*A[j+k+i];
                A[j+k+i]=x-w0*A[j+k+i];
                w0=w0*w;
            }
        }
    }
    if(f<0) for(int i=0;i<len;i++) A[i].x/=len;
}
void solve()
{
    memset(num,0,sizeof(num));
    sort(a + 1,a+n+1); int len1=a[n]+1;
    int len=1; while(len < len1*2)len<<=1;
    memset(b,0,sizeof(b)); memset(sum,0,sizeof(sum));
    for(int i=1;i<=n;i++) b[a[i]]++;
    for(int i=0;i<len;i++) A[i].x=b[i],A[i].y=0.0;
    FFT(A,len,1);
    for(int i=0;i<len;i++) A[i]=A[i]*A[i];
    FFT(A,len,-1);
    for(int i=0;i<len;i++) num[i]=A[i].x+0.5;
    ll tot=0;
    for(int i=1;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++) sum[i]=num[i] + sum[i-1];
    ans=0;
    for(int i=1;i<=n;i++){
        ans+=sum[len] - sum[a[i]];
        ans-=(ll)(n-i)*(n-i-1)/2;
        ans-=(ll)(n-i)*(i-1);
        ans-=(n-1);
    }
    tot=(ll)(n-1)*(n-2)*n/6;
    printf("%.7lf\n",(double)ans/tot);
}
int main()
{
   int t;
   scanf("%d",&t);
   while(t--){
    scanf("%d",&n);
    for(int i=1;i<=n;i++) scanf("%lld",&a[i]);
    solve();
   }
}