第4章 Matlab程序設(shè)計

2023/1/31MATLAB程序設(shè)計吳鵬（rocwoods）rocwoods@126.comMATLAB從零到進階2023/1/31主要內(nèi)容

M文件——腳本文件和函數(shù)文件

MATLAB程序流程控制與調(diào)試程序流程控制示例

匿名函數(shù)、子函數(shù)與嵌套函數(shù)

匿名函數(shù)、子函數(shù)與嵌套函數(shù)應(yīng)用案例

編寫高效的MATLAB程序

養(yǎng)成良好的編程風格

2023/1/31第一節(jié)M文件——腳本文件和函數(shù)文件

2023/1/31一、腳本文件1.什么時候適合用腳本文件2.什么是腳本文件

指令的增加并且需要頻繁重復計算時，直接從指令窗進行計算很麻煩。這時候使用腳本文件最為適宜

。

所謂腳本文件，就是一些按用戶意圖排列的MATLAB指令（包括流程控制指令在內(nèi)）組成的M文件。3.腳本文件和函數(shù)文件的區(qū)別

腳本文件運行后，產(chǎn)生的所有變量都駐留在MATLAB基本工作空間(baseworkspace)中。只要用戶不使用clear指令加以清除，或MATLAB不關(guān)閉，這些變量將一直保存在基本工作空間。2023/1/31二、函數(shù)文件其中out1,out2,…為輸出參數(shù)列表，in1,in2,…為輸入?yún)?shù)列表，funname為函數(shù)名。輸入?yún)?shù)和輸出參數(shù)個數(shù)根據(jù)問題的需要可以為0個到多個。函數(shù)文件也需要保存到MATLAB的搜索路徑下才能被調(diào)用。

格式：每個函數(shù)都由“function”引導，如下：function[out1,out2,…]=funname(in1,in2,…)注釋說明部分（%號引導的行）函數(shù)體

2023/1/31第二節(jié)MATLAB程序流程控制與調(diào)試

2023/1/31一、MATLAB程序流程控制if-else-end結(jié)構(gòu)switch-case-otherwise-end結(jié)構(gòu)for循環(huán)while循環(huán)try-catch結(jié)構(gòu)其他中斷、暫停語句2023/1/31二、MATLAB程序調(diào)試

2.設(shè)置斷點進行調(diào)試

1.語法錯誤和運行結(jié)果錯誤

語法錯誤即所書寫的代碼不符合MATLAB語法規(guī)范所造成的錯誤，比如由于粗心造成的拼寫錯誤，不十分了解某個函數(shù)的調(diào)用方法而造成的調(diào)用錯誤等

。運行結(jié)果錯誤是一類非常難以檢查的錯誤，程序能正常運行，只是運行結(jié)果與期望的不一樣，這類錯誤大多是由于算法錯誤引起的。2023/1/31在程序編輯窗口中編寫的M代碼的每一行的前面都標有行號，在可執(zhí)行的命令行的行號后面都有一個小短橫“-”，單擊某一行的“-”，就可在該行設(shè)置一個斷點，此時的“-”變成了紅色的圓點，表示設(shè)置斷點成功。將光標放到M代碼的某一行上，然后按快捷鍵F12，或者單擊程序編輯窗口工具欄上的圖標，也可以通過菜單項“Debug→Set/ClearBreakpoint”設(shè)置斷點。用函數(shù)dbstop設(shè)置斷點。2023/1/31第三節(jié)程序流程控制示例

2023/1/31一、概述

要通過一個完整的MATLAB程序來展示各種流程控制的用法，該程序來源于如下一個有趣的問題。一只失明的小貓不幸掉進山洞里，山洞有三個門，其中第一個門進去后走兩個小時后可以回到地面，第二個門進去后走四個小時又回到原始出發(fā)點，不幸的是第三個門進去后走六個小時還是回到原始出發(fā)點。小貓由于眼睛失明，每次都是隨機的選擇其中一個門走。問題是這只可憐的小貓走出山洞的期望時間是多少？2023/1/31二、問題分析

設(shè)小貓走出山洞的期望時間為t,如果小貓不幸進了第二個或第三個門，那么它過四個小時或是六個小時后又和進門之前面臨狀況一樣了，只不過這兩種不幸的情況發(fā)生的概率都為1/3，而萬幸一次性走出去的概率也是1/3。于是我們可以得到下面的方程：

