今天终于抽空把这个综（du）合（liu）知识点学了，心力交瘁……

多项式快速插值

给出 nnn 个点 (xi,yi)(x_i,y_i)(xi,yi) ，要求一个次数为 n−1n-1n−1 的多项式 F(x)F(x)F(x) 满足：F(xi)=yiF(x_i)=y_iF(xi)=yi
显然这个多项式是唯一确定的。
根据拉格朗日插值法，我们有：F(x)=∑i=1n∏j̸=i(x−xj)∏j̸=i(xi−xj)yiF(x)=\sum_{i=1}^{n}\frac{\prod_{j\not=i}(x-x_j)}{\prod_{j\not=i}(x_i-x_j)}y_iF(x)=i=1∑n∏j̸=i(xi−xj)∏j̸=i(x−xj)yi
这样我们是可以 O(n2)O(n^2)O(n2) 求的，考虑优化。
我们先考虑对于每个 iii ，如何快速得到 ∏j̸=i(xi−xj)\prod_{j\not=i}(x_i-x_j)∏j̸=i(xi−xj) 。
设 M(x)=∏i=1n(x−xi)M(x)=\prod_{i=1}^{n}(x-x_i)M(x)=∏i=1n(x−xi) ，我们即需求：M(x)x−xi\frac{M(x)}{x-x_i}x−xiM(x)
根据洛必达法则，当 xxx 取 xix_ixi 时，分子分母都等于0，一个0/0型，上下求导得：limx→xiM(x)x−xi=M′(x)lim_{x→x_i}\frac{M(x)}{x-x_i}=M'(x)limx→xix−xiM(x)=M′(x)
于是我们先分治NTT求出 M(x)M(x)M(x) ，再求导得到 M′(x)M'(x)M′(x) ，之后将 x1−nx_{1-n}x1−n 代入多点求值即可！
这一部分的复杂度是 O(nlog2n)O(n\ log^2n)O(n log2n) 的。
求导、多点求值这些前置知识我的博客里都有讲：

多项式的求逆、取模和多点求值学习小记

多项式的ln、exp和快速幂学习小记
又设 Vi=yi∏j̸=i(xi−xj)V_i=\frac{y_i}{\prod_{j\not=i}(x_i-x_j)}Vi=∏j̸=i(xi−xj)yi ，则此时 ViV_iVi 已知，我们要求：F(x)=∑i=1nVi∏j̸=i(x−xj)F(x)=\sum_{i=1}^{n}V_i\prod_{j\not=i}(x-x_j)F(x)=i=1∑nVij̸=i∏(x−xj)
还是分治NTT，设 L(x)=∑i=1n/2(x−xi)L(x)=\sum_{i=1}^{n/2}(x-x_i)L(x)=∑i=1n/2(x−xi) ，R(x)=∑i=n/2+1n(x−xi)R(x)=\sum_{i=n/2+1}^{n}(x-x_i)R(x)=∑i=n/2+1n(x−xi) ，则有：F(x)=∑i=1n/2Vi∏j̸=i,1≤j≤n/2(x−xj)R(x)+∑i=n/2+1nVi∏j̸=i,n/2+1≤j≤n(x−xj)L(x)F(x)=\sum_{i=1}^{n/2}V_i\prod_{j\not=i,1\leq j\leq n/2}(x-x_j)R(x)+\sum_{i=n/2+1}^{n}V_i\prod_{j\not=i,n/2+1\leq j\leq n}(x-x_j)L(x)F(x)=i=1∑n/2Vij̸=i,1≤j≤n/2∏(x−xj)R(x)+i=n/2+1∑nVij̸=i,n/2+1≤j≤n∏(x−xj)L(x)
递归即可求得，递归底层的 F(x)F(x)F(x) 就是 ViV_iVi 。
还有就是这里的 L(x)、R(x)L(x)、R(x)L(x)、R(x) 在多点求值中已经算过了，不用再算一遍啦。
这一部分的复杂度也是 O(nlog2n)O(n\ log^2n)O(n log2n) 。
故总时间复杂度即为 O(nlog2n)O(n\ log^2n)O(n log2n) ，常数很大很大。
模板题：洛谷 P5158 【模板】多项式快速插值

