问题
现在给出 n n n 次多项式 A ( x ) A(x) A(x) 和 m m m 次多项式 B ( x ) B(x) B(x),求出它们的商 C ( x ) C(x) C(x) 和余数 D ( x ) D(x) D(x)。也就是说,找出 n − m n-m n−m 次多项式 C ( x ) C(x) C(x) 和次数小于 m m m 的多项式 D ( x ) D(x) D(x) 使得 A ( x ) = B ( x ) C ( x ) + D ( x ) A(x)=B(x)C(x)+D(x) A(x)=B(x)C(x)+D(x)。系数对 998244353 998244353 998244353 取模。
思路
先考虑求商。一个很自然的思路就是求出
B
(
x
)
B(x)
B(x) 的乘法逆元再乘上
A
(
x
)
A(x)
A(x),但我们发现这样求出来的结果是错误的,原因是这样的结果求出来后根本没有余数我们考虑把余数通过某种方式处理掉。
考虑式子
A
(
x
)
=
B
(
x
)
C
(
x
)
+
D
(
x
)
A(x)=B(x)C(x)+D(x)
A(x)=B(x)C(x)+D(x),代入
1
x
\frac 1 x
x1 并两端乘上
x
n
x^n
xn,则可以得到
x
n
A
(
1
x
)
=
x
n
B
(
1
x
)
C
(
1
x
)
+
x
n
D
(
1
x
)
x^nA(\frac 1 x)=x^nB(\frac 1 x)C(\frac 1 x)+x^nD(\frac 1 x)
xnA(x1)=xnB(x1)C(x1)+xnD(x1)如果有
n
n
n 次多项式
F
(
x
)
F(x)
F(x)(不足
n
n
n 次就补
0
0
0),不难发现把
F
(
1
x
)
F(\frac 1 x)
F(x1) 乘上
x
n
x^n
xn 相当于把
F
(
x
)
F(x)
F(x) 的系数左右调转,得到的多项式记作
F
R
(
x
)
F_R(x)
FR(x)。显然,它可以
O
(
n
)
O(n)
O(n) 求出。
有了这个定义,上面的式子就可以简化为
A
R
(
x
)
=
B
R
(
x
)
C
R
(
x
)
+
x
n
−
m
+
1
D
R
(
x
)
A_R(x)=B_R(x)C_R(x)+x^{n-m+1}D_R(x)
AR(x)=BR(x)CR(x)+xn−m+1DR(x)
A
R
(
x
)
≡
B
R
(
x
)
C
R
(
x
)
(
mod
x
n
−
m
+
1
)
A_R(x)\equiv B_R(x)C_R(x)\ (\text{mod}\ x^{n-m+1})
AR(x)≡BR(x)CR(x) (mod xn−m+1)
这样余数的干扰就被排除,我们就可以用乘法逆元求出
C
R
(
x
)
≡
A
R
(
x
)
(
B
R
(
x
)
)
−
1
(
mod
x
n
−
m
+
1
)
C_R(x)\equiv A_R(x)(B_R(x))^{-1}\ (\text{mod}\ x^{n-m+1})
CR(x)≡AR(x)(BR(x))−1 (mod xn−m+1)这个结果恰巧为
n
−
m
n-m
n−m 次多项式,再把它反过来就可以得到
C
(
x
)
C(x)
C(x) 了。
现在求出了
C
(
x
)
C(x)
C(x),余数就可以用
D
(
x
)
=
A
(
x
)
−
B
(
x
)
C
(
x
)
D(x)=A(x)-B(x)C(x)
D(x)=A(x)−B(x)C(x) 求出来了。
代码
#include<iostream>
#include<cstdio>
#include<cstring>
#define re register
#define ll long long
using namespace std;
ll mod=998244353,G[2]={332748118,3},inv;
ll A[450001],B[450001],invB[450001],D[450001];
ll f[450001],g[450001];
int rev[450001],m,n,lenA,lenB;
inline ll qp(ll a,int b){
ll res=1ll;
while(b){
if(b&1) (res*=a)%=mod;
(a*=a)%=mod;
b>>=1;
}
return res;
}
ll w,w_n,x,y,z;
inline void NTT(ll *a,int flag){
for(re int i=1;i<n;++i){
if(rev[i]<i) swap(a[i],a[rev[i]]);
}
for(re int i=1;i<n;i<<=1){
w_n=qp(G[flag],(mod-1)/(i<<1));
for(re int j=0;j<n;j+=(i<<1)){
w=1ll; z=i+j;
for(re int k=j;k<z;++k){
x=a[k];
y=w*a[k+i]; if(y>=mod) y%=mod;
a[k]=x+y; if(a[k]>=mod) a[k]-=mod;
a[k+i]=x-y; if(a[k+i]<0) a[k+i]+=mod;
w*=w_n; if(w>=mod) w%=mod;
}
}
}
if(!flag){
inv=qp(n,mod-2);
for(re int i=0;i<n;++i) (a[i]*=inv)%=mod;
}
return;
}
int main(){
scanf("%d%d",&lenA,&lenB);
for(re int i=0;i<=lenA;++i) scanf("%lld",&A[i]);
for(re int i=lenB;i>=0;--i) scanf("%lld",&B[i]); //直接反过来输入
invB[0]=qp(B[0],mod-2);
for(re int i=1;i<=lenA-lenB;i<<=1){
n=(i<<2);
m+=1; memset(g,0,sizeof(g));
for(re int j=1;j<n;++j) rev[j]=((rev[j>>1]>>1)|((j&1)<<m));
for(re int j=i-1;j<(i<<1);++j) f[j]=B[j];
for(re int j=0;j<i;++j) g[j]=invB[j];
NTT(f,1); NTT(g,1);
for(re int j=0;j<n;++j) g[j]=(((g[j]*g[j])%mod)*f[j])%mod;
NTT(f,0); NTT(g,0);
for(re int j=i;j<(i<<1);++j) invB[j]=((invB[j]<<1)%mod-g[j]+mod)%mod;
} //求乘法逆元
memset(f,0,sizeof(f)); memset(g,0,sizeof(g));
for(re int i=0;i<=lenA-lenB;++i) f[i]=A[lenA-i];
for(re int i=0;i<=lenA-lenB;++i) g[i]=invB[i];
m=1; while((1<<m)<=((lenA-lenB)<<1)) m+=1; n=(1<<m);
for(re int i=1;i<n;++i) rev[i]=((rev[i>>1]>>1)|((i&1)<<(m-1)));
NTT(f,1); NTT(g,1);
for(re int i=0;i<n;++i) (f[i]*=g[i])%=mod; //得到 C_R(x)
NTT(f,0);
for(re int i=0;i<=lenA-lenB;++i) D[i]=f[lenA-lenB-i];
for(re int i=0;i<=(lenB>>1);++i) swap(B[i],B[lenB-i]);
m=1; while((1<<m)<=lenA) m+=1; n=(1<<m);
for(re int i=1;i<n;++i) rev[i]=((rev[i>>1]>>1)|((i&1)<<(m-1)));
NTT(B,1); NTT(D,1);
for(re int i=0;i<n;++i) (B[i]*=D[i])%=mod; //B(x)*D(x)
NTT(B,0); NTT(D,0);
for(re int i=0;i<=lenA-lenB;++i) printf("%lld ",D[i]);printf("\n");
for(int i=0;i<lenB;++i){
A[i]-=B[i];
if(A[i]<0) A[i]+=mod; //求解余数
}
for(re int i=0;i<lenB;++i) printf("%lld ",A[i]);printf("\n");
return 0;
}
谢谢观看,记得点赞