NTT学习笔记 & luogu P4245 【模板】任意模数NTT

背景:

NTT \text{NTT} NTT的学习咕了好久,反正有 FFT \text{FFT} FFT,直到多项式全家桶的出现。



题目传送门:

https://www.luogu.org/problemnew/show/P4245



正题:

想要学习 NTT \text{NTT} NTT请看:https://blog.csdn.net/linjiayang2016/article/details/80341958#_423

想知道更好的质数及其原根请看:http://blog.miskcoo.com/2014/07/fft-prime-table

这里只讲一下如何做到任意模数。
我们将给定的模数用三个模数( 469762049 , 998244353 , 1004535809 469762049,998244353,1004535809 469762049,998244353,1004535809)来代替,因为它们的原根均为 3 3 3
为什么要选这么大,我们考虑到 NTT \text{NTT} NTT时出现的数最大为 max ⁡ ( n , m ) ⋅ ( p − 1 ) 2 ≤ 1 0 23 \max(n,m)\cdot(p-1)^2≤10^{23} max(n,m)(p1)21023,而这三个数的乘积大于 1 0 23 10^{23} 1023,因此可以完美代替。
那么我们做 3 3 3次多项式乘法,其中 NTT \text{NTT} NTT时模数每次取一个。
这样我们得到 3 3 3个结果,我们相当于解一个同余方程组,因为模数互质,因而用中国剩余定理即可。
考虑到 m o d 0 ∗ m o d 1 ∗ m o d 2 mod_0*mod_1*mod_2 mod0mod1mod2这个数并不能用 long long \text{long long} long long存下,因此我们先做前面 2 2 2个,最后一个暴力合并即可。



代码:

#include<cstdio>
#include<cstring>
#include<algorithm>
typedef long long LL;
typedef long double LD;
const LL Mod[3]={469762049,998244353,1004535809},G=3;
using namespace std;
	LL d1[300010],d2[300010],a[300010],b[300010],ans[3][300010];
	int limit=1,n,m,k,l=0,mod,inv_G;
	int r[300010];
LL dg(LL x,LL k,LL mod)
{
	if(!k) return 1;
	LL op=dg(x,(k>>1),mod);
	if(k&1) return op*op%mod*x%mod; else return op*op%mod;
}
LL times(LL a,LL b,LL mod)
{
	a%=mod,b%=mod;
	return ((a*b-(LL)((LL)((LD)a/mod*b+1e-3)*mod))%mod+mod)%mod;
}
LL inv(LL x)
{
	return dg(x,mod-2,mod);
}
void NTT(LL *now,int op)
{
	for(int i=0;i<limit;i++)
		if(i<r[i]) swap(now[i],now[r[i]]);
	for(int mid=1;mid<limit;mid<<=1)
	{
		LL wn=dg(op==1?G:inv_G,(mod-1)/(mid<<1),mod);
		for(int j=0;j<limit;j+=(mid<<1))
		{
			LL w=1;
			for(int k=0;k<mid;k++,w=(w*wn)%mod)
			{
				LL x=now[j+k],y=w*now[j+k+mid]%mod;
				now[j+k]=(x+y)%mod;
				now[j+k+mid]=(x-y+mod)%mod;
			}
		}
	}
}
void work(int op)
{
	mod=Mod[op];
	inv_G=inv(G);
	memset(a,0,sizeof(a));
	memset(b,0,sizeof(b));
	for(int i=0;i<=n;i++)
		a[i]=d1[i]%mod;
	for(int i=0;i<=m;i++)
		b[i]=d2[i]%mod;
	NTT(a,1);
	NTT(b,1);
	for(int i=0;i<=limit;i++)
		a[i]=a[i]*b[i];
	NTT(a,-1);
	LL INV=inv(limit);
	for(int i=0;i<=n+m;i++)
		ans[op][i]=a[i]*INV%mod;
}
void merge()
{
	LL M=Mod[0]*Mod[1];
	for(int i=0;i<=n+m;i++)
	{
		LL A=(times(ans[0][i]*Mod[1]%M,dg(Mod[1]%Mod[0],Mod[0]-2,Mod[0]),M)+
			  times(ans[1][i]*Mod[0]%M,dg(Mod[0]%Mod[1],Mod[1]-2,Mod[1]),M))%M;
		LL K=((ans[2][i]-A)%Mod[2]+Mod[2])%Mod[2]*dg(M%Mod[2],Mod[2]-2,Mod[2])%Mod[2];
		printf("%lld ",((K%k)*(M%k)%k+A%k)%k);
	}
}
int main()
{
	scanf("%d %d %d",&n,&m,&k);
	while(limit<=n+m)
		limit<<=1,l++;
	for(int i=1;i<=limit;i++)
		r[i]=((r[i>>1]>>1)|((i&1)<<(l-1)));
	for(int i=0;i<=n;i++)
		scanf("%lld",&d1[i]);
	for(int i=0;i<=m;i++)
		scanf("%lld",&d2[i]);
	work(0),work(1),work(2);
	merge();
}
  • 1
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值