洛谷 P3338 [ZJOI2014]力

题意简述

读入\(n\)个数\(q_i\)

设\(F_j = \sum\limits_{i<j}\frac{q_i\times q_j}{(i-j)^2 }-\sum\limits_{i>j}\frac{q_i\times q_j}{(i-j)^2 }\)
令\(E_i=\frac{F_i}{q_i}\)，求\(E_i\)

题解思路

先推式子

\(E_j=\frac{F_j}{q_j}=\sum\limits_{i<j}\frac{q_i}{(i-j)^2 }-\sum\limits_{i>j}\frac{q_i}{(i-j)^2 }\)

设\(T_i=i^{-2}\)

则\(E_j=\sum\limits_{i=1}^{j-1}q_iT_{j-i}-\sum\limits_{i=j+1}^{n}q_iT_{i-j}\)

再设\(p_i=q_{n-i+1}\)

则\(E_j=\sum\limits_{i=1}^{j-1}q_iT_{j-i}-\sum\limits_{i=1}^{j-1}p_iT_{j-i}\)

可以发现，这两项都是两个多项式卷积的结果，所以可以用FFT来做

注意：如果用三次变两次优化，会掉精

代码

#include <cmath>
#include <iostream>
#include <algorithm>
const int N=400005;
const double Pi=acos(-1.0);
int n,nn,m,lim=1,l=-1,r[N];
struct Com {
    long double x,y;
    Com(long double xx=0,long double yy=0) { x=xx; y=yy; }
}p[N],q[N],t[N];
Com operator +(const Com& x,const Com& y) { return Com(x.x+y.x,x.y+y.y); }
Com operator -(const Com& x,const Com& y) { return Com(x.x-y.x,x.y-y.y); }
Com operator *(const Com& x,const Com& y) { return Com(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x); }
inline void FFT(Com *x,const int& lim,const int& type) {
    for (register int i=0;i<lim;++i) if (i<r[i]) std::swap(x[i],x[r[i]]);
    for (register int i=1;i<lim;i<<=1) {
        Com UR(cos(Pi/i),sin(Pi/i)*type);
        for (register int j=0,I=i<<1;j<lim;j+=I) {
            Com w(1,0);
            for (register int k=0;k<i;++k,w=w*UR) {
                Com t1=x[j+k],t2=w*x[j+k+i];
                x[j+k]=t1+t2; x[j+k+i]=t1-t2;
            }
        }
    }
}
int main() {
    std::ios::sync_with_stdio(0);
    std::cin.tie(0); std::cout.tie(0); std::cout<<std::fixed;
    std::cin>>n;
    for (register int i=1;i<=n;++i) std::cin>>p[i].x,t[i].x=1.0/i/i;
    for (register int i=1;i<=n;++i) q[i].x=p[n-i+1].x;
    for (nn=n<<1;lim<=nn;lim<<=1,++l);
    for (register int i=0;i<lim;++i) r[i]=(r[i>>1]>>1)|((i&1)<<l);
    FFT(p,lim,1); FFT(q,lim,1); FFT(t,lim,1);
    for (register int i=0;i<lim;++i) p[i]=p[i]*t[i],q[i]=q[i]*t[i];
    FFT(p,lim,-1); FFT(q,lim,-1);
    for (register int i=1;i<=n;++i)
        std::cout<<(p[i].x-q[n-i+1].x)/lim<<'\n';
}

转载于:https://www.cnblogs.com/xuyixuan/p/11326566.html