數學建模與matlab.ppt

1、數學實驗與Matlab,周曉陽 謝松法 梅正陽 華中科技大學出版社,周曉陽制作（第一版）,用計算機玩數學，獲得新的學習視野 學會Matlab軟件，對未來科研將大有幫助 數學實驗本身也是一種科研方法 努力培養(yǎng)自己的創(chuàng)新意識 最好與線性代數同步學習,實 驗 一,矩陣運算與Matlab命令,日常矩陣及其運算,矩陣應用實例: 欖球防護用品的生產管理,應用問題,一個工廠生產三種橄欖球用品 ： 防護帽、 墊肩、臀墊。 需要不同數量的：硬塑料 、 泡沫塑料 尼龍線 、 勞動力。 為監(jiān)控生產，管理者對它們之間的關系十分關心。 為把握這些量的關系，他列出下面的表,原料產品關系表,訂單,管理者接到四份訂單如上表所

2、示。 問應該如何計算每份訂單所需的原材料，以便組織生產？,將表格寫成矩陣形式,計 算,輸入下面Matlab指令 A=4 2 3;1 3 2;1 3 3;3 2 2, B=35 20 60 45;10 15 50 40;20 12 45 20 C=A*B 請自行計算觀看結果,Matlab基本指令,向量的創(chuàng)建和運算,1. 直接輸入向量,x1=1 2 4,x2=1,2,1,x3=x1 運行結果 x1 = 1 2 4 x2 = 1 2 1 x3 = 1 2 4,2.冒號創(chuàng)建向量,x1=3.4:6.7, x2=3.4:2:6.7 x3=2.6:-0.8:0 運算結果 x1 = 3.4000 4.4000

3、 5.4000 6.4000 x2 = 3.4000 5.4000 x3 = 2.6000 1.8000 1.0000 0.2000,3.生成線性等分向量,指令x=linspace(a,b,n) 在a,b區(qū)間產生 n 個等分點(包括端點) x=linspace(0,1,5) 結果 x = 0 0.2500 0.5000 0.7500 1.0000,工作空間,在Matlab窗口創(chuàng)建向量后并運行后，向量就存在于工作空間，可以被調用。,向量的運算,設x=x1 x2 x3; y=y1 y2 y3;為兩個三維向量，a,b為標量。 向量的數乘：a*x=a*x1 a*x2 a*x3 向量的平移： x+b=x

