1、Email: JPhone: 02985583997第五章 偏微分方程數值解 Numerical Methods for Partial Differential EquationsNumerical Methods for Partial Differential Equations&5.1 5.1 偏偏微分方程簡介微分方程簡介 &5.2 5.2 離散化公式離散化公式 &5.3 幾種常見偏微分方程的離散化計算幾種常見偏微分方程的離散化計算&5.4吸附床傳熱傳質模型中偏微分方程求解吸附床傳熱傳質模型中偏微分方程求解Email: JPhone: 02985583997本章要求本章要求教學目的教學

2、目的 講解： 偏微分方程離散格式及求解的一般過程教學要求教學要求E熟記 一階及二階偏微分方程的離散格式；C 精通 用EXCEL迭代對偏微分方程求解； 探索 用兩數組交替更新的辦法進行編程求解；F 延伸 對化學反應工程中物理場的模擬進行嘗試。教學重點教學重點各種偏微分方程的離散與求解EXCEL 循環迭代問題教學難點教學難點 特殊邊界條件的引入與應用Email: JPhone: 029855839975. 1 偏微分方程簡介&偏微分方程如果一個微分方程中出現多元函數的偏導數，或者說如果未知函數和幾個變量有關，而且方程中出現未知函數對幾個變量的導數，那么這種微分方程就是偏微分方程。在化工或化學動態模

3、擬方程中，常常有一個自變量是時間，其它的自變量為空間位置。如果只考慮一維空間，則只有兩個自變量；如果考慮兩維空間，則有3個自變量。 許多化工過程均是通過對偏微分方程的求解進行工藝參數的確定或數值模擬。Email: JPhone: 029855839975.1 偏微分方程簡介&偏微分方程的分類22222( )( )( )( )( )( )( )0uuuuuabcdefugx yxyxy 線性微分方程 Linear partial differencial equation擬線性微分方程 Quasilinear partial differencial equation非線性微分方程 Nonlin

4、ear partial differencial equation,x y 1, ,/,nx yuxy , ,/,nx yuxy Email: JPhone: 029855839975.1 偏微分方程簡介&數學上的分類：橢圓方程 Elliptic拋物線方程 Parabolic雙曲線方程 Hyperbolic&物理實際問題的歸類：波動方程(雙曲型）一維弦振動模型：熱傳導方程（拋物線型）一維線性熱傳導方程拉普拉斯方程（橢圓型）穩態靜電場或穩態溫度分布場）240bac 240bac 240bac 22222uutx 22uutx 22220uuxyEmail: JPhone: 02985583997

5、5.1 微分方程的求解思路&求微分方程數值解的一般步驟：Step1區域剖分：首先按一定規則將整個定義域分成若干小塊Step2微分方程離散：構造離散點或片的函數值遞推公式或方程Step3初始、邊界條件離散：根據遞推公式，將初值或邊界值離散化，補充方程，啟動遞推運算Step4 數值解計算：求解離散系統問題&微分方程的定解問題 離散系統的求解問題Email: JPhone: 029855839975.2 離散化公式&將自變量在時間和空間上以一定的間隔進行離散化，則應變量就變成了這些離散變量的函數。&一階偏導的離散化公式一般采用歐拉公式表示&有時為了保證系統和穩定性，對時間的差分往往采用向后公式, ,

6、( , , , )ni j kt n t x i x yj y z k zuu t x y z 1, , ,1, , ,1, , ,1, ,nni j ki j kt n t x i x yj y z k znnij ki j kt n t x i x yj y z k znni jki j kt n t x i x yj y z k znni j ki j kt n t x i x yj y z k zuuuttuuuxxuuuyyuuuzx 1, , ,(1),nni j ki j ktnt x i x yj y z k zuuutt Email: JPhone: 029855839975

7、.2 離散化公式P對于二階偏導，我們可以通過對泰勒展開式處理技術得到下面離散化計算公式：112, , , ,22,21, , ,1, ,22,2,1, ,1, ,22,222()2()2()nnni j ki j ki j kt n t x i x yj y z k znnnij ki j kij kt n t x i x yj y z k znnni jki j kij kt n t x i x yj y z k zt n tuuuuttuuuuxxuuuuyyuz , ,1, , ,12,2()nnni j ki j ki j kx i x yj y z k zuuuz Email: JP

8、hone: 029855839975.2 離散化公式推導&將uk+1在uk處按二階泰勒式展開：&將uk-1在uk處按二階泰勒式展開：&二式相加得：22312()2!kkkkuuhuuhO hxx 22312()2!kkkkuuhuuhO hxx 211222()kkkuuuuxx Email: JPhone: 029855839975.3幾種常見偏微分方程的離散化計算&1 1、 波動方程波動方程 其中： 為初值條件 為邊值條件 當該波動方程只提供初值條件時，稱此方程為波動方程的初值問題，二者均提供時稱為波動方程的混合問題。2222200012( , )( ),( )( ),( )ttxx lu

9、uaf x ttxuuxxtut ut 00012( ),( )( ),( )ttxx tuuxxtutut Email: JPhone: 029855839975.3.1 波動方程求解對于初值問題，是已知t=0時，u與 依賴于x的函數形式，求解不同位置，不同時刻的u值。而 u是定義在 的二元函數，即上半平面的函數。 對于混合問題除初值外，還有邊值。是已知初值及x=0及x=l 時u依賴于t的函數，求解不同位置x，不同時刻的u值。此時u是定義在 的帶形區域上的二元函數。ut 0,tx xt 0a)初值問題初值問題tx0lb)混合問題混合問題Email: JPhone: 029855839975.

