bzoj2194 快速傅立叶之二

最新推荐文章于 2020-09-21 22:30:57 发布

Euryale_

最新推荐文章于 2020-09-21 22:30:57 发布

阅读量451

点赞数

分类专栏： FFT\NTT 文章标签：快速fourier变换

本文链接：https://blog.csdn.net/u010336344/article/details/70160184

版权

FFT\NTT 专栏收录该内容

5 篇文章 0 订阅

订阅专栏

题目链接：bzoj2194
题目大意：
给出数组a，b，计算 $c[k]=\sum_{k≤i<n}a[i]\times b[i-k]$

题解：
FFT
按通常的卷积形式应该是 $c[k]=\sum_{i=1}^{k}a[i]\times b[k-i]$
要求和相等的
可是题目要的是 $c[k]=\sum_{i=k}^{n}a[i]\times b[i-k]$ 【把n-1一下】
是差相等的
于是可以考虑把其中一个数组翻转一下。即令b[i]=b[n-i+1]这样子的翻转。
那么式子就变成了

c [k] = \sum i = k n a [i] \times b [n - (i - k)] = \sum i = k n a [i] \times b [n + k - i]

$c[k]=\sum_{i=k}^na[i]\times b[n-(i-k)]=\sum_{i=k}^na[i]\times b[n+k-i]$
发现

i+(n+k−i)=n+k $i+(n+k-i)=n+k$ 符合通常的形式了。
令

d[n+k]=∑ni=ka[i]×b[n+k−i] $d[n+k]=\sum_{i=k}^na[i]\times b[n+k-i]$ ，这样就能直接上FFT了。
所以最后

c[0…n] $c[0\dots n]$ 对应的就是

d[n+1…2n] $d[n+1\dots 2n]$

总结：
形如 $c[k]=\sum a[i]\times b[i-k]$ 这样子的，可以把其中一个翻转过来做

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<cmath>
#include<iostream>
#include<algorithm>
using namespace std;
#define maxn 100010

const double pi=acos(-1);
struct Complex
{
    double x,y;
    Complex(){x=y=0;}
    Complex(double x,double y):x(x),y(y){}
    friend inline Complex operator + (Complex x,Complex y){return Complex(x.x+y.x,x.y+y.y);}
    friend inline Complex operator - (Complex x,Complex y){return Complex(x.x-y.x,x.y-y.y);}
    friend inline Complex operator * (Complex x,Complex y){return Complex(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);}
}a[maxn*4],b[maxn*4];
void DFT(Complex *c,int ln,int t)
{
    if (ln==1) return;
    Complex A0[ln>>1],A1[ln>>1];
    for (int i=0;i<=ln;i+=2) A0[i>>1]=c[i],A1[i>>1]=c[i+1];
    DFT(A0,ln>>1,t);DFT(A1,ln>>1,t);
    Complex wn(cos(2*pi/ln),t*sin(2*pi/ln)),w(1,0);
    for (int i=0;i<(ln>>1);i++,w=w*wn)
    {
        c[i]=A0[i]+w*A1[i];
        c[i+(ln>>1)]=A0[i]-w*A1[i];
    }
}
int ans[maxn];
void FFT(int ln)
{
    int fn=1;while (fn<=ln*2) fn<<=1;
    DFT(a,fn,1);DFT(b,fn,1);
    for (int i=0;i<=fn;i++) a[i]=a[i]*b[i];
    DFT(a,fn,-1);
    for (int i=0;i<=ln;i++) ans[i]=(int)(a[i+ln].x/fn+0.5);
}
int main()
{
    //freopen("a.in","r",stdin);
    //freopen("a.out","w",stdout);
    int i,n;
    scanf("%d",&n);n--;
    for (i=0;i<=n;i++)
     scanf("%lf%lf",&a[i].x,&b[i].x);
    for (i=0;i<=n/2;i++) swap(b[i],b[n-i]);
    FFT(n);
    for (i=0;i<=n;i++) printf("%d\n",ans[i]);
    return 0;
}