第2講MATLAB程序設(shè)計(jì)及其調(diào)試

M文件的功能和特點(diǎn)MATLAB是解釋性語言語法簡單調(diào)試容易人機(jī)交互性強(qiáng)開放性、可擴(kuò)展性強(qiáng)第2講MATLAB程序設(shè)計(jì)及其調(diào)試M文件輸入窗口?。?！例1：通過M腳本文件，畫出下列分段函數(shù)所表示的曲面.腳本文件first.m%first.m Thisismyfirstexample.a=2;b=2; clf;x=-a:0.2:a;y=-b:0.2:b;fori=1:length(y)forj=1:length(x)ifx(j)+y(i)>1z(i,j)=0.5457*exp(-0.75*y(i)^2-3.75*x(j)^2-1.5*x(j));elseifx(j)+y(i)<=-1z(i,j)=0.5457*exp(-0.75*y(i)^2-3.75*x(j)^2+1.5*x(j));elsez(i,j)=0.7575*exp(-y(i)^2-6.*x(j)^2);endendendaxis([-a,a,-b,b,min(min(z)),max(max(z))]);colormap(flipud(winter));surf(x,y,z);程序調(diào)試1程序調(diào)試概述一般來說，應(yīng)用程序的錯(cuò)誤有兩類，一類是語法錯(cuò)誤，另一類是運(yùn)行時(shí)的錯(cuò)誤。語法錯(cuò)誤包括詞法或文法的錯(cuò)誤，例如函數(shù)名的拼寫錯(cuò)、表達(dá)式書寫錯(cuò)等。程序運(yùn)行時(shí)的錯(cuò)誤是指程序的運(yùn)行結(jié)果有錯(cuò)誤，這類錯(cuò)誤也稱為程序邏輯錯(cuò)誤。1．控制單步運(yùn)行2．?dāng)帱c(diǎn)操作

3調(diào)試命令除了采用調(diào)試器調(diào)試程序外，MATLAB還提供了一些命令用于程序調(diào)試。命令的功能和調(diào)試器菜單命令類似，具體使用方法請讀者查詢MATLAB幫助文檔。例2：某實(shí)驗(yàn)對一根長10米的鋼軌進(jìn)行熱源的溫度傳播測試。用x表示測量點(diǎn)0:2.5:10(米)，用h表示測量時(shí)間0:30:60(秒)，用T表示測試所得各點(diǎn)的溫度(℃)。試用線性插值求出在一分鐘內(nèi)每隔20秒、鋼軌每隔1米處的溫度TI。T=[95,14,0,0,0;88,48,32,12,6;67,64,54,48,41];x=0:2.5:10;h=[0:30:60]';T=[95,14,0,0,0;88,48,32,12,6;67,64,54,48,41];xi=[0:10];hi=[0:20:60]';TI=interp2(x,h,T,xi,hi)例3數(shù)據(jù)的輸入輸出處理方法。Matlab環(huán)境下的文件與其它系統(tǒng)一樣，也有二類文件組成，一是文件，又稱M文件，另一類是數(shù)據(jù)文件.系統(tǒng)除提供了文件的一般管理功能外，還提供了對數(shù)據(jù)文件進(jìn)行操作的特殊功能函數(shù).3.1文件的打開與關(guān)閉Matlab提供了對數(shù)據(jù)文件建立、打開、讀、寫以及關(guān)閉等一系列函數(shù)，數(shù)據(jù)文件一般存放在磁盤等介質(zhì)上，用文件名標(biāo)識，系統(tǒng)對文件名沒有特殊要求。文件數(shù)據(jù)格式有二種形式，一是二進(jìn)制格式文件，二是文本文件，系統(tǒng)對這兩類文件提供了不同的讀寫功能函數(shù).文件的打開

fopen函數(shù)的調(diào)用格式為：

fid=fopen(文件名，打開方式)

其中文件名用字符串形式，表示待打開的數(shù)據(jù)文件.

常見的打開方式有：

‘r’

表示對打開的文件讀數(shù)據(jù)，

‘w’

表示對打開的文件寫數(shù)據(jù)，

‘a(chǎn)’

表示在打開的文件末尾添加數(shù)據(jù).

