MATLAB 實用教程：第5章MATLAB程序設計

《MATLAB 實用教程：第5章MATLAB程序設計》由會員分享，可在線閱讀，更多相關(guān)《MATLAB 實用教程：第5章MATLAB程序設計（31頁珍藏版）》請在裝配圖網(wǎng)上搜索。

1、第5章 MATLAB程序設計 5.1腳本文件和函數(shù)文件 M文件有兩種形式：M腳本文件和M函數(shù)文件。 5.1.1 M文本編輯器 MATLAB的M文件是通過M文件編輯／調(diào)試器窗口(Editor／Debugger)來創(chuàng)建的。 圖5.1 M文件編輯/調(diào)試器窗口 單擊MATLAB桌面上的圖標，或者單擊菜單“File”——“New”——“M-file”，可打開空白的M文件編輯器，也可以通過打開已有的M文件來打開M文件編輯器。如圖5.1所示為打開已創(chuàng)建的M文件。 5.1.2 M文件的基本格式 下面介紹繪制二階系統(tǒng)時域曲線的M文件，欠阻尼系統(tǒng)的時域輸出y與x的關(guān)系為，【例5.1】為M腳本文件

2、，【例5.2】為M函數(shù)文件。 【例5.1】用M腳本文件繪制二階系統(tǒng)時域曲線。 %EX0501 二階系統(tǒng)時域曲線 %畫阻尼系數(shù)為0.3的曲線 x=0:0.1:20; y1=1-1/sqrt(1-0.3^2)*exp(-0.3*x).*sin(sqrt(1-0.3^2)*x+acos(0.3)) plot(x,y1,'r') 【例5.2】創(chuàng)建一個畫二階系統(tǒng)時域曲線的函數(shù)，阻尼系數(shù)zeta為函數(shù)的輸入?yún)?shù)。 function y=Ex0502(zeta) % EX0502 Step response of quadratic system. % 二階系

3、統(tǒng)時域響應曲線 % zeta 阻尼系數(shù) % y 時域響應 % % copyright 2003-08-01 x=0:0.1:20; y=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin(sqrt(1-zeta^2)*x+acos(zeta)) plot(x,y) M函數(shù)文件的基本格式： 函數(shù)聲明行 H1行(用%開頭的注釋行) 在線幫助文本(用%開頭) 編寫和修改記錄(用%開頭) 函數(shù)體? 例如，在命令窗口輸入help和lookfor命令查看幫助信息： help Ex0502 EX0502

4、 Step response of quadratic system. 二階系統(tǒng)時域響應曲線 zeta 阻尼系數(shù) y 時域響應 lookfor '二階系統(tǒng)時域響應' Ex0502.m: %二階系統(tǒng)時域響應 5.1.3 M腳本文件 腳本文件的特點： (1) 腳本文件中的命令格式和前后位置，與在命令窗口中輸入的沒有任何區(qū)別。 (2) MATLAB在運行腳本文件時，只是簡單地按順序從文件中讀取一條條命令，送到MATLAB命令窗口中去執(zhí)行。 (3) 與在命令窗口中直接運行命令一樣，腳本文件運行產(chǎn)生的變量都是駐留在MATLAB的工作空間(workspace)

5、中，可以很方便地查看變量，除非用clear命令清除；腳本文件的命令也可以訪問工作空間的所有數(shù)據(jù)，因此要注意避免變量的覆蓋而造成程序出錯。 【例5.1續(xù)】在M文件編輯／調(diào)試器窗口中編寫M腳本文件繪制二階系統(tǒng)的多條時域曲線。 (1) 單擊MATLAB桌面上的圖標打開M文件編輯器。 (2) 將命令全部寫入M文件編輯器中，為了能標志該文件的名稱，在第一行寫入包含文件名的注釋。保存文件為Ex0501.m。 %EX0501 二階系統(tǒng)時域曲線 x=0:0.1:20; y1=1-1/sqrt(1-0.3^2)*exp(-0.3*x).*sin(sqrt(1-0.3^2)*x+acos(0.

6、3)) plot(x,y1,'r') %畫阻尼系數(shù)為0.3的曲線 hold on y2=1-1/sqrt(1-0.707^2)*exp(-0.707*x).*sin(sqrt(1-0.707^2)*x+acos(0.707)) plot(x,y2,'g') %畫阻尼系數(shù)為0.707的曲線 y3=1-exp(-x).*(1+x) plot(x,y3,'b') %畫阻尼系數(shù)為1的曲線 圖5.2 運行界面 (3) 選擇M文件編輯器菜單“Debug”——“Run”，就可以在圖形窗中看到如圖5.2所示的曲線。 查看工作空間的變量： whos

7、 Name Size Bytes Class x 1x201 1608 double array y1 1x201 1608 double array y2 1x201 1608 double array y3 1x201 1608 double array Grand total is 804 elem

8、ents using 6432 bytes 5.1.4 M函數(shù)文件 函數(shù)文件的特點： (1) 第一行總是以“function”引導的函數(shù)聲明行； ?函數(shù)聲明行的格式： function [輸出變量列表] = 函數(shù)名(輸入變量列表)? (2) 函數(shù)文件在運行過程中產(chǎn)生的變量都存放在函數(shù)本身的工作空間； (3) 當文件執(zhí)行完最后一條命令或遇到“return”命令時，就結(jié)束函數(shù)文件的運行，同時函數(shù)工作空間的變量就被清除； (4) 函數(shù)的工作空間隨具體的M函數(shù)文件調(diào)用而產(chǎn)生，隨調(diào)用結(jié)束而刪除，是獨立的、臨時的，在MATLAB運行過程中可以產(chǎn)生任意多個臨時的函數(shù)空間。 【例5

