HDU 4609 3-idiots
题意:
给出n条边,问选择3个边能构成三角形的概率
题解:
边权不大
考虑FFT的组合意义,得到的f[i]的值,是所有的j+k=i的a[j]*b[k]
如果把下标当做边权的话,FFT一次,就得到了选择两个边,能组成i的方案数。
三角形充要条件,最小的两个边之和大于最大边。
但是最大边可能之前在较小的用过了。
所以,考虑不合法的情况,最大边c>=a+b
这样,c一定就是一个>=a+b的后缀和了。
C(n,3)减去不合法的即可。
对于任意一个合法的,一定没有被减去,任意一个不合法的一定减去了。
正确。
注意:
1.多组数据清空,,用到的多项式是m>2*mx,buc的上界是mx
2.buc*buc可能爆int
3.FFT注意精度误差,绝对值可能在(0~0.5)之间,不能直接取整,采用四舍五入
a=floor(b+0.5)或者,a=round(b)
#include<cstdio> #include<cstdlib> #include<algorithm> #include<cstring> #include<iostream> #include<cmath> #define reg register int #define il inline #define numb (ch^'0') using namespace std; typedef long long ll; il void rd(int &x){ char ch;x=0;bool fl=false; while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true); for(x=numb;isdigit(ch=getchar());x=x*10+numb); (fl==true)&&(x=-x); } namespace Miracle{ const int N=1e5+5; const double Pi=acos(-1); ll n; ll buc[4*N],sum[4*N]; struct node{ double x,y; node(){} node(double xx,double yy){ x=xx;y=yy; } node friend operator *(node a,node b){ return node(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x); } node friend operator +(node a,node b){ return node(a.x+b.x,a.y+b.y); } node friend operator -(node a,node b){ return node(a.x-b.x,a.y-b.y); } void clear(){ x=0.0,y=0.0; } }a[4*N],b[4*N],c[4*N]; int rev[4*N]; void fft(node *f,int n,int c){ for(reg i=0;i<n;++i){ if(i<rev[i]) swap(f[i],f[rev[i]]); } for(reg p=2;p<=n;p<<=1){ int len=p/2; node gen=node(cos(Pi/len),c*sin(Pi/len)); for(reg l=0;l<n;l+=p){ node buf=node(1,0); for(reg k=l;k<l+len;++k){ node tmp=buf*f[k+len]; f[k+len]=f[k]-tmp; f[k]=f[k]+tmp; buf=buf*gen; } } } } ll mx; int m; void clear(){ for(reg i=0;i<m;++i) a[i].clear(),b[i].clear(),c[i].clear(); for(reg i=0;i<=mx;++i) buc[i]=0;mx=0; } int main(){ int T; rd(T); while(T--){ scanf("%lld",&n); clear(); int x; // for(reg i=0;i<16;++i){ // cout<<buc[i]<<" "; // }cout<<endl; for(reg i=1;i<=n;++i){ rd(x);mx=max(mx,(ll)x); buc[x]++; } m=1; for(;m<=mx*2;m<<=1); for(reg i=0;i<m;++i){ rev[i]=(rev[i>>1]>>1)|((i&1)?(m>>1):0); a[i].x=buc[i]; b[i].x=buc[i]; } //cout<<" mx "<<mx<<" "<<m<<endl; //cout<<" buc "<<endl; // for(reg i=0;i<m;++i){ // cout<<buc[i]<<" "; // }cout<<endl; fft(a,m,1);fft(b,m,1); for(reg i=0;i<m;++i) b[i]=a[i]*b[i]; fft(b,m,-1); //cout<<" dughrg "<<endl; for(reg i=0;i<m;++i) b[i].x=round(b[i].x/m); sum[mx]=buc[mx]; for(reg i=mx-1;i>=0;--i) sum[i]=sum[i+1]+buc[i]; double tot=n*(n-1)*(n-2)/6; double ans=tot; for(reg i=1;i<=mx;++i){ //cout<<" ii "<<i<<" "<<b[i].x<<endl; ll tmp=0; if(i%2==0){ tmp=buc[i/2]; } b[i].x-=tmp*tmp; b[i].x/=2; b[i].x+=tmp*(tmp-1)/2; //cout<<" b[i] "<<b[i].x<<endl; ans-=b[i].x*sum[i]; } printf("%.7lf\n",ans/tot); } return 0; } } signed main(){ Miracle::main(); return 0; } /* Author: *Miracle* Date: 2018/12/21 21:03:28 */