bzoj3527: [Zjoi2014]力(fft)

传送门
f f t fft fft菜题。
题意简述:给一个数列 a i a_i ai,对于 i = 1 → n i=1\rightarrow n i=1n求出 a n s i = ∑ i &lt; j a i ( i − j ) 2 − ∑ i &gt; j a i ( i − j ) 2 ans_i=\sum_{i&lt;j}\frac{a_i}{(i-j)^2}-\sum_{i&gt;j}\frac{a_i}{(i-j)^2} ansi=i<j(ij)2aii>j(ij)2ai


思路:
考虑分开求减号前后的两组和。
前面的直接是一个卷积的形式,后面的可以把 a a a数组翻转一下也是一个卷积的形式,然后上 f f t fft fft求即可。
代码:

#include<bits/stdc++.h>
#define ri register int
using namespace std;
const double pi=acos(-1.0);
const int N=1e5+5;
struct cp{
    double x,y;
    cp(double _x=0,double _y=0):x(_x),y(_y){}
    friend inline cp operator+(const cp&a,const cp&b){return cp(a.x+b.x,a.y+b.y);}
    friend inline cp operator-(const cp&a,const cp&b){return cp(a.x-b.x,a.y-b.y);}
    friend inline cp operator*(const cp&a,const cp&b){return cp(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
    friend inline cp operator/(const cp&a,const double&b){return cp(a.x/b,a.y/b);}
}a[N];
vector<int>pos;
vector<cp>A,B;
int lim,tim,n;
inline void init(const int&up){
    lim=1,tim=0;
    while(lim<=up)lim<<=1,++tim;
    pos.resize(lim),A.resize(lim),B.resize(lim),pos[0]=0;
    for(ri i=0;i<lim;++i)pos[i]=(pos[i>>1]>>1)|((i&1)<<(tim-1));
}
inline void fft(vector<cp>&a,const int&type){
    for(ri i=0;i<lim;++i)if(i<pos[i])swap(a[i],a[pos[i]]);
    cp wn,w,a0,a1;
    for(ri mid=1;mid<lim;mid<<=1){
        wn=cp(cos(pi/mid),sin(pi/mid)*type);
        for(ri j=0,len=mid<<1;j<lim;j+=len){
            w=cp(1,0);
            for(ri k=0;k<mid;++k,w=w*wn)a0=a[j+k],a1=w*a[j+k+mid],a[j+k]=a0+a1,a[j+k+mid]=a0-a1;
        }
    }
    if(type==-1)for(ri i=0;i<lim;++i)a[i]=a[i]/lim;
}
struct poly{
    vector<cp>a;
    poly(int k=0,cp x=cp(0,0)){a.resize(k+1),a[k]=x;}
    inline cp&operator[](const int&k){return a[k];}
    inline const cp&operator[](const int&k)const{return a[k];}
    inline int deg()const{return a.size()-1;}
    inline poly extend(const int&k){poly ret=*this;return ret.a.resize(k+1),ret;}
    friend inline poly operator*(const poly&a,const poly&b){
        int n=a.deg(),m=b.deg();
        init(n+m);
        for(ri i=0;i<=n;++i)A[i]=a[i];
        for(ri i=0;i<=m;++i)B[i]=b[i];
        for(ri i=n+1;i<lim;++i)A[i]=cp(0,0);
        for(ri i=m+1;i<lim;++i)B[i]=cp(0,0);
        fft(A,1),fft(B,1);
        for(ri i=0;i<lim;++i)A[i]=A[i]*B[i];
        poly ret;
        return fft(A,-1),ret.a=A,ret;
    }
}X,Y;
double ans[2][N];
inline void solve(int type){
    poly X(n),Y(n);
    X[0]=Y[0]=cp(0,0);
    for(ri i=1;i<=n;++i)X[i]=a[i],Y[i]=cp(1.0/(double)i/(double)i,0);
    X=(X*Y).extend(n);
    for(ri i=1;i<=n;++i)ans[type][i]=X[i].x;
}
int main(){
    scanf("%d",&n);
    for(ri i=1;i<=n;++i)scanf("%lf",&a[i].x);
    solve(1),reverse(a+1,a+n+1),solve(0),reverse(ans[0]+1,ans[0]+n+1);
    for(ri i=1;i<=n;++i)printf("%.3lf\n",ans[1][i]-ans[0][i]);
    return 0;
}
  • 1
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值