t=2*1/3+(4+t)*1/3+(6+t)*1/3

很簡單的一元一次方程，解得t=12。1.解析分析方法2.計算機模擬的思路2023/1/31輸入正整數(shù)n作為模擬小貓出洞的次數(shù)，生成一個1*n的數(shù)組T用來記錄每次小貓出洞時間，初始值為0。k從1循環(huán)到n，T(k)用來記錄每次實際出洞的時間。隨機等概率的生成{1,2,3}之間的一個數(shù)c，如果c=1，T(k)=T(k)+2，小貓走出山洞，開始下一次模擬。否則根據(jù)c=2或者c=3決定T(k)=T(k)+4或T(k)=T(k)+6，并繼續(xù)隨機生成c，直到c=1。模擬完n次后，計算T的均值得到小貓走出山洞的期望時間的近似值2023/1/31三、程序代碼

functionT=cat_in_holl(varargin)%varargin,使函數(shù)可以接受參數(shù)個數(shù)不定的輸入if~isempty(varargin)%輸入?yún)?shù)非空n=varargin{1};%varargin為cell型數(shù)組，取其第一個元素賦給nend%try-catch結(jié)構(gòu)用法示例try%如果n是正整數(shù)下面語句不會發(fā)生錯誤，進而執(zhí)行try-catch結(jié)構(gòu)之后的語句

%否則會發(fā)生錯誤,執(zhí)行由catch引導的語句。

ifn>0&&mod(n,1)==0;%n為正整數(shù)的判斷條件

%空語句，不會報錯

else%n不是正整數(shù)，報錯

error;end2023/1/31catchME%ME,用來記錄發(fā)生錯誤的一些信息

disp('函數(shù)沒有輸入?yún)?shù)或者輸入?yún)?shù)不是正整數(shù)標量');T=[];%給T賦空值

return;%函數(shù)返回，后面語句不再執(zhí)行，沒有return會接著執(zhí)行后面的語句end%switch-case-end結(jié)構(gòu)用法示例switchnargin%nargin，函數(shù)輸入?yún)?shù)的個數(shù)

case1%case1的情況是我們程序的核心部分，即整個模擬出洞的過程

T=zeros(1,n);fork=1:n%for循環(huán)用法示例

c=unidrnd(3,1);%等概率隨機生成{1，2，3}中某個數(shù)字2023/1/31whilec~=1%while循環(huán)用法示例

ifc==2T(k)=T(k)+4;elseT(k)=T(k)+6;endc=unidrnd(3,1);endT(k)=T(k)+2;endcase2T=[];disp('函數(shù)只能有一個輸入?yún)?shù)，且為正整數(shù)');otherwiseT=[];disp('函數(shù)輸入?yún)?shù)個數(shù)不能大于1，參數(shù)需為正整數(shù)');end2023/1/31四、程序驗證

在Editor中新建一個m文件，將上述程序復制過去，并以cat_in_holl.m文件名保存到MATLAB的搜索路徑下，對于n=10^4以及n=10^5，我們可以用如下指令得到近似期望時間T：

>>n=1e4; >>T=cat_in_holl(n); >>mean(T) ans= 11.8888 >>T=cat_in_holl(n); >>mean(T) ans= 12.02222023/1/31cat_in_holl.m是一個完整的m函數(shù)文件，里面用到了幾乎各種流程控制方法。讀者朋友可以試運行如下指令，觀看程序結(jié)果，體會各流程控制的功能：

T=cat_in_holl;%沒有輸入?yún)?shù)

T=cat_in_holl(10.1);%輸入?yún)?shù)為小數(shù)

T=cat_in_holl(-10000);%輸入?yún)?shù)為負數(shù)

T=cat_in_holl({1000});%輸入為元胞數(shù)組

T=cat_in_holl(1000,100);%輸入?yún)?shù)為2個

T=cat_in_holl(1000,100,100);%輸入?yún)?shù)大于2個2023/1/31第四節(jié)匿名函數(shù)、子函數(shù)與嵌套函數(shù)

2023/1/31一、匿名函數(shù)1. 匿名函數(shù)定義

