首页 > 其他分享 >【HNOI2017】礼物(FFT)

【HNOI2017】礼物(FFT)

时间:2022-10-28 20:57:00浏览次数:42  
标签:sum FFT widehat leq Complex HNOI2017 2n PN 礼物

显然,\(y_i\) 加上 \(c\) 可以看成是 \(x_i\) 减去 \(c\)。

所以就变成了 \(x_i\) 加上一个整数(可正可负)。

现将 \(x\) 环拆成一个长度为 \(2n\) 的序列 \(a\)(复制一遍),把 \(y\) 环拆成一个长度为 \(n\) 的序列 \(b\)。

那么旋转操作就可以看成是 \(b\) 序列与 \(a\) 序列中每一个长度为 \(n\) 的子串匹配求值。

也就是说,求这个东西的最小值:\(\sum_{j=0}^{n-1}(a_{i+j}-b_j+c)^2\)(\(0\leq i< n\))。

接下来推式子:(设 \(x\) 环的和是 \(suma\),\(y\) 环的和是 \(sumb\),\(x\) 环的平方和是 \(powa\),\(y\) 环的平方和是 \(powa\))

\[\begin{aligned} &\min_{i=0}^{n-1} \sum_{j=0}^{n-1}(a_{i+j}-b_j+c)^2\\ =&\min_{i=0}^{n-1}\sum_{j=0}^{n-1}(a_{i+j}-b_j)^2+c^2+2(a_{i+j}-b_j)c\\ =&\min_{i=0}^{n-1} [nc^2+\sum_{j=0}^{n-1}(a_{i+j}-b_j)^2+2c\sum_{j=0}^{n-1}(a_{i+j}-b_j)]\\ =&\min_{i=0}^{n-1} [nc^2+\sum_{j=0}^{n-1}(a_{i+j}^2+b_j^2-2a_{i+j}b_j)+2c(\sum_{j=0}^{n-1} a_{i+j}-\sum_{j=0}^{n-1} b_j)]\\ =&\min_{i=0}^{n-1} (nc^2+powa+powb-2\sum_{j=0}^{n-1}a_{i+j}b_j+2c(suma-sumb)) \end{aligned}\]

发现 \(powa+powb\) 是定值,而和 \(c\) 有关的 \(nc^2+2c(suma-sumb)\) 可以用二次函数最值 \(O(1)\) 求,也就是说我们只需要求 \(-2\sum_{j=0}^{n-1}a_{i+j}b_j\) 的最小值,即 \(\sum_{j=0}^{n-1}a_{i+j}b_j\) 的最大值。

我们可以把这个形式改变一下,把 \(\sum_{j=0}^{n-1}a_{i+j}b_j\) 改成 \(\sum_{i-j=k}a_ib_j\)(\(k\) 是旋转角度)。

令 \(S(k)=\sum_{i-j=k}a_ib_j\),设 \(\widehat{a_i}=a_{2n-i}\),则:

\[S(k)=\sum_{i-j=k}^{0\leq j< n}a_ib_j=\sum_{(2n-i)-j=k}^{0\leq j< n}\widehat{a_i}b_j=\sum_{i+j=2n-k}^{0\leq j< n}\widehat{a_i}b_j \]

把 \(n-k\) 代入:

\[S(2n-k)=\sum_{i+j=2n-(2n-k)}^{0\leq j< n}\widehat{a_i}b_j=\sum_{i+j=k}^{0\leq j< n}\widehat{a_i}b_j \]

有木有觉得很熟悉?

我们可以把 \(\widehat{a_i}\) 看成一个多项式 \(A\) 的系数,\(b_j\) 看成另一个多项式 \(B\) 的系数,那么 \(S(2n-k)\) 就是 \(A\times B\) 第 \(k\) 项的系数。

看会我们原来的问题:求 \(\sum_{i-j=k}a_ib_j\) 的最大值,也就是求 \(S(k)\) 的最大值(\(0\leq k<n\)),也就是 \(A\times B\) 第 \((n+1)\sim 2n\) 项系数的最大值。

那么现在就好处理了,我们先用 FFT 把 \(A\times B\) 算出来,再取最大值。

代码如下:

#include<bits/stdc++.h>

#define N 50010
#define PN 262200
#define INF 0x7fffffff

using namespace std;

struct Complex
{
	double x,y;
	Complex(){};
	Complex(double xx,double yy){x=xx,y=yy;}
}a[PN],b[PN],c[PN];

Complex operator + (Complex a,Complex b){return Complex(a.x+b.x,a.y+b.y);}
Complex operator - (Complex a,Complex b){return Complex(a.x-b.x,a.y-b.y);}
Complex operator * (Complex a,Complex b){return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}

const double pi=acos(-1);

int n,m,suma,sumb,powa,powb;
int limit=1,bit,rev[PN];