fid用于存儲文件句柄值，句柄值用來標(biāo)識該數(shù)據(jù)文件，其他函數(shù)可以利用它對該數(shù)據(jù)文件進(jìn)行操作.例如：打開一個(gè)名為std.dat的數(shù)據(jù)文件并進(jìn)行讀操作，其命令格式為：Fid=fopen(‘std.dat’,’r’)上述打開格式均為二進(jìn)制格式，如果想用ASCII文本格式，則必須在格式字符串中加上字符t，例如用’rt’表示以ASCII格式打開供讀操作的數(shù)據(jù)文件.

其他打開方式：‘r’

打開文件，讀數(shù)據(jù)，文件必須存在.‘w’

打開文件，寫數(shù)據(jù)，若文件不存在，系統(tǒng)會自動建立.‘a(chǎn)’

打開文件，在文件末尾添加數(shù)據(jù).‘r+’

打開文件，可以讀和寫數(shù)據(jù)，文件必須存在.‘w+’

打開文件，供讀與寫數(shù)據(jù)用.‘a(chǎn)+’

打開文件，供讀與添加數(shù)據(jù)用.2.文件的關(guān)閉fclose關(guān)閉文件文件在進(jìn)行完讀、寫等操作后，應(yīng)及時(shí)關(guān)閉，以保證文件的安全可靠.關(guān)閉文件命令格式為：Sta=fclose(Fid)

關(guān)閉Fid所表示的文件Sta表示關(guān)閉文件操作的返回代碼，若關(guān)閉成功，返回0，否則返回–1.二進(jìn)制文件的讀寫操作

(1)讀二進(jìn)制文件

fread函數(shù)可以讀取二進(jìn)制文件的數(shù)據(jù)，并將數(shù)據(jù)存入矩陣。其調(diào)用格式為：

[A,Count]=fread(fid,size,precision)

其中A用于存放讀取的數(shù)據(jù)，Count返回所讀取的數(shù)據(jù)元素個(gè)數(shù)，fid為文件句柄，size為可選項(xiàng)，若不選用則讀取整個(gè)文件內(nèi)容，若選用則它的值可以是下列值：3.2文件的讀寫操作

(1)N表示讀取N個(gè)元素到一個(gè)列向量.

(2)Inf表示讀取整個(gè)文件.

(3)[M,N]表示讀數(shù)據(jù)到M×N的矩陣中，數(shù)據(jù)按列存放.

precision代表讀寫數(shù)據(jù)的類型

如:'int8’:Integer;8bits'int16’:Integer;16bits(2)寫二進(jìn)制文件

fwrite

函數(shù)按照指定的數(shù)據(jù)類型將矩陣中的元素寫入到文件中.其調(diào)用格式為：

Count=fwrite(fid,A,precision)

其中Count返回所寫的數(shù)據(jù)元素個(gè)數(shù)，fid為文件句柄，A用來存放寫入文件的數(shù)據(jù)，precision用于控制所寫數(shù)據(jù)的類型，其形式與fread函數(shù)相同.

例.建立一數(shù)據(jù)文件magic5.dat，用于存放5階魔方陣.

程序如下：

fid=fopen(‘magic5.dat’,‘w’);

cnt=fwrite(fid,magic(5),‘int32’);

%int32:Integer;32bitsfclose(fid);

例.建立一數(shù)據(jù)文件test.dat，用于存放矩陣A的數(shù)據(jù).已知A=[-0.6515-0.2727-0.4354-0.3190-0.9047-0.7534-0.4567-0.3212-0.4132-0.3583-0.9264-0.8173-0.7823-0.3265-0.0631-0.1735-0.7373-0.0972-0.3267-0.6298-0.4768-0.6773-0.6574-0.1923-0.4389]Fid=fopen(‘test.dt’,‘w’)cnt=fwrite(Fid,A,‘float’)fclose(Fid)

程序段將矩陣A的數(shù)據(jù)以二進(jìn)制浮點(diǎn)數(shù)格式寫入文件test.dat中.Fid=fopen(‘test.dt’,‘r’)[B,cnt]=fread(Fid,[5,5],‘float’)fclose(Fid)

