MATLAB上机答案

1

一 熟悉Matlab工作环境

1、熟悉Matlab的5个基本窗口

思考题：

（1）变量如何声明，变量名须遵守什么规则、是否区分大小写。

答：变量一般不需事先对变量的数据类型进行声明，系统会依据变量被赋值的类型自动

进行类型识别，也就是说变量可以直接赋值而不用提前声明。变量名要遵守以下几条规则：

➢ 变量名必须以字母开头，只能由字母、数字或下划线组成。

➢ 变量名区分大小写。

➢ 变量名不能超过63个字符。

➢ 关键字不能作为变量名。

➢ 最好不要用特殊常量作为变量名。

（2）试说明分号、逗号、冒号的用法。

分号：分隔不想显示计算结果的各语句；矩阵行与行的分隔符。

逗号：分隔欲显示计算结果的各语句；变量分隔符；矩阵一行中各元素间的分隔符。

冒号：用于生成一维数值数组；表示一维数组的全部元素或多维数组某一维的全部元素。

（3）linspace（）称为“线性等分”函数，说明它的用法。

LINSPACE Linearly spaced vector. 线性等分函数

LINSPACE(X1, X2) generates a row vector of 100 linearly

equally spaced points between X1 and X2.

以X1为首元素，X2为末元素平均生成100个元素的行向量。

LINSPACE(X1, X2, N) generates N points between X1 and X2.

For N < 2, LINSPACE returns X2.

以X1为首元素，X2为末元素平均生成n个元素的行向量。如果n＜2，返回X2。

Class support for inputs X1,X2:

float: double, single

数据类型：单精度、双精度浮点型。

（4）说明函数ones（）、zeros（）、eye（）的用法。

ones（）生成全1矩阵。

zeros（）生成全0矩阵。

eye（）生成单位矩阵。

2、Matlab的数值显示格式

2

思考题：

（1）3次执行exist（’pi’）的结果一样吗？如果不一样，试解释为什么？

>> pi

ans =

3.1416

>> sin(pi);

>> exist('pi')

ans =

5

>> pi=0;

>> exist('pi')

ans =

1

>> pi

pi =

0

>> clear

>> exist('pi')

ans =

5

>> pi

ans =

3.1416

答：3次执行的结果不一样。exist（）函数是返回变量搜索顺序的一个函数。在第一次

执行时返回5代表变量pi是由Matlab构建的变量。在第二次执行时已经通过赋值语句定义

了变量pi，返回1代表pi是工作空间变量。第三次执行前清除了工作空间，此时pi为系统

默认常量，和第一次执行时性质一样，所以又返回5。

（2）圆周率pi是系统默认常量，为什么会被改变为0。

pi=0 为赋值语句，此时pi不再是系统默认常量，而是定义的变量了。

二 MATLAB语言基础

1、向量的生成和运算

练习：使用logspace（）创建1~4π的有10个元素的行向量。

>> A=logspace(0,1.0992,10)

A =

1.0000 1.3247 1.7550 2.3249 3.0799 4.0801 5.4051 7.1603

9.4856 12.5661

2、矩阵的创建、引用和运算

（1）矩阵的创建和引用

练习：创建以下矩阵：A为3×4的全1矩阵、B为3×3的0矩阵、C为3×3的单位矩

阵、D为3×3的魔方阵、E由C和D纵向拼接而成、F抽取E的2~5行元素生成、G由F经变

形为3×4的矩阵而得、以G为子矩阵用复制函数生成6×8的大矩阵H。

>> A=ones(3,4),B=zeros(3,3),C=eye(3,3),D=magic(3)

A =

1 1 1 1

1 1 1 1

1 1 1 1

3

B =

0 0 0

0 0 0

0 0 0

C =

1 0 0

>> E=[C;D], F=E(2:5,:), G=reshape(F,3,4)

E =

1 0 0

0 1 0

0 0 1

8 1 6

3 5 7

4 9 2

F =

0 1 0

0 0 1

8 1 6

3 5 7

G =

0 3 1 1

0 1 5 6

8 0 0 7

>> H=repmat(G,2)

H =

0 1 0

0 0 1

D =

8 1 6

3 5 7

4 9 2

0 3 1 1 0 3 1 1

0 1 5 6 0 1 5 6

8 0 0 7 8 0 0 7

0 3 1 1 0 3 1 1

0 1 5 6 0 1 5 6

4

8 0 0 7 8 0 0 7

2）矩阵运算

练习：1）用矩阵除法求下列方程组的解x=错误!未找到引用源。

>> A=[6 3 4;-2 5 7;8 -1 -3],B=[3;-4;-7]

A =

6 3 4

-2 5 7

8 -1 -3