4、1+b x2+b x3+b 向量和： x+y=x1+y1 x2+y2 x3+y3 向量差： x-y=x1-y1 x2-y2 x3-y3 數的乘冪： 如 a2,元素群運算(四則運算）,x.*y=x1*y1 x2*y2 x3*y3 (元素群乘積) x./y=x1/y1 x2/y2 x3/y3 (元素群右除，右邊的y做分母) x.y=y1/x1 y2/x2 y3/x3 (元素群左除，左邊的x做分母) x.5=x15 x25 x35 (元素群乘冪) 2.x=2x1 2x2 2x3 (元素群乘冪) x.y=x1y1 x2y2 x3y3 (元素群乘冪),元素群運算（函數計算）,Matlab有許多內部函數，

5、可直接作用于向量產生一個同維的函數向量。 x=linspace(0,4*pi,100);（產生100維向量x） y=sin(x); (y也自動為100維向量) y1=sin(x).2; y2=exp(-x).*sin(x); 觀察結果,創(chuàng)建矩陣（數值矩陣的創(chuàng)建）,直接輸入法創(chuàng)建簡單矩陣。 A=1 2 3 4; 5 6 7 8; 9 10 11 12 B=-1.3,sqrt(3);(1+2)*4/5,sin(5);exp(2),6 觀察運行結果,創(chuàng)建矩陣（符號矩陣的創(chuàng)建）,用指令“syms”說明符號變量。 syms a11 a12 a13 a14 a21 a22 a23 a24 a31 a32

6、a33 a34 b11 b12 b13 b14 b21 b22 b23 b24 b31 b32 b33 b34 A1=a11 a12 a13 a14 ;a21 a22 a23 a24; a31 a32 a33 a34, B1=b11 b12 b13 b14 ;b21 b22 b23 b24; b31 b32 b33 b34 運行,矩陣的運算(矩陣的加減、數乘、乘積),C=A1+B1 D=A1-B1 syms c, cA=c*A1 A2=A1(:,1:3), B1 G=A2*B1,矩陣的運算(矩陣的加減、數乘、乘積),A, A_trans=A H=1 2 3 ; 2 1 0 ; 1 2 3 ,

7、K=1 2 3 ; 2 1 0 ; 2 3 1 h_det=det(H), k_det=det(K), H_inv=inv(H),K_inv=K-1,矩陣的運算(左除和右除),左除“ ”： 求矩陣方程AX=B的解；（ A 、B的行要保持一致） 解為 X=AB； 當A為方陣且可逆時有X=AB=inv(A)*B； 右除“ / ”： 求矩陣方程XA=B的解 （A 、B的列要保持一致） 解為 X=B/A ， 當A為方陣且可逆時有X=B/A=B*inv(A),矩陣的運算(左除和右除),求矩陣方程： 設A、B滿足關系式：AB2B+A,求B。 其中A=3 0 1; 1 1 0; 0 1 4。 解：有(A-2

8、I)BA 程序 ： A=3 0 1; 1 1 0;0 1 4; B=inv(A-2*eye(3)*A, B=(A-2*eye(3)A 觀察結果：,分塊矩陣（矩陣的標識）,1.矩陣元素的標識 ： A(i,j)表示矩陣A 的第 i 行 j 列的元素； 2.向量標識方式 A(vr,vc)： vr=i1,i2,ik、vc=j1,j2,ju分別是含有矩陣A的行號和列號的單調向量。 A(vr,vc)是取出矩陣A的第i1,i2,ik行與j1,j2,ju列交叉處的元素所構成新矩陣。,分塊矩陣（矩陣的標識）,取出A的1、3行和1、3列的交叉處元素構成新矩陣A1 程序 A=1 0 1 1 2;0 1 -1 2 3

9、; 3 0 5 1 0;2 3 1 2 1, vr=1, 3;vc=1, 3; A1=A(vr, vc) 觀察結果,分塊矩陣（矩陣的標識）,將A分為四塊，并把它們賦值到矩陣B中，觀察運行后的結果。 程序 A11=A(1:2,1:2),A12=A(1:2,3:5), A21=A(3:4,1:2),A22=A(3:4,3:5) B=A11 A12;A21 A22 結果,分塊矩陣（矩陣的修改和提取）,修改矩陣A，將它的第1行變?yōu)?。 程序： A=1 0 1 1 2;0 1 -1 2 3; 3 0 5 1 0;2 3 1 2 1 A(1,:)=0 0 0 0 0; A 刪除上面矩陣A的第1、3行。 程

10、序： A(1,3,:)= 結果,生成特殊矩陣,全1陣 ones(n), ones(m,n), ones(size(A) 全零陣： zeros(n) ，zeros(m,n), zeros(size(A) 常常用于對某個矩陣或向量賦0初值 單位陣： eye(n)，eye(m,n) 隨機陣： rand(m,n)， rand(n)=rand(n,n)用于隨機模擬，常和rand(seed,k)配合使用。,生成特殊矩陣,將rand指令運行多次，觀察結果。 程序： y1=rand(1,5), y2=rand(1,5), rand(seed,3), x1=rand(1,5), rand(seed,3), x2

11、=rand(1,5) 結果,常用矩陣函數,det(A) ： 方陣的行列式； rank(A)： 矩陣的秩； eig(A)： 方陣的特征值和特征向量； trace(A)： 矩陣的跡； rref(A)： 初等變換階梯化矩陣A svd(A)： 矩陣奇異值分解。 cond(A)： 矩陣的條件數；,數據的簡單分析,1.當數據為行向量或列向量時，函數對整個向量進行計算. 2.當數據為矩陣時，命令對列進行計算，即把每一列數據當成同一變量的不同觀察值。 max(求最大)、min(求最小)、mean(求平均值)、sum(求和)、std(求標準差)、cumsum(求累積和)、median(求中值)、diff(差分)

12、、sort(升序排列)、sortrows(行升序排列)等等。,數據的簡單分析,觀察：生成一個36的隨機數矩陣，并將其各列排序、求各列的最大值與各列元素之和。 程序 rand(seed,1);A=rand(3,6), Asort=sort(A), Amax=max(A), Asum=sum(A) 結果,實驗二,函數可視化與Matlab作圖,函數的可視化,f (x), g (x)是周期函數嗎？觀察它們的圖象。 程序 clf, x=linspace(0,8*pi,100); F=inline(sin(x+cos(x+sin(x); y1=sin(x+cos(x+sin(x); y2=0.2*x+si

13、n(x+cos(x+sin(x); plot(x,y1,k:,x,y2,k-) legend(sin(x+cos(x+sin(x),0.2x+sin(x+cos(x+sin(x),2),令,繪制平面曲線(plot指令),plot(x,y)： 以x為橫坐標、y為縱坐標繪制二維圖形 x,y是同維數的向量； plot(y)： 相當于x=1,2,length(y)時情形。,繪制平面曲線（繪制多個圖形),1. plot(x,y1;y2;), x是橫坐標向量，y1;y2;是由若干函數的縱坐標拼成的矩陣 2. plot(x,y1), hold on, plot(x,y2), hold off 3. plot

14、(x,y1,x,y2,) 4.plotyy 兩個坐標系，用于繪制不同尺度的函數。,繪制平面曲線（線型、點形和顏色的控制）,plot(x,y,顏色線型點形) plot(x,y,顏色線型點形,x,y,顏色線型點形, ) 句柄圖形和set命令改變屬性值，可套用： h=plot(x,y), set(h,屬性,屬性值,屬性,屬性值,) 也可用plot(x,y,屬性,屬性值)設置圖形對象的屬性。,繪制平面曲線（屬性變量和屬性值）,線寬：LineWidth 點的大?。?MarkerSize 線型：LineStyle 顏色：color,繪制平面曲線（例）,觀察： 改變繪圖的線型和顏色。 用grid on 指令

15、為圖形窗口加上 網格線，并改變網格的線型和字體的大小。 程序 h=plot(0:0.1:2*pi,sin(0:0.1:2*pi); set(h,LineWidth,5,color,red); grid on set(gca,GridLineStyle,-,fontsize,16) 觀察結果,繪制平面曲線（坐標軸的控制）,axis指令 axis(xmin xmax ymin ymax)： 設定二維圖形的x和y坐標的范圍； axis(xmin xmax ymin ymax zmin ymax)： 設定三維圖形的坐標范圍 ； 其中xminxxmax， yminyymax ，zminzzmax。,繪制

16、平面曲線（gca屬性控制）,改變當前軸對象句柄gca屬性 用set(gca,屬性,屬性值,)可改變字體大小、坐標刻度等軸對象的內容。例如： set(gca,ytick,-1 -0.5 0 0.5 1) 將 y 坐標按向量-1 -0.5 0 0.5 1將刻度分成4格； set(gca,yticklabel,a|b|c|d|e) 改變y坐標刻度的說明。,繪制平面曲線（gca屬性控制,例）,設置y坐標的刻度并加以說明，并改變字體的大小。 程序 plot(0:0.1:2*pi,sin(0:0.1:2*pi),k.-,);grid on, axis(0 6.3 -1.1 1.1)， set(gca,yt

17、ick,-1 -0.5 0 0.5 1), set(gca,yticklabel,a|b|c|d|e), set(gca,fontsize,20) get(gca) 運行結果,繪制平面曲線（文字標注）,title(圖形標題)； xlabel（x軸名稱）；ylabel（y軸名稱）；zlabel（z軸名稱）； text(說明文字)：創(chuàng)建說明文字； gtext（說明文字）：用鼠標在特定位置輸入文字。 文字標注常用符號： pi （）；alpha （）；beta （）； leftarrow （左箭頭） rightarrow （右箭頭）； bullet （點號）,繪制平面曲線（程序講解，exp2_1.m）

18、,clf, t=0:0.1:3*pi;alpha=0:0.1:3*pi; plot(t,sin(t),r-);hold on; plot(alpha,3*exp(-0.5*alpha),k:); set(gca,fontsize,15,fontname,times New Roman), xlabel(itt(deg); ylabel(itmagnitude); title( itsine wave and itAe-alphaittwave);,繪制平面曲線（程序講解，exp2_1.m）,text(6,sin(6),fontsize15The Value itsin(t) at itt=6r

19、ightarrowbullet, HorizontalAlignment,right), text(2,3*exp(-0.5*2), fontsize15bulletleftarrow The Value of it3e-0.5 itt=, num2str(3*exp(-0.5*2), at itt =2 ); legend(itsin(t),itAe-alphat) 注1： num2str： string1 ,num2str,string2，用方括號 注2： legend 請結合圖形觀察此命令的使用,圖形窗口的創(chuàng)建和分割,subplot(m,n,k)命令。 在圖形區(qū)域中顯示多個圖形窗口。 m

20、為上下分割數，n為左右分割數，k為第k子圖編號。 例：將一個圖形分為9個子圖，在第k個子圖畫sin(kx) 的圖象. 程序： clf，b=2*pi;x=linspace(0,b,50); for k =1:9 y=sin(k * x); subplot(3,3,k),plot(x,y),axis(0,2*pi,-1,1) end,若干有用的指令,clf：清除圖形窗口已有的內容. shg：顯示圖形窗口。 clear、 clear x：清除工作空間的已有變量。 figure(n): 打開第n個圖形窗口 help: : 續(xù)行號,繪制二元函數,基本步驟： 1.生成二維網格點 2. 計算函數在網格點上的

21、值 3. 繪制函數圖形,三維繪圖（ meshgrid指令：生成網格點）,觀察meshgrid指令的效果。 程序： a=-0.98;b=0.98;c=-1;d=1;n=10; x=linspace(a,b,n); y=linspace(c,d,n); X,Y=meshgrid(x,y); plot(X,Y,+) 觀察結果,三維繪圖（計算函數值，定義域裁減）,程序： for i=1:n for j=1:n if (1-X(i,j)eps1|X(i,j)-Y(i,j)eps1 z(i,j)=NaN; else z(i,j)=1000*sqrt(1-X(i,j)-1.*log(X(i,j)-Y(i,j

22、); end end end,三維繪圖（繪圖指令）,mesh(X,Y,z) ： 在三維空間中繪出由(X,Y,z)表示的曲面； meshz(X,Y,z)： 除了具有mesh的功能外，還畫出上下高度線， meshc(X,Y,z)： 除了具有mesh的功能外，還在曲面的下方畫出函數z=f(x,y)的等值線圖， surf(X,Y,z)： 也是三維繪圖指令，與mesh的區(qū)別在于mesh繪出彩色的線，surf繪出彩色的面， 運行exp2_1,觀察效果,三維繪圖（等值線指令）,表現二維函數的圖形的另一種方式是繪制等值線圖。 contour(X,Y,z,n)： n條等高線，n可缺??； contourf(X,Y

23、,z,n)： 等值線間用不同的顏色填滿，有更好的視覺效果； contour3(X,Y,z,n)： 在三維空間畫出等值線圖 colorbar： 將顏色與函數值對應起來顯示在圖中。,三維繪圖（等值線指令，繼續(xù)exp2_2顯示效果）,clf,contour(X,Y,z,40),colorbar contourf(X,Y,z,40),colorbar contour3(X,Y,z,40),colormap(0,0,0) 為等值線標上函數值: 可套用下面程序的格式. cs,h=contour(X,Y,z,15); clabel(cs,h,labelspacing,244) labelspace是數值標記

24、之間相隔的寬度，默認值為144, 這里取了244，,空間曲線和運動方向的表現,一條空間曲線可以用矢量函數表示為,它的速度矢量表現為曲線的切矢量:,空間曲線和運動方向的表現,很顯然飛行曲線方程為：,繪制空間曲線（指令）,plot3(x,y,z)： 繪制三維空間曲線，用法和plot類似。 quiver(X,Y,u,v)：繪制二維矢量， 在坐標矩陣點X,Y處繪制矢量u,v, 其中u為矢量的x坐標，v為矢量的y 坐標，其維數不小于2。 quiver3（X,Y,Z,u,v,w）： 繪制三維矢量，用法與quiver類似。 Gradient： Fx,Fy,Fz=gradient(F)為函數F數值梯度,繪制空

25、間曲線（程序講解exp2_3）,exp2_3.m clf,t=linspace(0,1.5,20); x=t.2;y=(2/3)*t.3;z=(6/4)*t.4-(1/3)*t.3; plot3(x,y,z,r.-,linewidth,1,markersize,10),hold on Vx=gradient(x);Vy=gradient(y);Vz=gradient(z); h=quiver3(x,y,z,Vx,Vy,Vz),set(h,linewidth,1),grid on axis(0 1.5 0 1.5 0 40) xlabel(x),ylabel(y),zlabel(z)，box on 運行程序,應用、思考和練習,用平行截面法討論由曲面z=x2-y2構成的馬鞍面形狀。 對于二重積分，積分指示線方法是很有用的，當然你首先得了解一下什么是積分指示線法，請查閱高等數學相關的內容，然后設計一個數學實驗，然后用Matlab的繪圖工具表現這一方法。,應用、思考和練習,繪制微分方程 y/dx=xy, y(0)=0.4的斜率場， 并將解曲線畫在圖中，觀察斜率場和解曲線的關系。,應用、思考和練習,地球表面的氣溫差異很大，而且隨時間變化，