插值与拟合 (一) : 拉格朗日多项式插值 、Newton插值 、分段线性插值、Hermite插值 、样条插值、 B 样条函数插值、二维插值

插值:求过已知有限个数据点的近似函数。

拟合:已知有限个数据点,求近似函数,不要求过已知数据点,只要求在某种意义下它在这些点上的总偏差最小

插值和拟合都是要根据一组数据构造一个函数作为近似,由于近似的要求不同,二 者的数学方法上是完全不同的。而面对一个实际问题,究竟应该用插值还是拟合,有时容易确定,有时则并不明显。 

下面介绍几种基本的、常用的插值:拉格朗日多项式插值、牛顿插值、分段线性插 值、Hermite 插值和三次样条插值。

目录

1  拉格朗日多项式插值 

1.1  插值多项式 

范德蒙特(Vandermonde)行列式               截断误差 / 插值余项

1.2  拉格朗日插值多项式                                  1.3  用 Matlab 作 Lagrange 插值 

2  牛顿(Newton)插值 

 2.1 差商 : 定义与性质             2.2  Newton 插值公式 

Newton 插值的优点           差商与导数的关系         

2.3  差分 :向前差分、向后差分、中心差分              差分的两个性质

2.4  等距节点插值公式  、 Newton 向前插值公式

3  分段线性插值 

3.1  插值多项式的振荡              3.2  分段线性插值                         3.3  用 Matlab 实现分段线性插值 

4  埃尔米特(Hermite)插值 

4.1  Hermite 插值多项式                        4.2  用 Matlab 实现 Hermite 插值 

5  样条插值

5.1  样条函数的概念

    内节点 、边界点、k 次样条函数空间

二次样条函数                      三次样条函数

5.2  二次样条函数插值  

两类问题                                证明这两类插值问题都是唯一可解的

5.3  三次样条函数插值 

 3 种类型的边界条件:完备/Lagrange 、自然边界条件、周期条件 

5.4 三次样条插值在 Matlab 中的实现                    例 1  机床加工 

6   B 样条函数插值方法 

6.1  磨光函数                                          6.2  等距 B 样条函数 

6.3  一维等距 B 样条函数插值                   6.4  二维等距 B 样条函数插值 

7 二维插值 

7.1  插值节点为网格节点               7.2  插值节点为散乱节点                      习题


1  拉格朗日多项式插值 

1.1  插值多项式 

范德蒙特(Vandermonde)行列式

截断误差 / 插值余项

1.2  拉格朗日插值多项式 

1.3  用 Matlab 作 Lagrange 插值 

Matlab中没有现成的Lagrange插值函数,必须编写一个M文件实现Lagrange插值。 设n个节点数据以数组 x0 , y0  输入(注意 Matlat 的数组下标从 1 开始) ,m 个插值 点以数组 x输入,输出数组 y 为m 个插值。编写一个名为 lagrange.m 的 M 文件:

function y=lagrange(x0,y0,x); 
n=length(x0);m=length(x); 
for i=1:m    
    z=x(i);    
    s=0.0;    
    for k=1:n       
        p=1.0;       
        for j=1:n          
            if j~=k             
                p=p*(z-x0(j))/(x0(k)-x0(j));          
            end       
        end       
    s=p*y0(k)+s;    
    end    
y(i)=s; 
end 
 

2  牛顿(Newton)插值 

在导出 Newton 公式前,先介绍公式表示中所需要用到的差商、差分的概念及性质。

 2.1 差商 : 定义与性质

2.2  Newton 插值公式 

Newton 插值的优点

差商与导数的关系 

2.3  差分 :向前差分、向后差分、中心差分

节点等距时,即相邻两个节点之差(称为步长)为常数,Newton 插值公式的形 式会更简单。此时关于节点间函数的平均变化率(差商)可用函数值之差(差分)来表 示。

差分的两个性质

(i)各阶差分均可表成函数值的线性组合,例如 

(ii)各种差分之间可以互化。向后差分与中心差分化成向前差分的公式如下: 

2.4  等距节点插值公式  、 Newton 向前插值公式

3  分段线性插值 

3.1  插值多项式的振荡 

高次插值多项式的这些缺陷,促使人们转而寻求简单的低次多项式插值。 

3.2  分段线性插值 

 

用  I_{n}\left ( x \right ) 计算 x点的插值时,只用到 x左右的两个节点,计算量与节点个数n无关。 但n越大,分段越多,插值误差越小。实际上用函数表作插值计算时,分段线性插值就足够了,如数学、物理中用的特殊函数表,数理统计中用的概率分布表等。 

3.3  用 Matlab 实现分段线性插值 

用 Matlab 实现分段线性插值不需要编制函数程序,Matlab 中有现成的一维插值函 数 interp1

