写在前面

由于作者是一个初二蒟蒻，有一些地方可能存在问题，请多指教。喷轻点
你以为我会先讲插值吗？
当然了！

问题引入

所有数列都有规律——拉格朗日

众所周知，多项式插值定理告诉我们 n+1 个点唯一确定一个 n 次多项式。但是很多时候，我们需要根据这个多项式求值。
下面问题来了，已知 n 个在多项式 f(x)=\sum_{i=0}^{n}a_i x^i 上的点 \{(x_0,y_0),(x_1,y_1),\cdots,(x_n,y_n)\}，其中 f(x_i)=y_i,\forall_{i\in[1,n]\cap \mathbb{N}}，求 f(k)。
这就是祸害OIer万年的多项式插值问题。

初步思路

不会？？？[高斯消元](https://45475.blog.luogu.org/linear-equation-group)之前的日报讲得挺好的，至于秦九韶算法……就是将 $f(x)=\sum_{i=0}^{n}a_i x^i$ 转换为 $f(x)=a_0+x(a_1+x(a_2+x(a_3+\cdots+x(a_{n-1}+a_nx))))$。将 $n$ 次加法和 $\frac{n(n+1)}{2}$ 次乘法转换为对至多做 $n$ 次乘法和 $n$ 次加法的算法。 # 拉格朗日插值 先 $orz$ 拉格朗日一下： > 约瑟夫·拉格朗日（Joseph-Louis Lagrange，1736~1813）全名为约瑟夫·路易斯·拉格朗日，法国著名数学家、物理学家。1736年1月25日生于意大利都灵，1813年4月10日卒于巴黎。他在数学、力学和天文学三个学科领域中都有历史性的贡献，其中尤以数学方面的成就最为突出。 > 拉格朗日总结了18世纪的数学成果，同时又为19世纪的数学研究开辟了道路，堪称法国最杰出的数学大师。同时，他的关于月球运动（三体问题）、行星运动、轨道计算、两个不动中心问题、流体力学等方面的成果，在使天文学力学化、力学分析化上，也起到了历史性的作用，促进了力学和天体力学的进一步发展，成为这些领域的开创性或奠基性研究。 > ——摘自百度百科 然后让我们想想为什么慢？因为我们知道的太多了。 它要求 $f(k)$，而我们直接求出了 $f(x)$，再代入求 $f(k)$。 我们发现求出 $f(x)$ 的部分好像是多余的，考虑直接求出 $f(k)$。 于是拉格朗日就提出了拉格朗日插值法。先上公式： $$L_n(x)=\sum_{i=0}^{n}y_i\ell_i(x),\ell_i(x)=\prod_{i\neq j}^{n}\frac{x-x_j}{x_i-x_j}$$ $\ell_i(x)$ 叫做拉格朗日基本多项式（插值基函数），$L_n(x)$ 叫做拉格朗日插值多项式。 乍一看好像十分懵逼，我们一点一点地分析它。 首先显然 $L_n(x)$ 是一个 $n$ 次的多项式，然后再让我们分析一下 $\ell_i(x)$，我们可以发现它有一些神奇的性质： ------------ ### $\ell_i(x_i)=1

证明：

\ell_i(x_i)=\prod_{i\neq j}^{n}\frac{x_i-x_j}{x_i-x_j}=\prod_{i\neq j}^{n}1=1

\ell_i(x_k)=0,\forall_{k\neq i}

证明：

\ell_i(x_k)=\prod_{i\neq j}^{n}\frac{x_k-x_j}{x_i-x_j}=\frac{x_k-x_k}{x_i-x_k}\cdot\prod_{i\neq j\neq k}^{n}\frac{x_k-x_j}{x_i-x_j}=0\times\prod_{i\neq j\neq k}^{n}\frac{x_k-x_j}{x_i-x_j}=0

知道了这些性质之后，让我们想一想 L_n(x) 的正确性，证明如下：

L_n(x_j)=\sum_{i=0}^{n}y_i\ell_i(x_j)=\sum_{i\neq j}^{n}y_i\ell_i(x_j)+y_j\ell_j(x_j)=y_j

又因为 L_n(x) 是一个 n 次的多项式，n+1 个点唯一确定一个 n 次多项式。
所以拉格朗日插值是正确的。
我们不需要将 L_n(x) 求出，我们只需要求出 L_n(k) 就好了。
代码：

//求L_n(k) 
double L_n_k=0;
for (int i=1;i<=n;i++)
{
    double ell_i_k=0;
    for (int j=1;j<=n;j++)
        if (i!=j) ell_i_k*=(k-x[j])/(x[i]-x[j]);
    L_n_k+=y[i]*ell_i_k;
}

然后切一道板题——P4781 【模板】拉格朗日插值

当 x_i 取值连续时的做法

如果 x_i 取值连续的，即 x_i=i，我们可以做到 O(n)。
首先先看一下式子：

L_n(x)=\sum_{i=0}^{n}y_i\prod_{i\neq j}^{n}\frac{x-x_j}{x_i-x_j}=\sum_{i=0}^{n}y_i\prod_{i\neq j}^{n}\frac{x-j}{i-j}

要求 L_n(k)，先设 k 的前缀积 pre_i=\prod_{j=0}^ik-j 与后缀积 suf_i=\prod_{j=i}^{n}k-j。
那么式子变成：

L_n(k)=\sum_{i=0}^{n}y_i\cdot\frac{pre_{i-1}\cdot suf_{i+1}}{(-1)^{n-i}\cdot i!\cdot (n-i)!}

预处理阶乘、pre、suf，这个问题就是 O(n) 的啦。
代码：

//当x_i=i时求L_n(k) 
double L_n_k=0;
for (int i=1;i<=n;i++)
    if ((n-i)%2) L_n_k+=y[i]*((pre[i-1]*suf[i+1])/(-fac[i]*fac[n-i]));
    else L_n_k+=y[i]*((pre[i-1]*suf[i+1])/(fac[i]*fac[n-i]));

重心拉格朗日插值法（第一型）

我们发现拉格朗日插值有一点点小问题：当增加一个点时需要重新 O(n^2) 计算，如果要多次增加点时，拉格朗日插值就gg了。 设 \ell(x)=\prod_{i=0}^{n}x-x_i。 式子变成 L_n(x)=\sum_{i=0}^{n}y_i\frac{\ell(x)}{(x-x_i)\cdot\prod_{j\neq i}x_i-x_j}=\ell(x)\cdot\sum_{i=0}^{n}\frac{y_i}{(x-x_i)\cdot\prod_{j\neq i}x_i-x_j} 设重心权 w_i=\frac{y_i}{\prod_{j\neq i}x_i-x_j}，那么式子变成 L_n(x)=\ell(x)\cdot\sum_{i=0}^{n}\frac{w_i}{x-x_i}。 所以每加入一个点，我们用 O(n) 的时间复杂度计算出重心权 w_i，即可求出新的 L_n(x)。

于是就可以进行加入、撤销的操作了，它有向后稳定性。

重心拉格朗日插值法（第二型）

我们采用重心拉格朗日插值法（第一型）来插值多项式 g(x)\equiv 1，可得：

g(x)=\ell(x)\cdot\sum_{i=0}^{n}\frac{w_i}{x-x_i}

所以

\frac{L_n(x)}{g(x)}=\frac{\sum_{i=0}^{n}\frac{w_i}{x-x_i}\cdot y_i}{\sum_{i=0}^{n}\frac{w_i}{x-x_i}}

又因为 g(x)\equiv 1，所以

L_n(x)=\frac{\sum_{i=0}^{n}\frac{w_i}{x-x_i}\cdot y_i}{\sum_{i=0}^{n}\frac{w_i}{x-x_i}}

这就是重心拉格朗日插值法（第二型）或是叫做真正的重心拉格朗日插值公式。它比重心拉格朗日插值法（第一型）更好计算，常数更小，向前稳定，并且勒贝格常数很小，用 x_i=\cos(\frac{(2i+1)\pi}{2(n+1)}) 就是切比雪夫插值节点插值效果好，可以很好地模拟给定的函数，使得插值点个数趋于无穷时，最大偏差趋于零。

牛顿插值法

先 orz 牛顿一下：

艾萨克·牛顿（1643年1月4日—1727年3月31日）爵士，英国皇家学会会长，英国著名的物理学家，百科全书式的“全才”，著有《自然哲学的数学原理》、《光学》。

他在1687年发表的论文《自然定律》里，对万有引力和三大运动定律进行了描述。这些描述奠定了此后三个世纪里物理世界的科学观点，并成为了现代工程学的基础。
他通过论证开普勒行星运动定律与他的引力理论间的一致性，展示了地面物体与天体的运动都遵循着相同的自然定律；为太阳中心说提供了强有力的理论支持，并推动了科学革命。

在力学上，牛顿阐明了动量和角动量守恒的原理，提出牛顿运动定律。

在光学上，他发明了反射望远镜，并基于对三棱镜将白光发散成可见光谱的观察，发展出了颜色理论。他还系统地表述了冷却定律，并研究了音速。

在数学上，牛顿与戈特弗里德·威廉·莱布尼茨分享了发展出微积分学的荣誉。他也证明了广义二项式定理，提出了“牛顿法”以趋近函数的零点，并为幂级数的研究做出了贡献。

在经济学上，牛顿提出金本位制度。

                           ——摘自百度百科

反正就是大佬

然后让我们设被插函数为 f(x)，已知节点 \{(x_i,y_i)\},i\in [0,n]。
然后我们定义差商（均差）的概念：即导数的近似值，对等步长的离散函数 f(x) ，其 n 阶差商就是它的 n 阶差分与其步长的 n 次幂的比值。 ——摘自百度百科——差商

大概意思就是设 f[x_0,x_1,\cdots,x_k] 为 f(x) 的 k 阶差商，那么 f[x_0,x_1,\cdots,x_k]=\frac{f[x_1,x_2,\cdots,x_k]-f[x_0,x_1,\cdots,x_{k-1}]}{x_k-x_0}，f(x) 在 x_i 点处的零阶差商为 f(x_i)，即 f[x_0]=f(x_0)=y_0。

所以

0$ 阶差商 $f[x_0]=f(x_0) 1$ 阶差商 $f[x_0,x_1]=\frac{f[x_1]-f[x_0]}{x_1-x_0} 2$ 阶差商 $f[x_0,x_1,x_2]=\frac{f[x_1,x_2]-f[x_0,x_1]}{x_2-x_0}

……

n$ 阶差商 $f[x_0,x_1,\cdots,x_n]=\frac{f[x_1,x_2,\cdots,xn]-f[x_0,x_1,\cdots,x_{n-1}]}{x_n-x_0}

我们考虑存在未知点 x 的差商，移项可得

f[x]=f(x) f[x,x_0](x-x_0)+f[x_0]=f[x] f[x,x_0,x_1](x-x_1)+f[x_0,x_1]=f[x,x_0] \cdots\cdots f[x,x_0,x_1,\cdots,x_n](x-x_n)+f[x,x_0,x_1,\cdots,x_{n}]=f[x,x_0,x_1,\cdots,x_{n-1}]

将后面的项带入前面得

\footnotesize f(x)=f[x_0]+f[x_0,x_1](x-x_0)+f[x_0,x_1,x_2](x-x_0)(x-x_1)+\cdots+f[x_0,\cdots,x_n](x-x_0)\cdots(x-x_{n-1})+f[x,x_0,\cdots,x_n](x-x_0)\cdots(x-x_n)

因为最后一项无论 x 取哪一个都会是 0，可以去掉。
所以

f(x)=\sum_{i=0}^{n}f[x_0,x_1,\cdots,x_i]\cdot\prod_{j=0}^{i}(x-x_j)+R_n(x)=N_n(x)+R_n(x)

这就是牛顿插值公式。N_n(x) 叫做牛顿插值多项式，R_n(x) 叫做牛顿插值公式余项（截断误差）。在OI没什么意义。
实现的时候便可以增量构造，我们时刻维护 x_n 结尾的 n 个差商和多项式 \prod(x-x_i) 即可，差商用定义的公式递推，多项式因为每次乘二项式，直接递推，两个东西都可以 O(n) 维护。总复杂度为 O(n^2)。
不过它比拉格朗日插值更好的是，它可以插出多项式的系数。当给出的节点等距时，使用差分能达到更优的效果，其计算更加简便。
代码：

void insert(int xn, int yn)//插入点(xn,yn) 
{
    x[++n]=xn;
    y[n]=yn,now^=1,dt[now][0]=yn;//x^1,如果x是偶数,那么答案是x+1,如果x是奇数,那么答案是x-1
    for (int i=1;i<n;i++) dt[now][i]=dt[now][i-1]-dt[now^1][i-1]/x[n]-x[n-i];//维护差商 
    for (int i=n-1;i>=1;i--) ploy[i]=ploy[i-1]-x[n-1]*ploy[i];//维护多项式 
    ploy[0]*=-x[n-1];
    for(register int i=0;i<n;i++) ans[i]+=ploy[i]*dt[now][n-1];//累计答案 
}

泰勒插值

这是一个不太常用的办法，先看看泰勒多项式：

f(x)=\sum_{i=0}^{n}\frac{f^{(i)}(x_0)}{i!}(x-x_0)^i+R_n(x)

其中 f^{(i)}(x) 表示 f(x) 的 i 阶导。
因而，泰勒插值的条件就是已知 0-n 阶的导数 f^{(k)}(x_0),k\in[1,n]\cap\mathbb{N}，求 f(x)。
方法很简单，拿泰勒多项式倒过来用就好了，时间复杂度 O(n^2)。
又慢，限制又多。
不过有的时候还是比较方便的。

多项式快速插值

先给出板题——P5158 【模板】多项式快速插值

我们兴高采烈地打出牛顿插值/拉格朗日插值，然而……TLE……
于是我们考虑优化拉格朗日插值。
先给出了 n+1 个点的集合 X=\{(x_i,y_i):0\leq i\leq n\}，现在要求一个 n 次多项式 A(x) 满足 \forall(x,y)\in X,A(x)=y。
考虑分治，我们把要插值的点分成两半：

X_0=\left\{(x_i,y_i):1\leq i\leq \left \lfloor \frac{n}{2}\right \rfloor \right \} X_1=\left\{(x_i,y_i):\left \lfloor \frac{n}{2}\right \rfloor<i\leq n\right \}

假设已用 X_0 插值出多项式 A_0(x)，求 A(x)。
设 P(x)=\prod_{i=0}^{\left \lfloor \frac{n}{2}\right \rfloor}(x-x_i)。
那么构造 A(x)，使其满足

A(x)=A_1(x)P(x)+A_0(x)

有没有想到FFT？
其中 A_1(x) 是一个未知的多项式，由于 \forall(x,y)\in X_0,P(x)=0,A_0(x)=y，这样的话就满足 X_0 的点都在 A(x) 上，问题就变成要将 X_1 内的点插值，使得 \forall(x_i,y_i)\in X_1,y_i=A_1(x_i)P(x_i)+A_0(x_i) 化简之后得到

A_1(x_i)=\frac{A_0(x_i)-y_i}{P(x_i)}

这样就得到了新的待插值点，利用同样的方法求出插值出 A_1 然后合并就可以了，最终复杂度是 T(n)=2T(\frac{n}{2})+O(n\log^2 n)=O(n\log^3 n)。

更快的多项式快速插值

还是这个板题——P5158 【模板】多项式快速插值

如果 O(n\log^3 n)还是满足不了大佬AKIOI的需求时，我们就继续换方法。

注意：以下内容可能引起不适，建议学过极限的同学学习。

给出了 n+1 个点的集合 X=\{(x_i,y_i):0\leq i\leq n\}，要求一个 n 次多项式 A(x) 满足 \forall(x,y)\in X,A(x)=y。
先回顾一下拉格朗日插值公式：

L_n(x)=\sum_{i=0}^{n}y_i\ell_i(x)=\sum_{i=0}^{n}y_i\prod_{i\neq j}^{n}\frac{x-x_j}{x_i-x_j}

暴力求解是 O(n^3) 的，所以考虑优化。
首先求解分母的部分，设 g(x)=\prod_{i=0}^{n}(x-x_i)，可以用分治FFT，O(n\log^2n)
问题就变成了求当 x=x_i 时，\frac{g(x)}{x-x_i} 的值。
怎么样，有思路了吗？是不是一下就可以出解了？
当然不是
简单思考一下就可以发现，当 x=x_i 时，g(x)=x-x_i=0，所以\frac{g(x)}{x-x_i}=\frac{0}{0}……
那变一下问题，改成求 \lim_{x\to x_i}\frac{g(x)}{x-x_i} 的值，似乎可做。
这个时候，就要用洛必达法则啦！

洛必达法则（零比零型）

若函数 f(x) 和 g(x) 满足下列条件：

1. \lim_{x\to a}f(x)=\lim_{x\to a}g(x)=0
2. f(x)$、$g(x)$ 在点 $a$ 的某去心邻域内可导,且 $\lim_{x\to a}g'(x)\neq 0

那么

\lim_{x\to a}\frac{f(x)}{g(x)}=\lim_{x\to a}\frac{f'(x)}{g'(x)}

根据洛必达法则，\lim_{x\to x_i}\frac{g(x)}{x-x_i}=\lim_{x\to x_i}\frac{g'(x)}{(x-x_i)'}=g'(x)

那么我们就要求出 g'(x_0),g'(x_1),\cdots,g'(x_n)，使用多项式多点求值即可，成功优化到 O(n\cdot n\log^2n+n\log n)=O(n^2\log^2 n)！

【高能预警】下面的公式非常密集！！！

设 w_i=\frac{y_i}{\prod_{j\neq i}x_i-x_j}，分子已知，分母可以通过我们上边的方法求解。

原式变成

f(x)=\sum_{i=0}^{n}w_i\prod_{j=0,j\neq i}^{n}(x-x_j)

设 P_l(x)=\prod_{i=0}^{\frac{n}{2}}x-x_i,P_r(x)=\prod_{i=\frac{n}{2}+1}^{n}x-x_i。
再设 f_l(x)=\sum_{i=0}^{\frac{n}{2}}w_i\prod_{j=0,j\neq i}^{\frac{n}{2}}(x-x_j),f_r(x)=\sum_{i=\frac{n}{2}+1}^{n}w_i\prod_{j=\frac{n}{2}+1,j\neq i}^{n}(x-x_j)。
那么

f(x)=P_r(x)f_l(x)+P_l(x)f_r(x)

预处理出 P(x)，分治即可。最底层直接返回 w_i。
超长公式警告
总的来说，设 mid=\lfloor\frac{n}{2}\rfloor，f_{l,r}(x)表示 \{(x_i,y_i),i\in[l,r]\cap\mathbb{N}\}这些点插出来的多项式，则

f_{l,r}(x)=\sum_{i=l}^r \frac{y_i}{g'(x)}\prod_{j=l,j\neq i}^r(x-x_j)

所以

\small f_{l,r}(x)=\left[\prod_{j=mid+1}^r(x-x_j)\right]\left[\sum_{i=l}^{mid}\frac{y_i}{g'(x)}\prod_{j=l,j\neq i}^{mid}(x-x_j)\right]+\left[\prod_{j=l}^{mid}(x-x_j)\right]\left[\sum_{i=mid+1}^r\frac{y_i}{g'(x)}\prod_{j=mid+1,j\neq i}^r(x-x_j)\right]

所以

f_{l,r}=\left[\prod_{j=mid+1}^r(x-x_j)\right]f_{l,mid}(x)+\left[\prod_{j=l}^{mid}(x-x_j)\right]f_{mid+1,r}(x)

递归求解即可，时间复杂度 O(n\log^2n)，常数极大。

不过代码就不给了，因为我还没写出来
我太菜了

埃尔米特（\text{Hermite}）插值

有时候，我们不仅要求插值函数在给定节点上函数值重合，而且要求若干阶导数也重合，即，要求插值函数\varphi(x)满足：

\varphi^{(j)}(x_i)=f^{(j)}(x_i),\forall i\in[1,n]\cap\mathbb{N},j\in[1,m]\cap\mathbb{N}

以两点三次 \text{Hermite} 插值为例，讲解\text{Hermite}插值。

已知：x_0,x_1,f(x_0)=y_0,f(x_1)=y_1,f'(x_0)=y_0',f'(x_1)=y_1'，求三次多项式 H(x)，满足上述条件。

考虑模仿拉格朗日多项式的思想，设

H_3(x)=a_0f_0(x)+a_1f_1(x)+b_0g_0(x)+b_1g_1(x)

其中，f_0(x),f_1(x),g_0(x),g_1(x) 是三次多项式，且满足

f_i(x_j)=g_i'(x_j)=t_{ij},f_i'(x_j)=g_i(x_j)=0

带入插值条件得

H_3(x)=y_0f_0(x)+y_1f_1(x)+y_0'g_0(x)+y_1'g_1(x)

对于 f_0(x)，可得 f_0(x_0)=1,f_0(x_1)=f_0'(x_0)=f_0'(x_1)=0。
设 f_0(x)=(ax+b)(\frac{x-x_1}{x_0-x_1})^2
因为 f_0(x_0)=1,f_0'(x_0)=0，所以解得 a=-\frac{2}{x_0-x_1},b=\frac{3x_0-x_1}{x_0-x_1}=1+\frac{2x_0}{x_0-x_1}，所以 f_0(x)=(1-2\frac{x-x_0}{x_0-x_1})(\frac{x-x_1}{x_0-x_1})^2。
同理可得 f_1(x)=(1-2\frac{x-x_1}{x_1-x_0})(\frac{x-x_0}{x_1-x_0})^2。

对于 g_0(x)，可得 g_0'(x_0)=1,g_0(x_0)=g_0(x_1)=g_0'(x_1)=0。
设 g_0(x)=(cx+d)(\frac{x-x_1}{x_0-x_1})^2
因为 g_0(x_0)=0,g_0'(x_0)=1，所以解得 c=1,d=-x_0，所以 f_0(x)=(x-x_0)(\frac{x-x_1}{x_0-x_1})^2。
同理可得 g_1(x)=(x-x_1)(\frac{x-x_0}{x_1-x_0})^2。
带入得到插值公式：

H_3(x)=y_0(1-2\frac{x-x_0}{x_0-x_1})(\frac{x-x_1}{x_0-x_1})^2+y_1(1-2\frac{x-x_1}{x_1-x_0})(\frac{x-x_0}{x_1-x_0})^2+y_0'(x-x_0)(\frac{x-x_1}{x_0-x_1})^2+y_1'(x-x_1)(\frac{x-x_0}{x_1-x_0})^2

根据上述计算思想，n+1 个节点，唯一确定 2n+1 次埃尔米特插值多项式：

H_{2n+1}=\sum_{k=0}^{n}\left\{[1-2(x-x_k)L_k'(x_k)]f_kL_k(x)+(x-x_k)f_k'L_k^2(x)\right\}

这里的 L_k(x) 就是拉格朗日插值多项式。

不过其实我们日常生活不会用到IOI不考，所以实用性不强，但是它可以解决一些毒瘤的问题。

课后作业一道思考题

这是一道经典例题：给你两个正整数 k,n，请求出 \sum_{i=0}^{n}i^k 的通项公式，用各种插值法做一遍。
然后 A 了这道题——CF622F。
祝你们好运！

$$\Huge\color{purple}\text{完结撒花！！！}$$

参考资料

• 多项式的多点求值与快速插值——miskcoo

• @总结 - 4@ 多项式的多点求值与快速插值——Tiw_Air_OAO

• 题解 P5158 【【模板】多项式快速插值】——bztMinamoto

• 题解 P4781 【【模板】拉格朗日插值】——attack

• 从牛顿插值法到泰勒公式——马同学高等数学

• 【插值】插值方法原理详解——Y.D