2016 acm 北京网络赛 F题 Periodic Signal 循环卷积

#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<ctime>
#include<algorithm>
#define N 200005
#define pi acos(-1.0) // PI值
typedef long long ll;
#define LL long long
using namespace std;
struct complex
{
	double r,i;
	complex(double real=0.0,double image=0.0){
		r=real;	i=image;
	}
	// 以下为三种虚数运算的定义
	complex operator + (const complex o){
		return complex(r+o.r,i+o.i);
	}
	complex operator - (const complex o){
		return complex(r-o.r,i-o.i);
	}
	complex operator * (const complex o){
		return complex(r*o.r-i*o.i,r*o.i+i*o.r);
	}
}x1[N],x2[N];
ll a[N/2],b[N/2];
ll sum[N]; // 结果存在sum里
void brc(complex *y,int l) // 二进制平摊反转置换 O(logn)
{
	register int i,j,k;
	for(i=1,j=l/2;i<l-1;i++)
	{
		if(i<j)	swap(y[i],y[j]); // 交换互为下标反转的元素
								// i<j保证只交换一次
		k=l/2;
		while(j>=k) // 由最高位检索，遇1变0，遇0变1，跳出
		{
			j-=k;
			k/=2;
		}
		if(j<k)	j+=k;
	}
}
void fft(complex *y,int l,double on) // FFT O(nlogn)
				     		// 其中on==1时为DFT，on==-1为IDFT
{
	register int h,i,j,k;
	complex u,t;
	brc(y,l); // 调用反转置换
	for(h=2;h<=l;h<<=1) // 控制层数
	{
		// 初始化单位复根
		complex wn(cos(-on*2*pi/h),sin(-on*2*pi/h));
		for(j=0;j<l;j+=h) // 控制起始下标
		{
			complex w(1,0); // 初始化螺旋因子
			for(k=j;k<j+h/2;k++) // 配对
			{
				u=y[k];
				t=w*y[k+h/2];
				y[k]=u+t;
				y[k+h/2]=u-t;
				w=w*wn; // 更新螺旋因子
			} // 据说上面的操作叫蝴蝶操作…
		}
	}
	if(on==-1)	for(i=0;i<l;i++)	y[i].r/=l; // IDFT
}
void init(){
     memset(sum,0,sizeof(sum));
     memset(x1,0,sizeof(x1));
     memset(x2,0,sizeof(x2));
}
int main()
{
    int T;
    scanf("%d",&T);
    ll cnt1=0,cnt2=0;
    while(T--)
    {
        int n;
        cnt1=cnt2=0;
        init();
        scanf("%d",&n);
        for(int i = 0;i < n;i++) {scanf("%lld",&a[i]);cnt1+=a[i]*a[i];}
        for(int i = 0;i < n;i++) {scanf("%lld",&b[i]);cnt2+=b[i]*b[i];}
        int len = 1;
        while( len < 2*n ) len <<= 1;
        for(int i = 0;i < n;i++){
                x1[i] = complex(a[i],0);
            }
             for(int i = 0;i < n;i++){
                x2[i] = complex(b[n-i-1],0);
            }
      //      for(int i=n;i<len;i++) x1[i]=Complex(0,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 < len;i++){
                sum[i] = (LL)(x1[i].r+0.5);
            }
          //  for(int i = 0;i < len;i++) cout<<num[i]<<endl;
          LL ret=sum[n-1];
          int flag=0;
         // cout<<ret<<endl;
          for(int i=0;i<n-2;i++) {
            //    cout<<num[i]+num[i+n]<<endl;
            if(ret<sum[i]+sum[i+n])
                {ret=sum[i]+sum[i+n]; flag=n-1-i;}
                  //注意，此时得到的ret会有很小的浮点精度误差,
                //flag表示k,这个是正确的
          }
          ret=0;
          for(int i=0;i<n;i++){
            ret+=a[i]*b[(i+flag)%n]; //重新算一遍得到最后答案
          }
          printf("%lld\n",cnt1+cnt2-2*ret);
        }
    return 0;
}