Skip to content

拉格朗日插值法

2024-01-19

任意一个有穷数列都有通项公式。

给出 n 个点 (x1,y1),(x2,y2),,(xn,yn),怎么找到一个函数 f(x) 同时过这 n 个点?

拉格朗日的想法是这样的:构造一个开关函数

Lk(t)={1,t=xk0,t=other x

这个开关函数在学术上称为「插值基函数」。只有代入 xk 时开关才打开,代入其它 x 时开关关闭。

进而构造

f(x)=y1L1(x)+y2L2(x)++ynLn(x)

代入 xk 时,仅 yk 后面的开关是开着的,其余皆关闭,所以函数值为 yk,即 k,f(xk)=yk。这就是我们想要的函数。

问题就变成了如何构造开关函数。


容易发现 k,t=xk,

(txi)=(tx1)(tx2)(txn)=0

(txk) 这一项拿走,得到

ik(txi)={0,t=xk0,t=other x

为使 t=xk 时得到 1,再作以下改变

iktxixkxi={1,t=xk0,t=other x

于是我们成功地构造出了开关函数

Lk(t)=iktxixkxi

模板

cpp
#include <iostream>
#include <vector>

using namespace std;

double f(vector<double>& x, vector<double>& y, double t) {
    double result = 0.0;
    int n = x.size();
    
    for (int i = 0; i < n; i ++) {
        double term = y[i];
        for (int j = 0; j < n; j ++) {
            if (j != i) {
                term *= (t - x[j]) / (x[i] - x[j]);
            }
        }
        result += term;
    }
    
    return result;
}

int main() {
    vector<double> x = {1, 2, 3, 4, 5};
    vector<double> y = {1, 4, 9, 16, 25};

    double t;
    while (cin >> t) {
        cout << f(x, y, t) << endl;
    }

    return 0;
}