9、.2續(xù)】在M文件編輯／調(diào)試器窗口編寫計算二階系統(tǒng)時域響應的M函數(shù)文件，并在MATLAB命令窗口中調(diào)用該文件。 創(chuàng)建M函數(shù)文件并調(diào)用的步驟如下： (1) 編寫函數(shù)代碼 function y=Ex0502(zeta) %EX0502 畫二階系統(tǒng)時域曲線 x=0:0.1:20; y=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin(sqrt(1-zeta^2)*x+acos(zeta)) plot(x,y) (2) 將函數(shù)文件保存為“Ex0502.m”。 (3) 在MATLAB命令窗口輸入以下命令，則會出現(xiàn)f的計算值和繪制的曲線： f=Ex0502(0

10、.3) 程序分析： § 第一行指定該文件是函數(shù)文件，文件名為“Ex0502”，輸入?yún)?shù)為阻尼系數(shù)zeta，輸出參數(shù)為時域響應y。 § 當函數(shù)文件調(diào)用結(jié)束，查看x、y： x ??? Undefined function or variable 'x'. y ??? Undefined function or variable 'y'. 注意：M腳本文件和M函數(shù)文件的文件名及函數(shù)名的命名規(guī)則與MATLAB變量的命名規(guī)則相同。 5.2程序流程控制 5.2.1 for ... end循環(huán)結(jié)構(gòu) 語法： for 循環(huán)變量=array 循環(huán)體 end 說明

11、：循環(huán)體被循環(huán)執(zhí)行，執(zhí)行的次數(shù)就是array的列數(shù)，array可以是向量也可以是矩陣，循環(huán)變量依次取array的各列，每取一次循環(huán)體執(zhí)行一次。 【例5.3】使用for ... end循環(huán)的array向量編程求出 1+3+5...+100 的值。 % EX0503 使用向量for循環(huán) sum=0; for n=1:2:100 sum=sum+n; end sum sum = 2500 計算的結(jié)果為：sum =2500。 程序說明：循環(huán)變量為n，n對應為向量1:2:100，循環(huán)次數(shù)為向量的列數(shù)，每次循環(huán)n取一個元素。 【例5.4】

12、使用for ... end循環(huán)的array矩陣編程將單位陣轉(zhuǎn)換為列向量。 % EX0504 使用矩陣for循環(huán) sum=zeros(6,1); for n=eye(6,6) sum=sum+n; end sum sum = 1 1 1 1 1 1 程序分析：循環(huán)變量n對應為矩陣eye(6,6)的每一列，即第一次n為[1;0;0;0;0;0]，第一次n為[0;1;0;0;0;0]；循環(huán)次數(shù)為矩陣的列數(shù)6。 5.2.2 while ... end循環(huán)結(jié)構(gòu) 語法： whil

13、e 表達式 循環(huán)體 end 說明：只要表達式為邏輯真，就執(zhí)行循環(huán)體；一旦表達式為假，就結(jié)束循環(huán)。表達式可以是向量也可以是矩陣，如果表達式為矩陣則當所有的元素都為真才執(zhí)行循環(huán)體，如果表達式為nan，MATLAB認為是假，不執(zhí)行循環(huán)體。 【例5.5】與【例5.3】相同，計算1+3+5...+100 的值。 % EX0505 使用while循環(huán) sum=0; n=1; while n<=100 sum=sum+n; n=n+2 ; end sum n sum = 2500 n = 101 程序分析

14、：可以看出while ... end循環(huán)的循環(huán)次數(shù)由表達式來決定，當n=101就停止循環(huán)。 5.2.3 If…else…end條件轉(zhuǎn)移結(jié)構(gòu) 語法： if 條件式1 ?? 語句段1 elseif 條件式2 ? 語句段2 ... else ?? 語句段n+1? end 說明：當有多個條件時，條件式1為假再判斷elseif的條件式2，如果所有的條件式都不滿足，則執(zhí)行else的語句段n+1，當條件式為真則執(zhí)行相應的語句段；If…else…end結(jié)構(gòu)也可以是沒有elseif和else的簡單結(jié)構(gòu)。 【例5.6】用If結(jié)構(gòu)執(zhí)行二階系統(tǒng)時域響應，根據(jù)阻尼

15、系數(shù)00)&(zeta<1) y=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin(sqrt(1-zeta^2)*x+acos(zeta)); elseif zeta==1 y=1-exp(-x).*(1+x); end plot(x,y) 5.2.4 switch…case開關(guān)結(jié)構(gòu) 語法： switch 開關(guān)表達式

16、 case 表達式1 ?? 語句段1? case表達式2 ?? 語句段2 ... otherwise? ?? 語句段n? end 說明： (1) 將開關(guān)表達式依次與case后面的表達式進行比較，如果表達式1不滿足，則與下一個表達式2比較，如果都不滿足則執(zhí)行otherwise后面的語句段n；一旦開關(guān)表達式與某個表達式相等，則執(zhí)行其后面的語句段。 (2) 開關(guān)表達式只能是標量或字符串。 (3) case后面的表達式可以是標量、字符串或元胞數(shù)組，如果是元胞數(shù)組則將開關(guān)表達式與元胞數(shù)組的所有元素進行比較，只要某個元素與開關(guān)表達式相等，就執(zhí)行其后的

17、語句段。 【例5.7】用switch…case開關(guān)結(jié)構(gòu)得出各月份的季節(jié)。 % EX0507 使用switch結(jié)構(gòu) for month=1:12; switch month case{3,4,5} season='spring' case{6,7,8} season='summer' case{9,10,11} season='autumn' otherwise season='winter' end end season = winter