void FFT(Complex *a,int opt)
{
	for(int i=0;i<limit;i++)
		if(i<rev[i])
			swap(a[i],a[rev[i]]);
	for(int mid=1;mid<limit;mid<<=1)
	{
		Complex wn=Complex(cos(pi/mid),opt*sin(pi/mid));
		for(int i=0,len=(mid<<1);i<limit;i+=len)
		{
			Complex now=Complex(1,0);
			for(int j=0;j<mid;j++,now=now*wn)
			{
				Complex x=a[i+j],y=now*a[i+mid+j];
				a[i+j]=x+y,a[i+mid+j]=x-y;
			}
		}
	}
	if(opt==-1)
		for(int i=0;i<limit;i++)
			a[i].x/=limit;
}

int main()
{
	scanf("%d%d",&n,&m);
	for(int i=0;i<n;i++)
	{
		scanf("%lf",&a[i].x);
		a[n+i].x=a[i].x;
		suma+=a[i].x;
		powa+=a[i].x*a[i].x;
	}
	for(int i=0;i<n;i++)
	{
		scanf("%lf",&b[i].x);
		sumb+=b[i].x;
		powb+=b[i].x*b[i].x;
	}
	int d=round(-1.0*(suma-sumb)/n);
	reverse(a,a+2*n+1);//求a^
	while(limit<=3*n)
		limit<<=1,bit++;
	for(int i=0;i<limit;i++)
		rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
	FFT(a,1),FFT(b,1);
	for(int i=0;i<limit;i++)
		c[i]=a[i]*b[i];
	FFT(c,-1);
	int ans=-INF;
	for(int i=n+1;i<=2*n;i++)
		ans=max(ans,(int)round(c[i].x));//这里用round是怕被卡精度
	printf("%d\n",n*d*d+2*d*(suma-sumb)+powa+powb-2*ans);
	return 0;
}

标签:sum,FFT,widehat,leq,Complex,HNOI2017,2n,PN,礼物
From: https://www.cnblogs.com/ez-lcw/p/16837452.html

相关文章

  • 【CF553E】Kyoya and Train(期望dp,SPFA,FFT)
    考虑dp。发现正着dp好像不太好做,毕竟初值不太好设,而且时间一大于\(t\)费用就要加上\(x\),所以考虑倒着dp。设\(f_{u,j}\)表示现在已经到达\(u\)点,耗时\(j\),问......
  • fft
    重新改了一下板子抄起来比较短#include<bits/stdc++.h>usingnamespacestd;#definerep(i,h,t)for(inti=h;i<=t;i++)#definedep(i,t,h)for(inti=t;i>=h;i--......
  • (算法课)大整数相乘 |向量卷积|多项式相乘| 快速傅里叶变换FFT
    D(1021):很大的ABTimeLimit:1SecMemoryLimit:256MbSubmitted:6Solved:0Description如题,计算AB的值并输出。Input两行,分别代表A和B。A和B......
  • 天线阵列fft角度映射
    相位差:Δφ=2πdsin(θ)/λ=2π/N其中,d是天线间距,θ是入射角,λ是载频波长λ=f0/c,N是沿着天线阵列维度(例如:86根天线)的fft采样点数。将角度维fft映射到[0,2π]的相位,得到......
  • Python | 用Python制作送给女票的生日礼物
    HappyBirthday程序视频地址:​​https://www.bilibili.com/video/BV1R7411C7A1​​代码地址:​​https://github.com/borninfreedom/HappyBirthday​​截图: ......
  • EXCEL 计算FFT步骤
    原文链接:https://blog.csdn.net/WaliTool/article/details/1143696791,数据产生利用Excel模拟出一系列数据(本例子产生1024个数据)公式为:y=1.5sin(50∗2π102......
  • 多项式,FFT
    多项式,FFTPolynomialConvolution若$A(x)$与$B(x)$分别为两个多项式,则两个多项式的卷积为:$$A(x)\cdotB(x)=\sum_{i=0}n\sum_{j=0}ia_jb_{i-j}$$其中,$A(x)$......
  • P1194 买礼物
    P1194买礼物普及的题目,而且一眼就能看出该用什么做法。我主要是决定这道题建图的思想值得借鉴,每样东西原本的价格是a,所以新建一个节点0,0向i连边,边权为a,这样一共就有b......
  • 送爸妈最好的礼物,让爱与孝心一起归家
    还记得去年春节有一位儿子因为为父母手绘了一本“微信使用说明书”而感动了无数网友吗?他用自己的绘画特长,为父母专门手绘了一本“微信使用说明书”,几页图画配发详尽的文字,教......
  • 【STM32H7的DSP教程】第28章 FFT和IFFT的Matlab实现(幅频响应和相频响应)
    ​​​​第28章      FFT和IFFT的Matlab实现(幅频响应和相频响应)本章主要讲解fft,ifft和fftshift在matlab上的实现。28.1初学者重要提示28.2Matlab的FFT函数28.3Mat......