MATLAB(第四章和第八章)

一、实验内容：

题目1.

Jordan矩阵是矩阵分析中一类很实用的矩阵，其一般形式为

-0J= 010010000510051000,例如J1= 00510 0005100050试利用diag()函数给出构造J1的语句。

【分析】该题为对角矩阵的问题。对J要利用diag()能够构造对角矩阵和次对角矩阵的性

质。J1只需对角矩阵和次对角矩阵相加即可。这里需要对diag()函数的调用。如A=diag(V)---已知向量生成对角矩阵，A=diag(V,k)—生成主对角线上第k条对角线为V的矩阵（其中k可为正负）

【解答】：

输入如下语句：

>>J1=diag([-5 -5 -5 -5 -5])+diag([1 1 1 1],1) 按ENTER键，显示如下： J1=

-5 1 0 0 0 0 -5 1 0 0 0 0 -5 1 0 0 0 0 -5 1 0 0 0 0 -5

题目5.

a44b试求出Vandermonde矩阵A=c44de4a3b3c3d3e3a2b2c2d2e2abcde11并以最简的形式显示结果。 1,的行列式，

11【求解】该问题有两个知识点。一个构造是Vandermonde矩阵，另一个是求矩阵的行列式。

前者可以利用书中编写的面向符号矩阵的vander()函数构造出Vandermonde矩阵。需要用到V=vander(C)来调用。后者可以用MATLAB的det()函数来求解，他会自动采用解析解法求出其行列式的值。需要注意运用det()的前提是符号矩阵，本题中A已是符号矩阵，所以不用转换。最后，用simple()函数简化一下即可。

【解答】：

（1）构造矩阵：

输入如下语句：

>>syms a b c d e; A=vander([a b c d e])

按ENTER键，显示如下： A=

[ a^4, a^3, a^2, a, 1] [ b^4, b^3, b^2, b, 1] [ c^4, c^3, c^2, c, 1] [ d^4, d^3, d^2, d, 1] [ e^4, e^3, e^2, e, 1]

（2）以最简单的形式输出行列式： 输入如下语句：

>>det(A),simple(ans) 按ENTER键，显示如下： ans=

(c-d)*(b-d)*(b-c)*(a-d)*(a-c)*(a-b)*(-d+e)*(e-c)*(e-b)*(e-a)

77

15. 试求出线性代数方程组X 

26

61149352

7

21012=，并验证解的正确性 503126

-1，

【分析】:该题为线性方程的计算机求解问题。 需要考虑的是X=B*A在MATLAB中，需要

调用inv(A)*B函数，来得出方程的解。同时需要用到逆运算。

【解答】：

(1)输入如下语句：

>> A=[7,6,9,7;7,1,3,2;2,1,5,5;6,4,2,6];B=[2,1,0,1;0,3,1,2];A=A',B=B'; x=inv(A)*B,e1=norm(A*x-B),x1=inv(sym(A))*B,e2=norm(double(A*x1-B)) 语句运行后，显示如下： x =

-0.0057 0.4511 0.1034 -0.6207 -0.1609 -0.3678 0.2730 0.3204 e1 =

1.5879e-015 x1 =

[ -1/174, 157/348] [ 3/29, -18/29] [ -14/87, -32/87] [ 95/348, 223/696] e2 = 0

（2）对X进行逆运算，输入以下语句：

>> x1=x1'； x1*A'

语句运行后，显示如下：

ans =

[ 2, 1, 0, 1] [ 0, 3, 1, 2]

二、实验心得

这次是第三次上高等应用数学问题的MATLAB求解课程。通过老师上课的细心讲解与演示，我对MATLAB又有了更深的了解。原来MATLAB在线性代数问题矩阵问题中也可以用的如此灵活简便。

