MATLAB數(shù)學(xué)實(shí)驗(yàn)第二版答案(胡良劍).doc

《MATLAB數(shù)學(xué)實(shí)驗(yàn)第二版答案(胡良劍).doc》由會(huì)員分享，可在線(xiàn)閱讀，更多相關(guān)《MATLAB數(shù)學(xué)實(shí)驗(yàn)第二版答案(胡良劍).doc（24頁(yè)珍藏版）》請(qǐng)?jiān)谘b配圖網(wǎng)上搜索。

1、數(shù)學(xué)實(shí)驗(yàn)答案Chapter 1 Page20,ex1 (5) 等于exp(1),exp(2);exp(3),exp(4) (7) 3=1*3, 8=2*4 (8) a為各列最小值，b為最小值所在的行號(hào) (10) 1=4,false, 2=3,false, 3=2, ture, 4=1,ture (11) 答案表明：編址第2元素滿(mǎn)足不等式(30=20)和編址第4元素滿(mǎn)足不等式(40=10) (12) 答案表明：編址第2行第1列元素滿(mǎn)足不等式(30=20)和編址第2行第2列元素滿(mǎn)足不等式(40=10) Page20, ex2 (1)a, b, c的值盡管都是1， 但數(shù)據(jù)類(lèi)型分別為數(shù)值，字符， 邏輯

