多项式入门学习笔记

前言

这个傻狗啥都不会还敢来学多项式。

虽然这好像连入门都不算（（（

背板！！！背板！！！

前置的一些定义

点值表示法：将 \(n+1\) 个不同的值带入到 \(n\) 次多项式中形成 \(n+1\) 个点的坐标以代表唯一多项式，支持 \(O(n)\) 计算卷积。
单位根：复数意义下 \(x^{n}=1\) 即为 \(n\) 次单位根。(有 \(n\) 个)
本原单位根：\(e^{i\frac{2\pi}{n}}\) 称为 \(n\) 次单位根的本原单位根，记为 \(ω_n\) 。
原根：若 \(1\le g\le n-1\) 且 \(g\) 模 \(n\) 的阶为 \(\varphi(n)\) ，则称 \(g\) 为 \(n\) 的原根。

DFT

具体实现

暴力计算点值

\[b_k=\sum\limits_{i=0}^{n-1} a_{i}k \]

FFT

实现原理

在通过点值表示法来求卷积的过程中，利用转换和单位根的性质可以将求解的复杂度降为 \(O(n\log n)\) 。

具体实现

懒得码那么多又臭又长的证明了，直接给结论：(以下设 \(m=\frac{n}{2}\))

令

\[A(x)=\sum\limits_{i=0}^{n-1} a_ix^i \]

则

\[A(\omega_n^{k})=A_0(\omega_m^{k})+\omega_n^{k}A_1(\omega_m^{k})\ (k\le m) \]

\[A(\omega_n^{k+m})=A_0(\omega_m^{k})-\omega_n^{k}A_1(\omega_m^{k})\ (k\le m) \]

其中

\[A_0=\sum\limits_{i=0}^{m-1} a_{2i}x^i \]

\[A_1=\sum\limits_{i=0}^{m-1} a_{2i+1}x^i \]

就可以递归实现了。

优化(迭代实现)

建出递归树后我们可以发现，最后一层的树是最初的数的二进制下反转后排序的值。

看着很拗口，上例子：

0 1 2 3 4 5 6 7
0 2 4 6, 1 3 5 7
0 4, 2 6, 1 5, 3 7
0, 4, 2, 6, 1, 5, 3, 7
二进制：  
000, 100, 010, 110, 001, 101, 011, 111
翻转后：
000, 001, 010, 011, 100, 101, 110, 111

求出最后一行就可以向上合并了！

这种变换好像学名叫做 \(\text{Rader}\) 排序。

\(O(n)\) 预处理出 \(\text{Rader}\) 排序的映射关系，然后将数组交换就可以 \(O(n)\) 求出底层数了！！！

IDFT

具体实现

设

\[b_k=\sum\limits_{i=0}^{n-1} a_{i}\omega_n^{ik} \]

即 \(k\) 的点值。

则

\[a_k=\frac{1}{n}\sum\limits_{i=0}^{n-1} b_{i}\omega_n^{-ik} \]

反正就是 \(\text{FFT}\) 的逆变换。

不会证明，摆烂。

IFFT

快速地实现 \(\text{IDFT}\) 。

可以发现 \(\text{IDFT}\) 和 \(\text{DFT}\) 唯一相差的就是 \(\omega_m^{-ik}\) ，因此可以和 \(\text{FFT}\) 套用同种方法实现。

最后结果记得除以 \(n\) 。

NNT

这个 \(\text{FFT}\) 又要码又臭又长的复数，又要卡精度，我们考虑用别的东西代替单位根实现。

原根！好东西！！！

再次不加证明地扔出一个结论：

\[\omega_n\equiv g^{\frac{p-1}{n}}\pmod{p} \]

证明，狗都不会！！！

然后代换掉 \(\text{FFT}\) 中的复数单位根就好了。

Code

\(\text{FFT}\) (递归版)

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

inline int read(){
    int x=0,w=1;
    char ch=getchar();
    while(ch>'9'||ch<'0'){if(ch=='-')w=-1;ch=getchar();}
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
    return x*w;
}

const double phi=acos(-1);
int n,m,lim;
complex<double> f[N],g[N];

void F(complex<double>*a,int n,int opt){
    if(!n) return ;
    complex<double> a0[n],a1[n];
    for(int i=0;i<n;++i)
        a0[i]=a[i<<1],a1[i]=a[i<<1|1];
    F(a0,n>>1,opt);F(a1,n>>1,opt);
    complex<double> t(cos(phi/n),sin(phi/n)*opt),w(1,0);
    for(int i=0;i<n;++i,w*=t)
        a[i]=a0[i]+w*a1[i],a[i+n]=a0[i]-w*a1[i];
}