同时我还学到了很多MATLAB的应用。首先是矩阵的输入，我学会了如何用简单的函数语句直接输入零矩阵，幺矩阵，随机元素矩阵，对角元素矩阵，Hankel矩阵，伴随矩阵等。其调用的语句虽然看似简单，但还是要注意细节。就拿对角元素矩阵说，调用语句A=diag(V,k)是生成主对角线上第k条对角线为V的矩阵，这里要深刻理解k的含义，他可正可负，是对角线上第k条对角线。同时我们还要灵活掌握向量与矩阵的两两转化，不是只掌握一个就可以了。

我还学到了矩阵分析的基本概念及求解函数，比如如何做出行列式，迹，秩，范数，特征多项式，逆矩阵，广义逆矩阵，特征值与特征向量等。以行列式为例，我会运用简单的d=det(A)函数来调用，直接求行列式。但是要注意细节的是我们练习的题目中，是Hilbert矩阵，我们要用sym()函数把他符号化，再用det()调用。所以不符号化的话，就会遇到运行错误的问题。

我还会运用了矩阵的基本变换，知道了矩阵的正交分解，三角分解，对称矩阵的Cholesky分解，一半矩阵的伴随分解，Jordan变换等。我了解了线性代数方程对唯一解、无穷解、无解等问题的处理，矩阵方程的据算机求解等。就拿线性方程组的计算机求解来说，已AX=B为例，其实很简单，只要输入A,B矩阵，再运用X=inv(A)*B函数就可以编写出了。

第八章 数据插值、函数逼近问题的计算机求解 一、实验内容：

1 用y(t)te25tsint生成一组较稀疏的数据，并用一维数据插值的方法对给出的数据进

行曲线拟合，并将结果与理论曲线相比较。 【分析】：该题为一维差值问题。该题可以利用函数interp1(),其调用格式为

y=interp1(x,y,x1,f方法)。插值方法有很多，一般默认为可以选择’linear’’线性插值，‘nearest’对近点等值方法,’cubic’三次Hermite插值’spline’ 三次分段样条插值。本题运用最精确的三次分段样条插值进行求解，并用理论曲线与拟合的直线共同显示在一起进行比较。

【解答】：

(1) 根据所给函数用MATLAB生成较稀疏的数据，输入如下语句： >>t=0:0.2:2;

y=t.^2.*exp(-5*t).*sin(t); plot(t,y,'o') 显示图片如下：

(2) 调用一维插值函数interp1(),并用最精确的插样方法‘spline’(三次分段样条插值)。输

入如下语句:

>>ezplot('t.^2.*exp(-5*t).*sin(t)',[0,2]); hold on x1=0:0.01:2; y1=interp1(t,y,x1,'spline'); plot(x1,y1,'g')

如图所示，曲线拟合与原曲线完全一样。 还可以用四种方式求解，输入如下数据： >> t1=0:0.02:1;

y0=t1.^2.*exp(-5*t1).*sin(t1);

y1=interp1(t,y,t1);y2=interp1(t,y,t1,'nearest'); t1=0:0.02:1;

y0=t1.^2.*exp(-5*t1).*sin(t1); y1=interp1(t,y,t1);

y2=interp1(t,y,t1,'cubic');

y3=interp1(t,y,t1,'spline'); y4=interp1(t,y,t1,'nearest');

plot(t1,[y1',y2',y3',y4'],':',t,y,'o',t1,y0)

e1=max(abs(y0(1:49)-y2(1:49))),e2=max(abs(y0-y3)),e3=max(abs(y0-y4)) 按ENTER键，输出如下： e1 =

1.4195e-004 e2 =

1.0990e-004 e3 =

0.0018

并显示图片：

比较e,由结果可知三次分段样条插值的精确度最好。

2 用y(t)sin(10t3)在(0; 3) 区间内生成一组较稀疏的数据，并用一维数据插值的方法对给出的数据进行曲线拟合，并将结果与理论曲线相比较

【分析】：该题目和第一题一样，同样先生成数据，然后再用spline’(三次分段样条插值)函

数进行拟合。该问题需要考虑步长的取值，步长的大小取决了其精确度。

【解答】：

(1) 用MATLAB生成一组较稀疏的数据，输入如下语句： >>0:0.2:3;

y=sin(10*t.^2+3); plot(t,y,'o') 图片显示如下：

2

（2）画出理论曲线，拟合一维数据插值，进行比较,输入如下语句： >>ezplot('sin(10*t^2+3)',[0,3]) plot(t,y,'r')

ezplot('sin(10*t^2+3)',[0,3]); hold on x1=0:0.001:3; y1=interp1(t,y,x1,'spline'); plot(x1,y1,'r')

由图可知拟合的曲面一开始与原函数吻合的很好，但是却到后面越来越偏离，所以需要再减小样本点的步长，输入如下语句：

>> t=[0:0.1:1,1.1:0.01:3]; y=sin(10*t.^2+3); plot(t,y,'o') ezplot('sin(10*t^2+3)',[0,3]); hold on x1=0:0.001:3; y1=interp1(t,y,x1,'spline'); plot(x1,y1,'r') 显示图片如下：

即与原函数完全吻合。

该题还可以用四种插值方法进行求解，输入如下数据： >> t1=0:0.01:3;

>> y0=sin(10.*(t1.^2)+3); >> y1=interp1(t,y,t1);

>> y2=interp1(t,y,t1,'cubic'); >> y3=interp1(t,y,t1,'spline'); >> y4=interp1(t,y,t1,'nearest');

>> plot(t1,[y1',y2',y3',y4'],':',t,y,'o',t1,y0)

>> e1=max(abs(y0(1:49)-y2(1:49))),e2=max(abs(y0-y3)),e3=max(abs(y0-y4)) 按ENTER键，显示如下： e1 =

0.0023 e2 =

0.0282 e3 =

0.5748

由图可知比较并比较e,由结果可知三次分段样条插值的精确度最好。 题目三：用f(x,y)1x2y4esin(xy2x2y)原型函数生成一组网格数据或随机32xy数据，分别拟合出曲面，并和原曲面进行比较。

【分析】：该题考了两个方法。一个是生成一组网格数据，拟合曲面。另一个是生成随机

数据，拟合曲面。前者应先绘制已知数据的网格图，再采用'spline'拟合曲面，最后要对比插值结果和新网格下的函数值精确解，绘制出绝对插值误差曲面。后者应先选择随机分布的样本点，绘制出已知数据点的分布和已知数据点的三维分布图。然后对产生的三维曲面数据进行’v’算法差值，最后用abs()与原曲面进行比较，进行误差分析.

【解答】:

(1)形成网格数据：

绘制已知数据的网格图，输入如下语句： >>[x,y]=meshgrid(0.2:0.2:2);

z=exp(-x.^2-y.^4).*sin(x.*y.^2+x.^2.*y)./(3*x.^3+y); surf(x,y,z)

按ENTER键，图片显示如下：

采用三次分段样条插值方法，输入如下语句： >> [x1,y1]=meshgrid(0.2:0.02:2); z1=interp2(x,y,z,x1,y1,'spline'); surf(x1,y1,z1)

按ENTER键，图片显示如下：

对三次样条插值方法进行误差分析，对比插值结果和新网格下的函数值精确解，可以绘制出绝对插值误差曲面，如输入如下语句：

>> z0=exp(-x1.^2-y1.^4).*sin(x1.*y1.^2+x1.^2.*y1)./(3*x1.^3+y1); surf(x1,y1,abs(z1-z0)) 按ENTER键，图片显示如下：

由图可知，三次样条插值方法拟合的较为精确。 (2) 形成随即数据：

选择随机分布的样本点，输入如下语句：

>> x=0.2+1.8*rand(400,1); y=0.2+1.8*rand(400,1);

z=exp(-x.^2-y.^4).*sin(x.*y.^2+x.^2.*y)./(3*x.^3+y); plot(x,y,'.')

figure, plot3(x,y,z,'v')

显示已知数据点的分布，已知数据点的三维分布图片如下：

（数据点分布）

（数据点三维分布）

对产生的三维曲面数据进行’v’算法差值，输入如下语句： >>[x1,y1]=meshgrid(0.3:0.02:1.9); z1=griddata(x,y,z,x1,y1,'v4'); surf(x1,y1,z1)

按ENTER键，图片显示如下：

与原曲面进行比较，进行误差分析，输入如下语句：

>> z0=exp(-x1.^2-y1.^4).*sin(x1.*y1.^2+x1.^2.*y1)./(3*x1.^3+y1); surf(x1,y1,abs(z1-z0))

由图可知该方法的拟合结果较为精确。但是边界值还是存在较大误差。

题目四：假设已知一组数据，试用插值方法绘制出x 2 (¡2; 4:9) 区间内的光滑函数曲线，比较各种插值算法的优劣。

xi yi xi yi -2 1.6 -1.7 -1.4 0.1318 2.2 -1.1 0.14483 2.5 -0.8 0.15656 2.8 0.09768 -0.5 0.16622 3.1 0.08353 -0.2 0.17332 3.4 0.07015 0.1 0.1775 3.7 0.4 0.17853 4 0.7 1 1.3 0.16302 4.90 0.02236 0.10289 0.11741 1.9 0.17635 0.17109 4.3 4.6 0.15255 0.1402 0.12655 0.11219 0.05786 0.04687 0.03729 0.02914 【分析】：该题可以利用’linear’’线性插值，‘nearest’对近点等值方法,’cubic’三次Hermite插值’spline’ 三次分段样条插值这四种方法本题进行求解，并用理论曲线与拟合的直线共同显

示在一起进行比较。同时计算精确度，比较e的大小即可。 【解答】：

输入如下语句，进行四种方法的求解：

>> x=[-2,-1.7,-1.4,-1.1,-0.8,-0.5,-0.2,0.1,0.4,0.7,1,1.3,... 1.6,1.9,2.2,2.5,2.8,3.1,3.4,3.7,4,4.3,4.6,4.9];

y=[0.10289,0.11741,0.13158,0.14483,0.15656,0.16622,0.17332,... 0.1775,0.17853,0.17635,0.17109,0.16302,0.15255,0.1402,... 0.12655,0.11219,0.09768,0.08353,0.07019,0.05786,0.04687,... 0.03729,0.02914,0.02236]; >> x1=-2:0.02:4.9; >> plot(x,y,x,y,'o') y1=interp1(x,y,x1);

y2=interp1(x,y,x1,'cubic'); y3=interp1(x,y,x1,'spline'); y4=interp1(x,y,x1,'nearest');

plot(x1,[y1',y2',y3',y4'],':',x,y,'o',x,y) 显示图片如下：

输入如下数据，比较精确度：

>> e1=max(abs(y0(1:49)-y2(1:49))),e2=max(abs(y0-y3)),e3=max(abs(y0-y4)) 显示如下： e1 =

0.0039 e2 =

0.0274 e3 =

0.5479

由图可知三次分段样条插值的吻合度最好。比较e的大小，三次分段样条插值的e值最小。

题目五：

假设已知实测数据由下表给出，试对(x; y) 在(0:1; 0:1) ，(1:1; 1:1) 区域内的点进行插值，并用三维曲面的方式绘制出插值结果

【分析】：该题应采用最精确的三次分段样条插值方法‘spline‘进行拟合，然后可以利

用axis([0.1,1.1,0.1,1.1,min(),max()])函数进行制图。若该题不考虑插值数值，则可以直接采用MATLAB 下的shadinginterp 命令来实现。

【解答】：

（1)直接采用最精确的三次分段样条插值方法，得出插值曲面输入如下语句： >>[x,y]=meshgrid(0.1:0.1:1.1);

z=[0.83041,0.82727,0.82406,0.82098,0.81824,0.8161,0.81481,0.81463,0.81579,0.81853,0.82304;0.83172,0.83249,0.83584,0.84201,0.85125,0.86376,0.87975,0.89935,0.92263,0.94959,0.9801;0.83587,0.84345,0.85631,0.87466,0.89867,0.9284,0.96377,1.0045,1.0502,1.1,1.1529;0.84286,0.86013,0.88537,0.91865,0.95985,1.0086,1.0642,1.1253,1.1904,1.257,1.3222;0.85268,0.88251,0.92286,0.97346,1.0336,1.1019,1.1764,1.254,1.3308,1.4017,1.4605;0.86532,0.91049,0.96847,1.0383,1.118,1.2046,1.2937,1.3793,1.4539,1.5086,1.5335;0.88078,0.94396,1.0217,1.1118,1.2102,1.311,1.4063,1.4859,1.5377,1.5484,1.5052;0.89904,0.98276,1.082,1.1922,1.3061,1.4138,1.5021,1.5555,1.5573,1.4915,1.346;0.92006,1.0266,1.1482,1.2768,1.4005,1.5034,1.5661,1.5678,1.4889,1.3156,1.0454;0.94381,1.0752,1.2191,1.3624,1.4866,1.5684,1.5821,1.5032,1.315,1.0155,0.62477;0.97023,1.1279,1.2929,1.4448,1.5564,1.5964,1.5341,1.3473,1.0321,0.61268,0.14763]; [x1,y1]=meshgrid(0.1:0.02:1.1); z1=interp2(x,y,z,x1,y1,'spline'); surf(x1,y1,z1)

axis([0.1,1.1,0.1,1.1,min(z1(:)),max(z1(:))]) 绘制三维曲面如下显示：

（2）该题还可以用MATLAB 下的shadinginterp 命令来实现，输入如下语句： >> surf(x,y,z); shading interp 按ENTER键，显示图片如下：

由图可知，该插值曲面更光滑，这样的插值方法更好，但是该方法没有像样条差值一样能很好地吻合数据点。

题目九：已知如下的样本点(xi; yi) 数据，试对其进行分段三次多项式样条插值。 xi yi 1 2 3 4 5 6 7 8 9 10 244.0 221.0 208.0 208.0 211.5 216.0 219.0 221.0 221.5 220.0 【分析】： 该题为对已知数据点进行三次多项式样条差值。在MATLAB中可以利用csapi()函数来定义，其调用格式为S=csapi（x,y）。其中X=[X1，X2，…Xn],

Y=[Y1,Y2,…Yn].样条函数对象的插值结果可以用fnplt()绘制出来。其条用格式为fnplt(S), Y=fnplt(S，xp).

【解答】：

对样本点进行三次多项式样条差值，输入如下语句： >>x=1:10;

y=[244.0,221.0,208.0,208.0,211.5,216.0,219.0,221.0,221.5,220.0]; S=csapi(x,y); S.coefs

按ENTER键，显示差值多项式如下： ans =

1.0e+002 *

0.01144795602886 0.01565613191343 -0.25710408794229 2.44000000000000 0.01144795602886 0.05000000000000 -0.19144795602886 2.21000000000000 -0.02723978014428 0.08434386808657 -0.05710408794229 2.08000000000000 0.00251116454827 0.00262452765373 0.02986430779801 2.08000000000000 -0.00780487804878 0.01015802129852 0.04264685675026 2.11500000000000 0.00370834764686 -0.01325661284782 0.03954826520096 2.16000000000000 -0.00202851253865 -0.00213156990725 0.02416008244589 2.19000000000000 -0.00059429749227 -0.00821710752319 0.01381140501546 2.21000000000000 -0.00059429749227 -0.01000000000000 -0.00440570250773 2.21500000000000 绘制效果图，输入如下语句：

>> fnplt(S); hold on; plot(x,y,'o') 按ENTER键，显示如下：

题目十五：

假设习题5 中数据的原型函数为z(x,y)asin(xy)bcos(yx)cxdxye试用最小二乘方法识别出a， b;，c， d， e 的数值。

【分析】：该题为最小二乘曲线拟合问题。该题首先应输入一直参数，求解a,b,c,d,e的值。

222然后再绘制拟合的曲线图。

【解答】：

（1）用最小二乘的方法求出a,b,c,d,e的值，输入如下语句： >>[x,y]=meshgrid(0.1:0.1:1.1);

z=[0.83041,0.82727,0.82406,0.82098,0.81824,0.8161,0.81481,0.81463,0.81579,0.81853,0.82304;0.83172,0.83249,0.83584,0.84201,0.85125,0.86376,0.87975,0.89935,0.92263,0.94959,0.9801;0.83587,0.84345,0.85631,0.87466,0.89867,0.9284,0.96377,1.0045,1.0502,1.1,1.1529;0.84286,0.86013,0.88537,0.91865,0.95985,1.0086,1.0642,1.1253,1.1904,1.257,1.3222;0.85268,0.88251,0.92286,0.97346,1.0336,1.1019,1.1764,1.254,1.3308,1.4017,1.4605;0.86532,0.91049,0.96847,1.0383,1.118,1.2046,1.2937,1.3793,1.4539,1.5086,1.5335;0.88078,0.94396,1.0217,1.1118,1.2102,1.311,1.4063,1.4859,1.5377,1.5484,1.5052;0.89904,0.98276,1.082,1.1922,1.3061,1.4138,1.5021,1.5555,1.5573,1.4915,1.346;0.92006,1.0266,1.1482,1.2768,1.4005,1.5034,1.5661,1.5678,1.4889,1.3156,1.0454;0.94381,1.0752,1.2191,1.3624,1.4866,1.5684,1.5821,1.5032,1.315,1.0155,0.62477;0.97023,1.1279,1.2929,1.4448,1.5564,1.5964,1.5341,1.3473,1.0321,0.61268,0.14763]; x1=x(:); y1=y(:); z1=z(:);

A=[sin(x1.^2.*y1) cos(y1.^2.*x1) x1.^2 x1.*y1 ones(size(x1))]; theta=A\\z1

显示a,b,c,d,e的数据如下所示： theta = -0.8920 3.0938 -0.1220 2.7083 -2.4251

即a,b,c,d,e分别为-0.8920， 3.0938， -0.1220， 2.7083，-2.4251。 (2)对所得数据进行拟合，输入如下语句： >>[x,y]=meshgrid(0.1:0.02:1.1);

z=theta(1)*sin(x.^2.*y)+theta(2)*cos(y.^2.*x)+theta(3)*x.^2+... theta(4)*x.*y+theta(5); surf(x,y,z) 显示图片如下:

二、实验心得

这节课我深刻地认识到了MATLAB在数据插值，函数逼近问题的应用。他不但贴近生活，而且对我们的今后的工作有很大的帮助。

我学会了一维，二维甚至多维数据插值问题的求解，知道了‘spline‘三次分段样条插值拟合是最精确的。同时还了解了基于插值技术的求取数值积分的方法和离散数据的最优化问题求解方法。我还学会了利用S=csapi(x,y)来进行三次分段多项式的插值方法，利用S=spapi(k,,x,y)进行样条差值方法。同时我来学会了利用已知数据用P=polyfit(x,y,n)进行多项式的拟合。

总体来讲这节课的内容还是比较简单的，对于我们来说也只需按着书上，照猫画虎即可。但是实际操作中还是会遇到很多问题，比如错打了一个字等，都会使程序运行错误。可见实践是很重要的，同时还需要细心，脚踏实地的一步步来。这样才能使程序运行出来。 有些问题碰到了实在做不出来，我就和周围的同学讨论。在讨论与实际操作中，不但解决了问题，而且思维也拓宽了，学到了很多。