18、 season = winter season = spring season = spring season = spring season = summer season = summer season = summer season = autumn season = autumn season = autumn season = winter 程序分析：開關(guān)表達式為向量1:12，case后面的表達式為元胞數(shù)組，當元胞數(shù)組的某個元素與開關(guān)表達式相等，就執(zhí)行其后的語句段。 5.2.5 try... catch... end試探結(jié)構(gòu)

19、語法： try ?? 語句段1 catch?? 語句段2 end 說明：首先試探性地執(zhí)行語句段1，如果在此段語句執(zhí)行過程中出現(xiàn)錯誤，則將錯誤信息賦給保留的lasterr變量，并放棄這段語句，轉(zhuǎn)而執(zhí)行語句段2中的語句，當執(zhí)行語句段2又出現(xiàn)錯誤，則終止該結(jié)構(gòu)。 【例5.8】用try... catch... end結(jié)構(gòu)來進行矩陣相乘運算。 % EX0508 try結(jié)構(gòu) n=4; a=magic(n); m=3; b=eye(3); try c=a*b catch c=a(1:m,1:m)*b end lasterr

20、 c = 16 2 3 5 11 10 9 7 6 ans = Error using ==> * Inner matrix dimensions must agree. 程序分析：試探出矩陣的大小不匹配時，矩陣無法相乘，則再執(zhí)行catch后面的語句段，將a的子矩陣取出與b矩陣相乘?？梢酝ㄟ^這種結(jié)構(gòu)靈活地實現(xiàn)矩陣的乘法運算。 5.2.6流程控制語句 1. break命令 break命令可以使包含break的最內(nèi)層的for或while語句強制終止，立即跳出該結(jié)構(gòu)，執(zhí)行end后面的命令，brea

21、k命令一般和If結(jié)構(gòu)結(jié)合使用。 【例5.9】將【例5.5】增加條件用If與break命令結(jié)合，停止while循環(huán)。計算1+3+5...+100 的值，當和大于1000時終止計算。 % EX0509 用break終止while循環(huán) sum=0; n=1; while n<=100 if sum<1000 sum=sum+n; n=n+2; else break end end sum n sum = 1024 n = 65 程序分析：while…

22、end循環(huán)結(jié)構(gòu)嵌套If…else…end分支結(jié)構(gòu)，當sum為1024時跳出while循環(huán)結(jié)構(gòu)，終止循環(huán)。 2. continue命令 continue命令用于結(jié)束本次for或while循環(huán)，只結(jié)束本次循環(huán)而繼續(xù)進行下次循環(huán)。 【例5.10】將If命令與continue命令結(jié)合，計算的1～100中所有素數(shù)的和，判斷是否為素數(shù)是將100以內(nèi)的每個數(shù)都被2～整除，不能被整除的就是素數(shù)。 % EX0510 用continue終止while循環(huán) sum=2;ss=0; for n=3:100 for m=2:fix(sqrt(n)) if mod(n

23、,m)==0 ss=1; %能被整除就用ss為1表示 break; %能被整除就跳出內(nèi)循環(huán) else ss=0; %不能被整除就用ss為0表示 end end if ss==1 continue; %能被整除就跳出本次外循環(huán) end sum=sum+n; end sum sum = 1060 程序分析：fix(s

24、qrt(n))是將取整；本程序為雙重循環(huán)，兩個for循環(huán)嵌套還嵌套一個if結(jié)構(gòu)；當mod(n,m)==0時就用break跳出判斷是否為素數(shù)的內(nèi)循環(huán)，并繼續(xù)用continue跳出求素數(shù)和的外循環(huán)而繼續(xù)下次外循環(huán)。 3. return命令 return命令是終止當前命令的執(zhí)行，并且立即返回到上一級調(diào)用函數(shù)或等待鍵盤輸入命令，可以用來提前結(jié)束程序的運行。 注意：當程序進入死循環(huán)，則按Ctrl+break鍵來終止程序的運行。 4. pause命令 pause命令用來使程序運行暫停，等待用戶按任意鍵繼續(xù)。 語法： pause %暫停 pause(n) %暫停n

25、秒 5. keyboard命令 keyboard命令用來使程序暫停運行，等待鍵盤命令，執(zhí)行完自己的工作后，輸入return語句，程序就繼續(xù)運行。 6. input命令 input命令用來提示用戶應該從鍵盤輸入數(shù)值、字符串和表達式，并接受該輸入。 a=input('input a number:') %輸入數(shù)值給a input a number:45 a = 45 b=input('input a number:','s') %輸入字符串給b input a number:45 b = 45 input('input a number:')

26、 %將輸入值進行運算 input a number:2+3 ans = 5 5.3函數(shù)調(diào)用和參數(shù)傳遞 5.3.1子函數(shù)和私有函數(shù) 1. 子函數(shù) 在一個M函數(shù)文件中，可以包含一個以上的函數(shù)，其中只有一個是主函數(shù)，其它則為子函數(shù)。 (1) 在一個M文件中，主函數(shù)必須出現(xiàn)在最上方，其后是子函數(shù)，子函數(shù)的次序無任何限制； (2) 子函數(shù)不能被其它文件的函數(shù)調(diào)用，只能被同一文件中的函數(shù)(可以是主函數(shù)或子函數(shù))調(diào)用； (3) 同一文件的主函數(shù)和子函數(shù)變量的工作空間相互獨立； (4) 用help和lookfor命令不能提供子函數(shù)的幫助信息。 【例5.11】將【例5.2】