2、， 注意a與c相等， 但他們不等于b (2)double(fun)輸出的分別是字符a,b,s,(,x,)的ASCII碼 Page20,ex3 r=2;p=0.5;n=12; T=log(r)/n/log(1+0.01*p) Page20,ex4 x=-2:0.05:2;f=x.4-2.x; fmin,min_index=min(f) 最小值 最小值點(diǎn)編址 x(min_index) ans = 0.6500 最小值點(diǎn) f1,x1_index=min(abs(f) 求近似根-絕對(duì)值最小的點(diǎn) f1 = 0.0328 x1_index = 24 x(x1_index) ans = -0.8500 x(

3、x1_index)=;f=x.4-2.x; 刪去絕對(duì)值最小的點(diǎn)以求函數(shù)絕對(duì)值次小的點(diǎn) f2,x2_index=min(abs(f) 求另一近似根-函數(shù)絕對(duì)值次小的點(diǎn) f2 = 0.0630 x2_index = 65 x(x2_index) ans = 1.2500 Page20,ex5 z=magic(10) z = 92 99 1 8 15 67 74 51 58 40 98 80 7 14 16 73 55 57 64 41 4 81 88 20 22 54 56 63 70 47 85 87 19 21 3 60 62 69 71 28 86 93 25 2 9 61 68 75 52

4、 34 17 24 76 83 90 42 49 26 33 65 23 5 82 89 91 48 30 32 39 66 79 6 13 95 97 29 31 38 45 72 10 12 94 96 78 35 37 44 46 53 11 18 100 77 84 36 43 50 27 59 sum(z) sum(diag(z) z(:,2)/sqrt(3) z(8,:)=z(8,:)+z(3,:) Chapter 2 Page 45 ex1 先在編輯器窗口寫(xiě)下列M函數(shù)，保存為eg2_1.m function xbar,s=ex2_1(x) n=length(x); xbar=su

5、m(x)/n; s=sqrt(sum(x.2)-n*xbar2)/(n-1); 例如 x=81 70 65 51 76 66 90 87 61 77; xbar,s=ex2_1(x) Page 45 ex2 s=log(1);n=0; while se k=k+1;F(k)=F(k-1)+F(k-2); x=F(k)/F(k-1); end a,x,k 計(jì)算至k=21可滿(mǎn)足精度 Page 45 ex4 clear;tic;s=0; for i=1:1000000 s=s+sqrt(3)/2i; end s,toc tic;s=0;i=1; while i1.1)+x.*(x=-1.1)-1.1

6、*(x1); p=p+b*exp(-y.2-6*x.2).*(x+y-1).*(x+y=1); p=p+a*exp(-0.75*y.2-3.75*x.2+1.5*x).*(x+y A=1 2 3;4 5 6;7 8 0;C=2 -5 -22;-5 -24 -56;-22 -56 -16; X=lyap(A,C) X = 1.0000 -1.0000 -0.0000 -1.0000 2.0000 1.0000 -0.0000 1.0000 7.0000 Chapter 3 Page65 Ex1 a=1,2,3;b=2,4,3;a./b,a.b,a/b,ab ans = 0.5000 0.5000

7、 1.0000 ans = 2 2 1 ans = 0.6552 一元方程組x2,4,3=1,2,3的近似解 ans = 0 0 0 0 0 0 0.6667 1.3333 1.0000 矩陣方程1,2,3x11,x12,x13;x21,x22,x23;x31,x32,x33=2,4,3的特解 Page65 Ex 2(1) A=4 1 -1;3 2 -6;1 -5 3;b=9;-2;1; rank(A), rank(A,b) A,b為增廣矩陣 ans = 3 ans = 3 可見(jiàn)方程組唯一解 x=Ab x = 2.3830 1.4894 2.0213 (2) A=4 -3 3;3 2 -6;1

8、 -5 3;b=-1;-2;1; rank(A), rank(A,b) ans = 3 ans = 3 可見(jiàn)方程組唯一解 x=Ab x = -0.4706 -0.2941 0 (3) A=4 1;3 2;1 -5;b=1;1;1; rank(A), rank(A,b) ans = 2 ans = 3 可見(jiàn)方程組無(wú)解 x=Ab x = 0.3311 -0.1219 最小二乘近似解 (4) a=2,1,-1,1;1,2,1,-1;1,1,2,1;b=1 2 3;%注意b的寫(xiě)法 rank(a),rank(a,b) ans = 3 ans = 3 rank(a)=rank(a,b) ab ans =

9、1 0 1 0 一個(gè)特解 Page65 Ex3 a=2,1,-1,1;1,2,1,-1;1,1,2,1;b=1,2,3; x=null(a),x0=ab x = -0.6255 0.6255 -0.2085 0.4170 x0 = 1 0 1 0 通解kx+x0 Page65 Ex 4 x0=0.2 0.8;a=0.99 0.05;0.01 0.95; x1=a*x, x2=a2*x, x10=a10*x x=x0;for i=1:1000,x=a*x;end,x x = 0.8333 0.1667 x0=0.8 0.2; x=x0;for i=1:1000,x=a*x;end,x x = 0

10、.8333 0.1667 v,e=eig(a) v = 0.9806 -0.7071 0.1961 0.7071 e = 1.0000 0 0 0.9400 v(:,1)./x ans = 1.1767 1.1767 成比例，說(shuō)明x是最大特征值對(duì)應(yīng)的特征向量 Page65 Ex5 用到公式(3.11)(3.12) B=6,2,1;2.25,1,0.2;3,0.2,1.8;x=25 5 20; C=B/diag(x) C = 0.2400 0.4000 0.0500 0.0900 0.2000 0.0100 0.1200 0.0400 0.0900 A=eye(3,3)-C A = 0.7600

11、 -0.4000 -0.0500 -0.0900 0.8000 -0.0100 -0.1200 -0.0400 0.9100 D=17 17 17;x=AD x = 37.5696 25.7862 24.7690 Page65 Ex 6(1) a=4 1 -1;3 2 -6;1 -5 3;det(a),inv(a),v,d=eig(a) ans = -94 ans = 0.2553 -0.0213 0.0426 0.1596 -0.1383 -0.2234 0.1809 -0.2234 -0.0532 v = 0.0185 -0.9009 -0.3066 -0.7693 -0.1240 -0.

12、7248 -0.6386 -0.4158 0.6170 d = -3.0527 0 0 0 3.6760 0 0 0 8.3766 (2) a=1 1 -1;0 2 -1;-1 2 0;det(a),inv(a),v,d=eig(a) ans = 1 ans = 2.0000 -2.0000 1.0000 1.0000 -1.0000 1.0000 2.0000 -3.0000 2.0000 v = -0.5773 0.5774 + 0.0000i 0.5774 - 0.0000i -0.5773 0.5774 0.5774 -0.5774 0.5773 - 0.0000i 0.5773 +

13、0.0000i d = 1.0000 0 0 0 1.0000 + 0.0000i 0 0 0 1.0000 - 0.0000i (3) A=5 7 6 5;7 10 8 7;6 8 10 9;5 7 9 10 A = 5 7 6 5 7 10 8 7 6 8 10 9 5 7 9 10 det(A),inv(A), v,d=eig(A) ans = 1 ans = 68.0000 -41.0000 -17.0000 10.0000 -41.0000 25.0000 10.0000 -6.0000 -17.0000 10.0000 5.0000 -3.0000 10.0000 -6.0000

14、-3.0000 2.0000 v = 0.8304 0.0933 0.3963 0.3803 -0.5016 -0.3017 0.6149 0.5286 -0.2086 0.7603 -0.2716 0.5520 0.1237 -0.5676 -0.6254 0.5209 d = 0.0102 0 0 0 0 0.8431 0 0 0 0 3.8581 0 0 0 0 30.2887 (4)(以n=5為例) 方法一（三個(gè)for） n=5; for i=1:n, a(i,i)=5;end for i=1:(n-1),a(i,i+1)=6;end for i=1:(n-1),a(i+1,i)=1;

15、end a方法二（一個(gè)for） n=5;a=zeros(n,n); a(1,1:2)=5 6; for i=2:(n-1),a(i,i-1,i,i+1)=1 5 6;end a(n,n-1 n)=1 5; a方法三（不用for） n=5;a=diag(5*ones(n,1); b=diag(6*ones(n-1,1); c=diag(ones(n-1,1); a=a+zeros(n-1,1),b;zeros(1,n)+zeros(1,n);c,zeros(n-1,1) 下列計(jì)算 det(a) ans = 665 inv(a) ans = 0.3173 -0.5865 1.0286 -1.624

16、1 1.9489 -0.0977 0.4887 -0.8571 1.3534 -1.6241 0.0286 -0.1429 0.5429 -0.8571 1.0286 -0.0075 0.0376 -0.1429 0.4887 -0.5865 0.0015 -0.0075 0.0286 -0.0977 0.3173 v,d=eig(a) v = -0.7843 -0.7843 -0.9237 0.9860 -0.9237 0.5546 -0.5546 -0.3771 -0.0000 0.3771 -0.2614 -0.2614 0.0000 -0.1643 0.0000 0.0924 -0.0

17、924 0.0628 -0.0000 -0.0628 -0.0218 -0.0218 0.0257 0.0274 0.0257 d = 0.7574 0 0 0 0 0 9.2426 0 0 0 0 0 7.4495 0 0 0 0 0 5.0000 0 0 0 0 0 2.5505 Page65 Ex 7(1) a=4 1 -1;3 2 -6;1 -5 3;v,d=eig(a) v = 0.0185 -0.9009 -0.3066 -0.7693 -0.1240 -0.7248 -0.6386 -0.4158 0.6170 d = -3.0527 0 0 0 3.6760 0 0 0 8.3

18、766 det(v) ans = -0.9255 %v行列式正常, 特征向量線(xiàn)性相關(guān)，可對(duì)角化 inv(v)*a*v 驗(yàn)算 ans = -3.0527 0.0000 -0.0000 0.0000 3.6760 -0.0000 -0.0000 -0.0000 8.3766 v2,d2=jordan(a) 也可用jordan v2 = 0.0798 0.0076 0.9127 0.1886 -0.3141 0.1256 -0.1605 -0.2607 0.4213 特征向量不同 d2 = 8.3766 0 0 0 -3.0527 - 0.0000i 0 0 0 3.6760 + 0.0000i v

19、2a*v2 ans = 8.3766 0 0.0000 0.0000 -3.0527 0.0000 0.0000 0.0000 3.6760 v(:,1)./v2(:,2) 對(duì)應(yīng)相同特征值的特征向量成比例 ans = 2.4491 2.4491 2.4491 (2) a=1 1 -1;0 2 -1;-1 2 0;v,d=eig(a) v = -0.5773 0.5774 + 0.0000i 0.5774 - 0.0000i -0.5773 0.5774 0.5774 -0.5774 0.5773 - 0.0000i 0.5773 + 0.0000i d = 1.0000 0 0 0 1.000

20、0 + 0.0000i 0 0 0 1.0000 - 0.0000i det(v) ans = -5.0566e-028 -5.1918e-017i v的行列式接近0, 特征向量線(xiàn)性相關(guān)，不可對(duì)角化 v,d=jordan(a) v = 1 0 1 1 0 0 1 -1 0 d = 1 1 0 0 1 1 0 0 1 jordan標(biāo)準(zhǔn)形不是對(duì)角的，所以不可對(duì)角化 (3) A=5 7 6 5;7 10 8 7;6 8 10 9;5 7 9 10 A = 5 7 6 5 7 10 8 7 6 8 10 9 5 7 9 10 v,d=eig(A) v = 0.8304 0.0933 0.3963 0.

21、3803 -0.5016 -0.3017 0.6149 0.5286 -0.2086 0.7603 -0.2716 0.5520 0.1237 -0.5676 -0.6254 0.5209 d = 0.0102 0 0 0 0 0.8431 0 0 0 0 3.8581 0 0 0 0 30.2887 inv(v)*A*v ans = 0.0102 0.0000 -0.0000 0.0000 0.0000 0.8431 -0.0000 -0.0000 -0.0000 0.0000 3.8581 -0.0000 -0.0000 -0.0000 0 30.2887 本題用jordan不行, 原因未

22、知 (4)參考6(4)和7(1) Page65 Exercise 8 只有(3)對(duì)稱(chēng), 且特征值全部大于零, 所以是正定矩陣. Page65 Exercise 9(1) a=4 -3 1 3;2 -1 3 5;1 -1 -1 -1;3 -2 3 4;7 -6 -7 0 rank(a) ans = 3 rank(a(1:3,:) ans = 2 rank(a(1 2 4,:) 1,2,4行為最大無(wú)關(guān)組 ans = 3 b=a(1 2 4,:);c=a(3 5,:); bc 線(xiàn)性表示的系數(shù) ans = 0.5000 5.0000 -0.5000 1.0000 0 -5.0000 Page65 Ex

23、ercise 10 a=1 -2 2;-2 -2 4;2 4 -2 v,d=eig(a) v = 0.3333 0.9339 -0.1293 0.6667 -0.3304 -0.6681 -0.6667 0.1365 -0.7327 d = -7.0000 0 0 0 2.0000 0 0 0 2.0000 v*v ans = 1.0000 0.0000 0.0000 0.0000 1.0000 0 0.0000 0 1.0000 v確實(shí)是正交矩陣 Page65 Exercise 11 設(shè)經(jīng)過(guò)6個(gè)電阻的電流分別為i1, ., i6. 列方程組如下 20-2i1=a; 5-3i2=c; a-3i

24、3=c; a-4i4=b; c-5i5=b; b-3i6=0; i1=i3+i4;i5=i2+i3;i6=i4+i5; 計(jì)算如下 A=1 0 0 2 0 0 0 0 0;0 0 1 0 3 0 0 0 0;1 0 -1 0 0 -3 0 0 0; 1 -1 0 0 0 0 -4 0 0; 0 -1 1 0 0 0 0 -5 0;0 1 0 0 0 0 0 0 -3; 0 0 0 1 0 -1 -1 0 0;0 0 0 0 -1 -1 0 1 0; 0 0 0 0 0 0 -1 -1 1; b=20 5 0 0 0 0 0 0 0; Ab ans = 13.3453 6.4401 8.5420

25、3.3274 -1.1807 1.6011 1.7263 0.4204 2.1467 Page65 Exercise 12 A=1 2 3;4 5 6;7 8 0; left=sum(eig(A), right=sum(trace(A) left = 6.0000 right = 6 left=prod(eig(A), right=det(A) 原題有錯(cuò), (-1)n應(yīng)刪去 left = 27.0000 right = 27 fA=(A-p(1)*eye(3,3)*(A-p(2)*eye(3,3)*(A-p(3)*eye(3,3) fA = 1.0e-012 * 0.0853 0.1421 0

26、.0284 0.1421 0.1421 0 -0.0568 -0.1137 0.1705 norm(fA) f(A)范數(shù)接近0 ans = 2.9536e-013 Chapter 4Page84 Exercise 1(1) roots(1 1 1) (2) roots(3 0 -4 0 2 -1) (3) p=zeros(1,24); p(1 17 18 22)=5 -6 8 -5; roots(p) (4) p1=2 3; p2=conv(p1, p1); p3=conv(p1, p2); p3(end)=p3(end)-4; %原p3最后一個(gè)分量-4 roots(p3) Page84 Ex

27、ercise 2 fun=inline(x*log(sqrt(x2-1)+x)-sqrt(x2-1)-0.5*x); fzero(fun,2) Page84 Exercise 3 fun=inline(x4-2x); fplot(fun,-2 2);grid on; fzero(fun,-1),fzero(fun,1),fminbnd(fun,0.5,1.5) Page84 Exercise 4 fun=inline(x*sin(1/x),x); fplot(fun, -0.1 0.1); x=zeros(1,10);for i=1:10, x(i)=fzero(fun,(i-0.5)*0.0

28、1);end; x=x,-x Page84 Exercise 5 fun=inline(9*x(1)2+36*x(2)2+4*x(3)2-36;x(1)2-2*x(2)2-20*x(3);16*x(1)-x(1)3-2*x(2)2-16*x(3)2,x); a,b,c=fsolve(fun,0 0 0) Page84 Exercise 6 fun=(x)x(1)-0.7*sin(x(1)-0.2*cos(x(2),x(2)-0.7*cos(x(1)+0.2*sin(x(2); a,b,c=fsolve(fun,0.5 0.5) Page84 Exercise 7 clear; close; t

29、=0:pi/100:2*pi; x1=2+sqrt(5)*cos(t); y1=3-2*x1+sqrt(5)*sin(t); x2=3+sqrt(2)*cos(t); y2=6*sin(t); plot(x1,y1,x2,y2); grid on; 作圖發(fā)現(xiàn)4個(gè)解的大致位置，然后分別求解 y1=fsolve(x(1)-2)2+(x(2)-3+2*x(1)2-5,2*(x(1)-3)2+(x(2)/3)2-4,1.5,2) y2=fsolve(x(1)-2)2+(x(2)-3+2*x(1)2-5,2*(x(1)-3)2+(x(2)/3)2-4,1.8,-2) y3=fsolve(x(1)-2)2

30、+(x(2)-3+2*x(1)2-5,2*(x(1)-3)2+(x(2)/3)2-4,3.5,-5) y4=fsolve(x(1)-2)2+(x(2)-3+2*x(1)2-5,2*(x(1)-3)2+(x(2)/3)2-4,4,-4) Page84 Exercise 8(1) clear; fun=inline(x.2.*sin(x.2-x-2); fplot(fun,-2 2);grid on; 作圖觀察 x(1)=-2; x(3)=fminbnd(fun,-1,-0.5); x(5)=fminbnd(fun,1,2); fun2=inline(-x.2.*sin(x.2-x-2); x(2

31、)=fminbnd(fun2,-2,-1); x(4)=fminbnd(fun2,-0.5,0.5); x(6)=2 feval(fun,x) 答案: 以上x(chóng)(1)(3)(5)是局部極小，x(2)(4)(6)是局部極大，從最后一句知道x(1)全局最小， x(2)最大。 (2) clear; fun=inline(3*x.5-20*x.3+10); fplot(fun,-3 3);grid on; 作圖觀察 x(1)=-3; x(3)=fminsearch(fun,2.5); fun2=inline(-(3*x.5-20*x.3+10); x(2)=fminsearch(fun2,-2.5);

32、x(4)=3; feval(fun,x) (3) fun=inline(abs(x3-x2-x-2); fplot(fun,0 3);grid on; 作圖觀察 fminbnd(fun,1.5,2.5) fun2=inline(-abs(x3-x2-x-2); fminbnd(fun2,0.5,1.5) Page84 Exercise 9 close; x=-2:0.1:1;y=-7:0.1:1; x,y=meshgrid(x,y); z=y.3/9+3*x.2.*y+9*x.2+y.2+x.*y+9; mesh(x,y,z);grid on; 作圖觀察 fun=inline(x(2)3/9+

33、3*x(1)2*x(2)+9*x(1)2+x(2)2+x(1)*x(2)+9); x=fminsearch(fun,0 0) 求極小值 fun2=inline(-(x(2)3/9+3*x(1)2*x(2)+9*x(1)2+x(2)2+x(1)*x(2)+9); x=fminsearch(fun2,0 -5) 求極大值 Page84 Exercise 10 clear;t=0:24; c=15 14 14 14 14 15 16 18 20 22 23 25 28 . 31 32 31 29 27 25 24 22 20 18 17 16; p2=polyfit(t,c,2) p3=polyfi

34、t(t,c,3) fun=inline(a(1)*exp(a(2)*(t-14).2),a,t); a=lsqcurvefit(fun,0 0,t,c) 初值可以試探 f=feval(fun, a,t) norm(f-c) 擬合效果 plot(t,c,t,f) 作圖檢驗(yàn) fun2=inline(b(1)*sin(pi/12*t+b(2)+20,b,t); 原題修改f(x)+20 b=lsqcurvefit(fun2,0 0,t,c) figure f2=feval(fun2, b,t) norm(f2-c) 擬合效果 plot(t,c,t,f2) 作圖檢驗(yàn) Page84 Exercise 11

35、 fun=inline(1-x)*sqrt(10.52+x)-3.06*x*sqrt(1+x)*sqrt(5); x=fzero(fun, 0, 1) Page84 Exercise 12 r=5.04/12/100;N=20*12; x=7500*180 房屋總價(jià)格 y=x*0.3 首付款額 x0=x-y 貸款總額 a=(1+r)N*r*x0/(1+r)N-1) 月付還款額 r1=4.05/12/100;x1=10*10000; 公積金貸款 a1=(1+r1)N*r1*x1/(1+r1)N-1) x2=x0-x1 商業(yè)貸款 a2=(1+r)N*r*x2/(1+r)N-1) a=a1+a2 P

36、age84 Exercise 13 列方程th*R2+(pi-2*th)*r2-R*r*sin(th)=pi*r2/2 化簡(jiǎn)得sin(2*th)-2*th*cos(2*th)=pi/2 以下Matlab計(jì)算 clear;fun= inline(sin(2*th)-2*th*cos(2*th)-pi/2,th) th=fsolve(fun,pi/4) R=20*cos(th) Page84 Exercise 14 先在Editor窗口寫(xiě)M函數(shù)保存 function x=secant(fname,x0,x1,e) while abs(x0-x1)e, x=x1-(x1-x0)*feval(fnam

37、e,x1)/(feval(fname,x1)-feval(fname,x0); x0=x1;x1=x; end 再在指令窗口 fun=inline(x*log(sqrt(x2-1)+x)-sqrt(x2-1)-0.5*x); secant(fun,1,2,1e-8) Page84 Exercise 15 作系數(shù)為a,初值為xo,從第m步到第n步迭代過(guò)程的M函數(shù)： function f=ex4_15fun(a,x0,m,n) x(1)=x0; y(1)=a*x(1)+1;x(2)=y(1); if mm, plot(x(i),x(i),x(i+1),y(i-1),y(i),y(i); end e

38、nd hold off; M腳本文件 subplot(2,2,1);ex4_15fun(0.9,1,1,20); subplot(2,2,2);ex4_15fun(-0.9,1,1,20); subplot(2,2,3);ex4_15fun(1.1,1,1,20); subplot(2,2,4);ex4_15fun(-1.1,1,1,20); Page84 Exercise 16 設(shè)夾角t, 問(wèn)題轉(zhuǎn)化為 min f=5/sin(t)+10/cos(t) 取初始值pi/4, 計(jì)算如下 fun=(t)5/sin(t)+10/cos(t); t,f=fminsearch(fun, pi/4) t =

39、 0.6709 f = 20.8097 Page84 Exercise 17 提示：x(k+2)=f(x(k)=a2*x(k)*(1-x(k)*(1-a*x(k)*(1-x(k) 計(jì)算平衡點(diǎn)x |f(x)| ex4_18(0.9,1,20) ex4_18(-0.9,1,20) ex4_18(1.1,1,20) ex4_18(-1.1,1,20) Page84 Exercise 19 clear; close; x(1)=0; y(1)=0; for k=1:3000 x(k+1)=1+y(k)-1.4*x(k)2; y(k+1)=0.3*x(k); end plot(x(1000:1500),

40、y(1000:1500),+g);hold on plot(x(1501:2000),y(1501:2000),.b); plot(x(2001:2500),y(2001:2500),*y); plot(x(2501:3001),y(2501:3001),.r); Chapter 5Page101 Exercise 1 x=0 4 10 12 15 22 28 34 40; y=0 1 3 6 8 9 5 3 0; trapz(x,y) Page101 Exercise 2 x=0 4 10 12 15 22 28 34 40; y=0 1 3 6 8 9 5 3 0; diff(y)./di

41、ff(x) Page101 Exercise 3 xa=-1:0.1:1;ya=0:0.1:2; x,y=meshgrid(xa,ya); z=x.*exp(-x.2 -y.3); px,py = gradient(z,xa,ya); px Page101 Exercise 4 t=0:0.01:1.5; x=log(cos(t); y=cos(t)-t.*sin(t); dydx=gradient(y,x) x_1,id=min(abs(x-(-1);%找最接近x=-1的點(diǎn) dydx(id) Page101 Exercise 5(1)Fun=inline(1/(sqrt(2*pi).*exp

42、(-x.2./2);Quadl(fun,0,1)(2) fun=inline(exp(2*x).*cos(x).3); quadl(fun,0,2*pi) 或用trapz x=linspace(0,2*pi,100); y=exp(2*x).*cos(x).3; trapz(x,y) (3) fun=(x)x.*log(x.4).*asin(1./x.2); quadl(fun,1,3) 或用trapz x=1:0.01:3; y=feval(fun,x); trapz(x,y) (4) fun=(x)sin(x)./x; quadl(fun,1e-10,1) %注意由于下限為0，被積函數(shù)沒(méi)有

43、意義，用很小的1e-10代替 (5) %參考Exercise 5(4) (6) fun=inline(sqrt(1+r.2.*sin(th),r,th); dblquad(fun,0,1,0,2*pi) (7) 首先建立84頁(yè)函數(shù)dblquad2 clear; fun=(x,y)1+x+y.2; clo=(x)-sqrt(2*x-x.2); dup=(x)sqrt(2*x-x.2); dblquad2(fun,0,2,clo,dhi,100) Page101 Exercise 6 t=linspace(0,2*pi,100); x=2*cos(t);y=3*sin(t); dx=gradien

44、t(x,t);dy=gradient(y,t); f=sqrt(dx.2+dy.2); trapz(t,f) Page101 Exercise 7 xa=-1:0.1:1;ya=0:0.1:2; x,y=meshgrid(xa,ya); z=x.*exp(x.2+y.2); zx,zy=gradient(z,xa,ya); f=sqrt(1+zx.2+zy.2); s=0; for i=2:length(xa) for j=2:length(ya) s=s+(xa(i)-xa(i-1)*(ya(j)-ya(j-1)*(f(i,j)+f(i-1,j)+f(i,j-1)+f(i-1,j-1)/4;

45、 end end sPage101 Exercise 8 funl=inline(-(-x).0.2.*cos(x); funr=inline(x.0.2.*cos(x); quadl(funl,-1,0)+quadl(funr,0,1) Page101 Exercise 9 (以I32為例) fun=(x)abs(sin(x); h=0.1;x=0:h:32*pi;y=feval(fun,x);t1=trapz(x,y) h=pi;x=0:h:32*pi;y=feval(fun,x);t2=trapz(x,y)%步長(zhǎng)與周期一致，結(jié)果失真 q1=quad(fun,0,32*pi) q2=qua

46、dl(fun,0,32*pi) Page101 Exercise 10(2) 先在程序編輯器，寫(xiě)下列函數(shù)，保存為ex5_10_2f function d=ex5_10_2f(fname,a,h0,e) h=h0;d=(feval(fname,a+h)-2*feval(fname,a)+feval(fname,a-h)/(h*h); d0=d+2*e; while abs(d-d0)e d0=d;h0=h;h=h0/2; d=(feval(fname,a+h)-2*feval(fname,a)+feval(fname,a-h)/(h*h); end 再在指令窗口執(zhí)行 fun=inline(x.2

47、*sin(x.2-x-2),x); d=ex5_10_2f(fun,1.4,0.1,1e-3) Page101 Exercise 11 提示：f上升時(shí)，f0;f下降時(shí)，f0; f極值， f=0. Page101 Exercise 12 在程序編輯器，寫(xiě)下列函數(shù)，保存為ex5_12f function I=ex5_12(fname,a,b,n) x=linspace(a,b,n+1); y=feval(fname,x); I=(b-a)/n/3*(y(1)+y(n+1)+2*sum(y(3:2:n)+4*sum(y(2:2:n); 再在指令窗口執(zhí)行 ex5_12(inline(1/sqrt(2*

48、pi)*exp(-x.2/2),0,1,50) Page101 Exercise 13 fun=inline(5400*v./(8.276*v.2+2000),v); quadl(fun,15,30) Page101 Exercise 14 重心不超過(guò)凳邊沿。1/2, 2/3, 3/4, .,n/(n+1) Page101 Exercise15 利潤(rùn)函數(shù)fun=inline(p-c0+k*log(M*exp(-a*p)*M*exp(-a*p),p); 求p使fun最大 Page101 Exercise 16 clear; x=-3/4:0.01:3/4; y=(3/4+x)*2.*sqrt(1

49、-16/9.*x.2)*9.8; P=trapz(x,y) %單位：千牛 Page101 Exercise 17 clear; close; fplot(17-t(2/3)-5-2*t(2/3),0,20); grid; t=fzero(17-x(2/3)-5-2*x(2/3),7) t=0:0.1:8; y=17-t.(2/3)-5-2*t.(2/3); trapz(t,y)-20 %單位：百萬(wàn)元 Page101 Exercise 18 曲面面積計(jì)算 Chapter 6Page121 Exercise 1(1) fun=inline(x+y,x,y); t,y=ode45(fun,0 1 2

50、 3,1) %注意由于初值為y(0)=1,0 1 2 3中0不可缺 (3) 令y(1)=y,y(2)=y,化為方程組 y(1)=y(2),y(2)=0.01*y(2)2-2*y(1)+sin(t) 運(yùn)行下列指令 clear;close; fun=(t,y)y(2);0.01*y(2)2-2*y(1)+sin(t); t,y=ode45(fun,0 5,0;1); plot(t,y(:,1) (5) 令y(1)=y,y(2)=y,化為方程組 y(1)=y(2),y(2)=-mu*(y(1)2-1)*y(2)-y(1) 運(yùn)行下列指令,注意參數(shù)mu的處理 clear;close; fun=(t,y,

51、mu)y(2);-mu*(y(1)2-1)*y(2)-y(1); t,y=ode45(fun,0 20,2;0,1); plot(y(:,1),y(:,2);hold on; t,y=ode45(fun,0 20,2;0,2); plot(y(:,1),y(:,2),r);hold off; Page121 Excercise 2 roots(1 10 54 132 137 50) 通解A1*exp(-3*t)*cos(4*t)+A2*exp(-3*t)*sin(4*t)+A3*exp(-2*t)+A4*exp(-t)+A5*t*exp(-t) Page121 Excercise 3 dfun=inline(-1000.25*y(1)+999.75*y(2)+0.5;999.75*y(1)-1000.25*y(2)+0.5,x,y); x,y=ode45(dfun,0,50,1;-1);length(x) 所用節(jié)點(diǎn)很多 x,y=ode15s(dfun,0,50,1;-1);length(x) 所用節(jié)點(diǎn)很少 Page121 Excercise 4 clear; dfun=inline(x(2);2*x(3)+x(1)-(1-1/82.45)*(x(1)+1/82.45)/(