【BZOJ3527】【FFT】力

【问题描述】
给出n个数qi，给出Fj的定义如下：
令Ei=Fi/qi。试求Ei。
【输入格式】
输入文件force.in包含一个整数n，接下来n行每行输入一个数，第i行表示qi。
【输出格式】
输出文件force.out有n行，第i行输出Ei。与标准答案误差不超过1e-2即可。
【样例输入】
5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880
【样例输出】
-16838672.693
3439.793
7509018.566
4595686.886
10903040.872
【数据规模与约定】
对于30%的数据，n≤1000。
对于50%的数据，n≤60000。
对于100%的数据，n≤100000，0<qi<1000000000。

【分析】

这道题...在省选里面相当裸了。

自己把式子展开一下，发现跟卷积是类似的。

于是对公式的前半部分做一下FFT，后半部分再做一下FFT，减一下，然后就是公式的样子了。

感觉对FFT的理解更进一步了。

/*
宋代苏轼
《临江仙·夜饮东坡醒复醉》
夜饮东坡醒复醉，归来仿佛三更。家童鼻息已雷鸣。敲门都不应，倚杖听江声。
长恨此身非我有，何时忘却营营。夜阑风静縠纹平。小舟从此逝，江海寄余生。  
*/
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <queue>
#include <vector>
#include <iostream>
#include <string>
#include <ctime>
#define LOCAL
const double Pi = acos(-1.0);
const int MAXN = 100000 * 2 * 2 + 10;
using namespace std;
struct Num{
       double a, b;
       Num(double x = 0, double y = 0){a = x; b = y;}
       Num operator + (const Num &c){return Num(a + c.a, b + c.b);}
       Num operator - (const Num &c){return Num(a - c.a, b - c.b);}
       Num operator * (const Num &c){return Num(a * c.a - b * c.b, a * c.b + b * c.a);}
}x1[MAXN], x2[MAXN];
double data[MAXN], Ans[MAXN];
int n;
//交换成蝴蝶顺序 
void change(Num *t, int len, int loglen){
     for (int i = 0; i < len; i++){
         int k = 0, x = i, tmp = loglen;
         while (tmp--) {k = (k<<1) + (x & 1);x >>= 1;}
         if (k < i) swap(t[k], t[i]);
     } 
     return;
}
//0为逆向 
void FFT(Num *x, int len, int loglen, int type){
     if (type) change(x, len, loglen); 
     int t;//t代表长度
     t = (type ? 1 : (1<<loglen)); 
     for (int i = 0; i < loglen; i++){
         if (!type) t >>= 1;
         int l = 0, r = l + t;
         while (l < len){
               Num a, b;//临时变量
               Num tmp(1, 0), w(cos(Pi / t), (type ? 1 : -1) * sin(Pi / t)); 
               for (int j = l; j < l + t; j++){
                   if (type){
                      a = x[j];
                      b = x[j + t] * tmp;
                      x[j] = a + b;
                      x[j + t] = a - b;
                   }else{
                      a = x[j] + x[j + t];
                      b = (x[j] - x[j + t]) * tmp;
                      x[j] =  a;
                      x[j + t] = b;
                   }
                   tmp = tmp * w;
               }
               l = r + t;
               r = l + t; 
         }
         if (type) t <<= 1;
     }     
     if (!type){
        change(x, len, loglen);
        for (int i = 0; i < len; i++) x[i].a /= len;
     } 
}
void init(){
     memset(x1, 0, sizeof(x1));
     memset(x2, 0, sizeof(x2));
     int len = 0;
     while ((1 << len) < n) len++;
     len++;
     for (int i = 0; i < n; i++) x1[i] = Num(data[i], 0);
     for (int i = 1; i < n; i++) x2[i] = Num((double)1.0 / (double)(i * (double)i), 0);
     //for (int i = 1; i < n; i++) printf("%lf
", x2[i].a);
     
     FFT(x1, (1<<len), len, 1);
     FFT(x2, (1<<len), len, 1);
     for (int i = 0; i < (1 << len); i++) x1[i] = x1[i] * x2[i];
     FFT(x1, (1<<len), len, 0);
}
void debug(){
     int len = 0;
     scanf("%d", &n);
     while ((1<<len) <= (n << 1)) len++;
     for (int i = 0; i < n; i++) scanf("%lf", &x1[i].a);
     for (int i = 0; i < n; i++) scanf("%lf", &x2[i].a);
     FFT(x1, (1<<len), len, 1);
     FFT(x2, (1<<len), len, 1);
     for (int i = 0; i < (1 << len); i++) x1[i] = x1[i] * x2[i];
     FFT(x1, (1<<len), len, 0);
     for (int i = 0; i < n; i++) printf("%lf
", x1[i].a); 
}

int main() {
    
    scanf("%d", &n);
    for (int i = 0; i < n; i++) scanf("%lf", &data[i]);
    init();
    for (int i = 0; i < n; i++) Ans[i] = x1[i].a;
    reverse(data, data + n);
    init();
    for (int i = 0; i < n; i++) Ans[i] -= x1[n - 1 - i].a;
    for (int i = 0; i < n; i++) printf("%.3lf
", Ans[i]);
    //debug();
    return 0;
}