洛谷 P3338 [ZJOI2014]力解题报告

P3338 [ZJOI2014]力

题目描述

给出n个数qi，给出Fj的定义如下：

(F_j = sum_{i<j}frac{q_i q_j}{(i-j)^2 }-sum_{i>j}frac{q_i q_j}{(i-j)^2 })

令(E_i=frac{F_i}{q_i})，求(E_i).

输入输出格式

输入格式：

第一行一个整数(n)。

接下来(n)行每行输入一个数，第(i)行表示(q_i)。

输出格式：

(n)行，第(i)行输出(E_i)。

与标准答案误差不超过(10^{-2})即可。

说明

对于(30\%)的数据，(nle1000)。

对于(50\%)的数据，(nle60000)。

对于(100\%)的数据，(nle100000)，(0<qi<1000000000)。

[spj 0.01]

鉴于本傻子一开始啥也看不出来还是写一写好了。

把两边分开考虑。令(f(x)=frac{1}{x^2},g(x)=q_x)，然后发现就是(FFT)的形式，看成系数做就可以了。

右边转过来以后同理。

有一些细节比如

1/(i*i)会被卡精度，要写1/i/i
(f(0)=g(0)=0)

Code:

#include <cstdio>
#include <cmath>
#include <algorithm>
const int N=(1<<18)+10;
struct complex
{
    double x,y;
    complex(){}
    complex(double x,double y){this->x=x,this->y=y;}
    complex friend operator +(complex n1,complex n2){return complex(n1.x+n2.x,n1.y+n2.y);}
    complex friend operator -(complex n1,complex n2){return complex(n1.x-n2.x,n1.y-n2.y);}
    complex friend operator *(complex n1,complex n2){return complex(n1.x*n2.x-n1.y*n2.y,n1.x*n2.y+n1.y*n2.x);}
}a[N],b[N],tmpx,tmpy,wn,w;
const double pi=3.1415926535897632;
int n,turn[N],len=1,L=-1;
double out[N],q[N];
void FFT(complex *a,int typ)
{
    for(int i=0;i<len;i++)
        if(i<turn[i])
            std::swap(a[i],a[turn[i]]);
    for(int le=1;le<len;le<<=1)
    {
        wn=complex(cos(pi/le),typ*sin(pi/le));
        for(int p=0;p<len;p+=le<<1)
        {
            w=complex(1,0);
            for(int i=p;i<p+le;i++,w=w*wn)
            {
                tmpx=a[i],tmpy=w*a[i+le];
                a[i]=tmpx+tmpy,a[i+le]=tmpx-tmpy;
            }
        }
    }
}
int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;i++) scanf("%lf",q+i);
    for(int i=1;i<=n;i++) a[i].x=q[i],b[i].x=1.0/i/i;
    while(len<=n+1<<1) len<<=1,++L;
    for(int i=0;i<len;i++) turn[i]=turn[i>>1]>>1|(i&1)<<L;
    FFT(a,1),FFT(b,1);
    for(int i=0;i<len;i++) a[i]=a[i]*b[i];
    FFT(a,-1);
    for(int i=1;i<=n;i++) out[i]+=a[i].x/len;
    for(int i=0;i<len;i++) a[i]=complex(0,0);
    for(int i=1;i<=n;i++)
        a[n+1-i]=complex(q[i],0);
    FFT(a,1);
    for(int i=0;i<len;i++) a[i]=a[i]*b[i];
    FFT(a,-1);
    for(int i=1;i<=n;i++) out[n+1-i]-=a[i].x/len;
    for(int i=1;i<=n;i++) printf("%lf
",out[i]);
    return 0;
}

2018.12.3