y=interp1(x0,y0,x,'method') 

method 指定插值的方法,默认为线性插值。其值可为:

  • 'nearest'   最近项插值

  • 'linear'    线性插值

  • 'spline'    逐段 3 次样条插值

  • 'cubic'    保凹凸性 3 次插值

 所有的插值方法要求 x0 是单调的。 当 x0 为等距时可以用快速插值法,使用快速插值法的格式为'*nearest'、'*linear'、 '*spline'、'*cubic'。

4  埃尔米特(Hermite)插值 

4.1  Hermite 插值多项式 

如果对插值函数,不仅要求它在节点处与函数同值,而且要求它与函数有相同的一 阶、二阶甚至更高阶的导数值,这就是 Hermite 插值问题。本节主要讨论在节点处插值 函数与函数的值及一阶导数值均相等的 Hermite 插值。 

4.2  用 Matlab 实现 Hermite 插值 

Matlab 中没有现成的 Hermite 插值函数,必须编写一个 M 文件实现插值。 

function y=hermite(x0,y0,y1,x); 
n=length(x0);m=length(x); 
for k=1:m    
    yy=0.0;    
    for i=1:n       
        h=1.0;       
        a=0.0;       
        for j=1:n          
            if j~=i             
                h=h*((x(k)-x0(j))/(x0(i)-x0(j)))^2;             
                a=1/(x0(i)-x0(j))+a;          
            end       
        end       
        yy=yy+h*((x0(i)-x(k))*(2*a*y0(i)-y1(i))+y0(i));    
    end    
    y(k)=yy; 
end 

 

5  样条插值

许多工程技术中提出的计算问题对插值函数的光滑性有较高要求,如飞机的机翼外 形,内燃机的进、排气门的凸轮曲线,都要求曲线具有较高的光滑程度,不仅要连续, 而且要有连续的曲率,这就导致了样条插值的产生。

5.1  样条函数的概念

 

所谓样条(Spline)本来是工程设计中使用的一种绘图工具,它是富有弹性的细木 条或细金属条。绘图员利用它把一些已知点连接成一条光滑曲线(称为样条曲线),并使连接点处有连续的曲率。 

    内节点 、边界点、k 次样条函数空间

二次样条函数

三次样条函数

利用样条函数进行插值,即取插值函数为样条函数,称为样条插值。例如分段线性插值 是一次样条插值。下面我们介绍二次、三次样条插值。  

5.2  二次样条函数插值  

两类问题

证明这两类插值问题都是唯一可解的

5.3  三次样条函数插值 

 3 种类型的边界条件:完备/Lagrange 、自然边界条件、周期条件 

5.4 三次样条插值在 Matlab 中的实现 

在 Matlab 中数据点称之为断点。如果三次样条插值没有边界条件,最常用的方法, 就是采用非扭结(not-a-knot)条件。这个条件强迫第 1 个和第 2 个三次多项式的三阶 导数相等。对最后一个和倒数第 2 个三次多项式也做同样地处理。

Matlab 中三次样条插值也有现成的函数:

y=interp1(x0,y0,x,'spline'); 

y=spline(x0,y0,x); 

pp=csape(x0,y0,conds),y=ppval(pp,x)

 其中 x0,y0 是已知数据点,x 是插值点,y 是插值点的函数值。 对于三次样条插值,我们提倡使用函数 csape,csape 的返回值是 pp 形式,要求出插值点的函数值,必须调用函数 ppval

pp=csape(x0,y0):使用默认的边界条件,即 Lagrange 边界条件。

pp=csape(x0,y0,conds)中的 conds 指定插值的边界条件,其值可为:

  • 'complete'    边界为一阶导数,即默认的边界条件
  • 'not-a-knot'   非扭结条件  
  • 'periodic'     周期条件
  • 'second'      边界为二阶导数,二阶导数的值[0, 0]。
  • 'variational'   设置边界的二阶导数值为[0,0]。

对于一些特殊的边界条件,可以通过 conds 的一个 1× 2 矩阵来表示,conds 元素的 取值为 1,2。此时,使用命令

pp=csape(x0,y0_ext,conds) 

其中 y0_ext=[left, y0, right],这里 left 表示左边界的取值,right 表示右边界的取值。

conds(i)=j 的含义是给定端点i的 j 阶导数,即 conds 的第一个元素表示左边界的条 件,第二个元素表示右边界的条件;

conds=[2,1]表示左边界是二阶导数,右边界是一阶 导数,对应的值由 left 和 right 给出。

详细情况请使用帮助 help csape。 

例 1  机床加工 

解  编写以下程序: 

