快速数论变换

分类

缩写全称作用时间复杂度
DFT离散傅立叶变换时频域转换 O ( n 2 ) O(n^2) O(n2)
FFT快速傅立叶变换时频域转换 ( ( (有精度误差 ) ) ) O ( 大常数 + n l o g 2 n ) O(大常数+nlog_2n) O(大常数+nlog2n)
NTT/FNTT快速数论变换模意义下的时频域转换 O ( 小常数 + n l o g 2 n ) O(小常数+nlog_2n) O(小常数+nlog2n)
MTT任意模数的NTT任意模意义下的时频域转换 O ( n l o g 2 n ) O(nlog_2n) O(nlog2n)
FWT快速沃尔什变换快速集合卷积 O ( 不定 ) O(不定) O(不定)
FMT快速莫比乌斯变换逆莫比乌斯反演? O ( 不定 ) O(不定) O(不定)

快速数论变换(坑)

前置技能

%   本文不包含但必不可少的前置技能:

  1. 快速傅里叶变换
  2. 基础数论(互质,同余运算,欧拉定理,欧拉函数)

%   由于 FFT 涉及到复数运算,难免有精度问题,而且有的时候精度还不小,这便让我们考虑是否有在模意义下快速计算的方法,这就是快速数论变换
( Fast Number-Theoretic Transform,FNT ) (\text{Fast Number-Theoretic Transform,FNT}) (Fast Number-Theoretic Transform,FNT) %   理解FFT后,可以发现,FFT用的是单位根的五大性质。

  1. ω n 0 = ω n n = 1 \omega_n^0=\omega_n^n=1 ωn0=ωnn=1
  2. ∀   i ≠ j , ω n i ≠ ω n j \forall\ i\ne j,\omega_n^i\ne \omega_n^j  i=j,ωni=ωnj (用于还原系数)
  3. ω 2 n 2 k = ω n k \omega^{2k}_{2n}=\omega^k_n ω2n2k=ωnk
  4. ω n k + n 2 = − ω n k \omega_n^{k+\frac{n}{2}}=-\omega_n^k ωnk+2n=ωnk
  5. ∑ i = 0 n − 1 ( ω n j − k ) i = { 0 , k ≠ j n , k = j \begin{aligned}\sum\limits_{i=0}^{n-1}(\omega_n^{j-k})^i =\begin{cases}0,k\ne j\\n,k=j\\\end{cases}\end{aligned} i=0n1(ωnjk)i={0,k=jn,k=j (用于逆FFT)
原根

(数论)阶 对于 ( a , n ) = 1 (a,n)=1 (a,n)=1 的整数,满足 a r ≡ 1 ( m o d n ) a^r≡1 \pmod n ar1(modn)最小整数 r r r,称为 a a a n n n(数论)阶
原根 对于正整数 n n n,整数 a a a,若 a a a n n n 的阶 r r r 等于 φ ( n ) φ(n) φ(n),则称 a a a 为模 n n n 的一个原根。

%   感性地理解一下,一个整数 n n n 的数论阶为其最小开始循环的次方数(且与 n n n 互质)。若一个正整数 a a a p p p 的数论阶恰好为 φ ( n ) \varphi(n) φ(n) ,则 a a a p p p 的原根。

正题

%   根据原根的定义可得,对于质数 p p p,假如 g g g p p p 的原根,则 g 0 ≡ g φ ( p ) ≡ 1 ( m o d p ) g^0≡g^{\varphi(p)}≡1\pmod p g0gφ(p)1(modp),进一步可以得到,对于任意 a ≢ b ( m o d φ ( p ) ) a\not≡ b\pmod {\varphi(p)} ab(modφ(p)),有 g a ≢ g b ( m o d p ) g^a\not ≡g^b\pmod p gagb(modp)

%   对于质数 p = k ⋅ 2 N + 1 p=k\cdot2^N+1 p=k2N+1,设其原根为 g g g。我们令 g n ≡ g p − 1 n ( m o d p ) {g_n≡g^{\frac{p-1}{n}}\pmod p} gngnp1(modp)。注意,这里的 g n g_n gn 是定义出来的,只是在数值上等于 g p − 1 n g^{\frac {p-1}{n}} gnp1,与 g g g 没有关系。容易发现,这样的定义使得 g n g_n gn 满足了第一条、第二条和第三条性质。

性质1证明:
%   由于 g n n ≡ ( g p − 1 n ) n ≡ g p − 1 ≡ 1 ≡ g n 0 ( m o d p ) {g^n_n≡\Big(g^{^{p-1\over n}}\Big)^n≡g^{p-1}≡1≡g_n^0\pmod p} gnn(gnp1)ngp11gn0(modp) %   因此其满足性质1。

性质2证明:
%   由于 p ∈ prime p\in \text{prime} pprime,因而 φ ( p ) = p − 1 \varphi(p)=p-1 φ(p)=p1,因此对于所有 a ≢ b ( m o d ( p − 1 ) ) {a\not≡ b\pmod {(p-1)}} ab(mod(p1)),有 g a ≢ g b ( m o d p ) {g^a\not≡ g^b\pmod p} gagb(modp)
%   因此对于所有 a ≢ b ( m o d ( p − 1 ) ) {a\not≡ b\pmod {(p-1)}} ab(mod(p1)),有 g n a ≡ ( g p − 1 n ) a ≢ ( g p − 1 n ) b ≡ g n b ( m o d p ) {g_n^a≡(g^{\frac{p-1}{n}})^a\not≡ (g^{\frac{p-1}{n}})^b≡g_n^b\pmod p} gna(gnp1)a(gnp1)bgnb(modp)  因此其满足性质2。

性质3证明:
%   由于 g 2 n 2 k ≡ ( g p − 1 2 n ) 2 k ≡ ( g p − 1 n ) k ≡ g n k ( m o d p ) {g^{2k}_{2n}≡(g^{p-1\over 2n})^{2k}≡(g^{p-1\over n})^{k}}≡g_n^k\pmod p g2n2k(g2np1)2k(gnp1)kgnk(modp) %   因此其满足性质3。

性质4证明:
%   由于 p p p 是质数,并且 g n n ≡ 1 ( m o d p ) {g_n^n \equiv 1 \pmod p} gnn1(modp),因而有
g n n 2 ≡ 1 ( m o d p ) or g n n 2 ≡ − 1 ( m o d p ) {g^{\frac n2}_n≡1\pmod p}\quad \text{or}\quad {g^{\frac n2}_n≡-1\pmod p} gn2n1(modp)orgn2n1(modp) %   由第二条性质可得 g n n 2 ≢ g n n ≡ 1 ( m o d p ) {g^{\frac n2}_n\not≡ g_n^{n} ≡1\pmod p} gn2ngnn1(modp),因而舍去前者,取后者,因而有 g n k + n 2 ≡ g n k ⋅ g n n 2 = − g n k ( m o d p ) g_n^{k+\frac n2}≡g_n^k\cdot g_n^{\frac n2}=-g_n^k\pmod p gnk+2ngnkgn2n=gnk(modp) %   因此其满足性质4。

性质5证明:
%   当 k ≠ j k\ne j k=j 时,根据等比数列的求和公式,可得:
∑ i = 0 n − 1 ( g n j − k ) i = ( g n k ) n − 1 g n k − 1 = ( g n n ) k − 1 g n k − 1 = 1 − 1 g n k − 1 =   0 \begin{aligned}\sum_{i=0}^{n-1}(g_n^{j-k})^i =&\dfrac{(g_n^k)^{n}-1}{g_n^k-1}\\ =&\dfrac{(g_n^n)^{k}-1}{g_n^k-1}\\ =&\dfrac{1-1}{g_n^k-1}\\ =&\ 0\\ \end{aligned} i=0n1(gnjk)i====gnk1(gnk)n1gnk1(gnn)k1gnk111 0

%   当 k = j k=j k=j 时,可得 ∑ i = 0 n − 1 ( g n j − k ) i = ∑ i = 0 n − 1   1 = n \begin{aligned}\sum\limits_{i=0}^{n-1}(g_n^{j-k})^i =\sum\limits_{i=0}^{n-1}\ 1=n\end{aligned} i=0n1(gnjk)i=i=0n1 1=n
%   综上,有 ∑ i = 0 n − 1 ( g n j − k ) i = { 0 , k ≠ j n , k = j \begin{aligned}\sum\limits_{i=0}^{n-1}(g_n^{j-k})^i =\begin{cases}0,k\ne j\\n,k=j\\\end{cases}\end{aligned} i=0n1(gnjk)i={0,k=jn,k=j

%   至此,我们发现 g n g_n gn 满足 w n w_n wn 的五条性质,因而我们可以大方地用 g n g_n gn 代替 w n w_n wn

细节

%   为什么模数必须取 k ⋅ 2 N + 1 k\cdot2^N+1 k2N+1?注意到我们取 g n ≡ g p − 1 n ( m o d p ) {g_n≡g^{\frac{p-1}{n}}\pmod p} gngnp1(modp),因而 p − 1 n \dfrac {p-1}{n} np1 必须是整数。考虑到之前FFT的过程是不断地二分,因而我们可以保证 n n n 总是 2 2 2 的倍数,因而模数必须有足够多的质因子 2 2 2,因而取 k ⋅ 2 N + 1 k\cdot2^N+1 k2N+1,且 N N N 必须足够大。

%   但这也注定了FNT的一个弊端,对于如果 r ⋅ 2 N + 1 r⋅2^N+1 r2N+1 是个素数,那么在模 k ⋅ 2 N + 1 k⋅2^N+1 k2N+1 意义下,只可以处理 2 N 2^N 2N 以内规模的数据。对于一些好的质数及其原根见这里

%   在逆 F N T FNT FNT 的时候,我们需要求 g n − 1 g_n^{-1} gn1,在模意义下,也就是 g n g_n gn 的逆元,预处理即可。

代码

Accepted -O2 \text{Accepted -O2} Accepted -O2 / 用时: 1925 m s 1925ms 1925ms / 内存: 243752 KB 243752\text{KB} 243752KB

#include<bits/stdc++.h>
using namespace std;
#define LL long long
const int maxn=3*1e6+10,P=998244353,G=3,inv_G=332748118;
int N,M,limit,L,r[maxn];
LL a[maxn],b[maxn];
inline LL fastpow(LL a,LL k) {
    LL base=1;
    while(k) {
        if(k&1) base=(base*a)%P;
        a=(a*a)%P;
        k>>=1;
    } return base%P;
}
#define inv(x) fastpow(x,P-2)
inline void NTT(LL *A,int type) {
    for(int i=0;i<limit;i++)
        if(i<r[i]) swap(A[i],A[r[i]]);
    for(int mid=1;mid<limit; mid <<= 1) {
        LL Wn=fastpow(type==1?G:inv_G,(P-1)/(mid<<1));
        for(int j=0;j<limit;j+=(mid<<1)) {
            LL w=1;
            for(int k=0;k<mid;k++,w=(w*Wn)%P) {
                int x=A[j+k],y=w*A[j+k+mid]%P;
                A[j+k]=(x+y)%P,
                A[j+k+mid]=(x-y+P)%P;
            }
        }
    }
}
int main() {
    scanf("%d%d",&N,&M);
    for(int i=0;i<=N;i++)
        scanf("%d",&a[i]),a[i]=(a[i]+P)%P;
    for(int i=0;i<=M;i++)
        scanf("%d",&b[i]),b[i]=(b[i]+P)%P;
    for(limit=1,L=0;limit<=N+M;limit<<=1,L++);
    for(int i=0;i<limit;i++)
        r[i]=(r[i>>1]>>1)|((i&1)<<(L-1));
    NTT(a,1);NTT(b,1);
    for(int i=0;i<limit;i++)
        a[i]=(a[i]*b[i])%P;
    NTT(a,-1);
    LL inv=inv(limit);
    for(int i=0;i<=N+M;i++)
        printf("%d ",(a[i]*inv)%P);
    return 0;
}

可以发现,FNT和FFT的效率差别还是很大的,FNT完美碾压FFT。

任意模数

%   FNT虽快,但使用局限太多,当模数不为 k ⋅ 2 N + 1 k⋅2^N+1 k2N+1 时,FNT便无法发挥其效果,有没有解决办法呢?当然有。(坑)
%https://www.hrwhisper.me/introduction-to-simplex-algorithm/

  • 4
    点赞
  • 12
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值