HDU_4609_3-idiots

来源：互联网发布：mac模拟器windows 编辑：程序博客网时间：2024/05/29 15:26

题型：数论

题意：n条线段，取3个组成三角形的概率

分析：

gg bin对于这题的blog写的最棒了~

http://www.cnblogs.com/kuangbin/archive/2013/07/24/3210565.html

我just贴代码。

代码：

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>#include<algorithm>#define LL long long#define mt(a,b) memset(a,b,sizeof(a))using namespace std;const int M = 1e5+100;struct Complex {///复数结构体    double x,y;///实部和虚部 x+yi    Complex(double _x=0,double _y=0) {        x=_x;        y=_y;    }    friend Complex operator -(const Complex &a,const Complex &b) {        return Complex(a.x-b.x,a.y-b.y);    }    friend Complex operator +(const Complex &a,const Complex &b) {        return Complex(a.x+b.x,a.y+b.y);    }    friend Complex operator *(const Complex &a,const Complex &b) {        return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);    }};class FFT {///快速傅里叶变换    Complex u,t;    void change(Complex y[],int len) { ///len必须是2的幂        for(int i=1,j=len>>1,k; i<len-1; i++) {            if(i<j) swap(y[i],y[j]);            k=len>>1;            while(j>=k) {                j-=k;                k>>=1;            }            if(j<k) j+=k;        }    }public:    void fft(Complex y[],int len,int on) { ///on==1时是DFT，on==-1时是IDFT        change(y,len);        double tmp=-on*2*acos(-1.0);        for(int h=2; h<=len; h<<=1) {            Complex wn(cos(tmp/h),sin(tmp/h));            for(int j=0; j<len; j+=h) {                Complex w(1,0);                int h2=h>>1;                for(int k=j; k<j+h2; k++) {                    u=y[k];                    t=w*y[k+h2];                    y[k]=u+t;                    y[k+h2]=u-t;                    w=w*wn;                }            }        }        if(on==-1) {            for(int i=0; i<len; i++) {                y[i].x/=len;            }        }    }} g;Complex x[M*4];int a[M*2];LL flag[M*4];LL sum[M*4];int main() {    int _;    int n;    scanf("%d",&_);    while(_--) {        mt(flag,0);        scanf("%d",&n);        for(int i=0; i<n; i++) {            scanf("%d",&a[i]);            flag[a[i]]++;        }        sort(a,a+n);        int len = 1;        while(len<(a[n-1]+1)*2) len <<= 1;        for(int i=0; i<len; i++) {            x[i] = Complex(flag[i],0);        }        g.fft(x,len,1);        for(int i=0; i<len; i++) {            x[i] = x[i] * x[i];        }        g.fft(x,len,-1);        for(int i=0; i<len; i++) {            flag[i] = (LL)(x[i].x+0.5);        }        len = a[n-1]*2;        for(int i=0; i<n; i++) {            flag[a[i]+a[i]] -- ;        }        ///注意这边是1~len        for(int i=1; i<=len; i++) {            flag[i] /= 2;        }        sum[0] = 0;        for(int i=0; i<len; i++) {            sum[i+1] = sum[i] + flag[i+1];        }        LL cnt = 0;        for(int i = 0; i < n; i++) {            cnt += sum[len]-sum[a[i]];            cnt -= (LL)(n-1-i)*i;            cnt -= (n-1);            cnt -= (LL)(n-1-i)*(n-i-2)/2;        }//            printf("%d\n",cnt);        LL all = (LL)n*(n-1)*(n-2)/6;        printf("%.7f\n",cnt*1.0/all);    }    return 0;}/**241 3 3 442 3 3 4*/

0 0