fhandle=@(arglist)expr

其中，expr是具體的函數(shù)表達式，arglist是指定的函數(shù)的自變量。 匿名函數(shù)按照不同的分類方法可以分為不同的種類。按照自變量的個數(shù)以及層數(shù)可以分為以下幾種：單變量匿名函數(shù)，多變量匿名函數(shù)，單重匿名函數(shù)，多重匿名函數(shù)。2. 單變量匿名函數(shù)

f=@(x)x.^2;3. 多變量匿名函數(shù)

g=@(x,y)x.^2+y.^2;2023/1/314. 單重匿名函數(shù) 無論單變量還是多變量匿名函數(shù)，都屬于單重匿名函數(shù)。這類匿名函數(shù)的特點是，只有一個“@”符號引導，“@”之后就是具體的函數(shù)表達式。自變量輸入單重匿名函數(shù)后，得到的是具體的數(shù)值。除了單重匿名函數(shù)外，還有雙層、乃至多重匿名函數(shù)。這些多重匿名函數(shù)在參數(shù)傳遞方面非常方便。5. 多重匿名函數(shù)

f=@(a,b)@(x)a*x+b f= @(a,b)@(x)a*x+b

其中“a,b”是外層變量，“x”是內(nèi)層變量?？梢赃@樣理解這個表達式：每個“@”符號后面括號里的變量的作用域一直到表達式的結(jié)尾。這樣，“a,b”的作用域就是“@(x)a*x+b”，而“x”的作用域就是“a*x+b”。因此，對于給定的“a,b”，“gab=f(a,b)”是一個單層以x為變量的單變量匿名函數(shù)。2023/1/31二、子函數(shù)函數(shù)M文件中，我們稱第一個“function”引導的函數(shù)為主函數(shù)，外部調(diào)用函數(shù)文件時總是從主函數(shù)開始執(zhí)行。如果整個函數(shù)文件只有一個函數(shù)，那這個函數(shù)就是主函數(shù)。在設(shè)計比較復雜的程序的時候，為了能夠使程序之間邏輯關(guān)系清楚，易于閱讀和維護，通常會采用模塊化的設(shè)計方式，這時就會涉及子函數(shù)的使用。從函數(shù)編寫格式以及形態(tài)上來說，子函數(shù)和主函數(shù)并無區(qū)別，區(qū)別僅在于在函數(shù)文件中的位置以及調(diào)用關(guān)系。2023/1/31子函數(shù)是在同一函數(shù)文件中，主函數(shù)之后的由“function”引導的函數(shù)，一個函數(shù)文件可以有多個子函數(shù)。這些子函數(shù)在函數(shù)文件中的排列順序可以隨意，前提是都要位于主函數(shù)后。子函數(shù)可以被其所在的函數(shù)文件中的主函數(shù)和其他子函數(shù)調(diào)用，子函數(shù)也可以調(diào)用主函數(shù)，在設(shè)計子函數(shù)調(diào)用主函數(shù)的程序時一定要小心，避免形成“死”的調(diào)用循環(huán)。

2023/1/31三、嵌套函數(shù)在程序編寫中，有的數(shù)學表達式比較長，不方便用一行來表示，而且這個表達式里還有其他的參數(shù)。這時候我們就需要將其寫成函數(shù)，用嵌套函數(shù)可以方便的解決參數(shù)共享問題。當求解微分方程或者進行優(yōu)化問題求解時，一些復雜的表達式往往涉及參數(shù)傳遞過程，這個時候，比較適合將其表達式寫成嵌套函數(shù)。在編寫GUI用戶界面的回調(diào)函數(shù)時，參數(shù)傳遞往往是個棘手的問題。利用嵌套函數(shù)，可以方便解決這個問題。2023/1/311. 嵌套函數(shù)定義

嵌套函數(shù)，即nestedfunction，顧名思義，是嵌套在函數(shù)體內(nèi)部的函數(shù)。嵌套函數(shù)以function聲明，結(jié)束的時候加上end。需要說明的是，包含有嵌套函數(shù)的函數(shù)，無論他是主函數(shù)、子函數(shù)或者嵌套函數(shù)，都應(yīng)該在末尾的地方加上end。下面代碼是一個簡單的嵌套函數(shù)的例子：functionr=MyTestNestedFun(input)a=5;c=sin(input)+tan(input);functiony=nestedfun(b)y=a*c+b;endr=nestedfun(5);end

