数学实验练习题参考答案

第一次练习

教学要求：熟练掌握Matlab 软件的基本命令和操作，会作二维、三维几何图形，能够用Matlab 软件解决微积分、线性代数与解析几何中的计算问题。

补充命令

vpa(x,n) 显示x 的n 位有效数字，教材102页

fplot(‘f(x)’,[a,b]) 函数作图命令，画出f(x)在区间[a,b]上的图形

在下面的题目中m 为你的学号的后3位（1-9班）或4位（10班以上） 1.1 计算30sin lim

x mx mx x →-与3sin lim x mx mx

x

→∞- 程序：

syms x

limit((1001*x-sin(1001*x))/x^3,x,0) 结果：

1003003001/6

程序： syms x

limit((1001*x-sin(1001*x))/x^3,x,inf) 结果： 0

1.2 cos

1000

x

mx

y e =，求''y 程序： syms x

diff(exp(x)*cos(1001*x/1000),2) 结果：

-2001/1000000*exp(x)*cos(1001/1000*x)-1001/500*exp(x)*sin(1001/1000*x)

1.3 计算

2

2

11

00

x

y e dxdy +??

程序：

dblquad(@(x,y) exp(x.^2+y.^2),0,1,0,1) 结果：

2.13935019514228

1.4 计算4

2

2

4x dx m x

+? 程序： syms x

int(x^4/(1000^2+4*x^2)) 结果：

1/12*x^3-1002001/16*x+1003003001/32*atan(2/1001*x)

1.5 (10)cos ,x y e mx y =求

程序： syms x

diff(exp(x)*cos(1000*x),10) 结果：

-1009999759158992000960720160000*exp(x)*cos(1001*x)-10090239998990319040000160032*exp(x)*sin(1001*x)

1.6 0x =的泰勒展式(最高次幂为4).

程序： syms x

taylor(sqrt(1001/1000+x),5) 结果：

1/100*10010^(1/2)+5/1001*10010^(1/2)*x-1250/1002001*10010^(1/2)*x ^2+625000/1003003001*10010^(1/2)*x^3-390625000/1004006004001*10010^(1/2)*x^4

1.7 Fibonacci 数列{}n x 的定义是121,1x x ==，

12,(3,4,)n n n x x x n --=+= 用循环语句编程给出该数列的前20项（要求

将结果用向量的形式给出）。 程序： x=[1,1]; for n=3:20

x(n)=x(n-1)+x(n-2); end x

结果：

Columns 1 through 10

1 1

2

3 5 8 13 21 3

4 55

Columns 11 through 20

89 144 233 377 610 987 1597 2584 4181 6765

1.8对矩阵

211

020

41

1000

A

m

??

?

-

?

= ?

?

- ?

??

，求该矩阵的逆矩阵，特征值，特

征向量，行列式，计算6

A，并求矩阵,P D（D是对角矩阵），使得1

A PDP-

=。

程序与结果：

a=[-2,1,1;0,2,0;-4,1,1001/1000];

inv(a)

0.50100100100100 -0.00025025025025 -0.50050050050050

0 0.50000000000000 0

2.00200200200200 -0.50050050050050 -1.00100100100100

eig(a)

-0.49950000000000 + 1.32230849275046i

-0.49950000000000 - 1.32230849275046i

2.00000000000000

[p,d]=eig(a)

p =

0.3355 - 0.2957i 0.3355 + 0.2957i 0.2425

0 0 0.9701

0.8944 0.8944 0.0000

注：p的列向量为特征向量

d =

-0.4995 + 1.3223i 0 0

0 -0.4995 - 1.3223i 0

0 0 2.0000

a^6

11.9680 13.0080 -4.9910

0 64.0000 0

19.9640 -4.9910 -3.0100

1.9 作出如下函数的图形（注：先用M 文件定义函数，再用fplot 进行函数作图）：

1202

()12(1)12

x x f x x x ?

≤≤??=??-<≤??

函数文件f.m: function y=f(x) if 0<=x&x<=1/2

y=2.0*x;

else 1/2

y=2.0*(1-x);

end

程序：fplot(@f,[0,1])

1.10 在同一坐标系下作出下面两条空间曲线（要求两条曲线用不同的颜色表示）

（1）cos sin x t y t z t =??=??=? （2）2cos 2sin x t y t z t =??

=??=?

程序：

t=-10:0.01:10; x1=cos(t); y1=sin(t); z1=t;

plot3(x1,y1,z1,'k');hold on x2=cos(2*t); y2=sin(2*t); z2=t;

plot3(x2,y2,z2,'r');hold off

1.11 已知422134305,203153211A B m -???? ? ?

=-=-- ? ? ? ?-????

，在MA TLAB 命令窗口中

建立A 、B 矩阵并对其进行以下操作：

(1) 计算矩阵A 的行列式的值det()A

(2) 分别计算下列各式：1122,*,.*,,,,T A B A B A B AB A B A A --- 解：（1）程序：

