数学实验.docx
《数学实验.docx》由会员分享，可在线阅读，更多相关《数学实验.docx（35页珍藏版）》请在冰点文库上搜索。
数学实验
第一次练习题
1、求
的所有根。
（先画图后求解）（要求贴图）
>>solve（'exp（x）-3*x^2',0）
ans=
-2*lambertw（-1/6*3^（1/2））
-2*lambertw（-1,-1/6*3^（1/2））
-2*lambertw（1/6*3^（1/2））
2、求下列方程的根。
1）
a=solve（'x^5+5*x+1',0）;a=vpa（a,6）
a=
1.10447+1.05983*i
-1.00450+1.06095*i
-.199936
-1.00450-1.06095*i
1.10447-1.05983*i
2）
至少三个根
>>fzero（'x*sin（x）-1/2',3）
ans=
2.9726
>>fzero（'x*sin（x）-1/2',-3）
ans=
-2.9726
>>fzero（'x*sin（x）-1/2',0）
ans=
-0.7408
3）
所有根
>>fzero（'sin（x）*cos（x）-x^2',0）
ans=
0
>>fzero（'sin（x）*cos（x）-x^2',0.6）
ans=
0.7022
3、求解下列各题：
1）
>>symx;
>>limit（（x-sin（x））/x^3）
ans=
1/6
2）
>>symx;
>>diff（exp（x）*cos（x）,10）
ans=
（-32）*exp（x）*sin（x）
3）
>>symx;
>>vpa（（int（exp（x^2）,x,0,1/2））,17）
ans=
0.54498710418362222
4）
>>symx;
>>int（x^4/（25+x^2）,x）
ans=
125*atan（x/5）-25*x+x^3/3
5）
>>symsx
f=sqrt（1+x）;
taylor（f,0,9）
ans=
-（429*x^8）/32768+（33*x^7）/2048-（21*x^6）/1024+（7*x^5）/256-（5*x^4）/128+x^3/16-x^2/8+x/2+1
6）
>>subs（diff（exp（sin（1/x））,x,3）,2）
ans=
-0.5826
4.
（1）
>>A=[-211;020;-413];inv（A）,eig（A）,[P,D]=eig（A）
ans=
-1.50000.50000.5000
00.50000
-2.00000.50001.0000
ans=
-1
2
2
P=
-0.7071-0.24250.3015
000.9045
-0.7071-0.97010.3015
D=
-100
020
002
4、求点（1,1,4）到直线L:
的距离
>>M0=[1,1,4];M1=[3,0,1];M0M1=M1-M0;
v=[-1,0,2];
d=norm（cross（M0M1,v））/norm（v）
d=
1.0954
5、已知
分别在下列条件下画出
的图形：
（要求贴图）
，在同一坐标系里作图
>>symsx;
>>fplot（'（1/sqrt（2*pi））*exp（-（（x）^2）/2）',[-3,3],'r'）
>>holdon
>>fplot（'（1/sqrt（2*pi））*exp（-（（x-1）^2）/2）',[-3,3],'y'）
>>holdon
>>fplot（'（1/sqrt（2*pi））*exp（-（（x+1）^2）/2）',[-3,3],'g'）
>>holdoff
，在同一坐标系里作图。
>>symsx;
fplot（'（1/sqrt（2*pi））*exp（-（（x）^2）/2）',[-3,3],'r'）
holdon
fplot（'（1/（sqrt（2*pi）*2））*exp（-（（x）^2）/（2*2^2））',[-3,3],'y'）
holdon
fplot（'（1/（sqrt（2*pi）*4））*exp（-（（x）^2）/（2*4^2））',[-3,3],'g'）
holdoff
6、画下列函数的图形：
（要求贴图）
（1）
>>ezmesh（'u*sin（t）','u*cos（t）','t/4',[0,20,0,2]）
（2）
>>x=0:
0.1:
3;y=x;
[XY]=meshgrid（x,y）;
Z=sin（X*Y）;
>>mesh（X,Y,Z）
（3）
ezmesh（'sin（t）*（3+cos（u））','cos（t）*（3+cos（u））','sin（u）',[0,2*pi,0,2*pi]）
第二次练习题
1、设
，数列
是否收敛？
若收敛，其值为多少？
精确到6位有效数字。
f=inline（'（x+7/x）/2'）;
>>x0=3;
>>fori=1:
20
x0=f（x0）;
fprintf（'%g,%g\n',i,x0）;
end
1,2.66667
2,2.64583
3,2.64575
4,2.64575
5,2.64575
6,2.64575
7,2.64575
8,2.64575
9,2.64575
10,2.64575
11,2.64575
12,2.64575
13,2.64575
14,2.64575
15,2.64575
16,2.64575
17,2.64575
18,2.64575
19,2.64575
20,2.64575
该数列收敛，收敛于2.64575。
2、设
是否收敛？
若收敛，其值为多少？
精确到17位有效数字。
（注：
学号为单号的取
，学号为双号的取
）
学号为10，所以取p=8
>>f=inline（'1/（x^8）'）;
x0=0;
fori=1:
100
x0=（x0+f（i））;
fprintf（'%g,%.16f\n',i,x0）;
end
1,1.0000000000000000
2,1.0039062500000000
3,1.0040586657902759
4,1.0040739245793384
5,1.0040764845793384
6,1.0040770799535192
7,1.0040772534200448
8,1.0040773130246896
9,1.0040773362552626
10,1.0040773462552626
11,1.0040773509203365
12,1.0040773532460168
13,1.0040773544719115
14,1.0040773551495150
15,1.0040773555396993
16,1.0040773557725300
17,1.0040773559158835
18,1.0040773560066281
19,1.0040773560655085
20,1.0040773561045711
21,1.0040773561310101
22,1.0040773561492331
23,1.0040773561620027
24,1.0040773561710874
25,1.0040773561776410
26,1.0040773561824297
27,1.0040773561859704
28,1.0040773561886174
29,1.0040773561906164
30,1.0040773561921406
31,1.0040773561933130
32,1.0040773561942224
33,1.0040773561949334
34,1.0040773561954934
35,1.0040773561959375
36,1.0040773561962919
37,1.0040773561965766
38,1.0040773561968066
39,1.0040773561969933
40,1.0040773561971459
41,1.0040773561972711
42,1.0040773561973744
43,1.0040773561974599
44,1.0040773561975311
45,1.0040773561975906
46,1.0040773561976406
47,1.0040773561976826
48,1.0040773561977181
49,1.0040773561977483
50,1.0040773561977738
51,1.0040773561977956
52,1.0040773561978142
53,1.0040773561978302
54,1.0040773561978440
55,1.0040773561978560
56,1.0040773561978664
57,1.0040773561978753
58,1.0040773561978831
59,1.0040773561978900
60,1.0040773561978960
61,1.0040773561979011
62,1.0040773561979057
63,1.0040773561979097
64,1.0040773561979133
65,1.0040773561979164
66,1.0040773561979193
67,1.0040773561979217
68,1.0040773561979239
69,1.0040773561979259
70,1.0040773561979277
71,1.0040773561979293
72,1.0040773561979306
73,1.0040773561979319
74,1.0040773561979330
75,1.0040773561979339
76,1.0040773561979348
77,1.0040773561979357
78,1.0040773561979364
79,1.0040773561979370
80,1.0040773561979377
81,1.0040773561979381
82,1.0040773561979386
83,1.0040773561979390
84,1.0040773561979395
85,1.0040773561979399
86,1.0040773561979404
87,1.0040773561979406
88,1.0040773561979408
89,1.0040773561979410
90,1.0040773561979413
91,1.0040773561979415
92,1.0040773561979417
93,1.0040773561979419
94,1.0040773561979421
95,1.0040773561979424
96,1.0040773561979426
97,1.0040773561979428
98,1.0040773561979430
99,1.0040773561979430
100,1.0040773561979430
所以原数列收敛，收敛于1.0040773561979430
书上习题
1.
>>f=inline（'（x-1）/（x+1）'）;
x0=5.5;
fori=1:
20
x0=f（x0）;
fprintf（'%g,%g\n',i,x0）;
end
1,0.692308
2,-0.181818
3,-1.44444
4,5.5
5,0.692308
6,-0.181818
7,-1.44444
8,5.5
9,0.692308
10,-0.181818
11,-1.44444
12,5.5
13,0.692308
14,-0.181818
15,-1.44444
16,5.5
17,0.692308
18,-0.181818
19,-1.44444
20,5.5
>>
所以该函数不收敛，迭代数列不是总收敛的。
2先分别求出分式线性函数
、
的不动点，再编程判断它们的迭代序列是否收敛．
（1）解方程
得到x=-1,是函数f1（x）不动点。
>>f=inline（'（x-1）/（x+3）'）;
x0=-5;
fori=1:
1000
x0=f（x0）;
fprintf（'%g,%g\n',i,x0）;
end
1,3
2,0.333333
3,-0.2
4,-0.428571
5,-0.555556
6,-0.636364
7,-0.692308
8,-0.733333
9,-0.764706
10,-0.789474
11,-0.809524
12,-0.826087
13,-0.84
14,-0.851852
15,-0.862069
16,-0.870968
17,-0.878788
18,-0.885714
19,-0.891892
20,-0.897436
21,-0.902439
22,-0.906977……….
…………….
981,-0.99796
982,-0.997962
983,-0.997964
984,-0.997966
985,-0.997969
986,-0.997971
987,-0.997973
988,-0.997975
989,-0.997977
990,-0.997979
991,-0.997981
992,-0.997983
993,-0.997985
994,-0.997987
995,-0.997989
996,-0.997991
997,-0.997993
998,-0.997995
999,-0.997997
1000,-0.997999
可以看出迭代序列逐渐接近-1，所以收敛于-1
（2）解方程
，得到x=－5和3，是函数f2（x）的不动点。
>>f=inline（'（-x+15）/（x+1）'）;
x0=5;
fori=1:
20
x0=f（x0）;
fprintf（'%g,%g\n',i,x0）;
end
1,1.66667
2,5
3,1.66667
4,5
5,1.66667
6,5
7,1.66667
8,5
9,1.66667
10,5
11,1.66667
12,5
13,1.66667
14,5
15,1.66667
16,5
17,1.66667
18,5
19,1.66667
20,5
由此可见迭代序列不收敛
练习4能否找到一个分式线性函数
，使它产生的迭代序列收敛到给定的数？
用这种办法近似计算
．
由
得
只要有解，才有可能收敛，并且收敛的数为其中的根。
可以令a=1,b=2,c=1,d=1;
f=inline（'（x+2）/（x+1）'）;%先定义函数
symsx;
x0=5.5;
fori=1:
1:
10
x0=f（x0）;
fprintf（'%g,%g\n',i,x0）;
end
1,1.15385
2,1.46429
3,1.4058
4,1.41566
5,1.41397
6,1.41426
7,1.41421
8,1.41421
9,1.41421
10,1.41421
7
（1）、函数
的迭代是否会产生混沌？
>>x1=0:
0.05:
0.5;
y1=2*x1;
x2=0.5:
0.05:
1;
y2=2*（1-x2）;
figure
plot（x1,y1,x2,y2）
gtext（'2*x'）
gtext（'2*（1-x）'）
练习8函数
=
（0
1）称为Logistic映射，试从“蜘蛛网”图观察它取初值为
=0.5产生的迭代序列的收敛性，将观察记录填人下表，若出现循环，请指出它的周期．
表4.3Logistic迭代的收敛性
3.3
3.5
3.56
3.568
3.6
3.84
序列收敛情况
no
no
no
no
no
no
0.3418
0.4942
0.5390
0.5465
0.5767
0.8193
f=inline（'a*x*（1-x）'）;
x=[0,1];
y=[];
x
（1）=0.5;
y
（1）=0;x
（2）=x
（1）;y
（2）=x
（1）;
fori=1:
100
x（1+2*i）=x（2*i）;
x（2+2*i）=f（x（1+2*i））;
y（1+2*i）=x（2+2*i）;
y（2+2*i）=y（1+2*i）;
end
plot（x,y,'r'）;
holdon;
symsxy;
y=x;
ezplot（x,[0,2]）;
ezplot（f（x）,[0,2]）;
axis（[0,2,0,2]）;
holdoff
练习12对例2，试着提高迭代次数至26000、28000、100000、500000等观察图形有什么变化．
functionMartin（a,b,c,N）
f=@（x,y）（y-sign（x）*sqrt（abs（b*x-c）））;
g=@（x）（a-x）;
m=[0;0];
forn=1:
N
m（:
n+1）=[f（m（1,n）,m（2,n））,g（m（1,n））];
end
plot（m（1,:
）,m（2,:
）,'kx'）;
axisequal
Martin（45.2,-300,5000）
Martin（45,2,-300,26000）
Martin（45,2,-300,28000）
Martin（45,2,-300,500000）
练习13取参数
为其他的值会得到什么图形？
参考表4.4.
表4.4Martin迭代参数参考表
-1000
0.1
-10
0.4
1
0
90
30
10
10
-10
100
-200
-4
-80
-137
17
4
10
100
-10
Martin（-1000,0.1,-10,5000）
Martin（0.4,1,0,5000）
Martin（90,30,10,5000）
Martin（10,-10,100,5000）
Martin（-137,17,4,5000）
Martin（10,100,-10,5000）
练习14 设A，B，C为某三角形的顶点，现作这样的迭代：
计算两个点的中点，这两个点分别是A,B,C中随机取得的一点，与前一步求得的中点（初始点任取）．当迭代次数大于10000时，试观察所得的散点图．
第三次练习题
练习2、对
，，求出平面映射
的通项，并画出这些点的散点图。
A=[4,2;1,3];
t=[];
fori=1:
20
x=2*rand（2,1）-1;
t（length（t）+1,1:
2）=x;
forj=1:
40
x=A*x;
t（length（t）+1,1:
2）=x;
end
end
plot（t（:
1）,t（:
2）,'*'）
grid（'on'）
练习3对于练习1中的
，
，求出
的通项.
B=[0.4,0.2;0.1,0.3];
t=[];
fori=1:
20
x=2*rand（2,1）-1;
t（length（t）+1,1:
2）=x;
forj=1:
40
x=A*x;
t（length（t）+1,1:
2）=x;
end
end
plot（t（:
1）,t（:
2）,'*'）
grid（'on'）
练习4对随机给出的
，观察数列
.该数列有极限吗?
第四次练习题
1.
编程找出
的所有勾股数，并问：
能否利用通项表示
?
>>forb=1:
995
a=sqrt（（b+5）^2-b^2）;
if（a==floor（a））
fprintf（'a=%i,b=%i,c=%i\n',a,b,b+5）
end
end
a=15,b=20,c=25
a=25,b=60,c=65
a=35,b=120,c=125
a=45,b=200,c=205
a=55,b=300,c=305
a=65,b=420,c=425
a=75,b=560,c=565
a=85,b=720,c=725
a=95,b=900,c=905
4、用MonteCarlo方法计算圆周率
。
clear
a=10000;
fori=1:
a
x=rand（）;
f=inline（'1/（1+（x^2））'）;
F（1,i）=f（x）;
end
jifen=mean（F）
jifen=
0.7870
>>0.7807*4
ans=
3.1228
2、12个篮球队A、B、C、D、E、F、G、H、I、J、K、L进行单循环比赛，其比赛结果如下：
B
C
D
E
F
G
H
I
J
K
L
A
A胜
C胜
A胜
A胜
F胜
G胜
A胜
I胜
A胜
K胜
L胜
B
B胜
B胜
B胜
F胜
G胜
H胜
B胜
J胜
B胜
B胜
C
D胜
E胜
C胜
C胜
C胜
I胜
C胜
K胜
L胜
D
E胜
D胜
G胜
D胜
D胜
J胜
D胜
L胜
E
F胜
E胜
H胜
E胜
J胜
K胜
E胜
F
G胜
F胜
I胜
J胜
F胜
F胜
G
H胜
G胜
G胜
K胜
L胜
H
H胜
J胜
H胜
L胜
I
J胜
I胜
L胜
J
升级会员 