链接: http://acm.hdu.edu.cn/showproblem.php?pid=4609
题意: 给定 N 个正整数, 表示 N 条线段的长度, 问任取 3 条, 可以构成三角形的概率为多少~
数据范围: N<=10^5 ~~
思路:设三边分别为 x, y, z (x<=y<=z) 枚举 z ,统计 x+y 大于 z 的数目 .
比赛时能想到的只有 O(n^2) 的算法,无力 AC~
赛后才知道有种东西叫 FFT ~
以下为官方解题报告:
/*
记录 A_i 为长度为 i 的树枝的数量，并让 A 对它本身做 FFT，得到任意选两个树枝能得到的各个和的数量。枚举第三边，
计算出所有两边之和大于第三条边的方案数，并把前两条边包含最长边的情况减掉就是答案。
*/
设长度为 a1,a2, ....an, 统计每种长度个数为 cnt[ ai ] ; 可表示多项式
多项式自乘后， 其指数即为 ai + aj 的值， 系数即为 方案数， 减去不合法的即为所求~
#include <cstdio>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <cstring>
using namespace std;
typedef __int64 LL;
const double pi=acos(-);
const int MAXN=3e5;
const double eps=1e-;
struct C
{
double r, i;
C (){}
C(double _r, double _i):r(_r),i(_i){}
inline C operator +(const C a)const{
return C(r+a.r, i+a.i);
}
inline C operator - (const C a)const {
return C(r-a.r, i-a.i);
}
inline C operator * (const C a)const{
return C(r*a.r-i*a.i, r*a.i+i*a.r);
}
}a[MAXN], b[MAXN];
int num[MAXN], cnt[MAXN];
LL res[MAXN], sum[MAXN]; void brc(C *y,int l) // 二进制平摊反转置换 O(logn)
{
register int i,j,k;
for(i=,j=l>>;i<l-;i++){
if(i<j) swap(y[i],y[j]); // 交换互为下标反转的元素
// i<j保证只交换一次
k=l>>;
while(j>=k) {// 由最高位检索，遇1变0，遇0变1，跳出
j-=k;
k>>=;
}
if(j<k) j+=k;
}
}
void FFT(C *y,int l,int on) // FFT O(nlogn)
// 其中on==1时为DFT，on==-1为IDFT
{
register int h,i,j,k;
C u,t;
brc(y,l); // 调用反转置换
for(h=;h<=l;h<<=) // 控制层数
{
// 初始化单位复根
C wn(cos(on**pi/h),sin(on**pi/h));
for(j=;j<l;j+=h) // 控制起始下标
{
C w(,); // 初始化螺旋因子
for(k=j;k<j+h/;k++) // 配对
{
u=y[k];
t=w*y[k+h/];
y[k]=u+t;
y[k+h/]=u-t;
w=w*wn; // 更新螺旋因子
} // 据说上面的操作叫蝴蝶操作…
}
}
if(on==-) for(i=;i<l;i++) y[i].r/=l; // IDFT
} int n, L , Max, T;
int main() {
scanf("%d", &T);
while (T--) {
Max = ;
memset(cnt, , sizeof(cnt));
scanf("%d", &n);
for (int i = ; i < n; ++i) {
scanf("%d", num + i);
cnt[num[i]]++;
Max = max(Max, num[i]);
}
++Max;L = ;
while (L < Max <<) {
L <<= ;
}
for (int i = ; i < Max; ++i) {
a[i] = C(cnt[i], );
} for (int i = Max; i < L; ++i) {
a[i] = C(, );
}
FFT(a, L, );
for (int i = ; i < L; ++i)
a[i] = a[i] * a[i]; // 多项式自乘
FFT(a, L, -); for (int i = ; i < L; ++i) {
res[i] = (LL) (a[i].r + 0.5);
} for (int i = ; i <= Max; ++i)
res[i << ] -= cnt[i];// 自己和自己乘， 即 枚举的是 x+x 的 for (int i = ; i < L; ++i)
res[i] >>= ; // x+y 与 y+x 算一次
for (int i = ; i < L; ++i) { //前缀和
sum[i] = sum[i - ] + res[i];
}
double tot = , den = 1.0*n * (n - ) * (n - )/; for (int i = ; i < n; ++i) {
tot += sum[num[i]] / den;// x+y<=z 的个数 }
double ans = - tot ;
printf("%.7f\n", ans);
} return ;
}