27、畫二階系統(tǒng)時域曲線的函數(shù)作為子函數(shù)，編寫畫多條曲線的程序。 function Ex0511() % EX0511 使用函數(shù)調(diào)用繪制二階系統(tǒng)時域響應 z1=0.3; Ex0502(z1); %調(diào)用Ex0502 hold on z1=0.5 Ex0502(z1) %調(diào)用Ex0502 z1=0.707; Ex0502(z1) %調(diào)用Ex0502 function y=Ex0502(zeta) %子函數(shù)，畫二階系統(tǒng)時域曲線 x=0:0.1:20; y=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin(sqrt(1-ze

28、ta^2)*x+acos(zeta)) plot(x,y) 程序分析：主函數(shù)是Ex0511，子函數(shù)是Ex0502，在主函數(shù)中三次調(diào)用子函數(shù)。程序保存為Ex0511.m文件。 2. 私有函數(shù) 私有函數(shù)是指存放在private子目錄中的M函數(shù)文件，具有以下性質(zhì)： (1) 在private目錄下的私有函數(shù)，只能被其父目錄的M函數(shù)文件所調(diào)用，而不能被其它目錄的函數(shù)調(diào)用，對其它目錄的文件私有函數(shù)是不可見的，私有函數(shù)可以和其它目錄下的函數(shù)重名； (2) 私有函數(shù)父目錄的M腳本文件也不可調(diào)用私有函數(shù)； (3) 在函數(shù)調(diào)用搜索時，私有函數(shù)優(yōu)先于其它MATLAB路徑上的函數(shù)。 3. 調(diào)用

29、函數(shù)的搜索順序 在MATLAB中調(diào)用一個函數(shù)，搜索的順序如下： § 查找是否子函數(shù)； § 查找是否私有函數(shù)； § 從當前路徑中搜索此函數(shù)； § 從搜索路徑中搜索此函數(shù)。 5.3.2局部變量和全局變量 1. 局部變量 局部變量(Local Variables)是在函數(shù)體內(nèi)部使用的變量，其影響范圍只能在本函數(shù)內(nèi)，只在函數(shù)執(zhí)行期間存在。 2. 全局變量 全局變量(Global Variables)是可以在不同的函數(shù)工作空間和MATALB工作空間中共享使用的變量。 【例5.12】修改【例5.11】在主函數(shù)和子函數(shù)中使用全局變量。 function Ex0512() % E

30、X0512 使用全局變量繪制二階系統(tǒng)時域響應 global X X=0:0.1:20; z1=0.3; Ex0502(z1); hold on z1=0.5; Ex0502(z1); z1=0.707; Ex0502(z1); function Ex0502(zeta) %子函數(shù)，畫二階系統(tǒng)時域曲線 global X y=1-1/sqrt(1-zeta^2)*exp(-zeta*X).*sin(sqrt(1-zeta^2)*X+acos(zeta)); plot(X,y); 程序分析：X變量為全局變量，在需要使用的主函數(shù)和子函數(shù)中都需要用global定義

31、；同樣如果在工作空間中定義X為全局變量后也可以使用： global X who Your variables are: X 注意：由于全局變量在任何定義過的函數(shù)中都可以修改，因此不提倡使用全局變量；必須使用時應十分小心，建議把全局變量的定義放在函數(shù)體的開始，全局變量用大寫字符命名。 5.3.3函數(shù)的參數(shù) 1. 參數(shù)傳遞規(guī)則 function Ex0513() % EX0513 參數(shù)傳遞繪制二階系統(tǒng)時域響應 z1=0.3; [x1,y1]=Ex0502(z1); plot(x1,y1) hold on z1=0.5; [x2,y2]=Ex0502(z1);

32、plot(x2,y2) z1=0.707; [x3,y3]=Ex0502(z1); plot(x3,y3); function [x,y]=Ex0502(zeta) %子函數(shù)，畫二階系統(tǒng)時域曲線 x=0:0.1:20; y=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin(sqrt(1-zeta^2)*x+acos(zeta)); function [x,y]=Ex0502(zeta) %子函數(shù)，畫二階系統(tǒng)時域曲線 x=0:0.1:20; y=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin... (sqrt(

33、1-zeta^2)*x+acos(zeta)); 圖5.3 參數(shù)傳遞 【例5.13】將【例5.11】畫二階系統(tǒng)時域的函數(shù)修改，使用輸入輸出參數(shù)來實現(xiàn)參數(shù)傳遞，如圖5.3所示。 程序分析：主函數(shù)Ex0513調(diào)用子函數(shù)Ex0502，子函數(shù)中的zeta為輸入?yún)?shù)，函數(shù)調(diào)用時將z1傳遞給子函數(shù)zeta，子函數(shù)計算后將輸出參數(shù)x和y傳回給主函數(shù)的x1、y1；主函數(shù)調(diào)用子函數(shù)三次，后面兩次參數(shù)的傳遞也是同樣。 2. 函數(shù)參數(shù)的個數(shù) (1) nargin和nargout變量 函數(shù)的輸入輸出參數(shù)的個數(shù)可以通過變量nargin和nargout獲得，nargin用于獲得輸入?yún)?shù)的個數(shù)，nargou

34、t用于獲得輸出參數(shù)的個數(shù)。 語法： nargin %在函數(shù)體內(nèi)獲取實際輸入變量的個數(shù) nargout %在函數(shù)體內(nèi)獲取實際輸出變量的個數(shù) nargin(’fun’) %在函數(shù)體外獲取定義的輸入?yún)?shù)個數(shù) nargout(‘fun’) %在函數(shù)體外獲取定義的輸出參數(shù)個數(shù) 【例5.14】計算兩個數(shù)的和，根據(jù)輸入的參數(shù)個數(shù)不同使用不同的運算表達式。 function [sum,n]=Ex0514(x,y) % EX0514 參數(shù)個數(shù)可變，計算x和y的和 if nargin==1 sum=x+0; %輸入一個參數(shù)就計算與0的和 e

