zoukankan      html  css  js  c++  java
  • [Luogu P3723] [AH2017/HNOI2017]礼物 (FFT 卷积)

    题面

    传送门:洛咕


    Solution

    调得我头大,我好菜啊

    好吧,我们来颓柿子吧:
    我们可以只旋转其中一个手环。对于亮度的问题,因为可以在两个串上增加亮度,我们也可以看做是可以为负数的。

    所以说,我们可以假设我们旋转(B)串,上下要加上的亮度差为(p),可以直接拍出一个最暴力的柿子:
    (f(x))表示(B)串以(x)为开头的差异值,有:
    (f(x)=sum_{i=0}^{x-1}(B[i]-A[i+n-x]+p)^2+sum_{i=x}^{n-1}(B[i]-A[i-x]+p)^2)

    大力展开化简后有:
    (f(x)=sum_{i=0}^{n-1}A[i]^2+sum_{i=0}^{n-1}B[i]^2+n*p^2-2psum_{i=0}^{n-1}(A[i]-B[i])-2sum_{i=0}^{x-1}(B[i]*A[i+n-x])-2sum_{i=x}^{n-1}(B[i]*A[i-x]))

    前两项(sum_{i=0}^{n-1}A[i]^2+sum_{i=0}^{n-1}B[i]^2)显然(O(n)预处理出来)
    中间两项(n*p^2-2psum_{i=0}^{n-1}(A[i]-B[i]))是一个关于(p)的二次函数,我们找最小值就好。(因为这题(m)非常小,我们也可以暴力枚举),复杂度(O(1))(O(m))
    最后两项(-2sum_{i=0}^{x-1}(B[i]*A[i+n-x])-2sum_{i=x}^{n-1}(B[i]*A[i-x]))看起来非常像卷积,但是并不是,因此我们得做点处♂理。

    蒟蒻本人是这样处理的:
    首先,后面那个循环范围是肯定没法卷的,因此我们先把后面的循环处理一下得:
    (-2sum_{i=0}^{x-1}(B[i]*A[i+n-x])-2sum_{i=0}^{n-x-1}(A[i]*B[i+x]))
    然后我们可以考虑把前面那项的(A)反转(这样可以处理掉(n)来方便卷积),把后面那项的(B)反转(这样可以制造(n)(sum)对应)
    (-2sum_{i=0}^{x-1}(B[i]*A'[x-1-i])-2sum_{i=0}^{n-x-1}(A[i]*B'[n-1-i-x]))

    哦豁,卷积,搞定。
    时间复杂度(O(n*logn))


    Code

    我什么时候才能一次性写对FFT啊

    //Luogu P3723 [AH2017/HNOI2017]礼物
    //Jan,20th,2019
    //颓柿子+FFT加速计算
    #include<iostream>
    #include<cstdio>
    #include<cmath>
    #include<complex>
    #include<algorithm>
    using namespace std;
    long long read()
    {
    	long long x=0,f=1; char c=getchar();
    	while(!isdigit(c)){if(c=='-') f=-1;c=getchar();}
    	while(isdigit(c)){x=x*10+c-'0';c=getchar();}
    	return x*f;
    }
    const int M=50000+100;
    const int N=M*4;
    const double PI=acos(-1);
    typedef complex <double> cp;
    inline cp omega(int K,int n)
    {
    	return cp(cos(2*PI*K/n),sin(2*PI*K/n));
    }
    void FFT(cp a[],int n,bool type)
    {
    	static int tmp[N],num=n-1,len=0;
    	while(num!=0) num/=2,len++;
    	for(int i=0,j;i<=n;i++)
    	{
    		for(j=0,num=i;j<len;j++)
    			tmp[j]=num%2,num/=2;
    		reverse(tmp,tmp+len);
    		for(j=0,num=0;j<len;j++)
    			num+=tmp[j]*(1<<j);
    		if(i<num) swap(a[i],a[num]);
    	}
    	for(int l=2;l<=n;l*=2)
    	{
    		cp x0=omega(1,l);
    		if(type==true) x0=conj(x0);
    		int m=l/2;
    		for(int j=0;j<n;j+=l)
    		{
    			cp x(1,0);
    			for(int k=0;k<m;k++,x*=x0)
    			{
    				cp temp=x*a[j+k+m];
    				a[j+k+m]=a[j+k]-temp;
    				a[j+k]=a[j+k]+temp;
    			}
    		}
    	}
    }
    int n,m,a[N],b[N];
    cp B1[N],A1[N],B2[N],A2[N];
    long long ans;
    int main()
    {
    	freopen("3723.in","r",stdin);
    	
    	n=read(),m=read();
    	for(int i=0;i<n;i++)
    		a[i]=read();
    	for(int i=0;i<n;i++)
    		b[i]=read();
    	
    	long long dif=0;
    	for(int i=0;i<n;i++)
    		ans+=a[i]*a[i]+b[i]*b[i],dif+=a[i]-b[i];
    	long long t_ans=0x3f3f3f3f*0x3f3f3f3f;
    	for(int i=-m;i<=m;i++)
    		t_ans=min(t_ans,n*i*i-2*i*dif);
    	ans+=t_ans;
    	
    	int t=1;
    	while(t<2*n) t*=2;
    	reverse(a,a+n);
    	for(int i=0;i<n;i++)
    		A1[i]=a[i],B1[i]=b[i];
    	FFT(A1,t,false);
    	FFT(B1,t,false);
    	for(int i=0;i<t;i++)
    		A1[i]*=B1[i];
    	FFT(A1,t,true);
    	for(int i=0;i<t;i++)
    		A1[i].real()/=t;
    		
    	reverse(a,a+n);
    	reverse(b,b+n);
    	for(int i=0;i<n;i++)
    		A2[i]=a[i],B2[i]=b[i];
    	FFT(A2,t,false);
    	FFT(B2,t,false);
    	for(int i=0;i<t;i++)
    		A2[i]*=B2[i];
    	FFT(A2,t,true);
    	for(int i=0;i<t;i++)
    		A2[i].real()/=t;
    	
    	t_ans=(long long)(2*floor(A2[n-1].real()+0.5));
    	for(int i=1;i<n;i++)
    		t_ans=max(t_ans,(long long)(2*floor(A1[i-1].real()+A2[n-i-1].real()+0.5)));
    	ans-=t_ans;
    	printf("%lld",ans);
    	return 0;
    }
    
    
  • 相关阅读:
    OZ Report 오즈 리포트 개발
    日期 英文 英语 韩文 韩语
    行合并
    R内存扩展 win7内存扩展
    计算日期函数
    C#颜色 Color.FromArgb ColorTranslator 16进制
    Phone
    浏览器 UserAgent
    context.Request.Files post 上传问题件
    The connection to adb is down, and a severe error has occured.
  • 原文地址:https://www.cnblogs.com/GoldenPotato/p/10296961.html
Copyright © 2011-2022 走看看