matlab实现插值法和曲线拟合电子教案
- 格式:docx
- 大小:187.24 KB
- 文档页数:14
Matlab多项式插值与曲线拟合
一、基本统计处理
1、查取最大值
MAX函数的命令格式有:
[Y,I]= max (X):将max(X)返回矩阵X的各列中的最大元素值及其该元素的位置赋予行向量Y与I;当X为向量时,则Y与I为单变量。
[Y,I]=max(X,[],DIM):当DIM=1时按数组X的各列查取其最大的元素值及其该元素的位置赋予向量Y与I;当DIM=2时按数组X的各行查取其最大的元素值及其该元素的位置赋予向量Y 与I.
max(A,B):返回一个与A,B同维的数组,其每一个元素是由A,B同位置上的元素的最大值组成。
【例1】查找下面数列x的最大值。
x=[3 5 9 6 1 8] % 产生数列x
x = 3 5 9 6 1 8
y=max(x) % 查出数列x中的最大值赋予y
y = 9
[y,l]=max(x) % 查出数列x中的最大值及其该元素的位置赋予y,l
y = 9
l = 3
【例2】分别查找下面3×4的二维数组x中各列和各行元素中的最大值。
x=[1 8 4 2;9 6 2 5;3 6 7 1] % 产生二维数组x
x = 1 8 4 2。
MATLAB软件中软件拟合与插值运算的方法内容目录
1MATLAB中软件拟合与插值运算的方法1
1.1拟合函数的选择1
1.1.1线性拟合1
1.1.2非线性拟合2
1.2拟合函数的求解2
1.2.1直接法2
1.2.2迭代法3
1.3MATLAB插值函数4
1.3.1样条插值函数4
1.3.2拉格朗日插值函数5
1.3.3指数插值函数5
结论6
近来,随着科学技术的进步,数据采集技术的发展,大量的实验数据和实验结果越来越多,如何合理地分析处理数据,描绘实际趋势,就变得十分重要,MATLAB中的软件拟合与插值是目前应用最多的数据处理技术之一、本文介绍了MATLAB中软件拟合与插值运算的方法及其具体实现。
1.1拟合函数的选择
1.1.1线性拟合
线性拟合是指拟合函数可以用一元线性方程描述,MATLAB中的拟合
函数有polyfit、polyval和 polyconf等。
其中,polyfit函数用来根据
输入的拟合数据拟出一元多项式,polyval函数用来求出拟合后的拟合值,polyconf函数用来计算拟合的参数的置信范围。
例如,用polyfit函数
拟合下面的数据,输入x = [1 2 3 4 5]和y = [4.3 7.3 11.1 14.1
18.4],拟出的拟合函数为y = 4.1 + 2.3x,即拟合函数为y = 4.1 +
2.3x。
1.1.2非线性拟合。
插值法和曲线拟合电子科技大学摘要:理解拉格朗日多项式插值、分段线性插值、牛顿前插,曲线拟合,用matlab 编程求解函数,用插值法和分段线性插值求解同一函数,比较插值余项;用牛顿前插公式计算函数,计算函数值;对于曲线拟合,用不同曲线拟合数据。
关键字:拉格朗日插值多项式;分段线性插值;牛顿前插;曲线拟合引言:在数学物理方程中,当给定数据是不同散点时,无法确定函数表达式,求解函数就需要很大的计算量,我们有多种方法对给定的表格函数进行求解,我们这里,利用插值法和曲线拟合对函数进行求解,进一步了解函数性质,两种方法各有利弊,适合我们进行不同的散点函数求解。
正文:一、插值法和分段线性插值 1拉格朗日多项式原理对某个多项式函数,已知有给定的k + 1个取值点:其中对应着自变量的位置,而对应着函数在这个位置的取值。
假设任意两个不同的x j 都互不相同,那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为:其中每个为拉格朗日基本多项式(或称插值基函数),其表达式为:[3]拉格朗日基本多项式的特点是在 上取值为1,在其它的点 上取值为0.2分段线性插值原理给定区间[a ,b ], 将其分割成a=x 0 <x 1 〈…〈x n =b, 已知函数y= f(x ) 在这些插值结点的函数值为 y k =f(x k )(k=0,1,…,n)求一个分段函数I h (x), 使其满足:(1) I h (x k )=y k ,(k=0,1,…,n) ;(2) 在每个区间[x k ,x k+1 ] 上,I h (x )是个一次函数。
易知,I h (x)是个折线函数, 在每个区间[x k ,x k+1 ]上,(k=0,1,…,n)k 1k k1k 1k k 1k k k ,1)()()(x x x x x f x x x x x f x L --+--=++++,于是, I h (x )在[a,b ]上是连续的,但其一阶导数是不连续的.3拉格朗日插值多项式算法 ○1输入,(0,1,2,,)i i x y i n =,令0)(=x L n .错误!对0,1,2,,i n =,计算0,()()/()ni j i j j j il x x x x x -≠=--∏()()()n n i i L x L x l x y ←−−+4分段线性插值算法错误!输入(x k ,y k ),k=0,1,…,n;○2计算k 1k k1k 1k k 1k k k ,1)()()(x x x x x f x x x x x f x L --+--=++++5插值法和分段线性插值程序按下列数据分别作五次插值和分段线性插值,画出两条插值曲线以及给定数据点。
第五次课Matlab曲线拟合与插值运算一、本次课学习要点1、Matlab曲线拟合与插值运算2、符号表达式二、本次课教学重点利用M文件的Matlab曲线拟合与插值运算三、教学基本内容1、曲线拟合在许多应用领域中,人们经常需要从一系列已知离散点上的数据集[(x1,y1),(x2,y2)],…(x n,y n)]得到一个解析函数y=f (x)。
得到的解析函数f(x)应当在原离散点x i上尽可能接近给定的y i的值。
这一过程称为曲线拟合。
最常用的曲线拟合是最小二乘法曲线拟合。
似合结果可使误差的平方和最小,MATLAB提供的函数polyfit,根据给定的自变量数组x和函数数组y,按照拟合的阶数要求自动求解满足最小二乘意义的一阶或高阶解析函数f(x),使用很方便。
为了说明这个问题,我们取以下函数为例:2-=y⨯25.0xx=0:0.1:1;for i=1:length(x);y(i)=0.5-2*x(i)^2;end显示为:将y值进行一定的修改,输入如下的程序y=[0.52 0.45 0.4 0.35 0.18 0.02 -0.25 -0.4 -0.81 -1.1 -1.5]m=1;fxy1=polyfit(x,y,m)m=2;fxy2=polyfit(x,y,m)y1=polyval(fxy1,x) %多项式求值,x为输入值,fx1为一次拟合出来的多项式y2=polyval(fxy2,x)plot(x,y,'o',x,y1,'k:',x,y2,'k')显示为:2、插值运算与曲线拟合不同,插值运算不是试图找出适合于所有自变量数组x的全局最优拟合函数Y=f(x〕,而是要找到一个解析函数连接自变量相邻的两个点(xi,xi十1),由此还可以找到两点间的数值。
根据自变量的维数不同,插值方法可以分为一维插值y=f(x)和二维插值z=f(x,y)等。
在许多工程问题上,我们只能获得无规律的离散点上的数值,插值可以帮助我们得到近似的连续过程,便于用数学的解析方法对已有数据进行处理和运算,因此是一个很有用的工具。
目录【一维插值】interp1 (1)yi = interp1(x,y,xi,method) (1)例1 (1)例2 (2)【二维插值】interp2 (4)ZI = interp2(X,Y,Z,XI,YI,method) (4)插值方式比较示例 (4)例3 (8)例4 (9)【三角测量和分散数据插值】 (13)【数据拟合】 (17)例5 (17)例6 (18)【一维插值】interp1yi = interp1(x,y,xi,method)例1在1-12 的11 小时内,每隔1 小时测量一次温度,测得的温度依次为:5,8,9,15,25,29,31,30,22,25,27,24。
试估计每隔1/10 小时的温度值。
建立M文件temp.mhours=1:12;temps=[5 8 9 15 25 29 31 30 22 25 27 24];h=1:0.1:12;t=interp1(hours,temps,h,'spline');plot(hours,temps,'kp',h,t,'b');35302520151050 2 4 6 8 10 12例2已知飞机下轮廓线上数据如下,求x每改变0.1时的y值。
X0357911 12131415 Y0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6建立M文件plane.mx0=[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=interp1(x0,y0,x,'nearest');y2=interp1(x0,y0,x);y3=interp1(x0,y0,x,'spline');plot(x0,y0,'kp',x,y1,'r')2.521.510.50 5 10 15 plot(x0,y0,'kp',x,y2,'r')2.521.510.50 5 10 15 plot(x0,y0,'kp',x,y3,'r')2.521.510.50 5 10 15 【二维插值】interp2ZI = interp2(X,Y,Z,XI,YI,method)插值方式比较示例用较大间隔产生peaks 函数数据点[x,y] = meshgrid(-3:1:3);z = peaks(x,y);surf(x,y,z)642-2-4-642 40 2-2 -2-4 -4●产生一个较好的网格[xi,yi] = meshgrid(-3:0.25:3);●利用最近邻方式插值zi1 = interp2(x,y,z,xi,yi,'nearest');surf(xi,yi,zi1)●双线性插值方式zi2 = interp2(x,y,z,xi,yi,'bilinear');surf(xi,yi,zi2)●双立方插值方式zi3 = interp2(x,y,z,xi,yi,'bicubic');surf(xi,yi,zi3)●不同插值方式构造的等高线图对比contour(xi,yi,zi1)321-1-2-3-3 -2 -1 0 1 2 3 contour(xi,yi,zi2)321-1-2-3-3 -2 -1 0 1 2 3 contour(xi,yi,zi3)321-1-2-3-3 -2 -1 0 1 2 3例3测得平板表面3*5 网格点处的温度分别为:82 81 80 82 8479 63 61 65 8184 84 82 85 86试作出平板表面的温度分布曲面z=f(x,y)的图形。
在Matlab中如何进行数据插值与拟合引言:数据处理是科学研究与工程开发中不可或缺的环节之一。
而数据插值和拟合则是数据处理中常用的技术手段。
在Matlab这一强大的数值分析工具中,提供了丰富的函数与工具箱,使得数据插值与拟合变得更加便捷高效。
本文将详细阐述在Matlab中如何进行数据插值与拟合,并介绍几个常用的插值与拟合方法。
一、数据插值数据插值是通过已知的有限个数据点,推导出数据点之间未知位置上的数值。
在Matlab中,可以利用interp1函数进行数据插值。
假设我们有一组离散的数据点,存储为两个向量x和y。
那么,可以通过以下步骤进行数据插值:1. 调用interp1函数,并传入x和y作为输入参数。
```matlabxi = linspace(min(x), max(x), n);yi = interp1(x, y, xi, '方法');```其中,xi是插值点的位置,min和max分别是x向量的最小值和最大值,n是插值点的数量。
'方法'是要使用的插值方法,可以选择线性插值(method='linear')、样条插值(method='spline')等。
2. 绘制插值结果曲线。
```matlabplot(x, y, 'o', xi, yi)legend('原始数据','插值结果')```使用plot函数可以绘制原始数据点和插值结果的曲线。
通过设置不同的插值方法和插值点的数量,可以探索不同的插值效果。
二、数据拟合数据拟合是通过已知的一组数据点,找到一个符合数据趋势的函数模型。
在Matlab中,可以利用polyfit函数进行多项式拟合。
假设我们有一组离散的数据点,存储为两个向量x和y。
那么,可以通过以下步骤进行数据拟合:1. 调用polyfit函数,并传入x和y作为输入参数。
```matlabp = polyfit(x, y, n);```其中,n是多项式的次数,p是拟合多项式的系数。
Matlab中的曲线拟合与插值技巧在数据科学和工程领域中,曲线拟合和插值技术是常用的数学方法。
在Matlab 中,有许多工具和函数可用于处理这些技术。
本文将讨论Matlab中的曲线拟合和插值技巧,并介绍一些实际应用案例。
一、曲线拟合技术曲线拟合是根据已知数据点来构造一个与这些点最匹配的曲线模型。
在Matlab 中,常用的曲线拟合函数包括polyfit和lsqcurvefit。
1. polyfit函数polyfit函数是Matlab中一个功能强大的多项式拟合函数。
它可以拟合多项式曲线模型,并通过最小二乘法找到最佳拟合系数。
例如,我们有一组数据点(x,y),我们想要拟合一个二次多项式曲线来描述这些数据。
可以使用polyfit函数:```matlabx = [1, 2, 3, 4, 5];y = [2, 4, 6, 8, 10];degree = 2;coefficients = polyfit(x, y, degree);```在上述例子中,degree参数设置为2,表示拟合一个二次多项式曲线。
polyfit 函数将返回一个包含拟合系数的向量,可以用来构造拟合曲线。
2. lsqcurvefit函数lsqcurvefit函数是Matlab中一个用于非线性最小二乘拟合的函数。
与polyfit函数不同,lsqcurvefit函数可以用于拟合任意曲线模型,不局限于多项式。
例如,我们想要拟合一个指数函数曲线来拟合数据:```matlabx = [1, 2, 3, 4, 5];y = [1.1, 2.2, 3.7, 6.5, 12.3];model = @(params, x) params(1)*exp(params(2)*x);params0 = [1, 0];estimated_params = lsqcurvefit(model, params0, x, y);```在上述例子中,model是一个函数句柄,表示要拟合的曲线模型。
matlab学习——05插值和拟合(⼀维⼆维插值,拟合)05插值和拟合1.⼀维插值(1) 机床加⼯零件,试⽤分段线性和三次样条两种插值⽅法计算。
并求x=0处的曲线斜率和13<=x<=15范围内y的最⼩值。
x0=[0 3 5 7 9 11 12 13 14 15];y0=[0 1.2 1.7 2 2.1 2.0 1.8 1.2 1.0 1.6];x=0:0.1:15;% interp1现有插值函数,要求x0单调,'method'有% nearest 最近项插值 linear 线性插值% spline ⽴⽅样条插值 cubic ⽴⽅插值y1=interp1(x0,y0,x);y2=interp1(x0,y0,x,'spline');pp1=csape(x0,y0);y3=fnval(pp1,x);pp2=csape(x0,y0,'second');y4=fnval(pp2,x);[x',y1',y2',y3',y4']subplot(1,4,1)plot(x0,y0,'+',x,y1)title('Piecewise linear 分段线性')subplot(1,4,2)plot(x0,y0,'+',x,y2)title('spline1')subplot(1,4,3)plot(x0,y0,'+',x,y3)title('spline2')subplot(1,4,4)plot(x0,y0,'+',x,y4)title('second')dx=diff(x);dy=diff(y3);dy_dx=dy./dx;dy_dx0=dy_dx(1);ytemp=y3(131:151);ymin=min(ytemp);index=find(y3==ymin);xmin=x(index);[xmin,ymin](2) 已知速度的四个观测值,⽤三次样条求位移S=0.15到0.18上的vd(t)积分t 0.15 0.16 0.17 0.18vt 3.5 1.5 2.5 2.8format compact;% 已知速度的四个观测值,⽤三次样条求位移S=0.15到0.18上的vd(t)积分% t 0.15 0.16 0.17 0.18% vt 3.5 1.5 2.5 2.8clc,clearx0=0.15:0.01:0.18;y0=[3.5 1.5 2.5 2.8];% csape 三次样条插值,返回要求插值的的函数值pp=csape(x0,y0) % 默认的边界条件,Lagrange边界条件format long gxishu = pp.coefs % 显⽰每个区间上三次多项式的系数s=quadl(@(t)ppval(pp,t),0.15,0.18) % 求积分format % 恢复短⼩数的显⽰格式% 画图t=0.15:0.001:0.18;y=fnval(pp,t);plot(x0,y0,'+',t,y)pp =包含以下字段的 struct:form: 'pp'breaks: [0.1500 0.1600 0.1700 0.1800]coefs: [3×4 double]pieces: 3order: 4dim: 1xishu =1 ⾄2 列-616666.666666667 33500-616666.666666667 15000-616666.666666668 -3499.999999999993 ⾄4 列-473.333333333334 3.511.6666666666671 1.5126.666666666667 2.5s =0.0686252.⼆维插值(1) 丘陵测量⾼度。
插值法和曲线拟合电子科技大学摘要:理解拉格朗日多项式插值、分段线性插值、牛顿前插,曲线拟合,用matlab编程求解函数,用插值法和分段线性插值求解同一函数,比较插值余项;用牛顿前插公式计算函数,计算函数值;对于曲线拟合,用不同曲线拟合数据。
关键字:拉格朗日插值多项式;分段线性插值;牛顿前插;曲线拟合引言:在数学物理方程中,当给定数据是不同散点时,无法确定函数表达式,求解函数就需要很大的计算量,我们有多种方法对给定的表格函数进行求解,我们这里,利用插值法和曲线拟合对函数进行求解,进一步了解函数性质,两种方法各有利弊,适合我们进行不同的散点函数求解。
正文:一、插值法和分段线性插值1拉格朗日多项式原理对某个多项式函数,已知有给定的k + 1个取值点:其中对应着自变量的位置,而对应着函数在这个位置的取值。
假设任意两个不同的x j都互不相同,那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为:其中每个为拉格朗日基本多项式(或称插值基函数),其表达式为:[3]拉格朗日基本多项式的特点是在上取值为1,在其它的点上取值为0。
2分段线性插值原理给定区间[a,b], 将其分割成a=x0<x1<…<xn=b, 已知函数y= f(x) 在这些插值结点的函数值为y k =f(xk)(k=0,1,…,n)求一个分段函数Ih(x), 使其满足:(1) I h (x k )=y k ,(k=0,1,…,n) ;(2) 在每个区间[x k ,x k+1 ] 上,I h (x)是个一次函数。
易知,I h (x)是个折线函数, 在每个区间[x k ,x k+1 ]上,(k=0,1,…,n) k1k k 1k 1k k 1k k k ,1)()()(x x x x x f x x x x x f x L --+--=++++,于是, I h (x)在[a,b]上是连续的,但其一阶导数是不连续的。
3拉格朗日插值多项式算法○1输入,(0,1,2,,)i i x y i n = ,令0)(=x L n 。
四、实验结果(包括程序或图表、结论陈述、数据记录及分析等,可附页)1.数据插值:(1)轮船的甲板成近似半椭圆面形为了得到甲板的面积,首先测量得到横向最大相间8.534米;然后等间距地测得纵向高度,自左向右分别为:0.914, 5.060, 7.772, 8.717, 9.083, 9.144, 9.083, 8.992, 8.687, 7.376, 2.073,计算甲板的面积。
解:输入程序:x=0:0.8534:8.534;y=[0.914,5.060,7.772,8.717,9.083,9.144,9.083,8.992,8.687,7.376,2.073];xi=0:0.001:8.534;yi=interp1(x,y,xi,'spline');plot(x,y,'o',xi,yi)S=trapz(xi,yi)则运行结果如右图,即面积S为65.2824平方米。
(2) 山区地貌图在某山区(平面区域(0,2800)(0,2400)内,单位:米)测得一些地点的高程(单位:米)如下表所示,试作出该山区的地貌图和等高线图。
○2输入程序:x=0:400:2800;y=0:400:2400;z=[1430,1450,1470,1320,1280,1200,1080,940;1450,1480,1500,1550,1510,1430,1300,12 00;1460,1500,1550,1600,1550,1600,1600,1600;1370,1500,1200,1100,1550,1600,1550,1 380;1270,1500,1200,1100,1350,1450,1200,1150;1230,1390,1500,1500,1400,900,1100,1060;1180,1320,1450,1420,1400,1300,700,900];[xi,yi]=meshgrid(0:10:2800,0:10:2400);zi=interp2(x,y,z,xi,yi,'spline');subplot(2,2,1)%在二维上标注等高线[C,H]=contour(xi,yi,zi);clabel(C,H)xlabel('x')ylabel('y')title('在二维上标注等高线')%在三维上标注等高线subplot(2,2,2)[C,H]=contour3(xi,yi,zi);clabel(C,H)xlabel('x')ylabel('y')zlabel('z')title('在三维上标注等高线')%带有基准平面的网格图subplot(2,2,3)meshz(xi,yi,zi)%axis off tight;xlabel('x')ylabel('y')zlabel('z')title('带有基准平面的网格图线')画出的等高线图为:2.数据拟合Malthus人口指数增长模型中参数从1790—1980年间美国每隔10年的人口记录如下表:。
插值与拟合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=-10hold offx=[-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 onplot(x0,y1,'b:')legend('插值曲线','原数据曲线')因此插值多项式一般不要超过四次为宜。
m a t l a b实现插值法和曲线拟合插值法和曲线拟合电子科技大学摘要:理解拉格朗日多项式插值、分段线性插值、牛顿前插,曲线拟合,用matlab编程求解函数,用插值法和分段线性插值求解同一函数,比较插值余项;用牛顿前插公式计算函数,计算函数值;对于曲线拟合,用不同曲线拟合数据。
关键字:拉格朗日插值多项式;分段线性插值;牛顿前插;曲线拟合引言:在数学物理方程中,当给定数据是不同散点时,无法确定函数表达式,求解函数就需要很大的计算量,我们有多种方法对给定的表格函数进行求解,我们这里,利用插值法和曲线拟合对函数进行求解,进一步了解函数性质,两种方法各有利弊,适合我们进行不同的散点函数求解。
正文:一、插值法和分段线性插值1拉格朗日多项式原理对某个多项式函数,已知有给定的k + 1个取值点:其中对应着自变量的位置,而对应着函数在这个位置的取值。
假设任意两个不同的x j都互不相同,那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为:其中每个为拉格朗日基本多项式(或称插值基函数),其表达式为:[3]拉格朗日基本多项式的特点是在上取值为1,在其它的点上取值为0。
2分段线性插值原理给定区间[a,b], 将其分割成a=x0 <x1<…<xn=b, 已知函数y= f(x) 在这些插值结点的函数值为y k =f(xk)(k=0,1,…,n)求一个分段函数Ih(x), 使其满足:(1) Ih (xk)=yk,(k=0,1,…,n) ;(2) 在每个区间[xk ,xk+1] 上,Ih(x)是个一次函数。
易知,I h (x)是个折线函数, 在每个区间[x k ,x k+1 ]上,(k=0,1,…,n)k 1k k1k 1k k 1k k k ,1)()()(x x x x x f x x x x x f x L --+--=++++,于是, I h (x)在[a,b]上是连续的,但其一阶导数是不连续的。
3拉格朗日插值多项式算法○1输入,(0,1,2,,)i i x y i n =L ,令0)(=x L n 。
○2对0,1,2,,i n =L ,计算 0,()()/()ni j i j j j il x x x x x -≠=--∏()()()n n i iL x L x l x y ←−−+4分段线性插值算法○1输入(x k ,y k ),k=0,1,…,n;○2计算k 1k k1k 1k k 1k k k ,1)()()(x x x x x f x x x x x f x L --+--=++++5插值法和分段线性插值程序按下列数据分别作五次插值和分段线性插值,画出两条插值曲线以及给定数据点。
求x 1=0.32, x 2=0.55, x 3=0.68 时的函数近似值,并比较两种方法的插值余function lagrintxi=[0.32,0.55,0.68]; %xi=[0.2:0.001:0.8];x=[0.3,0.42,0.50,0.58,0.66,0.72];y=[1.04403,1.08462,1.11803,1.15603,1.19817,1.23223]; L=zeros(size(y)); m=length(xi); for i=1:mdxi=xi(i)-x;L(1)=prod(dxi(2:6))/prod(x(1)-x(2:6)); L(6)=prod(dxi(1:6-1))/prod(x(6)-x(1:6-1)); for j=2:6-1num=prod(dxi(1:j-1))*prod(dxi(j+1:6));den=prod(x(j)-x(1:j-1))*prod(x(j)-x(j+1:6));L(j)=num/den;endyi(i)=sum(y.*L);fprintf('x=%f,y=%f\n',xi(i),yi(i));endplot(xi,yi,'r');axis([0.2 0.8 1.03 1.24]);hold onplot(x,y,'b.','markersize',20)grid on分段线性插值算法程序:function [y]=div%xi=[0.3:0.001:0.72];x0=[0.3,0.42,0.50,0.58,0.66,0.72];y0=[1.04403,1.08462,1.11803,1.15603,1.19817,1.23223]; k=1;xi=[0.32,0.55,0.68];for j=1:3for i=1:5if xi(j)>=x0(i) && xi(j)<=x0(i+1) && k<=3lx(1)=(xi(j)-x0(i+1))/(x0(i)-x0(i+1));lx(2)=(xi(j)-x0(i))/(x0(i+1)-x0(i));y(k)=lx(1)*y0(i)+lx(2)*y0(i+1);k=k+1;endendendplot(xi,y,'r');axis([0.2 0.8 1.03 1.24]);hold onplot(x0,y0,'b.','markersize',20)grid on6运算结果拉格朗日插值结果x=0.320000,y=1.049958 x=0.550000,y=1.141271 x=0.680000,y=1.209300 拉格朗日插值余项:)72.0)(66.0)(58.0)(5.0)(42.0)(3.0(!6)()()()(655------=-=x x x x x x f x L x f x R ξ分段插值结果ans =1.0508 1.1418 1.2095 分段线性插值余项:)(max 8))((2)(max)()(max 11121,1x f h x x x x f x L x f R i i i i i i x x x i i i i x x x i x x x ''≤--''≤-=+++≤≤+≤≤≤≤ξ由于拉格朗日插值的余项比分段线性插值的余项要求更为严格,点少、区间小的时候,拉格朗日插值要更好。
但在区间较大、节点较多的时候,分段线性插值要更好。
二、牛顿前插 1牛顿前插原理n 次牛顿前插公式:)1()1(!)1(!2!1)(!)(0020010000+--∆++-∆+∆+=-∏∆=+-==∑n t t t n f t t f t f f j t k f th x N n k j nk k n ΛΛ插值余项:)()!1()()1()()1(10ξ+++--=+n n n f h n n t t t th x R Λ,),(0n x x ∈ξm 阶差分记作i m i m i m f f f 111-+-∆-∆=∆。
[]i f x k 阶差商是[]∑∏=≠=-=km kmi i i mm k x xx f x x x f 0010)()(,,,Λ差分和差商之间的关系是[]k ik k i i i h k f x x x f !,,,1∆=++Λ2牛顿前插算法○1输入),2,1,0(,,n i y x n ii Λ=。
○2对n k ,,3,2,1Λ=,k i ,,2,1Λ=计算各阶差分i m i m i m f f f 111-+-∆-∆=∆ ○3计算函数值)1()1(!)1(!2!1)(!)(0020010000+--∆++-∆+∆+=-∏∆=+-==∑n t t t n f t t f t f f j t k f th x N n k j nk k n ΛΛ3牛顿前插程序:编写一个用牛顿前插公式计算函数值的程序,要求先输出差分表,再计算x 点的函数值分别求x =0.158和x =0.636的三次插值的值,并比较二者的插值余项。
这里以x=0.636为例 function [P]=newtoncha x0=0.636;X=[0.125 0.250 0.375 0.500 0.625 0.750]; Y=[0.796 0.773 0.744 0.704 0.656 0.602]; h=abs(X(2)-X(1)); n=find(abs(x0-X)<3*h);X=X(n(1):n(end));Y=Y(n(1):n(end));w=length(X);R=zeros(w,w);R(:,1)=Y(:);for k=2:wfor j=k:wR(j,k)=R(j,k-1)-R(j-1,k-1); endendt=(x0-X(1))/h;T=1;for m=1:w-1T=T*(t-m+1);N(m)=R(m+1,m+1)*T/factorial(m); endP=R(1,1)+sum(N);4运行结果:差分表X=0.636时ans =0.651459661824000 x=0,158时ans =0.790229818880000 三、曲线拟合 1曲线拟合原理: 给定数据(,),1,2,,j j x y j n=L 。
记拟合函数的形式为0011()()()()m m p x a x a x a x ϕϕϕ=+++L (1.1),其中{}0()mk k x ϕ=为已知的线性无关函数。
求系数***01,,,ma a a L 使得2201110(,,,)[()][()]nnmm j j k k j j j j k a a a p x y a x y ϕϕ====-=-∑∑∑L (1.2)取最小值。
称**0()()mk k k p x a x ϕ==∑(1.3)为拟合函数或经验公式。
如果),,1,0()(m k x x k k Λ==ϕ,则(1.3)为m 次最小二乘拟合多项式2曲线拟合算法:已知数据对),,2,1)(,(n j y x i i Λ=,求多项式)()(0n m x a x P mi i i <=∑=,使得∑∑==-=Φnj j mi i j i m y x a a a a 1210)(),,,(Λ为最小。
注意到此时kk x x =)(ϕ,多项式系数ma a a ,,,10Λ满足下面的线性方程组:⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡++m m m m m m m T T T a a a S S S S S SS S S M M ΛM M M ΛΛ10102112110其中)2,,2,1,0(1m k xS nj k jk Λ==∑=,),,2,1,0(1m k x y T nj k jj k Λ==∑=然后只要调用线性方程组的函数程序即可 3曲线拟合程序:试分别用抛物线y =a +bx 2和指数曲线y =ae bx 拟合下列数据画出数据点和两条拟合曲线,并通过计算2个拟合函数残差向量的2范数来比较拟合优劣。