35、lseif nargin==0 sum=0; %無輸入?yún)?shù)就輸出0 else sum=x+y; %輸入的是兩個數(shù)則就計算和 end 在命令窗口調(diào)用Ex0514函數(shù)，分別使用2個、1個和無輸入?yún)?shù)結(jié)果如下： [y,n]=Ex0514(2,3) y = 5 n = 2 [y,n]=Ex0514(2) y = 2 n = 1 [y,n]=Ex0514 y = 0 n = 0 注意：如果輸入的參數(shù)多于輸入?yún)?shù)個數(shù)，則會出錯。 [y,

36、n]=Ex0514(1,2,3) ??? Error using ==> ex0514 Too many input arguments. 也可以在工作空間查看函數(shù)體定義的輸入?yún)?shù)個數(shù)： nargin('Ex0514') ans = 2 【例5.14續(xù)】添加以下程序，查看用nargout變量獲取輸出參數(shù)個數(shù)。 if nargout==0 %當輸出參數(shù)個數(shù)為0時，運算結(jié)果為0 sum=0; end 在命令窗口調(diào)用Ex0514函數(shù)，當輸出參數(shù)格式不同時，結(jié)果如下： Ex0514(2,3) %當輸出參數(shù)個數(shù)為0時 ans = 0

37、 y=Ex0514(2,3) %當輸出參數(shù)個數(shù)為1時 y = 5 [y,n,x]=Ex0514 %當輸出參數(shù)個數(shù)為太多時 ??? Error using ==> ex0514 Too many output arguments. 程序分析：當輸出參數(shù)個數(shù)為0時，即使有兩個輸入?yún)?shù)，運算結(jié)果也為0，結(jié)果送給ans變量；當輸出的參數(shù)個數(shù)太多，也會出錯。 (2) varargin和varargout變量 varargin和varargout可以獲得輸入輸出變量的各元素內(nèi)容。 【例5.15】計算所有輸入變量的和。 function [y

38、,n]=Ex0515(varargin) % EX0515 使用可變參數(shù)varargin if nargin==0 %當沒有輸入變量時輸出0 disp('No Input variables.') y=0; elseif nargin==1 %當一個輸入變量時，輸出該數(shù) y=varargin{1}; else n=nargin; y=0; for m=1:n y=varargin{m}+y; %當有多個輸入變量時，取輸入變量循環(huán)相加 end end n=nargin;

39、 在MATLAB的命令窗口中輸入不同個數(shù)的變量調(diào)用函數(shù)Ex0515，結(jié)果如下： [y,n]=Ex0515(1,2,3,4) %輸入4個參數(shù) y = 10 n = 4 [y,n]=Ex0515(1) %輸入1個參數(shù) y = 1 n = 1 [y,n]=Ex0515 %無輸入?yún)?shù) No Input variables. y = 0 n = 0 程序分析：n為輸入?yún)?shù)的個數(shù)；y為求和運算的結(jié)果。 5.3.4程序舉例 【例5.16】根據(jù)阻尼系數(shù)繪制

40、不同二階系統(tǒng)的時域響應，當欠阻尼時，當臨界阻尼時，當過阻尼時。 M文件的程序代碼如下： function y=Ex0516(z1) % EX0516 主函數(shù)調(diào)用子函數(shù)，根據(jù)阻尼系數(shù)繪制二階系統(tǒng)時域曲線 t=0:0.1:20; if (z1>=0)&(z1<1) y=plotxy1(z1,t); elseif z1==1 y=plotxy2(z1,t); else y=plotxy3(z1,t); end function y1=plotxy1(zeta,x) %畫欠阻尼二階系統(tǒng)時域曲線 y1=1-1/sqrt(1-zeta^2)*e

41、xp(-zeta*x).*sin(sqrt(1-zeta^2)*x+acos(zeta)); plot(x,y1) function y2=plotxy2(zeta,x) %畫臨界阻尼二階系統(tǒng)時域曲線 y2=1-exp(-x).*(1+x); plot(x,y2) function y3=plotxy3(zeta,x) %畫過阻尼二階系統(tǒng)時域曲線 y3=1-1/(2*sqrt(zeta^2-1))*(exp(-((zeta-sqrt(zeta^2-1))*x))./(zeta-sqrt(zeta^2-1))... -exp(-((zeta+sqrt(zeta^

42、2-1))*x))./(zeta+sqrt(zeta^2-1))); plot(x,y3) 程序分析：主函數(shù)名為Ex0516，三個子函數(shù)名分別為plotxy1、plotxy2、plotxy3，文件保存為Ex0516.m。當在命令窗口中輸入以下命令： y=Ex0516(0.3); hold on y=Ex0516(0.707); y=Ex0516(1); y=Ex0516(2); 則產(chǎn)生如圖5.4所示時域響應曲線。 圖5.4 不同阻尼系數(shù)的時域響應曲線 【例5.17】編寫M函數(shù)文件，通過流程控制語句，建立如下的矩陣。 functio

43、n y=Ex0517(m) % EX0517 用循環(huán)流程控制語句創(chuàng)建矩陣 y=0; m=m-1; for n=1:m y=[0,y]; %創(chuàng)建全0行 end for n=1:m a=[1:1:n]; b=a; for k=m:-1:n b=[0,b]; end y=[b;y]; n=n+1; end 程序分析：將矩陣的行列數(shù)用輸入?yún)?shù)m來確定，輸出參數(shù)為矩陣y。使用雙重循環(huán)來創(chuàng)建矩陣，將文件保存為Ex0517.m。 在命令窗口中調(diào)用Ex0517函數(shù)： y=Ex0517(5

44、) y = 0 1 2 3 4 0 0 1 2 3 0 0 0 1 2 0 0 0 0 1 0 0 0 0 0 5.4 M文件性能的優(yōu)化和加速 5.4.1 P碼文件 1. P碼文件的生成 P碼文件使用pcode命令生成，生成的P碼文件與原M文件名相同，其擴展名為“.p”。 語法： pcode Filename.m %在當前目錄生成Filen