>> x=AB

x =

1.0200

-14.0000

9.7200

2）求矩阵的秩；

>> r=rank(A)

r =

3

3）求矩阵的特征值与特征向量

>> [X,Lamda]=eig(A)

X =

0.8013 -0.1094 -0.1606

0.3638 -0.6564 0.8669

0.4749 0.7464 -0.4719

4）矩阵的乘幂（平方）与开方

>> A^2

ans =

62 29 33

34 12 6

26 22 34

Lamda =

9.7326 0 0

0 -3.2928 0

0 0 1.5602

B =

3

-4

-7

5

>> A1=sqrtm(A)

A1 =

2.2447 + 0.2706i 0.6974 - 0.1400i 0.9422 - 0.3494i

-0.5815 + 1.6244i 2.1005 - 0.8405i 1.7620 - 2.0970i

1.9719 - 1.8471i -0.3017 + 0.9557i 0.0236 + 2.3845i

5）矩阵的指数与对数（以e为底）

>> Ae=expm(A)

Ae =

1.0e+004 *

1.0653 0.5415 0.6323

0.4830 0.2465 0.2876

0.6316 0.3206 0.3745

>> Ael=logm(A)

Ael =

1.7129 + 0.4686i 0.5305 - 0.2425i 0.5429 - 0.6049i

1.1938 + 2.8123i 0.3658 - 1.4552i -0.5514 - 3.6305i

-0.0748 - 3.1978i 0.7419 + 1.6546i 1.8333 + 4.1282i

6）矩阵的提取（取右上三角）与翻转（逆时针转90度）

>> a=triu(A)

a =

6 3 4

0 5 7

0 0 -3

3、多维数组的创建及运算

练习：创建三维数组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) =

1 3

4 2

A(:,:,2) =

1 2

2 1

>> a1=rot90(A)

a1 =

4 7 -3

3 5 -1

6 -2 8

6

A(:,:,3) =

>> B=reshape(A,3,2,2)

B(:,:,1) =

1 2

4 1

3 2

3 5 7 1

B(:,:,2) =

2 7

1 5

3 1

三 Matlab数值运算

1、多项式运算

练习：求错误!未找到引用源。的商及余多项式。

>> p1=conv([1 0 1],conv([1 3],[1 1]))

p1 =

1 4 4 4 3

>> [q r]=deconv(p1,[1 0 2 1])

q =

1 4

r =

0 0 2 -5 -1

2、多形式插值和拟合

有一组实验数据如附表1-1所示。请分别用拟合（二阶至三阶）和插值（线性和三次样条）

的方法来估测X=9.5时Y的值

X

Y

1

16

2

32

3

70

4

142

5

260

6

436

7

682

8

1010

9

1432

10

1960

>> x=1:10;y=[16 32 70 142 260 436 682 1010 1432 1960];

>> p1=polyfit(x,y,1)

p1 =

204.8000 -522.4000

>> y1=polyval(p1,9.5)

y1 =

1.4232e+003

>>

p2=polyfit(x,y,2),y2=polyval(p2,9.5)

p2 =

32.0000 -147.2000 181.6000

y2 =

1.6712e+003

7

>>

p3=polyfit(x,y,3),y3=polyval(p3,9.5)

p3 =

2.0000 -1.0000 5.0000

10.0000

y3 =

3、习题

1.6820e+003

>> y4=interp1(x,y,9.5)

y4 =

1696

>> y5=spline(x,y,9.5)

y5 =

1682

（1）用函数roots求方程错误!未找到引用源。的根

>> roots([1 -1 -1])

ans =

-0.6180

1.6180

（2）错误!未找到引用源。，在n个节点（n不要太大，如取5~11）上用分段线性和三次

样条插值方法，计算m个插值点（m可取50~100）的函数值。通过数值和图形输出，将两种

插值结果与精度进行比较。适当增加n，再作比较。

>> x=linspace(0,2*pi,8),y=sin(x)

x =

0 0.8976 1.7952 2.6928 3.5904 4.4880 5.3856

6.2832

y =

0 0.7818 0.9749 0.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)

8

1

0.8

0.6

0.4

0.2

0

-0.2

-0.4

-0.6

-0.8

-1

01234567

>> e1=y1-y0;e2=y2-y0;

>> plot(xi,e1)

0.1

0.08

0.06

0.04

0.02

0

-0.02

-0.04

-0.06

-0.08

-0.1

01234567

>> plot(xi,e2)

9

0.015

0.01

0.005

0

-0.005

-0.01

-0.015

01234567

（3）大气压强p随高度x变化的理论公式为错误!未找到引用源。，为验证这一公式，测

