[BJOI2020] 主音识别

最新推荐文章于 2021-03-18 18:58:07 发布

原创最新推荐文章于 2021-03-18 18:58:07 发布 · 287 阅读

0 ·

CC 4.0 BY-SA版权

多项式专栏收录该内容

1 篇文章

订阅专栏

本文介绍了一种使用快速傅立叶变换(FFT)优化的算法，解决了一个关于求解最小g(x)的问题。该问题涉及到集合和数列的模运算，并通过建立多项式进行卷积计算来找到最优解。

~~咕咕咕的北京省选题~~

description

给定参数 $k$ ,一个大小为 $m$ 的集合 ${b_1,b_2,b_3...b_m\}$ ，一个长度为 $n$ 的数列 ${a_i\}$

定义 $k)]f(x,a_i)=[a_i\equiv x(\bmod k)]+[\exist1\leq j\leq m,a_i\equiv x+b_j(\bmod k)]$

通俗来讲：
若 $x$ 与 $a_i$ 模 $k$ 同余，则 $f(x,a_i)$ 为 $2$ ；
否则，若 $B$ 中存在一个元素 $b_j$ ，使得 $x+b_j$ 与 $a_i$ 模 $k$ 同余，则 $f(x,a_i)$ 为 $1$ ；
否则， $f(x,a_i)$ 为 0。

定义 $g(x)=∑i=1nf(x,ai)g(x)=\sum_{i=1}^n f(x,a_i)$

求在 $[0, k)$ 上使得 $g (x)$ 取得最小值的 $x$ ，如有多解输出多个

$1≤m≤k≤1.2×105,1≤n≤1.2×105,0=b1<b2<b3⋯<bm<k,0≤ai≤1.2×1051\leq m\leq k\leq 1.2\times 10^5,1\leq n\leq 1.2\times 10^5,0=b_1<b_2<b_3\cdots<b_m<k,0\leq a_i\leq 1.2\times 10^5$

solution

首先可以想到 $70$ 分的暴力做法

对于每个 $a_i$ 模上 $k$ ，开一个桶去维护

枚举 $x$ ， $O (m)$ 计算答案为 $cnt[i]+∑j=1mcnt[bj+i]cnt[i]+\sum_{j=1}^m cnt[b_j+i]$

复杂度 $O (k m)$

考虑优化，枚举的复杂度基本无法优化，考虑计算答案的优化

观察式子

$\begin{aligned} a_i&\equiv x+b_j(\bmod k)\\ a_i-b_j+k&\equiv x(\bmod k) \end{aligned}$

我们考虑建立两个多项式 $A (x), B (x)$

其中

$A(x)=∑i=0k−1cntaixiB(x)=2+∑i=1k−1cntbk−ixiA(x)=\sum_{i=0}^{k-1}cnta_ix^i \\ B(x)=2+\sum_{i=1}^{k-1}cntb_{k-i}x^i$

那么当我们枚举到 $x$ 的时候，其实就是 $A (x), B (x)$ 的卷积的 $x$ 次项系数和 $x + k$ 次项系数之和

注意忽略掉 $b_1=0$ 的情况，要不然会出问题

利用 $FFT\text{FFT}$ 优化即可

复杂度 $O(klog⁡k)O(k\log k)$

#include <bits/stdc++.h>
using namespace std;

# define Rep(i,a,b) for(int i=a;i<=b;i++)
# define _Rep(i,a,b) for(int i=a;i>=b;i--)
# define RepG(i,u) for(int i=head[u];~i;i=e[i].next)

typedef long long ll;

const int N=1e6+5;
const double pi=acos(-1.0);

template<typename T> void read(T &x){
   x=0;int f=1;
   char c=getchar();
   for(;!isdigit(c);c=getchar())if(c=='-')f=-1;
   for(;isdigit(c);c=getchar())x=(x<<1)+(x<<3)+c-'0';
    x*=f;
}

int k,m,n;
int A[N],B[N],C[N];
int pos[N];
int cnt[N];
int ans,out[N];

struct cp{
	double real,imag;
	cp(double a=0.0,double b=0.0){real=a,imag=b;}	
}a[N],b[N],c[N];

cp operator + (cp a,cp b){return cp(a.real+b.real,a.imag+b.imag);}
cp operator - (cp a,cp b){return cp(a.real-b.real,a.imag-b.imag);}
cp operator * (cp a,cp b){return cp(a.real*b.real-a.imag*b.imag,a.real*b.imag+a.imag*b.real);}

void fft(int limit,cp *a,int opt){
	for(int i=0;i<limit;i++)c[i]=a[pos[i]];
	for(int i=0;i<limit;i++)a[i]=c[i];
	for(int i=1;i<limit;i<<=1){
		int len=i<<1;
		for(int j=0;j<limit;j+=len){
			cp z0(cos(pi/i),opt*sin(pi/i)),z(1.0,0.0);
			for(int k=0;k<i;k++,z=z*z0){
				cp x=a[j+k],y=z*a[j+k+i];
				c[j+k]=x+y;
				c[j+k+i]=x-y;	
			}
		}
		for(int j=0;j<limit;j++)a[j]=c[j];
	}
}

int main()
{
	freopen("tonic.in","r",stdin);
	freopen("tonic.out","w",stdout);
	read(k),read(m),read(n);
	Rep(i,1,m)read(B[i]),b[k-B[i]].real++;
	Rep(i,1,n)read(A[i]),A[i]%=k,a[A[i]].real++;
	b[0].real=2,b[k].real=0;
	int limit=1;
	while(limit<=2*k)limit<<=1;
	for(int i=1;i<limit;i<<=1)for(int j=0;j<i;j++)pos[i+j]=pos[j]+limit/i/2;
	fft(limit,a,1),fft(limit,b,1);
	for(int i=0;i<limit;i++)a[i]=a[i]*b[i];
	fft(limit,a,-1);
	for(int i=0;i<limit;i++)C[i]=(int)(a[i].real/limit+0.5);
	for(int i=0;i<k;i++){
		int tot=C[i]+C[i+k];
		if(tot>ans){
			ans=tot;
			out[0]=0;
			out[++out[0]]=i;
		}
		else if(tot==ans)out[++out[0]]=i;
	}
	printf("%d\n",out[0]);
	Rep(i,1,out[0])printf("%d ",out[i]);
	puts("");
	return 0;
}
/*
12 7 8
0 2 4 5 7 9 11
1 3 5 6 8 10 12 1
*/