BZOJ 3527 [Zjoi2014]力(FFT)

题意:中文题

思路:我们对于这个公式拆成利两部分看,一部分是i<j,一部分是剩下的,我们在i<j的部分可以直接进行化简成卷积的形式,而对于大于的部分比较麻烦,相当于对所有的都替换成n-j,这样你两部分加起来的和还是定值,然后在做一次卷积即可

代码:(抄胡小兔爷的--->传送门

 
#include <bits/stdc++.h>
using namespace std;
 
const int N = 1000005;
const double PI=acos(-1.0);
struct cp
{
    double r,i;
    cp(){}
    cp(double _r,double _i)
    {
        r=_r;i=_i;
    }
    cp operator +(const cp a)const
    {
        return cp(a.r+r,a.i+i);
    }
    cp operator -(const cp a)const
    {
        return cp(r-a.r,i-a.i);
    }
    cp operator *(const cp a)const
    {
        return cp(r*a.r-i*a.i,r*a.i+i*a.r);
    }
    cp conj()
    {
        return cp(r,-i);
    }
};
 
int n=1,m;
double q[N],res[N];
cp f[N],f1[N],g[N],omg[N],inv[N];
void FFT_init(){
    for(int i=0;i<n;i++){
        omg[i]=cp(cos(2*PI*i/n),sin(2*PI*i/n));
        inv[i]=omg[i].conj();
    }
}
void fft(cp *a,cp *omg){
    int lim=0;
    while((1 << lim) < n) lim++;
    for(int i=0;i<n;i++){
        int t=0;
        for(int j=0;j<lim;j++)
            if(i >> j&1) t |= 1 << (lim - j - 1);
        if(i<t) swap(a[i],a[t]);
    }
    for(int l=2;l<=n;l*=2){
        int m=l/2;
        for(cp *p=a;p!=a+n;p+=l){
            for(int i=0;i<m;i++){
                cp t=omg[n/l*i]*p[m+i];
                p[m+i]=p[i]-t;
                p[i]=p[i]+t;
            }
        }
    }
}
 
int main(){
    scanf("%d",&m);
    for(int i=0;i<m;i++){
        scanf("%lf",&q[i]);
        f[i].r=q[i];
        f1[m - 1 - i].r=q[i];
        if(i) g[i].r=(1.0 / i / i);
    }
    while(n<2*m)n*=2;
    FFT_init();
    fft(f,omg), fft(f1,omg), fft(g,omg);
    for(int i=0;i<n;i++) f[i]=f[i]*g[i],f1[i]=f1[i]*g[i];
    fft(f,inv), fft(f1,inv);
    for(int i=0;i<m;i++){
        res[i] = f[i].r/n - f1[m-1-i].r/n;
    }
    for(int i = 0; i < m; i++){
        printf("%lf\n", res[i]);
    }
    return 0;
}

 

posted @ 2018-08-03 20:50  啦啦啦天啦噜  阅读(137)  评论(0编辑  收藏  举报