题目描述:
Profess X is an expert in signal processing. He has a device which can send a particular 1 second signal repeatedly. The signal is A0 ...An-1 under n Hz sampling.
One day, the device fell on the ground accidentally. Profess X wanted to check whether the device can still work properly. So he ran another n Hz sampling to the fallen device and got B0 ...Bn-1.
To compare two periodic signals, Profess X define the DIFFERENCE of signal A and B as follow:
You may assume that two signals are the same if their DIFFERENCE is small enough.
Profess X is too busy to calculate this value. So the calculation is on you.
输入:
The first line contains a single integer T, indicating the number of test cases.
In each test case, the first line contains an integer n. The second line contains n integers, A0 ...An-1. The third line contains n integers, B0 ... Bn-1.
T≤40 including several small test cases and no more than 4 large test cases.
For small test cases, 0<n≤6⋅103.
For large test cases, 0<n≤6⋅104.
For all test cases, 0≤Ai,Bi<220.
输出:
For each test case, print the answer in a single line.
样例输入:
2
9
3 0 1 4 1 5 9 2 6
5 3 5 8 9 7 9 3 2
5
1 2 3 4 5
2 3 4 5 1
样例输出:
80
0
题目大意, 呃.....
给你举个例子吧, 比如给你两组数列A和b, 都是n个数, 这里我们简单点描述, n取4, 编号从0开始, 到n-1结束.
然后有以下几种组合, 问你这四种组合里面最小的组合的值是多少.
组合0: ( A0 - B0 )^2 + ( A1 - B1 )^2 + ( A2 - B2 )^2 + ( A3 - B3 )^2 (对应原题面的k=0)
组合1: ( A0 - B1 )^2 + ( A1 - B2 )^2 + ( A2 - B3 )^2 + ( A3 - B0 )^2 (对应原题面的k=1)
组合2: ( A0 - B2 )^2 + ( A1 - B3 )^2 + ( A2 - B0 )^2 + ( A3 - B1 )^2 (对应原题面的k=2)
组合3: ( A0 - B3 )^2 + ( A1 - B0 )^2 + ( A2 - B1 )^2 + ( A3 - B2 )^2 (对应原题面的k=3, 也就是最大的k=n-1=3)
如果n取任意数的话, 用形象一点的描述的话, 就是, 固定住A的各项, 然后不断地 "滚动" B的各项.
就像上面举的例子, 一开始B的各项的顺序是0123, 滚动一下变成1230, 在滚动一下变成2301, 再滚动一下变成3012, 再滚就滚回0123了.
现在我们以5为例讲一下这个题是怎么跟FFT扯上关系的吧.
n=5, 可以 "滚" 5次,那就是5个组合.
组合0: ( A0 - B0 )^2 + ( A1 - B1 )^2 + ( A2 - B2 )^2 + ( A3 - B3 )^2 + ( A4 - B4 )^2 (对应原题面的k=0)
组合1: ( A0 - B1 )^2 + ( A1 - B2 )^2 + ( A2 - B3 )^2 + ( A3 - B4 )^2 + ( A3 - B0 )^2 (对应原题面的k=1)
组合2: ( A0 - B2 )^2 + ( A1 - B3 )^2 + ( A2 - B4 )^2 + ( A3 - B0 )^2 + ( A3 - B1 )^2 (对应原题面的k=2)
组合3: ( A0 - B3 )^2 + ( A1 - B4 )^2 + ( A2 - B0 )^2 + ( A3 - B1 )^2 + ( A3 - B2 )^2 (对应原题面的k=3)
组合4: ( A0 - B4 )^2 + ( A1 - B0 )^2 + ( A2 - B1 )^2 + ( A3 - B2 )^2 + ( A3 - B3 )^2 (对应原题面的k=4, 也就是最大的k=n-1=4)
如果打开的的话, 那么会发现这五种组合都有共同的10项, 那就是
A0^2 + A1^2 + A2^2 + A3^2 + A4^2 + B0^2 + B1^2 + B2^2 + B3^2 + B4^2
那这个时候就相当于比较这些东西谁比较大了(原本的式子减去公共的10项再取个负就是下面的东西)
组合0' : A0*B0 *2 + A1*B1 *2 + A2*B2 *2 + A3*B3 *2 + A4*B4 *2
组合1' : A0*B1 *2 + A1*B2 *2 + A2*B3 *2 + A3*B4 *2 + A4*B0 *2
组合2' : A0*B2 *2 + A1*B3 *2 + A2*B4 *2 + A3*B0 *2 + A4*B1 *2
组合3' : A0*B3 *2 + A1*B4 *2 + A2*B0 *2 + A3*B1 *2 + A4*B2 *2
组合4' : A0*B4 *2 + A1*B0 *2 + A2*B1 *2 + A3*B2 *2 + A4*B3 *2
这些东西怎么算呢?
上图!(主要是电脑的公式编辑器我不怎么会用, 所以就干脆手写了)(加了滤镜把一些有的没的的东西去掉了)
啊对了, 我字丑, 下面的十六个字是 项的次数 对应系数 项的次数 对应系数
有没有看出点门道?
组合0' 正好是x^4的系数和x^9的系数之和 恰巧9=5+4
组合1' 正好是x^3的系数和x^8的系数之和 恰巧8=5+3
组合2' 正好是x^2的系数和x^7的系数之和 恰巧7=5+2
组合3' 正好是x^1的系数和x^6的系数之和 恰巧6=5+1
组合4' 正好是x^0的系数和x^5的系数之和 恰巧5=5+0
还不会找规律我觉得你不适合写这个.....
而图里面的最上面的多项式相乘, 可以用FFT进行优化, 找到之后生成的生成函数系数数组就是图上面的东西, 挨个儿看过去哪个x^i的系数和x^(n+i)的系数之和最大就行, 那我们也找到了最大的组合(n-1-i)', 也就找到了最小的组合(n-1-i).
比如我发现x^2的系数和x^(n+2)的系数之和最大, 那就是组合(n-1-2)'是最大的, 那么组合(n-1-2)就是最小的.
行了 , 比如我找到了某个组合是最大的, 但是, 我要告诫你, 由于某些原因, 如果你直接把找到的最大的组合某撇的值取负然后加上A B两个数组各项平方之和的当结果的话, 你一定会WA, 因为精度不够.
那重算就行啦, 不知道你有没有意识到我刻意设置组合标号的含义---"滚"B数组的次数, 这个时候不难了吧, 比如我想求滚了num次之后的组合对应的值, A[i]要减掉的是B[(i+num)%n] 这个模可以理解成从末尾又滚回去了.
然后求n下和, 就是结果
代码如下:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <algorithm>
#include <iostream>
using namespace std;
const double PI=acos(-1);//这个是一定要有的, 因为fft有应用到三角函数,因此头文件math.h也要包含进去
const int maxx=300000;
struct complex
{
double a,b;//a表示实部,b表示虚部
complex(double aa=0.0,double bb=0.0)
{
a=aa;
b=bb;
}
complex operator +(const complex &e)
{//重载运算符+用以直接进行虚数相加,下同
return complex(a+e.a,b+e.b);
}
complex operator -(const complex &e)
{//重载运算符-
return complex(a-e.a,b-e.b);
}
complex operator *(const complex &e)
{//重载运算符*
return complex(a*e.a-b*e.b,a*e.b+b*e.a);
}
};
void change(complex *y,long long int len)
{//change这个词还挺常用于给自己写的功能函数重名的,要注意一下,不过有可能以后更新模板就换掉了
long long int i,j,k;
for(i=1,j=len/2;i<len-1; i++)
{
if(i<j)
swap(y[i],y[j]);
k=len/2;
while(j>=k)
{
j-=k;
k>>=1;
}
if(j<k)
{
j+=k;
}
}
}
void fft(complex *y,long long int len,int on)
{
change(y, len);
int i,j,k;
for(i=2;i<=len;i<<=1)
{
complex wn(cos(-on*2*PI/i),sin(-on*2*PI/i));
for(j=0;j<len;j+=i)
{
complex w(1,0);
for(k=j;k<j+i/2;k++)
{
complex u=y[k],t=w*y[k+i/2];
y[k]=u+t;
y[k+i/2]=u-t;
w=w*wn;
}
}
}
if(on==-1)//还原回原来的多项式
for(i=0;i<len;i++)
y[i].a/=len;
}
complex a[maxx],b[maxx],ans[maxx];
long long int A[100000],B[100000],ab[maxx];
int main()
{
long long int T,i,n,len_a,len,num,res;
scanf("%lld",&T);
while(T--)
{
memset(a,0,sizeof(a));
memset(b,0,sizeof(b));
memset(A,0,sizeof(A));
memset(B,0,sizeof(B));
memset(ab,0,sizeof(ab));
memset(ans,0,sizeof(ans));
scanf("%lld",&n);
for(i=0;i<n;i++)
scanf("%lld",&A[i]);
for(i=0;i<n;i++)
scanf("%lld",&B[n-1-i]);//注意这里是倒置哦!
len=1;
len_a=n;
while(len<2*len_a)
len=len<<1;
for(i=0;i<n;i++)
{
a[i].a=A[i];
b[i].a=B[i];
}
fft(a,len,1);
fft(b,len,1);
for(i=0;i<len;i++)
ans[i]=a[i]*b[i];
fft(ans,len,-1);
for(i=0;i<len;i++)
ab[i]=(long long int)(ans[i].a+0.5);
res=0;
num=0;
for(i=0;i<n;i++)
{
if(res<(ab[i]+ab[i+n]))//为什么这么写参见上面的解释
{
res=ab[i]+ab[i+n];
num=n-1-i;
}
}
res=0;
for(i=0;i<n;i++)
res+=(A[i]-B[n-1-((i+num)%n)])*(A[i]-B[n-1-((i+num)%n)]);//因为是倒置输入的,取B[(i+num)%n]实际上在代码里面得写B[n-1-((i+num)%n)]
printf("%lld\n",res);
}
return 0;
}