[FFT]luogu 3803 模板多项式乘法

Posted 2020-12-27 mastervan

题目背景

这是一道FFT模板题

注意：虽然本题开到3s，但是建议程序在1s内可以跑完，本题需要一定程度的常数优化。

题目描述

给定一个n次多项式F(x)，和一个m次多项式G(x)。

请求出F(x)和G(x)的卷积。

输入输出格式

输入格式：

 

第一行2个正整数n,m。

接下来一行n+1个数字，从低到高表示F(x)的系数。

接下来一行m+1个数字，从低到高表示G(x))的系数。

 

输出格式：

 

一行n+m+1个数字，从低到高表示F(x)?G(x)的系数。

 

输入输出样例

输入样例

1 2
1 2
1 2 1

输出样例

1 4 5 2

说明

保证输入中的系数大于等于 0 且小于等于9。

对于100%的数据： n, m leq {10}^6n,m≤106 , 共计20个数据点，2s。

数据有一定梯度。

空间限制：256MB

 

分析

FFT的算法详解作为一个初中OIer写不出来，但是大体还是能讲讲的。

首先多项式有两种表达方式：系数法和点值法。

系数法就是我们常用的方法，如：f(x)=a+bx+cx²

点值法只能有点抽象的说，如根据两点确定一条直线（一次函数），三点确定一条抛物线（二次函数），……，n+1个点确定一个n次函数，如下：

f(x)=(x[0],f(x[0])),(x[1],f(x[1])),...,(x[n],f(x[n]))

离散傅里叶变换（DFT）就是把一个用系数法表达的多项式转化为点值法的形式

离散傅里叶反变换（IDFT）则是把点值法转化为系数法表达

（已经说的尽可能不抽象了）

那么回到我们要解决的问题：多项式乘法，显然朴素算法是O（n²）的（系数法）

显然我们可以先把所需求积的两个多项式变换为点值法，然后有F(x)=f(x)*g(x)

可是F(x)这个n+1元方程组不太好解，高斯消元？好像甚至不如朴素算法。

因为我们计算x[0],x[0]²,...,x[0]ⁿ时，需要大量时间计算，那么有没有什么数简单一些呢？有：1和-1

可是两个数哪里够？这时候，我们引入虚数（虚数单位为i）。

显然，i²=-1,那么i²*i²=(-1)*(-1)=1=i⁴

这样看起来有点多了，可还是不够。

我们可以在复平面直角坐标系的原点做一个半径为1的圆，把这个圆n等分，相当于产生了n个弧，弧的端点则是我们可以利用的数为ω_n^k，这个数我们称为n的k次复根，意思显然，这个数自乘k次就会归1

（憋不出来啦，原谅博主数学能力有限，如对屏幕前的dalao您有所帮助本人非常荣幸！下面附上博主认为有帮助的一篇博文ovo）

感谢GGN大佬的博客指导

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <complex>
using namespace std;
const int N=1e7+1;
const double Pi=3.14159265358979;
typedef complex<double> cn;
cn a[N],b[N];
int r[N],mxb,bit;

void Get_Rev(int bit) {
    for (int i=0;i<(1<<bit);i++)
    r[i]=(r[i>>1]>>1)|((i&1)<<(bit-1));
}

void FFT(cn *a,int mx,int inot) {
    for (int i=0;i<mx;i++)
    if (i<r[i]) swap(a[i],a[r[i]]);
    for (int mlen=1;mlen<mx;mlen<<=1) {
        cn we=exp(cn(0,inot*Pi/mlen));
        for (int len=mlen<<1,l=0;l<mx;l+=len) {
            cn wk=cn(1,0);
            for (int k=l;k<l+mlen;k++) {
                cn x=a[k],y=wk*a[k+mlen];
                a[k]=x+y;a[k+mlen]=x-y;
                wk*=we;
            }
        }
    }
    if (inot==-1)
    for (int i=0;i<mx;i++) a[i]/=mx;
}

int main() {
    int n,m;
    scanf("%d%d",&n,&m);
    for (int i=0;i<=n;i++) {
        int x;
        scanf("%d",&x);
        a[i]=(double)x;
    }
    for (int i=0;i<=m;i++) {
        int x;
        scanf("%d",&x);
        b[i]=(double)x;
    }
    mxb=1;
    while (mxb<=n+m) mxb<<=1,bit++;
    Get_Rev(bit);
    FFT(a,mxb,1);FFT(b,mxb,1);
    for (int i=0;i<=mxb;i++) a[i]*=b[i];
    FFT(a,mxb,-1);
    for (int i=0;i<=n+m;i++) printf("%d ",(int)(a[i].real()+0.5));
}

View Code