HDU4609 计数问题+FFT

Posted 2020-10-10

题目大意：

给出n(1e5)条线段(长度为整数，<=1e5)，求任取3条线段能组成一个三角形的概率。

用cnt[i]记录长度为i的线段的个数，通过卷积可以计算出两条线段长度之和为i的方案数sum[i]：先用FFT计算出cnt[i]的卷积sum[i]，为取两条线段长度和为i的排列数（包括自己和自己），去掉自己和自己的方案数，再对所有sum[i]除以2即为所求方案数。

之后对所有线段a[i]有大到小排列，考虑第i条线段是三角形最长边的情况（长度相同则将编号大的视为更长，就没有长度相同的情况了。实际计算时无影响）。首先是sum[i + 1] + sum[i + 2] + ... + sum[len - 1]，即其他两边之和要大于他的方案数。然后去掉：1.另两边有一条比他大，有一条比他小。 2.两条都比他长。 3.有一条是他自己。这三种情况。

要注意有些地方爆int。还有len的计算要注意边界（一直wa，wa了半天）。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;

const double pi = acos(-1);

struct Complex{
    double real, imag;
    Complex (){}
    Complex (double _real, double _imag){
        real = _real;
        imag = _imag;
    }
};

Complex operator + (Complex x, Complex y){
    return Complex(x.real + y.real, x.imag + y.imag);
}

Complex operator - (Complex x, Complex y){
    return Complex(x.real - y.real, x.imag - y.imag);
}

Complex operator * (Complex x, Complex y){
    return Complex(x.real * y.real - x.imag * y.imag, x.real * y.imag + x.imag * y.real);
}

Complex operator / (Complex x, double y){
    return Complex(x.real / y, x.imag / y);
}

int reverse(int x, int len){
    int t = 0;
    for (int i = 1; i < len; i <<= 1){
        t <<= 1;
        if (x & i) t |= 1;
    }
    return t;
}

Complex A[400010];
void FFT(Complex *a, int n, int DFT){
    for (int i = 0; i < n; ++i) A[reverse(i, n)] = a[i];
    for (int i = 2; i <= n; i <<= 1){
        Complex wn = Complex(cos(2 * pi / i), DFT * sin(2 * pi / i));
        for (int j = 0; j < n; j += i){
            Complex w = Complex(1, 0);
            for (int k = 0; k < (i >> 1); ++k){
                Complex x = A[j + k];
                Complex y = w * A[j + k + (i >> 1)];
                A[j + k] = x + y;
                A[j + k + (i >> 1)] = x - y;
                w = w * wn;
            }
        }
    }
    if (DFT == -1) for (int i = 0; i < n; ++i) A[i] = A[i] / n;
    for (int i = 0; i < n; ++i) a[i] = A[i];
}

int T, n;
int a[400010];
int cnt[400010];
Complex B[400010];
long long sum[400010];

int main(){

    scanf("%d", &T);
    while (T--){
        scanf("%d", &n);
        int maxL = 0;
        memset(cnt, 0, sizeof(cnt));
        memset(sum, 0, sizeof(sum));
        for (int i = 0; i < n; ++i){
            scanf("%d", a + i);
            if (a[i] > maxL) maxL = a[i];
            ++cnt[a[i]];
        }
        int len = 1;
        while (len <= maxL * 2) len <<= 1;
        for (int i = 0; i <= maxL; ++i) B[i] = Complex(cnt[i], 0);
        for (int i = maxL + 1; i < len; ++i) B[i] = Complex(0, 0);
        FFT(B, len, 1);
        for (int i = 0; i < len; ++i) B[i] = B[i] * B[i];
        FFT(B, len, -1);
        for (int i = 0; i < len; ++i) sum[i] = (long long)(B[i].real + 0.5);
        for (int i = 0; i < n; ++i) --sum[a[i] + a[i]];
        for (int i = 0; i < len; ++i) sum[i] >>= 1;
        for (int i = len - 2; i >= 0; --i) sum[i] += sum[i + 1];
        sort(a, a + n);
        long long ans = 0;
        for (int i = 0; i < n; ++i){
            long long tmp = sum[a[i] + 1];
            tmp -= ((long long)n - i - 1) * i;
            tmp -= ((long long)n - i - 1) * (n - i - 2) / 2LL;
            tmp -= n - 1;
            ans += tmp;
        }
        double Ans = (double)ans * 6.0 / n / (n - 1) / (n - 2);
        printf("%.7f\n", Ans);
    }

    return 0;
}