實(shí)用-matlab教學(xué)資料-25-26常微分方程

實(shí)用MATLAB1精選ppt例1生產(chǎn)決策問(wèn)題——如何收入最高？某廠甲乙兩種產(chǎn)品,每種產(chǎn)品所需原料量如表。若1kg產(chǎn)品甲和乙的售價(jià)分別為6萬(wàn)元和5萬(wàn)元,原料ABC的限量分別為100kg,160kg,180kg。試確定生產(chǎn)這兩種產(chǎn)品各多少kg才能使總銷(xiāo)售收入最高？原料產(chǎn)品甲產(chǎn)品乙A,kg53B,kg67C,kg052精選pptlinprog[x,fval,exitflag,output,lamda]=

linprog(f,A,B,Aeq,Beq,lb,ub,x0,options)f:目標(biāo)函數(shù)系數(shù),列向量標(biāo)準(zhǔn)模型3精選ppt例2,3解極小值問(wèn)題funx=@(x)exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);funx=@(x)56*(x(2)-x(1)^2)^2+(1-x(1))^2;4精選ppt

例6求最大營(yíng)業(yè)額及方營(yíng)銷(xiāo)案

設(shè)兩種設(shè)備的銷(xiāo)售量分別為x1,x2,售價(jià)及其售出所需營(yíng)業(yè)時(shí)間如下表。求在總營(yíng)業(yè)時(shí)間800h內(nèi)最大營(yíng)業(yè)額以及銷(xiāo)售計(jì)劃。設(shè)備售價(jià),元售出所需營(yíng)業(yè)時(shí)間,h1300.524502+0.25x25精選pptmaxf(x)min-f(x)minf(x)s.t.約束條件programming極大值問(wèn)題極小值問(wèn)題6精選ppt9.6MATLAB優(yōu)化工具箱7精選ppt一、優(yōu)化工具的啟動(dòng)在命令窗口輸入optimtoolMATLAB主界面APPS→OptimizationTool

8精選ppt9精選ppt問(wèn)題描述及結(jié)果顯示ProblemSetupandResults優(yōu)化參數(shù)設(shè)置Options幫助QuickReference右上角<<，隱藏二、界面簡(jiǎn)介10精選ppt1、問(wèn)題描述SolverAlgorithmProblemConstraints11精選ppt結(jié)果顯示RunsolverandviewresultsFinalpoint將結(jié)果輸出到workspace12精選ppt2．優(yōu)化參數(shù)設(shè)置Stoppingcriteria:停止準(zhǔn)則Functionvaluecheck:函數(shù)值檢查User-suppliedderivatives:用戶自定義導(dǎo)數(shù)Approximatedderivatives:數(shù)值微分Hessian矩陣Algorithmsettings:算法設(shè)置Inneriterationstoppingcriteria:內(nèi)部迭代停止準(zhǔn)則Plotfunctions:用戶自定義繪圖函數(shù)Outputfunctions:用戶自定義輸出函數(shù)Displaytocommandwindow:輸出到命令行窗口類(lèi)似options=optimoptions(……)功能13精選ppt三、使用步驟1、選擇求解器和優(yōu)化算法；

2、定義目標(biāo)函數(shù)和相關(guān)參數(shù)；

3、設(shè)置優(yōu)化選項(xiàng)；

4、單擊“Start”求解；

5、查看求解狀態(tài)和求解結(jié)果；

6、導(dǎo)出目標(biāo)函數(shù)、選項(xiàng)和結(jié)果。提醒：可隨時(shí)查看第三列QuickReference，了解到參數(shù)意義。在該界面中選中文字，右鍵點(diǎn)擊find，可進(jìn)一步查詢。14精選ppt例Optimizationtool15精選ppt例

求表面積為300m2的最大圓柱體體積。Optimizationtool16精選ppt例