2023/1/312. 嵌套函數(shù)種類嵌套函數(shù)可以分為單重嵌套函數(shù)和多重嵌套函數(shù)。單重嵌套函數(shù)，這樣的函數(shù)嵌套在別的函數(shù)體內(nèi)，自己內(nèi)部不再有嵌套的函數(shù)，例如上面介紹的是單重嵌套函數(shù)。多重嵌套函數(shù)嵌套在別的函數(shù)體內(nèi)，同時自己內(nèi)部又嵌套著別的另一層或幾層函數(shù)。如下面類型的函數(shù)：2023/1/31functionx=A(p1,p2)...functiony1=B1(p3)...functionz1=C1(p4)...end...end...functiony2=B2(p5)...functionz2=C2(p6)...functionw=D(p7)...endendendend

2023/1/313. 嵌套函數(shù)的變量作用域變量的作用域指的變量能夠被程序訪問、修改、設(shè)置等的代碼范圍。

父函數(shù)和嵌套在其內(nèi)的函數(shù)，他們各自的變量是可以互相訪問的。但是必須注意的是，嵌套函數(shù)訪問父函數(shù)的變量，可以在函數(shù)定義里直接拿過來用，而父函數(shù)訪問嵌套在其內(nèi)的函數(shù)里的變量則必須要經(jīng)過調(diào)用之后。彼此沒有嵌套關(guān)系的嵌套函數(shù)間是不能簡單得共享變量的。如果非要共享，則只能通過它們所在的父函數(shù)來進行。2023/1/31四、嵌套函數(shù)彼此調(diào)用關(guān)系將主函數(shù)看成“父親”，嵌套函數(shù)依據(jù)嵌套深度可以看成“兒子”，“孫子”，“重孫”等等。這樣最外層的主函數(shù)可以認為是所有內(nèi)部函數(shù)的“祖先”，位于不同嵌套函數(shù)里的函數(shù)之間的關(guān)系類似于一個家族之間的“親屬關(guān)系”。嵌套函數(shù)之間的調(diào)用關(guān)系可以總結(jié)成下面幾句話：父親可以求助兒子，兒子可以求助父親，也即父子可以互相求助。一個人不能求助孫子、重孫等后代，但是可以求助自己的直系祖宗（祖父、曾祖父等等）以及和其直系祖宗是親兄弟的先人。一個人可以求助自己的親兄弟，或者親叔叔、親伯伯但不能求助侄兒。

2023/1/31第五節(jié)匿名函數(shù)、子函數(shù)與嵌套函數(shù)應(yīng)用案例2023/1/31一、匿名函數(shù)應(yīng)用實例1. 匿名函數(shù)在求解方程中的應(yīng)用【例4.5-1】求下列方程的根：求解代碼如下：>>f=@(x)exp(x)+x^2+x^(sqrt(x))-100%構(gòu)造方程的匿名函數(shù)形式f=@(x)exp(x)+x^2+x^(sqrt(x))-100>>formatlong>>x0=fzero(f,3)%求方程的根，初始值為3x0=4.163549956946139%求出的解>>f(x0)%代入原方程驗證ans=2.842170943040401e-0142023/1/31一、匿名函數(shù)應(yīng)用實例【例4.5-1序】對于 ，求下列方程相應(yīng)的值，并畫出和相應(yīng)的的圖像。

求解代碼如下：f=@(a)@(x)exp(x)+x^a+x^(sqrt(x))-100;%構(gòu)造函數(shù)句柄formatlongaa=0:0.01:2;plot(aa,arrayfun(@(a)fzero(f(a),4),aa),'*-')%利用arrayfun求解不同的a對應(yīng)的xxlabel('$a$','interpreter','latex','fontsize',15)%標注x、y坐標軸，按照latex語法ylabel('$x$','interpreter','latex','fontsize',15)title('$\mathrm{e}^{x}+x^{\sqrt{x}}+x^a-100$','interpreter','latex',...'fontsize',15)

