MATLAB偏微分方程求解（課堂PPT）

1、.1題目：用MATLAB求解偏微分方程 主講人： 班級(jí)： 時(shí)間：.2基礎(chǔ)知識(shí)預(yù)習(xí) 微分方程的求解包含 ：常微分方程的求解（上節(jié)課已經(jīng)講過(guò)）這里不再贅述。 ：偏微分方程的求解（本次教學(xué)內(nèi)容）.3偏微分方程概念 偏微分方程（Partial Differential Equation，簡(jiǎn)稱(chēng)PDE）指含有未知函數(shù)及其偏導(dǎo)數(shù)的方程。描述自變量、未知函數(shù)及其偏導(dǎo)數(shù)之間的關(guān)系。 偏微分方程分為線(xiàn)性偏微分方程式與非線(xiàn)性偏微分方程式，常常有幾個(gè)解而且涉及額外的邊界條件。常微分方程：在微分方程中，若自變量的個(gè)數(shù)只有一個(gè)的微分方程。偏微分方程：自變量的個(gè)數(shù)有兩個(gè)或兩個(gè)以上的微分方程。.4求解偏微分方程的方法 求解偏

2、微分方程的數(shù)值方法： 1. 有限元法（Finite Element Method, FEM）- hp-FEM 2. 有限體積法（Finite Volume Method, FVM） 3. 有限差分法（Finite Difference Method, FDM）。 其它：廣義有限元法（Generalized Finite Element Method, FFEM）、擴(kuò)展有限元法（eXtended Finite Element Method, XFEM）、無(wú)網(wǎng)格有限元法（Meshfree Finite Element Method）、離散迦遼金有限元法（Discontinuous Galerkin

3、 Finite Element Method, DGFEM）等。.5MATLAB解偏微分方程 MATLAB提供了兩種方法解決PDE 問(wèn)題：pdepe()函數(shù)，它可以求解一般的PDEs，具有較大的通用性，但只支持命令行形式調(diào)用。 PDE 工具箱，可以求解特殊PDE 問(wèn)題，PDEtool 有較大的局限性，比如只能求解二階PDE 問(wèn)題，并且不能解決偏微分方程組，但是它提供了GUI界面，從繁雜的編程中解脫出來(lái)了，同時(shí)還可以通過(guò)File-Save As直接生成M代碼 使用pdeval()直接計(jì)算某個(gè)點(diǎn)的函數(shù)值？.6一般偏微分方程組(PDEs)的MATLAB求解 直接求解一般偏微分方程(組)，它的調(diào)用格式

4、為sol=pdepe(m,pdefun,pdeic,pdebc,x,t)(1) ) u,t,s(x,),()u,( cxuxuutxfxxxtuxutxmm，問(wèn)題描述函數(shù)初值條件邊界條件輸出參數(shù)自變量參數(shù).7【輸入?yún)?shù)】（1） pdefun：是PDE 的問(wèn)題描述函數(shù)，它必須換成下面的標(biāo)準(zhǔn)形式 PDE 就可以編寫(xiě)下面的入口函數(shù) c,f,s=pdefun(x,t,u,du) m,x,t就是對(duì)應(yīng)于(式1)中相關(guān)參數(shù)和自變量，du是u的一階導(dǎo)數(shù)，由給定的輸入變量即可表示出出c,f,s這三個(gè)函數(shù)(1) ) u,t,s(x,),()u,( cxuxuutxfxxxtuxutxmm，.8【輸入?yún)?shù)】（2）

5、pdeic：是PDE 的初值條件，必須化為下面的形式 我們使用下面的簡(jiǎn)單的函數(shù)來(lái)描述為u0=pdeic(x)00),(uutx.9【輸入?yún)?shù)】（3） pdebc：是PDE的邊界條件描述函數(shù)，必須先化為下面的形式 于是邊值條件可以編寫(xiě)下面函數(shù)描述為pa,qa,pb,qb=pdebc(x,t,u,du)其中a 表示下邊界，b 表示下邊界0 ) xu u,t,(x, f *u).t,q(x, u) t,p(x,.10【輸入?yún)?shù)】（4） m：就是對(duì)應(yīng)于(式1)中相關(guān)參數(shù) x,t：就是對(duì)應(yīng)于(式1)中自變量(1) ) u,t,s(x,),()u,( cxuxuutxfxxxtuxutxmm，.11【輸出

6、參數(shù)】 sol：是一個(gè)三維數(shù)組，sol(:,:,i)表示ui的解，換句話(huà)說(shuō)uk對(duì)應(yīng)x(i)和t(j)時(shí)的解為sol(i,j,k).12實(shí)例講解（題目） 例：.13初值條件邊界條件.14實(shí)例講解（解法） 【解】第一步根據(jù)（1）對(duì)照給出的偏微分方程，則原方程可以改寫(xiě)為.15輸入?yún)?shù)（1）目標(biāo)PDE函數(shù) % 目標(biāo)PDE函數(shù) function c,f,s=pdefun (x,t,u,du) c=1;1; f=0.024*du(1);0.17*du(2); temp=u(1)-u(2); s=-1;1.*(exp(5.73*temp)-exp(-11.46*temp);.16輸入?yún)?shù)（2）初值條件 初值

7、條件改寫(xiě)為 % 初值條件函數(shù) function u0=pdeic(x) u0=1;0;.17輸入?yún)?shù)（3）邊界條件 邊界條件改寫(xiě)為 % 邊界條件函數(shù) function pa,qa,pb,qb=pdebc(xa,ua,xb,ub,t) %a表示左邊界，b表示右邊界 pa=0;ua(2);qa=1;0; pb=ub(1)-1;0;qb=0;1;.18（4）主調(diào)函數(shù)clcx=0:0.05:1;t=0:0.05:2;m=0;sol=pdepe(m,pdefun,pdeic,pdebc,x,t);figure(numbertitle,off,name,PDE Demoby Matlabsky)%創(chuàng)建個(gè)窗

8、口，窗口名字是name后邊的名字NumberTitle,off是關(guān)掉默認(rèn)顯示名字。subplot(211)surf(x,t,sol(:,:,1)%sol(:,:,i)表示ui的解title(The Solution of u_1)xlabel(X)ylabel(T)zlabel(U)subplot(212)surf(x,t,sol(:,:,2)%sol(:,:,i)表示ui的解title(The Solution of u_2)xlabel(X)ylabel(T)zlabel(U).19.20PDEtool求解特殊PDE問(wèn)題 MATLAB的偏微分工具箱(PDE toolbox)可以比較規(guī)范的求

9、解各種常見(jiàn)的二階偏微分方程(特殊二階的PDE).21典型偏微分方程的描述 .22 .23 （3）雙曲線(xiàn)型偏微分方程的一般形式 .24 （4）特征值型偏微分方程的一般形式，注 意它是（1）的變形，不能算獨(dú)立的一類(lèi) .25MATLAB 采用有限元的方法求解各種PDE MATLAB 為我們提供一個(gè)pdetool (在command window 中鍵輸pdetool打開(kāi))的交互界面，可以求解二元偏微分u(x1,x2)(注意只能求解二元)。方程的參數(shù)由a、c、d和f確定，求解域由圖形確定，求解域確定好后，需要對(duì)求解域進(jìn)行柵格化(這個(gè)是自動(dòng))。.26 .27 .28偏微分方程邊界條件的描述 Dirich

10、let(狄利克萊)條件 Neumann(紐曼)條件.29 .30 .31求解實(shí)例 .32 【解】由給定的PDE，可以得出d=1,c=1,a=2,f=10.33step1:點(diǎn)擊工具欄的【PDE】按鈕，如下輸入PDE的參數(shù)，注意選擇Hyperbolic .34step2:繪制求解域?qū)ψ鴺?biāo)軸的操作可以在【Options】主菜單中操作，包括設(shè)置網(wǎng)格、坐標(biāo)系范圍等 (1)【Options】-Axis Limits設(shè)置如下.35 .36 (2)點(diǎn)擊工具欄上的第三個(gè)按鈕【繪制橢圓】，任意繪制一個(gè)橢圓，雙擊橢圓，設(shè)置如下 .37重復(fù)上面的操作，參數(shù)如下 .38得到 .39 (3)在set formula 中如

11、下輸入，“+”表示求并集，“-”表示求差集，注意沒(méi)有直接求交接的操作符.40step3:邊界條件和初值條件 初值條件可以通過(guò)【Solve】-【Parameters】設(shè)置 邊值條件設(shè)置如下 (1)點(diǎn)擊工具欄的第6 個(gè)按鈕【區(qū)域邊界】，顯示如下.41 (2)【Boundary】-【Remove All Subdomain Borders】移除所有子域的邊界，將得到所有子域合并成一個(gè)求解域 (3) 【Boundary】-【Secify Boundary Conditons】設(shè)置邊界如下，注意我們這里只有Dirichlet條件.42step4:生成使用有限元方法求解方程所需的柵格 點(diǎn)擊工具欄的第8/9 個(gè)按鈕，對(duì)求解域生成柵格，多次點(diǎn)擊可以在原來(lái)基礎(chǔ)上繼續(xù)細(xì)化柵格，直到自己覺(jué)得滿(mǎn)意 為止，當(dāng)然可以通過(guò)【Mesh】主菜單進(jìn)行精確控制.43 step5:求解方程點(diǎn)解工具欄的第10 個(gè)按鈕“=”【求解方程】 step6:求解結(jié)果繪圖 點(diǎn)擊第11 個(gè)按鈕【繪制圖形】，里面的選項(xiàng)很豐富，可以繪制等高線(xiàn)等好多，甚至播放動(dòng)畫(huà)，具體大家可以自己慢慢摸索.44動(dòng)畫(huà)播放設(shè)置： (1)【Solve】-【Parameters】設(shè)置合適的時(shí)間向量Time (2)【Plot】-【Parameters】選