您的位置:首页 > 大数据 > 物联网

HDU 4609 3-idiots(FFT学习)

2015-11-15 23:54 501 查看
附上kuangbin大神博客:讲解的非常清楚

分析:暴力枚举aka_k,那么问题变成了找满足ai+aj>ak(i,j,k互不相等)a_i+a_j>a_k(i,j,k互不相等)的(ai,aj)(a_i,a_j)的对数。FFT求出所有可能的和的对数,然后去掉坐标相同的就可以了。

PS:目前就写了两道FFT,基本上都是卷积的应用,FFT神奇的地就是能在nlognnlogn的时间内算出所有的X(k)=∑a[i]∗a[k−i]X(k)=\sum{a[i]*a[k-i]}

附上代码:

#include <bits/stdc++.h>
#define LL long long
#define FOR(i,x,y)  for(int i = x;i < y;++ i)
#define IFOR(i,x,y) for(int i = x;i > y;-- i)

using namespace std;

//FFT copy from kuangbin
const double pi = acos (-1.0);
// Complex  z = a + b * i
struct Complex {
double a, b;

Complex(double _a=0.0,double _b=0.0):a(_a),b(_b){}

Complex operator + (const Complex &rhs) const {
return Complex(a + rhs.a , b + rhs.b);
}

Complex operator - (const Complex &rhs) const {
return Complex(a - rhs.a , b - rhs.b);
}

Complex operator * (const Complex &rhs) const {
return Complex(a * rhs.a - b * rhs.b , a * rhs.b + b * rhs.a);
}
};

//len = 2 ^ k
void change (Complex y[] , int len) {
for (int i = 1 , j = len / 2 ; i < len -1 ; i ++) {
if (i < j) swap(y[i] , y[j]);
int k = len / 2;
while (j >= k) {
j -= k;
k /= 2;
}
if(j < k) j += k;
}
}

// FFT
// len = 2 ^ k
// on = 1  DFT    on = -1 IDFT
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].a /= len;
}
}
}

const int maxn = 100010;
const int maxm = 100010;

LL cnt[maxm<<2],sum[maxm<<2];
int a[maxn],n;

Complex x[maxm<<2];

void work(){
memset(cnt,0,sizeof(cnt));
FOR(i,0,n)  scanf("%d",&a[i]),++ cnt[a[i]];
sort(a,a+n);
int mx_len = (a[n-1] + 1)<<1;
int len = 1;
while(len < mx_len)    len <<= 1;
FOR(i,0,len)    x[i] = Complex(cnt[i],0.0);
FFT(x,len,1);
FOR(i,0,len)    x[i] = x[i]*x[i];
FFT(x,len,-1);
FOR(i,0,len)    cnt[i] = (LL)(x[i].a+0.5);
FOR(i,0,n)      -- cnt[a[i]<<1];
FOR(i,0,len)    cnt[i] >>= 1;
sum[0] = cnt[0];
FOR(i,1,len)    sum[i] = sum[i-1]+cnt[i];
LL ans = 0;
FOR(i,0,n){
ans += sum[len-1] - sum[a[i]];
ans -= (LL)(n-i-1)*i;
ans -= (LL)(n-1);
ans -= (LL)(n-i-1)*(n-i-2)>>1;
}
LL tot = (LL)n*(n-1)*(n-2)/6;
printf("%.7f\n",(double)ans/tot);
}

int main()
{
//freopen("test.in","r",stdin);
int T;  scanf("%d",&T);
while(T--){
scanf("%d",&n);
work();
}
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: