《项目二 数值计算》由会员分享，可在线阅读，更多相关《项目二 数值计算（9页珍藏版）》请在金锄头文库上搜索。
1、项目二 数值计算实验类型 验证型 实验学时 8实验目的 学习Matlab软件的数值计算功能及操作实验原理 Matlab软件的数值计算基本功能及操作实验内容1. 数据插值插值法是实用的数值方法，是函数逼近的重要方法。在生产和科学实验中，自变量x与因变量y的函数y = f(x)的关系式有时不能直接写出表达式，而只能得到函数在若干个点的函数值或导数值。当要求知道观测点之外的函数值时，需要估计函数值在该点的值。如何根据观测点的值，构造一个比较简单的函数y=(x)，使函数在观测点的值等于已知的数值或导数值。用简单函数y=(x)在点x处的值来估计未知函数y=f(x)在x点的值。寻找这样的函数(x)，办法是
2、很多的。(x)可以是一个代数多项式，或是三角多项式，也可以是有理分式；(x)可以是任意光滑（任意阶导数连续）的函数或是分段函数。函数类的不同，自然地有不同的逼近效果。在许多应用中，通常要用一个解析函数（一、二元函数）来描述观测数据。根据测量数据的类型：1测量值是准确的，没有误差。2测量值与真实值有误差。这时对应地有两种处理观测数据方法：1插值或曲线拟合。2回归分析（假定数据测量是精确时，一般用插值法，否则用曲线拟合）。MATLAB中提供了众多的数据处理命令。有插值命令，有拟合命令，有查表命令。1.1 插值命令命令1 interp1功能 一维数据插值（表格查找）。该命令对数据点之间计算内插值。它
3、找出一元函数f(x)在中间点的数值。格式 yi = interp1(x,Y,xi) %返回插值向量yi，每一元素对应于参量xi，同时由向量x与Y的内插值决定。参量x指定数据Y的点。若Y为一矩阵，则按Y的每列计算。yi是阶数为length(xi)*size(Y,2)的输出矩阵。yi = interp1(Y,xi) %假定x=1:N，其中N为向量Y的长度，或者为矩阵Y的行数。yi = interp1(x,Y,xi,method) %用指定的算法计算插值：nearest：最近邻点插值，直接完成计算；linear：线性插值（缺省方式），直接完成计算；spline：三次样条函数插值。对于该方法，命令in
4、terp1调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline用它们执行三次样条函数插值；pchip：分段三次Hermite插值。对于该方法，命令interp1调用函数pchip，用于对向量x与y执行分段三次内插值。该方法保留单调性与数据的外形；cubic：与pchip操作相同；v5cubic：在MATLAB 5.0中的三次插值。例1.1 已知某产品从1900至2010每间隔10年的产量，试用不同的插值方法计算其在1995年的产量，并对不同的插值方法得到的结果进行比较。year = 1900:10:2010;product = 7
5、5.995 91.972 105.711 123.203 131.669 150.697 179.323 203.212 226.505 249.633 256.344 267.893 ;p19951 = interp1(year,product,1995,pchip);%分段三次插值p19952 = interp1(year,product,1995,nearest);%最邻近插值p19953 = interp1(year,product,1995,linear)；%线形插值p19954 = interp1(year,product,1995,spline)；%三次样条插值x = 1900:
6、1:2010;y1 = interp1(year,product,x,pchip); y2 = interp1(year,product,x, nearest); y3 = interp1(year,product,x, linear); y4 = interp1(year,product,x, spline);subplot(2,2,1);plot(year,product,o,x,y1)subplot(2,2,2);plot(year,product,o,x,y2)subplot(2,2,3);plot(year,product,o,x,y3)subplot(2,2,4);plot(yea
7、r,product,o,x,y4)结果：p19951 = 75.9950； p19952 = 75.9950； p19953 = 75.9950；p19954 = 75.9950命令2 interp2功能 二维数据内插值（表格查找）格式 ZI = interp2(X,Y,Z,XI,YI) %返回矩阵ZI，其元素包含对应于参量XI与YI（可以是向量、或同型矩阵）的元素，即Zi(i,j)Xi(i,j),yi(i,j)。用户可以输入行向量和列向量Xi与Yi，此时，输出向量Zi与矩阵meshgrid(xi,yi)是同型的。同时取决于由输入矩阵X、Y与Z确定的二维函数Z=f(X,Y)。参量X与Y必须是单
8、调的，且相同的划分格式，就像由命令meshgrid生成的一样。若Xi与Yi中有在X与Y范围之外的点，则相应地返回nan（Not a Number）。ZI = interp2(Z,XI,YI) %缺省地，X=1:n、Y=1:m，其中m,n=size(Z)。再按第一种情形进行计算。ZI = interp2(Z,n) %作n次递归计算，在Z的每两个元素之间插入它们的二维插值，这样，Z的阶数将不断增加。interp2(Z)等价于interp2(z,1)。ZI = interp2(X,Y,Z,XI,YI,method) %用指定的算法method计算二维插值：linear：双线性插值算法（缺省算法）；n
9、earest：最临近插值；spline：三次样条插值；cubic：双三次插值。例2 已知1950到1990每间隔10年，工龄10年，20年，30年的劳动报酬，试用不同的二维插值方法计算1975年，工龄15年的平均工资，并对结果进行比较。years = 1950:10:1990;service = 10:10:30;wage = 150.697 199.592 187.625 179.323 195.072 250.287 203.212 179.092 322.767 226.505 153.706 426.730 249.633 120.281 598.243;w = interp2(ser
10、vice,years,wage,15,1975)2. 数据拟合问题：给定一批数据点，需确定满足特定要求的曲线或曲面若要求所求曲线（面）通过所给所有数据点，就是插值问题；若不要求曲线（面）通过所有数据点，而是要求它反映对象整体的变化趋势，这就是数据拟合，又称曲线拟合或曲面拟合函数插值与曲线拟合都是要根据一组数据构造一个函数作为近似，由于近似的要求不同，二者在数学方法上是完全不同的2.1多项式拟合：命令：p=polyfit(X,Y,n)给出已知数据（X，Y）的n次多项式拟合，输出为多项式系数（按降幂排列）命令：polyval(p,x0)利用上面拟合出的多项式p,计算在x0处的值，x0可以是单个点，
11、也可以是数组。例2.1 求如下数据的拟合曲线：X=0.5:0.5:3;Y=1.75 2.45 3.81 4.80, 7.00 8.60;P=polyfit(X,Y,2)Y1=polyval(P,X);f=inline(c(1)*exp(c(2)*x),c,x);c=lsqcurvefit(f,0,0,X,Y);sum(Y-Y1).2)sum(Y-f(c,X).2)P =0.5614 0.8287 1.1560c = 1.4632 0.60030.1781 0.40432.2 plot(X,Y,+,X,Y1,X,f(c,X),*)非线性拟合命令：c=lsqcurvefit(fun,c0,x,y)
12、说明：fun是预先定义的一个函数文件（M文件），c0是初值。例xdata=-0.2 0 0.1 0.2 0.15 0.3;ydata=1.50 1.06 0.95 0.84 0.86 0.72;f=inline(c(1)*exp(c(2)*x),c,x);c=lsqcurvefit(f,0,0,xdata,ydata);plot(xdata,ydata,*,xdata,f(c,xdata)sum(f(c,xdata)-ydata).2)例3根据美国人口从1790年到1990年间的人口数据（如下表），确定人口指数增长模型（Logistic模型）中的待定参数，估计出美国2010年的人口，同时画出拟
13、合效果的图形。表1 美国人口统计数据年 份1790180018101820183018401850人口(106)3.95.37.29.612.917.123.2年 份1860187018801890190019101920人口(106)31.438.650.262.976.092.0106.5年 份193019401950196019701980人口(106)123.2131.7150.7179.3204.0226.5提示：可采用以下两种模型指数增长模型：T=1790:10:1980;X=3.9 5.3 7.2 9.6 12.9 17.1 23.2 31.4 38.6 50.2 62.9 76
14、.0 92.0 106.5 123.2 131.7 150.7 179.3 204.0 226.5;f=inline(r(1)*exp(r(2)*T),r,T);r=lsqcurvefit(f,0,0,T,X)plot(T,X,*,T,f(r,T)norm(feval(f,r,T)-X)2Logistic模型：T=1790:10:1980;X=3.9 5.3 7.2 9.6 12.9 17.1 23.2 31.4 38.6 50.2 62.9 76.0 92.0 106.5 123.2 131.7 150.7 179.3 204.0 226.5;f=inline(r(1)/(1+(r(1)/r(2)-1)*exp(-r(3)*T),r,T);r=lsqcurvefit(f,227,3.9,0,T,X)plot(T,X,*,T,f(r,T)norm(feval(f,r,T)-X)2可参考拟合函数：a=lsqcurvefit(example_curvefit_fun,a0,x,y);3.方程及方程组的求解命令：solve，fsolve格式：solve的基本功能是求解符合方程（组）SOLVE(eqn1,eqn2
