hdu4906：3-idiots【FFT】

题意

给出n个正整数，求从中任选3个数(下标不重复)，求以这三个数为长度的边能构成三角形的概率。

分析与解

可以用FFT求出两个数的和（注意去重），然后枚举最长边，利用前缀和寻找比他短的和，总和即为不能形成三角形的数量。用总数减去这个数就是构成三角形的数量。

#include <bits/stdc++.h>
using namespace std;

int n, a[(1<<18)+5];
int rev[(1<<18)+5];
int lg(int n)
{
    int i = 0;
    for (int j = 1; j < n; j <<= 1, i++);
    return i;
}
const double pi = 3.1415926536;

typedef complex<double> Complex;
Complex A[(1<<18)+5];

void fft(Complex a[], int n, int flag)
{
    int lgn = lg(n); rev[0] = 0;
    for (int i = 0; i < n; i++)
        rev[i] = (rev[i>>1]>>1)|((i&1)<<(lgn-1));
    for (int i = 0; i < n; i++)
        A[rev[i]] = a[i];
    Complex u, t;
    for (int k = 2; k <= n; k <<= 1) {
        Complex dw = Complex(cos(2*pi/k), sin(flag*2*pi/k));
        for (int i = 0; i < n; i += k) {
            Complex w = 1;
            for (int j = 0; j < k>>1; j++) {
                t = w*A[i+j+(k>>1)], u = A[i+j];
                A[i+j] = u+t, A[i+j+(k>>1)] = u-t;
                w *= dw;
            }
        }
    }
    if (flag == 1)
        for (int i = 0; i < n; i++)
            a[i] = A[i];
    else
        for (int i = 0; i < n; i++)
            a[i] = A[i]/(double)n; 
}

Complex d[(1<<18)+5];
long long sum[(1<<18)+5];

int main()
{
    int T;
    scanf("%d", &T);
    while (T--) {
        int mx = INT_MIN;
        scanf("%d", &n);
        for (int i = 1; i <= n; i++) scanf("%d", &a[i]), d[a[i]] = Complex(d[a[i]].real()+1, 0), mx = max(mx, a[i]);
        int len = 1;
        for (; len <= mx*2; len <<= 1);
        fft(d, len, 1);
        for (int i = 0; i < len; i++) d[i] = d[i]*d[i];
        fft(d, len, -1);
        for (int i = 1; i <= n; i++) d[a[i]+a[i]] = Complex(d[a[i]+a[i]].real()-1, 0);
        memset(sum, 0, sizeof sum);
        sum[0] = 0; d[0] = Complex(0, 0);
        for (int i = 1; i < len; i++) sum[i] = sum[i-1] + (int)(d[i].real()+0.001)/2, d[i] = Complex(0, 0);
        long long cnt = 0;
        for (int i = 1; i <= n; i++) cnt += sum[a[i]];
        printf("%.7lf\n", (double)((long long)n*(n-1)*(n-2)/6-cnt)/((long long)n*(n-1)*(n-2)/6));
    }
    return 0;
}