河南城建学院MATLAB上机实验内容答案解析

1、-/一 熟悉Matlab工作环境1、熟悉Mat lab的5个基本窗口思考题：变量如何声明，变量名须遵守什么规则、是否区分大小写。答：变量一般不需事先对变量的数据类型进行声明，系统会依据变量被赋值的类型自动 进行类型识别，也就是说变量可以直接赋值而不用提前声明。变量名要遵守以下几条规则：变量名必须以字母开头，只能由字母、数字或下划线组成。变量名区分大小写。变量名不能超过63个字符。关键字不能作为变量名。最好不要用特殊常量作为变量名。试说明分号、逗号、冒号的用法。分号：分隔不想显示计算结果的各语句；矩阵行与行的分隔符。逗号：分隔欲显示计算结果的各语句；变量分隔符；矩阵一行中各元素间的分隔符。冒号：

2、用于生成一维数值数组；表示一维数组的全部元素或多维数组某一维的全部元素。linspace()称为“线性等分”函数，说明它的用法。LINSPACE Linearly spaced vector. 线性等分函数LINSPACE(X1, X2) generates a row vector of 100 linearlyequally spaced points between X1 and X2.以X1为首元素，X2为末元素平均生成100个元素的行向量。LINSPACE(X1, X2, N) generates N points between X1 and X2.For N 2, LINSPAC

3、E returns X2.以X1为首元素，X2为末元素平均生成n个元素的行向量。如果nV2,返回X2。Class support for inputs X1,X2:float: double, single数据类型：单精度、双精度浮点型。说明函数 ones ()、zeros ()、eye ()的用法。ones ()生成全1矩阵。zeros ()生成全0矩阵。eye ()生成单位矩阵。2、Matlab的数值显示格式-/-/思考题:(1) 3次执行exist ( pi)的结果一样吗？如果不一样，试解释为什么？ pi pi=0; clearans exist(pi) exist(pi) pi pi=

4、0; clearans exist(pi) exist(pi)3.1416ansans sin(pi); exist(pi) pi pianspians exist(pi) pi pianspians3.1416答：3次执行的结果不一样。exist()函数是返回变量搜索顺序的一个函数。在第一次 执行时返回5代表变量pi是由Matlab构建的变量。在第二次执行时已经通过赋值语句定义 了变量pi，返回1代表pi是工作空间变量。第三次执行前清除了工作空间，此时pi为系统 默认常量，和第一次执行时性质一样，所以又返回5。(2)圆周率pi是系统默认常量，为什么会被改变为0。pi=0为赋值语句，此时pi不

5、再是系统默认常量，而是定义的变量了。二MATLAB语言基础1、向量的生成和运算练习：使用logspace()创建14n的有10个元素的行向量。 A=logspace(0,1.0992,10)A =1.00001.32471.75502.32493.07994.08015.40517.16039.4856 12.56612、矩阵的创建、引用和运算(1)矩阵的创建和引用练习：创建以下矩阵：A为3X4的全1矩阵、B为3X3的0矩阵、C为3X3的单位矩 阵、D为3X3的魔方阵、E由C和D纵向拼接而成、F抽取E的25行元素生成、G由F经变 形为3X4的矩阵而得、以G为子矩阵用复制函数生成6X8的大矩阵H

6、。 A=ones(3,4),B=zeros(3,3),C=eye(3,3),D=magic(3)A=11111111B =010000001000D =000816C =357100492 E=C;D, F=E(2:5,:), G=reshape(F,3,4)E =100010001816357492010001816357031101568007 H=repmat(G,2)H =0311031101560156800780070311031101560156800780072)矩阵运算练习：1）用矩阵除法求下列方程组的解x二：一二：二匚6x1 3x- 4* = 3 5x- 7x3 = 48x

