NTT学习笔记

粗略的学习笔记~

NTT(快速数论变换)

• 解决多项式卷积在取模意义下的值，并且不使用复数（浮点数）而是使用整数(在取模的意义下)来减少精度误差。

其实FFT预处理单位根和将较大数据拆成两部分也可以减少精度误差，但是还是没有NTT的精度高，而且不方便取模。

---

这里假设你会了FFT。
不会的话看这里FFT

FFT之所以精度误差大，就是在于单位根为无理数$e$的多少次方，所以必须用浮点数，那么由于计算机的特性，误差会很大。

而这里是在取模的一个前提下，那么如果我们找到了一个新的且为整数的单位根（这里是在取模意义下），我们记其为$G$，那么用这个根去替代原来的$e$，进行运算。由于所有运算就会变成整数运算，所以精度误差大大降低。

如何找到这么一个神奇$G$满足原来的单位根的性质呢？

---

原根

• 百度百科

原根：对于一个数字$m$，当$gcd(a,m)=1$（$a$和$m$的最大公约数为1，也就是互质）时，我们记满足$a^xequiv 1( m mod m)$式子的最小$x$为$ord_m(a)$，然后当$ord_m(a)=varphi(m)$时称$a$为$m$的原根。（不一定所有的数字都有原根）。

扩展:
费马小定理:$a^{p-1}equiv 1( m mod p)$，$p$为素数。（这个可以用来求逆元，也可写为$a^{varphi(p)}equiv 1( m mod p)$因为$varphi(p)=p-1$当$p$为质数）
费马小定理还有个用途就是当幂的指数特别大时，需要对$p$取模，我们可以将指数对$p$取模然后求（因为根据费马小定理，将模出来等于一的可以先模掉）。
离散对数:满足$a^xequiv k( m mod m)$我们将$x$称为$a$在模$m$的意义下的离散对数，我们可以简记为$x=ind_a(k)$。

性质（我们将原根一般用$g$表示）:

1. $g^i m mod m(iin[0,m-1])$可以将$[0,m-1]$内的数字全部取到。
2. $(g^{frac{m-1}{p_i}} m mod m)≠1$，其中$p_i$为$m$的素数因子。

---

• 原根的求法:
  因为一般原根很小，我们采用暴力枚举的方法。
  对于一个数字$m$我们求它的原根，根据性质2，我们先线性筛出$m$的所有素数因子，然后从$2$开始枚举原根$g$，如果对于一个枚举的$g$，满足对于$m$的所有质因数$p_i$，满足性质2，那么这个数就是原根，然后break掉就好了。

要求原根的一个题目BZOJ3992

---

转换

通过离散对数的定义，我们可以将原来FFT单位根的方程$x^n=1$写为$x^nequiv 1( m mod m)$，那么我们求的新的整数$x$不就是原根了吗？

所以我们用$g$来代替了单位根，那么$我们仍旧学FFT一样记omega_n^k=g^{frac{k(m-1)}{n}}(kin Z,)$，$m$为模数，这样就代替了单位复数根，但是只有当$m$为素数时，原根才一定存在，所以我们取模数一般取素数，而且这种质数一般取费马质数(NTT质数)。

一些常用的取模素数以及原根:IN There

---

但是模数并不总是素数，而是一些普通的数字，那么我们怎么办呢?

那么对于模数，我们可以用几个乘积大于最大卷积的系数的模数，分别求取取模后的卷积，然后用中国剩余定理合并即可。

不会中国剩余定理的->百度

但是有时卷积的系数会非常之大，导致你不能直接合并，否则爆long long，所以我们要使用一个小技巧:

假如我们有三个模数$m_1,m_2,m_3$，其中$m_1 imes m_2$不会爆long long，但是$m_1 imes m_2 imes m_3$会爆，那么我们先直接合并前两个得到$M=m_1 imes m_2$的卷积，然后我们就会剩下以下这个方程式:

$egin{cases} x equiv a_1 ({mod} M) \ xequiv a_2 ({mod} m_3) end{cases}$

那么转换一下就可以写成：

$egin{cases} x = a_1+k_1 imes M \ x= a_3+k_3 imes m_3 end{cases}$

然后可以写成：

$x=k_1 imes M+a_1=k_3 imes m_3+a_3$

那么：

$k_1 imes Mequiv a_3-x( m mod m_3)$

$k_3 imes m_3$被$m_3$模掉了。

所以我们解出$k_1$，然后带入原式即可算出答案$x$，那么这个问题就解决了。