讀取文件test.dat的內(nèi)容.[5,5]改成[5,4]，有何不同？

2文本文件fscanf

讀ASCII文本文件

[A,COUNT]=fscanf(Fid,format,size)其中A為數(shù)據(jù)矩陣，用以存放讀取的數(shù)據(jù)，COUNT返回所讀取的數(shù)據(jù)元素個(gè)數(shù)。format用以控制讀取的數(shù)據(jù)格式，由%加上格式符組成，格式符為：

d,i,o,u,x,e,f,g,s,c與[...]

例如：s=fscanf(fid,‘%s’)讀取一個(gè)字符串a(chǎn)=fscanf(fid,‘%5d’)讀取5位數(shù)的整數(shù)b=fscanf(fid,‘%6.2d’)讀取浮點(diǎn)數(shù)fprintf

寫ASCII數(shù)據(jù)文件,其格式為：

COUNT=fprintf(Fid,format,A,…)其中A為要寫入文件的數(shù)據(jù)矩陣，先按format格式化數(shù)據(jù)矩陣A，后寫入到Fid所指定的文件.例如：x=0:0.1:1;y=[x;exp(x)];Fid=fopen('exp.dt','w');fprintf(Fid,'%6.2f%12.8f\n',y);fclose(Fid);1.求極限函數(shù)limit用于求符號函數(shù)f的極限.系統(tǒng)可以根據(jù)用戶要求，計(jì)算變量從不同方向趨近于指定值的極限值.該函數(shù)的格式及功能：

limit(f,x,a)：求符號函數(shù)f（x）的極限值.即計(jì)算當(dāng)變量x趨近于常數(shù)a時(shí)，f（x）函數(shù)的極限值.3.3微積分函數(shù)

limit(f,a):求符號函數(shù)f（x）的極限值.由于沒有指定符號函數(shù)f（x）的自變量，既變量x趨近于a.

limit(f):求符號函數(shù)f（x）的極限值.沒有指定變量的目標(biāo)值時(shí)，系統(tǒng)默認(rèn)變量趨近于0，即a=0的情況.

limit(f,x,a,'right'):求符號函數(shù)f的極限值.'right'表示變量x從右邊趨近于a.

limit(f,x,a,'left'):求符號函數(shù)f的極限值.'left'表示變量x從左邊趨近于a.2.符號變量的定義參與符號運(yùn)算的對象可以是符號變量、符號表達(dá)式或符號矩陣.

符號變量要先定義，后引用.可以用sym函數(shù)、syms函數(shù)將運(yùn)算量定義為符號型數(shù)據(jù).

引用符號運(yùn)算函數(shù)時(shí)，用戶可以指定函數(shù)執(zhí)行過程中的變量參數(shù)；若用戶沒有指定變量參數(shù)，則使用findsym函數(shù)默認(rèn)的變量作為函數(shù)的變量參數(shù).2.1sym函數(shù)

sym函數(shù)的主要功能是創(chuàng)建符號變量，以便進(jìn)行符號運(yùn)算，也可以用于創(chuàng)建符號表達(dá)式或符號矩陣.用sym函數(shù)創(chuàng)建符號變量的一般格式為：

x=sym(‘x’)

其目的是將’x’創(chuàng)建為符號變量，以x作為輸出變量名.每次調(diào)用該函數(shù),可以定義一個(gè)符號變量.例作符號計(jì)算：

a,b,x,y均為符號運(yùn)算量.在符號運(yùn)算前，應(yīng)先將a,b,x,y定義為符號運(yùn)算量.a=sym(‘a(chǎn)’);

%定義‘a(chǎn)’為符號運(yùn)算量，輸出變量名為ab=sym(‘b’);x=sym(‘x’);y=sym(‘y’);[x,y]=solve(a*x-b*y-1,a*x+b*y-5,x,y)%以a,b為符號常數(shù)，x,y為符號變量即可得到方程組的解：

x=3/ay=2/b2.2syms函數(shù)

syms函數(shù)的功能與sym函數(shù)類似.syms函數(shù)可以在一個(gè)語句中同時(shí)定義多個(gè)符號變量，其一般格式為：

symsarg1arg2…argN

例求極限symsx；

