• 快速数论变换NTT


    前置知识

    FFT

    我们知道,\(FFT\)通过代入单位根这样特殊的数,实现了系数表示法与点值表示法之间的快速转化,从而快速进行多项式乘法。

    但是,单位根性质虽然优秀,但它们毕竟是复数,对精度要求过高,因此,我们需要一种特殊的整数,让它们替代单位根:


    原根

    \(a,p\)互素,且\(p>1\)

    对于\(a^n \equiv 1 \pmod{p}\)最小的n,我们称之为\(a\)\(p\)的阶,记做\(\delta_p(a)\)

    原根

    \(p\)是正整数,\(a\)是整数,若\(\delta_p(a)\)等于\(\phi(p)\),则称\(a\)为模\(p\)的一个原根。

    一个数的原根是不唯一的,

    如果\(p\)有原根,那么它一定有\(\phi(\phi(p))\)个原根。

    也不是所有的数都有原根,原根存在的重要条件为\(m = 2,4,p^a,2p^{a}\),其中\(p\)为奇素数,\(a \ge 1\)

    原根有一个重要的性质:

    • \(P\)为素数,假设一个数\(g\)\(P\)的原根,那么\(g^i \pmod P (1<g<P,0<i<P)\)的结果两两不同

    快速数论变换NTT

    我们先给出一个结论:我们可以将在\(FFT\)中使用的\(w_n\)替换为\(g^{\frac {p-1}{n}}\)

    为什么可以呢?

    考虑我们之前用到了\(w_n\)的那些性质:

    1.\(w_n^0=w_n^n=1\)

    证明:显然\(g^0=1\)\(w_n^n=g^{p-1} \equiv 1 \pmod p\),后者是因为费马小定理。

    2.\(w_{2n}^{2k}=w_{n}^{k}\)

    证明:

    \[(g^{\frac {p-1}{2n}})^{2k}=g^{\frac {2k(p-1)}{2n}}=g^{\frac {k(p-1)}{n}}=(g^{\frac {p-1}{n}})^{k}​ \]

    3.\(w_{n}^{k+\frac n2}=-w_n^k\)

    证明:

    \[(w_n^{\frac n2})^2=g^{p-1} \equiv 1\pmod p \]

    \[w_n^{\frac n2}\equiv 1或-1\pmod p \]

    根据上面提到的性质\(w_n^{\frac n2}\not=w_n^n\),所以\(w_n^{\frac n2}\equiv -1\pmod p\)

    所以这条性质自然成立。

    综上,又因为\(w_n^k\)互不相同,所以我们可以用\(g^{\frac {p-1}{n}}\)作为\(w_n\)

    这就是快速数论变换\(NTT\)的原理了,它比\(FFT\)更优,因为它通过取模避免了精度上的问题,但它也有自身的问题,即对模数有限制——必须是有原根的数。

    常用的模数有\(998244353,1004535809,469762049\),它们的原根都是\(3\)

    对于其他的模数,我们有方法————任意模数\(NTT\),这个我们之后再说。

    至于实现,就是将\(FFT\)中算\(w_n\)的时候改一下即可,这里给大家提供一个有封装的代码,可以通过模板题

    	#include<bits/stdc++.h>
    using namespace std;
    const int N=(1<<21)+10;
    const int mod=998244353;
    const int gg=3;
    const int giv=(mod+1)/3;
    typedef long long ll;
    
    namespace Math{
    	int inv[N],base[N],jc[N],tp=2;
    	inline void init(int n){
    		if(tp==2) inv[0]=inv[1]=base[0]=base[1]=jc[0]=jc[1]=1;
    		for(tp;tp<=n;++tp)
    			inv[tp]=1ll*(mod-mod/tp)*inv[mod%tp]%mod; 
    	}
    	inline int ksm(int a,ll b){
    		int ret=1;
    		for(;b;a=1ll*a*a%mod,b>>=1) (b&1)&&(ret=1ll*ret*a%mod);
    		return ret;
    	}
    	inline int add(int x,int y){return (x+y>=mod)?x+y-mod:x+y;}
    	inline int dec(int x,int y){return (x-y<0)?x-y+mod:x-y;}
    }
    using namespace Math;
    
    namespace Container{
    	struct poly{
    		vector<int>v;
    		inline int& operator[](int x){while(x>=v.size())v.push_back(0);return v[x];}
    		inline poly(int x=0):v(1){v[0]=x;}
    		inline int size(){return v.size();}
    		inline void resize(int x){v.resize(x);}
    		inline void mem(int l,int lim){
    			int t=min(lim,(int)v.size());
    			for(int i=l;i<t;++i) v[i]=0;
    		}
    	};
    }
    using namespace Container;
    
    namespace basic{
    	int r[N],Wn[N];
    	inline void NTT(int lim,poly& f,int tp){
    		for(int i=0;i<lim;++i) if(i<r[i]) swap(f[i],f[r[i]]);
    		for(int mid=1;mid<lim;mid<<=1){
    			int len=mid<<1,wn=ksm(tp==1?gg:giv,(mod-1)/len);
    			Wn[0]=1;for(int i=1;i<mid;++i) Wn[i]=1ll*Wn[i-1]*wn%mod;
    			for(int l=0;l+len-1<lim;l+=len){
    				for(int k=l;k<=l+mid-1;++k){
    					int w1=f[k],w2=1ll*Wn[k-l]*f[k+mid]%mod;
    					f[k]=add(w1,w2);f[k+mid]=dec(w1,w2);
    				}
    			}
    		}
    	}
    	inline poly poly_mul(int n,int m,poly f,poly g){
    		int lim=1,len=0;
    		while(lim<(n+m)) lim<<=1,len++;
    		for(int i=0;i<lim;++i) r[i]=(r[i>>1]>>1)|((i&1)<<(len-1));
    		f.mem(n,lim);g.mem(m,lim);
    		NTT(lim,f,1);NTT(lim,g,1);
    		for(int i=0;i<lim;++i) f[i]=1ll*f[i]*g[i]%mod;
    		NTT(lim,f,-1);
    		int iv=ksm(lim,mod-2);
    		for(int i=0;i<lim;++i) f[i]=1ll*f[i]*iv%mod;
    		return f;
    	}
    }
    using namespace basic;
    int n,m;poly f,g;
    int main(){
    	scanf("%d%d",&n,&m);++n;++m;
    	for(int i=0;i<n;++i) scanf("%d",&f[i]);
    	for(int i=0;i<m;++i) scanf("%d",&g[i]);
    	f=poly_mul(n,m,f,g);
    	for(int i=0;i<n+m-1;++i) printf("%d ",f[i]);
    	return 0;
    } 
    

    至此,我们已经能够快速完成特殊模数下的多项式乘法,但是多项式还有更多的运算等待我们去探究,多项式基础运算

  • 相关阅读:
    ConcurrentHashMap之实现细节
    Java 开发 2.0: 用 Hadoop MapReduce 进行大数据分析
    mapreduce从wordcount开始
    centos 5.5 安装mysql 5.5 全程详细记录 RPM方式安装
    使用GDAL工具对OrbView3数据进行正射校正
    centos 5.5 mysql5.5 乱码
    netty vs mina netty和mina的区别
    VC欣赏、家人是阻力,极客化、国际化——90后创业生态
    悲惨而又丢人的创业经历:草根创业者含恨倾诉为什么失败
    悲惨而又丢人的创业经历:草根创业者含恨倾诉为什么失败
  • 原文地址:https://www.cnblogs.com/tqxboomzero/p/14164745.html
Copyright © 2020-2023  润新知