a=[4,-2,2;-3,0,5;1,5*1001,3]; b=[1,3,4;-2,0,3;2,-1,1];det(a) -130158

(2) 2*a-b 7 -7 0 -4 0 7

0 10011 5

a*b 12 10 12 7 -14 -7

-10003 0 15022

a.*b 4 -6 8

6 0 15 2 -5005 3

a*inv(b) 1.0e+003 *

-0.0000 0 0.0020 0.0000 0.0016 0.0001 1.1443 -1.0006 -1.5722

inv(a)*b 0.3463 0.5767 0.5383

0.0004 -0.0005 -0.0005 -0.1922 0.3460 0.9230

a^2 24 10002 4

-7 25031 9 -15008 15013 25036 A' 4 -3 1

-2 0 5005 2 5 3

1.12

已知22

()2()x f x μσ--=分别在下列条件下画出)(x f 的图形：

(1)/600m σ=,μ分别为0,1,1-(在同一坐标系上作图); (2)0μ=,σ分别为1,2,4,/100m (在同一坐标系上作图). (1)程序： x=-5:0.1:5;

h=inline('1/sqrt(2*pi)/s*exp(-(x-mu).^2/(2*s^2))');

y1=h(0,1001/600,x);y2=h(-1,1001/600,x);y3=h(1,1001/600,x); plot(x,y1,'r+',x,y2,'k-',x,y3,'b*')

(2)程序：

z1=h(0,1,x);z2=h(0,2,x);z3=h(0,4,x); z4=h(0,1001/100,x); plot(x,z1,'r+',x,z2,'k-',x,z3,'b*',x,z4, 'y:')

1.13 作出24

z mx y =+的函数图形。 程序：x=-5:0.1:5;y=-10:0.1:10;

[X Y]=meshgrid(x,y);Z=1001*X.^2+Y.^4; mesh(X,Y,Z);

1.14对于方程5

0.10200

m

x x -

-=，先画出左边的函数在合适的区间上的图形，借助于软件中的方程求根的命令求出所有的实根,找出函数的单调区间,结合高等数学的知识说明函数为什么在这些区间上是单调的,以及该方程确实只有你求出的这些实根。最后写出你做此题的体会。 解：作图程序：（注：x 范围的选择是经过试探而得到的）

x=-1.7:0.02:1.7;y=x.^5-1001/200*x-0.1; plot(x,y)；grid on;

-2

-1.5-1-0.500.51 1.52

由图形观察，在x=-1.5,x=0,x=1.5附近各有一个实根 求根程序：solve('x^5-1001/200*x-0.1') 结果：

-1.4906852047544424910680160298802 -.19980020616193485540810824654811e-1

.49944480891598282491814739731534e-2-1.4957641717395114847435704202656*i

.49944480891598282491814739731534e-2+1.4957641717395114847435704202656*i 1.5006763291923163201104639065887 三个实根的近似值分别为：

-1.490685，-0.019980，1.500676 由图形可以看出，函数在区间(,1)-∞-单调上升，在区间(1,1)-单调下降，

在区间(1,)∞单调上升。

diff('x^5-1001/200*x-0.1',x) 结果为5*x^4-1001/200

solve('5*x^4-1001/200.')得到两个实根：-1.0002499与1.0002499 可以验证导函数在(, 1.0002499)-∞-内为正，函数单调上升 导函数在( 1.0002499,1.0002499)-内为负，函数单调下降 导函数在(1.0002499,)∞内为正，函数单调上升 根据函数的单调性，最多有3个实根。

1.15 求23m 0x e x -=的所有根。（先画图后求解）（要求贴图） 作图命令：（注：x 范围的选择是经过试探而得到的） x=-5:0.001:15;y=exp(x)-3*1001*x.^2; plot(x,y);grid on;

6

可以看出，在(-5,5)内可能有根，在(10,15)内有1个根

将(-5,5)内图形加细，最终发现在(-0.03,0.03)内有两个根。 用solve('exp(x)-3*1001.0*x^2',x)可以求出3个根为： .18417113274368129311145677478702e-1 13.162041092091149185726742857195 -.18084038990284796648194134222365e-1 即：-0.018417,0.018084,13.16204

第二次练习

教学要求：要求学生掌握迭代、混沌的判断方法，以及利用迭代思想解决实际问题。

2.1 设11

()/23n n

n m x x x x +?

=+???=?，数列{}n x 是否收敛？若收敛，其值为多少？

精确到8位有效数字。 解：程序代码如下（m=1000）： >> f=inline('(x+1000/x)/2'); x0=3; for i=1:20; x0=f(x0);

fprintf('%g,%g\n',i,x0); end 运行结果：

1,168.167 11,31.6228 2,87.0566 12,31.6228 3,49.2717 13,31.6228 4,34.7837 14,31.6228 5,31.7664 15,31.6228 6,31.6231 16,31.6228 7,31.6228 17,31.6228 8,31.6228 18,31.6228 9,31.6228 19,31.6228 10,31.6228 20,31.6228

由运行结果可以看出，，数列{}n x 收敛，其值为31.6228。

2.2求出分式线性函数

2

12

1

(),()

x x m

f x f x

x m x m

-+

==

++

的不动点，再编程判

断它们的迭代序列是否收敛。

解：取m=1000.

（1）程序如下：

f=inline('(x-1)/(x+1000)');

x0=2;

for i=1:20;

x0=f(x0);

fprintf('%g,%g\n',i,x0);

end

运行结果：

1,0.000998004 11,-0.001001

2,-0.000999001 12,-0.001001

3,-0.001001 13,-0.001001

4,-0.001001 14,-0.001001

5,-0.001001 15,-0.001001

6,-0.001001 16,-0.001001

7,-0.001001 17,-0.001001

8,-0.001001 18,-0.001001

9,-0.001001 19,-0.001001

10,-0.001001 20,-0.001001

由运行结果可以看出，，分式线性函数收敛，其值为-0.001001。易见函数的不动点为-0.001001（吸引点）。

（2）程序如下：

f=inline('(x+1000000)/(x+1000)');

x0=2;

for i=1:20; x0=f(x0);

fprintf('%g,%g\n',i,x0); end 运行结果：

1,998.006 11,618.332 2,500.999 12,618.302 3,666.557 13,618.314 4,600.439 14,618.309 5,625.204 15,618.311 6,615.692 16,618.31 7,619.311 17,618.311 8,617.929 18,618.31 9,618.456 19,618.31 10,618.255 20,618.31

由运行结果可以看出，，分式线性函数收敛，其值为618.31。易见函数的不动点为618.31（吸引点）。 2.3 下面函数的迭代是否会产生混沌？（56页练习7（1））

1202

()12(1)1

2

x x f x x x ?

≤≤??=??-<≤??

解：程序如下：

f=inline('1-2*abs(x-1/2)'); x=[]; y=[]; x(1)=rand();

y(1)=0;x(2)=x(1);y(2)=f(x(1));

for i=1:100;

x(1+2*i)=y(2*i);

x(2+2*i)=x(1+2*i);

y(2+2*i)=f(x(2+2*i));

end

plot(x,y,'r');

hold on;

syms x;

ezplot(x,[0,1/2]);

ezplot(f(x),[0,1]);

axis([0,1/2,0,1]);

>> hold off

运行结果：

1 -

2 abs(x - 1/2)

00.050.10.150.20.250.30.350.40.450.5

x

2.4 函数()(1)(01)f x x x x α=-≤≤称为Logistic 映射，试从“蜘蛛网”图观察它取初值为00.5x =产生的迭代序列的收敛性，将观察记录填人下表，若出现循环，请指出它的周期．

解：当α=3.3时，程序代码如下： f=inline('3.3*x*(1-x)'); x=[]; y=[]; x(1)=0.5;

y(1)=0;x(2)=x(1);y(2)=f(x(1)); for i=1:1000; x(1+2*i)=y(2*i); x(2+2*i)=x(1+2*i); y(1+2*i)=x(1+2*i); y(2+2*i)=f(x(2+2*i)); end

plot (x,y,'r'); hold on; syms x;

ezplot(x,[0,1]); ezplot(f(x),[0,1]); axis([0,1,0,1]); hold off 运行结果：

-(33 x (x - 1))/10

00.10.20.30.40.50.60.70.80.91

x

当α=3.5时，上述程序稍加修改，得：

-(7 x (x - 1))/2

00.10.20.30.40.50.60.70.80.91

x

当α=3.56时，得：

-(89 x (x - 1))/25

00.10.20.30.40.50.60.70.80.91

x

当 =3.568时，得：

-(446 x (x - 1))/125

00.10.20.30.40.50.60.70.80.91

x

当α=3.6时，得：

-(18 x (x - 1))/5

00.10.20.30.40.50.60.70.80.91

x

当α=3.84时，得：

-(96 x (x - 1))/25

00.10.20.30.40.50.60.70.80.91

x

a b c为其它的值会得到什么图形？参考下

2.5对于Martin迭代，取参数,,

表（取自63页练习13）

解：取m=10000；迭代次数N=20000；

在M-文件里面输入代码：

function Martin(a,b,c,N)

f=@(x,y)(y-sign(x)*sqrt(abs(b*x-c)));

g=@(x)(a-x);

m=[0;0];

for n=1:N

m(:,n+1)=[f(m(1,n),m(2,n)),g(m(1,n))];

end

plot(m(1,:),m(2,:),'kx');

axis equal

在命令窗口中执行Martin（10000，10000，10000，20000），得：

x 104执行Martin（-10000，-10000，10000，20000），得：

x 104

执行Martin（-10000，10，-10000，20000），得：

执行Martin（10，10，0.5，20000），得：