%定義符號變量f=(x*(exp(sin(x))+1)-2*(exp(tan(x))-1))/sin(x)^3；

%確定符號表達(dá)式w=limit(f)

%求函數(shù)的極限w=-1/23.微分函數(shù)diff函數(shù)用于對符號表達(dá)式s求微分.該函數(shù)的一般引用格式為：

diff(s,’v’,n)

應(yīng)用diff（s）沒有指定微分變量和微分階數(shù)，則系統(tǒng)按findsym函數(shù)指示的默認(rèn)變量對符號表達(dá)式s求一階微分.

應(yīng)用diff（s，‘v’）或diff（s，sym（‘v’））格式，表示以v為自變量，對符號表達(dá)式s求一階微分.

應(yīng)用diff（s，n）格式，表示對符號表達(dá)式s求n階微分，n為正整數(shù).

應(yīng)用diff（s，‘v’，n）格式，表示以v為自變量，對符號表達(dá)式s求n階微分.例求導(dǎo)數(shù)：x=sym('x');%定義符號變量diff(sin(x^2))%求導(dǎo)運(yùn)算ans=2*cos(x^2)*x4．積分函數(shù)積分函數(shù)int（s，v，a，b)可以對被積函數(shù)或符號表達(dá)式s求積分。其引用格式為：

int（s，v，a，b)說明：應(yīng)用int（s）格式，表示沒有指定積分變量和積分階數(shù)時(shí)，系統(tǒng)按findsym函數(shù)指示的默認(rèn)變量對被積函數(shù)或符號表達(dá)式s求一階積分.應(yīng)用int（s，v）格式，表示以v為自變量，對被積函數(shù)或符號表達(dá)式s求一階不定積分.應(yīng)用積分函數(shù)時(shí)，如果給定a、b兩項(xiàng)，表示是進(jìn)行定積分運(yùn)算.a、b分別表示定積分的下限和上限。不指定積分的下限和上限表示求不定積分.例求下述積分.求積分：symsxint(1/(1+x^2))ans=atan(x)解方程函數(shù)的格式為：

solve(expr1,expr2,...,exprN,var1,var2,...varN)或

solve(expr1,expr2,...,exprN)其功能為：求解代數(shù)方程組expr1,expr2,...,exprN的根，未知數(shù)為var1,var2,...varN.說明：若不指明符號表達(dá)式expr1,expr2,...,exprN的值，系統(tǒng)默認(rèn)為0.例如給出一個(gè)表達(dá)式x^2-3*x-8,則系統(tǒng)將按x^2-3*x-8=0進(jìn)行運(yùn)算.5.解方程例解代數(shù)方程：a*x2-b*x-6=0symsabxsolve(a*x^2-b*x-6)ans=[1/2/a*(b+(b^2+24*a)^(1/2))][1/2/a*(b-(b^2+24*a)^(1/2))]即該方程有兩個(gè)根:x1=1/2/a*(b+(b^2+24*a)^(1/2))；

x2=1/2/a*(b-(b^2+24*a)^(1/2))6.級數(shù)(級數(shù)求和)級數(shù)求和運(yùn)算是數(shù)學(xué)中常見的一種運(yùn)算。例如：

f(x)=a0+a1x+a2x2+a3x3+…+anxn函數(shù)symsum可以用于此類對符號函數(shù)f的求和運(yùn)算.該函數(shù)的引用時(shí)，應(yīng)確定級數(shù)的通項(xiàng)式s，變量的變化范圍a和b.該函數(shù)的引用格式為：

symsum(s,a,b)例求級數(shù)的和:

鍵入：1/12+1/22+1/32+1/42+……symsksymsum(1/k^2,1,Inf)%k值為1到無窮大ans=1/6*pi^2其結(jié)果為：1/12+1/22+1/32+1/42+……=π2/6MATLAB求解優(yōu)化問題的主要函數(shù)優(yōu)化函數(shù)的輸入變量

優(yōu)化函數(shù)的輸出變量

非線性優(yōu)化有約束的一元函數(shù)的最小值單變量函數(shù)求最小值的標(biāo)準(zhǔn)形式為