2023/1/31【例4.5-1序】得到的圖形2023/1/31一、匿名函數(shù)應(yīng)用實例2. 匿名函數(shù)在顯式表示隱函數(shù)方面的應(yīng)用 隱函數(shù)一般無法在數(shù)學上顯式的表達。這里說的顯式表示指的是構(gòu)造一個MATLAB函數(shù)來表達隱函數(shù)，具體思路是對于給定的隱函數(shù)的自變量x，通過數(shù)值方法求解出因變量y，這樣就相當于顯式表達隱函數(shù)。請看下面的例子：

【例4.5-2】顯式表示下列y關(guān)于x的隱函數(shù)： 利用匿名函數(shù)，我們可以在MATLAB中顯式得寫出y和x的關(guān)系式如下：

y=@(x)fzero(@(y)(exp(y)+x^y)^(1/y)-x^2*y,1);2023/1/31一、匿名函數(shù)應(yīng)用實例 【例4.5-3】顯式表示下列z關(guān)于x，y的隱函數(shù)：

z=@(x,y)fzero(@(z)z-sin((z*x-0.5)^2+x*2*y^2-z/10)*...exp(-((x-0.5-exp(-y+z))^2+y^2-z/5+3)),rand);z(1,2)ans=7.369600885222808e-004z(0,0.5)ans=0.0028075073342922023/1/31一、匿名函數(shù)應(yīng)用實例3. 匿名函數(shù)在求積分區(qū)域方面的應(yīng)用

有的時候，我們需要根據(jù)已知的積分值和被積函數(shù)求對應(yīng)的積分區(qū)域，當被積函數(shù)表達式不是很復雜的時候，用匿名函數(shù)比較適合。請看下例：

【例4.5-4】要使sin(x)^2./(x^2)這個式子的積分值為0.99*pi，求其關(guān)于0對稱的積分區(qū)域。

u0=fzero(@(u)0.99*pi/2-quadl(@(x)sin(x).^2./(x.^2),0,u),1) u0= 32.3138103582174102023/1/31一、匿名函數(shù)應(yīng)用實例4. 匿名函數(shù)和符號計算的結(jié)合

【例4.5-5】求下面函數(shù)三階導數(shù)在[0,1]區(qū)間的圖像：

本例如果用手算的話，比較繁瑣。如果用符號運算得到三階導數(shù)的解析表達式，然后再轉(zhuǎn)化成匿名函數(shù)，則比較方便。本例展示了符號計算和數(shù)值計算結(jié)算結(jié)合的一種途徑。代碼如下：symsxf=(x+tan(x))^(sin(x));c=diff(f,3);f3=eval(['@(x)'vectorize(c)]);x=linspace(0,1,100);plot(x,f3(x))title('(x+tan(x))\^(sin(x))三階導數(shù)圖像')2023/1/312023/1/31一、匿名函數(shù)應(yīng)用實例5. 匿名函數(shù)在優(yōu)化中的應(yīng)用 匿名函數(shù)在優(yōu)化中的應(yīng)用主要是以表示目標函數(shù)的形式出現(xiàn)的。

【例4.5-6】求下面函數(shù)的最小值 用匿名函數(shù)來表示目標函數(shù)如下：

f=@(x)3*x(1)^2+2*x(1)*x(2)+x(2)^2;

進一步求解:

x0=[1,1];%初始值

[x,fval]=fminunc(f,x0) x= 1.0e-006* 0.2541-0.2029 fval= 1.3173e-0132023/1/31二、子函數(shù)和嵌套函數(shù)應(yīng)用實例在求解積分上限中的應(yīng)用【例4.5-7】如下述積分表達式，如何已知和，求得？2023/1/31functionsol=exampleIntLimit1(a,e,l)%用嵌套函數(shù)表示被積表達式

