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;
}