第5章Matlab数值计算20141030

第五章MATLAB数值计算（2）1插值与逼近2常微分方程数值求解3非线性方程数值求解插值与逼近•插值：在实际问题中常遇到这样的函数y=f(x)其在某个区间[a,b]上是存在的，但通过观测或测量或实验只能的到在这个区间[a,b]上有限个离散点x0,x1,…,xn上的函数值y=f(xi)（i=0,1,…,n）或f(x)的函数表达式是已知的，但确很复杂而不便于计算，希望用一个简单的函数来描述它。3•假设函数y=f(x)是[a,b]上的实值函数，x0,x1,…,xn是[a,b]上n+1个互异的点，f(x)在这些点上的取值分别为y0,y1,…,yn,求一个确定的函数P(x),使之满足：P(xi)=yi(i=0,1,2,…,n)(10-5)称x0,x1,…,xn为插值节点，关系式(1)称为插值原则，函数P(x)称为函数y=f(x)的插值函数，区间[a,b]称为插值区间。插值法的概念:要求该多项式的函数曲线要经过y=f(x)上已知的这n+1个点，(x0,y0),(x1,y1),…(xn,yn)。同时在其它上要估计误差R(x)=f(x)-P(x)。],[bax45•解：设，其中为待定系数。•由给定的条件有:••联立上式得：•于是，6线性插值线性插值也叫两点插值，已知函数y=f(x)在给定互异点x0,x1上的值为y0=f(x0)，y1=f(x1)线性插值就是构造一个一次多项式P1(x)=ax+b使它满足条件P1(x0)=y0P1(x1)=y1（10-6）其几何解释就是一条直线，通过已知点A(x0,y0)，B(x1,y1)。)()(0010101xxxxyyyxP（点斜式）101001011)(yxxxxyxxxxxP（对称式）P1(x)就是所求的一次多项多，称为f(x)的线性插值多项式。7例：用线性插值求（x*=10.723805）•解：设，取x0=100，x1=121•则y0=10y1=11•从而•Matlab(实现)yi=interp1([100121],[1011],115)8逼近•用简单的函数p(x)近似地代替函数f(x)，是计算数学中最基本的概念和方法之一。近似代替又称为逼近，函数f(x)称为被逼近的函数，p(x)称为逼近函数，两者之差称为逼近的误差或余项。9•一致逼近以函数f(x)和p(x)的最大误差作为度量误差f(x)-p(x)的“大小”的标准，在这种意义下的函数逼近称为一致逼近或均匀逼近，讲得更具体一点，也即对于任意给定的一个小正数0，如果存在函数p(x)，使不等式成立，则称该函数p(x)在区间[a,b]上一致逼近或均匀逼近于函数f(x)。10•平方逼近：如果采用作为度量误差的“大小”的标准，在这种意义下的函数逼近称为平方逼近或均方逼近。这种方法要比一致逼近的相应问题简单得多。11matlab数值插值1.一维数值插值interp1函数调用格式为：Y1=interp1(X,Y,X1,'method')函数根据X、Y的值，计算函数在X1处的值。X、Y是两个等长的已知向量，分别描述采样点和样本值，X1是一个向量或标量，描述欲插值的点，Y1是一个与X1等长的插值结果。method是插值方法，允许的取值有'linear'(线性插值)、'nearest'(最近插值)、'spline'(三次样条插值)、'cubic'（三次多项式插值），缺省值是'linear'。例用不同的插值方法计算sin(x)在π/2点的值。这是一个一维插值问题。在MATLAB命令窗口，输入命令：X=0:0.2:pi;Y=sin(X);%给出X、Yinterp1(X,Y,pi/2)%用缺省方法(即线性插值方法)计算sin(π/2)interp1(X,Y,pi/2,'nearest')%用最近方法计算sin(π/2)interp1(X,Y,pi/2,'linear')%用线性方法计算sin(π/2)interp1(X,Y,pi/2,'spline')%用三次样条方法计算sin(π/2)interp1(X,Y,pi/2,'cubic')%用三次多项式方法计算sin(π/2)MATLAB中有一个专门的三次样条插值函数Y1=spline(X,Y,X1)，其功能及使用方法与函数Y1=interp1(X,Y,X1,'spline')完全相同。11.522.533.540.99750.9980.99850.9990.999511-nearest,2-linear,3-spline,4-cubic例已知检测参数f随时间t的采样结果，用数值插值法计算t=2，7，12，17，22，17，32，37，42，47，52，57时f的值。这是一个一维数值插值问题，命令如下：T=0:5:65;X=2:5:57;F=[3.2015,2.2560,879.5,1835.9,2968.8,4136.2,5237.9,6152.7,...6725.3,6848.3,6403.5,6824.7,7328.5,7857.6];F1=interp1(T,F,X)%用线性方法插值F1=interp1(T,F,X,'nearest')%用最近方法插值F1=interp1(T,F,X,'spline')%用三次样条方法插值F1=interp1(T,F,X,'cubic')%用三次多项式方法插值2.二维数值插值MATLAB中，提供了解决二维插值问题的函数。其调用格式为：Z1=interp2(X,Y,Z,X1,Y1,'method')其中X、Y是两个向量，分别描述两个参数的采样点，Z是与参数采样点对应的采样变量的样本值，X1、Y1是两个向量或标量，描述欲插值的点。method的取值与一维插值函数相同。例设Z=x2+y2，对Z函数在(0，1)×(0，2)区域内进行插值。命令如下：x=0:0.1:10;y=0:0.2:20;[X,Y]=meshgrid(x,y);Z=X.^2+Y.^2;interp2(x,y,Z,0.5,0.5)%对函数在(0.5,0.5)点进行插值interp2(x,y,Z,[0.50.6],0.4)%对函数在(0.5,0.4)点和(0.6，0.4)点进行插值interp2(x,y,Z,[0.50.6],[0.40.5])%对函数在(0.5,0.4)点和(0.6，0.5)点进行插值interp2(x,y,Z,[0.50.6]‘,[0.40.5])%对函数在(0.5,0.4),(0.6,0.4),(0.5,0.5)和(0.6,0.5)点进行插值3.三维数值插值对三维函数插值的函数是interp3，其使用方法和interp2相同。其调用格式为：W1=interp3(X,Y,Z,W,X1,Y1,Z1,'method')函数返回三维插值结果。其中X、Y、Z是三个向量，分别描述三个参数的采样点，W是与参数采样点对应的采样变量的样本值，X1、Y1、Z1是三个向量或标量，描述欲插值的点。method是插值方法，可选，其缺省值是‘line'。method的取值与一、二维插值函数相同。数据插值•例：已知某处山区地形选点测量坐标数据为：–x=00.511.522.533.544.55–y=00.511.522.533.544.555.56•海拔高度数据为：–z=8990878592919693908782–9296989995918986848284–9698959290888584838185–8081828995969392898686–8285879899969788858283–8285899495939291868488–8892939495898786838192–9296979896939584828184–8585818280808185909395–8486819899989796958487–8081858283848790958688–8082818485868382818082–8788899899979698949287•求对数据插值加密后形成的地貌图数据插值•绘原始数据图–x=0:.5:5;–y=0:.5:6;–z=[8990878592919693908782–9296989995918986848284–9698959290888584838185–8081828995969392898686–8285879899969788858283–8285899495939291868488–8892939495898786838192–9296979896939584828184–8585818280808185909395–8486819899989796958487–8081858283848790958688–8082818485868382818082–8788899899979698949287];–mesh(x,y,z)012345024680859095100数据插值•绘插值后的地貌图–xi=linspace(0,5,50);%加密横坐标数据到50个–yi=linspace(0,6,80);%加密纵坐标数据到60个–[xii,yii]=meshgrid(xi,yi);%生成网格数据–zii=interp2(x,y,z,xii,yii,'cubic');%插值–mesh(xii,yii,zii)%加密后的地貌图01234502467580859095100曲线拟合(一)多项式曲线拟合MATLAB中，提供了解决使用最小二乘法进行曲线拟合的函数。调用格式为：[P,S]=polyfit(X,Y,m)函数根据采样点X和采样点函数值Y，产生一个m次多项式P及其在采样点的误差向量S。其中X、Y是两个等长的向量，P是一个长度为m+1的向量。P(1)*x^m+P(2)*x^(m-1)+……+P(m)*x+p(m+1)例用一个5次多项式在区间[0，2π]内逼近函数sin(x)。命令如下：X=linspace(0,2*pi,50);Y=sin(X);[P,S]=polyfit(X,Y,5)%得到5次多项式的系数和误差plot(X,Y,'k-*',X,polyval(P,X),'r-')01234567-1.5-1-0.500.511.5真值拟合值使用指定函数进行曲线拟合x=lsqnonlin(fun,x0)x返回拟合参数，fun为拟合曲线对应的函数，x0为初始值例：根据所给试验数据，用双曲线模型拟合，确定模型参数0xyyabx0,0.05,0.1,0.2,0.4,0.6,0.8,1.21.335,1.253,1.180,1.058,0.887,0.803,0.752,0.685xyfunction[]=curvefit()data=[01.335;0.051.253;0.11.180;0.21.058;0.40.887;0.60.803;0.80.752;1.20.685];x=data(:,1);y=data(:,2);plot(x,y,'g-');holdon;x0=[10]';p=lsqnonlin(@curve_model,x0);y1=y(1);y=y1-x./(p(1)+p(2)*x);plot(x,y,'bo');%定义curve_model函数functionf=curve_model(p,data)data=[01.335;0.051.253;0.11.180;0.21.058;0.40.887;0.60.803;0.80.752;1.20.685];x=data(:,1);y=data(:,2);y1=y(1);z=y1-x./(p(1)+p(2)*x);f=z-y;函数零点[x,feval,exitflag,output]=fzero(fun,x0,options)x为函数fun在x0附近的零点，x0为x的初始值，其他符号意义与一元函数极值时相同。例求函数的极值，初始值为。2()56fxxx00.0x[xx,fx,flag]=fzero(‘x.*x-5*x+6’,0.0)%一个零点[xx,fx,flag]=fzero(‘x.*