Code

#include<cstdio>
#include<algorithm>
#include<cctype>
using namespace std;
typedef long long LL;
const int N=1e5+5,M=18,G=3,mo=998244353;
int tot;
int xx[N],yy[N],val[N];
int a[N],b[N],c[N],rr[N];//a=b*c+rr
int ra[N],rb[N<<2],irb[N<<2];
int f[N*M<<1],stf[N<<2],enf[N<<2];
int g[N*M<<1],stg[N<<2],eng[N<<2];
int h[N],sth[N],enh[N],f3[N];
int f1[N<<2],f2[N<<2],wn[N<<2],rev[N<<2];
inline int read()
{int X=0,w=0; char ch=0;while(!isdigit(ch)) w|=ch=='-',ch=getchar();while(isdigit(ch)) X=(X<<3)+(X<<1)+(ch^48),ch=getchar();return w?-X:X;
}
void write(int x)
{if(x>9) write(x/10);putchar(x%10+'0');
}
inline int ksm(int x,int y)
{int s=1;while(y){if(y&1) s=(LL)s*x%mo;x=(LL)x*x%mo;y>>=1;}return s;
}
inline void NTT(int *y,int len,int ff)
{for(int i=0;i<len;i++)if(i<rev[i]) swap(y[i],y[rev[i]]);for(int h=2,d=len>>1;h<=len;h<<=1,d>>=1)for(int i=0,k=h>>1;i<len;i+=h)for(int j=0,cnt=0;j<k;j++,cnt+=d){int u=y[i+j],t=(LL)wn[cnt]*y[i+j+k]%mo;y[i+j]=u+t>=mo?u+t-mo:u+t;y[i+j+k]=u-t<0?u-t+mo:u-t;}if(ff==-1){reverse(y+1,y+len);int inv=ksm(len,mo-2);for(int i=0;i<len;i++) y[i]=(LL)y[i]*inv%mo;}
}
void make(int v,int l,int r)
{if(l==r){g[stg[v]=++tot]=mo-xx[l];g[eng[v]=++tot]=1;return;}int mid=l+r>>1,ls=v<<1,rs=ls|1;make(ls,l,mid),make(rs,mid+1,r);int na=eng[ls]-stg[ls]+1,nb=eng[rs]-stg[rs]+1;int len=1,ll=0;while(len<na+nb) len<<=1,ll++;for(int i=0;i<len;i++) rev[i]=rev[i>>1]>>1|(i&1)<<ll-1;int w0=ksm(G,(mo-1)/len);for(int i=wn[0]=1;i<=len;i++) wn[i]=(LL)wn[i-1]*w0%mo;for(int i=0;i<na;i++) f1[i]=g[stg[ls]+i];for(int i=na;i<len;i++) f1[i]=0;for(int i=0;i<nb;i++) f2[i]=g[stg[rs]+i];for(int i=nb;i<len;i++) f2[i]=0;NTT(f1,len,1),NTT(f2,len,1);for(int i=0;i<len;i++) f1[i]=(LL)f1[i]*f2[i]%mo;NTT(f1,len,-1);stg[v]=tot+1;na+=nb-1;for(int i=0;i<na;i++) g[++tot]=f1[i];eng[v]=tot;
}
void getinv(int len,int ll)
{if(len==1){irb[0]=ksm(rb[0],mo-2);return;}getinv(len>>1,ll-1);for(int i=0;i<len;i++) rev[i]=rev[i>>1]>>1|(i&1)<<ll-1;int w0=ksm(G,(mo-1)/len);for(int i=wn[0]=1;i<=len;i++) wn[i]=(LL)wn[i-1]*w0%mo;for(int i=0;i<len>>1;i++) f1[i]=rb[i];for(int i=len>>1;i<len;i++) f1[i]=0;NTT(f1,len,1),NTT(irb,len,1);for(int i=0;i<len;i++) irb[i]=(2-(LL)f1[i]*irb[i]%mo+mo)*irb[i]%mo;NTT(irb,len,-1);for(int i=len>>1;i<len;i++) irb[i]=0;
}
void solve(int v,int l,int r,int fa)
{int na=enf[fa]-stf[fa],nb=eng[v]-stg[v];if(na>=nb){int nc=na-nb;for(int i=0;i<=na;i++) a[i]=f[stf[fa]+i];for(int i=0;i<=nb;i++) b[i]=g[stg[v]+i];for(int i=0;i<=nc;i++) ra[i]=a[na-i];for(int i=0;i<=nb;i++) rb[i]=b[nb-i];for(int i=nc+1;i<=nb;i++) rb[i]=0;int len=1,ll=0;while(len<=nc*2+1) len<<=1,ll++;for(int i=nb+1;i<len;i++) rb[i]=0;for(int i=0;i<len;i++) irb[i]=0,f1[i]=0;getinv(len,ll);for(int i=0;i<=nc;i++) f1[i]=ra[i],f2[i]=irb[i];for(int i=nc+1;i<len;i++) f1[i]=f2[i]=0;NTT(f1,len,1),NTT(f2,len,1);for(int i=0;i<len;i++) f1[i]=(LL)f1[i]*f2[i]%mo;NTT(f1,len,-1);for(int i=0;i<=nc;i++) c[nc-i]=f1[i];for(int i=nc+1;i<nb;i++) c[i]=0;len=1,ll=0;while(len<nb<<1) len<<=1,ll++;for(int i=0;i<len;i++) rev[i]=rev[i>>1]>>1|(i&1)<<ll-1;int w0=ksm(G,(mo-1)/len);for(int i=wn[0]=1;i<=len;i++) wn[i]=(LL)wn[i-1]*w0%mo;for(int i=0;i<nb;i++) f1[i]=b[i],f2[i]=c[i];for(int i=nb;i<len;i++) f1[i]=0,f2[i]=0;NTT(f1,len,1),NTT(f2,len,1);for(int i=0;i<len;i++) f1[i]=(LL)f1[i]*f2[i]%mo;NTT(f1,len,-1);for(int i=0;i<nb;i++) rr[i]=(a[i]-f1[i]+mo)%mo;while(nb>1 && !rr[nb-1]) nb--;stf[v]=tot+1;for(int i=0;i<nb;i++) f[++tot]=rr[i];enf[v]=tot;}else{stf[v]=tot+1;for(int i=stf[fa];i<=enf[fa];i++) f[++tot]=f[i];enf[v]=tot;}if(l==r){val[l]=f[stf[v]];return;}int mid=l+r>>1;solve(v<<1,l,mid,v);solve(v<<1|1,mid+1,r,v);
}
void work(int v,int l,int r)
{if(l==r) return;int mid=l+r>>1,ls=v<<1,rs=ls|1;work(ls,l,mid);work(rs,mid+1,r);int na=enh[l]-sth[l]+1,nb=eng[rs]-stg[rs]+1;int len=1,ll=0;while(len<na+nb) len<<=1,ll++;for(int i=0;i<len;i++) rev[i]=rev[i>>1]>>1|(i&1)<<ll-1;int w0=ksm(G,(mo-1)/len);for(int i=wn[0]=1;i<=len;i++) wn[i]=(LL)wn[i-1]*w0%mo;for(int i=0;i<na;i++) f1[i]=h[sth[l]+i];for(int i=na;i<len;i++) f1[i]=0;for(int i=0;i<nb;i++) f2[i]=g[stg[rs]+i];for(int i=nb;i<len;i++) f2[i]=0;NTT(f1,len,1),NTT(f2,len,1);for(int i=0;i<len;i++) f1[i]=(LL)f1[i]*f2[i]%mo;NTT(f1,len,-1);na+=nb-1;for(int i=0;i<na;i++) f3[i]=f1[i];for(int i=na;i<len;i++) f3[i]=0;na=enh[mid+1]-sth[mid+1]+1,nb=eng[ls]-stg[ls]+1;for(int i=0;i<na;i++) f1[i]=h[sth[mid+1]+i];for(int i=na;i<len;i++) f1[i]=0;for(int i=0;i<nb;i++) f2[i]=g[stg[ls]+i];for(int i=nb;i<len;i++) f2[i]=0;NTT(f1,len,1),NTT(f2,len,1);for(int i=0;i<len;i++) f1[i]=(LL)f1[i]*f2[i]%mo;NTT(f1,len,-1);na+=nb-1;for(int i=0;i<na;i++) h[sth[l]+i]=f3[i]+f1[i]>=mo?f3[i]+f1[i]-mo:f3[i]+f1[i];enh[l]=sth[l]+na-1;
}
int main()
{int n=read();for(int i=1;i<=n;i++) xx[i]=read(),yy[i]=read();make(1,1,n);int m=eng[1]-stg[1];for(int i=0;i<=m;i++) f1[i]=g[stg[1]+i];for(int i=0;i<m;i++) f1[i]=(LL)f1[i+1]*(i+1)%mo;f1[m--]=0;stf[tot=0]=1;for(int i=0;i<=m;i++) f[enf[0]=++tot]=f1[i];solve(1,1,n,0);for(int i=1;i<=n;i++) val[i]=(LL)yy[i]*ksm(val[i],mo-2)%mo;tot=0;for(int i=1;i<=n;i++){sth[i]=enh[i]=++tot;h[tot]=val[i];}work(1,1,n);for(int i=0;i<n;i++) write(h[sth[1]+i]),putchar(' ');return 0;
}

