算法简介：同袍的义理

拉格朗日插值法，是一种已知多项式函数 $f(x)$ 上的 $k$ 个点 $(x_i,y_i)$ 求这个多项式的系数的东西。是以法国十八世纪数学家约瑟夫·拉格朗日命名的一种多项式插值方法。

算法演示：僚佐的才巧

算法推导 ：御敌的执定

多项式的系数表示法

众所周知一个 $k$ 次多项式可以被它的 $k+1$ 个系数唯一地确定，即：

这种确定多项式的方法叫做多项式的系数表示法。

多项式的点值表示法

众所周知，一个 $k$ 次多项式同样的可以被图上的 $k+1$ 个点唯一确定。我们设这个多项式的形式为：

那么我们把这几个点带入到这个式子里，就可以得到 $k+1$ 个 $k+1$ 元一次方程组（小学生都知道吧）

这样的一个方程组肯定不会存在一样的，也就是说，我们绝对可以把系数 $a_0,a_1,\dots,a_k$ 全部解出来，也就确定了这个多项式，这种表示法叫做多项式的点值表示法。

系数转换点值表示法

这个就很简单，想要从系数转化为点值，我们只需要对于这个多项式代入 $k+1$ 个不同的 $x$ 值，他们得到的 $y$ 值就是系数表示法转换出来的点至表示法。

点值转换系数表示法

我们看前面点值表示法的推导过程，这里其实就是一个求解 $k+1$ 元方程组的过程。我们就可以通过高斯消元，在 $O(k^3)$ 的时间复杂度内求解出方程。

但是这样的复杂度显而易见是不够优秀的，我们可以尝试考虑用拉格朗日插值法在 $O(k^2)$ 的复杂度内确定这个多项式的系数。我们首先从一个结论入手：

这个结论的正确性也显然，当我们取到 $x = a$ 时，就有 $f(x) -f(a) = 0$ 。

我们把这里的 $f(a)$ 中的 $a$ 全部取成 $x_1,x_2,\dots,x_{k+1}$ 的值，那么得到的结果就是 $y_1,y_2,\dots,y_{k+1}$ 。那么我们就可以把 $f(x)$ 整体带入，得到一个关于 $f(x)$ 的同余方程组：

此时方程组的形式就变得特别像中国剩余定理了。我们想一下中国剩余定理的求解过程：

  1. 计算模数的乘积 $M = \prod\limits_{i=1}^{k+1} (x-x_i)$ 。

  2. 对于每个模数，我们都可以得到一个 $m_i = \cfrac M {x-x_i}=\prod\limits_{i\neq j}(x-x_j)$ 。

  3. 那么现在我们就要计算这 $m_i$ 在第 $i$ 个模数 $(x-x_i)$ 意义下的逆元，$m_i^{-1} = \prod\limits_{i\neq j}\cfrac 1 {x_i - x_j}$ 。

  4. 然后就计算一个 $c_i = \prod\limits_{i\neq j} \cfrac{x-x_j}{x_i-x_j}$ 。

  5. 然后最后的解就是所有 $y_i \times c_i$ 的和，$f(x)$ 的表达式也就是：

不难发现计算这个式子的复杂度就是 $O(k^2)$ 的，搞定。

例题选讲：用臣的久计

大秦为你打开题目传送门

这道题目就是很简单的一个模板题了。所有的值都知道了，只需要求 $f(k)$ 即可。

这里直接写代码。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
// #pragma GCC optimize(2)
// #pragma GCC optimize(3, "Ofast", "inline")
#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
typedef __int128 int128;
typedef pair<int , int > pii;
typedef unsigned long long ull;

namespace FastIO 
{
    // char buf[1 << 20], *p1, *p2;
    // #define getchar() (p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 1 << 20, stdin), p1 == p2) ? EOF : *p1++)
    template<typename T> inline T read(T& x) {
        x = 0;
        int f = 1;
        char ch;
        while (!isdigit(ch = getchar())) if (ch == '-') f = -1;
        while (isdigit(ch)) x = (x << 1) + (x << 3) + (ch ^ 48), ch = getchar();
        x *= f;
        return x;
    }
    template<typename T, typename... Args> inline void read(T& x, Args &...x_) {
        read(x);
        read(x_...);
        return;
    }
    inline ll read() {
        ll x;
        read(x);
        return x;
    }
};
using namespace FastIO;

const int N = 2010;
const int P = 998244353;

ll n, k;
ll x[N], y[N];

inline void Input() {
    read(n, k);
    for(int i = 1; i <= n; i++) {
        read(x[i], y[i]);
    }
}

inline ll Pow(ll a, ll b) {
    ll ans = 1;
    while(b) {
        if(b & 1) ans = ans * a % P;
        b >>= 1;
        a = a * a % P;
    }
    return ans;
}

inline ll Lagrange(int t) {
    ll tot = 0;
    for(int i = 1; i <= n; i++) {
        ll p = y[i], q = 1;
        for(int j = 1; j <= n; j++) {
            if(i != j) {
                p = p * (t - x[j]) % P;
                q = q * (x[i] - x[j]) % P;
            }
        }
        tot = (tot + (p * Pow(q, P - 2) % P + P) % P) % P;
    }
    return (tot + P) % P;
}

inline void Work() {
    printf("%lld\n", Lagrange(k));
}

int main() {
    int T = 1;
    while(T--) {
        Input();
        Work();
    }
    return 0;
}

算法实战：野火的豪烈

查看相关标签喵！