MATLAB偏微分方程求解

題目：用MATLAB求解偏微分方程主講人：班級：時間：1整理課件根底知識預習微分方程的ＭＡＴＬＡＢ求解包含１：常微分方程的ＭＡＴＬＡＢ求解〔上 節(jié)課已經講過〕這里不再贅述。２：偏微分方程的ＭＡＴＬＡＢ求解〔本 次教學內容〕2整理課件偏微分方程概念偏微分方程〔PartialDifferentialEquation，簡稱PDE〕指含有未知函數(shù)及其偏導數(shù)的方程。描述自變量、未知函數(shù)及其偏導數(shù)之間的關系。偏微分方程分為①線性偏微分方程式與②非線性偏微分方程式，常常有幾個解而且涉及額外的邊界條件。常微分方程：在微分方程中，假設自變量的個數(shù)只有一個的微分方程。偏微分方程：自變量的個數(shù)有兩個或兩個以上的微分方程。3整理課件求解偏微分方程的方法求解偏微分方程的數(shù)值方法：1.有限元法〔FiniteElementMethod,FEM〕---hp-FEM2.有限體積法〔FiniteVolumeMethod,FVM〕3.有限差分法〔FiniteDifferenceMethod,FDM〕。其它：廣義有限元法〔GeneralizedFiniteElementMethod,FFEM〕、擴展有限元法〔eXtendedFiniteElementMethod,XFEM〕、無網(wǎng)格有限元法〔MeshfreeFiniteElementMethod〕、離散迦遼金有限元法〔DiscontinuousGalerkinFiniteElementMethod,DGFEM〕等。4整理課件MATLAB解偏微分方程MATLAB提供了兩種方法解決PDE問題：①pdepe()函數(shù)，它可以求解一般的PDEs，具有較大的通用性，但只支持命令行形式調用。②PDE工具箱，可以求解特殊PDE問題，PDEtool有較大的局限性，比方只能求解二階PDE問題，并且不能解決偏微分方程組，但是它提供了GUI界面，從繁雜的編程中解脫出來了，同時還可以通過File->SaveAs直接生成M代碼使用pdeval()直接計算某個點的函數(shù)值？？？5整理課件一般偏微分方程組(PDEs)的MATLAB求解直接求解一般偏微分方程(組)，它的調用格式為sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t)問題描述函數(shù)初值條件邊界條件輸出參數(shù)自變量參數(shù)6整理課件【輸入?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ù)7整理課件【輸入?yún)?shù)】〔2〕@pdeic：是PDE的初值條件，必須化為下面的形式我們使用下面的簡單的函數(shù)來描述為u0=pdeic(x)8整理課件【輸入?yún)?shù)】〔3〕@pdebc：是PDE的邊界條件描述函數(shù)，必須先化為下面的形式于是邊值條件可以編寫下面函數(shù)描述為[pa,qa,pb,qb]=pdebc(x,t,u,du)其中a表示下邊界，b表示下邊界9整理課件【輸入?yún)?shù)】〔4〕m：就是對應于(式1)中相關參數(shù)x,t：就是對應于(式1)中自變量10整理課件【輸出參數(shù)】sol：是一個三維數(shù)組，sol(:,:,i)表示ui的解，換句話說uk對應x(i)和t(j)時的解為sol(i,j,k)11整理課件實例講解〔題目〕例：12整理課件初值條件邊界條件13整理課件實例講解〔解法〕【解】第一步根據(jù)〔1〕對照給出的偏微分方程，那么原方程可以改寫為14整理課件輸入?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(1)-u(2);s=[-1;1].*(exp(5.73*temp)-exp(-11.46*temp));15整理課件輸入?yún)?shù)〔2’〕初值條件初值條件改寫為%%初值條件函數(shù)functionu0=pdeic(x)u0=[1;0];16整理課件輸入?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];17整理課件〔4'〕主調函數(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','PDEDemo——byMatlabsky')%創(chuàng)立個窗口，窗口名字是name后邊的名字'NumberTitle','off'是關掉默認顯示名字。subplot(211)surf(x,t,sol(:,:,1))%sol(:,:,i)表示ui的解title('TheSolutionofu_1')xlabel('X')ylabel('T')zlabel('U')subplot(212)surf(x,t,sol(:,:,2))%sol(:,:,i)表示ui的解title('TheSolutionofu_2')xlabel('X')ylabel('T')zlabel('U')18整理課件19整理課件PDEtool求解特殊PDE問題MATLAB的偏微分工具箱(PDEtoolbox)可以比較標準的求解各種常見的二階偏微分方程(特殊二階的PDE)20整理課件典型偏微分方程的描述

21整理課件

22整理課件

〔3〕雙曲線型偏微分方程的一般形式

23整理課件

〔4〕特征值型偏微分方程的一般形式，注意它是〔1〕的變形，不能算獨立的一類24整理課件MATLAB采用有限元的方法求解各種PDEMATLAB為我們提供一個pdetool(在commandwindow中鍵輸pdetool翻開)的交互界面，可以求解二元偏微分u(x1,x2)(注意只能求解二元)。方程的參數(shù)由a、c、d和f確定，求解域由圖形確定，求解域確定好后，需要對求解域進行柵格化(這個是自動)。25整理課件

26整理課件

27整理課件偏微分方程邊界條件的描述Dirichlet(狄利克萊)條件Neumann(紐曼)條件28整理課件

29整理課件

30整理課件求解實例

31整理課件【解】由給定的PDE，可以得出d=1,c=1,a=2,f=1032整理課件step1:點擊工具欄的【PDE】按鈕，如下輸入PDE的參數(shù)，注意選擇Hyperbolic

33整理課件step2:繪制求解域

對坐標軸的操作可以在【Options】主菜單中操作，包括設置網(wǎng)格、坐標系范圍等(1)【Options】->AxisLimits設置如下34整理課件

35整理課件

(2)點擊工具欄上的第三個按鈕【繪制橢圓】，任意繪制一個橢圓，雙擊橢圓，設置如下36整理課件重復上面的操作，參數(shù)如下

37整理課件得到

38整理課件

(3)在setformula中如下輸入，“+〞表示求并集，“-〞表示求差集，注意沒有直接求交接的操作符39整理課件step3:邊界條件和初值條件初值條件可以通過【Solve】->【Parameters…】設置邊值條件設置如下(1)點擊工具欄的第6個按鈕【區(qū)域邊界】，顯示如下40整理課件

(2)【Boundary】->【RemoveAllSubdomainBorders】移除所有子域的邊界，將得到所有子域合并成一個求解域(3)【Boundary】->【SecifyBoundaryConditons…】設置邊界如下，注意我們這里只有Dirichlet條件41整理課件step4:生成使用有限元方法求解方程所需的柵格點擊工具欄的第8/9個按鈕，對求解域生成柵格，屢次點擊可以在原來根底上繼續(xù)細化柵格，直到自己覺得滿意為止，當然可以通過【Mesh】主菜單進行精確控制42整理課件

step5:求解方程點解工具欄的第10個按鈕“=〞【求解方程】

step6:求解結果繪圖點擊第11個按鈕【繪制圖形】，里面的選項很豐富，可以繪制等高線等好多，甚至播放動畫，具體大家可以自己慢慢摸索43整理課件動畫播放設置：(1)【Solve】->【Parameters】設置適宜的時間向量Time(2)【Plot】->【Parameters】選中【Animation】，點擊后面的【Options