题解 P5050 【模板】多项式多点求值

· · 题解

整了一天，没救了（

使用了非（多次）取模做法，但是一个点跑 1.2s，难过。

理论推导

约定 n 为 F(x) 的项数，其最高次为 n-1 次。

由余式定理有 F(x_i)=F(z) \mod (z-x_i)，余式 R_i 只有常数项。

设 G_i=x_i+z，Q_i,R_i 的意义与多项式取模中的相同。

由于 G_i 是一次式，则 Q_i 最高次为 n-2 次。

多项式取模告诉我们 Q_{iR}\equiv F_R \cdot G_{iR}^{-1}\ \ \ (\mod z^{n-m+1})。

设 Q_{R(l,r)}=F_R\cdot (\prod_{i=l}^rG_{iR})^{-1}。因为需要用于分治计算，不再考虑原来的模意义。

可以发现有 Q_{R(l,mid)}\equiv Q_{R(l,r)} \cdot \prod_{i=mid+1}^rG_{iR}，Q_{R(mid+1,r)} 同理。于是可以分治解决。

又因为 R_i 只有常数项，我们只需要求出 Q_i 的常数项，即 [z^{n-2}]Q_{iR}，即 Q_{iR} 的最高位。可以发现当分治进入子区间时，由于还能乘上的 G_i 只有 len-1 个，最终能贡献到最高位的只有后 len 个位置，只保存有贡献的位置即可（此处 len 为子区间长度）。

注意此处的后 len 个位置指的是相对于 z^{n-2} 在本区间对应位置而言的，这是令人迷惑的，可以自己画个图感受一下。具体实现见下文。

于是复杂度为 O(n\log^2n)。

具体实现

仅供参考。

首先预处理一下每个分治段的 \prod_{i=l}^{r}G_{iR}。

以便处理，我们钦定 z^{n-2} 在当前区间 (l,r) 的对应位置是 Q_{R(l,r)} 的最高位，即 z^{len-1}，len=r-l+1。

举个例子，分治到左区间时：

\prod_{j=mid+1}^{r}G_{jR})

那我们大概需要将第一层分治的 Q_R 右移一位。

搞清楚这个，大概就可以实现了。

其他两篇题解的代码是真的没看懂，不知道为啥要翻转 G_{R}^{-1}，好像用了一些神秘操作。

Code

给点核心代码。

因为 \prod G_{iR} 常数项必然为 1，所以可以去掉，那么项数就与区间长度相同，于是可以用静态数组保存，而不必动态开点。

此处将 m 补齐至二的次幂，同时满足 n>m，省去大量细节。

EVA->新世纪福音战士

void eva_solve(ull *F,int l,int r,int dep,int id){
    static ull H[maxn],A[maxn],B[maxn],x,y;
    if(l==r) return;
    int len=r-l+1;
    int lim=len<<1;
    int mid=(l+r)>>1; 
    memcpy(H,EVA[dep]+l,8*len);
    memcpy(A+1,EVA[dep+1]+l,4*len),A[0]=1;
    memcpy(B+1,EVA[dep+1]+mid+1,4*len),B[0]=1;
    NTT(H,1,lim,id),NTT(A,1,lim,id),NTT(B,1,lim,id);
    for(int i=0;i<lim;++i){
        x=A[i],y=B[i];
        A[i]=H[i]*y%p[id];
        B[i]=H[i]*x%p[id];
    }
    NTT(A,-1,lim,id),NTT(B,-1,lim,id);
    memcpy(EVA[dep+1]+l,A+(len>>1),4*len);//后len>>1个
    memcpy(EVA[dep+1]+mid+1,B+(len>>1),4*len);
    memset(H,0,8*lim),memset(A,0,8*lim),memset(B,0,8*lim);
    eva_solve(F,l,mid,dep+1,id);
    eva_solve(F,mid+1,r,dep+1,id);
}
void eva_init(ull *X,int l,int r,int dep,int id){
    static ull A[maxn],B[maxn];
    if(l==r){
        EVA[dep][l]=X[l]>0?p[id]-X[l]:-X[l];
        return;
    }
    int len=r-l+1;
    int mid=(l+r)>>1;
    eva_init(X,l,mid,dep+1,id);
    eva_init(X,mid+1,r,dep+1,id);
    int lim=len<<1;
    memcpy(A+1,EVA[dep+1]+l,4*len),A[0]=1;
    memcpy(B+1,EVA[dep+1]+mid+1,4*len),B[0]=1;
    NTT(A,1,lim,id),NTT(B,1,lim,id);
    for(int i=0;i<lim;++i) 
    A[i]=A[i]*B[i]%p[id];
    NTT(A,-1,lim,id);
    memcpy(EVA[dep]+l,A+1,8*len);
    memset(A,0,8*lim);
    memset(B,0,8*lim);
}
void eva(ull *F,ull *G,ull *X,int len,int id){//X即xi
    static ull H[maxn],T[maxn];
    int lim=len<<1;
    eva_init(X,0,len-1,0,id);
    memcpy(H+1,EVA[0],8*len),H[0]=1;
    invf(H,T,len+1,id);
    memcpy(H,F,8*len),H[len]=0;
    reverse(H,H+len);
    NTT(H,1,lim,id),NTT(T,1,lim,id);
    for(int i=0;i<lim;++i)
    H[i]=H[i]*T[i]%p[id];
    NTT(H,-1,lim,id);
    memcpy(EVA[0]+1,H,8*(len-1)),EVA[0][0]=0;//右移一位
    eva_solve(F,0,len-1,0,id);
    for(int i=0;i<len;++i)
    G[i]=add(F[0],EVA[to[len]][i]*X[i]%p[id],p[id]);
    memset(T,0,8*(len+1));
    memset(H,0,8*len);
}