最后我们求出的是在$m mod(m_1 imes m_2 imes m_3)$的意义下的卷积，如果原来让你$m mod P(P<(m_1 imes m_2 imes m_3))$，那么你再将答案对$P$取模即可。

---

模板代码：

#include<cstdio>
#include<cstring>
#include<algorithm>
#define ll long long
using namespace std;

const int M=3e6+10;
const ll p1=469762049ll,p2=998244353ll,p3=1004535809ll,g=3;
const ll Mod=468937312667959297ll;
//三个乘积大于模数且大于最大卷积的系数P*(n-1)^2，原根均为3
int n,m,p;
ll t1[M],t2[M],t3[M],t4[M],ans[3][M];
int R[M];
int fpow(int a,int b,int mod){
	int ans=1;
	for(;b;b>>=1,a=(1ll*a*a)%mod){
		if(b&1) ans=(1ll*ans*a)%mod;
	}
	return ans;
}
ll mul(ll a,ll b,ll mod){
    a%=mod,b%=mod;
    return ((a*b-(ll)((ll)((long double)a/mod*b+1e-3)*mod))%mod+mod)%mod;
}

void DNT(ll *a,int n,int f,int mod){
	for(int i=0;i<n;i++) if(i<R[i]) swap(a[i],a[R[i]]);
	for(int i=2;i<=n;i<<=1){
		int now=i>>1;
		int wn=fpow(g,(mod-1)/i,mod);
		if(f==-1) wn=fpow(wn,mod-2,mod);//逆变换时，逆元，或者采用如下的交换法
		for(int j=0;j<n;j+=i){
			int w=1,x,y;
			for(int k=j;k<j+now;k++,w=1ll*w*wn%mod){
				x=a[k],y=1ll*w*a[k+now]%mod;
				a[k]=(x+y)%mod;a[k+now]=(x-y+mod)%mod;
			}
		}
	}
	if(f==-1){
		int inv=fpow(n,mod-2,mod);
		for(int i=0;i<=n;i++){
			a[i]=1ll*a[i]*inv%mod;
		}
//		for(int i=1;i<=n/2;i++) swap(a[i],a[n-i]);//交换法，因为g^(-1)=g^(n-1)，所以可以先正变换，翻转就为逆变换。
	}
}

void NTT(ll *a,ll *b,int k,int mod){
	DNT(a,k,1,mod),DNT(b,k,1,mod);
	for(int i=0;i<=k;i++) a[i]=1ll*a[i]*b[i]%mod;
//	DNT(a,k,-1,mod);
}
int k;
void mcpy(int d){
	for(int i=0;i<=k;i++) ans[d][i]=t3[i];
	if(d==2) return;
	memset(t3,0,sizeof(t3));memset(t4,0,sizeof(t4));
	for(int i=0;i<=n;i++) t3[i]=t1[i];
	for(int i=0;i<=m;i++) t4[i]=t2[i];
}

void REV(){
	DNT(ans[0],k,-1,p1);
	DNT(ans[1],k,-1,p2);
	DNT(ans[2],k,-1,p3);
}

int lg;
int main(){
	scanf("%d%d%d",&n,&m,&p);
	for(int i=0;i<=n;i++) scanf("%lld",&t1[i]),t3[i]=t1[i];
	for(int i=0;i<=m;i++) scanf("%lld",&t2[i]),t4[i]=t2[i];
	k=1;while(k<=n+m)k<<=1,++lg;
	for(int i=0;i<k;i++) R[i]=(R[i>>1]>>1)|((i&1)<<(lg-1));
	NTT(t3,t4,k,p1);
	mcpy(0);
	NTT(t3,t4,k,p2);
	mcpy(1);
	NTT(t3,t4,k,p3);
	mcpy(2);
	REV();
	for(int i=0;i<=n+m;i++){//中国剩余定理合并+后面推导出的公式
		ll x=((mul(1ll*ans[0][i]*p2%Mod,fpow(p2%p1,p1-2,p1),Mod))+
			  (mul(1ll*ans[1][i]*p1%Mod,fpow(p1%p2,p2-2,p2),Mod)))%Mod;
		ll y=((((ans[2][i]-x)%p3+p3)%p3)*fpow(Mod%p3,p3-2,p3))%p3;
//		printf("%lld %lld %lld ",x,y,ans[2][i]);
		printf("%lld ",(1ll*(y%p)*(Mod%p)%p+x%p)%p);
	}
	return 0;
}

---

一些例题

【模板】任意模数NTT

HDU5829 Rikka with Subset - 题解

[SDOI2015]序列统计

---

End

目前初学，还有很多不懂，望大佬指点错误！

---