7、t 3x3 = 7 A=6 3 4;-2 5 7;8 -1 -3,B=3;-4;-7B =6343-257-48-1-3-7 x=ABx =1.0200-14.00009.72002)求矩阵的秩; r=rank(A)r =33）求矩阵的特征值与特征向量 X,Lamda=eig（A）X =Lamda =0.8013-0.1094-0.16069.7326000.3638-0.65640.86690-3.292800.47490.7464-0.4719001.56024）矩阵的乘幂（平方）与开方 A2ans =62293334126262234262234 A1=sqrtm(A)A1 =2.244

8、7 + 0.2706i0.6974 -0.1400i0.9422 - 0.3494i-0.5815 + 1.6244i2.1005 -0.8405i1.7620 - 2.0970i1.9719 - 1.8471i-0.3017 +0.9557i0.0236 + 2.3845i5）矩阵的指数与对数（以e为底） Ae=expm(A)Ae =1.0e+004 *1.06530.54150.63230.48300.24650.28760.63160.32060.3745 Ael=logm(A)Ael =1.7129+0.4686i0.5305-0.2425i0.5429-0.6049i1.1938+2

9、.8123i0.3658-1.4552i-0.5514-3.6305i-0.0748-3.1978i0.7419+1.6546i1.8333+4.1282i6）矩阵的提取（取右上三角）与翻转（逆时针转90度） a=triu(A) a1=rot90(A)a =a1 =63447-305735-100-36-283、多维数组的创建及运算练习：创建三维数组a，第一页为;，第二页为;了，第三页为;；。然后用 reshape函数重排为数组B， B为3行、2列、2页。 a=1 3;4 2,b=1 2;2 1,c=3 5;7 1 A=cat(3,a,b,c)A(:,:,1)=A(:,:,2)=A(:,:,3

10、)=131235422171 B=reshape(A,3,2,2)B(:,:,1)=B(:,:,2)=122741153231三Matlab数值运算1、多项式运算练习：求厂一-.的商及余多项式。3& + 2E! +1 p1=conv(1 0 1,conv(1 3,1 1)pl = TOC o 1-5 h z 14443 q r=deconv(p1,1 0 2 1)q =14r =002-5-12、多形式插值和拟合有一组实验数据如附表1-1所示。请分别用拟合（二阶至三阶）和插值（线性和三次样条） 的方法来估测X=9.5时Y的值X12345678910Y1632701422604366821010

11、14321960 x=1:10;y=16 32 70 142 260 436 682 1010 1432 1960;y11.4232e+003 p1=polyfit(x,y,1) y1=polyval(p1,9.5)y11.4232e+003p1 =204.8000-522.40002.0000-1.00005.000010.00002.0000-1.00005.000010.0000y3 =1.6820e+003 y4=interp1(x,y,9.5)y4 =1696 y5=spline(x,y,9.5)y5 =1682p2二polyfit(x,y,2),y2=polyval(p2,9.5)

12、p2 =32.0000-147.2000181.6000y2 =1.6712e+003p3=polyfit(x,y,3),y3=polyval(p3,9.5)p3 =3、习题用函数roots求方程二-二-：=J的根 roots(1 -1 -1)ans =-0.61801.6180=二n: ：二2-，在n个节点(n不要太大，如取511)上用分段线性和三 次样条插值方法，计算m个插值点(m可取50100)的函数值。通过数值和图形输出，将两 种插值结果与精度进行比较。适当增加n，再作比较。 x=linspace(0,2*pi,8),y=sin(x)x =00.89761.79522.69283.59

13、044.48805.38566.2832y =00.78180.97490.4339-0.4339-0.9749-0.7818-0.0000 xi=linspace(0,2*pi,100);y0=sin(xi);y1=interp1(x,y,xi);y2=interp1(x,y,xi,spline ); plot(xi,y0,*,xi,y1,-.,xi,y2)0.80.60.40.20叩-0.2-0.4-0.6-0.8-101234512345671234567 e1=y1-y0;e2=y2-y0; plot(xi,e1)0.10.080.060.040.020-0.02-0.04-0.06-

14、0.08-0.1 plot(xi,e2)-0.0150.0150.010.0050-0.005-0.0101234567大气压强p随高度x变化的理论公式为：=,为验证这一公式， 测得某地大气压强随高度变化的一组数据如表所示。试用插值法和拟合法进行计算并绘图， 看那种方法较为合理，且总误差最小。高度/m0300600100015002000压强/Pa0.96890.93220.89690.85190.79890.7491插值法： x=0 300 600 1000 1500 2000; p=0.9689 0.9322 0.8969 0.8519 0.7989 0.7491; xi=linspace

15、(0,2000);p0=1.0332*exp(-(xi+500)/7756); p1=interp1(x,p,xi,spline); plot(xi,p0,*,xi,p1) e1=p1-p0; e=sum(e1.”2)e =1.8652e-005拟合法： x=0 300 600 1000 1500 2000; p=0.9689 0.9322 0.8969 0.8519 0.7989 0.7491; P=log10(p)-0.0137-0.0305-0.0473-0.0696-0.0975-0.1255 p1=polyfit(x,P,1) pl =-0.0001-0.0137 b=p1(1)/0

16、.4343,a=10.p1(2) b =-1.2863e-004 a =0.9689 xi=linspace(0,2000);p0=1.0332*exp(-(xi+500)/7756); p2=polyval(p1,xi);P2=10.”p2; e2=P2p0;e=sum(e2.”2) e =1.8116e-005四Matlab数值运算1、数值微积分练习：瑞士地图如图所示，为了算出其国土面积，首先对地图作如下测量：以由西向东 方向为X轴，由南到北方向为Y轴，选择方便的原点，并将从最西边界点到最东边界点但轴上 的区间适当划分为若干段，在每个分点的Y方向测出南边界点和北边界点的Y坐标Y1和Y2,根

17、 据地图比例尺知道18mm相当于40km，试由测量数据计算瑞士国土近似面积，与其精确值 41228km2 比较。X710.51317.53440.544.548566168.576.580.591Y14445475050383030343634414546Y24459707293100110110110117118116118118X96101104106.5111.5118123.5136.5142146150157158Y143373328326555545250666668Y2121124121121121116122838182868568x=7,10.5,13,17.5,34,40.

18、5,44.5,48,56,61,68.5,76.5,80.5,91,96,101,104,106.5,111.5,118,123.5,136.5,142,146,150,157,158;y1=44,45,47,50,50,38,30,30,34,36,34,41,45,46,43,37,33,28,32,65,55,54,52,50,66,6;y2=44,59,70,72,93,100,110,110,110,117,118,116,118,118,121,124,121,121,121,116,1; X=x./18*40;Y1=y1./18*40;Y2=y2./18*40; t1=trapz

19、(X,Y1),t2=trapz(X,Y2), t=t2-t1t1 =3.3819e+004t2 =7.6328e+004t =4.2510e+004 expt=t-41228expt =1.2819e+0032、习题利用梯形法和辛普森法求定积分的值，并对结果进行比较。如果积分区间改为55 结果有何不同？梯形积分中改变自变量x的维数，结果有何不同？ x=linspace(-3,3);y=exp(-x.”2/2); t=(1/2*pi)*trapz(x,y)q =t =3.92683.9267q=(1/2*pi)*quad(exp(-x.2/2),-3,3) x=linspace(-5,5);y=

20、exp(-x.”2/2); t=(1/2*pi)*trapz(x,y)q =t =3.93743.9374q=(1/2*pi)*quad(exp(-x.2/2),-5,5) x=linspace(-3,3,150);y=exp(-x.”2/2); t=(1/2*pi)*trapz(x,y) t =3.9268分别用矩形法、梯形法、辛普森法和牛顿-科茨4种方法近似计算定积分-三 取n=4，保留4位有效数字。矩形法： t=cumsum(y)*1/99;T=t(100) x=linspace(0,1);y=x./(x.”2+4);0.1126梯形法：0.1126梯形法： x=linspace(0,1

21、);y=x./(x.”2+4);辛普森法： q=quad(x./(x.”2+4),0,1)q =0.1116-/ t=trapz(x,y)t =0.1116牛顿-科茨法： q=quadl(x./(x.”2+4),0,1) q =0.1116五Matlab符号运算1、符号矩阵创建练习：分别用sym和syms创建符号表达式：_ 二:二二- ,-二n:， f1=sym(cos(x) + (sin(x)2)(1/2)f1 二cos(x) + (sin(x)2)(1/2) syms y e t f2=y/exp(-2*t)f2 =y/exp(-2*t)2、习题试创建以下2个矩阵：sinlsin.2sin

22、.3eee:lA =sin4sin.5sin6B =e-eD.sin7sin.8sin.9 -e6、符号表达式的变量替换练习：(1)已知工一二一：一三 :一，按照自变量x和自变量a， 表达式f分别进行降幕排列。 f=sym(a*x2+b*x+c-3)”3-a*(c*x2+4*b*x-1)f =(a*x2+b*x+c-3)”3-a*(c*x2+4*b*xT) fl=collect(f),f2=collect(f,a)fl =a3*x6+3*b*a2*x5+(c-3)*a2+2*b2*a+a*(2*(c-3)*a+b2)*x4+(4*(c-3)*b*a+b*(2*(c-3)*a+b2)*x3+(c

23、-3)*(2*(c-3)*a+b2)+2*b2*(c-3)+a*(c-3)2-a*c)*x2+(3*(c-3)2*b-4*b*a)*x+(c-3)”3+af2 =a3*x6+3*(b*x+c-3)*x4*a”2+(3*(b*x+c-3)”2*x2-c*x2-4*b*x+1)*a+(b*x+c-3)”38、符号方程的求解练习：(1) f=sym(x”2-1)/(x2-3*x+2); limit(f,x,2)ans =NaN 求函数f (x) =cos2x-sin2x的积分；求函数二=:二、n:的导数。 f=sym(cos(2*x)-sin(2*x); int(f)ans =1/2*sin(2*x

24、)+1/2*cos(2*x) g=sym(exp(x)+x*sin(x)”(1/2); diff(g)ans =1/2/(exp(x)+x*sin(x)”(1/2)*(exp(x)+sin(x)+x*cos(x)计算定积分匚二、n二一:二： f=sym(sin(x)+2); int(f,x,0,pi/6)ans =-1/2*3(1/2)+1/3*pi+1求下列线性方程组的解(n - y - z = 103x- 2y- z = 14I - 3y - z = 1 f1=sym(x+y+z=10); f2=sym(3*x+2*y+z=14); f3=sym(2*x+3*y-z=1); g=solve

25、(f1,f2,f3,x,y,z)g =x: 1x1 symy： 1x1 symz: 1x1 sym g.xans =1 g.yans =2 g.zans =7求解当y (0) =2, z (0) 二7时，微分方程组的解(dy _ .云一笠smx1 dz原 +y=l + X g_y,g_z=dsolve(Dy-z=sin(x),Dz+y=1+x,y(0)=2,z(0)=7,x)g_y =cos(x)+6*sin(x)+1/2*sin(x)*x+1+xg_z =-3/2*sin(x)+6*cos(x)+1+1/2*cos(x)*x六Matlab程序设计1、程序流程控制结构练习：(1)请把exp2.

26、函数文件用while循环改写。function s=exp3(x)n=1;s=0;while n=10e-6k=k+1;if rem(k,2)=0jspi=jspi+1/i;elsejspi二jspi-1/i;endi=i+2;endp=4*jspi ,k2、子函数和参数传递练习：编写求矩形面积函数rect，当没有输入参数时，显示提示信息；当只输入一个参 数时，则以该参数作为正方形的边长计算其面积；当有两个参数时，则以这两个参数为长和 宽计算其面积。function s=mianji(a,b)switch nargincase 0error(没有输入参数)case 1s=a*a;case 2s=a*b;end3、习题编写一个函数project1.m，其功能是判断某一年是否为闰年。function ryear(year)s=0;if rem(year,4)=0s=s+1;endif rem(year,100)=0s=s-1;endif rem(year,400)=0s=s