【bzoj3513】[MUTC2013]idiots【FFT/NTT】

题目传送门
题意：给定n个长度分别为ai$a_i$的木棒，问随机选择3个木棒能够拼成三角形的概率。
题解：正面可能很不好考虑，所以我们可以反向考虑，考虑用总方案数减去无法组成三角形的方案总数。我们先构造出木棒长度的生成函数A。则选两根木棒，使得它们长度之和为x的方案总数为A2$A^2$的次数为x的项的系数再减掉重复的。这个A2$A^2$可以FFT求出。设a[i]为长度为i的木棍的数量之和，sum[i]为长度≥i的木棍数量之和，f[i]为选两根木棍长度为i的方案总数，maxl为选两根木棍的最大长度之和，B=A2$B=A^2$。则显然不能形成三角形的方案总数为∑maxli=1f[i]sum[i]$\sum _{i=1}^{maxl}f[i]sum[i]$。
如何求f[i]？
f[i]=i≡0(mod2)?B[i]−a[i2]:B[i]2$f[i]=\frac{i\equiv 0(mod2)?B[i]-a[\frac{i}{2}]:B[i]}{2}$
i模2为0时会出现选两根长度相同的情况，共a[i/2]$a[i/2]$种长度相同的情况，所以要减去。由于B=A2$B=A^2$，所以B$B$的每项的系数都算了两次(比如说B次数为n的项，A的1和n-1乘起来以及A的n-1和1乘起来，这两种重复了)，所以要除以二。
最后用总方案数减去不行的方案数，再求个概率，就解决了！
代码

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
using namespace std;
typedef long long ll;
const int N=270005;
const double pi=3.141592653589793;
int t,n,m,tmp,x,rev[N];
ll tot,ans,sum[N];
struct complex{
    double x,y;
    complex(){
        x=y=0;
    }
    complex(double x,double y):x(x),y(y){}
    friend complex operator + (const complex &a,const complex &b){
        return complex(a.x+b.x,a.y+b.y);
    }
    friend complex operator - (const complex &a,const complex &b){
        return complex(a.x-b.x,a.y-b.y);
    }
    friend complex operator * (const complex &a,const complex &b){
        return complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);
    }
}a[N];
void fft(complex *a,int dft){
    for(int i=0;i<n;i++){
        if(i<rev[i]){
            swap(a[i],a[rev[i]]);
        }
    }
    for(int i=1;i<n;i<<=1){
        complex wn=complex(cos(pi/i),dft*sin(pi/i));
        for(int j=0;j<n;j+=i<<1){
            complex w=complex(1,0),x,y;
            for(int k=j;k<j+i;k++,w=w*wn){
                x=a[k];
                y=w*a[k+i];
                a[k]=x+y;
                a[k+i]=x-y;
            }
        }
    }
    if(dft==-1){
        for(int i=0;i<n;i++){
            a[i].x/=n;
        }
    }
}
int main(){
    scanf("%d",&t);
    while(t--){
        tmp=0;
        memset(a,0,sizeof(a));
        memset(sum,0,sizeof(sum));
        scanf("%d",&m);
        for(int i=1;i<=m;i++){
            scanf("%d",&x);
            tmp=max(tmp,x);
            a[x].x++;
        }
        for(int i=tmp;i>=1;i--){
            sum[i]=sum[i+1]+a[i].x;
        }
        for(n=1;n<=2*tmp+1;n<<=1);
        for(int i=0;i<n;i++){
            rev[i]=(rev[i>>1]>>1)|((i&1)*(n>>1));
        }
        fft(a,1);
        for(int i=0;i<n;i++){
            a[i]=a[i]*a[i];
        }
        fft(a,-1);
        tot=ans=1LL*m*(m-1)*(m-2)/6;
        for(int i=1;i<n;i++){
            ll tmp=a[i].x+0.5;
            if(i%2==0){
                tmp-=sum[i/2]-sum[i/2+1];
            }
            tmp>>=1;
            ans-=tmp*sum[i];
        }
        printf("%.7lf\n",1.0*ans/tot);
    }
    return 0;
}