【数值计算方法】曲线拟合与插值：Lagrange插值、Newton插值及其python/C实现

2024-07-30 09:23:47 浏览数 (2)

一、近似表达方式

插值、拟合和投影都是常用的近似表达方式，用于对数据或函数进行估计、预测或表示。

  • 插值（Interpolation）
    • 指通过已知数据点之间的插值方法，来估计或推算出在这些数据点之间的数值。插值可以用于构建平滑的曲线或曲面，以便在数据点之间进行预测或补充缺失的数据。
  • 拟合（Fitting）
    • 指通过选择合适的函数形式和参数，将一个数学模型与已知数据点拟合得最好的过程。拟合的目标是找到一个函数，使其在数据点附近的值与实际观测值尽可能接近。拟合可以用于数据分析、曲线拟合、回归分析等领域。
  • 投影（Projection）
    • 指将一个向量或一组向量映射到另一个向量空间或子空间上的过程。在线性代数中，投影可以用来找到一个向量在另一个向量或向量空间上的投影或投影分量。投影可以用于降维、数据压缩、特征提取等领域，以及计算机图形学中的投影变换。

二、插值

Lagrange插值和Newton插值都是常见的多项式插值方法，用于通过给定的一组数据点来估计在其他点上的函数值。它们之间的主要区别在于插值多项式的构建方法。

  • Lagrange插值使用基于Lagrange多项式的方法来构建插值多项式。
    • Lagrange多项式是通过将每个数据点与一个基函数相乘，并使得在其他数据点上该基函数为零来构造的。最终的插值多项式是将所有这些基函数相加得到的。
    • Lagrange插值的优点是易于理解和实现，但在数据点较多时可能会导致计算复杂度较高的问题。
  • Newton插值使用差商的概念来构建插值多项式。
    • 差商是一个递归定义的概念，用于计算插值多项式中的系数。差商的计算可以通过表格形式进行，其中每一列都表示不同阶数的差商。通过计算差商，可以逐步构建插值多项式。
    • Newton插值的优点是在计算差商时可以重复使用已计算的差商值，从而减少计算量。

1. Lagrange插值

Lagrange插值是一种用于通过已知数据点构造一个多项式函数的方法。它是基于拉格朗日插值多项式的原理，该多项式通过每个数据点并满足相应的条件。拉格朗日插值可用于估计数据点之间的值，而不仅仅是在给定数据点上进行插值。

使用Lagrange插值的基本步骤如下：

  • 给定一组已知的数据点，包括横坐标和纵坐标的值。
  • 根据数据点的数量，构造相应次数的拉格朗日插值多项式。
  • 将每个数据点的函数值乘以对应的拉格朗日插值多项式，并将它们相加，得到最终的插值函数。

通过这种方法，可以在给定的数据点上获得一个平滑的插值函数，使得在这些数据点之间的任何位置上都可以估计函数的值。Lagrange插值在数据点较少或数据点之间存在较大间隔时可能会出现一些问题，例如插值多项式可能会产生振荡现象，这被称为Runge现象。

Lagrange插值公式
线性插值(n=1)
抛物插值(n=2)

范德蒙行列式 - 知乎 (zhihu.com)

icon-default.png?t=N7T8icon-default.png?t=N7T8

https://zhuanlan.zhihu.com/p/161300510

python实现
代码语言：javascript复制
import numpy as np


# 定义Lagrange插值函数
def lagrange_interpolation(x, y, xi):
    n = len(x)
    yi = 0.0

    for i in range(n):
        # 计算拉格朗日插值多项式的每一项
        term = y[i]
        for j in range(n):
            if j != i:
                term *= (xi - x[j]) / (x[i] - x[j])
        yi  = term

    return yi


# 示例数据点
x = np.array([0.32, 0.34, 0.36])
y = np.array([0.314567, 0.333487, 0.352274])

# 要进行插值的点
xi = 0.3367

# 进行插值
yi = lagrange_interpolation(x, y, xi)

print("插值结果：", yi)
print("真实结果：", np.sin(xi))

输出：

代码语言：javascript复制
插值结果： 0.3303743620374999
真实结果： 0.330374191555628
C语言实现
代码语言：javascript复制
#include <stdio.h>

// 计算Lagrange插值多项式的值
double lagrange_interpolation(double x[], double y[], int n, double xi) {
    double yi = 0.0;

    for (int i = 0; i < n; i  ) {
        double term = y[i];
        for (int j = 0; j < n; j  ) {
            if (j != i) {
                term *= (xi - x[j]) / (x[i] - x[j]);
            }
        }
        yi  = term;
    }

    return yi;
}

