HDU-4609(FFT/NTT)

Posted 2021-03-20 chasedeath

HDU-4609(FFT/NTT)

题意: 给出n个木棒，现从中不重复地选出3根来，求能拼出三角形的概率。

计算合法概率容易出现重复，所以建议计算不合法方案数

枚举选出的最大边是哪条，然后考虑剩下两条边之和小于等于它

两条边之和为(x)的方案数可以(FFT/NTT)得到，是一个简单的构造

即(f(x)=sum x^{length_i})，求出(f(x)^2)，就能得到和的方案数，但是会重复，包括自己和自己算，一对算两次

处理一下前缀和即可

#include<bits/stdc++.h>
using namespace std;

#define reg register
typedef long long ll;
#define rep(i,a,b) for(reg int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(reg int i=a,i##end=b;i>=i##end;--i)

template <class T> inline void cmin(T &a,T b){ ((a>b)&&(a=b)); } 
template <class T> inline void cmax(T &a,T b){ ((a<b)&&(a=b)); } 

char IO;
int rd(){
    int s=0,f=0;
    while(!isdigit(IO=getchar())) if(IO=='-') f=1;
    do s=(s<<1)+(s<<3)+(IO^'0');
    while(isdigit(IO=getchar()));
    return f?-s:s;
}

const double PI=acos(-1);
const int N=262150,P=998244353;
const int g=3;

int n,m;
int b[100010];
int rev[N];
typedef complex <double> Cp;
Cp d[N];
ll a[N];

void FFT(int n,Cp *a,int f){
    rep(i,0,n-1) if(rev[i]>i) swap(a[i],a[rev[i]]);
    for(reg int i=1;i<n;i<<=1) {
        Cp w(cos(PI/i),f*sin(PI/i));
        for(reg int l=0;l<n;l+=i*2) {
            Cp e(1,0);
            for(reg int j=l;j<l+i;j++,e=e*w) {
                Cp t=a[j+i]*e;
                a[j+i]=a[j]-t;
                a[j]=a[j]+t;
            }
        }
    }
}

int main(){
    rep(kase,1,rd()) {
        n=rd();
        int ma=0;
        rep(i,1,n) b[i]=rd(),d[b[i]]+=Cp(1,0),ma=max(ma,b[i]);
        int R=1,c=0;
        while(R<=ma*2) R<<=1,c++;
        rep(i,1,R) rev[i]=(rev[i>>1]>>1)|((i&1)<<(c-1));
        FFT(R,d,1);
        rep(i,0,R-1) d[i]=d[i]*d[i];
        FFT(R,d,-1);
        rep(i,0,R-1) a[i]=(ll)(d[i].real()/R+0.5);
        rep(i,1,n) a[b[i]*2]--;
        rep(i,1,R) a[i]/=2;    //处理重复情况
        rep(i,1,R) a[i]+=a[i-1];
        ll ans=0,sum=0;
        rep(i,1,n) {
            ans+=a[b[i]]; // 计算不合法方案数
            sum+=1ll*(i-1)*(i-2)/2;
        }
        rep(i,0,R+3) a[i]=0,d[i]=Cp(0,0);
        rep(i,0,ma+3) b[i]=0;
        ans=sum-ans;
        printf("%.7lf
",1.0*ans/sum);
    }
}