10、3.1 波動方程求解2222200012( , )( ),( )( ),( )ttxx luuaf x ttxuuxxtut ut 1121122(i1,2,m-1)22( , ) ()()(1,2,)nnnnnniiiiiiuuuuuuaf x ttxn 22212221211222()()()(22)()( , )()()()nnnnniiiiitttuauauauuxf x txxx 100012(),() (i1,2,m)(),() (n1,2,)iiinnmuuuj xixtun tun t 方程離散化整理可得：邊界條件初始條件離散化xxinniu1niu 1niu Email: J

11、Phone: 029855839975.3.1 波動方程求解&例例5.1:5.1: 用數值法求解下面偏微分方程。 0012()3150,30,2501,0WnWjxttTtxTtttxx 110.01,0.1nnjjnnjjttttttxxx 110.020.680.3nnnjWjjtTtt 112()3nnnnjjjjjWntttttTx 此微分方程，是在不考慮流體本身熱傳導時的套管傳熱微分方程.由計算結果可知，當計算的時間序列進行到7272時，傳熱過程已達到穩態，各點上的溫度已不隨時間的增加而改變。如果改變套管長度或傳熱系數，則達到穩態的時間亦會改變。Email: JPhone: 0298

12、55839975.3.2 一維流動熱傳導方程 與波動方程的情形類似，用差商近似代替偏商，可以得到一維流動傳熱傳導方程的混合問題的差分方程，以其解作為流動傳熱傳導方程的近似解。222001( , ) 00 00 0 tx lxuuuabf u t (xl , t)txxu(x), (xl)u(t)xu (t) 0 (t) 12111201012(,)()() (i1 2)0 (0,1,2,)(nnnnnnniiiiiiiinnmmnuuuuuuuabf i x n ttxxui x, ,muunxun ) (0,1,2,)tn 1, , ,21, , ,1, ,22,2,1, ,1, ,22,2

13、2,2()2()nni j ki j kt n t x i x yj y z k znnnij ki j kij kt n t x i x yj y z k znnni jki j kij kt n t x i x yj y z k zt n t x i x yj y z kuuuttuuuuxxuuuuyyuz , ,1, , ,122()nnni j ki j ki j kzuuuz 2 2、一維流動熱傳導方程的混合問題、一維流動熱傳導方程的混合問題Email: JPhone: 02985583997將上式進行處理得到： 該式是顯式格式。只要保證式中各項系數大于零，一般情況下是穩定的，可以

