【算法笔记】Lagrange 插值法

合集 - OI算法笔记(13)

1.【算法笔记】K.M.P.字符串匹配2022-01-182.【算法笔记】Tarjan 算法 · 上2022-01-213.【算法笔记】Manacher 算法2022-01-314.【算法笔记】二维凸包2022-01-305.【算法笔记】扩展 Euclid 算法2022-01-276.【算法笔记】Tarjan 算法 · 下2022-01-24

7.【算法笔记】Lagrange 插值法2022-04-14

8.【算法笔记】中国剩余定理（CRT）与扩展中国剩余定理（exCRT）2022-04-129.【算法笔记】Baby Step Giant Step（BSGS）及其扩展2022-04-0910.【算法笔记】Gauss 消元 · 行列式2022-04-0311.【算法笔记】Euler-Fermat 定理2022-03-1912.【算法笔记】快速 Fourier 变换（FFT）2022-03-1313.【算法笔记】乘法逆元2022-02-01

收起

• 本文总计约 5000 字，阅读大约需要 30 分钟。
• 警告！警告！警告！本文有大量的 $\LaTeX$ 公式渲染，可能会导致加载异常缓慢！

前言

又一次没有前言，以后要不就不写前言了吧，怪麻烦的 QwQ。

题目引入

根据代数基本定理，我们知道 $(n+1)$ 个横坐标互不相同的点可以唯一确定一个最高次项为 $n$ 的多项式。例如：

点 $(1,2)$ 和 $(2,3)$ 就可以唯一确定多项式 $f(x)=x+1$；点 $(1,1),(2,4),(3,9)$ 可以唯一确定多项式 $g(x)=x^2$。

那么，给定任意 $n$ 个点，你能计算出其所唯一确定的多项式 $f(x)$ 吗？其中 $x\le 10^4$。

暴力怎么做

当然，我们依旧是要想一想，暴力怎么做？

假如这 $n$ 个点分别为 $(x_1,y_1),(x_2,y_2),\cdots,(x_n,y_n)$，那么我们可以定义其确定的多项式为 $f(x)=a_0+a_1x+a_2x^2+\cdots+a_{n-1}x^{n-1}$。

分别代入这 $n$ 个点，得到一个方程：

$$\begin{cases} y_1=a_0+a_1x_1+a_2x_1^2+\cdots+a_{n-1}x_1^{n-1},\\ y_2=a_0+a_1x_2+a_2x_2^2+\cdots+a_{n-1}x_2^{n-1},\\ \cdots\\ y_n=a_0+a_1x_n+a_2x_n^2+\cdots+a_{n-1}x_n^{n-1}.\\ \end{cases}$$

当然，这个方程，我们就可以用 Gauss 消元法来解决，复杂度为 $\mathcal{O}(n^3)$。

『哇啊，那这个复杂度太高了吧？能过 $10^4$ 的数据吗？』

显然是不能的，所以，我们要做时间复杂度的优化。

Lagrange 插值

为了优化时间复杂度，著名的数学家 $\mathcal{J-L.Lagrange}$ 给出了一个更加高效的算法，即 Lagrange 插值法。

我们知道，对任意多项式 $f(x)$，若 $x=a$，那么 $f(x)\equiv f(a)\pmod{(x-a)}$。

因此，我们分别代入 $(x_1,y_1),(x_2,y_2),\cdots,(x_n,y_n)$，就有了：

$$\begin{cases} f(x)\equiv y_1\pmod{(x-x_1)},\\ f(x)\equiv y_2\pmod{(x-x_2)},\\ \cdots\\ f(x)\equiv y_n\pmod{(x-x_n)}.\\ \end{cases}$$

我们使用中国剩余定理，得到：

$$M=\prod_{i=1}^n (x-x_i),\ M_i=\dfrac{M}{x-a_i},\ t_i=\prod_{i\neq j}\dfrac{1}{x-x_j},\\ f(x)=\sum_{i=1}^n y_iM_it_i。$$

因此，有：

$$f(x)=\sum_{i=1}^n y_i\prod_{i\neq j}\dfrac{x-x_j}{x_i-x_j}.$$

上述多项式就是 Lagrange 插值法的公式，通过这个公式，我们可以通过 $n$ 个点值求出其唯一确定的多项式。

通过这个公式，求多项式的时间复杂度为 $\Theta(n^2)$。

代码