clc,clear 
x0=[0   3   5   7   9   11   12   13   14  15]; 
y0=[0  1.2  1.7  2.0  2.1  2.0  1.8  1.2   1.0  1.6]; 
x=0:0.1:15; 
y1=lagrange(x0,y0,x);  %调用前面编写的Lagrange插值函数 
y2=interp1(x0,y0,x); 
y3=interp1(x0,y0,x,'spline'); 
pp1=csape(x0,y0); 
y4=ppval(pp1,x); 
pp2=csape(x0,y0,'second'); 
y5=ppval(pp2,x); 
fprintf('比较一下不同插值方法和边界条件的结果:\n') 
fprintf('x     y1      y2      y3      y4     y5\n') 
xianshi=[x',y1',y2',y3',y4',y5']; 
fprintf('%f\t%f\t%f\t%f\t%f\t%f\n',xianshi') 
subplot(2,2,1), plot(x0,y0,'+',x,y1), title('Lagrange') 
subplot(2,2,2), plot(x0,y0,'+',x,y2), title('Piecewise linear') 
subplot(2,2,3), plot(x0,y0,'+',x,y3), title('Spline1') 
subplot(2,2,4), plot(x0,y0,'+',x,y4), title('Spline2') 
dyx0=ppval(fnder(pp1),x0(1))  %求x=0处的导数 
ytemp=y3(131:151); 
index=find(ytemp==min(ytemp)); 
xymin=[x(130+index),ytemp(index)] 

计算结果略。 可以看出,拉格朗日插值的结果根本不能应用,分段线性插值的光滑性较差(特别 是在x =14 附近弯曲处),建议选用三次样条插值的结果。 

6   B 样条函数插值方法 

6.1  磨光函数 

实际中的许多问题,往往是既要求近似函数(曲线或曲面)有足够的光滑性,又要 求与实际函数有相同的凹凸性,一般插值函数和样条函数都不具有这种性质。如果对于 一个特殊函数进行磨光处理生成磨光函数(多项式),则用磨光函数构造出样条函数作 为插值函数,既有足够的光滑性,而且也具有较好的保凹凸性,因此磨光函数在一维插 值(曲线)和二维插值(曲面)问题中有着广泛的应用。 由积分理论可知,对于可积函数通过积分会提高函数的光滑度,因此,我们可以利 用积分方法对函数进行磨光处理。 

6.2  等距 B 样条函数 

6.3  一维等距 B 样条函数插值 

等距 B 样条函数与通常的样条有如下的关系: 

6.4  二维等距 B 样条函数插值 

7 二维插值 

前面讲述的都是一维插值,即节点为一维变量,插值函数是一元函数(曲线)。若 节点是二维的,插值函数就是二元函数,即曲面。如在某区域测量了若干点(节点)的 高程(节点值),为了画出较精确的等高线图,就要先插入更多的点(插值点),计算这些点的高程(插值)。 

7.1  插值节点为网格节点 

Matlab 中有一些计算二维插值的程序。如    

z=interp2(x0,y0,z0,x,y,'method') 

如果是三次样条插值,可以使用命令

pp=csape({x0,y0},z0,conds,valconds),z=fnval(pp,{x,y}) 

解  编写程序如下: 

clear,clc 
x=100:100:500; 
y=100:100:400; 
z=[636    697    624    478   450      
   698    712    630    478   420 
   680    674    598    412   400    
   662    626    552    334   310]; 
pp=csape({x,y},z') 
xi=100:10:500; yi=100:10:400 
cz1=fnval(pp,{xi,yi}) 
cz2=interp2(x,y,z,xi,yi','spline') 
[i,j]=find(cz1==max(max(cz1))) 
x=xi(i),y=yi(j),zmax=cz1(i,j) 

 

7.2  插值节点为散乱节点 

对上述问题,Matlab 中提供了插值函数 griddata,其格式为:

ZI = GRIDDATA(X,Y,Z,XI,YI) 

例 3  在某海域测得一些点(x,y)处的水深 z 由下表给出,在矩形区域(75,200) ×(-50,150) 内画出海底曲面的图形。 

解  编写程序如下: 

x=[129  140  103.5  88  185.5  195  105  157.5  107.5  77  81  162  162  117.5]; 
y=[7.5  141.5  23   147  22.5  137.5  85.5  -6.5  -81   3  56.5  -66.5  84 -33.5]; 
z=-[4     8    6     8    6     8     8     9     9   8    8    9    4    9]; 
xi=75:1:200; 
yi=-50:1:150; 
zi=griddata(x,y,z,xi,yi','cubic') 
subplot(1,2,1), plot(x,y,'*') 
subplot(1,2,2), mesh(xi,yi,zi) 

习题

 

  • 45
    点赞
  • 311
    收藏
    觉得还不错? 一键收藏
  • 4
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 4
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值