BZOJ3527 [ZJOI2014]力

2014.12.07 17:50 Sun| 17 visits oi_2015| 2015_刷题日常| Text

Problem

Description

给出n个数qi,给出Fj的定义如下:

令Ei=Fi/qi。试求Ei。

Input

输入文件force.in包含一个整数n,接下来n行每行输入一个数,第i行表示qi。

Output

输出文件force.out有n行,第i行输出Ei。与标准答案误差不超过 1e-2 即可。

Sample Input

5

4006373.885184

15375036.435759

1717456.469144

8514941.004912

1410681.345880

Sample Output

-16838672.693

3439.793

7509018.566

4595686.886

10903040.872

HINT

对于 30%的数据,n≤1000。

对于 50%的数据,n≤60000。

对于 100%的数据,n≤100000,0<qi<1000000000。

Solution

关于FFT,参见 Fast Fourier Transformation

由题意,可得出如下表格(以 n=5 为例)。其中第 i 行第 j 列的数字表示 $q_j$ 对答案 $E_i$ 的贡献。

可见,答案以“0”为界,被划分成了相似的两部分,其中左下的部分对答案的贡献为正,右上的部分贡献为负。先考虑左下部分,显然可看做是卷积的形式,即 $( q_1, q_2, \dots, q_n ) * ( 0, {1 \over 1^2}, \dots, {1 \over n^2} )$ 。

同理,右上部分也是卷积的形式 $( q_n, q_{n-1}, \dots, q_1 ) * ( 0, {1 \over 1^2}, \dots, {1 \over n^2} )$ 。求出卷积后和左下部分做一次减法,即为答案。总时间复杂度 $\mathcal O(n \log n)$ 。快速傅里叶变换还真是强大呢!

快乐地敲一遍FFT的短♂小♂精♂悍的模板。搞定!

// 注意求 ${1 / i^2}$ 不能写做 $1/(i*i)$ ……

Code

#include <bits/stdc++.h>
using namespace std;
#define N 1000000

int n,m; double q[N];
complex<double> a[N],b[N],p1[N],p2[N];
const double pi=3.14159265358979323;

void FFT(complex<double> x[],int n,int type)
{
    for(int i=0,t=0;i<n;i++)
    {
        if(i>t) swap(x[i],x[t]);
        for(int j=n>>1;(t^=j)<j;j>>=1);
    }
    for(int s=2;s<=n;s<<=1)
    {
        complex<double> wn(cos(type*2*pi/s),sin(type*2*pi/s));
        for(int i=0;i<n;i+=s)
        {
            complex<double> w(1,0),t;
            for(int j=0;j<s>>1;w*=wn,j++)
                t=w*x[i+j+(s>>1)],
                x[i+j+(s>>1)]=x[i+j]-t,x[i+j]+=t;
        }
    }
}

int main()
{
    cin>>m;
    for(int i=0;i<m;i++) scanf("%lf",&q[i]);
    for(int t=m,i=1;i>>2<t;i<<=1) n=i;
    for(int i=1;i<m;i++) b[i]=1.0/i/i; FFT(b,n,1);
    for(int i=0;i<m;i++) a[i]=q[i];
    FFT(a,n,1); for(int i=0;i<n;i++) p1[i]=a[i]*b[i]; FFT(p1,n,-1);
    memset(a,0,sizeof a);
    for(int i=0;i<m;i++) a[i]=q[m-i-1];
    FFT(a,n,1); for(int i=0;i<n;i++) p2[i]=a[i]*b[i]; FFT(p2,n,-1);
    for(int i=0;i<m;i++)
        printf("%.3lf\n",(p1[i].real()-p2[m-i-1].real())/n);
    return 0;
}