HDU4609 NTT||FFT

先用NTT求出两两组合的方案数
然后就能o(n) 求出能组成三角形的方案数
NTT:

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=4e5+10;
const ll mod=( 1ll << 47 ) * 7 * 4451 + 1 ;
const ll g=3;
ll mul( ll x, ll y )
{
    return ( x * y - ( long long ) ( x / ( long double ) mod * y + 1e-3 ) * mod + mod ) % mod ;
}
ll power ( ll a, ll b )
{
    ll res = 1, tmp = a ;
    while ( b )
    {
        if ( b & 1 ) res = mul ( res, tmp ) ;
        tmp = mul ( tmp, tmp ) ;
        b >>= 1 ;
    }
    return res ;
}
void DFT ( ll y[], int n, bool rev )
{
    for ( int i = 1, j, t, k ; i < n ; ++ i )
    {
        for ( k = n >> 1, t = i, j = 0 ; k ; k >>= 1, t >>= 1 )
        {
            j = j << 1 | t & 1 ;
        }
        if ( i < j ) swap ( y[i], y[j] ) ;
    }
    for ( int s = 2, ds = 1 ; s <= n ; ds = s, s <<= 1 )
    {
        ll wn = power ( g, ( mod - 1 ) / s ) ;
        if ( !rev ) wn = power ( wn, mod - 2 ) ;
        for ( int k = 0 ; k < n ; k += s )
        {
            ll w = 1, t ;
            for ( int i = k ; i < k + ds ; ++ i, w = mul ( w, wn ) )
            {
                y[i + ds] = ( y[i] - ( t = mul ( y[i + ds], w ) ) + mod ) % mod ;
                y[i] = ( y[i] + t ) % mod ;
            }
        }
    }
}
void NTT ( ll x1[], ll x2[], int n )
{
    DFT ( x1, n, 1 ) ;
    DFT ( x2, n, 1 ) ;
    for ( int i = 0 ; i < n ; ++ i ) x1[i] = mul ( x1[i], x2[i] ) ;
    DFT ( x1, n, 0 ) ;
    ll vn = power ( n, mod - 2 ) ;
    for ( int i = 0 ; i < n ; ++ i ) x1[i] = mul ( x1[i], vn ) ;
}
ll s[N],vis[N],sum[N];
ll x1[N],x2[N];
int main()
{
    int t;
    scanf("%d",&t);
    while(t--)
    {
        ll n;
        scanf("%lld",&n);
        ll mx;
        memset(vis,0,sizeof(vis));
        for(int i=0;i<n;i++)
        {
            scanf("%lld",&s[i]);
            vis[s[i]]++;
        }
        sort(s,s+n);
        mx=s[n-1]+1;
        int len=1;
        while(len<2*mx) len<<=1;
        for(int i=0; i<mx; i++)
            x2[i]=x1[i]=vis[i];
        for(int i=mx; i<len; i++)
            x2[i]=x1[i]=0;
        NTT(x1,x2,len);
        len=2*mx;
        ll res=0;
        for(int i=0;i<n;i++)
            x1[s[i]+s[i]]--;
        for(int i=1;i<=len;i++)
            x1[i]>>=1;
        for(int i=1;i<=len;i++)
           sum[i]=sum[i-1]+x1[i];
        for(int i=0;i<n;i++)
        {
            res+=sum[len]-sum[s[i]];
            res-=(n-1);
            res-=(n-1-i)*i;
            res-=(n-1-i)*(n-i-2)/2;
        }
        printf("%.7lf\n",6.0*res/(n*(n-1)*(n-2)));
    }
}