HDU4609 3-idiots

传送门


最近比赛遇到了要用FFT的题，就趁机加强一下。

一句话题意：在(n)个数中随便拿出三个数，问能组成三角形的概率。（(1 leqslant n leqslant 10^ 5, 1 leqslant a_i leqslant10 ^ 5)）


看到(a_i)的限制，就能想到开一个桶(num[i])表示大小为(i)的数的个数。

那么(f = num ^ 2)，(f(i))就表示拿两个数，组成大小为(i)的数的方案数。但这里面有不合法的和重复的，所以要将(f(a_i+a_i))减(1)，以及将所有(f(i))除以(2)，因为(a_x+a_y)和(a_y+a_x)是同一种选法。

上述的(f(i))直接用FFT即可计算，接下来考虑怎么组成三角形。


我们可以枚举长边。
那么对于每一个长边(a_i)，可以和(sumlimits_{j>a_i} f(j))组成三角形，但这样算会多。

一是可能会有一个比(a_i)长的边，一个比(a_i)短的边，那么方案数就要减去((i-1)*(n-i)).
二是可能一个刚好也选(a_i)（注意，是选(a_i)本身，而选和(a_i)相同的数是合法的，而且两个短边不可能同时选到(a_i)，因为这种情况再刚求完(f)的时候已经去掉了），另一个随便选，因此要再减去(n-1).
最后一种就是选了两个都比(a_i)大的数，再减去(C_{n - i}^2).


总的来说这题难点不在于卷积，而是后面组合数的推导，要不重不漏的考虑到所有情况。

#include<cstdio>
#include<iostream>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<cstdlib>
#include<cctype>
#include<vector>
#include<queue>
#include<assert.h>
#include<ctime>
using namespace std;
#define enter puts("") 
#define space putchar(' ')
#define Mem(a, x) memset(a, x, sizeof(a))
#define In inline
#define forE(i, x, y) for(int i = head[x], y; ~i && (y = e[i].to); i = e[i].nxt)
typedef long long ll;
typedef double db;
const int INF = 0x3f3f3f3f;
const db eps = 1e-8;
const int maxn = 3e5 + 5;
const db PI = acos(-1);
In ll read()
{
	ll ans = 0;
	char ch = getchar(), las = ' ';
	while(!isdigit(ch)) las = ch, ch = getchar();
	while(isdigit(ch)) ans = (ans << 1) + (ans << 3) + ch - '0', ch = getchar();
	if(las == '-') ans = -ans;
	return ans;
}
In void write(ll x)
{
	if(x < 0) x = -x, putchar('-');
	if(x >= 10) write(x / 10);
	putchar(x % 10 + '0');
}

int n, a[maxn];

int len = 1, lim = 0, rev[maxn];
struct Comp
{
	db x, y;
	In Comp operator + (const Comp& oth)const {return (Comp){x + oth.x, y + oth.y};}
	In Comp operator - (const Comp& oth)const {return (Comp){x - oth.x, y - oth.y};}
	In Comp operator * (const Comp& oth)const {return (Comp){x * oth.x - y * oth.y, x * oth.y + y * oth.x};}
	friend In void swap(Comp &a, Comp& b) {swap(a.x, b.x), swap(a.y, b.y);}
}A[maxn];
In void fft(Comp* a, int flg)
{
	for(int i = 0; i < len; ++i) if(i < rev[i]) swap(a[i], a[rev[i]]);
	for(int i = 1; i < len; i <<= 1)
	{
		Comp omg = (Comp){cos(PI / i), sin(PI / i) * flg};
		for(int j = 0; j < len; j += (i << 1))
		{
			Comp o = (Comp){1, 0};
			for(int k = 0; k < i; ++k, o = o * omg)
			{
				Comp tp1 = a[k + j], tp2 = o * a[k + j + i];
				a[k + j] = tp1 + tp2, a[k + j + i] = tp1 - tp2;
			}
		}
	}
}

ll num[maxn], sum[maxn];
In ll solve()
{
	sort(a + 1, a + n + 1);
	for(int i = 1; i <= n; ++i) num[a[i]]++;
	int N = a[n] + a[n];
	while(len <= N) len <<= 1, ++lim;
	for(int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (lim - 1));
	for(int i = 0; i <= a[n]; ++i) A[i] = (Comp){num[i], 0};
	for(int i = a[n] + 1; i < len; ++i) A[i] = (Comp){0, 0};
	fft(A, 1);
	for(int i = 0; i < len; ++i) A[i] = A[i] * A[i];
	fft(A, -1);
	for(int i = 0; i <= N; ++i) num[i] = A[i].x / len + 0.5;
	for(int i = 1; i <= n; ++i) num[a[i] + a[i]]--;
	for(int i = 1; i <= N; ++i) num[i] /= 2, sum[i] = sum[i - 1] + num[i];
	ll ret = -1LL * n * (n - 1);
	for(int i = 1; i <= n; ++i)
	{
		ret += sum[N] - sum[a[i]];
		ret -= 1LL * (i - 1) * (n - i);
		ret -= (1LL * (n - i) * (n - i - 1) >> 1);
	}
	return ret;
}

In void init()
{
	Mem(num, 0);
	len = 1, lim = 0;
}

int main()
{
	int T = read();
	while(T--)
	{
		init();
		n = read();
		for(int i = 1; i <= n; ++i) a[i] = read();
		printf("%.7lf
", 1.0 * solve() / (n - 2) / (n - 1) / n * 6);
	}
	return 0;
}