插值与拟合matlab实现

插值与拟合matlab实现

插值与拟合的Matlab实现

王正盛编写

在科技工程中，除了要进行一定的理论分析外，通过实验、观测数据，做分析、处理也是必不可少的一种途径。由于实验测定实际系统的数据具有一定的代表性，因此在处理时必须充分利用这些信息；又由于测定过程中不可避免会产生误差，故在分析经验公式时又必须考虑这些误差的影响。两者相互制约。据此合理建立实际系统数学模型的方法成为数值逼近法。

一、插值法

1、数学原理

工程实践和科学实验中，常常需要从一组实验观测数据中，求自变量与因变量的一个近似的函数关系式。

例如：观测行星的运动，只能得到某时刻所对应的行星位置（用经纬度表示），想知道任何时刻的行星位置。

例如：大气压测定问题；导弹发射问题；程序控制铣床加工

精密工件问题；飞机船舶制造问题等等。都属于此类问题。因为考虑到代数多项式既简单又便于计算，所以人们就用代数多项式近似地表示满足个点的函数关系式——插值法建模。

（1）计算方法课程中学习了两种多项式插值：Lagrange插值和Newton均差插值：

已知n+1个数据点：

n次Lagrange插值公式：

特别地，当n=1时，————线性插值

当n=2时，

———————抛物线插值或二次插值

Newton均差插值公式：，其中是k阶均差，可由均差表方便计算得到。

Lagrange插值和Newton均差插值本质上是一样的，只是形式不同而已，因为插值多项式是唯一的。

（2）Runge现象和分段低次插值：

如在[-5，5]上各阶导数存在，但在此区间取n个节点构造的Lagrange插值多项式在区间并非都收敛，而且分散得很厉害。（matlab\bin\ Lagrange.m是自己编写的M文件）[例]取n=-10

hold off

x=[-5:1:5];

y=1./(1+x.^2);

x0=[-5:0.1:5];

y0=lagrange(x,y,x0);

y1=1./(1+x0.^2);

plot(x0,y0)

hold on

plot(x0,y1,'b:')

legend('插值曲线','原数据曲线')

因此插值多项式一般不要超过四次为宜。为避免Runge现象提出分段插值，所谓分段插值就是首先把插值点分开，在每一段上用低次插值，再连接起来。

如分段线性插值、分段二次插值、分段三次插值等等。（3）三次样条插值：

三次样条插值算法的插值精度较高，所构造的曲线比较光滑（即在节点处连续可导）。因此，在许多工程设计或制造业中，例如飞机、导弹、高速机车、万吨级轮船等外

形设计常常利用本算法进行插值计算。

三次样条插值函数的定义：

设已知n+1个数据，如果函数满足如下性质：

（1）在每个子区间上，是次数不超过三次的多项式；（2）

（3）在插值区间上二次连续可微，记为

则称为三次样条插值函数，简称三次样条。

三次样条函数在每个子区间上可由4个系数唯一确定，因此，在上有4n个待定系数。由于，

故有

这里给出了3n-3个条件，加上插值条件（2），共有4n-2个条件。为了确定，通常还需要补充两个边界条件。

常用的边界条件有三类：

第一种边界条件：给定两边界节点处的一阶导数，并要求满足：

第二种边界条件：给定两边界节点处的二阶导数，并要求满足：

特别地，若，则所得的样条称为自然样条（MATLAB中即是自然样条）。

第三种边界条件：被插函数是以为周期的周期函数，要求满足：

定理：三次样条插值问题的解存在且唯一。（证明略）

求三次样条函数有多种方法，这里介绍一种著名的三弯矩法。

1、表达式：设三次样条函数为，且，在上，由定义中的（2）及处的二阶导数为，则有

(1)

(2)

由于是三次多项式，所以是线性函数，于是线性插值，得

令，得

（3）

对（3）式积分两次，得

（4）

其中为积分常数，利用（1），有

解得

代入（4）并化简，得

——————————————————（5）

这就是系数用二阶导数表示的三次样条的表达式。因此，只要确定，即得三次样条函数。

2、结点的关系式

为了求，利用的一阶导数的连续性：；

将（5）式对求导数，得

类似地，可得由，得

上式两边同乘以，得

令

则得到方程组

（6）

称为结点的M关系式（力学上称为三弯矩方程组）

3、边界条件

（6）式中是n+1个未知数的n-1个方程的方程组，不能唯一的确定，要唯一确定，必须附加两个条件，可由实际情况在两端点处给出的条件，称为边界条件（端点条件）。常见两种边界条件：

（1）给定两端点的二阶导数，

（7）

这时由结点的M关系式（6）可求出，特别地当时为自然样条。

（2）给出两端点的一阶导数（斜率）：，此时由一阶导数的连续性，可得到两个方程：

（8）

边界条件（7）和（8）可以写成统一形式

（9）

其中

当时，即为（7）式；当时，即为（8）式。

将（6）式与（9）式合在一起得到以下三对角方程组：

此方程组对角占优，用追赶法一定可解出。

4、解题步骤

第一步：确定边界条件；

第二步：建立结点的M关系式，求出各点的；

第三步：将所得的代回样条函数表达式（5），即得在个子区间上的表达式。

例题：某飞机头部外形曲线数据如下：x