https://blog.csdn.net/jwg2732/article/details/123179632?spm=1001.2014.3001.5501
任意模数NTT
已知多项式\(F(x),G(x)\),求\(H(x)=F(x)G(x)(mod\;P)\)
其中\(P\)是给定模数,不保证是NTT模数
我们先取三个NTT模数
常用的是(998244353,469762049,1004535809)
因为他们的原根都是3
我们分别用这三个模数做NTT,然后再用中国剩余定理合并
不过这三个模数乘起来爆longlong
所以我们必须一步一合并
我们可以对于每一位分别考虑
设\(A=998244353,B=469762049,C=1004535809\)
设这一位的真实值是\(x\)
我们已知
\[\left\{\begin{matrix}x\equiv x_1(mod\;A)\\x\equiv x_2(mod\;B)\\x\equiv x_3(mod\;C)\end{matrix}\right.
\]
也就是
第一个方程的通解是\(x_1+k_1A\)
也就是
\[x_1+k_1A\equiv x_2(mod \;B)
\]
移项可得
\[k_1\equiv\frac{x_2-x_1}{A}
\]
我们代入可以的到\(X=x_1+k_1A\)
那么前两个方程的通解是
\[X+K(AB)
\]
类似的代入第三个方程
\[X+K(AB)\equiv x_3(mod\;C)
\]
移项
\[K\equiv \frac{x_3-X}{AB}
\]
那么我们可以算出一个特解\(x_0\)
总的通解式就是
\[x_0+k_0ABC(mod\;ABC)
\]
因为\(ABC\)一定大于\(x_0\)所以答案就是\(x_0\)
因为要做\(9\)次NTT,常数很大
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const LL A = 998244353;
const LL B = 1004535809;
const LL C = 469762049;
#define Poly vector<LL>
#define len(x) (int)x.size()
#define turn(x,v) x.resize(v)
const int N = 4e5+7;
LL Pow(LL a,LL b,LL mod)
{
LL res=1;
while(b)
{
if(b&1) res=res*a%mod;
a=1ll*a*a%mod;
b>>=1;
}
return res;
}
int tr[N];
void Retry(int n)
{
for(int i=0;i<n;i++)
tr[i]=((tr[i>>1]>>1)|((i&1)?(n>>1):0));
}
void NTT(Poly &f,int opt,LL mod)
{
int n=len(f);
LL G=3,IG=Pow(G,mod-2,mod);
for(int i=0;i<n;i++)
if(i<tr[i]) swap(f[i],f[tr[i]]);
for(int k=2;k<=n;k<<=1)
{
int l=(k>>1);
LL w=Pow(opt==1?G:IG,(mod-1)/k,mod);
for(int i=0;i<n;i+=k)
{
LL g=1;
for(int j=i;j<i+l;j++)
{
LL v=f[j+l]*g%mod;
f[j+l]=(f[j]-v+mod)%mod;
f[j]=(f[j]+v)%mod;
g=g*w%mod;
}
}
}
}
Poly Times(Poly f,Poly g,int len,LL mod)
{
NTT(f,1,mod);NTT(g,1,mod);
for(int i=0;i<len;i++)
f[i]=f[i]*g[i]%mod;
NTT(f,-1,mod);
LL invn=Pow(len,mod-2,mod);
for(int i=0;i<len;i++)
f[i]=f[i]*invn%mod;
return f;
}
LL mul(LL a,LL b,LL mod)
{
LL res=0;
while(b)
{
if(b&1) res=(res+a)%mod;
a=(a+a)%mod;
b>>=1;
}
return res;
}
LL Inv(LL x,LL mod)
{
LL a=x,b=mod-2,res=1;
while(b)
{
if(b&1) res=res*a%mod;
a=a*a%mod;
b>>=1;
}
return res;
}
LL Merge(LL x1,LL x2,LL x3,LL mod)
{
LL k1=((x2-x1)%B+B)%B*Inv(A,B)%B;
LL x0=(x1+k1*A%(A*B))%(A*B);
LL k0=((x3-x0)%C+C)%C*Inv(A*B%C,C)%C;
return (x0+k0*A%mod*B%mod)%mod;
}
inline LL read() {
LL x = 0, w = 1;
char ch = 0;
while (ch < '0' || ch > '9') { // ch 不是数字时
if (ch == '-') w = -1; // 判断是否为负
ch = getchar(); // 继续读入
}
while (ch >= '0' && ch <= '9') { // ch 是数字时
x = x * 10 + (ch - '0'); // 将新读入的数字’加’在 x 的后面
// x 是 int 类型,char 类型的 ch 和 ’0’ 会被自动转为其对应的
// ASCII 码,相当于将 ch 转化为对应数字
// 此处也可以使用 (x<<3)+(x<<1) 的写法来代替 x*10
ch = getchar(); // 继续读入
}
return x * w; // 数字 * 正负号 = 实际数值
}
void Put(Poly f,int len)
{
for(int i=0;i<len;i++)
printf("%lld ",f[i]);
printf("\n");
}
Poly FreeMul(Poly f,Poly g,int len,LL mod)
{
Retry(len);
Poly X1=Times(f,g,len,A);
Poly X2=Times(f,g,len,B);
Poly X3=Times(f,g,len,C);
// Put(X1,len);
Poly res;
turn(res,len);
for(int i=0;i<len;i++)
res[i]=Merge(X1[i],X2[i],X3[i],mod);
return res;
}
int n,m,p;
int main()
{
cin>>n>>m>>p;
n++;m++;
Poly f,g;
int k=n+m;
int T;
for(T=1;T<k;T<<=1);
turn(f,T);turn(g,T);
for(int i=0;i<n;i++)
f[i]=read();
for(int i=0;i<m;i++)
g[i]=read();
Poly res=FreeMul(f,g,T,p);
for(int i=0;i<n+m-1;i++)
printf("%lld ",res[i]);
return 0;
}

浙公网安备 33010602011771号