14、獲得穩定的解。 分析上式可以發現，當為了提高數值精度取適當小的x 時，最有可能小于零的系數是 uin的系數，若要保證此項系數大于零，此時t必須相應地更小,會導致計算量將大大增加，這是顯式格式的缺點，為了克服此缺點，下面提出一種隱式格式： 偏微分方程在 點上進行離散化，且對時間的偏微分采用向后歐拉公式得到原偏微分方程的離散化公式：111111211122(,(1)()nnnnnnniiiiiiiuuuuuuuabf i x nttxx 122211222(,)()(12)()()()nnnniiiitttttutf i x n tabuabuauxxxxx (,(1)i x nt 5.3.2 一

15、維流動熱傳導方程Email: JPhone: 02985583997 從圖5-3中可見要由初值及邊界條件一排一排推上去是不行的，需解線性方程組,同時添上二邊界條件： 正好共有m+2個方程，同時有m+2個變量，就能解出n+1排上各點值。這樣，每解一個線性方程組，就可以往上推算一排點的u值，雖然引入了方程組的求解，有可能增加計算量，但由于隱式格式無條件穩定，t的取法與x 無關，可以少計算許多排節點上的u 值，相應于顯式格式來說，最終反而節省了計算量。 111011(1), nnnmmuntuu 5.3.2 一維流動熱傳導方程Email: JPhone: 02985583997&例5.2 考慮縱向導

16、熱的套管換熱器內管各點溫度分布微分方程： 解：首先根據前面的知識，將所求 的方程離散化：代入微分方程并化簡得：分析上式可知，如果知道了某一時刻的各點t，（，（j=0,1,2.10,11），），就可以求下一時刻的各點溫度值t(j=1,2.10),現在已經知道了零時刻管內各點的溫度分布及入口處在任何時刻的溫度，如想求下一時刻的溫度值，根據上面的離散化計算公式，還需知道在j=11處的溫度，這個溫度可利用給定的邊界條件離散化求得： 有了以上各式，上面的微分方程就可以求解了。220012()0.0013150,30,3001,0WnWjxtttTtxxTtttxx 11211222()0.01,0.1n

17、njjnnjjnnnjjjttttttxxttttxxx 1110.020.0010.6780.301nnnnjWjjjtTttt 110nnjjxtttxx 5.3.2 一維流動熱傳導方程Email: JPhone: 029855839975.3.3 穩態導熱/擴散方程&3 3、穩態導熱、穩態導熱/ /擴散方程擴散方程 在化工導熱及擴散過程中，沒有物流的流動，僅靠導熱及擴散進行熱量及質量的傳遞。如果此時系統達到穩定狀態，也就是說系統中每一個控制單元的各項性質如溫度、濃度等不再隨時間的改變而改變，系統中的各種性質只與其所處的位置有關，利用化工知識，我們可以得到下面二維、三維的穩態導熱或擴散偏微

18、分方程： 二維： 三維：二維的穩態導熱或擴散偏微分方程又稱調和方程。22220uuxy2222220uuuxyzuu un ()uun 常見有三種邊界條件： 第一類邊界條件： 第二類邊界條件： 第三類邊界條件：Email: JPhone: 02985583997&離散化公式：&取 ，經化簡得：&外節點（邊界節點）和內節點&求解方法劃分網格建立節點離散方程迭代求解（或解稀疏方程組）1,1,1,122220()()iji jiji ji ji juuuuuuxyxy ,1,11,11()4i jiji jiji juuuuu xy求解區域N 節點邊界五點格式示意圖5.3.3穩態導熱/擴散方程求解E

19、mail: JPhone: 029855839975.3.3穩態導熱/擴散方程求解 常用的3種迭代格式： (1)同步迭代： (2)異步迭代： (3)超松弛迭代：當計算范圍R 為 矩陣區域，x方向m等分，y方向n等分，最佳松弛因子為：由數學知識可知，用這些迭代法求解上面的偏微分方程均收斂。 (1)( )( )( )( ),1,11,11()4kkkkki jiji jiji juuuuu (1)( )( )(1)(1),1,11,11()4kkkkki jiji jiji juuuuu ( )( )(1)(1)1,11,1(1)( ),1()4(1)kkkkiji jiji jkki ji ju

20、uuuuuwuw u (,)axb cyd 22coscos112wmn 緊湊迭代Email: JPhone: 029855839975.3.3 穩態導熱/擴散方程求解例例5.35.3 ：處于傳熱平衡狀態的某保溫，假設其形狀為長方體，在x，y兩個方向上存在熱傳導，且導熱系數相等，已知邊界溫度分布如下圖所示： 解：取某一微元進行能量衡算，由于已達傳熱平衡狀態，故可得： 傳導入熱量-傳導出熱量=0 1xy10(1,1)24( , ) 30 27001,0t x yyyx 1( , ) 30 27001,0t x yxxy 2( , )30001,0tx yyx3( , ) 30001,1t x y

21、xy ()yttyyx zy ()xttxxy zx ytx zy xty zx xyz22220 xyttxyEmail: JPhone: 029855839975.3.3 穩態導熱/擴散方程求解&Microsoft Excel 迭代計算公式中的循環引用在“工具工具”菜單上，單擊“選項”，再單擊“重新計算”選項卡。 選中“迭代”復選框。 若要設置 Microsoft Excel 進行重新計算的最大次數，請在“最多迭代次數”框中鍵入迭代次數。迭代次數越高，Excel 用于計算工作表的時間越多。 若要設置兩次迭代結果之間可以接受的最大誤差，請在“最大誤差”框中鍵入所需的數值。數值越小，結果越精確

22、，Excel 用于計算工作表的時間也越多。Email: JPhone: 029855839975.4 吸附床傳熱傳質模型中偏微分方程求解實例&5.4.1 5.4.1 基本設定及假設基本設定及假設 1.吸附器結構參數的設定 上圖所示的是套筒式吸附器，該吸附器的有效長度為L，其有效內徑為D,環隙寬度為，吸附器壁厚為b。導熱流體通過環隙將熱量傳入或傳出吸附器，吸附質通過吸附器上端的小管進入或離開吸附器。吸附器結構示意圖吸附器結構示意圖 DL熱流體Email: JPhone: 029855839975.4.1 基本設定及假設2.吸附床外流體傳熱的一些基本假設： 1). 忽略流體在環隙寬度上的溫度梯度；

23、 2). 忽略熱損失； 3). 忽略吸附器壁厚b上的溫度梯度，用集中參數法求取吸附器壁面溫度。.吸附床內傳熱傳質的一些基本假設： 1). 吸附床內的吸附質氣體處于氣滯狀態； 2). 忽略蒸發器、冷凝器和吸附床之間的壓力差； 3). 吸附床內各計算微元內達到吸附平衡。吸附量可利用回歸方程計算； 4. 吸附熱利用微分吸附熱，隨吸附量和吸附溫度的改變而改變；比熱采用有效比熱，亦隨溫度改變，但在計算微元內，可認為是常數；5. 床層活性炭導熱系數采用當量導熱系數，可由實驗測量得到。 Email: JPhone: 029855839975.4.2流體傳熱模型的建立在軸方向上取一環隙微元，作能量分析如下：1