int main() {
    // 示例数据点

    double x[] = {0.32, 0.34, 0.36};
    double y[] = {0.314567, 0.333487, 0.352274};

    // 要进行插值的点
    double xi = 0.3367;

    // 数据点的个数
    int n = sizeof(x) / sizeof(x[0]);

    // 进行插值
    double yi = lagrange_interpolation(x, y, n, xi);

    printf("插值结果：%fn", yi);

    return 0;
}

输出：

代码语言：javascript复制
插值结果：0.330374

2. Newton插值

Newton插值基于差商的概念。通过给定的一组数据点，Newton插值可以生成一个通过这些点的多项式，从而在给定的数据范围内进行插值和外推。

Newton插值的基本思想是使用差商来递归地构建一个多项式。差商是通过递归地计算数据点之间的差分来定义的。具体而言，对于给定的数据点 (x0, y0), (x1, y1), ..., (xn, yn)，差商可以表示为：

f[x_{0}] = y_{0}f[x_{0}] = y_{0}
f[x_{1}, x_{0}] =frac{ (f[x_{1}] - f[x_{0}]) }{ (x_{1} - x_{0})}f[x_{1}, x_{0}] =frac{ (f[x_{1}] - f[x_{0}]) }{ (x_{1} - x_{0})}
f[x_{2}, x_{1}, x_{0}] =frac{ (f[x_{2}, x_{1}] - f[x_{1}, x_{0}]) }{ (x_{2} - x_{0})}f[x_{2}, x_{1}, x_{0}] =frac{ (f[x_{2}, x_{1}] - f[x_{1}, x_{0}]) }{ (x_{2} - x_{0})}

………… f[xn, xn-1, ..., x0] = (f[xn, xn-1, ..., x1] - f[xn-1, ..., x0]) / (xn - x0)

然后，通过将这些差分商逐步添加到多项式中，可以得到一个多项式，表示为：

P(x) = f[x_{0}] f[x_{1}, x_{0}](x - x_{0}) f[x_{2}, x_{1}, x_{0}](x - x_{0})(x - _{x1}) ...P(x) = f[x_{0}] f[x_{1}, x_{0}](x - x_{0}) f[x_{2}, x_{1}, x_{0}](x - x_{0})(x - _{x1}) ...

这个多项式可以用于在给定的数据范围内进行插值，即通过已知的数据点来估计其他点的函数值。Newton插值的优点之一是它可以通过添加更多的数据点来逐步改进插值结果。然而，它也存在一些问题，比如所得到的多项式可能会出现龙格现象（Runge's phenomenon），导致在边界处产生振荡。

python实现
代码语言：javascript复制
def newton_interpolation(x, y, xi):
    # 计算差分商
    n = len(x)
    f = [[0] * n for _ in range(n)]
    for i in range(n):
        f[i][0] = y[i]

    for j in range(1, n):
        for i in range(n - j):
            f[i][j] = (f[i   1][j - 1] - f[i][j - 1]) / (x[i   j] - x[i])

    # 构建插值多项式
    result = f[0][0]
    for j in range(1, n):
        term = f[0][j]
        for i in range(j):
            term *= (xi - x[i])
        result  = term

    return result


# 示例数据
x = [0.32, 0.34, 0.36]
y = [0.314567, 0.333487, 0.352274]
xi = 0.3367

# 进行插值
interpolated_value = newton_interpolation(x, y, xi)
print("插值结果:", interpolated_value)

输出：

代码语言：javascript复制
插值结果: 0.3303743620375
C语言实现
代码语言：javascript复制
#include <stdio.h>

double newton_interpolation(double x[], double y[], int n, double xi) {
    // 计算差分商
    double f[n][n];
    for (int i = 0; i < n; i  ) {
        f[i][0] = y[i];
    }

    for (int j = 1; j < n; j  ) {
        for (int i = 0; i < n - j; i  ) {
            f[i][j] = (f[i 1][j-1] - f[i][j-1]) / (x[i j] - x[i]);
        }
    }

    // 构建插值多项式
    double result = f[0][0];
    for (int j = 1; j < n; j  ) {
        double term = f[0][j];
        for (int i = 0; i < j; i  ) {
            term *= (xi - x[i]);
        }
        result  = term;
    }

    return result;
}

int main() {
    // 示例数据
    double x[] = {0.32, 0.34, 0.36};
    double y[] = {0.314567, 0.333487, 0.352274};
    int n = sizeof(x) / sizeof(x[0]);
    double xi = 0.3367;

    // 进行插值
    double interpolated_value = newton_interpolation(x, y, n, xi);
    printf("插值结果: %fn", interpolated_value);

    return 0;
}

输出：

代码语言：javascript复制
插值结果: 0.330374
double 递归 函数 数据 python

0 人点赞