MATLAB偏微分方程求解課件

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

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

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

4、組)，它的調(diào)用格式 為 sol=pdepe(m,pdefun,pdeic,pdebc,x,t) (1) ) u,t,s(x,),()u,( c x u x u utxfx x x t u x u tx mm ， 問題描 述函數(shù) 初值 條件 邊界 條件 輸出 參數(shù) 自變 量 參數(shù) 【輸入?yún)?shù)】（1） pdefun：是PDE 的問題描述函數(shù)，它必 須換成下面的標準形式 PDE 就可以編寫下面的入口函數(shù) c,f,s=pdefun(x,t,u,du) m,x,t就是對應于(式1)中相關參數(shù)和自變量， du是u的一階導數(shù)，由給定的輸入變量即可 表示出出c,f,s這三個函數(shù) (1) ) u,t,s(x,)

5、,()u,( c x u x u utxfx x x t u x u tx mm ， 【輸入?yún)?shù)】（2） pdeic：是PDE 的初值條件，必須化為下 面的形式 我們使用下面的簡單的函數(shù)來描述為 u0=pdeic(x) 00 ),(uutx 【輸入?yún)?shù)】（3） pdebc：是PDE的邊界條件描述函數(shù)，必 須先化為下面的形式 于是邊值條件可以編寫下面函數(shù)描述為 pa,qa,pb,qb=pdebc(x,t,u,du)其中a 表示下 邊界，b 表示下邊界 0 ) x u u,t,(x, f *u).t,q(x, u) t,p(x, 【輸入?yún)?shù)】（4） m：就是對應于(式1)中相關參數(shù) x,t：就是對

6、應于(式1)中自變量 (1) ) u,t,s(x,),()u,( c x u x u utxfx x x t u x u tx mm ， 【輸出參數(shù)】 sol：是一個三維數(shù)組，sol(:,:,i)表示ui的解， 換句話說uk對應x(i)和t(j)時的解為sol(i,j,k) 實例講解（題目） 例： 初值條件 邊界條件 實例講解（解法） 【解】第一步根據(jù)（1）對照給出的偏微分 方程，則原方程可以改寫為 輸入?yún)?shù)（1）目標PDE函數(shù) % 目標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(

7、1)-u(2); s=-1;1.*(exp(5.73*temp)-exp(- 11.46*temp); 輸入?yún)?shù)（2）初值條件 初值條件改寫為 % 初值條件函數(shù) function u0=pdeic(x) u0=1;0; 輸入?yún)?shù)（3）邊界條件 邊界條件改寫為 % 邊界條件函數(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; （4）主調(diào)函數(shù) clc x=0:0.05:1; t=0:0.05:2; m=0; sol=pdepe(m,pdefun,pdei

8、c,pdebc,x,t); figure(numbertitle,off,name,PDE Demoby Matlabsky)%創(chuàng)建個窗口，窗口名字 是name后邊的名字NumberTitle,off是關掉默認顯示名字。 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) xlabe

9、l(X) ylabel(T) zlabel(U) PDEtool求解特殊PDE問題 MATLAB的偏微分工具箱(PDE toolbox)可 以比較規(guī)范的求解各種常見的二階偏微分 方程(特殊二階的PDE) 典型偏微分方程的描述 （3）雙曲線型偏微分方程的一般形式 （4）特征值型偏微分方程的一般形式，注 意它是（1）的變形，不能算獨立的一類 MATLAB 采用有限元的方法求解各 種PDE MATLAB 為我們提供一個pdetool (在 command window 中鍵輸pdetool打開)的交 互界面，可以求解二元偏微分u(x1,x2)(注意 只能求解二元)。方程的參數(shù)由a、c、d和f 確定，

10、求解域由圖形確定，求解域確定好 后，需要對求解域進行柵格化(這個是自動)。 偏微分方程邊界條件的描述 Dirichlet(狄利克萊)條件 Neumann(紐曼)條件 求解實例 【解】由給定的PDE，可以得出 d=1,c=1,a=2,f=10 step1:點擊工具欄的【PDE】按鈕，如下輸入PDE 的參數(shù)，注意選擇Hyperbolic step2:繪制求解域 對坐標軸的操作可以在【Options】主菜單中操作，包括設 置網(wǎng)格、坐標系范圍等 (1)【Options】-Axis Limits設置如下 (2)點擊工具欄上的第三個按鈕【繪制橢 圓】，任意繪制一個橢圓，雙擊橢圓，設 置如下 重復上面的操作

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