sub.to函數(shù)fminbnd格式x=fminbnd(fun,x1,x2)%返回自變量x在區(qū)間上函數(shù)fun取最小值時(shí)x值，fun為目標(biāo)函數(shù)的表達(dá)式字符串或MATLAB自定義函數(shù)的函數(shù)柄.x=fminbnd(fun,x1,x2,options)[x,fval]=fminbnd(…)[x,fval,exitflag]=fminbnd(…)[x,fval,exitflag,output]=fminbnd(…)例4計(jì)算下面函數(shù)在區(qū)間(0，1)內(nèi)的最小值.

解：>>[x,fval,exitflag,output]=fminbnd('(x^3+cos(x)+x*log(x))/exp(x)',0,1)x=0.5223fval=0.3974exitflag=1output=iterations:9funcCount:9algorithm:'goldensectionsearch,parabolicinterpolation'例5在[0，5]上求下面函數(shù)的最小值解：先自定義函數(shù)：在MATLAB編輯器中建立M文件為：functionf=myfun(x)f=(x-3).^2-1;保存為myfun.m，然后在命令窗口鍵入命令：>>x=fminbnd(@myfun,0,5)則結(jié)果顯示為：x=3無約束多元函數(shù)最小值多元函數(shù)最小值的標(biāo)準(zhǔn)形式為其中：x為向量.命令利用函數(shù)fminsearch求無約束多元函數(shù)最小值.函數(shù)fminsearch格式x=fminsearch(fun,x0)%x0為初始點(diǎn)，fun為目標(biāo)函數(shù)的表達(dá)式字符串或MATLAB自定義函數(shù)的函數(shù)柄.x=fminsearch(fun,x0,options)%options查optimset.[x,fval]=fminsearch(…)%最優(yōu)點(diǎn)的函數(shù)值.[x,fval,exitflag]=fminsearch(…)%exitflag與單變量情形一致.[x,fval,exitflag,output]=fminsearch(…)%output與單變量情形一致.例6求的最小值點(diǎn).解：>>X=fminsearch('2*x(1)^3+4*x(1)*x(2)^3-10*x(1)*x(2)+x(2)^2',[0,0])結(jié)果為

X=1.00160.8335或在MATLAB編輯器中建立函數(shù)文件.functionf=myfun(x)f=2*x(1)^3+4*x(1)*x(2)^3-10*x(1)*x(2)+x(2)^2;保存為myfun.m，在命令窗口鍵入>>X=fminsearch('myfun',[0,0])或>>X=fminsearch(@myfun,[0,0])結(jié)果為：

X=1.00160.8335有約束的多元函數(shù)最小值非線性有約束的多元函數(shù)的標(biāo)準(zhǔn)形式為：sub.to其中：x、b、beq、lb、ub是向量，A、Aeq為矩陣，C(x)、Ceq(x)是返回向量的函數(shù)，f(x)為目標(biāo)函數(shù)，f(x)、C(x)、Ceq(x)可以是非線性函數(shù).在MATLAB5.x中，它的求解由函數(shù)constr實(shí)現(xiàn).函數(shù)

fmincon格式x=fmincon(fun,x0,A,b)x=fmincon(fun,x0,A,b,Aeq,beq)x=fmincon(fun,x0,A,b,Aeq,beq,lb,ub)x=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon)x=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)[x,fval]=fmincon(…)[x,fval,exitflag]=fmincon(…)[x,fval,exitflag,output]=fmincon(…)[x,fval,exitflag,output,lambda]=fmincon(…)[x,fval,exitflag,output,lambda,grad]=fmincon(…)[x,fval,exitflag,output,lambda,grad,hessian]=fmincon(…)參數(shù)說明：fun為目標(biāo)函數(shù)，它可用前面的方法定義；nonlcon的作用是通過接受的向量x來計(jì)算非線性不等約束和等式約束分別在x處的估計(jì)C和Ceq，通過指定函數(shù)柄來使用，如：>>x=fmincon(@myfun,x0,A,b,Aeq,beq,lb,ub,@mycon)，先建立非線性約束函數(shù)，并保存為mycon.m：function[C,Ceq]=mycon(x)C=…%計(jì)算x處的非線性不等約束的函數(shù)值.Ceq=…%計(jì)算x處的非線性等式約束的函數(shù)值.lambda是Lagrange乘子，它體現(xiàn)哪一個(gè)約束有效.output輸出優(yōu)化信息；grad表示目標(biāo)函數(shù)在x處的梯度；hessian表示目標(biāo)函數(shù)在x處的Hessian值.例7求下面問題在初始點(diǎn)（0，1）處的最優(yōu)解minsub.to解：約束條件的標(biāo)準(zhǔn)形式為：sub.to先在MATLAB編輯器中建立非線性約束函數(shù)文件：

