码迷,mamicode.com
首页 > 其他好文 > 详细

hdu 4609 FFT

时间:2015-04-28 20:29:10      阅读:184      评论:0      收藏:0      [点我收藏+]

标签:

题意:给出一堆数,问从这些数中取3个能组成三角形的概率?

 

sol:其实就是问从这些数里取3个组成三角形有多少种取法

脑洞大开的解法:用FFT

设一开始的数是1 3 3 4

作一个向量x,其中x[i]=边长为i的边的个数

那么就有x=[0 1 0 2 1 0 0 0 0]

令y=x,对x和y作DFT,得到dx和dy。令dn=dx*dy,再对dn作IDFT得到n

那么就得到n=[0 0 1 0 4 2 4 4 1 0 ]

其中n[i]=在x和y中各选一条边,使得两条边之和为i有几种方案

这时得到的n并不好,包含了各种重复的方案。还得减一下

最后得到的n=[0 0 0 0 2 1 1 2 0]

 

然后对n做奇怪的操作。。。

先求前缀和,得到S

技术分享
对于a[i].  我们假设a[i]是形成的三角形中最长的。这样就是在其余中选择两个和>a[i],而且长度不能大于a[i]的。(注意这里所谓的大于小于,不是说长度的大于小于,其实是排好序以后的,位置关系,这样就可以不用管长度相等的情况,排在a[i]前的就是小于的,后面的就是大于的)。

根据前面求得的结果。
长度和大于a[i]的取两个的取法是sum[len]-sum[a[i]].

但是这里面有不符合的。

一个是包含了取一大一小的
cnt -= (long long)(n-1-i)*i;

一个是包含了取一个本身i,然后取其它的
cnt -= (n-1);

还有就是取两个都大于的了
cnt -= (long long)(n-1-i)*(n-i-2)/2;
View Code

 

 

Code:

技术分享
  1 #include  <stdio.h>
  2 #include  <string.h>
  3 #include  <iostream>
  4 #include  <algorithm>
  5 #include  <math.h>
  6 using namespace  std;
  7 const double PI = acos(-1.0);
  8 #define LL long long
  9 #define MAXL 524230
 10 #define MAXN 100010
 11 
 12 //复数结构体
 13 struct  Complex
 14 {
 15     double x,y;//实部和虚部  x+yi
 16     Complex(double _x = 0,double _y = 0)
 17     {
 18         x = _x;
 19         y = _y;
 20     }
 21     Complex operator -(const Complex &b)const
 22     {
 23         return  Complex(x-b.x,y-b.y);
 24     }
 25     Complex operator +(const Complex &b)const
 26     {
 27         return  Complex(x+b.x,y+b.y);
 28     }
 29     Complex operator *(const Complex &b)const
 30     {
 31         return  Complex(x*b.x-y*b.y,x*b.y+y*b.x);
 32     }
 33 };
 34 /*
 35 *进行FFT和IFFT前的反转变换。
 36 *位置i和(i二进制反转后位置)互换
 37 * len必须去2的幂
 38 */
 39 void change(Complex y[],int  len)
 40 {
 41     int  i,j,k;
 42     for(i = 1, j = len/2; i <len-1; i++)
 43     {
 44         if(i < j)swap(y[i],y[j]);
 45         //交换互为小标反转的元素,i<j保证交换一次
 46         //i做正常的+1,j左反转类型的+1,始终保持i和j是反转的
 47         k = len/2;
 48         while(j >= k)
 49         {
 50             j -= k;
 51             k /= 2;
 52         }
 53         if(j < k)j += k;
 54     }
 55 }
 56 /*
 57 *做FFT
 58 * len必须为2^k形式,
 59 * on==1时是DFT,on==-1时是IDFT
 60 */
 61 //fft(x,len,1):对向量x做DFT(时域->频域),向量长度为1--len
 62 //fft(x,len,-1):做IDFT(频域->时域)
 63 void fft(Complex y[],int len,int  on)
 64 {
 65     change(y,len);
 66     for(int h = 2; h <= len; h <<= 1)
 67     {
 68         Complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
 69         for(int j = 0; j < len; j+=h)
 70         {
 71             Complex  w(1,0);
 72             for(int k = j; k < j+h/2; k++)
 73             {
 74                 Complex u = y[k];
 75                 Complex t = w*y[k+h/2];
 76                 y[k] = u+t;
 77                 y[k+h/2] = u-t;
 78                 w = w*wn;
 79             }
 80         }
 81     }
 82     if(on == -1)
 83         for(int i = 0; i < len; i++)
 84             y[i].x /= len;
 85 }
 86 
 87 Complex x[MAXL];
 88 LL nx[MAXL],tri[MAXN],sx[MAXL],NUM[MAXL];
 89 int TIMES,N,tm;
 90 
 91 int main()
 92 {
 93     scanf("%d",&TIMES);
 94     for (int Times=1;Times<=TIMES;Times++)
 95     {
 96         scanf("%d",&N);
 97 
 98         memset(nx,0,sizeof(nx));
 99         memset(sx,0,sizeof(sx));
100         memset(x,0,sizeof(x));
101         memset(NUM,0,sizeof(NUM));
102 
103         for (int i=0;i<N;i++)
104         {
105             scanf("%d",&tri[i]);
106             NUM[tri[i]]++;
107         }
108 
109         sort(tri,tri+N);
110         int mx=tri[N-1];
111         int l1=mx+1;
112         int len=1;
113         while (len<2*l1) len<<=1;
114         for (int i=0;i<l1;i++)
115             x[i]=Complex(NUM[i],0);
116         for (int i=l1;i<len;i++)
117             x[i]=Complex(0,0);
118 
119         fft(x,len,1);
120 
121         for (int i=0;i<len;i++)
122             x[i]=x[i]*x[i];
123 
124         fft(x,len,-1);
125         for (int i=0;i<len;i++)
126             nx[i]=(LL)(x[i].x+0.5);
127 
128         len=2*mx;
129         for (int i=0;i<N;i++)
130             nx[tri[i]*2]--;
131         for (int i=1;i<=len;i++)
132             nx[i]=nx[i]/2;
133 /*
134         for (int i=0;i<=len;i++)
135             cout<<nx[i]<<" ";
136         cout<<endl;
137 */
138         sx[0]=0;
139         for (int i=1;i<=len;i++)
140             sx[i]=sx[i-1]+nx[i];
141 /*
142         for (int i=0;i<=len;i++)
143             cout<<sx[i]<<" ";
144         cout<<endl;
145 */
146         LL tot=0;
147         for (int i=0;i<N;i++)
148         {
149             tot+=sx[len]-sx[tri[i]];
150             tot-=(LL)(N-i-1)*i;
151             tot-=(N-1);
152             tot-=(LL)(N-1-i)*(N-i-2)/2;
153         }
154         //cout<<tot<<endl;
155         LL kk=(LL)(N*(N-1)*(N-2)/6);
156         //cout<<tot<<" "<<kk<<endl;
157         printf("%.7lf\n",(double)tot/kk);
158     }
159 
160     return 0;
161 }
View Code

 

 

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

 

hdu 4609 FFT

标签:

原文地址:http://www.cnblogs.com/pdev/p/4463823.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!