某化學(xué)反應(yīng)實(shí)驗(yàn)所得生成物的濃度隨時(shí)間的變化數(shù)據(jù)，（1）擬合模型參數(shù)（2）繪出擬合曲線和數(shù)據(jù)點(diǎn)t12345678y×1034.006.408.008.809.229.509.709.86t910111213141516y×10310.0010.2010.3210.4210.5010.5510.5810.60Optimizationtool17精選ppt第10章常微分方程數(shù)值解OrdinaryDifferentialEquations18精選pptODE數(shù)值求解思路首先將常微分方程(組)及其邊界條件離散化，即轉(zhuǎn)化為差分方程；然后求得常微分方程(組)在離散點(diǎn)上的函數(shù)近似值，這些近似值即為ODE數(shù)值解。離散化19精選ppt幾何意義——折線近似原函數(shù)曲線x1x2x0xnP0P1PnP2xyy=y(x)20精選ppt1、初值問(wèn)題的描述在自變量的一端給定邊界條件I.C.一階常微分方程初值問(wèn)題10.1初值問(wèn)題一、初值問(wèn)題數(shù)值解法21精選ppt2、數(shù)值求解方法Euler法Runge-Kutta法線性多步法22精選pptCarlRunge（卡爾·龍格）(1856－1927）德國(guó)數(shù)學(xué)家，物理學(xué)家，光譜學(xué)家柏林大學(xué)數(shù)學(xué)博士師從德國(guó)著名數(shù)學(xué)家、被譽(yù)為“現(xiàn)代分析之父”的卡爾·魏爾施特拉斯(ErnstKummer)Runge-Kutta法23精選ppt主要貢獻(xiàn)1、龍格現(xiàn)象2、龍格－庫(kù)塔法3、拉普拉斯-龍格-楞次矢量（LRL）CarlRunge（卡爾·龍格）(1856－1927）24精選ppt二、MATLAB功能函數(shù)求解器Solver算法ODE類(lèi)型精度適用場(chǎng)合ode45單步顯式4,5階R-K法非剛性中等大部分場(chǎng)合的首選算法ode23單步顯式2，3階R-K法非剛性低適用于精度要求較低或中等剛度的場(chǎng)合ode113多步法：Adams法非剛性低~高適用于精度要求較嚴(yán)格的場(chǎng)合或計(jì)算量較大的問(wèn)題。odesolver功能函數(shù)25精選ppt1、一階常微分方程O(píng)DE26精選ppt[x,y]=ode45(odefun,tspan,y0,options)1、一階常微分方程ode4527精選ppt[x,y]=ode45(odefun,tspan,y0,options)(1)odefun：定義微分方程

函數(shù)式m文件匿名函數(shù)28精選ppt(2)tspan:給出求解區(qū)間或節(jié)點(diǎn)，

可以指定區(qū)間[a,b]可以指定節(jié)點(diǎn)[x0,x1,…]等距節(jié)點(diǎn)[a:h:b]注意：h是輸出步長(zhǎng)，不是計(jì)算步長(zhǎng)[x,y]=ode45(odefun,tspan,y0,options)(3)y0:初值29精選ppt(4)optionsoptions=odeset('name1',value1,...)[x,y]=ode45(odefun,tspan,y0,options)誤差控制輸出控制步長(zhǎng)控制……30精選ppt(4)options誤差控制namevaluedescriptionRelToldefault,1e-3相對(duì)誤差A(yù)bsToldefault,1e-6絕對(duì)誤差options=odeset('name1',value1,...)[x,y]=ode45(odefun,tspan,y0,options)例options=odeset('RelTol',1e-5)31精選ppt[x,y]=ode45(odefun,tspan,y0,options)或Sol=ode45(odefun,tspan,y0,options)(5)x:自變量序列(6)y:因變量序列(7)Sol:解的結(jié)構(gòu)數(shù)組，包括Sol.x,Sol.y,Sol.solver32精選ppt例1解常微分方程在[0,1]上的解，輸出步長(zhǎng)h=0.1

[x,y]=ode45(odefun,tspan,y0,options)33精選ppt2、一階常微分方程組ODE34精選ppt[x,y]=ode45(odefun,tspan,y0,options)(1)odefun：定義微分方程組

函數(shù)式m文件或匿名函數(shù)odefun返回的必須是列向量2、一階常微分方程組35精選ppt例函數(shù)式m文件定義微分方程functiondydx=odefun(x,y)dydx(1)=y(2);dydx(2)=1+2*x+4*x^2*y(1);dydx=[dydx(1);dydx(2)];[x,y]=ode45(@odefun,tspan,y0,options)36精選ppt例匿名函數(shù)定義微分方程odefun=@(x,y)[y(2);1+x^2*y(1)];[x,y]=ode45(odefun,tspan,y0,options)37精選ppt(2)y0:初值，列向量[y1(a);y2(a);y3(a);....][x,y]=ode45(odefun,tspan,y0,options)38精選ppt例2求常微分方程組在[0,1]上的解[x,y]=ode45(odefun,tspan,y0,options)39精選ppt例2求常微分方程組在[0,1]上的解，輸出步長(zhǎng)0.1[x,y]=ode113(odefun,tspan,y0,options)40精選pptsave

使用文件保存數(shù)據(jù)

格式：save(filename,va1,…,format)功能：采用format規(guī)定的格式存儲(chǔ)變量va1,…

的值到名為filename的文件中41精選pptsave(filename,va1,…,format)（1）filename:文件名在程序運(yùn)行時(shí)，自動(dòng)在當(dāng)前文件夾創(chuàng)建，并寫(xiě)入數(shù)據(jù)。僅當(dāng)該文件關(guān)閉時(shí)，數(shù)據(jù)才能寫(xiě)入其中。42精選pptsave(filename,va1,…,format)filename:文件名文件類(lèi)型后綴word*.docexcel*.xlstxt*.txtmat*.mat(存儲(chǔ)為matlab的變量格式)43精選pptsave(filename,va1,…,format)format:存儲(chǔ)格式如果是excel文件，最好使用-tabs格式'-ascii'8位數(shù)'-ascii','-tabs'8位數(shù)，定位輸出'-ascii','-double'16位數(shù)'-ascii','-double','-tabs'16位數(shù)，定位輸出'-mat'存儲(chǔ)為matlab中的格式，配合.mat文件使用。44精選pptload功能：將文件中的數(shù)據(jù)導(dǎo)入工作空間格式：load文件名45精選ppt(1)變量替換，將高階方程降為一階常微分方程組(2)再使用ode系列求解器3、高階常微分方程ode系列求解器只能求解一階常微分方程(組),解高階常微分方程時(shí)，46精選ppt例3求二階常微分方程在[1,6]上的解，輸出步長(zhǎng)0.1，要求ε<1×10-5

[x,y]=ode45(odefun,tspan,y0,options)47精選ppt二階常微分方程一階常微分方程組48精選ppt

常微分方程中同時(shí)包含有快變分量和慢變分量，其變化速度相差非常大的量級(jí)，這種特點(diǎn)在數(shù)學(xué)上稱(chēng)為剛性；

描述這種過(guò)程的常微分方程組稱(chēng)為剛性方程組。4、剛性方程(Stiffproblem)49精選ppt剛性問(wèn)題兩個(gè)子系統(tǒng)效應(yīng)速度相差懸殊，構(gòu)成剛性系統(tǒng)?？臻g飛行器自動(dòng)控制系統(tǒng)，包含兩個(gè)子系統(tǒng)：控制飛行器質(zhì)心運(yùn)動(dòng)的系統(tǒng)，相對(duì)變化緩慢；控制飛行器運(yùn)動(dòng)姿態(tài)的系統(tǒng)，