24、.流體通過流動流入環隙微元的能量為 2.流體通過流動流出環隙微元的能量3.流體熱傳導在x 處的熱量導入 7 總能量平衡方程fffTTTxxxxx 流體傳熱微元模流體傳熱微元模型型 infffpffquSCT ()fou tfffpffTquSCTxx outqinqxqfxffTqSx 其中： f 流體的密度 uf 環隙的流體速度， Sf 環隙的橫截面積，Cpf流體的比熱。4. 流體熱傳導在x+ x處的熱量導入 5. 微元體傳遞給吸附床的熱量 qt 6. 微元體內的能量變化率 為流體的橫截面積xxq(/)()fFxxffTTxqxSxx ()tffwqhxDTT cq()fpfffcCSx T

25、qt 221()ffffffwffffTuThDTTTtxSx fffpfC fSEmail: JPhone: 029855839975.4.3 吸附床內吸附劑傳熱傳質模型的建立&吸附床內發生著熱量和質量的傳遞，但質量的傳遞是建立在熱量傳遞基礎上的，故只要建立熱量傳遞方程，就可以根據平衡吸附量方程求出各處的吸附量。吸附床內的熱量傳遞主要以熱傳導為主，既有經向的熱傳導，也有軸向的熱傳導，為了便于建模分析，選取如圖所示的吸附床微元體，進行衡算： 2bxaTqrrx x+ x x xrr+ r吸附床內傳熱傳質微元體1.軸向導入熱量 :2.軸向導出熱量3.徑向導入熱量 4.徑向導出熱量 5.微元體內的

26、能量變化率 其中 為吸附床層內的有效比熱。6.總能量平衡方程 222()bbxxaTTqrrxxx 2braTqrxr 222 ()()bbdraTTqrrxrrr 2bceffaTqrrx Ct ()effpapbPmCCmCHT 222211bbbbeffTTTTtrrxraeffaeffC Email: JPhone: 029855839975 .4 .4 吸附器內/外無量綱化方程吸附器內/外無量綱化方程2,021()()fwwwfwwbswswshDThDTTTTTxSSx ,sswpsssC 221()ffffffwffffTuThDTTTtxSx fffpfC 222211bbbb

27、effTTTTtrrxraeffaeffC *2000000,ffwbininintxrxrtLRLTTTTTTTTTTTT 無量綱化處理無量綱化處理Email: JPhone: 029855839975 .4 .4 吸附器內吸附器內/ /外無量綱化方程外無量綱化方程2*2()PeKA Bitxx 20*2()sswsfswsfDrBiBiBiBitx 22222,fffffbbefffbbffffwsswsfsswswuC LhDLSDrPeBiKASShDLhDLLDrBiBiCSSR 22*2*2*bDrCCtxrrr整理可得：其中：Email: JPhone: 029855839975 .4 .4 吸附器內吸附器內/ /外無量綱化方程外無量綱化方程*(,0)0(,0)0(,0)0 x tx tx rt*0*1(0,)1 0(1,)(1 ,) 0 xxxtxxtxtx *1000 0