45、ame.p pcode Filename.m -inplace %在Filename.m所在目錄生成Filename.p pcode Ex0517.m 則在當前目錄就生成了P碼文件Ex0517.p。 2. P碼文件的特點 (1) P碼文件的運行速度比原M文件速度快 (2) 存在同名的M文件和P碼文件時則P碼文件被調(diào)用 (3) P碼文件保密性好 用字處理軟件打開Ex0517.p文件，看到的是亂碼。 5.4.2 M文件性能優(yōu)化 1. 使用循環(huán)時提高速度的措施 循環(huán)語句及循環(huán)體是MATLAB編程的瓶頸問題，改進這種狀況有三種方法： (1) 盡量用向量的運算來

46、代替循環(huán)操作。 (2) 在必須使用多重循環(huán)的情況下，如果兩個循環(huán)執(zhí)行的次數(shù)不同，則建議在循環(huán)的外環(huán)執(zhí)行循環(huán)次數(shù)少的，內(nèi)環(huán)執(zhí)行循環(huán)次數(shù)多的，也可以顯著提高速度。 (3) 應用Mex技術(shù) 如果耗時的循環(huán)不可避免，就應該考慮用其他語言，如C或Fortran語言，按照Mex技術(shù)要求的格式編寫相應部分的程序，然后通過編譯聯(lián)接，形成在MATLAB可以直接調(diào)用的動態(tài)鏈接庫(DLL)文件，這樣就可以顯著地加快運算速度(在8.1.1小節(jié)介紹)。 2. 大型矩陣的預先定維 給大型矩陣動態(tài)地定維是個很費時間的事。 【例5.18】將【例5.17】中的雙重循環(huán)改為單循環(huán)，并先用zeros函數(shù)定維來提高

47、運行速度，創(chuàng)建矩陣。 function y=Ex0518(m) % EX0518 先定維再創(chuàng)建矩陣 m=m-1; y=zeros(m); for n=1:m-1 a=1:m-n; y(n,n+1:m)=a; end y 在命令窗口中分別調(diào)用Ex0517和Ex0518函數(shù)，比較其運行速度的不同： y=Ex0517(200) y=Ex0518(200) 3. 優(yōu)先考慮內(nèi)在函數(shù) 矩陣運算應該盡量采用MATLAB的內(nèi)在函數(shù)，因為內(nèi)在函數(shù)是由更底層的C語言構(gòu)造的，其執(zhí)行速度顯然很快。 4. 采用高效的算法 在實際應用中，解決同樣的數(shù)學問題經(jīng)

48、常有各種各樣的算法。因此，應尋求更高效的算法。 5. 盡量使用M函數(shù)文件代替M腳本文件 由于M腳本文件每次運行時，都必須把程序裝入內(nèi)存，然后逐句解釋執(zhí)行，十分費時。 5.4.3 JIT和加速器 1. JIT和加速器的加速范圍 JIT和加速器的應用范圍如下： § 只對維數(shù)不超過3的“非稀疏”數(shù)組起作用； § 只對“雙精度”、“整數(shù)”、“字符串”和“邏輯”等四種數(shù)據(jù)類型起作用； § 只對內(nèi)部函數(shù)的調(diào)用起作用，對用戶的M函數(shù)或MEX文件不起作用； § 只對控制語句for、while、if、elseif、switch中標量運算的條件表達式起作用； § 當一行程序中含有“不可加速

49、的”命令或變量時，整行都不被加速； § 當變量改變數(shù)據(jù)類型或維數(shù)，則該語句不被加速； § 如果i和j不以虛數(shù)單位形式使用，則該語句不被加速； § 在Intel系列CPU硬件上，Windows和Linux系統(tǒng)加速能力最強。 2. 使用程序性能剖析窗口 (1) 打開程序性能剖析窗口 指示燈 命令輸入欄 時間指示器 開始剖析 圖5.5 程序性能剖析窗口 選擇菜單“View”——“Profiler”；或使用在命令窗口輸入“profile viewer”命令都可以打開程序性能剖析窗口，如圖5.5所示。 (2) 對MATLAB的命令進行剖析 對MATLAB的命令進行剖析有兩

50、種方式： § 在上圖中的“命令輸入欄”中輸入需要剖析的命令，然后單擊“Start Profiling”按鈕，則指示燈變?yōu)榫G色，“Start Profiling”按鈕變灰色，時間指示器的計時在累加；當命令運行結(jié)束時，狀態(tài)恢復，并出現(xiàn)剖析報告。 § 將“命令輸入欄”清空，然后單擊“Start Profiling”按鈕開始剖析，則“Start Profiling”按鈕變?yōu)椤癝top Profiling”，指示燈變?yōu)榫G色表示啟動剖析，時間指示器的計時在累加；然后在MATLAB命令窗口中輸入需要剖析的命令，進行剖析；當命令運行結(jié)束，單擊“Stop Profiling”按鈕停止剖析，則狀態(tài)恢復，并出現(xiàn)

51、剖析報告。這種方式的時間指示器顯示的只是單擊按鈕“Stop Profiling”開始到單擊“Stop Profiling”按鈕的時間，并不表示命令的運行時間。 (3) 查看剖析報告 在程序性能分析窗口中以表格顯示“剖析分析匯總表”(Profile Summary)，從上到下按占用時間的多少排列。 將第三章的例題進行剖析，【例3.24】求微分方程組的解，保存為“Ex0324.m”文件，內(nèi)容如下： %符號微分方程組求解 [x,y]=dsolve('Dx=y,Dy=-x') [x,y]=dsolve('Dx=y,Dy=-x','t') 圖5.6 剖析報告 將該文件進行剖析，在“命令輸

