Loading...

洛谷P3338 ZJOI2014 力

check评论:0 条 remove_red_eye浏览量:1538 change_historyTags:编程学习笔记
作者 : deco date_range日期 : 2019-03-09

推式子吧

$E_j=\sum\limits_{i=0}^{j-1} \frac{q_i}{(j-i)^2}-\sum\limits_{i=j+1}^{n-1} \frac{q_i}{(j-i)^2}$

令$f[i]=\frac{1}{i^2}$

则$E_j=\sum\limits_{i=0}^{j-1} q[i]f[i]-\sum\limits_{i=j+1}^{n-1} q[i]f[j-i]$

后面一部分,翻转$f$数组,就和前面的一样是个卷积形式,$FFT$求解即可

#include <bits/stdc++.h>
using namespace std;
#define pi acos(-1)
#define cp complex<double>
cp w[400010],g[400010];
int limit=1,l,r[400010],n;
double t[400010],a[400010],b[400010],q[400010],p[400010];
void FFT(cp x[],int len,int flag)
{
    for(int i=0;i<len;i++)
    {
        if(i<r[i])
        {
            swap(x[i],x[r[i]]);
        }
    }
    for(int mid=1;mid<len;mid<<=1)
    {
        cp wn=(cp){cos(pi/mid),flag*sin(pi/mid)};
        for(int i=0,r=(mid<<1);i<len;i+=r)
        {
            cp ww=(cp){1,0};
            for(int k=0;k<mid;k++,ww=ww*wn)
            {
                cp m1=x[i+k],m2=ww*x[i+mid+k];
                x[i+k]=(m1+m2),x[i+mid+k]=(m1-m2);
            }
        }
    }
}
void work(double *x,double *y,double *z)
{
    for(int i=0;i<=limit;i++)
    {
        w[i]={x[i],0},g[i]={y[i],0};
    }
    FFT(w,limit,1),FFT(g,limit,1);
    for(int i=0;i<=limit;i++)
    {
        w[i]*=g[i];
    }
    FFT(w,limit,-1);
    for(int i=1;i<=n;i++)
    {
        z[i]=(w[i].real()/limit+0.5);
    }
}
int main()
{
    cin>>n;
    for(int i=1;i<=n;i++)
    {
        scanf("%lf",&p[i]);
        t[i]=(double)1.0/i/i;
        q[i]=p[i];
    }
    while(limit<(n<<1))
    {
        limit<<=1,++l;
    }
    for(int i=1;i<=limit;i++)
    {
        r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    }
    reverse(q+1,q+n+1);
    work(p,t,a);
    work(q,t,b);
    for(int i=1;i<=n;i++)
    {
        printf("%.3lf\n",a[i]-b[n-i+1]);
    }
}

暂无评论

正在回复给  
去登陆?

   点击刷新验证码

标签云

文章留名