数值分析·学习 | 拉格朗日插值法matlab实现

不变因子

已于 2022-10-20 16:09:11 修改

阅读量2w

点赞数 28

文章标签：学习 matlab 数据分析

于 2022-06-29 13:12:34 首次发布

本文链接：https://blog.csdn.net/qq_59235603/article/details/125478171

版权

前言

本篇内容为个人所学知识分享

一、拉格朗日（Lagrange）插值是什么？

对于构造通过n+1个不同的节点 $x_0<x_1<\cdot \cdot \cdot <x_n$ 的n次插值多项式 $L_n(x)$ ,假定它满足条件

$L_n(x_j)=y_j,j=0,1,\cdot \cdot \cdot ,n.$

为了构造 $L_n(x)$ ，我们先定义n次插值基函数。

若n次多项式 $l_j(x)(j=0,1,...,n)$ 在n+1个节点 $x_0<x_1<\cdot \cdot \cdot <x_n$ 上满足条件

$l_j(x_k)=\left\{\begin{matrix} 1,k=j, & & \\ 0,k\neq j, & & \end{matrix}\right.j,k=0,1,\cdot \cdot \cdot ,n,$

则称这n+1个n次多项式 $l_0(x),l_1(x),\cdot \cdot \cdot ,l_n(x)$ 为节点 $x_0,x_1,\cdot \cdot \cdot ,x_n$ 上的n次插值基函数。

且 $l_k(x)=\frac{(x-x_0)\cdot \cdot \cdot(x-x_{k-1})(x-x_{k+1})\cdot \cdot \cdot (x-x_n)}{(x_k-x_0)\cdot \cdot \cdot(x_k-x_{k-1})(x_k-x_{k+1})\cdot \cdot \cdot (x_k-x_n)},k=0,1,\cdot \cdot \cdot ,n.$

使得插值多项式可以表示为

$L_n(x)=\sum_{k=0}^{n}y_kl_k(x),$

由 $l_k(x)$ 的定义，知

$L_n(x_j)=\sum_{k=0}^{n}y_kl_k(x_j)=y_j,j=0,1,\cdot \cdot \cdot ,n.$

形如上式的插值多项式 $L_n(x)$ 便称为拉格朗日（Lagrange）插值多项式。线性插值和抛物线插值只是拉格朗日插值的特殊情况。

拉格朗日（Lagrange）插值多项式的另外一种表示形式如下：

$L_n(x)=\sum_{k=0}^{n}y_k\frac{w_{n+1}(x)}{(x-x_k)w'_{n+1}(x_k)}.$

其中

$w_{n+1}(x)=(x-x_0)(x-x_1)\cdot \cdot \cdot (x-x_n),$

$w'_{n+1}(x_k)=(x_k-x_0)\cdot \cdot \cdot (x_k-x_{k-1})(x_k-x_{k+1})\cdot \cdot \cdot (x_k-x_n).$

二、matlab实现代码

1.线性插值：

即n=1的时候,一次的插值函数，即已知条件为两个插值节点及其值

function y0=Linear_interpolation(x,y,x0)
%功能：线性插值
%输入：x为插值节点，y为插值节点对应的值，x0为计算点
%输出：x0处的值
y0=y(1)+(y(2)-y(1))/(x(2)-x(1))*(x0-x(1));%点斜式
y0=y(1)*(x(2)-x0)/(x(2)-x(1))+y(2)*(x0-x(1))/(x(2)-x(1));%两点式
end

2.抛物线插值：

即n=2的时候，二次的插值函数，即已知条件为三个插值节点及其值

function y0=Parabolic_interpolation(x,y,x0)
%功能：抛物线插值
%输入：x为插值节点，y为插值节点对应的值，x0为计算点
%输出：x0处的值
y0=y(1)*(x0-x(2))*(x-x(3))/(x(1)-x(2))/(x(1)-x(3))+y(2)*(x0-x(1))*(x0-x(3))...
    /(x(2)-x(1))/(x(2)-x(3))+y(3)*(x0-x(1))*(x0-x(2))/(x(3)-x(1))/(x(3)-x(2));
end

3.拉格朗日（Lagrange）插值：

（一般情况）即已知条件为一串插值节点及其值

function y0 = Lagrange_interpolation(x,y,x0)
%功能：拉格朗日插值多项式求解
%输入：x为插值节点，y为插值节点对应的值，x0为计算点集
%输出：x0处值的集合y0
y0 = zeros(1,length(x0));
for i = 1:1:length(x0)  %循环作用：遍历x0
   X = ones(1,length(x));  
    for k = 1:1:length(x)%循环作用：(x0(i)-x(j))/(x(k)-x(j))的加和;
        for j= 1:1:length(x)%循环作用：(x0(i)-x(j))/(x(k)-x(j))的乘积
            if j ~= k  
                X(k) = X(k) * (x0(i)-x(j)) / (x(k)-x(j));
            end
        end 
        y0(i) = y0(i) + X(k)*y(k);  
    end
end
end