函数逼近——(Lagrange)拉格朗日插值法 | 北太天元 or Matlab

函数逼近——(Lagrange)拉格朗日插值法 | 北太天元 or Matlab

引用

CSDN

1.

https://m.blog.csdn.net/Math_Boy_W/article/details/137800068

拉格朗日插值法是数值分析中的重要方法，广泛应用于数学、物理、工程等领域。本文将详细介绍拉格朗日插值法的理论基础、算法实现及其在北太天元或Matlab中的具体应用。

一、Lagrange插值法

对于 (n+1) 个样本点 ((x_i,y_i),i=0,1,2,\dots ,n)，拉格朗日插值多项式 (L_n(x)) 可以表示为：

[L_n(x)=\sum_{i=0}^{n}y_iI_i(x)]

其中，(I_i(x)) 是拉格朗日插值基函数，定义为：

[I_i(x)=\frac{(x-x_0)(x-x_1)\dots(x-x_{i-1})(x-x_{i+1})\dots(x-x_n)}{(x_i-x_0)(x_i-x_1)\dots(x_i-x_{i-1})(x_i-x_{i+1})\dots(x_i-x_n)}]

也可以表示为：

[I_i(x)=\frac{\omega_n(x)}{(x-x_i)\omega_n'(x_i)}]

其中，(\omega_n(x)=\prod_{j=0}^{n}(x-x_j))。

二、算法实现步骤

1. 输入已知插值点的横坐标向量 (x0) 和纵坐标向量 (y0)。
2. 计算每个 (x) 对应的拉格朗日基函数 (I_i(x))。
3. 结合 (I_i(x)) 与 (y0) 计算每个 (x) 对应的 (L_n(x))。
4. 输出对应的 (L_n(x))。

三、北太天元 or Matlab 实现

以下是拉格朗日插值法的Matlab实现代码：

function  [y] = Lag_interp_v1(x0,y0,x)
    n1 = length(x0);   % n1表示样本点的个数 
    I = zeros(1,n1);   % 预留出要用的空间
    n2 = length(x);    % n2表示所求点的个数
    Ln = zeros(1,n2);
    for j=1:1:n2     %  依次代入 自变量 x(j) 
        omega_x = x(j)-x0;    % 数 - 矩阵 ，表示 [x(j) - x0(1),x(j)-x0(2),...,x(j)-x0(n1)]
        for i = 1:1:n1   %  对于x(j)求对应的基函数 I（i）
            w = x0(i) - x0;  % 同样是 数 - 矩阵
            % 这里使用 if 和 内置的 prod 函数代替了 for 循环
            % prod 表示矩阵内所有元素的乘积
            if i == 1    
                % omega_x(i+1:n1)表示向量的节选，第i+1个到第n1个元素
                I(i) = prod(omega_x(i+1:n1))/prod(w(i+1:n1));
            elseif i == n1
                I(i) = prod(omega_x(1:i-1))/prod(w(1:i-1));
            else
                I(i) = prod(omega_x(1:i-1))/prod(w(1:i-1));
                I(i) = I(i) * prod(omega_x(i+1:n1))/prod(w(i+1:n1));
            end 
        end
        %使用矩阵的乘积，行向量 × 列向量 得到一个值
        Ln(j) = y0*I'; 
    end    
    y = Ln;
end

四、数值算例

利用 (f(x) = \ln x) 的如下数据：

进行Lagrange插值:

1. 计算 (x = [0.412, 0.511, 0.666]) 处的近似值
2. 计算 (x_i = 0.3 + ih, h=0.01, i = 0, 1, 2, \dots, 50) 处的近似值，并作图

调用函数 [y] = Lag_interp_v1(x0,y0,x)，相应的实现代码为：

clc, clear all, format long;
x0 = linspace(0.4, 0.7, 4); % 输入样本点的横坐标
y0 = [-0.916291, -0.693147, -0.510826, -0.35765]; % 输入样本点的纵坐标
% 简单的算几个点
x1 = [0.412, 0.511, 0.666];
y1 = Lag_interp_v1(x0, y0, x1);
for i = 1:length(x1)
    fprintf('f(%f) = %f \n', x1(i), y1(i));
end
% 利用很多的点来画图
x2 = linspace(0.3, 0.8, 51); % 共51个要求的点
% 利用写好的 Lag_interp_v1 函数计算要求点的纵坐标
y2 = Lag_interp_v1(x0, y0, x2);
delta = abs(y2 - log(x2));
%作图
figure(1); %画出第一个图像
plot(x2, y2, 'b');

figure(2);
plot(x2, y2, 'b');
hold on
plot(x2, log(x2), 'r'); % y = lnx 的图像
hold off
% Ln(x) 与 lnx 的 误差
figure(3); % 画出第二个图像
plot(x2, delta, 'g');
% 文字形式表示出来
for j = 1:51
    fprintf('f(%f) = %f \n', x2(j), y2(j));
end

运行后得到的结果：

f(0.412000) = -0.886972 
f(0.511000) = -0.671305 
f(0.666000) = -0.407185 
f(0.300000) = -1.191936 
f(0.310000) = -1.161676 
f(0.320000) = -1.132046 
f(0.330000) = -1.103035 
f(0.340000) = -1.074630 

此处省略

f(0.770000) = -0.261515 
f(0.780000) = -0.248246 
f(0.790000) = -0.235059 
f(0.800000) = -0.221941 

插值结果可视化

下面是 (L_n(x)) 的图像：

下面是 (L_n(x)) 与 (\ln(x)) 的图像对比：

下面是 (L_n(x)) 与 (\ln(x)) 的误差图像：

本文详细介绍了拉格朗日插值法的理论基础、算法实现及其在北太天元或Matlab中的具体应用。通过理论与实践的结合，帮助读者更好地理解拉格朗日插值法的应用场景和实现方法。

本文原文来自CSDN