首页 > 其他 > 详细

[洛谷P3338] ZJOI2014 力

时间:2019-06-29 11:23:26      阅读:88      评论:0      收藏:0      [点我收藏+]

问题描述

给出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 } \]
令Ei=Fi/qi,求Ei.

输入格式

第一行一个整数n。

接下来n行每行输入一个数,第i行表示qi。

输出格式

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

解析

首先,我们先来化简一下式子:
\[ \begin{align} E[i]=&\frac{F[i]}{p[i]}\\=&\sum_{j<i}{\frac{p[j]}{(i-j)^2}}-\sum_{j>i}{\frac{p[j]}{(i-j)^2}} \end{align} \]
那么,设\(f1[i]=p[i],g[i]=\frac{1}{i^2}\),则原式可以化简得:
\[ E[i]=\sum_{j=1}^{i-1}f1[j]*g[i-j]-\sum_{j=i+1}^{n}f1[j]*g[j-i] \]
为了将等式右边的两部分变得形式一样,不妨设\(f2[i]=f1[n-i+1]\),即将\(f1\)翻转得到\(f2\),就可以达到我们的目的:
\[ E[i]=\sum_{j=1}^{i-1}f1[j]*g[i-j]-\sum_{j=1}^{i-1}f2[j]*g[i-j] \]
观察到等式的形式与多项式卷积的形式十分的相似,就是卷积中第\(i-1\)次方项的系数表达式。因此,我们可以利用FFT分别求出\(f1\)\(g\)\(f2\)\(g\)的卷积,然后答案即为对应次数项的差。

代码

#include <iostream>
#include <cstdio>
#include <cmath>
#define N 400002
using namespace std;
const double PI=acos(-1);
struct Complex{
    double r,i;
}a[N],b[N];
Complex operator + (Complex a,Complex b) {return (Complex){a.r+b.r,a.i+b.i};}
Complex operator - (Complex a,Complex b) {return (Complex){a.r-b.r,a.i-b.i};}
Complex operator * (Complex a,Complex b) {return (Complex){a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r};}
int n=1,m,i,lim=1,r[N];
double f1[N],f2[N],g[N],ansa[N],ansb[N];
void trans()
{
    double tmp[N];
    for(int i=m;i>=1;i--) tmp[m-i+1]=f2[i];
    for(int i=1;i<=m;i++) f2[i]=tmp[i];
}
void FFT(Complex *a,int inv)
{
    for(int i=0;i<n;i++){
        if(i<r[i]) swap(a[i],a[r[i]]);
    }
    for(int l=2;l<=n;l*=2){
        int mid=l/2;
        Complex omg=(Complex){cos(2*PI/l),inv*sin(2*PI/l)};
        for(int i=0;i<n;i+=l){
            Complex w=(Complex){1,0};
            for(int j=0;j<mid;j++,w=w*omg){
                Complex tmp=a[i+j+mid]*w;
                a[i+j+mid]=a[i+j]-tmp;
                a[i+j]=a[i+j]+tmp;
            }
        }
    }
}
void solve(double *f,double *g,double *ans)
{
    for(int i=0;i<n;i++) a[i].r=f[i],a[i].i=0;
    for(int i=0;i<n;i++) b[i].r=g[i],b[i].i=0;
    FFT(a,1);FFT(b,1);
    for(int i=0;i<n;i++) a[i]=a[i]*b[i];
    FFT(a,-1);
    for(int i=1;i<=m;i++) ans[i]=a[i].r/n;
}
int main()
{
    cin>>m;
    while(n<2*m) n*=2;
    while((1<<lim)<n) lim++;
    for(i=0;i<n;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(lim-1));
    for(i=1;i<=m;i++){
        cin>>f1[i];
        g[i]=1.0/i/i;
        f2[i]=f1[i];
    }
    trans();
    solve(f1,g,ansa);
    solve(f2,g,ansb);
    for(i=1;i<=m;i++) printf("%.3lf\n",ansa[i]-ansb[m-i+1]);
    return 0;
}

[洛谷P3338] ZJOI2014 力

原文:https://www.cnblogs.com/LSlzf/p/11105283.html

(0)
(0)
   
举报
评论 一句话评论(0
关于我们 - 联系我们 - 留言反馈 - 联系我们:wmxa8@hotmail.com
© 2014 bubuko.com 版权所有
打开技术之扣,分享程序人生!