BZOJ 3527[Zjoi2014]力 FFT

版权声明:本文为博主原创文章,未经博主允许不得转载。 https://blog.csdn.net/qq_33929112/article/details/54590319

题目链接:BZOJ3527
第一次学会如何写数学公式,虽然只是简单的入门,但还是有点激动。。。
首先这个题很明显是多项式乘法,但是强迫症的我过于纠结下标,以至于困惑了好久,简直SB。

注:下标均从0开始。
现在进入正题。

Fj=i<jqiqj(ij)2i>jqiqj(ij)2

现在我们转化一下:记住下标从0开始
Fj=i=0j1qiqj(ij)2i=j+1n1qiqj(ij)2

Ej=Fjqj=i=0j1qi(ij)2i=j+1n1qi(ij)2

f(i)=qig(i)=1i2,其中g(0)=0

Ej=i=0j1f(i)g(ji)i=j+1n1f(i)g(ji)

等号右边第一个式子j1i=0f(i)g(ji)=ji=0f(i)g(ji),就是卷积的形式,直接FFT算即可。

等号右边第二个式子n1i=j+1f(i)g(ji)=n1i=jf(i)g(ji)=nj1i=0f(i+j)g(i)
倒序处理,令h(n1ij)=f(i+j),有nj1i=0f(i+j)g(i)=nj1i=0h(n1ij)g(i),记为Xi=nj1i=0h(n1ij)g(i),则Xnj1=ji=0h(ji)g(i),也是卷积的形式,直接FFT即可。
最后Ej=ji=0f(i)g(ji)Xnj1,其中ji=0f(i)g(ji)Xnj1均可由FFT直接算的,题目即可解决。
代码如下:

/**************************************************************
    Problem: 3527
    User: sfailsthy
    Language: C++
    Result: Accepted
    Time:3764 ms
    Memory:17680 kb
****************************************************************/

#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>
#include <cmath>
using namespace std;
const int maxn=262200+10;
const double PI =acos(-1.0);

struct Complex{
    double real,image;
    Complex(){
        real=image=0.0;
    }

    Complex(double a,double b){
        real=a;
        image=b;
    }

    Complex operator + (const Complex &s) const{
        return Complex(real+s.real,image+s.image);
    }

    Complex operator - (const Complex &s) const{
        return Complex(real-s.real,image-s.image);
    }

    Complex operator * (const Complex &s) const{
        return Complex(real*s.real-image*s.image,real*s.image+image*s.real);
    }
}x1[maxn],x2[maxn],A[maxn];

int n,rev[maxn];
double a[maxn/2],b[maxn/2],c[maxn/2];

void init(int &len,int len1,int len2){
    int k=1,L=0,r,t;
    while(k<2*len1||k<2*len2){
        k<<=1;
        L++;
    }

    len=k;
    for(int i=0;i<len;i++){
        t=i,r=0,k=L;
        while(k--){
            r<<=1;
            r|=t&1;
            t>>=1;
        }
        rev[i]=r;
    }
}

void FFT(Complex x[],int len,int op){
    Complex u,t;
    for(int i=0;i<len;i++){
        A[rev[i]]=x[i];
    }
    for(int i=0;i<len;i++){
        x[i]=A[i];
    }

    for(int k=2;k<=len;k<<=1){
        Complex wn(cos(2*PI/k*op),sin(2*PI/k*op));
        for(int i=0;i<len;i+=k){
            Complex w(1,0);
            for(int j=0;j<k/2;j++){
                u=x[i+j];
                t=w*x[i+j+k/2];
                x[i+j]=u+t;
                x[i+j+k/2]=u-t;
                w=w*wn;
            }
        }
    }

    if(op==-1){
        for(int i=0;i<len;i++){
            x[i].real/=len;
        }
    }
}

int main(){
    scanf("%d",&n);
    for(int i=0;i<n;i++){
        scanf("%lf",&a[i]);
    }
    for(int i=1;i<n;i++){
        b[i]=1.0/i/i;
    }

    int len,len1=n,len2=n;
    init(len,len1,len2);
    for(int i=0;i<len;i++){
        x1[i]=Complex();
        x2[i]=Complex();
    }
    for(int i=0;i<len;i++){
       x1[i]=Complex(a[i],0);
       x2[i]=Complex(b[i],0);
    }

    FFT(x1,len,1);
    FFT(x2,len,1);
    for(int i=0;i<len;i++){
        x1[i]=x1[i]*x2[i];
    }
    FFT(x1,len,-1);
    for(int i=0;i<n;i++){
        c[i]=x1[i].real;
    }

    for(int i=0;i<len;i++){
        x1[i]=Complex();
        x2[i]=Complex();
    }
    for(int i=0;i<n;i++){
        x1[i]=Complex(a[n-1-i],0);
        if(i) x2[i]=Complex(b[i],0);
    }
    FFT(x1,len,1);
    FFT(x2,len,1);
    for(int i=0;i<len;i++){
        x1[i]=x1[i]*x2[i];
    }
    FFT(x1,len,-1);
    for(int i=0;i<n;i++){
        c[i]-=x1[n-1-i].real;
    }

    for(int i=0;i<n;i++){
        printf("%.3lf\n",c[i]);
    }
    return 0;
}
展开阅读全文

没有更多推荐了,返回首页