function[c,ceq]=mycon(x)c=(x(1)-1)^2-x(2);ceq=[];%無等式約束.然后，在命令窗口鍵入如下命令或建立M文件：>>fun='x(1)^2+x(2)^2-x(1)*x(2)-2*x(1)-5*x(2)';%目標(biāo)函數(shù).>>x0=[01];>>A=[-23];%線性不等式約束.>>b=6;>>Aeq=[];%無線性等式約束.>>beq=[];>>lb=[];%x沒有下、上界.>>ub=[];>>[x,fval,exitflag,output,lambda,grad,hessian]=fmincon(fun,x0,A,b,Aeq,beq,lb,ub,@mycon)則結(jié)果為x=34fval=-13exitflag=1%解收斂.

output=iterations:2funcCount:9stepsize:1algorithm:'medium-scale:SQP,Quasi-Newton,line-search'firstorderopt:[]cgiterations:[]lambda=lower:[2x1double]%x下界有效情況，通過lambda.lower可查看.upper:[2x1double]%x上界有效情況，為0表示約束無效.

eqlin:[0x1double]%線性等式約束有效情況，不為0表示約束有效.eqnonlin:[0x1double]%非線性等式約束有效情況.ineqlin:2.5081e-008%線性不等式約束有效情況.neqnonlin:6.1938e-008%非線性不等式約束有效情況.grad=%目標(biāo)函數(shù)在最小值點(diǎn)的梯度.1.0e-006*-0.1776hessian=%目標(biāo)函數(shù)在最小值點(diǎn)的Hessian值.1.0000-0.0000-0.00001.0000二次規(guī)劃問題二次規(guī)劃問題（quadraticprogramming）的標(biāo)準(zhǔn)形式為：

sub.to

其中，H、A、Aeq為矩陣，f、b、beq、lb、ub、x為向量其它形式的二次規(guī)劃問題都可轉(zhuǎn)化為標(biāo)準(zhǔn)形式.MATLAB5.x版中的qp函數(shù)已被6.0版中的函數(shù)quadprog取代。函數(shù)

quadprog格式x=quadprog(H,f,A,b)%其中H,f,A,b為標(biāo)準(zhǔn)形中的參數(shù)，x為目標(biāo)函數(shù)的最小值.x=quadprog(H,f,A,b,Aeq,beq)%Aeq,beq滿足等約束條件.x=quadprog(H,f,A,b,Aeq,beq,lb,ub)%lb,ub分別為解x的下界與上界.x=quadprog(H,f,A,b,Aeq,beq,lb,ub,x0)%x0為設(shè)置的初值x=quadprog(H,f,A,b,Aeq,beq,lb,ub,x0,options)%options為指定的優(yōu)化參數(shù).[x,fval]=quadprog(…)%fval為目標(biāo)函數(shù)最優(yōu)值.[x,fval,exitflag]=quadprog(…)%exitflag與線性規(guī)劃中參數(shù)意義相同.[x,fval,exitflag,output]=quadprog(…)%output與線性規(guī)劃中參數(shù)意義相同.[x,fval,exitflag,output,lambda]=quadprog(…)%lambda與線性規(guī)劃中參數(shù)意義相同.例8求二次規(guī)劃的最優(yōu)解

maxf(x1,x2)=x1x2+3sub.tox1+x2-2=0解：化成標(biāo)準(zhǔn)形式：

sub.tox1+x2=2在Matlab中實(shí)現(xiàn)如下：H=[0,-1;-1,0];f=[0;0];Aeq=[11];b=2;[x,fval,exitflag,output,lambda]=quadprog