【hihocoder#1388】Periodic Signal NTT

题目链接:http://hihocoder.com/problemset/problem/1388?sid=974337

题目大意:找出一个$k$,使得$\sum_{i=0}^{n-1}(A_{i}-B_{i+k \quad mod \quad n})^{2}$最小 

 

把那个式子拆开..得到:

$\sum _{i=0}^{n-1} A_{i}^{2}+\sum_{i=0}^{n-1}B_{i}^{2}-\sum_{i=0}^{n-1}2*A_{i}*B_{i+k}$

前面的两个和式是定值,所以最小化后面的值,后面这个是个卷积的形式,所以把$B$翻转硬上就好了...

 

然后问题就来了..FFT了之后貌似精度爆炸了..然后就写NTT..然后又炸了..发现模数规模不够..于是查表换了个巨大的模数,然后又炸了..

然后去看题解,意识到一个严重的问题..相乘爆longlong...然后又学习了TA爷的$O(1)$快速乘..

附上一张神表:折跃

 

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
using namespace std;
#define LL long long
inline LL read()
{
	LL x=0,f=1; char ch=getchar();
	while (ch<'0' || ch>'9') {if (ch=='-') f=-1; ch=getchar();}
	while (ch>='0' && ch<='9') {x=x*10+ch-'0'; ch=getchar();}
	return x*f;
}
#define P (LL)((29LL<<57)+1)
#define G 3
#define MAXN 800100
int T,N,len;
LL Wn[30],A[MAXN],B[MAXN],C[MAXN],wn[MAXN],a[MAXN],b[MAXN];
inline LL Mul(LL x,LL y) {return (x*y-(LL)(x/(long double)P*y+1e-3)*P+P)%P;}
inline LL Pow(LL x,LL y)
{
    LL re=1;
    for (LL i=y; i; i>>=1,x=Mul(x,x))
        if (i&1) re=Mul(re,x);
    return re;
}
inline LL Inv(LL x) {return Pow(x,P-2);}
inline void Prework()
{
    len=1;
    while (len<(N<<1)) len<<=1;
    for (int i=0; i<N; i++) A[i]=a[i];
    for (int i=0; i<N; i++) B[i]=b[N-1-i];
    for (int i=N; i<len; i++) A[i]=0;
    for (int i=N; i<len; i++) B[i]=0;
//    for (int i=0; i<len; i++) printf("%I64d ",A[i]); puts("");
//    for (int i=0; i<len; i++) printf("%I64d ",B[i]); puts("");
    for (int i=0; i<=28; i++) wn[i]=Pow(G,(P-1)/(1<<i));
}
inline void Rader(LL *x)
{
    for (int i=1,j=len>>1,k; i<len-1; i++)
        {
            if (i<j) swap(x[i],x[j]);
            k=len>>1;
            while (j>=k) j-=k,k>>=1;
            if (j<k) j+=k;
        }
}
inline void DFT(LL *x,int opt)
{
    Rader(x);
    for (int h=2,id=0; h<=len; h<<=1)
        {
            LL Wn=wn[++id];
            for (int i=0; i<len; i+=h)
                {
                    LL W=1;
                    for (int j=i; j<i+h/2; j++)
                        {
                            LL u=x[j]%P,t=Mul(W,x[j+h/2]);
                            x[j]=(u+t)%P; x[j+h/2]=(u-t+P)%P;
                            W=Mul(W,Wn);
                        }
                }
        }
    if (opt==-1)
        {
            for (int i=1; i<len/2; i++) swap(x[i],x[len-i]);
            for (int i=0; i<len; i++)
                x[i]=Mul(x[i],Inv(len));
        }
}
inline void NTT(LL *A,LL *B)
{
    DFT(A,1); DFT(B,1);
    for (int i=0; i<len; i++) C[i]=Mul(A[i],B[i]);
    DFT(C,-1);
}
int main()
{
	T=read();
	while (T--)
		{
			N=read();
			LL ans=0;
			for (int i=0; i<=N-1; i++) a[i]=read();
			for (int i=0; i<=N-1; i++) b[i]=read();			
			for (int i=0; i<=N-1; i++) ans+=(LL)a[i]*a[i];
			for (int i=0; i<=N-1; i++) ans+=(LL)b[i]*b[i];
			Prework();
			NTT(A,B);
			LL mx=C[N-1];
			for (int i=0; i<N-1; i++) mx=max(mx,C[i]+C[i+N]);
//			for (int i=0; i<len; i++) printf("%I64d ",C[i]); puts("");
			printf("%lld\n",ans-2*mx); 
		}
	return 0;
} 

  

posted @ 2017-01-09 21:30  DaD3zZ  阅读(211)  评论(0编辑  收藏  举报