signed main(){
    n=read();m=read();
    for(int i=0;i<=n;++i) f[i]=read();
    for(int i=0;i<=m;++i) g[i]=read();
    for(lim=1;lim<=n+m;lim<<=1);
    F(f,lim>>1,1);F(g,lim>>1,1);
    for(int i=0;i<lim;++i) f[i]*=g[i];
    F(f,lim>>1,-1);
    for(int i=0;i<=n+m;++i)printf("%.0f ",fabs(f[i].real())/lim);
    return 0;
}

\(\text{FFT}\) (迭代版)

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

inline int read(){
    int x=0,w=1;
    char ch=getchar();
    while(ch>'9'||ch<'0'){if(ch=='-')w=-1;ch=getchar();}
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
    return x*w;
}

const double phi=acos(-1);
int n,m,lim,l;
int r[N];
complex<double> f[N],g[N];

void F(complex<double>*a,int n,int opt){
    for(int i=0;i<lim;++i)
        if(i<r[i]) swap(a[i],a[r[i]]);
    for(int i=1;i<lim;i<<=1){
        complex<double> W(cos(phi/i),sin(phi/i)*opt);
        for(int j=0;j<lim;j+=(i<<1)){
            complex<double> w(1,0);
            for(int k=0;k<i;++k,w*=W){
                complex<double> x=a[j+k],y=w*a[j+i+k];
                a[j+k]=x+y; a[j+k+i]=x-y;
            }
        }
    }
}

signed main(){
    n=read();m=read();
    for(int i=0;i<=n;++i) f[i]=read();
    for(int i=0;i<=m;++i) g[i]=read();
    for(lim=1;lim<=n+m;lim<<=1,++l);
    for(int i=0;i<lim;++i) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    F(f,lim>>1,1);F(g,lim>>1,1);
    for(int i=0;i<lim;++i) f[i]*=g[i];
    F(f,lim>>1,-1);
    for(int i=0;i<=n+m;++i)printf("%.0f ",fabs(f[i].real())/lim);
    return 0;
}

\(\text{NNT}\) (只有迭代版)

#include<bits/stdc++.h>
#define N 4000005
#define int long long
using namespace std;

inline int read(){
    int x=0,w=1;
    char ch=getchar();
    while(ch>'9'||ch<'0'){if(ch=='-')w=-1;ch=getchar();}
    while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
    return x*w;
}

const int mod = 998244353;//good moshu !!!
const int g1 = 3,g2= 332748118;
int n,m,lim,l;
int r[N];
int f[N],g[N];

int qpow(int x,int k){
    int res=1;
    while(k){
        if(k&1) res=res*x%mod;
        x=x*x%mod;k>>=1;
    }
    return res;
}

void NTT(int *a,int opt){
    for(int i=0;i<lim;++i)
        if(i>r[i]) swap(a[i],a[r[i]]);
    for(int i=1;i<lim;i<<=1){
        int W=qpow(opt==1?g1:g2,(mod-1)/(i<<1));
        for(int j=0;j<lim;j+=(i<<1)){
            int w=1;
            for(int k=0;k<i;++k,w=w*W%mod){
                int x=a[j+k],y=w*a[j+i+k]%mod;
                a[j+k]=(x+y)%mod; a[j+k+i]=(x-y+mod)%mod;
            }
        }
    }
}

signed main(){
    n=read();m=read();
    for(int i=0;i<=n;++i) f[i]=read();
    for(int i=0;i<=m;++i) g[i]=read();
    for(lim=1;lim<=n+m;lim<<=1,++l);
    for(int i=0;i<lim;++i) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    NTT(f,1);NTT(g,1);
    for(int i=0;i<lim;++i) f[i]*=g[i];
    NTT(f,-1);
    int inv =qpow(lim,mod-2);
    for(int i=0;i<=n+m;++i)printf("%lld ",f[i]*inv%mod);
    return 0;
}

posted @ 2022-05-10 20:24 lxg_swrz 阅读(95) 评论(0) 收藏举报

刷新页面返回顶部

ChengJY_的狗窝

多项式入门学习笔记

前言

前置的一些定义

DFT

FFT

IDFT

IFFT

NNT

Code

公告