HDU 4609 3-idiots FFT+计数问题

题目大意：

就是现在给出n个长度不超过100000的树枝 的长度, n <= 100000, 长度都是整数, 求在这n根树枝当中任意取其中3根, 能组成三角形的概率

大致思路：

这个题还是很好的一道计数问题, 首先由于n <= 100000不难想到考虑O(nlogn)的算法, 如果用一个多项式A(x), 其x^k的系数ak表示长度为k的树枝的数量, 那么A(x)的平方当中x^k的系数就是从这些树枝中取两根可重复的和为k的排列数, 稍作处理即可得到A[i]数组其中i表示两根树枝的长度和为i, A[i]表示这样的组合数量是A[i]种

那么对于一个三角形, 考虑其最大的边, 首先将n根树枝排序, 然后从小到大一次考虑其作为组成的三角形中的最长边即可

对于第i根树枝最为最长边, 另外两条边的和需要大于第i根树枝的长度, 也就是之前处理出的A[i]~A[maxLength]的和, 然后去掉用了与第i根之后树枝的组合, 在减去用了第i根树枝作为宁外两个的组合的组数即可, 细节可参见代码注释

代码如下：

Result  :  Accepted     Memory  :  14112 KB     Time  :  764 ms

/*
 * Author: Gatevin
 * Created Time:  2015/7/15 10:33:43
 * File Name: HDU4609.cpp
 */
#include<iostream>
#include<sstream>
#include<fstream>
#include<vector>
#include<list>
#include<deque>
#include<queue>
#include<stack>
#include<map>
#include<set>
#include<bitset>
#include<algorithm>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cctype>
#include<cmath>
#include<ctime>
#include<iomanip>
using namespace std;
const double eps(1e-8);
typedef long long lint;

const double PI = acos(-1.0);

struct Complex
{
    double real, image;
    Complex(double _real, double _image)
    {
        real = _real;
        image = _image;
    }
    Complex(){}
};

Complex operator + (const Complex &c1, const Complex &c2)
{
    return Complex(c1.real + c2.real, c1.image + c2.image);
}

Complex operator - (const Complex &c1, const Complex &c2)
{
    return Complex(c1.real - c2.real, c1.image - c2.image);
}

Complex operator * (const Complex &c1, const Complex &c2)
{
    return Complex(c1.real*c2.real - c1.image*c2.image, c1.real*c2.image + c1.image*c2.real);
}

int rev(int id, int len)
{
    int ret = 0;
    for(int i = 0; (1 << i) < len; i++)
    {
        ret <<= 1;
        if(id & (1 << i)) ret |= 1;
    }
    return ret;
}

Complex tA[264000];
void FFT(Complex* a, int len, int DFT)
{
    for(int i = 0; i < len; i++)
        tA[rev(i, len)] = a[i];
    for(int s = 1; (1 << s) <= len; s++)
    {
        int m = (1 << s);
        Complex wm = Complex(cos(DFT*2*PI/m), sin(DFT*2*PI/m));
        for(int k = 0; k < len; k += m)
        {
            Complex w = Complex(1, 0);
            for(int j = 0; j < (m >> 1); j++)
            {
                Complex u = w*tA[k + j + (m >> 1)];
                Complex t = tA[k + j];
                tA[k + j] = t + u;
                tA[k + j + (m >> 1)] = t - u;
                w = w*wm;
            }
        }
    }
    if(DFT == -1) for(int i = 0; i < len; i++) tA[i].real /= len, tA[i].image /= len;
    for(int i = 0; i < len; i++) a[i] = tA[i];
    return;
}

int branch[100010];
int num[200010];
Complex a[264000];//(1 << 17 = 131072, 1 << 18 = 262144)
lint A[200010];
lint sumA[200010];//表示A[i]的前缀和

int main()
{
    int T, n;
    scanf("%d", &T);
    while(T--)
    {
        scanf("%d", &n);
        int maxBranch = 0;
        for(int i = 0; i < n; i++)
        {
            scanf("%d", branch + i);
            maxBranch = max(maxBranch, branch[i]);
        }
        memset(num, 0, sizeof(num));
        for(int i = 0; i < n; i++)
            num[branch[i]]++;
        for(int i = 0; i <= maxBranch; i++)
            a[i] = Complex(num[i], 0);
        int len = 1;
        while(len <= maxBranch) len <<= 1;
        len <<= 1;
        for(int i = maxBranch + 1; i < len; i++)
            a[i] = Complex(0, 0);
        FFT(a, len, 1);
        for(int i = 0; i < len; i++)
            a[i] = a[i]*a[i];
        FFT(a, len, -1);
        for(int i = 0; i <= 2*maxBranch; i++)
            A[i] = (lint)(a[i].real + 0.5);
        for(int i = 0; i <= 2*maxBranch; i += 2)
            A[i] -= num[i >> 1];
        for(int i = 0; i <= 2*maxBranch; i++)
            A[i] /= 2;
        //到现在为止A[i]表示的是去两根不同的branch的长度和为i的组合种数
        sumA[0] = 0;
        for(int i = 1; i <= 2*maxBranch; i++)
            sumA[i] = sumA[i - 1] + A[i];
        lint ans = 0;
        sort(branch, branch + n);
        for(int i = 1; i <= n; i++)//以第i根作为边最长的
        {
            lint tmp = sumA[2*maxBranch] - sumA[branch[i]];//另外两条边长度和要大于branch[i]
            tmp -= (lint)(n - i)*(i - 1);//另外两条一条比branch[i]长, 一条不比它长
            tmp -= (lint)(n - i)*(n - i - 1) / 2;//两条都比他长
            tmp -= n - 1;//另外两条的组合中包括它自己的组合
            ans += tmp;
        }
        double p = ans*6./n/(n - 1)/(n - 2);
        printf("%.7f\n", p);
    }
    return 0;
}