洛谷 P3723 [AH/HNOI2017]礼物(bzoj P4827 [Hnoi2017]礼物)

来源:互联网 发布:2017年python饱和了 编辑:程序博客网 时间:2024/05/29 08:39

传送门

当初省选时看这题一脸懵逼,打了70分的代码事实满分还是很容易的(啪啪啪,自己wa了十次还说容易??

好了说下题解吧,这题就是求循环卷积,我的代码是0~n-1,式子是1~n(主要是博主FFT习惯写0~n-1),照样发图片(博主太蠢,不会输数学公式。。)

注意:比较大小时的精度!!!在这里wa了10次!!!

Tips:round()四舍五入函数

代码:

#include<iostream>#include<stdio.h>#include<math.h>#include<stdlib.h>using namespace std;const int Maxn=1000005;const int oo=2147483647;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;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);} int p[Maxn],q[Maxn],rev[Maxn],dig[Maxn];int ans1=0,sumx,sumy,N,L,ans2=-oo,ans=oo,ansx,ansy;int 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=a[i].real/N;}int main(){freopen("s(2).in","r",stdin);int n,m;scanf("%d%d",&n,&m);for(int i=0;i<n;i++)scanf("%d",&p[i]);for(int i=0;i<n;i++){scanf("%d",&q[i]);q[i+n]=q[i];}for(int i=0;i<n;i++)ans1=ans1+(p[i]*p[i]+q[i]*q[i]);for(L=0,N=1;N<2*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[n-i-1],0);for(int i=0;i<2*n;i++)b[i]=complex(q[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++){ans2=max(ans2,(int)round(c[n+i-1].real));ansx+=p[i];ansy+=q[i];} for(int i=-m-1;i<=m+1;i++){int ans3=ans1-2*ans2;ans=min(ans,ans3+2*i*ansx-2*i*ansy+n*i*i);}printf("%d",ans);return 0;} 


原创粉丝点击