functionf=fun1(beta)f=a.*(1-e.^2)./(1-e.^2.*sin(beta).^2).^(3/2);end%調(diào)用fzero求滿足條件的beta0值sol=fzero(@(beta0)fun2(beta0,l,@fun1),3);end%用子函數(shù)表示積分functiong=fun2(beta0,l,fhdle)g=quadl(fhdle,0,beta0)-l;end2023/1/31二、子函數(shù)和嵌套函數(shù)應(yīng)用實例2. 在GUI中的應(yīng)用（三角界面，見書中代碼）在3D作圖中的應(yīng)用【例4.5-9】畫出下列函數(shù)的圖像：其中2023/1/31function[m,n,TT]=plot3dnmT(N,L)%N:inf的近似，L:[0,2]區(qū)間的剖分個數(shù)C=zeros(N,1);%nested-function：Tmn=calcT(mm,nn)中用來存儲計算結(jié)果m=linspace(0,2,L);[m,n]=meshgrid(m,m);TT=zeros(size(n));%和網(wǎng)格數(shù)據(jù)m,n對應(yīng)的計算出來的T(m,n)網(wǎng)格數(shù)據(jù)forii=1:Lforjj=1:LTT(ii,jj)=calcT(m(ii,jj),n(ii,jj));endend%=====計算T(m,n)的nest-functionfunctionTmn=calcT(mm,nn)forN1=1:NC(N1)=(mm^N1/gamma(N1+1))*sum(nn.^(0:N1-1)./gamma(1:N1));Tmn=1.0-exp(-mm-nn)*sum(C);endendmesh(n,m,TT);end

2023/1/31嵌套函數(shù)表示待優(yōu)化的目標函數(shù)【例4.5-10】已知求下面表達式的最小值。其中，m在[0,2]這個區(qū)間范圍內(nèi)

2023/1/31functionm=Findmw=[pi/2,pi,pi*1.5];N=[pi/2-1,-2,-1.5*pi-1];functiony=ObjecFun(m)y=(quadl(@(t)t.^m.*cos(t),0,w(1))-N(1))^2+(quadl(@(t)t.^m.*cos(t),0,...w(2))-N(2))^2+(quadl(@(t)t.^m.*cos(t),0,w(3))-N(3))^2;endm=fminbnd(@ObjecFun,0,2);end2023/1/31二、子函數(shù)和嵌套函數(shù)應(yīng)用實例5. 嵌套函數(shù)在表示微分方程方面的應(yīng)用（求解微分方程，見書中代碼）2023/1/31第六節(jié)編寫高效的MATLAB程序2023/1/31一、重新認識循環(huán)1. 高版本MATLAB對循環(huán)結(jié)構(gòu)的優(yōu)化

從MATLAB6.5版開始，MATLAB引入了JIT（justintime）技術(shù)和加速器（accelerator），并在后續(xù)版本中不斷優(yōu)化。到了R2009a，很多情況下，循環(huán)體本身已經(jīng)不是程序性能提高的瓶頸了，瓶頸更多的來源于循環(huán)體內(nèi)部的代碼實現(xiàn)方式，以及使用循環(huán)的方式。(詳見例【例4.6-1】)選擇循環(huán)還是向量化凡是涉及到矩陣運算的時候，則盡量用向量化如果向量化可能導致超大型矩陣的產(chǎn)生，使用前要慎重向量化的使用要靈活，多分析其運行機制，并與循環(huán)做對比MATLAB初學者要盡量多用向量化思路編程高效的編程包括高效的開發(fā)和高效的運行

2023/1/31二、提高代碼效率的方法

1. 預分配內(nèi)存 當矩陣尺寸較大時，預分配內(nèi)存（preallocation），是MATLAB高效編程中最應(yīng)該注意，同時也是很容易把握的一個技巧(詳見例【例4.6-5】)選用恰當?shù)暮瘮?shù)類型

MATLAB有多種函數(shù)類型，不同的函數(shù)類型調(diào)用效率差別相當大，為了寫出高效的MATLAB代碼，我們就要選用恰當?shù)暮瘮?shù)類型來實現(xiàn)我們的任務(wù)。(詳見例【例4.6-6】【例4.6-7】)選用恰當?shù)臄?shù)據(jù)類型

MATLAB同樣有多種數(shù)據(jù)類型，不同的數(shù)據(jù)類型在存儲和訪問效率上也不盡相同，選用恰當?shù)臄?shù)據(jù)類型對寫出高效的MATLAB代碼同樣具有重要的意義。(詳見例【例4.6-8】)2023/1/31二、提高代碼效率的方法

4. 減少減少無謂損耗——給一些函數(shù)“瘦身”

MATLAB自身以及工具箱里的很多函數(shù)具有較