多项式快速插值学习小记相关推荐

多项式相关操作学习笔记
多项式相关操作学习笔记标签: 多项式说在前边记录一下相关的多项式操作,顺便存个模板.(多点求值之后的部分,有点写不动了...留坑留坑多项式定义给定一个环\(R\)(\(R\)通常是交换环, ...
多项式的ln、exp、快速幂和开根学习小记
不妨又学习了一下多项式的求ln.exp.快速幂和开根操作. 这些操作比之前的求逆更上了一层台阶,应用同样很广. 多项式求逆等知识在我的博客里有讲:多项式的求逆.取模和多点求值学习小记多项式ln 给出 ...
@总结 - 4@ 多项式的多点求值与快速插值
目录 @0 - 参考资料@ @1 - 多点求值@ @理论推导@ @参考代码@ @例题与应用@ @2 - 快速插值@ @理论推导@ @(不建议参考的)代码@ @例题与应用@(暂无) @0 - 参考资料@ ...
[多项式算法](Part 4)FWT 快速沃尔什变换学习笔记
其他多项式算法传送门: [多项式算法](Part 1)FFT 快速傅里叶变换学习笔记 [多项式算法](Part 2)NTT 快速数论变换学习笔记 [多项式算法](Part 3)MTT 任意模数FF ...
annoy(快速近邻向量搜索包)学习小记 - pip命令学习与annoy基础使用
1. 写在前面在写fun-rec新闻推荐系统的YouTubeDNN召回的时候, 得到用户向量和新闻向量,基于用户向量,需要从海量新闻里面得到最相似的TopK个新闻, 此时需要用到快速向量检索技术,之 ...
【学习笔记】超简单的多项式快速幂
整理的算法模板合集: ACM模板目录 P5245 [模板]多项式快速幂普通版(a0=1a_0=1a0=1) vector版本AC代码加强版(a0≠1a_0 \neq 1a0=1) 点我看 ...
学习matlab（五）——多项式、插值、极限
针对数据分析和处理,MATL AB提供了大量的函数,非常方便和灵活.本章将详细的介绍利用MATLAB进行一些基本的数据分析,主要包括多项式及其函数,插值,以及函数的极限.MATLAB能够很好的解决多项 ...
多项式全家桶学习笔记【持续更新】
本文完成的时间跨度较长,文风变化可能较大-- 最近更新于2020年2月17日 Part 1.主线乘法前面讲过FFT 然而FFT常数感人,适用范围还窄,比如不能取模于是有了NTT 其实就是取模的F ...
（2012.01.12－2012.04.01）八十二天的学习小记
(2012.01.12-2012.04.01)八十二天的学习小记哈哈,原来又是过了八十二个日子了,真快啊~这次发的学习小记日期记录时间有点长,回看1月份的东西,原来已经隔了八十多个日子了,对于这 ...

多项式快速插值学习小记

多项式快速插值

Code

多项式快速插值学习小记相关推荐

最新文章

热门文章