HDU 4609 3-idiots（FFT学习）

附上kuangbin大神博客：讲解的非常清楚

分析：暴力枚举$a_k$，那么问题变成了找满足$a_i+a_j>a_k(i,j,k互不相等)$的$(a_i,a_j)$的对数。FFT求出所有可能的和的对数，然后去掉坐标相同的就可以了。

PS：目前就写了两道FFT，基本上都是卷积的应用，FFT神奇的地就是能在$nlogn$的时间内算出所有的

$$X(k)=\sum{a[i]*a[k-i]}$$

附上代码：

#include <bits/stdc++.h>
#define LL long long
#define FOR(i,x,y)  for(int i = x;i < y;++ i)
#define IFOR(i,x,y) for(int i = x;i > y;-- i)

using namespace std;

//FFT copy from kuangbin
const double pi = acos (-1.0);
// Complex  z = a + b * i
struct Complex {
    double a, b;

    Complex(double _a=0.0,double _b=0.0):a(_a),b(_b){}

    Complex operator + (const Complex &rhs) const {
        return Complex(a + rhs.a , b + rhs.b);
    }

    Complex operator - (const Complex &rhs) const {
        return Complex(a - rhs.a , b - rhs.b);
    }

    Complex operator * (const Complex &rhs) const {
        return Complex(a * rhs.a - b * rhs.b , a * rhs.b + b * rhs.a);
    }
};

//len = 2 ^ k
void change (Complex y[] , int len) {
    for (int i = 1 , j = len / 2 ; i < len -1 ; i ++) {
        if (i < j) swap(y[i] , y[j]);
        int k = len / 2;
        while (j >= k) {
            j -= k;
            k /= 2;
        }
        if(j < k) j += k;
    }
}

// FFT
// len = 2 ^ k
// on = 1  DFT    on = -1 IDFT
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].a /= len;
        }
    }
}

const int maxn = 100010;
const int maxm = 100010;

LL cnt[maxm<<2],sum[maxm<<2];
int a[maxn],n;

Complex x[maxm<<2];

void work(){
    memset(cnt,0,sizeof(cnt));
    FOR(i,0,n)  scanf("%d",&a[i]),++ cnt[a[i]];
    sort(a,a+n);
    int mx_len = (a[n-1] + 1)<<1;
    int len = 1;
    while(len < mx_len)    len <<= 1;
    FOR(i,0,len)    x[i] = Complex(cnt[i],0.0);
    FFT(x,len,1);
    FOR(i,0,len)    x[i] = x[i]*x[i];
    FFT(x,len,-1);
    FOR(i,0,len)    cnt[i] = (LL)(x[i].a+0.5);
    FOR(i,0,n)      -- cnt[a[i]<<1];
    FOR(i,0,len)    cnt[i] >>= 1;
    sum[0] = cnt[0];
    FOR(i,1,len)    sum[i] = sum[i-1]+cnt[i];
    LL ans = 0;
    FOR(i,0,n){
        ans += sum[len-1] - sum[a[i]];
        ans -= (LL)(n-i-1)*i;
        ans -= (LL)(n-1);
        ans -= (LL)(n-i-1)*(n-i-2)>>1;
    }
    LL tot = (LL)n*(n-1)*(n-2)/6;
    printf("%.7f\n",(double)ans/tot);
}

int main()
{
    //freopen("test.in","r",stdin);
    int T;  scanf("%d",&T);
    while(T--){
        scanf("%d",&n);
        work();
    }
    return 0;
}