52、入欄”中輸入“Ex0324”，然后單擊“Start Profiling”按鈕開始剖析，則“剖析分析匯總表”如圖5.6所示，單擊“Ex0324”超鏈接則顯示詳細的列表： 3. JIT和加速器的開關(guān)函數(shù) 在MATLAB6.5版總是把JIT和加速器置于開啟狀態(tài)，可以使用命令來控制JIT和加速器的開啟和關(guān)閉。 語法： feature JIT on %開啟JIT feature JIT off %關(guān)閉JIT feature accel on %開啟加速器 feature accel off %關(guān)閉加速器 將【例5.18】創(chuàng)建矩陣的程序在開啟或關(guān)閉JIT

53、和加速器時，查看其運行速度。 在圖5.6的“命令輸入欄”中輸入“y=Ex0518(200)” 當默認狀態(tài)時JIT和加速器置于開啟狀態(tài)，在程序性能剖析窗口中查看運行時間為0.43；當使用“feature JIT off”命令關(guān)閉JIT時，運行時間為0.451；當將JIT和加速器都關(guān)閉時則運行時間為0.481。 注意：JIT和加速器對于不同的M文件作用不同，而不同的運行環(huán)境、不同軟件、不同機器都會得出不同的運行時間，上面的運行時間只作為參考。 5.5內(nèi)聯(lián)函數(shù) 1. 內(nèi)聯(lián)函數(shù)的創(chuàng)建 創(chuàng)建內(nèi)聯(lián)函數(shù)可以使用inline命令實現(xiàn)。 語法： inline(‘string’,arg1,ar

54、g2,…) %創(chuàng)建內(nèi)聯(lián)函數(shù) 說明：’string’必須是不帶賦值號(“=”)的字符串；arg1和arg2是函數(shù)的輸入變量。 【例5.19】創(chuàng)建內(nèi)聯(lián)函數(shù)實現(xiàn)。 f=inline('sin(x)*exp(-z*x)','x','z') %創(chuàng)建內(nèi)聯(lián)函數(shù) f = Inline function: f(x,z) = sin(x)*exp(-z*x) y=f(5,0.3) %調(diào)用函數(shù)f y = -0.2140 2. 查看內(nèi)聯(lián)函數(shù) MATLAB可以用char、class和argnames命令方便地查看內(nèi)聯(lián)

55、函數(shù)的信息。 語法： char(inline_fun) %查看內(nèi)聯(lián)函數(shù)的內(nèi)容 class(inline_fun) %查看內(nèi)聯(lián)函數(shù)的類型 argnames(inline_fun) %查看內(nèi)聯(lián)函數(shù)的變量 【例5.19續(xù)】查看內(nèi)聯(lián)函數(shù)的信息。 char(f) ans = sin(x)*exp(-z*x) class(f) ans = inline argnames(f) ans = 'x' 'z' 3. 使內(nèi)聯(lián)函數(shù)適用于數(shù)組運算 內(nèi)聯(lián)函數(shù)的輸入變量不能是數(shù)組，但可以使用v

56、ectorize命令將內(nèi)聯(lián)函數(shù)適用于數(shù)組運算。 語法： vectorize(inline_fun) %使內(nèi)聯(lián)函數(shù)適用于數(shù)組運算 【例5.19續(xù)】使內(nèi)聯(lián)函數(shù)適用于數(shù)組運算。 ff=vectorize(f) %使內(nèi)聯(lián)函數(shù)f轉(zhuǎn)換為適合數(shù)組運算 ff = Inline function: ff(x,z) = sin(x).*exp(-z.*x) x=0:0.1:20; y=ff(x,0.3); 4. 執(zhí)行內(nèi)聯(lián)函數(shù) 內(nèi)聯(lián)函數(shù)還可以直接使用feval命令執(zhí)行。 語法： [y1，y2,…]=feval(inline_

57、fun,arg1,arg2…) 【例5.19續(xù)】執(zhí)行內(nèi)聯(lián)函數(shù)。 x=0:0.1:20; z=0:0.05:10; y=feval(ff,x,z) 5.6利用函數(shù)句柄執(zhí)行函數(shù) 5.6.1函數(shù)句柄的創(chuàng)建 1. 函數(shù)句柄的創(chuàng)建 語法： h_fun=@fun %創(chuàng)建函數(shù)句柄 h_fun=str2func(‘fun’) %創(chuàng)建函數(shù)句柄 h_array=str2func({‘fun1’,’fun2’,…}) %創(chuàng)建函數(shù)句柄數(shù)組 說明：fun是函數(shù)名，h_fun是函數(shù)句柄，h_array是函數(shù)句柄數(shù)組。 【例5.20】創(chuàng)建MATLA

58、B內(nèi)部函數(shù)的句柄。 h_sin=@sin; %創(chuàng)建函數(shù)句柄 h_cos=str2func('cos'); %創(chuàng)建函數(shù)句柄數(shù)組 h_array=str2func({'sin','cos','tan'}) h_array = @sin @cos @tan 2. 使用函數(shù)句柄的優(yōu)點 (1) 在更大范圍調(diào)用函數(shù) 函數(shù)句柄包含了函數(shù)文件的路徑和函數(shù)類型，即函數(shù)是否為內(nèi)部函數(shù)、M或P文件、子函數(shù)、私有函數(shù)等，因此無論函數(shù)所在的文件是否在搜索路徑上，是否是當前路徑，是否是子函數(shù)或私有函數(shù)，只要函數(shù)句柄存在，函

