三次样条插值多项式matlab
- 格式:doc
- 大小:60.00 KB
- 文档页数:6
三次样条插值函数MATLAB 编程实现三次样条插值函数为()()[)()[]1011,,,,n n n S x x x x S x S x x x x-⎧∈⎪=⎨⎪∈⎩ 利用三次埃尔米特插值函数表示三次样条插值函数,即()()()()())111111,,j j j j j j j j j j j S x y x y x m x m x x x x ααββ++++++⎡=+++∈⎣(0,1,,1j n =-)基函数满足()()()()()()21112111121121111212jj j j j j j j j j j j j j j j j j j jj j j j x x x x x x x x xx xx x x x x x xx xx x x x xx x x x x x xααββ++++++++++++⎛⎫⎛⎫--=+ ⎪⎪ ⎪⎪--⎝⎭⎝⎭⎛⎫⎛⎫--=+ ⎪⎪ ⎪⎪--⎝⎭⎝⎭⎛⎫-=-⎪ ⎪-⎝⎭⎛⎫-=-⎪ ⎪-⎝⎭由上式易得()()()()()()()()()()()()()()1331111331112211112211612612246246j j j j j j j j j j j j j j j j j j j j j j jj j j j j x x x x x x xx x x x x x x x x x x x x xx xx x x x x x x x x ααββ+++++++++++++++''=---+''=-+--+''=---+''=---则有()()()()()()()()()()()111111113333111111122221111661212242466j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j S x y x y x m x m x x x x x y x y x x x x x x x x x x x x x m x m x x x x x x x x x ααββ+++++++++++++++++++''''''''''=+++⎡⎤⎡⎤++⎢⎥⎢⎥=-+-+⎢⎥⎢⎥----⎣⎦⎣⎦⎡⎤++⎢⎥+-+-⎢⎥----⎣⎦)1,j j x x x +⎡⎤⎢⎥⎡∈⎣⎢⎥⎣⎦(0,1,,1j n =-)同理有()()()()()()()()()()()()()()()11111113333111111122221111661212242466j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j j S x y x y x m x m x x x x x y x y x x x x x x x x x x x x x m x m x x x x x x x x x ααββ------------------''''''''''=+++⎡⎤⎡⎤++⎢⎥⎢⎥=-+-+⎢⎥⎢⎥----⎣⎦⎣⎦⎡⎤⎡++⎢⎥+-+-⎢⎥----⎣⎦⎣)1,j j x x x -⎤⎢⎥⎡∈⎣⎢⎥⎦(1,,j n =)根据样条函数二阶导数连续性,即()()100j j j j S x S x +''''+=-(1,,1j n =-)即()()()()()()()()()()()()()()()()111111332211111111113322111166426624j j jj j j j j j j jj j j j j j j j j j j j j j jj j j j j jj j j j j j j j x x y x x y x x x x m m x x x x x x x x x x y x x y x x x x m m x x x x x x x x ++++++++++--------------+++--------=+++----(1,,1j n =-)化简得()()()()()111111111111233j j j j j j j j j j j j j j j j j j jj j xx m x x m x x m x x x x y y y y x x x x +-+--+-++-+--+-+---=-+---(1,,1j n =-)可得线性方程组()()()()()()()()()()0121201023231213121221111110212110211032213221322122233333n n n n n n n n n n n n m m x x x x x x m x x x x x x m x x x x x x m m m x x x x y y y y x x x x x x x x y y y y x x x x y ------⨯+-+⨯⎛⎫ ⎪ ⎪ ⎪---⎛⎫⎪ ⎪--- ⎪ ⎪⎪ ⎪⎪ ⎪ ⎪--- ⎪⎝⎭ ⎪ ⎪ ⎪⎝⎭---+------+---=()()()121112112113n n n n n n n n n n n n n x x x x y y y x x x x ----------⨯⎛⎫⎪ ⎪ ⎪ ⎪⎪⎪ ⎪-- ⎪-+- ⎪--⎝⎭为了使样条插值问题有惟一解,我们在原有方程基础上增加两个边界条件。
样条插值函数 matlab样条插值函数是一种广泛应用于数值计算中的插值技术。
它的主要功能是使用已经给定的数据点来构建一个连续且平滑的函数,该函数可以通过插值来预测在给定数据点之间的值。
在 MATLAB 中,我们可以使用 spline() 函数来实现样条插值函数的操作。
样条插值函数的主要思想是在给定的数据点上构建一组多项式,以实现对在数据点之间位置的函数的插值。
这些多项式通常在相邻数据点之间的区间上定义,并满足一系列的平滑条件。
在实践中,可以使用不同的插值多项式来构建样条插值函数,例如线性样条,二次样条或三次样条。
为了使用 MATLAB 中的 spline() 函数进行样条插值的操作,我们需要先确定要插值的数据点以及在这些数据点上的值。
通常我们会使用一维数组来存储这些数据点的值,例如:x = [0 0.5 1 1.5 2];接下来,我们需要为这些数据点上的值生成一个一维数组,例如:y = [0 0.4794 0.8415 0.9975 0.9093];为了使用 spline() 函数进行样条插值,我们需要将这些数据点作为输入参数传递给该函数,例如:yy = spline(x,y);这将返回一个样条插值函数,可以使用该函数来计算在数据点之间位置的函数值。
例如,要计算在 0.25 处的函数值,我们可以使用以下代码:val = ppval(yy,0.25);上述代码将返回插值函数在 0.25 处的函数值。
我们还可以使用 plot() 函数来可视化样条插值函数的结果,例如:plot(x,y,'o',0:0.01:2,ppval(yy,0:0.01:2));该代码将绘制原始数据点以及样条插值函数在给定范围内的函数值。
样条插值函数是一种强大的数值计算技术,可以用于许多科学和工程应用。
在 MATLAB 中,我们可以使用 spline() 函数轻松实现样条插值函数的操作并可视化结果。
例5.6.1 给定以下数据, 求出三次样条函数,并计算函数分别在-0.15,-0.05,0.05,0.18,0.25处的近似值,并作图。
解:编程如下:clearx=[0.1,0.2,0.15,0,-0.2,0.3];y=[0.95,0.84,0.86,1.06,1.50,0.72];pp=spline(x,y);pp.coefsxx=[-0.15,-0.05,0.05,0.18,0.25];yy=ppval(pp,xx) %or:yy=spline(x,y,xx)fnplt(pp,'k')hold onplot(x,y,'o')hold onplot(xx,yy,'r*')运行结果:ans =-36.3850 21.8592 -5.1164 1.5000-36.3850 0.0282 -0.7390 1.0600227.6995 -10.8873 -1.8249 0.9500-143.0047 23.2676 -1.2059 0.8600-143.0047 1.8169 0.0484 0.8400yy =1.2943 1.1016 1.0186 0.8409 0.8291写出三次样条函数为:⎪⎪⎪⎩⎪⎪⎪⎨⎧≤≤+-+-+--≤≤+---+--≤≤+-----≤≤+-+-≤≤-++-+++-=3.02.0,84.0)2.0(0484.0)2.0(8169.1)2.0(0047.1432.015.0,86.0)15.0(2059.1)15.0(2676.23)15.0(0047.14315.01.0,95.0)1.0(8249.1)1.0(8873.10)1.0(6995.2271.00,06.1739.00282.0385.3602.0,5.1)2.0(1164.5)2.0(859.21)2.0(385.36)(2323232323x x x x x x x x x x x x x x x x x x x x x s练习: 给定函数11)(2+=x x f ,55≤≤-x ,节点)10,,1,0(5 =+-=k k x k ,用三次样条插值求)(x S 及在某些点处的近似值.MATLAB 实现:yi=interp1(x,y,xi) %根据数据(x,y)给出分段线性插值函数的xi 的值yi.yi=interp1(x,y,xi,’spline’) %使用三次样条插值.yi=interp1(x,y,xi,’cubic’) %使用分段三次插值.Yi=spline(x,y,xi) %同yi=interp1(x,y,xi,’spline’)pp=spline(x,y) %返回三次样条插值的分段多项式(pp)形式结构(使用“非扭结”端点条件).MATLAB 中的spline 函数可以对数据点进行三次插值,默认的边界条件强制使得插值第一段的三次项系数与第二段的三次项系数相同,最后一段与倒数第二段的多项式系数相同,这就所谓的非扭结(not-a-knot )条件. 1一维插值函数Interp1()命令格式: yi=interp1(x,y,xi,’method’)x 为插值节点构成的向量,y 为插值节点函数值构成的向量,yi 是被插值点xi 的插值结果,‘method ‘是采用的插值方法,缺省时表示分线段性插值,’nearest ‘为最邻近插值;’linear ‘为分线段性插值;’spline’为三次样条插值;’pchip’为分段Hermite 插值;’cubic’为分段Hermite 插值例子:画出y=sin(x)在区间[0 10]的曲线,并在曲线上插值节点xk=k,k=0,1…10及函数值,画出分段线性插值折线图x=0:10;y=sin(x);xi=0:0.25:10;yi1=interp1(x,y,xi,'nearest');yi2=interp1(x,y,xi,'linear');yi3=interp1(x,y,xi,'spline');yi4=interp1(x,y,xi,'pchip');yi5=interp1(x,y,xi,'cubic');subplot(1,5,1)plot(x,y,'o',xi,yi1,'k--',xi,sin(xi),'k:');title('\bfNearest');subplot(1,5,2)plot(x,y,'o',xi,yi2,'k--',xi,sin(xi),'k:');title('\bfLinear');subplot(1,5,3)plot(x,y,'o',xi,yi3,'k--',xi,sin(xi),'k:');title('\bfSpline');subplot(1,5,4)plot(x,y,'o',xi,yi4,'k--',xi,sin(xi),'k:');title('\bfPchip');subplot(1,5,1)plot(x,y,'o',xi,yi5,'k--',xi,sin(xi),'k:');title('\bfCubic');spline()为三次样条函数命令格式1:yi=spline(x,y,xi),意义等同于yi=interp1(x,y,xi,'spline')命令格式2:pp=spline(x,y) ,输出三次样条函数分段表示的结构pchip()命令格式与spline()完全相同csape()为可输入边界条件的三次样条函数命令格式:pp=csape(x,y,conds,valconds),x为插值节点构成的向量,y为插值节点函数值构成的向量;conds为边界类型,缺省为非扭结边界条件;valconds表示边界值。
MATLAB 中的曲线拟合和插值在大量的使用领域中,人们经常面临用一个分析函数描述数据(通常是测量值)的任务。
对这个问题有两种方法。
在插值法里,数据假定是正确的,要求以某种方法描述数据点之间所发生的情况。
这种方法在下一节讨论。
这里讨论的方法是曲线拟合或回归。
人们设法找出某条光滑曲线,它最佳地拟合数据,但不必要经过任何数据点。
图11.1说明了这两种方法。
标有'o'的是数据点;连接数据点的实线描绘了线性内插,虚线是数据的最佳拟合。
11.1 曲线拟合曲线拟合涉及回答两个基本问题:最佳拟合意味着什么?应该用什么样的曲线?可用许多不同的方法定义最佳拟合,并存在无穷数目的曲线。
所以,从这里开始,我们走向何方?正如它证实的那样,当最佳拟合被解释为在数据点的最小误差平方和,且所用的曲线限定为多项式时,那么曲线拟合是相当简捷的。
数学上,称为多项式的最小二乘曲线拟合。
如果这种描述使你混淆,再研究图11.1。
虚线和标志的数据点之间的垂直距离是在该点的误差。
对各数据点距离求平方,并把平方距离全加起来,就是误差平方和。
这条虚线是使误差平方和尽可能小的曲线,即是最佳拟合。
最小二乘这个术语仅仅是使误差平方和最小00.20.40.60.81-2024681012xy =f (x )Second O rder C urv e Fitting图11.1 2阶曲线拟合在MATLAB 中,函数polyfit 求解最小二乘曲线拟合问题。
为了阐述这个函数的用法,让我们以上面图11.1中的数据开始。
» x=[0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1]; » y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];为了用polyfit ,我们必须给函数赋予上面的数据和我们希望最佳拟合数据的多项式的阶次或度。
如果我们选择n=1作为阶次,得到最简单的线性近似。
三次样条众所周知,使用高阶多项式的插值常常产生病态的结果。
目前,有多种消除病态的方法。
在这些方法中,三次样条是最常用的一种。
在MATLAB中,实现基本的三次样条插值的函数有spline,ppval,mkpp和unmkpp。
在这些函数中,仅spline 在《MATLAB参考指南》中有说明。
下面几节,将展示在M文件函数中实现三次样条的基本特征基本特征在三次样条中,要寻找三次多项式,以逼近每对数据点间的曲线。
在样条术语中,这些数据点称之为断点。
因为,两点只能决定一条直线,而在两点间的曲线可用无限多的三次多项式近似。
因此,为使结果具有唯一性。
在三次样条中,增加了三次多项式的约束条件。
通过限定每个三次多项式的一阶和二阶导数,使其在断点处相等,就可以较好地确定所有内部三次多项式。
此外,近似多项式通过这些断点的斜率和曲率是连续的。
然而,第一个和最后一个三次多项式在第一个和最后一个断点以外,没有伴随多项式。
因此必须通过其它方法确定其余的约束。
最常用的方法,也是函数spline所采用的方法,就是采用非扭结(not-a-knot)条件。
这个条件强迫第一个和第二个三次多项式的三阶导数相等。
对最后一个和倒数第二个三次多项式也做同样地处理。
基于上述描述,人们可能猜想到,寻找三次样条多项式需要求解大量的线性方程。
实际上,给定N个断点,就要寻找N-1个三次多项式,每个多项式有4个未知系数。
这样,所求解的方程组包含有4*(N-1)个未知数。
把每个三次多项式列成特殊形式,并且运用各种约束,通过求解N个具有N个未知系数的方程组,就能确定三次多项式。
这样,如果有50个断点,就有50个具有50个未知系数的方程组。
幸好,用稀疏矩阵,这些方程式能够简明地列出并求解,这就是函数spline 所使用的计算未知系数的方法。
《数值分析》课外课堂大作业论文题目:基于多项式插值与三次样条插值曲线拟合的比较姓名:学号:学院:专业方向:联系方式:(QQ号)(手机号)导师姓名:完成人(亲笔)签字基于多项式插值与三次样条插值曲线拟合的比较摘要:在数值计算中经常要计算函数,当函数只在有限点集上给定函数值要包含改点集的区间上用公式给出函数的简单表达式,这就涉及在已知区间上用简单函数逼近已知复杂函数问题。
本文为了解决这类问题就采用多项式插值与三次样条插值两种插值法并利用MATLAB数值分析软件进行编程,实现相应数据的曲线拟合以获得最佳曲线模型与相应数据的曲线拟合,选出最优的插值法以解决所给数据的曲线拟合问题。
关键词:函数;多项式插值;三次样条插值;曲线拟合;MATLABAbstract:In numerical analysis ,the function value is often calculated .when the function is only given a function point set ,the simple expression of the function is given by the interval .which involves the use of a simple function to approximate the known complex function .in order to solve this problem ,we use polynomial interpolation and cubic spline interpolation tow kind of interpolation method and use MATLAB numerical analysis software to program ,to achieve the curve fitting of the corresponding date to obtain the best cure fitting ,and to choose the best interpolation method to solve the problem of curve fitting to the date.Keyword: Function ; Polynomial interpolation ; Cubic spline interpolation ; Fitting of a curve ; MATLAB前言现代科学研究中,物理量之间的相互关系通量是用函数来描述的,许多实际问题都用函数y=f(x)来表示某种内在规律的数量关系其中相当一部分函数是通过试验或观测得到的也有少量函数关系是由经典物理分析推导得到的,但许多实际问题很难用经典理论分析得出,因为虽然f(x)在某个区间[a,b]上是存在的,有的还是连续的,但往往这个f(x)并不包含我们所得函数表的所有值因此我们希望根据给定的函数表做一个即能反应函数f(x)的特行,又便于计算的简单函数p(x),用p(x)近似f(x),这样确定的p(x)就是我们希望得得到的插值函数。
Matlab中的插值和平滑方法1. 引言在数值分析和数据处理中,插值和平滑是常用的技术手段,可以用于填补数据的空缺以及降低数据中的噪声。
Matlab作为一种强大的数值计算和数据处理软件,提供了丰富的插值和平滑方法,本文将介绍其中的一些常用方法及其应用。
2. 插值方法2.1 线性插值线性插值是最简单的一种插值方法,它假设待插值函数在相邻数据点之间是线性变化的。
Matlab中提供了interp1函数实现线性插值,可以通过设定插值点的横坐标向量和已知数据点的横坐标向量,以及对应的纵坐标向量,得到插值结果。
2.2 分段插值分段插值是一种更精确的插值方法,它假设待插值函数在相邻数据点之间是分段线性变化的。
Matlab中的interp1函数也可以实现分段插值,通过指定'linear'插值方法和 'pchip'插值方法,可以得到不同的插值结果,前者得到的结果比较平滑,而后者更接近原始数据的形状。
2.3 样条插值样条插值是一种更高阶的插值方法,它假设待插值函数在相邻数据点之间是多项式变化的。
Matlab中的spline函数可以实现三次样条插值,它通过计算每个数据点处的二阶导数,得到一个以每个数据点为节点的三次多项式函数。
样条插值可以更加精确地还原数据,但也容易受到离群点的干扰。
3. 平滑方法3.1 移动平均移动平均是一种常用的平滑方法,它通过计算数据点周围一定范围内的平均值,得到平滑后的结果。
Matlab中的smoothdata函数提供了不同的平滑方法,包括简单移动平均、指数移动平均和加权移动平均等,可以根据具体需求选择适当的方法。
3.2 Savitzky-Golay滤波Savitzky-Golay滤波是一种基于最小二乘法的平滑方法,它通过拟合多项式曲线来实现数据的平滑。
Matlab中的sgolay函数可以实现Savitzky-Golay滤波,通过指定不同的拟合阶数和窗口大小,可以得到不同程度的平滑结果。
实验五插值法5.1实验目的掌握插值的基本思想与方法,会借助数学软件Matlab求解并讨论其收敛性.5.2实验内容1、Lagrange插值法、Newton插值法的Matlab求解方法,在对Runge现象的观察基础上,了解高次插值的不稳定性及其改进方法;2、熟悉Matlab中的插值求解函数,掌握三次样条插值的Matlab求解;3、会求解某些简单的实际问题.5.3实验步骤5.5.1 Lagrange插值法和Newton插值法教师示范:通过计算实例,学习Lagrange插值法和Newton插值法的Matlab 程序编制及其应用.实例1. 拉格朗日插值法计算插值.已知:x:0 1 2 3y:-5 -6 -1 16,求x从0到3间隔0.1的函数值.实例2. 拉格朗日插值法求插值多项式.程序见interpEg3.m.Lagrange插值:自编程序,interpH.m的M文件,yi=interpH(x,y,xi).Newton插值:自编程序, newinter.m的M文件,yi=newinter(x,y,xi).5.5.2 Runge现象教师示范:观察Rung现象,了解高次插值的不稳定性.程序参见rungeinterp.m.5.5.3 分段低次插值和三次样条插值学习Matlab的插值求解命令。
分段线性插值: yi=interp1(x,y,xi,’linear’,’pp’)三次样条插值:yi=interp1(x,y,xi,’spline’,’pp’)或yi=spline(x,y,xi)二维插值: interp2(x,y,z,xi,yi,’spline’)griddata(x,y,z,xi,yi)教师示范:机翼下轮廓线,见PPT文件。
学生练习1:对5.5.2中的问题分别采用分段线性插值和三次样条插值求解,了解消除Rung现象的基本思路和低次插值的优点.学生练习2:画手练习.在Matlab中输入命令:figure('position',get(0,'screensize'))axes('position',[0 0 1 1])[x,y] = ginput;将你的手放在屏幕上,沿着手的边界,用鼠标点击选取一些点,按回车键结束选取。
3.3 插值与拟合的MATLAB实现简单的插值与拟合可以通过手工计算得出,但复杂的只能求助于计算机了。
3.3.1 线性插值在MATLAB 中,一维的线性插值可以用函数interpl 来实现。
函数interpl 的调用格式如下:yi = interpl ( x , y , xi ) ,其中yi 表示在插值向量xi 处的函数值,x 与y 是数据点。
这个函数还有如下两种形式:yi = interpl(y , xi),省略x,x 此时为l : N,其中N 为向量y 的长度。
yi = interpl(x , y , xi , method ) ,其中method 为指定的插值方法,可取以下凡种:nearest :最近插值。
linear :线性插值。
spline :三次样条插值。
cubic :三次插值。
注意:对于上述的所有的调用格式,都要求向量x 为单调。
例如:对以下数据点:( 2 * pi , 2 ) , ( 4 * pi , 3 ) , ( 6 * pi , 5 ) , ( 8 * pi , 7 ) , ( 10 * pi , 11 ) , ( 12 * pi , 13 ) , ( 14 * pi , 17) 进行插值,求x = pi , 6 的函数值。
>> x=linspace(0, 2 * pi, 8 );>> y=[2, 3, 5, 7, 11, 13, 17, 19 ];>> xl=[pi , 6 ];>> yl=interpl(x, y, xl)yl =90000 1836903.3.2 Lagrange 插值Lagrange 插值比较常用,是MATLAB 中相应的函数,但根据Lagrange 插值函数公式,可以用M 文件实现:Lagrange.mfunctions = Larange(x, y, x0 )% Lagrange 插值,x 与y 为已知的插值点及其函数值,x0 为需要求的插值点的值nx = length( x );ny = length( y );if nx ~=nywaming( ‘向量x 与y 的长度应该相同’)return;endm = length ( x0 ) ;%按照公式,对需要求的插值点向量x0 的元素进行计算for i = l: mt =0.0;for j = l : nxu = 1.0;for k = l : nxif k~=ju=j * ( x0( i )-x ( k ) ) / ( x( j )-( k ) ) ;endendt = t + u * y( j );ends( i ) = t ;endreturn例如:对(l , 2 ) , ( 2 , 4 ) , ( 3 , 6 ) , ( 4 , 8 ) , ( 5 , 10 ) 进行Lagrange 插值,求x = 23 , 3.7 的函数值。
MATLAB插值法引言MATLAB是一种高级编程语言和环境,特别适用于数值计算和数据可视化。
插值法是一种在给定有限的数据点的情况下,通过构造插值函数来估计其他数据点的方法。
在MATLAB中,有多种插值方法可供选择,例如拉格朗日插值、牛顿插值和样条插值等。
本文将详细介绍MATLAB中常用的插值方法及其应用。
一、拉格朗日插值法拉格朗日插值法是一种多项式插值方法,通过构造一个满足给定数据点要求的多项式函数,来估计其他数据点的函数值。
其基本思想是通过一个多项式函数对已知数据点进行拟合,以实现函数值的估计。
以下是使用MATLAB实现拉格朗日插值法的步骤:1.确定待插值的数据点集合,假设有n个数据点。
2.构造拉格朗日插值多项式。
拉格朗日插值多项式的表达式为:其中,为拉格朗日基函数,其表达式为:3.利用构造的拉格朗日插值多项式求解其他点的函数值。
二、牛顿插值法牛顿插值法是一种基于差商的插值方法,通过构造一个n次多项式函数来拟合已知数据点,并利用差商的性质来求解其他点的函数值。
使用MATLAB实现牛顿插值法的步骤如下:1.确定待插值的数据点集合,假设有n个数据点。
2.计算差商表。
差商表的计算公式为:3.构造牛顿插值多项式。
牛顿插值多项式的表达式为:4.利用构造的牛顿插值多项式求解其他点的函数值。
三、样条插值法样条插值法是一种通过多段低次多项式来逼近原始数据,以实现光滑插值的方法。
它在相邻数据点处保持一定的连续性,并通过边界条件来确定插值函数的特性。
以下是使用MATLAB实现样条插值法的步骤:1.确定待插值的数据点集合,假设有n个数据点。
2.根据数据点的个数确定样条插值的次数。
一般情况下,插值多项式的次数小于或等于n-1。
3.利用边界条件构造样条插值函数。
常用的边界条件有:自然边界、固定边界和周期边界。
4.利用MATLAB中的插值函数csape或interp1等进行样条插值。
5.利用样条插值函数求解其他点的函数值。
插值插值法是实用的数值方法,是函数逼近的重要方法。
在生产和科学实验中,自变量x与因变量y的函数y = f(x)的关系式有时不能直接写出表达式,而只能得到函数在若干个点的函数值或导数值。
当要求知道观测点之外的函数值时,需要估计函数值在该点的值。
如何根据观测点的值,构造一个比较简单的函数y=φ(x),使函数在观测点的值等于已知的数值或导数值。
用简单函数y=φ(x)在点x处的值来估计未知函数y=f(x)在x点的值。
寻找这样的函数φ(x),办法是很多的。
φ(x)可以是一个代数多项式,或是三角多项式,也可以是有理分式;φ(x)可以是任意光滑(任意阶导数连续)的函数或是分段函数。
函数类的不同,自然地有不同的逼近效果。
在许多应用中,通常要用一个解析函数(一、二元函数)来描述观测数据。
根据测量数据的类型:1.测量值是准确的,没有误差。
2.测量值与真实值有误差。
这时对应地有两种处理观测数据方法:1.插值或曲线拟合。
2.回归分析(假定数据测量是精确时,一般用插值法,否则用曲线拟合)。
MATLAB中提供了众多的数据处理命令。
有插值命令,有拟合命令,有查表命令。
2.2.1 插值命令命令1 interp1功能一维数据插值(表格查找)。
该命令对数据点之间计算内插值。
它找出一元函数f(x)在中间点的数值。
其中函数f(x)由所给数据决定。
各个参量之间的关系示意图为图2-14。
Yi:插值点图2-14 数据点与插值点关系示意图格式yi = interp1(x,Y,xi) %返回插值向量yi,每一元素对应于参量xi,同时由向量x与Y的内插值决定。
参量x指定数据Y的点。
若Y为一矩阵,则按Y的每列计算。
yi是阶数为length(xi)*size(Y,2)的输出矩阵。
注意,这个2是列。
yi = interp1(Y,xi) %假定x=1:N,其中N为向量Y的长度,或者为矩阵Y的行数。
yi = interp1(x,Y,xi,method) %用指定的算法计算插值:’nearest’:最近邻点插值,直接完成计算;’linear’:线性插值(缺省方式),直接完成计算;’spline’:三次样条函数插值。
三次样条插值多项式
——计算物理实验作业四
陈万物理学2013级
主程序:
clear,clc;
format rat
x = [1,4,9,16,25,36,49,64];
y = [1,2,3,4,5,6,7,8];
f1 = ;
fn = 1/16;
[a,b,c,d,M,S] = spline(x,y,f1,fn);
子程序1:
function [a,b,c,d,M,S]=spline(x,y,f1,fn)
% 三次样条插值函数
% x是插值节点的横坐标
% y是插值节点的纵坐标
% u是插值点的横坐标
% f1是左端点的一阶导数
% fn是右端点的一阶导数
% a是三对角矩阵对角线下边一行
% b是三对角矩阵对角线
% c是三对角矩阵对角线上边一行
% S是插值点的纵坐标
n = length(x);
h = zeros(1,n-1);
deltay = zeros(1,n);
miu = zeros(1,n-1);
lamda = zeros(1,n-1);
d = zeros(1,n-1);
for j = 1:n-1
h(j) = x(j+1)-x(j);
deltay(j) = y(j+1)-y(j);
end % 得到h矩阵
for j = 2:n-1
sumh = h(j-1) + h(j);
miu(j) = h(j-1) / sumh;
lamda(j) = h(j) / sumh;
d(j) = 6*( deltay(j)/h(j)-(deltay(j-1)/h(j-1)))/sumh; end
% 根据第一类边界条件,作如下规定
lamda(1) = 1;
d(1) = 6*(deltay(1)/h(1)-f1)/h(1);
miu(1) = 1;
d(n) = 6*(fn-deltay(n-1)/h(n-1))/h(n-1);
% 输出三对角矩阵的a,b,c
a = miu;
b = 2*ones(1,n);
c = lamda;
M = chase(a,b,c,d); %调用chase函数得到M
sym u;
for j = 1:n-1
u = x(j)::x(j+1);
v = ones(size(u));
S = (M(j)*(x(j+1)*v-u).^3/(6*h(j))+M(j+1)*(u-x(j)*v).^3/(6*h(j))... +(y(j)-M(j)*h(j)^2/6)*(x(j+1)*v-u)/h(j)+(y(j+1)...
-M(j+1)*h(j)^2/6)*(u-x(j)*v)/h(j));
plot(u,S,'-k');
hold on
end
plot(x,y,'-.*r');
xlabel('x'),ylabel('y'),title('cubic spline interp');
end
子程序2:
function M = chase(a,b,c,f)
% 追赶法求解三对角矩阵方程,Ax=f
% a是对角线下边一行的元素
% b是对角线元素
% c是对角线上边一行的元素
n = length(b);
beta = ones(1,n-1);
y = ones(1,n);
M = ones(1,n);
for i = (n-1):(-1):1
a(i+1) = a(i);
end
% 将a矩阵和n对应
beta(1) = c(1)/b(1);
for i = 2:(n-1)
beta(i) = c(i)/( b(i)-a(i)*beta(i-1) );
end
y(1) = f(1)/b(1);
for i = 2:n
y(i) = (f(i)-a(i)*y(i-1))/(b(i)-a(i)*beta(i-1)); end
M(n) = y(n);
for i = (n-1):(-1):1
M(i) = y(i)-beta(i)*M(i+1);
end
end
三次样条插值结果:
与拉格朗日插值作对比:
分析图一知,三次样条插值结果与预期结果吻合得很好,曲线平滑连续性好,在左右短点处的
小区间也吻合得很好,可以延伸到区间[a,b]外的一小段,用最邻近的小区间插值函数可以近似求得[a,b]区间外一小段范围内的函数值。
分析图一和图二,三次样条插值在x(i)和x(i+1)间隔很小的地方插值效果稍差,间隔稍大时插值效果非常好。
拉格朗日插值则不同,x(i)和x(i+1)间隔越大插值效果越差。