洛谷 P3338 [ZJOI2014]力（bzoj P3527 [ZJOI2014]力）

来源：互联网发布：python 股票数据库编辑：程序博客网时间：2024/05/29 02:26

传送门

很难受，不会输公式，于是手写了一份，比较丑，希望不要介意。fft的套路果然深。。照着网上模板学习了一下。。

代码：

#include<iostream>#include<stdio.h>#include<string.h>#include<math.h> using namespace std;const int Maxn=1000005;int n,N,L;int rev[Maxn];int dig[Maxn];double p[Maxn];struct complex{double real,image;complex(){real=image=0.0;}complex(double a,double b){real=a;image=b;}}a[Maxn],b[Maxn],c[Maxn],A[Maxn],x,y,ans[Maxn];complex operator +(complex a,complex b){return complex(a.real+b.real,a.image+b.image);}complex operator -(complex a,complex b){return complex(a.real-b.real,a.image-b.image);}complex operator *(complex a,complex b){return complex(a.real*b.real-a.image*b.image,a.real*b.image+a.image*b.real);}void FFT(complex a[],int flag){      for(int i=0;i<N;i++)A[i]=a[rev[i]];      for(int i=0;i<N;i++)a[i]=A[i];      for(int i=2;i<=N;i<<=1){          complex wn(cos(2*M_PI/i),flag*sin(2*M_PI/i));          for(int k=0;k<N;k+=i){              complex w(1,0);              for(int j=0;j<i/2;j++){                x=a[k+j];                  y=w*a[k+j+i/2];                  a[k+j]=x+y;                  a[k+j+i/2]=x-y;                  w=w*wn;                }          }      }      if(flag==-1)for(int i=0;i<N;i++)a[i].real/=N;  }  int main(){scanf("%d",&n);for(int i=0;i<n;i++)scanf("%lf",&p[i]);for(L=0,N=1;N<n;N<<=1,L++);N<<=1;L++;for(int i=0;i<N;i++){        int len=0;          for(int t=i;t;t>>=1)dig[len++]=t&1;          for(int j=0;j<L;j++)rev[i]=rev[i]*2+dig[j];      }    for(int i=0;i<n;i++)a[i]=complex(p[i],0);    for(int i=1;i<n;i++)b[i]=complex(1.0/i/i,0);    FFT(a,1);FFT(b,1);for(int i=0;i<N;i++)c[i]=a[i]*b[i];FFT(c,-1);for(int i=0;i<n;i++)ans[i].real=c[i].real;memset(a,0,sizeof(a));      memset(b,0,sizeof(b));  for(int i=0;i<n;i++)a[i]=complex(p[n-i-1],0);    for(int i=1;i<n;i++)b[i]=complex(1.0/i/i,0);FFT(a,1);FFT(b,1);for(int i=0;i<N;i++)c[i]=a[i]*b[i];FFT(c,-1);for(int i=0;i<n;i++)ans[i].real-=c[n-i-1].real;for(int i=0;i<n;i++)printf("%f\n",ans[i].real);return 0;}

阅读全文

0 0