目录

线性插值 

原理

流程图

 代码

 抛物插值

原理

流程图

 代码

 拉格朗日插值

代码

牛顿插值

原理

代码

分段线性插值

代码

线性插值 

原理

流程图

 单个点的线性插值代码

X=[0.2 0.4];
Y=[21 25];
x=0.7;
x0=X(1)
y0=Y(1);
x1=X(2);
y1=Y(2);


L0=(x-x1)/(x0-x1);
L1=(x-x0)/(x1-x0);
y=y0*L0+y1*L1;

多个点的线性插值代码

time = [1 3 8 12 15 20 24];
tem = [8 9 16 23 22 18 10];
time_i = 1:0.01:24;
tem_i = self_interp1(time,tem,time_i,'linear');
plot(time,tem,'o',time_i,tem_i);

%自己写一个interp1类似功能的接口
%在这个接口中参数x需要从大到小排列，y随意
function [yi]=self_interp1(x,y,xi,method)
% 初始化yi，给它xi对应的列
col_xi = size(xi,2);
yi = zeros(1,col_xi);
% 检测使用的插值方法 这里期望的是'linear'
if strcmp(method,'linear')
    % 找到每个xi在x序列中的位置
    col_x = size(x,2);
    for i = 1:col_xi,
        for j = 1:col_x-1, 
            % 假如需要计算插值公式
            if x(j+1) > xi(i), 
                yi(i) = y(j)+(y(j+1)-y(j))/(x(j+1)-x(j))*(xi(i)-x(j));
                break;
            end
            % 假如插值处的数据已经测得了，就直接把值给它，节约计算资源
            if x(j) == xi(i),
                yi(i) = y(j);
                break;
            end
        end
        % 以上没有把最后一个数据点考虑进去,需要加上
        yi(col_xi) = y(col_x);
    end
else
    error('插值方法请选择(linear)\n');
end
end

 抛物插值

原理

流程图

 代码

X=[0.2 0.4 0.6];
Y=[21 25 23];
x0=X(1);
x1=X(2);
x2=X(3);
y0=Y(1);
y1=Y(2);
y2=Y(3);
x=0.7;

L0=(x-x1)*(x-x2)/(x0-x1)/(x0-x2);
L1=(x-x0)*(x-x2)/(x1-x0)/(x1-x2);
L2=(x-x0)*(x-x1)/(x2-x0)/(x2-x1);

y=y0*L0+y1*L1+y2*L2

 拉格朗日插值

代码

x=[0.2 0.4 0.6 0.8];
y=[21 25 23 20];
yh=lagrange(x,y,0.7)

function yh=lagrange (x,y,xh)
    n = length(x);
    m = length(xh);
    yh = zeros(1,m); 
    c1 = ones(n-1,1);
    c2 = ones(1,m);
    for i=1:n
      xp = x([1:i-1 i+1:n]);
      yh = yh + y(i)*prod((c1*xh-xp'*c2)./(x(i)-xp'*c2));
    end
end

牛顿插值

原理

代码

xi=[1 4 9];
yi=[1 2 3];
x=7;
p= Newton_fun(x,xi,yi)


function p= Newton_fun(x,xi,yi)
n=length(xi);
f=zeros(n,n);

% 对差商表第一列赋值
for k=1:n      
    f(k)=yi(k);
end
% 求差商表
for i=2:n       % 差商表从0阶开始；但是矩阵是从1维开始存储！！！！！！
    for k=i:n
        f(k,i)=(f(k,i-1)-f(k-1,i-1))/(xi(k)-xi(k+1-i));  
    end
end
disp('差商表如下：');
disp(f);

%求插值多项式
p=0;          
for k=2:n
    t=1;
    for j=1:k-1
        t=t*(x-xi(j));
        disp(t)
    end
    p=f(k,k)*t+p;
    disp(p)
end
p=f(1,1)+p;

end

分段线性插值

原理

代码

x = [1 3 8 12 15 20 24];
y = [8 9 16 23 22 18 10];

yy=fdxx(x,y,7)

function yy=fdxx(x,y,xx)
    n=size(x,2);
    for i=1:n-1
        if x(i)<xx&&xx<x(i+1)
            L1=(xx-x(i+1))/(x(i)-x(i+1));
            L2=(xx-x(i))/(x(i+1)-x(i));
            yy=L1*y(i)+L2*y(i+1);
            break;
        elseif x(i)==xx
             yy=y(i);      
        end
    end
    
end

分段抛物插值

原理

代码

x = [1 3 8 12 15 20 24];
y = [8 9 16 23 22 18 10];
y=fenduanpaowu(x,y,7)

function y=fenduanpaowu(xi,yi,x)
    n=size(xi,2);
    if x<xi(2)
        L1=(x-xi(1))*(x-xi(3))/(xi(1)-xi(2))/(xi(1)-xi(3));
        L2=(x-xi(1))*(x-xi(3))/(xi(2)-xi(1))/(xi(2)-xi(3));
        L3=(x-xi(1))*(x-xi(2))/(xi(3)-xi(1))/(xi(3)-xi(2));
        y=L1*yi(1)+L2*yi(2)+L3*yi(3);
    elseif x>xi(end-1)
        L1=(x-xi(end-1))*(x-xi(end))/(xi(end-2)-xi(end-1))/(xi(end-2)-xi(end));
        L2=(x-xi(end-2))*(x-xi(end))/(xi(end-1)-xi(end-2))/(xi(end-1)-xi(end));
        L3=(x-xi(end-2))*(x-xi(end-1))/(xi(end)-xi(end-2))/(xi(end)-xi(end-1));
        y=L1*yi(1)+L2*yi(2)+L3*yi(3);
   else
        for k=2:n-1
            if xi(k+1)>x  
                if abs(x-xi(k))<abs(x-xi(k+1))
                    i=k-1;
                    
                else
                    i=k;
                end
                L1=(x-xi(i+1))*(x-xi(i+2))/(xi(i)-xi(i+1))/(xi(i)-xi(i+2));
                L2=(x-xi(i))*(x-xi(i+2))/(xi(i+1)-xi(i))/(xi(i+1)-xi(i+2));
                L3=(x-xi(i))*(x-xi(i+1))/(xi(i+2)-xi(i))/(xi(i+2)-xi(i+1));
                y=L1*yi(i)+L2*yi(i+1)+L3*yi(i+2);
            end
        end
   end
    
