问题
给出(n)次多项式(A(x)),(m)次多项式(B(x)),求多项式(D(x)),(R(x))使得$$A(x)=B(x)D(x)+R(x)$$,满足(degle n-m,deg R<m)。
即求多项式(A(x))对(B(x))的带余除法。
做法
首先,对于任意多项式$$P(x)=sum _{i=0}^na_ix^i$$,显然有$$P^R(x)=x^nP(frac{1}{x})=sum {i=0}^na{n-i}x^i$$。
这就是说,通过这个操作,我们可以把多项式的系数翻转。例如:
[P(x)=x^2+2x+3 \
x^2P(frac{1}{x})=3x^2+2x+1
]
对(A(x)=B(x)D(x)+R(x))进行这个操作:
[A(x)=B(x)D(x)+R(x) \
A(frac{1}{x})=B(frac{1}{x})D(frac{1}{x})+R(frac{1}{x}) \
x^nA(frac{1}{x})=x^nB(frac{1}{x})D(frac{1}{x})+x^nR(frac{1}{x}) \
x^nA(frac{1}{x})=x^mB(frac{1}{x})x^{n-m}D(frac{1}{x})+x^{n-m+1}x^{m-1}R(frac{1}{x}) \
A^R(x)=B^R(x)D^R(x)+x^{n-m+1}R^R(x)
]
注意到(x^{n-m+1}R^R(x))的最低次项至少为(n-m+1)次,我们对整个式子取模:
[A^R(x)=B^R(x)D^R(x)mod x^{n-m+1}
]
可以这样做是因为(A)和(B)是直接给出的,不会出问题,而(deg Dle n-m),(deg D^Rle n-m),所以也不会出问题。
这样就可以轻松地算出$$D^R=A_R(B^R)^{-1}$$,运用多项式求逆元的方法。之后把(D^R)转换回(D),带入减一减就能得到(R)的值。
整个做法的关键就是消去余数。
代码
#include<cstdio>
#include<cctype>
#include<cstring>
#include<algorithm>
using namespace std;
typedef long long giant;
const int q=998244353;
const int basic=3;
const int maxn=1<<17|1;
const int maxj=18;
int read() {
int x=0,f=1;
char c=getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=-1;
for (;isdigit(c);c=getchar()) x=x*10+c-'0';
return x*f;
}
int A[maxn],B[maxn],C[maxn],D[maxn],R[maxn];
int AR[maxn],BR[maxn],IBR[maxn],DR[maxn],T[maxn];
int rev[maxn],Wn[maxj][2];
void print(int a[],int len) {
len=max(len,0);
for (int i=0;i<=len;++i) printf("%d ",a[i]);
puts("");
}
int mi(int x,int y) {
int ret=1;
for (;y;y>>=1,x=(giant)x*x%q) if (y&1) ret=(giant)ret*x%q;
return ret;
}
int inv(int x) {
return mi(x,q-2);
}
void ntt(int a[],int len,bool op=0) {
for (int i=0;i<len;++i) if (i<rev[i]) swap(a[i],a[rev[i]]);
for (int i=2,the=1;i<=len;i<<=1,++the) {
int wn=Wn[the][op];
for (int j=0;j<len;j+=i) {
int w=1;
for (int k=j;k<j+i/2;++k,w=(giant)w*wn%q) {
int u=a[k]%q,v=(giant)a[k+i/2]*w%q;
a[k]=(u+v)%q,a[k+i/2]=((u-v)%q+q)%q;
}
}
}
if (op) {
int iv=inv(len);
for (int i=0;i<len;++i) a[i]=(giant)a[i]*iv%q;
}
}
void TIME(int a[],int la,int b[],int lb,int c[]) {
int len=max(la,lb),M,xj;
for (M=1,xj=0;M<=(len<<1);M<<=1,++xj);
for (int i=0;i<M;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(xj-1));
ntt(a,M);
ntt(b,M);
for (int i=0;i<M;++i) c[i]=(giant)a[i]*b[i]%q;
ntt(c,M,1);
ntt(a,M,1);
ntt(b,M,1);
}
void INV(int a[],int c[],int len,int mod) {
if (mod==1) {
c[0]=inv(a[0]);
return;
}
INV(a,c,len,(mod+1)>>1);
int M,xj,chang=max(len,mod);
for (M=1,xj=0;M<=(chang<<1);M<<=1,++xj);
for (int i=0;i<M;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(xj-1));
memset(T,0,sizeof T);
copy(c,c+chang,T);
ntt(a,M);
ntt(T,M);
for (int i=0;i<M;++i) {
c[i]=(2ll-(giant)a[i]*T[i]%q)%q;
if (c[i]<0) c[i]+=q;
c[i]=(giant)c[i]*T[i]%q;
if (c[i]<0) c[i]+=q;
}
ntt(c,M,1);
ntt(a,M,1);
for (int i=mod;i<M;++i) c[i]=0;
}
int main() {
#ifndef ONLINE_JUDGE
freopen("7.in","r",stdin);
freopen("my.out","w",stdout);
#endif
Wn[0][0]=Wn[0][1]=1;
for (int i=1;i<maxj;++i) {
Wn[i][0]=mi(basic,(q-1)>>i);
Wn[i][1]=inv(Wn[i][0]);
}
int n=read();
for (int i=0;i<=n;++i) A[i]=read();
int m=read();
for (int i=0;i<=m;++i) B[i]=read();
if (n<m) {
puts("0");
print(A,n);
return 0;
}
for (int i=0;i<=n-m;++i) AR[i]=A[n-i],BR[i]=B[m-i];
INV(BR,IBR,m,n-m+1);
TIME(AR,n-m,IBR,n-m,DR);
for (int i=0;i<=n-m;++i) D[i]=DR[n-m-i];
TIME(B,m,D,n-m,C);
for (int i=0;i<m;++i) R[i]=((giant)(A[i]-C[i])%q+q)%q;
print(D,n-m);
print(R,m-1);
return 0;
}