得某地大气压强随高度变化的一组数据如表所示。试用插值法和拟合法进行计算并绘图，看

那种方法较为合理，且总误差最小。

高度/m

压强/Pa

0

0.9689

300

0.9322

600

0.8969

1000

0.8519

1500

0.7989

2000

0.7491

插值法：

>> x=[0 300 600 1000 1500 2000]; p=[0.9689 0.9322 0.8969 0.8519 0.7989 0.7491];

>> xi=linspace(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)

10

P =

-0.0137 -0.0305 -0.0473 -0.0696 -0.0975 -0.1255

>> p1=polyfit(x,P,1)

p1 =

-0.0001 -0.0137

>> b=p1(1)/0.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=P2-p0;e=sum(e2.^2)

e =

1.8116e-005

11

四 Matlab数值运算

1、数值微积分

练习：瑞士地图如图所示，为了算出其国土面积，首先对地图作如下测量：以由西向东

方向为X轴，由南到北方向为Y轴，选择方便的原点，并将从最西边界点到最东边界点在X轴上

的区间适当划分为若干段，在每个分点的Y方向测出南边界点和北边界点的Y坐标Y1和Y2，根

据地图比例尺知道18mm相当于40km，试由测量数据计算瑞士国土近似面积，与其精确值

41228km比较。

X

Y1

Y2

X

Y1

Y2

>>

x=[7,10.5,13,17.5,34,40.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

6,68];

>>

y2=[44,59,70,72,93,100,110,110,110,117,118,116,118,118,121,124,121,121,121,116,1

22,83,81,82,86,85,68];

>> X=x./18*40;Y1=y1./18*40;Y2=y2./18*40;

>> t1=trapz(X,Y1),t2=trapz(X,Y2), t=t2-t1

t1 =

3.3819e+004

t2 =

7.6328e+004

t =

96

43

101 104 106.5 111.5 118 123.5 136.5 142 146 150 157 158

37 33 28 32

121

65 55 54

83

52

81

50

82

66

86

66

85

68

68

7

44

44

10.5 13

45

59

47

70

17.5 34

50

72

50

93

40.5 44.5 48

38 30 30

56

34

61

36

68.5 76.5 80.5 91

34 41 45 46

2

100 110 110 110 117 118 116 118 118

121 124 121 121 116 122

12

4.2510e+004

>> expt=t-41228

expt =

1.2819e+003

2、习题

（4）利用梯形法和辛普森法求定积分的值，并对结果进行比较。如果积分区间改为-5~5

结果有何不同？梯形积分中改变自变量x的维数，结果有何不同？

>> x=linspace(-3,3);y=exp(-x.^2/2);

>> t=(1/2*pi)*trapz(x,y)

t =

3.9267

>>

q=(1/2*pi)*quad('exp(-x.^2/2)',-3,3)

>> x=linspace(-5,5);y=exp(-x.^2/2);

>> t=(1/2*pi)*trapz(x,y)

t =

3.9374

>>

q=(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

（5）分别用矩形法、梯形法、辛普森法和牛顿-科茨4种方法近似计算定积分错误!未找

到引用源。，取n=4，保留4位有效数字。

矩形法：

>> x=linspace(0,1);y=x./(x.^2+4);

>> t=cumsum(y)*1/99;T=t(100)

T =

q =

3.9374

q =

3.9268

13

0.1126

梯形法：

>> x=linspace(0,1);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创建符号表达式：错误!未找到引用源。，错误!未找到引用源。。

>> 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）试创建以下2个矩阵：

6、符号表达式的变量替换

练习：（1）已知错误!未找到引用源。,按照自变量x和自变量a，对表达式f分别进行降幂

排列。

>> f=sym('(a*x^2+b*x+c-3)^3-a*(c*x^2+4*b*x-1)')

f =(a*x^2+b*x+c-3)^3-a*(c*x^2+4*b*x-1)

>> f1=collect(f),f2=collect(f,'a')

14

f1 =

a^3*x^6+3*b*a^2*x^5+((c-3)*a^2+2*b^2*a+a*(2*(c-3)*a+b^2))*x^4+(4*(c-3)*b*a+b

*(2*(c-3)*a+b^2))*x^3+((c-3)*(2*(c-3)*a+b^2)+2*b^2*(c-3)+a*(c-3)^2-a*c)*x^2+(3*(

c-3)^2*b-4*b*a)*x+(c-3)^3+a

f2 =

a^3*x^6+3*(b*x+c-3)*x^4*a^2+(3*(b*x+c-3)^2*x^2-c*x^2-4*b*x+1)*a+(b*x+c-3)^3

8、符号方程的求解

练习：（1）求错误!未找到引用源。

>> f=sym('(x^2-1)/(x^2-3*x+2)');

>> limit(f,'x',2)

ans =NaN

（2）求函数f（x）=cos2x-sin2x的积分；求函数错误!未找到引用源。的导数。

>> f=sym('cos(2*x)-sin(2*x)');

>> int(f)

ans =1/2*sin(2*x)+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))