end

拉格朗日插值即对所要插值的函数进行拉格朗日多项式拟合

这是matlab插值系列的第二期，第一期：[数值分析拟合]Matlab三次样条插值拟合数据

（以后会有时间的时候再更新更多的插值方法）

这篇文章我推导过程参考过了一些其他文章，代码是自己写的，如有不对或者公式打错了欢迎批评指正

首先，对于所需要插值的自变量x和所需插值的数据点y:

一、我们先来了解它的插值原理：

        对于在一组数值散点中的任意一点进行插值，找到一个满足相应条件的n次多项式，我们希望能用所有点的函数值去表示它，并且每一点的函数值都与原来的函数值相符合。

        因此，设原数据的每一个点的函数值为,为了组成插值所得到的,前面配凑的系数是，那么则有：

        记从0开始推导，则
        

        那么由两个点的插值开始，所得的应当是线性插值

则有  那么显然有

对于中间的每一个点，都应当有

整理得到

那么从0开始，则其插值的系数为

, 有

这时，规律不够明显， 因此我们采用更多的点来进行相关的推导

当点的个数为n+1时,为了保证依然成立

首先，为了保证下方条件，分子为

 （图片来自参考文章

拉格朗日（Lagrange）插值多项式的基函数构造法(详细推导)）

此时，根据,可以有n+1个方程，那么由于有n+1个c(c0,c1,c2....cn)，则可通过解方程唯一确定对应的c0,解出c0后，即可确定

,

通过以上的公式推导，下面使用MATLAB来实现上述的计算内容

先加上我在里面调用的len.m

%返回一维数组长度(行数或列数其中大的一个，必须是一位数组)
function [length]=len(A)

    s=size(A);
    if s(2)==1
        length=s(1); %取行
    end
    if s(1)==1
        length=s(2);
    end
    if s(1)~=1 && s(2)~=1
        disp('必须是单行向量或单列向量')
        return
    end
            
end

编写函数     Lagrange.m计算拉格朗日插值:

其中许多的注释我用英语写的看着别嫌麻烦哈，这段时间为了学英语我也在天天头皮发麻）   

% *** Lagrange  interpolation ***
function L = Lagrange(x,y,x_2)
% x --origional x vector
% y --origional y vector need to be interpolated
% x_2 --the vector x that you want to interpolate to 
% author:FriedParrot   --2022-2-4
if len(x) ~= len(y)
    error('The length of x and y should be correspond');
end

xi = x_2;  %  generate the vector that need to be interpolated
% define the initial vectors
L = zeros(1,len(xi));

for i = 1:1:len(xi)  % This is the index of the xi(to be interpolated)
    l = ones(1,len(x));  % it is used every time
    for k = 1:1:len(x)
        for j= 1:1:len(x)
            if j ~= k  % similar to prod but the denominator shouldn't be zero
                l(k) = l(k) * ( xi(i)-x(j)) / (x(k)-x(j));
                % the index should be k --the first index of the iteration
            end
        end % cacultte the l_j
        L(i) = L(i) + l(k)*y(k);  % just for every single loop of s
    end
end

if nargout == 0
    figure('name','Lagrange Interpolation');
    plot(xi,L);
end

end

下面加上测试代码以及效果：

（代码使用外插值，即公式中的范围比数据的范围更广）

x = [1,2,4,6,8,9];
y =  cos(x);
x_2 = 0:0.05:10;
Lagrange(x,y,x_2);

插值效果： 

参考文章：

拉格朗日（Lagrange）插值多项式的基函数构造法(详细推导)

mathworld里关于Lagrange插值的简介（用翻墙软件可以进）

《专题一：插值法（1）拉格朗日插值法》

目录

前言

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

二、matlab实现代码

1.线性插值：

2.抛物线插值：

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

总结：

前言

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

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

对于构造通过n+1个不同的节点的n次插值多项式,假定它满足条件

为了构造，我们先定义n次插值基函数。

若n次多项式在n+1个节点上满足条件

则称这n+1个n次多项式为节点上的n次插值基函数。

且

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

由的定义，知

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

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

其中

二、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=0;
for m=length(x)%循环作用：不同项作和
    X=1;
    for n=length(x)%循环作用：构成每一项的插值基函数的形式
        if x(m)~=x(n)%判断作用：排除分子为零的循环
           X=X.*(x0-x(n))./(x(m)-x(n));
        else
            continue;
        end
    end
    y0=y(m)*X./(x0-x(m))+y0;
end
end

总结

拉格朗日插值的插值对象为：一串不同的插值节点，且知道插值节点及其值。

代码部分由于线性插值和抛物线插值是拉格朗日的特殊情况，所以小编在编写的时候，为了让看起来没有重复，选择了直接按照运算形式编写代码。