1、第一次练习教学要求:熟练掌握 Matlab 软件的基本命令和操作,会作二维、三维几何图形,能够用Matlab 软件解决微积分、线性代数与解析几何中的计算问题。补充命令vpa(x,n) 显示 x 的 n 位有效数字,教材 102 页fplot(f(x),a,b) 函数作图命令,画出 f(x)在区间a,b上的图形在下面的题目中 为你的学号的后 3 位(1-9 班)或 4 位(10 班以上)m1.1 计算 与30sinlixsilixmsyms xlimit(216*x-sin(216*x)/x3)ans =366935404/3limit(216*x-sin(216*x)/x3,inf)ans =
2、01.2 ,求 cos1xmyeysyms xdiff(exp(x)*cos(216*x/1000),2)ans =(46599*cos(451*x)/500)*exp(x)/250000 - (451*sin(451*x)/500)*exp(x)/2501.3 计算 210xyeddblquad(x,y) exp(x.2+y.2),0,1,0,1)ans =2.13941.4 计算42xdmsyms xint(x4/(2162+4*x2)ans =(91733851*atan(x/451)/4 - (203401*x)/4 + x3/121.5 (10)cos,xyemy求syms xdif
3、f(exp(x)*cos(216*x),10)ans =-356485076957717053044344387763*cos(216*x)*exp(x)-3952323024277642494822005884*sin(216*x)*exp(x)1.6 给出 在 的泰勒展式(最高次幂为 4). 10.mx0syms xtaylor(216/1000+x)(1/3),5,x)ans =-(9765625*451(1/2)*500(1/2)*x4)/82743933602 +(15625*451(1/2)*500(1/2)*x3)/91733851 -(125*451(1/2)*500(1/2)
4、*x2)/406802 + (451(1/2)*500(1/2)*x)/216 +(451(1/2)*500(1/2)/5001.7 Fibonacci 数列 的定义是 用循环语句编nx12,x12,(3,4)nnx程给出该数列的前 20 项(要求将结果用向量的形式给出) 。x=1,1;for n=3:20x(n)=x(n-1)+x(n-2);endxx=Columns 1 through 10 1 1 2 3 5 8 13 21 34 55Columns 11 through 20 89 144 233 377 610 987 1597 2584 4181 67651.8 对矩阵 ,求该矩阵
5、的逆矩阵,特征值,特征向量,行列式,210410Am计算 ,并求矩阵 ( 是对角矩阵) ,使得 。6,PD1APDA=-2,1,1;0,2,0;-4,1,216/1000;inv(A)ans =0.4107 0.0223 -0.45540 0.5000 01.8215 -0.4554 -0.9107eig(A)ans =-0.5490 + 1.3764i-0.5490 - 1.3764i2.0000 det(A)ans =4.3920P,D=eig(A)P = %特征向量0.3245 - 0.3078i 0.3245 + 0.3078i 0.2425 0 0 0.9701 0.8944 0.8
6、944 0.0000 D =-0.5490 + 1.3764i 0 0 0 -0.5490 - 1.3764i 0 0 0 2.0000 P*D6*inv(P) %A6 的值ans =15.3661 12.1585 + 0.0000i -5.8531 0 64.0000 0 23.4124 -5.8531 + 0.0000i -1.6196 1.9 作出如下函数的图形(注:先用 M 文件定义函数,再用 fplot 进行函数作图):122()1)xxfm 文件:function y=fenduan(x) if x syms n A=sym(4,2;1,3);x=1;2;P,D=eig(A) %没
7、有 sym 下面的矩阵就会显示为小数P = -1, 2 1, 1D = 2, 0 0, 5 An=P*Dn*inv(P)An = 2n/3 + (2*5n)/3, (2*5n)/3 - (2*2n)/3 5n/3 - 2n/3, (2*2n)/3 + 5n/3 xn=An*xxn =2*5n - 2n2n + 5n3.2 对于练习 1 中的 , ,求出 的通项. B3.0124ABTx)2,1(,()0(21nx syms n A=sym(2/5,1/5;1/10,3/10); x=1;2;P,D=eig(A) P = -1, 2 1, 1D = 1/5, 0 0, 1/2 An=P*Dn*i
8、nv(P)An = (2*(1/2)n)/3 + (1/5)n/3, (2*(1/2)n)/3 - (2*(1/5)n)/3 (1/2)n/3 - (1/5)n/3, (1/2)n/3 + (2*(1/5)n)/3xn =2*(1/2)n - (1/5)n(1/2)n + (1/5)n3.3 对随机给出的 ,观察数列 .该数列有极限吗?Tx),(0(2)1 )(12nx A=4,2;1,3;a=;x=2*rand(2,1)-1;for i=1:20a(i,1:2)=x;x=A*x; end for i=1:20if a(i,1)=0else t=a(i,2)/a(i,1);fprintf(%g
9、,%gn,i,t);endend 结论:在迭代 17 次后,发现数列 存在极限为 0.5)(12nx3.4 对 120 页中的例子,继续计算 .观察 及 的极限是否),21(,nyx,nyx)(nxm存在. (120 页练习 9) A=2.1,3.4,-1.2,2.3;0.8,-0.3,4.1,2.8;2.3,7.9,-1.5,1.4;3.5,7.2,1.7,-9.0;x0=1;2;3;4;x=A*x0;for i=1:1:100a=max(x);b=min(x);m=a*(abs(a)abs(b)+b*(abs(a) A=2.1,3.4,-1.2,2.3;0.8,-0.3,4.1,2.8;2
10、.3,7.9,-1.5,1.4;3.5,7.2,1.7,-9.0;P,D=eig(A)P =-0.3779 -0.8848 -0.0832 -0.3908-0.5367 0.3575 -0.2786 0.4777-0.6473 0.2988 0.1092 -0.7442-0.3874 -0.0015 0.9505 0.2555D =7.2300 0 0 00 1.1352 0 00 0 -11.2213 00 0 0 -5.8439结论:A 的绝对值最大特征值等于上面的 的极限相等,为什么呢?)(nxm还有,P 的第三列也就是-11.2213 对应的特征向量和上题求解到的 y 也有系数关系,两
11、者都是-11.2213 的特征向量。3.6 设 ,对问题 2 求出若干天之后的天气状态,并找出其特点(取Tp)25.0,.()04 位有效数字). (122 页练习 12) A2=3/4,1/2,1/4;1/8,1/4,1/2;1/8,1/4,1/4;P=0.5;0.25;0.25;for i=1:1:20P(:,i+1)=A2*P(:,i);endPP =Columns 1 through 100.5000 0.5625 0.5938 0.6035 0.6069 0.6081 0.6085 0.6086 0.6087 0.60870.2500 0.2500 0.2266 0.2207 0.2
12、185 0.2178 0.2175 0.2174 0.2174 0.21740.2500 0.1875 0.1797 0.1758 0.1746 0.1741 0.1740 0.1739 0.1739 0.1739Columns 11 through 200.6087 0.6087 0.6087 0.6087 0.6087 0.6087 0.6087 0.6087 0.6087 0.60870.2174 0.2174 0.2174 0.2174 0.2174 0.2174 0.2174 0.2174 0.2174 0.21740.1739 0.1739 0.1739 0.1739 0.1739
13、 0.1739 0.1739 0.1739 0.1739 0.1739Column 210.60870.21740.1739结论:9 天后,天气状态趋于稳定 P*=(0.6087,0.2174,0.1739) T3.7 对于问题 2,求出矩阵 的特征值与特征向量,并将特征向量与上一题中的结论作对2A比. (122 页练习 14) A2=3/4,1/2,1/4;1/8,1/4,1/2;1/8,1/4,1/4; P,D=eig(A2)P =-0.9094 -0.8069 0.3437-0.3248 0.5116 -0.8133-0.2598 0.2953 0.4695D =1.0000 0 00
14、0.3415 00 0 -0.0915分析:事实上,q=k(-0.9094, -0.3248, -0.2598)T 均为特征向量,而上题中 P*的 3 个分量之和为 1,可令 k(-0.9094, -0.3248, -0.2598)T=1,得 k=-0.6696.有 q=(0.6087, 0.2174, 0.1739),与 P*一致。3.8 对问题 1,设 为 的两个线性无关的特征向量,若2,p1A,具体求出上述的 ,将 表示成 的线性组合,求 的Tp)()0vu,)0(p21,p)(kp具体表达式,并求 时 的极限,与已知结论作比较. (123 页练习 16)k)(kp A=3/4,7/18
15、;1/4,11/18;P,D=eig(A);syms k pk;a=solve(u*P(1,1)+v*P(1,2)-1/2,u*P(2,1)+v*P(2,2)-1/2,u,v);pk=a.u*D(1,1).k*P(:,1)+a.v*D(2,2).k*P(:,2)pk =-5/46*(13/36)k+14/235/46*(13/36)k+9/23或者:p0=1/2;1/2;P,D=eig(sym(A);B=inv(sym(P)*p0B =5/469/23syms kpk=B(1,1)*D(1,1).k*P(:,1)+B(2,1)*D(2,2).k*P(:,2)pk =-5/46*(13/36)k
16、+14/235/46*(13/36)k+9/23 vpa(limit(pk,k,100),10)ans =.6086956522.3913043478结论:和用练习 12 中用迭代的方法求得的结果是一样的。第四次练习教学要求:会利用软件求勾股数,并且能够分析勾股数之间的关系。会解简单的近似计算问题。4.1 求满足 , 的所有勾股数,能否类似于(11.8) ,把它们用一个公式2bc10表示出来?解法程序 1:for b=1:998a=sqrt(b+2)2-b2);if(a=floor(a)fprintf(a=%i,b=%i,c=%in,a,b,b+2)endend运行结果:a=4,b=3,c=5
17、a=6,b=8,c=10a=8,b=15,c=17a=10,b=24,c=26a=12,b=35,c=37a=14,b=48,c=50a=16,b=63,c=65a=18,b=80,c=82a=20,b=99,c=101a=22,b=120,c=122a=24,b=143,c=145a=26,b=168,c=170a=28,b=195,c=197a=30,b=224,c=226a=32,b=255,c=257a=34,b=288,c=290a=36,b=323,c=325a=38,b=360,c=362a=40,b=399,c=401a=42,b=440,c=442a=44,b=483,c=4
18、85a=46,b=528,c=530a=48,b=575,c=577a=50,b=624,c=626a=52,b=675,c=677a=54,b=728,c=730a=56,b=783,c=785a=58,b=840,c=842a=60,b=899,c=901a=62,b=960,c=962解法程序 2: n=0;m=;for a=1:100for c=a+1:1000b=sqrt(c2-a2);if (b=floor(b) m(:,l)=a,b,c;endendendm勾股数 , 的解是:2bc102,),(, 2uua以下是推导过程:由 ,有22)(b4ba显然 , ,从而 是 2 的倍数
19、.设 ,代入上式得到:4 )1(2uau因为 ,从而 .2cc4.2 将上一题中 改为 , , , ,分别找出所有的勾股数.将它们与bb4567时的结果进行比较,然后用公式表达其结果。,1b(1) 时通项:4 )2(),2(),1, uuaa=8,b=6,c=10a=12,b=16,c=20a=16,b=30,c=34a=20,b=48,c=52a=24,b=70,c=74a=28,b=96,c=100a=32,b=126,c=130a=36,b=160,c=164a=40,b=198,c=202a=44,b=240,c=244a=48,b=286,c=290a=52,b=336,c=340a
20、=56,b=390,c=394a=60,b=448,c=452a=64,b=510,c=514a=68,b=576,c=580a=72,b=646,c=650a=76,b=720,c=724a=80,b=798,c=802a=84,b=880,c=884a=88,b=966,c=970(2) 5 时通项: bc 22,5(1),(),5(1)abcuua=15,b=20,c=25a=25,b=60,c=65a=35,b=120,c=125a=45,b=200,c=205a=55,b=300,c=305a=65,b=420,c=425a=75,b=560,c=565a=85,b=720,c=72
21、5a=95,b=900,c=905(3) 6 时通项bc )2(3),2(),16, uucbaa=12,b=9,c=15a=18,b=24,c=30a=24,b=45,c=51a=30,b=72,c=78a=36,b=105,c=111a=42,b=144,c=150a=48,b=189,c=195a=54,b=240,c=246a=60,b=297,c=303a=66,b=360,c=366a=72,b=429,c=435a=78,b=504,c=510a=84,b=585,c=591a=90,b=672,c=678a=96,b=765,c=771a=102,b=864,c=870a=10
22、8,b=969,c=975(4) 7 时通项bc )12(7),2(),1(7, uucbaa=21,b=28,c=35a=35,b=84,c=91a=49,b=168,c=175a=63,b=280,c=287a=77,b=420,c=427a=91,b=588,c=595a=105,b=784,c=791综上:当 c-b=k 为奇数时,通项 )12(),2(),1(, ukukcba当 c-b=k 为偶数时,通项 /,uk4.3 对 , ( ) ,对哪些 存在本原勾股数?(140 页练习 12)10ckb20程序:for k=1:200for b=1:999a=sqrt(b+k)2-b2)
23、;if(a=floor(a)break;endendend运行结果:1,2,8,9,18,25,32,49,50,72,81,98,121,128,162,169,200,4.4 设方程(11.15)的解构成数列 ,观察数列 , ,nqpnpq, , .你能得到哪些等式?试根据这些等式推导出关于nqp2n的递推关系式. (142 页练习 20),解:1000 以内解构成的数列 , , , , 如npqnp2nqpnqp下:n 1 2 3 4 5 62 7 26 97 362 1351p1 4 15 56 209 780nq3 11 41 153 571 21314 15 56 209 780
24、2911n1 3 11 41 153 571p我们发现这些解的关系似乎是:=1nqn= 2因为 = ,所以 。n1npnnqpq112np有以下结论:(4.1)123nn可以看成一个线性映射,令,TnnyxX)(2A(4.1)可写成: 1nA4.5 选取 对随机的 ,根据 的概率求出 的近似值。 (取自 130 页练习10m,ab(,)7)提示:(1)最大公约数的命令:gcd(a,b)(2)randint(1,1,u,v)产生一个在u,v区间上的随机整数程序: m=21600;s=0;for i=1:ma=randint(1,2,1,109);if gcd(a(1),a(2)=1;s=s+1;
25、endendpi=sqrt(6*m/s)pi =3.14314.6 用求定积分的 Monte Carlo 法近似计算 。 (102 页练习 16)提示:Monte Carlo 法近似计算 的一个例子。对于第一象限的正方形 ,内画出四分之一个圆01,xy向该正方形区域内随即投点,则点落在扇形区域内的概率为 .4投 次点,落在扇形内的次数为 ,则 ,因此 .nnc4nc程序如下n=100000;nc=0;for i=1:nx=rand;y=rand;if(x2+y2=1)nc=nc+1;endendpi=4*nc/n解:程序:a=0;b=1;m=1000;H=1;s=0;for i=1:mxi=rand();yi=H*rand();if yisqrt(1-xi2);s=s+1;endendpi=4*H*(b-a)*s/m运行结果:pi =3.1480综合题一、方程求根探究设方程 420x1.用 matlab 命令求该方程的所有根;0.2 0.4 0.6 0.8 10.20.40.60.81