59、數(shù)就能執(zhí)行。 (2) 提高函數(shù)調(diào)用的速度 不使用函數(shù)句柄時，對函數(shù)的每次調(diào)用都要為該函數(shù)進行全面的路徑搜索，直接影響了速度。 (3) 使函數(shù)調(diào)用象使用變量一樣方便、簡單。 (4) 可迅速獲得同名重載函數(shù)的位置、類型信息。 5.6.2用feval命令執(zhí)行函數(shù) 函數(shù)也可以使用feval命令直接執(zhí)行，feval命令可以使用函數(shù)句柄或函數(shù)名。 語法： [y1,y2,…]=feval(h_fun,arg1,arg2…) [y1,y2,…]=feval(‘funname’,arg1,arg2…) 說明：h_fun是函數(shù)句柄，’funname’是函數(shù)名，arg1、arg2…是輸入

60、參數(shù)，y1、y2…是輸出參數(shù)。 【例5.21】將【例5.17】編寫的繪制二階系統(tǒng)時域響應曲線中的調(diào)用各子函數(shù)改為利用函數(shù)句柄實現(xiàn)。 function y=Ex0521(z1) % EX0521 利用函數(shù)句柄執(zhí)行函數(shù)，二階系統(tǒng)時域響應 t=0:0.1:20; h_plotxy1=str2func('plotxy1') %創(chuàng)建函數(shù)句柄 h_plotxy2=str2func('plotxy2') %創(chuàng)建函數(shù)句柄 h_plotxy3=str2func('plotxy3') %創(chuàng)建函數(shù)句柄 if (z1>=0)&(z1<1) y=feval(h_plotxy1,z1

61、,t); %執(zhí)行函數(shù) elseif z1==1 y=feval(h_plotxy2,z1,t); %執(zhí)行函數(shù) else y=feval(h_plotxy3,z1,t); %執(zhí)行函數(shù) end function y1=plotxy1(zeta,x) %畫欠阻尼二階系統(tǒng)時域曲線 y1=1-1/sqrt(1-zeta^2)*exp(-zeta*x).*sin(sqrt(1-zeta^2)*x+acos(zeta)); plot(x,y1) function y2=plotxy2(zeta,x) %畫臨界阻尼二階系統(tǒng)時域曲線 y2=1-exp(-x)

62、.*(1+x); plot(x,y2) function y3=plotxy3(zeta,x) %畫過阻尼二階系統(tǒng)時域曲線 y3=1-1/(2*sqrt(zeta^2-1))*(exp(-((zeta-sqrt(zeta^2-1))*x))./(zeta-sqrt(zeta^2-1))... -exp(-((zeta+sqrt(zeta^2-1))*x))./(zeta+sqrt(zeta^2-1))); plot(x,y3) 程序分析：在主函數(shù)中使用函數(shù)句柄執(zhí)行各子函數(shù)，運行的結(jié)果與【例5.17】相同。 在MATLAB的命令窗口調(diào)用該Ex0521函數(shù)有三種格式：

63、 (1) 用feval命令利用函數(shù)句柄執(zhí)行 h_Ex0521=str2func('Ex0521') h_Ex0521 = @Ex0521 y=feval(h_Ex0521,1); h_plotxy1 = @plotxy1 h_plotxy2 = @plotxy2 h_plotxy3 = @plotxy3 (2) 用feval命令利用函數(shù)名執(zhí)行 y=feval('Ex0521',1); (3) 直接調(diào)用函數(shù) y=Ex0521(1); 5.7利用泛函命令進行數(shù)值分析 在MATLA

64、B中，所有以函數(shù)為輸入變量的命令，都稱為泛函命令。 常見語法： [輸出變量列表]=函數(shù)名(h_fun,輸入變量列表) [輸出變量列表]=函數(shù)名(‘funname’,輸入變量列表) 說明：h_fun是要被執(zhí)行的M函數(shù)文件的句柄，或者是內(nèi)聯(lián)函數(shù)和字符串；‘funname’是M函數(shù)文件名。 5.7.1求極小值 1. fminbnd函數(shù) fminbnd函數(shù)用來計算單變量非線性函數(shù)的最小值。 語法： [x,y]=fminbnd(h_fun,x1,x2,options) [x,y]=fminbnd(‘funname’,x1,x2,options) 說明：h_fun是函

65、數(shù)句柄，’funname’是函數(shù)名，必須是單值非線性函數(shù)；options是用來控制算法的參數(shù)向量，默認值為0可省略；x是fun函數(shù)在區(qū)間x1

66、為圖中的圓點(0.6370，11.2528)，誤差小于10-4。 圖5.7 求humps函數(shù)最小值 2. fminsearch函數(shù) fminsearch函數(shù)是求多變量無束縛非線性最小值。 語法： x=fminsearch(h_fun,x0) x=fminsearch(‘funname’,x0) 說明：x0是最小值點的初始猜測值。 【例5.23】求著名的Banana測試函數(shù)f(x,y)=100(y-x2)2+(1-x)2的最小值，它的理論最小值是x=1，y=1。該測試函數(shù)有一片淺谷，很多算法都難以逾越。 fn=inline('100*(x(2)-x(1)^2)^2+(1-x(1))^2' ,'x') %用inline產(chǎn)生內(nèi)聯(lián)函數(shù),x和y用二元數(shù)組表示 fn = Inline function: fn(x) = 100*(x(2)-x(1)^2)^2+(1-x(1))^2 y=fminsearch(fn,[0.5,-1]) %從(0.5,-1)為初始值開始搜索求最小值 y = 1.