任意模数NTT（学习笔记）

最新推荐文章于 2023-07-04 20:41:35 发布

QAQQQQQQQQQQQ

最新推荐文章于 2023-07-04 20:41:35 发布

阅读量1.9k

点赞数 3

CC 4.0 BY-SA版权

分类专栏： ——多项式—— 笔记 NTT 文章标签： NTT

本文链接：https://blog.youkuaiyun.com/sizeof_you/article/details/84557728

笔记同时被 3 个专栏收录

11 篇文章

订阅专栏

——多项式——

6 篇文章

订阅专栏

NTT

5 篇文章

订阅专栏

本文深入探讨了快速傅立叶变换(FFT)和数论变换(NTT)在处理大整数乘法和多项式乘法中的作用。详细介绍了NTT的几种实现方式，包括特殊模数法、三模数法、拆系数法和MTT法，以及如何通过中国剩余定理进行结果合并，以避免精度损失。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

$F F T$ 有时候会被卡精度？所以可能会有模数，有了模数以后就需要模数的原根。

原根是什么？（留坑待填）

$N T T$ 有很多种解决方法
$1 .$ 特殊模数
$(2 k + 1) ∣ (p - 1) ， (p - 1) > D F T 的长度$ ，可以直接暴力求原根 $g$ ，用 $g$ 代替单位复数根

$2 .$ 一般模数

三模数法（9次DFT）
如果这个模数的原根不好求，而模数又很大，可以用三个模数
要求 $mod1\times mod2\times mod3\ge n\times p^2$
常用 $998244353, 1004535809, 469762049$ ，因为他们的原根都是 $3$
对这几个模数分别做 $D F T$ ，然后用中国剩余定理合并，然后对原模数取模即可
中国剩余定理部分：
$ans\equiv c1\ mod\ m1$
$ans\equiv c2\ mod\ m2$
$ans\equiv c3\ mod\ m3$
如果直接合并的话会爆 $long\ long$ ，可以先合并两个，再用奇技淫巧合并第三个
合并前两个：
$ans\equiv (c1\times m2\times inv(m2,m1)+c2\times m1\times inv(m1,m2)) mod\ (m1\times m2)$
其中 $i n v (x, y)$ 表示 $x$ 对 $y$ 取模的逆元。
把上式化简为 $ans\equiv C\ mod\ M$
设 $ans=x\times M+C=y\times m3+c3$
接下来很重要：求出 $x\ mod\ m3$ 意义下的值：
$x\equiv (c3-C)\times M^{-1}\ mod\ m3$
这样就可以算出右半部分的值 $q$ ，令 $x=k\times m3+q$ ，代入 $a n s$ 得：
$ans=k\times m1\times m2\times m3+q\times M+C$
因为 $ans\in [0,m1\times m2\times m3)$ ，所以 $k = 0$ ，于是 $a n s$ 就可直接计算。
做9次DFT，常数极大
有一道模板题luogu4245
直接用三模数法，代码如下：

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#define LL long long
#define maxn 400005
using namespace std;
const LL mod1=998244353,mod2=1004535809,mod3=469762049,g=3;
const LL M=1LL*mod1*mod2;

inline int rd(){
    int x=0,f=1;char c=' ';
    while(c<'0' || c>'9') f=c=='-'?-1:1,c=getchar();
    while(c<='9' && c>='0') x=x*10+c-'0',c=getchar();
    return x*f; 
}

int n,m,p,rev[maxn],limit=1,l;
LL a[3][maxn],b[3][maxn],ans[maxn];

inline LL mul(LL x,int k,LL MOD){
    LL ret=0;
    while(k){
        if(k&1) (ret+=x)%=MOD;
        (x+=x)%=MOD; k>>=1;
    } return ret%MOD;
}

inline LL qpow(LL x,int k,int MOD){
    LL ret=1;
    while(k){
        if(k&1) ret=ret*x%MOD;
        x=x*x%MOD; k>>=1;
    } return ret%MOD;
}

inline void NTT(LL *F,int type,int MOD){
    for(int i=0;i<limit;i++){
        F[i]%=MOD;
        if(i<rev[i]) swap(F[i],F[rev[i]]);
    }
    for(int mid=1;mid<limit;mid<<=1){
        LL Wn=qpow(g,type==1?(MOD-1)/(mid<<1):(MOD-1-(MOD-1)/(mid<<1)),MOD);//!
        for(int r=mid<<1,j=0;j<limit;j+=r){
            LL w=1;
            for(int k=0;k<mid;k++,w=w*Wn%MOD){
                LL x=F[j+k],y=w*F[j+mid+k]%MOD;
                F[j+k]=(x+y)%MOD; F[j+mid+k]=(x-y+MOD)%MOD;
            }
        }
    }
    if(type==-1){
        LL INV=qpow(limit,MOD-2,MOD);//除以limit
        for(int i=0;i<limit;i++) F[i]=F[i]*INV%MOD; 
    }
}

inline void CRT(){
    for(int i=0;i<limit;i++){
        LL tmp=0;
        (tmp+=mul(a[0][i]*mod2%M,qpow(mod2,mod1-2,mod1),M))%=M;
        (tmp+=mul(a[1][i]*mod1%M,qpow(mod1,mod2-2,mod2),M))%=M;
        a[1][i]=tmp;
    }
    for(int i=0;i<limit;i++){//x=(c3-C)/M mod mod3
        LL tmp=(a[2][i]-a[1][i]%mod3+mod3)%mod3*qpow(M%mod3,mod3-2,mod3)%mod3;
        ans[i]=(M%p*tmp%p+a[1][i]%p)%p;
    }
}

inline void solve(int x,int MOD){
    NTT(a[x],1,MOD); NTT(b[x],1,MOD);
    for(int i=0;i<limit;i++) a[x][i]=a[x][i]*b[x][i]%MOD;
    NTT(a[x],-1,MOD);
}

int main(){
    n=rd(); m=rd(); p=rd();
    for(int i=0;i<=n;i++){
        int x=rd();
        for(int j=0;j<3;j++) a[j][i]=x%p;
    }
    for(int i=0;i<=m;i++){
        int x=rd();
        for(int j=0;j<3;j++) b[j][i]=x%p;
    }
    while(limit<=n+m) limit<<=1,++l;
    for(int i=0;i<limit;i++)
        rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    solve(0,mod1); solve(1,mod2); solve(2,mod3);
    CRT();
    for(int i=0;i<=n+m;i++) printf("%lld ",ans[i]);
    return 0;
}