hdu4609(fft)

题目链接: http://acm.hdu.edu.cn/showproblem.php?pid=4609

题意: 给出 n 根木棒求从中任取三根可以组成三角形的概率.

思路: fft

将 a 数组转化为 num 数组, 其中 num[i] 为长度为 i 的木棒的数目, 再做 num 和 num 的卷积, 这步可以用 fft 模板完成. 将卷积结果存入 num 中, 则 num[i] 为从 n 根木棒中任选两根且其长度和为 i 的组合数. 注意: 里面包含了自己和自己组合以及正反序的情况. 所以要先对 num 处理一下.

若两边的长和为 x (x > 0), 再有一条边长度小于 x 且大于 0, 那么这三条边一定可以组成三角形. 前面已经得到了 num 数组, 那么接下来可以枚举第三条边, 然后累计对答案的贡献. 为了避这里免重复计算, 先给 a 数组排序, 然后枚举 a[i] 作为当前三角形中的最长边, 累计对答案的贡献即可.

得到了能组成三角形的方案数 sol, 则概率为 sol / (n * (n - 1) * (n - 2) / 6).

代码:

#include <iostream>
#include <algorithm>
#include <math.h>
#include <string.h>
#include <stdio.h>
#define ll long long
using namespace std;

const double PI = acos(-1.0);

struct Complex{//复数结构体
    double x, y;//实部,虚部
    Complex(double _x = 0.0, double _y = 0.0){
        x = _x;
        y = _y;
    }
    Complex operator -(const Complex &b) const{
        return Complex(x - b.x, y - b.y);
    }
    Complex operator +(const Complex &b) const{
        return Complex(x + b.x, y + b.y);
    }
    Complex operator *(const Complex &b) const{
        return Complex(x * b.x - y * b.y, x * b.y + y * b.x);
    }
};

//进行FFT和IFFT反转变化
//位置i和(i二进制反转后位置)互换
void change(Complex y[], int len){//len必须为2的幂
    for(int i = 1, j = len / 2; i < len - 1; i++){
        if(i < j) swap(y[i], y[j]); //交换互为下标反转的元素,i<j保证只交换一次
        int k = len >> 1;
        while(j >= k){
            j -= k;
            k /= 2;
        }
        if(j < k) j += k;
    }
}

//做FFT,len必须为2的幂,on=1是DFT,on=-1是IDTF
void fft(Complex y[], int len, int on){
    change(y, len);//调用反转置换
    for(int h = 2; h <= len; h <<= 1){
        Complex wn(cos(-on * 2 * PI / h), sin(-on * 2 * PI / h));
        for(int j = 0; j < len; j += h){
            Complex w(1, 0);//初始化螺旋因子
            for(int k = j; k < j + h / 2; k++){//配对
                Complex u = y[k];
                Complex t = w * y[k + h / 2];
                y[k] = u + t;
                y[k + h / 2] = u - t;
                w = w * wn;//更新螺旋因子
            }
        }
    }
    if(on == -1){
        for(int i = 0; i < len; i++){
            y[i].x /= len;//IDTF
        }
    }
}

const int MAXN = 4e5 + 10;
ll sum[MAXN], num[MAXN];
Complex x1[MAXN];
int a[MAXN];

int main(void){
    int t, n;
    scanf("%d", &t);
    while(t--){
        scanf("%d", &n);
        memset(num, 0, sizeof(num));
        memset(sum, 0, sizeof(sum));
        for(int i = 0; i < n; i++){
            scanf("%d", &a[i]);
            num[a[i]]++;
        }
        sort(a, a + n);
        int len = 1;
        int len1 = a[n - 1] + 1;
        while(len < (len1 << 1)) len <<= 1;
        for(int i = 0; i < len1; i++){
            x1[i] = Complex(num[i], 0);
        }
        for(int i = len1; i < len; i++){
            x1[i] = Complex(0, 0);
        }
        fft(x1, len, 1);
        for(int i = 0; i < len; i++){
            x1[i] = x1[i] * x1[i];
        }
        fft(x1, len, -1);//IDFT(x1*x1)
        for(int i = 0; i < len; i++){
            num[i] = (ll)(x1[i].x + 0.5);//四舍五入
        }
        //此时的num[i]为从a中选中任意两根木棒的和等于i的组合数,包含了自己和自己组合的情况以及正反序的组合情况
        ll sol = 0;
        len = a[n - 1] << 1;
        for(int i = 0; i < n; i++){//减去自己和自己组合的
            num[a[i] + a[i]]--;
        }
        for(int i = 0; i <= len; i++){//得到的组合数是无序的,所以要除2
            num[i] /= 2;
        }
        for(int i = 1; i <= len; i++){
            sum[i] += sum[i - 1] + num[i];//sum[i]即从n根木棍中选出两根长度和大于等于i的方案数
        }
        for(int i = 0; i < n; i++){//当a[i]为三角形最长边时对答案的贡献
            sol += sum[len] - sum[a[i]];
            sol -= (ll)(n - i - 1) * i;//减去有一条边大于a[i]的情况
            sol -= (ll)(n - i - 1) * (n - i - 2) / 2;//减去两条边都大于a[i]的情况
            sol -= n - 1;//减去包括a[i]的情况
        }
        ll cnt = (ll)n * (n - 1) * (n - 2) / 6;
        double gel = sol * 1.0 / cnt;
        printf("%.7lf
", gel);
    }
    return 0;
}