（3）计算定积分错误!未找到引用源。

>> f=sym('sin(x)+2');

>> int(f,'x',0,pi/6)

ans =-1/2*3^(1/2)+1/3*pi+1

（4）求下列线性方程组的解

>> f1=sym('x+y+z=10');

>> f2=sym('3*x+2*y+z=14');

>> f3=sym('2*x+3*y-z=1');

>> g=solve(f1,f2,f3,'x','y','z')

g =

15

x: [1x1 sym]

y: [1x1 sym]

z: [1x1 sym]

>> g.x

ans =1

>> g.y

ans =2

>> g.z

ans =7

（5）求解当y（0）=2，z（0）=7时，微分方程组的解

>> [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+x

g_z =-3/2*sin(x)+6*cos(x)+1+1/2*cos(x)*x

六 Matlab程序设计

1、程序流程控制结构

练习：（1）请把exp2.函数文件用while循环改写。

function s=exp3(x)

n=1;s=0;

while n<=x

s=s+n;

n=n+1;

end

s

（2）用错误!未找到引用源。公式求pi的近似值，直到最后一项的绝对值小于10为止，

试编写其M脚本文件。

k=0;jspi=1;i=3;

while (1/i)>=10e-6

k=k+1;

-6

16

if rem(k,2)==0

jspi=jspi+1/i;

else

jspi=jspi-1/i;

end

i=i+2;

end

p=4*jspi ,k

2、子函数和参数传递

练习：编写求矩形面积函数rect，当没有输入参数时，显示提示信息；当只输入一个参

数时，则以该参数作为正方形的边长计算其面积；当有两个参数时，则以这两个参数为长和

宽计算其面积。

function s=mianji(a,b)

switch nargin

case 0

error('没有输入参数')

case 1

s=a*a;

case 2

s=a*b;

end

3、习题

（3）编写一个函数project1.m，其功能是判断某一年是否为闰年。

function ryear(year)

s=0;

if rem(year,4)==0

s=s+1;

end

if rem(year,100)==0

s=s-1;

end

17

if rem(year,400)==0

s=s+1;

end

if s==1

fprintf('%4d 是闰年.n',year)

else

fprintf('%4d 不是闰年.n',year)

end

（4）编制一个函数，使得该函数能对输入的两个数值进行比较并返回其中的最小值。

function c = bijiao(a,b)

if nargin==2

if a < b

c=a;

else

c=b;

end

else

error('输入参数不正确')

end

（6）观察以下循环语句，计算每个循环的循环次数和循环结束之后var的值。

var=1;

while mod(var,10)~=0

var=var+1

end

循环次数10，var=10。

var=2;

>> while var<=100

var=var^2;

end

循环次数4，var=256。

>> var=3;

18

>> while var>100

var=var^2;

end

循环次数0次，var=3。

七 Matlab数据可视化

1、二维图形绘制

练习：写出图A2的绘制方法。

y1=sin(x);y2=cos(x);

plot(x,y1,'r -',x,y2,'m --')

x=linspace(0,4*pi);

y1=sin(x);y2=cos(x);

plot(x,y1,'r --',x,y2,'m -')

ylabel('·ù¶È'),xlabel('ʱ¼ä')

legend('sinx','cosx')

gtext('leftarrowsinx')

gtext('leftarrowcosx')

axis([0 16 -1 1])

xdate=0:0.5:16;ydate=0

line(xdate,ydate,'Color','k','Marker','.')

2、三维曲线和三维曲面绘制

练习：（1）绘制以上空间螺旋线的俯视图、左侧视图和前视图。

z=0:0.1:6*pi;

x=cos(z);y=sin(z);

subplot(2,2,1);plot3(cos(z),sin(z),z);

title('ÈýάÇúÏß')

subplot(2,2,2);plot(cos(z),sin(z));

title('¸©ÊÓͼ')

subplot(2,2,3);plot(cos(z),z);

title('×óÊÓͼ')

subplot(2,2,4);plot(sin(z),z);

19

title('Ç°ÊÓͼ')

（2）设错误!未找到引用源。，求定义域x=[-2，2]，y=[-2，2]内的z值（网格取0.1）。

请把z的值用网面图形象的表示出来，如图A3所示。

x=-2:0.1:2;y=x;

[X,Y]=meshgrid(x,y);

Z=X.^2.*exp(-(X.^2+Y.^2));

surf(X,Y,Z)