hihocoder 1388 Periodic Signal

ACM-ICPC国际大学生程序设计竞赛北京赛区(2016)网络赛 1006 Periodic Signal

传送门

题目大意

给你两个序列，求一个差值平方和最小
具体求法就是
min{(a1-b1)²+...+(an-bn)²,(a1-b2)²+...+(an-b1)²，...，(a1-bn)²+...+(an-b1)²}

知道题意后化简的
$\sum_{i=1}^{n}a_{i}^{2}+\sum_{i=1}^{n}b_{i}^{2}-2(a_{i}*b{j})(i\epsilon n,j\epsilon n)$
现在问题就变成了怎么快速求出$a_{i}*b{j}$ 的最大值
普通方法的时间复杂度O($n^{2}$),肯定要超时。
考虑下面的循环矩阵

$$\begin{bmatrix} a_{1}&a_{2}&a_{3}&...&a_{n}\\ a_{n}&a_{1}&a_{2}&...&a_{n-1}\\ a_{n-1}&a_{n}&a_{1}&...&a_{n-2}\\ .&.&.&...&.\\ a_{2} &a_{3}&a_{4}&...&a_{1} \end{bmatrix}*\begin{bmatrix} b_{1}\\ b_{2}\\ b_{3}\\ .\\ b_{n} \end{bmatrix}= \begin{bmatrix} c_{1}\\ c_{2}\\ c_{3}\\ .\\ c_{n} \end{bmatrix}$$
只需要求 $c_{n}$ 的最大值即可。
可以参考下面一个链接
http://www.zybang.com/question/dd8b336ad690b3e2f9aa4dc0b596e1ea.html
将循环矩阵*向量转换为求卷积
所以用快速傅里叶变换将时间复杂度降为了O( $nlog_{2}n$)
具体做法是
令 $\overrightarrow{a}$=( $a_{1},a_{2}...,a_{n}$)
令 $\overrightarrow{b}$=( $b_{n},b_{n-1},...b_{2},b_{1},b_{n},b_{n-1}...b_{2}$)
之后对 $\overrightarrow{a}$ 求 $\overrightarrow{c}$=DFT( $\overrightarrow{a}$)
对 $\overrightarrow{b}$ 求 $\overrightarrow{d}$= DFT( $\overrightarrow{b}$)
时间复杂度为O( $nlog_{2}n$)
然后求 $\overrightarrow{S}=\overrightarrow{c}*\overrightarrow{d}$
时间复杂度为O(n)
再对 $\overrightarrow{S}$ 求 $\overrightarrow{ans}$=IDFT( $\overrightarrow{S}$)
时间复杂度为O( $nlog_{2}n$)
最后求 $\overrightarrow{ans}$ 中最大的元素即可
不过这道题的坑点在于精度问题
WA了若干次
发现精度怎么都够不了
于是发现在 $\overrightarrow{ans}$ 中答案的位置就是移动的位置
于是先求出移动的位置
然后暴力跑一边就求出答案了
代码如下(缩进懒得改了套的模板，将就看)

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#define N 150005
#define pi acos(-1.0)
using namespace std;
int l,T,n,num,wz;
long long ans,tot,sum[N],A[N],B[N];
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];
void brc(complex *y,int l)
{  
    register int i,j,k;
    for(i=1,j=l/2;i<l-1;i++)  
    {
        if(i<j)swap(y[i],y[j]);
        k=l/2;
        while(j>=k)j-=k,k/=2;
        if(j<k)j+=k;
    }
}
void fft(complex *y,int l,double on)
{
    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;
}
long long Getans(int k)
{
    ans=0;
    for(int i=0;i<n;i++)ans+=(A[i]-B[(i+k)%n])*(A[i]-B[(i+k)%n]);
    return ans;
}
int main()  
{
    register int i,j;
    scanf("%d",&T);
    while(T--)
    {
        l=1;
        tot=0;
        scanf("%d",&n);
        while(l<n*2)l<<=1;
        for(i=n-1;i>=0;i--)
        {
            scanf("%d",&num);
            x1[i].r=num;
            x1[i].i=0.0;
            tot+=num*num;
            A[n-i-1]=num;
        }
        for(i=n;i<2*n-1;i++)x1[i]=x1[i-n];
        for(i=2*n-1;i<l;i++)x1[i].r=x1[i].i=0.0;
        for(i=0;i<n;i++)
        {
            scanf("%d",&num);
            x2[i].r=num;
            x2[i].i=0.0;
            tot+=num*num;
            B[i]=num;
        }
        for(i=n;i<l;i++)x2[i].r=x2[i].i=0.0;  
        fft(x1,l,1);//DFT
        fft(x2,l,1);//DFT
        for(i=0;i<l;i++)x1[i]=x1[i]*x2[i];
        fft(x1,l,-1);//IDFT
        double m=-1;
        for(i=n-1;i<2*n;i++)if(m<x1[i].r)m=x1[i].r,wz=i;//求位置
        printf("%lld\n",Getans((wz+1)%n));
    }
    return 0;
}