这个算法在洛谷上的模板为 P4781【模板】拉格朗日插值，因为我们只需要求 $f(k)$ 的值，所以直接把 $x=k$ 代入到上式中即可知：

$$f(k)=\sum_{i=1}^n y_i\prod_{j\neq i}\dfrac{k-x_j}{x_i-x_j}.$$

于是就可以直接模拟上述公式了，代码如下：

#include <cstdio>
using namespace std;

typedef long long ll;
const ll MOD = 998244353;
const int MAXN = 5001;

int n;
ll x[MAXN], y[MAXN], k, ans;

ll inv(ll val) {  //Fermat 小定理求逆元 
    val = (val % MOD + MOD) % MOD;

    ll p = MOD - 2, res = 1;
    
    while(p) {
        if(p & 1) res = (res * val) % MOD;
        val = (val * val) % MOD; p >>= 1;
    }
    return res;
}

int main(void) {
    scanf("%d%lld", &n, &k);
    
    for(int i = 1; i <= n; ++i) 
        scanf("%lld%lld", &x[i], &y[i]);
        
    for(int i = 1; i <= n; ++i) {  //插值主体 
        ll tmp = y[i];
        
        for(int j = 1; j <= n; ++j) {
            if(i == j) continue;
            tmp = (tmp * ((k - x[j]) % MOD + MOD) % MOD * inv(x[i] - x[j])) % MOD; 
        }
        
        ans = (ans + tmp) % MOD;
    }
    
    printf("%lld", ans);
    return 0;
}
//by CaO

重心 Lagrange 插值

Lagrange 插值固然效率高，但我们发现其有一个缺点：不可维护。

所谓不可维护，就是如果我们在已知的点的基础上新增一些未知的点，那么我们就需要整个重新套用一次公式。

那么，单次重构的时间复杂度依旧为 $\Theta(n^2)$。

『如果我要进行 $m$ 次重构，那么总时间复杂度就是 $\Theta(mn^2)$ 了，如果 $n\le 10^3$，$m\le 10^3$，那么时间复杂度就不够了！』

所以，单次重构的时间复杂度太高，那么我们就要想办法降低重构的复杂度。

我们手玩一下上面的那个式子：

$$f(x)=\sum_{i=1}^n y_i\prod_{i\neq j}\dfrac{x-x_j}{x_i-x_j}.$$

我们可以把求积符号拆分：

$$\begin{aligned} f(x) = \sum_{i=1}^n \left(y_i\prod_{i\neq j}(x-x_j)\prod_{i\neq j}\dfrac{1}{x_i-x_j}\right)\\ \end{aligned}$$

在求和号里面同时乘 $(x-x_i)$ 并除以 $(x-x_i)$。就有：

$$f(x) = \sum_{i=1}^n \left(y_i\prod_{j=1}^n (x-x_j)\cdot\dfrac{1}{x-x_i}\prod_{i\neq j}\dfrac{y_i}{x_i-x_j}\right).$$

注意到 $\prod_{j=1}^n x-x_j$ 是一个常数，提出来，则有：

$$\begin{aligned} f(x) &= \prod_{j=1}^n(x-x_j)\left(\sum_{i=1}^n y_i\dfrac{1}{x-x_i}\prod_{i\neq j}\dfrac{1}{x_i-x_j}\right)\\ &= \prod_{i=1}^n(x-x_i)\left(\sum_{i=1}^n \dfrac{y_i}{\prod_{i\neq j}(x_i-x_j)}\cdot \dfrac{1}{x-x_i}\right)\\ \end{aligned}$$

定义：

$$A(x)=\prod_{i=1}^n (x-x_i),\\ W_i= \dfrac{y_i}{\prod_{i\neq j}(x_i-x_j)}.$$

原式即变为：

$$f(x)=A(x)\sum_{i=1}^n\dfrac{W_i}{x-x_i}.$$

我们会发现，每次新插值时，维护 $A(x)$ 的复杂度为 $\mathcal{O}(1)$ 的，维护 $W(1),W(2),\cdots,W(n+1)$ 的总复杂度为 $\mathcal{O}(n)$ 的。所以单次重构的时间复杂度为 $\mathcal{O}(n)$。

这种算法被称为重心 Lagrange 插值法，代码留给读者作为练习（其实是我太菜了 QwQ）。

例题

本题目列表会持续更新。

• [P4187]【模板】拉格朗日插值
• [CF622F]The Sum of the k-th Powers