多项式多点求值
重点写了我认为的疑难点,若其他部分由疑问或含糊不清,欢迎提出,积极改正
题面:给一n次多项式\(F(x)\),求m个数代入的值
构造关于\(x_0\)的函数,使得代入\(x_0\)后,值为\(0\),则有\(G(x)=x-x_0\)。做多项式取模\(F(x)=Q(x)G(x)+R(x)\),\(F(x_0)=Q(x_0)G(x_0)+R(x_0)\),\(R(x_0)\)只有常数项,即答案
先分治算出\(\prod_{i=1}^{m/2}(x-x_i)\)和\(\prod_{i=m/2+1}^{m}(x-x_i)\),再分治取模。复杂度\(O(nlogn+mlog^2m)\)
进行优化,多项式除法中
\[Q_r(x)=F_r(x){G_r}^{-1}(x) \pmod {x^{n-m+1}}
\]
设\(G(x)=\prod_{i=1}^m(x-x_i),G_0(x)=\prod_{i=1}^{m/2}(x-x_i),G1(x)=\prod_{i=m/2+1}^{m}(x-x_i)\)
由\(G(x)=G_0(x)G_1(x)\),得\({G_{0r}}^{-1}={G_r}^{-1}(x)G_1(x)\)
\[\begin{aligned}Q_{0r}(x)&=F_r(x){G_{0r}}^{-1}(x)\\
&=F_r(x){G_r}^{-1}(x)G_1(x)\\
&=Q(x)G_1(x)
\end{aligned}\]
在分治最底层计算答案时,我们只需要常数项,系数翻转后的最高项。又由于运算在\(\pmod {x^n}\)意义下进行,我们考虑的是对n-1项贡献。分治过程中只保留最高的区间长度项,这些可能在最底层贡献最高项。
具体表现为目前区间长度为n的多项式,这一层递归会乘项数为n/2的多项式。至少n/2项才会对n-1项有贡献,现在因为砍掉前n/2项,原n-1项变为n/2-1项,则下一层长度为n/2的递归区间,考虑对n/2-1项的贡献。
流程:
- 分治预处理\(G_r(x)\)
- 求最外层逆元,算出最外层\(Q_r(x)=F_r(x){G_r}^{-1}(x)\),并保留后n项
- 分治计算\(Q_r(x)\),最底层最高项即常数项记为\(c_i\),这一位答案为\([x^0]F(x)-c_ix_i\)
点击查看代码
#include<bits/stdc++.h>
using namespace std;
const int MAX=4e5+10;
#define int long long
int n,m,f[MAX],x[MAX],ans[MAX],a[MAX],b[MAX],A[MAX],bl,bc,rev[MAX],inv;
const int mod=998244353;int iG[MAX];
inline int read(){
int x=0,f=1;char c=getchar();
while(c>'9'||c<'0'){if(c=='-')f=-1;c=getchar();}
while(c>='0'&&c<='9'){x=(x<<3)+(x<<1)+(c^'0');c=getchar();}
return x*f;
}
inline void eva(int*,int*,int*,int);
signed main(){
n=read();m=read();
for(int i=0;i<=n;++i) f[i]=read();
for(int i=1;i<=m;++i) x[i]=read();
eva(f,x,ans,max(n,m));
for(int i=1;i<=m;++i) printf("%lld\n",ans[i]);
}
int *g[MAX<<2],*h[MAX<<2],bin[MAX<<6],*np(bin);
inline int power(int a,int b){
int res=1;
while(b){
if(b&1) res=res*a%mod;
a=a*a%mod;b>>=1;
}return res;
}inline void work(int n){
bl=1;bc=0;
while(bl<=n) bl<<=1,++bc;
for(int i=0;i<bl;++i)
rev[i]=(rev[i>>1]>>1)|((i&1)<<bc-1);
}inline void NTT(int *a,int n,int op){
for(int i=0;i<n;++i)
if(i<rev[i]) swap(a[i],a[rev[i]]);
for(int i=1;i<n;i<<=1){
int wn=power(3,(op*(mod-1)/(i<<1)+mod-1)%(mod-1));
for(int j=0;j<n;j+=i<<1){
int w=1;
for(int k=0;k<i;++k){
int x=a[j+k],y=w*a[j+k+i]%mod;
a[j+k]=(x+y)%mod;a[j+k+i]=(x-y+mod)%mod;
w=w*wn%mod;
}
}
}if(op==-1){
inv=power(n,mod-2);
for(int i=0;i<n;++i) a[i]=a[i]*inv%mod;
}
}inline void mul(int *f,int *g,int *c,int n){
work(n<<1);
for(int i=0;i<n;++i) a[i]=f[i],b[i]=g[i];
for(int i=n;i<bl;++i) a[i]=b[i]=0;
NTT(a,bl,1);NTT(b,bl,1);
for(int i=0;i<bl;++i) a[i]=a[i]*b[i]%mod;
NTT(a,bl,-1);
for(int i=0;i<n;++i) c[i]=a[i];
}void solve(int *f,int *g,int n){
if(n==1){g[0]=power(f[0],mod-2);return;}
solve(f,g,n>>1);work(n);
for(int i=0;i<n;++i) a[i]=f[i],b[i]=g[i];
for(int i=n;i<bl;++i) a[i]=b[i]=0;
NTT(a,bl,1);NTT(b,bl,1);
for(int i=0;i<bl;++i) a[i]=a[i]*b[i]%mod*b[i]%mod;
NTT(a,bl,-1);
for(int i=0;i<n;++i) g[i]=(2*g[i]-a[i]+mod)%mod;
}
void solve1(int pos,int l,int r){
g[pos]=np;np+=(r-l+2)*2;h[pos]=np;np+=(r-l+2)*2;
if(l==r){g[pos][0]=1;g[pos][1]=(mod-x[l])%mod;return;}
int mid=l+r>>1;solve1(pos<<1,l,mid);solve1(pos<<1|1,mid+1,r);
mul(g[pos<<1],g[pos<<1|1],g[pos],r-l+2);
}void solve2(int pos,int l,int r,int *ans){
if(l==r){ans[l]=h[pos][0];return;}
int mid=l+r>>1;
reverse(g[pos<<1|1],g[pos<<1|1]+r-mid+1);
mul(h[pos],g[pos<<1|1],A,r-l+1);
for(int i=0;i<mid-l+1;++i) h[pos<<1][i]=A[(r-mid)+i];
solve2(pos<<1,l,mid,ans);
reverse(g[pos<<1],g[pos<<1]+mid-l+2);
mul(h[pos],g[pos<<1],A,r-l+1);
for(int i=0;i<r-mid;++i) h[pos<<1|1][i]=A[(mid-l+1)+i];
solve2(pos<<1|1,mid+1,r,ans);
}inline void eva(int *a,int *b,int *c,int n){
solve1(1,1,n);int M=1;while(M<=n) M<<=1;
solve(g[1],iG,M);
for(int i=0;i<=n;++i) A[i]=f[n-i];
mul(A,iG,A,n+n+1);reverse(A,A+n+n+1);
for(int i=0;i<n;++i) h[1][i]=A[n+1+i];
solve2(1,1,n,c);
for(int i=1;i<=n;++i) c[i]=(a[0]+c[i]*b[i])%mod;
}