【HNOI2017/BZOJ4827】礼物 FFT

最新推荐文章于 2020-09-24 22:18:03 发布

CRTorlonia

最新推荐文章于 2020-09-24 22:18:03 发布

阅读量240

点赞数

本文链接：https://blog.csdn.net/crtorlonia/article/details/79947587

版权

原题走这里

原题相当于求 $\sum(a[i]-b[i]+C)^2$ ，其中b可以任意旋转，C为整数，由于原题中有a,b小于m的限制条件，所以C范围最大也大不过 $(-m,m)$ 。
那么我们可以把式子拆开：

\sum i = 0 n - 1 (a [i] 2 + b [i] 2 - 2 * a [i] * b [i] + 2 * C (a [i] - b [i]) + C 2)

$\sum_{i=0}^{n-1}(a[i]^2+b[i]^2-2*a[i]*b[i]+2*C(a[i]-b[i])+C^2)$
由于除了

a[i]∗b[i] a [ i ] ∗ b [ i ] $a[i]*b[i]$ 这一项以外，其他所有项都可以预处理或枚举，所以问题就转化为了：如何最大化：

∑n−1i=0a[i]∗b[i] ∑ i = 0 n − 1 a [ i ] ∗ b [ i ] $\sum_{i=0}^{n-1}a[i]*b[i]$
把b反转，就变成了

∑n−1i=0a[i]∗b[n−1−i] ∑ i = 0 n − 1 a [ i ] ∗ b [ n − 1 − i ] $\sum_{i=0}^{n-1}a[i]*b[n-1-i]$ ，就可以FFT了。

然而旋转我们怎么处理呢？我们可以把上式拆开成：
$\sum_{i=0}^{n-1}a[i]*b[n-j-1-i]+\sum_{i=n-j}^{n-1}a[i]+b[2*n-j-1-i]$
FFT求出 $D[i]=A[i]*B[i]$ ，上式就等于 $D[n-j-1]+D[2*n-j-1]$
枚举C求出最小值就可以了。

具体实现见代码如下：

#include <bits/stdc++.h>
using namespace std;
int n,m,N,sum,mx,mn=2e9,summ;
complex<double> a[200000],b[200000],aa[200000],w[200000];
void init() {
    const double pi=acos(-1);
    for(int i=0; i<N; i++) {
        w[i]=complex<double>(cos(2.0*pi*i/N),sin(2.0*pi*i/N));
    }
}
void FFT(complex<double> *a) {
    for(int i=0,j=0; i<N; i++) {
        if(i>j)swap(a[i],a[j]);
        for(int k=N>>1; (j^=k)<k; k>>=1);
    }
    for(int i=2; i<=N; i<<=1) {
        int m=i>>1,l=N/i;
        for(int j=0; j<N; j+=i) {
            for(int k=0; k!=m; k++) {
                complex<double> t=a[j+m+k]*w[l*k];
                a[j+m+k]=a[j+k]-t;
                a[j+k]+=t;
            }
        }
    }
}
int main() {
    cin>>n>>m;
    for(int i=0; i<n; i++) {
        cin>>a[i];
    }
    for(int i=0; i<n; i++) {
        cin>>b[n-i-1];
    }
    for(int i=0; i<n; i++) {
        sum+=round((a[i]*a[i]+b[n-i-1]*b[n-i-1]).real());
        summ+=round((b[i]-a[i]).real());
    }
    N=n<<1;
    for(int i=1; i<N; N=max(i<<=1,N));
    init();
    FFT(a);
    FFT(b);
    reverse(w+1,w+N);
    for(int i=0; i<N; i++) {
        a[i]*=b[i];
    }
    FFT(a);
    for(int i=0; i<N; i++) {
        a[i]/=N;
    }
    mx=a[n-1].real();
    for(int i=0; i<n-1; i++) {
        mx=max(mx,int(round((a[i]+a[i+n]).real())));
    }
    for(int i=-m; i<=m; i++) {
        mn=min(mn,n*i*i+2*i*summ);
    }
    cout<<sum-2*mx+mn<<endl;
    return 0;
}