【模板】任意模数多项式乘法

题目描述

给定 \(2\) 个多项式 \(F(x), G(x)\) ，请求出 \(F(x) \times G(x)\)。

系数对 \(p\) 取模，且不保证 \(p\) 可以分解成 \(p = a \cdot 2^k + 1\) 之形式。

题解

可以用快速数论变换NTT算法，关键在于取的那个素数。

由于系数最大为 \(10^5 \times 10^{9+9} = 10^{23}\) 所以可以用 __int128 来存储，取一个可以表示为 \(2^a \times b+1\) 的大素数。

用 Python 来找大素数，取两个随机数 \(a\) 和 \(b\) ，然后用 rabinMiller 来判断是否为素数。

找到了一个大素数：\(2^{73} \times 39+1\)

将这个数的欧拉函数分解质因数，从 \(2\) 开始判断，发现 \(17\) 是它的一个原根。

然后就是 __int128 乘法，先转换成 long double 求出两个数相乘除以模数，再转换回 __int128。

代码：

#include<bits/stdc++.h>
#define int long long
#define ll __int128
#define mxn 300000
#define ld long double
#define pb push_back
#define mkp make_pair
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define rept(i,a,b) for(int i=a;i<b;++i)
using namespace std;
inline int read(){
	int x=0;char ch=getchar();
	while(!isdigit(ch))ch=getchar();
	while(isdigit(ch))x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
	return x;
}
const ll md=((ll)1<<73)*39+1;
ll mul(ll x,ll y){
	ll a=(ld)x*y/md;
	ll ans=x*y-a*md;
	ans%=md;
	if(ans<0)ans+=md;
	return ans;
}
ll power(ll x,ll y){
	ll ans=1;
	for(;y;y>>=1){
		if(y&1)ans=mul(ans,x);
		x=mul(x,x);
	}
	return ans;
}
inline ll mod(ll a,ll b){
	return mul(a,power(b,md-2));
}
int n,m,k,s,p,rev[mxn];
ll f[mxn],g[mxn];
void ntt(ll *a,int n,int flag){
	rept(i,0,n)if(i<rev[i])swap(a[i],a[rev[i]]);
	for(int h=1;h<n;h<<=1){
		ll x,y,s=power(17,(md-1)/2/h);
		for(int j=0;j<n;j+=h<<1){
			ll w=1;
			for(int k=j;k<j+h;++k){
				x=a[k],y=mul(w,a[k+h]);
				a[k]=(x+y)%md;
				a[k+h]=(x-y+md)%md;
				w=mul(w,s);
			}
		}
	}
	if(flag==-1){
		ll p=power(n,md-2);
		reverse(a+1,a+n);
		rept(i,0,n)a[i]=mul(a[i],p);
	}
}
signed main(){
	n=read(),m=read(),p=read();
	rep(i,0,n)f[i]=read();
	rep(i,0,m)g[i]=read();
	for(k=0,s=1;s<n+m+1;s<<=1,++k);
	rept(i,0,s)rev[i]=(rev[i>>1]>>1)|((i&1)<<(k-1));
	ntt(f,s,1);ntt(g,s,1);
	rept(i,0,s)f[i]=mul(f[i],g[i]);
	ntt(f,s,-1);
	rep(i,0,n+m)printf("%lld ",(int)(f[i]%p));
	return 0;
}