1、numerical solutionnumerical solution數值解數值解 第四章第四章 導熱數值解法基礎導熱數值解法基礎 Principles of Numerical Principles of Numerical Methods of Conduction Methods of Conduction 主要內容主要內容核心知識點3第四章導熱問題的數值解法45Fluent Phoenics Parabolic Hyperbolic Or Elliptic Numerical Integration Code Series 64-0 4-0 概述概述一、求解導熱問題的三種基本方法一、

2、求解導熱問題的三種基本方法（一）理論分析法（一）理論分析法 TheoreticalTheoretical1. 1. 基本求解過程基本求解過程 在理論分析的基礎上，直接對微分方程在給定的定解條件下進行求解，在理論分析的基礎上，直接對微分方程在給定的定解條件下進行求解，得得到溫度與空間變量和時間變量之間的函數關系式，到溫度與空間變量和時間變量之間的函數關系式，通過這種關系式，獲得通過這種關系式，獲得物體內任意時刻的溫度值，稱為物體內任意時刻的溫度值，稱為分析解分析解，或，或理論解理論解、精確解精確解；2. 2. 主要特點主要特點（1 1） 能獲得問題的精確解，可為實驗、數值計算提供比較依據；能獲得

3、問題的精確解，可為實驗、數值計算提供比較依據；（2 2） 局限性很大，對復雜的問題無法求解；局限性很大，對復雜的問題無法求解；（3 3） 分析解分析解具有普遍性，各種因素的影響清晰可見。具有普遍性，各種因素的影響清晰可見。 7（二）實驗法（二）實驗法 ExperimentalExperimental1. 1. 基本求解過程基本求解過程 在傳熱學基本理論的指導下，采用對所研究對象的在傳熱學基本理論的指導下，采用對所研究對象的傳熱過程利用相似理論搭建試驗臺，進行測試計算。傳熱過程利用相似理論搭建試驗臺，進行測試計算。2.2. 主要特點主要特點 是傳熱學的基本研究方法是傳熱學的基本研究方法 （1 1

4、） 適應性不好適應性不好 （2 2） 費用昂貴費用昂貴8（三）數值法（三）數值法 NumericalNumerical1.1. 基本求解過程基本求解過程 把原來在時間、空間上連續的物理量的場，用把原來在時間、空間上連續的物理量的場，用有有限個離散點上的值的集合限個離散點上的值的集合來代替。來代替。 把原來在空間和時間上連續的物理量的場轉變為把原來在空間和時間上連續的物理量的場轉變為有限個離散的網格單元節點上的物理量的集合，通過有限個離散的網格單元節點上的物理量的集合，通過求解按一定方法建立起來的關于這些值的代數方程，求解按一定方法建立起來的關于這些值的代數方程，從而獲得離散點上被求物理量的值，

5、稱之為從而獲得離散點上被求物理量的值，稱之為數值解數值解，或或離散解離散解。建立控制方程、定解條件建立控制方程、定解條件 確定節點（區域離散化）確定節點（區域離散化）建立節點離散方程建立節點離散方程設立設立t 場的迭代初值場的迭代初值求解代數方程求解代數方程改進初場改進初場是否收斂是否收斂是是否否解的分析解的分析 連續連續t t10 2. 2. 主要特點主要特點（1 1）已成為傳熱學一個分支）已成為傳熱學一個分支 計算傳熱學（數值傳熱學）計算傳熱學（數值傳熱學）；（2 2）在很大程度上彌補了分析法的缺點，適應性強，）在很大程度上彌補了分析法的缺點，適應性強， 特別對于復雜問題更顯其優越性；特別

6、對于復雜問題更顯其優越性；（3 3）與實驗法相比成本低）與實驗法相比成本低（如航天飛機）（如航天飛機）。 11有限差分法（有限差分法（Finite Difference ModelFinite Difference Model）有限元法（有限元法（Finite Element ModelFinite Element Model） 邊界元法（邊界元法（Boundary Element ModelBoundary Element Model） 分子動力學模擬（分子動力學模擬（Molecular DynamicsMolecular Dynamics）二、數值計算方法二、數值計算方法(i,j-1)(i

7、,j)(i,j+1)(i+1,j)(i-1,j)網格線：網格線：沿沿x和和y y方向分別按間距方向分別按間距 x和和 y y，用一系列與坐標軸平行的，用一系列與坐標軸平行的 網格線，把求解區域分割成許多網格線，把求解區域分割成許多 小的矩形網絡，稱為小的矩形網絡，稱為子區域子區域。節點：節點：網格線的交點（內、邊界網格線的交點（內、邊界節點節點）邊界節點：邊界節點：網格線與物體邊界的交點網格線與物體邊界的交點步長：步長：相鄰兩節點的距離，相鄰兩節點的距離， x x， y y均勻網格：均勻網格： x= x= y y微元體：微元體：以每一個節點為中心的小區域以每一個節點為中心的小區域（元體、控制容

8、積）（元體、控制容積）一、區域、時間的離散化一、區域、時間的離散化1. 1. 描述：描述：二維、穩態、常物性、無二維、穩態、常物性、無q qv v 、矩形域、矩形域4-1 4-1 建立離散方程的方法建立離散方程的方法二維二維13xyxyji(i,j)MN2. 2. 步驟步驟1 1）劃定網格線）劃定網格線2 2）確定節點）確定節點x x0 0+i+i x, yx, y0 0+j+j y y = (i , j)= (i , j)3 3）以每個節點為）以每個節點為中心劃分單元格。中心劃分單元格。節點溫度為所在微節點溫度為所在微元體的平均溫度元體的平均溫度。圖圖2 2 導熱問題數值求解示例導熱問題數值

9、求解示例(i-1,j)(i+1,j)(i,j+1)(i,j-1)144 4）劃分單元網格后，物體內溫度由原來連續函數成）劃分單元網格后，物體內溫度由原來連續函數成為有限個離散數值，溫度曲線成為階梯狀變化。為有限個離散數值，溫度曲線成為階梯狀變化。5 5）列節點代數方程，求解有限個離散溫度值）列節點代數方程，求解有限個離散溫度值i,j111114i, ji, ji , ji , jt( tttt) x x = = y y，( (i i, ,j j) )節點節點: :15二、建立離散方程的常用方法二、建立離散方程的常用方法（2 2）熱平衡法（重點）熱平衡法（重點） 將相鄰節點間的將相鄰節點間的溫度

10、分布視為線性溫度分布視為線性，利用熱，利用熱平衡原理，認為平衡原理，認為某一節點某一節點與周圍各節點區域之間與周圍各節點區域之間的總換熱量的總換熱量=0=0，進而建立離散方程。，進而建立離散方程。 適用于內節點、邊界節點。適用于內節點、邊界節點。（1 1）Taylor SeriesTaylor Series（泰勒級數）展開法（泰勒級數）展開法 應用泰勒級數展開式，把導熱微分方程中的應用泰勒級數展開式，把導熱微分方程中的各階各階導數（微分）導數（微分）用相應的用相應的差分差分表達式來代替，表達式來代替，進而建立離散方程。進而建立離散方程。 適用于內節點。適用于內節點。161. 1. Taylor

11、 SeriesTaylor Series（泰勒級數）展開法泰勒級數）展開法1 1）根據泰勒級數）根據泰勒級數(1715)(1715)展開式，用節點展開式，用節點( (i,ji,j）的溫度的溫度t ti,ji,j 表示表示 節點節點( (i+1,ji+1,j) )的溫度的溫度t ti+1,ji+1,j，展開式為，展開式為22334412341234i,ji,ji,ji,ji,ji,jttxtxtxttx( )xx!x!x! 2 2）移項整理并歸并可得（一階截差公式）移項整理并歸并可得（一階截差公式）141i,ji ,ji ,jtttO(x ) ()xx截斷誤差截斷誤差:未明確寫出的級數未明確寫出

12、的級數余項中的余項中的x的最低階數為的最低階數為1一階導數一階導數向前差分向前差分表達式表達式173 3）同理，用節點）同理，用節點( (i,ji,j) ) 的溫度的溫度t ti,ji,j 表示節點表示節點(i-1,j)(i-1,j)的的溫度溫度t ti-1,j i-1,j 22334412342234i,ji,ji,ji,ji,ji,jttxtxtxttx()xx!x!x! 式式(1) -(1) -式式(2)(2)：4 4）移項整理并歸并可得（一階截差公式）移項整理并歸并可得（一階截差公式）142i ,ji,ji ,jtttO(x ) ()xx1122i,ji,ji ,jtttO(x )xx

13、一階導數一階導數中心差分表達式中心差分表達式一階導數一階導數向后差分表達式向后差分表達式二階截差公式二階截差公式185 5）式）式(1)+(1)+式式(2)(2)，移項整理，可得，移項整理，可得，二階導數二階導數的中心差分：的中心差分：2112222i, ji , ji, ji , jttttO (x)xx2112222i,ji,ji,ji,jttttO( y )yy截斷誤差截斷誤差：未明確寫出的級數：未明確寫出的級數 余項中的余項中的x的最低階數為的最低階數為26 6）同樣可得：）同樣可得：0tt()()xxyy1 19vttttc()()()q()xxyyzz 舉例：舉例：二維、穩態、常物

14、性、無二維、穩態、常物性、無q qv v 導熱問題導熱問題A: A: 首先，寫出通用微分表達式：首先，寫出通用微分表達式：22220ttxy111122220i,ji,ji,ji,ji,ji,jttttttxy化化簡簡20小小結結一階導數的各種差分形式一階導數的各種差分形式1. 1. 向前差分（以向前差分（以x向為例）向為例） 2. 2. 向后差分向后差分3 3. 中心差分中心差分1. 1. 向前差分向前差分2. 2. 向后差分向后差分3. 3. 中心差分中心差分1i,ji,ji,jttt()xx112i,ji,ji,jttt()xx1i,ji,ji,jttt()xx1i ,jkki,ji,j

15、ttt()112i ,ji ,jkki,jttt()1i ,jkki,ji,jttt()t t 對對坐坐標標t t 對對 2. 2. 熱平衡法熱平衡法0LPRPTPBP11111111i,ji,jLPi,ji,jRPi,ji,jTPi,ji,jBPttyxttyxttxyttxy xyxyjiMNLPTRB仍以二維、穩態、常物性、無仍以二維、穩態、常物性、無q qv v 導熱問題為例導熱問題為例由由傅里葉定律傅里葉定律：由由能量守恒定律能量守恒定律：111122220i,ji,ji,ji,ji,ji,jttttttxy22小結小結 物理問題的數值求解過程物理問題的數值求解過程建立控制方程及定解

16、條件建立控制方程及定解條件確定節點（區域離散化）確定節點（區域離散化）建立節點物理量的代數方程建立節點物理量的代數方程設立溫度場的迭代初值設立溫度場的迭代初值求解代數方程求解代數方程 解的分析解的分析改進初場改進初場是是否否相對誤差相對誤差1010-6 -6 10 10-3-3111kkiikkiikikkiikmaxmax ttttmaxtttmaxt23 仍以二維、穩態、常物性、無仍以二維、穩態、常物性、無q qv v導熱問題為例。導熱問題為例。1 1）在直角坐標中，其導熱微分方程為：）在直角坐標中，其導熱微分方程為：22220ttxy1i,j11122220i,ji,ji,ji,ji,j

17、ttttttxy代入上式，得節點方程為：代入上式，得節點方程為：1. 1. 內節點離散方程的建立內節點離散方程的建立2222i,ji,jttxy和2 2）將二階差分式）將二階差分式4-2 4-2 穩態導熱的數值計算穩態導熱的數值計算2D2D2411114i, ji, ji , j, jjii ,ttttt3 3）若）若 x= x= y y，化簡得，化簡得 ：重要說明：重要說明： 所求節點溫度前的系數所求節點溫度前的系數所有相鄰節點溫度前的系數。所有相鄰節點溫度前的系數。但不包括熱流但不包括熱流( (或熱流密度或熱流密度) )前的系數前的系數。 這一結論也適用于邊界節點。這一結論也適用于邊界節點

18、。 P P8989表表4-14-1( (內部節點的導出在本書中簡單，在內部節點的導出在本書中簡單，在 數值計算中內部節點數值計算中內部節點的導出、具體方程還有許多的導出、具體方程還有許多) )與左鄰右與左鄰右舍關系如舍關系如何何?25 第一類邊界條件：處理比較簡單，因為已知邊界溫第一類邊界條件：處理比較簡單，因為已知邊界溫度，可將其以數值的形式加入到內節點的離散方程中，度，可將其以數值的形式加入到內節點的離散方程中，組成封閉的代數方程組，直接求解。組成封閉的代數方程組，直接求解。 第二類邊界條件或第三類邊界條件第二類邊界條件或第三類邊界條件：必須用熱平衡法，：必須用熱平衡法，建立邊界節點的建立

19、邊界節點的離散方程離散方程，邊界節點與內節點的離散方，邊界節點與內節點的離散方程一起組成封閉的代數方程組，才能求解。程一起組成封閉的代數方程組，才能求解。 為求解方便，將第二類、第三類邊界條件合并考慮，為求解方便，將第二類、第三類邊界條件合并考慮，用用q qww表示邊界上的熱流密度值或熱流密度表達式。表示邊界上的熱流密度值或熱流密度表達式。用用q qv v表示內熱源強度。表示內熱源強度。2. 2. 邊界節點離散方程的建立邊界節點離散方程的建立26qwxyqw(1) (1) 平直邊界上的節點平直邊界上的節點 P P8989 序號序號2 2111222 20i,ji,ji,ji,jfh xh xt

20、tt()tt111202i,ji,ji,ji,ji,ji,jwttttxyxyttxqyy wfi,jxy,qh(tt) 27(2) (2) 外部角點外部角點 序號序號3 321122022i,ji,ji,jv,i,jfh xtt()th xxtq112220222v,i,i,ji,jwi,ji,wjjttyyqxtxytxxqqywfi,jxy,qh(tt) xyqw課后作業：課后作業：將將q qww= =h h ( (t tf f - -t ti,ji,j) ) 帶入外部角點的溫帶入外部角點的溫度離散方程，化簡度離散方程，化簡28(3) (3) 內部角點內部角點29q qww的情況：的情況

21、：（1 1）第二類邊界條件：將）第二類邊界條件：將q qww=Const=Const，帶入上面各式即可，帶入上面各式即可 絕熱或對稱邊界條件？絕熱或對稱邊界條件？（2 2）第三類邊界條件：將）第三類邊界條件：將q qww= =h h ( (t tf f - -t ti,ji,j) )，帶入上面各式即可，帶入上面各式即可 （3 3）輻射邊界條件：）輻射邊界條件：44,()wfi jqTTconstqw或其他或其他30設有一矩形平板，如圖示。設有一矩形平板，如圖示。a=2ba=2b。在邊界在邊界x=0 x=0和和y=0y=0是是絕熱的，在絕熱的，在x=ax=a處給出處給出第三第三類類邊界條件，即給

22、定邊界條件，即給定h h和和t tf f 。而而y=y=b b處，是第一類邊界條處，是第一類邊界條件，件，即溫度為已知即溫度為已知，t=ct=c1111，c c1212，C C1515 。試寫出各節點的離散方程。試寫出各節點的離散方程。例例1 1解：解：采用均勻網格采用均勻網格y y = =x=b/2x=b/2，給各未知節點給各未知節點編號編號t t1 1，t t2 2， t t10 10 。 按節點所在位置和題目所示邊界條件，寫按節點所在位置和題目所示邊界條件，寫出各節點的離散方程。出各節點的離散方程。0battttttttt1234567810ccccc1111123451t9310bat

23、tttttttt1234567810ccccc1111123451t9123241113415423612563764581378598671014971010896242542442442422224123745891022fttttttcttttttttcttttttttcttttttttchbh b()tt.tth.b)tt.(t序序1522fh bct2bx 3212324111341542361256376458137859867101497101089156.7.8.9.101 2424244244242(2)2(4 2)2.2.3.4.2225.ffttttttctttttttt

24、cttttttttcttttttttchbhbtttth bh btttct (2)(4)444444421 10 00 0 0 0 010 2 00 0 0 0 01 02 10 0 0 0 00 1 101 0 0 0 00 01 02 1 0 0 00 00 1 10 1 0 00 00 0 102 1 00 00 0 01 10 10 00 0 00 1 010 00 0 00 0 2 1hbhb=12345678910tttttttttt111213140000ffhbhbttcccc系數系數矩陣矩陣變量變量矩陣矩陣常數常數矩陣矩陣 Abt332.2.節點方程組的求解節點方程組的求解

25、111 112 211221 122 2221 12 221 12 2.n nn niiiin nnnnta ta ta tcta ta ta tcta ta ta tcta ta tnn nna tc前已寫出所有內節點和邊界節點的溫度差分方程前已寫出所有內節點和邊界節點的溫度差分方程n n個未知節點溫度個未知節點溫度t t1 1、t t2 2、t t3 3.t tn n，n n個代數方程式：個代數方程式：代數方程組的求解方法：代數方程組的求解方法：直接解法、迭代解法直接解法、迭代解法34直接解法：直接解法：通過有限次運算獲得代數方程精確解通過有限次運算獲得代數方程精確解; ; 矩陣求逆矩陣求

26、逆、高斯消元法高斯消元法迭代解法迭代解法：先對要計算的場作出假設、在迭代計算過程中不：先對要計算的場作出假設、在迭代計算過程中不斷予以改進、直到計算結果與假定值的結果相差小于允許值。斷予以改進、直到計算結果與假定值的結果相差小于允許值。稱迭代計算已經收斂。稱迭代計算已經收斂。缺點：缺點：所需內存較大、方程數目多時不便、不適用于非線性問所需內存較大、方程數目多時不便、不適用于非線性問題（若物性為溫度的函數，節點溫度差分方程中的系數不再是題（若物性為溫度的函數，節點溫度差分方程中的系數不再是常數，而是溫度的函數。這些系數在計算過程中不斷更新）常數，而是溫度的函數。這些系數在計算過程中不斷更新）迭代

27、解法有多種：迭代解法有多種：簡單迭代（簡單迭代（JacobiJacobi 雅克比迭代）、雅克比迭代）、 高斯高斯- -賽德爾迭代賽德爾迭代、塊迭代、交替方向迭代等、塊迭代、交替方向迭代等35解解：二維穩態導熱問題。：二維穩態導熱問題。對內部節點，應用對內部節點，應用P P8989表表4-14-1序號序號1 1的公式的公式一各向同性材料的方形物體，其一各向同性材料的方形物體，其 為常量。已知各邊界的溫度為常量。已知各邊界的溫度如圖所示，且如圖所示，且xx= =yy。試用。試用高斯高斯- -賽德爾賽德爾迭代求其內部網格迭代求其內部網格節點節點1, 2, 3, 41, 2, 3, 4的溫度。的溫度。

28、例例2 21234t =240t =60t =60t =60 1001231240604ttt 1102141240604ttt 1103141240604ttt 1114231240604ttt36開始，假設開始，假設t t1 1(0)(0)= = t t2 2(0)(0) =120 =120; ; t t3 3(0)(0) = =t t4 4(0)(0) =80=80。帶入上式進行迭。帶入上式進行迭代，第代，第1515次迭代值匯總如下：次迭代值匯總如下：例例2 21234t =240t =60t =60t =60迭代次數迭代次數t1/t2/t3/t4/012012080801125126.

29、2581.2581.8752126.875127.1982.1982.3453127.345127.4282.4282.464127.46127.4882.4882.495127.49127.495 82.495 82.50374-3 4-3 非穩態導熱的數值計算非穩態導熱的數值計算-1D-1D一、顯式差分格式一、顯式差分格式2. 2. 物體溫度沿物體溫度沿x方向變化方向變化，物體在物體在x方向被分割為方向被分割為n段，段，所取步長所取步長 x ，其中第，其中第i個節個節點坐標為點坐標為i i x ，簡記為，簡記為i。3. 3. 物體中溫度隨時間變化物體中溫度隨時間變化，所取時間間隔步長，所取

30、時間間隔步長，其中，其中第第k個節點坐標為個節點坐標為k，簡記為，簡記為k。 表示物體中表示物體中k k時刻，時刻，i i x x位置處溫度：位置處溫度：kit22ttax 1. 1. 導熱微分方程式導熱微分方程式一維、非穩態、常物性、無一維、非穩態、常物性、無q qv v384.4. 內節點（內節點（i i, , k k）溫度對）溫度對x x的二階導數，的二階導數， 采用中心差分表達式應為：采用中心差分表達式應為：5.5. 溫度對時間的一階導數，取溫度對時間的一階導數，取向前差分向前差分：6.6. 將將4 4與與5 5帶入帶入1 1的非穩態導熱方程的非穩態導熱方程7.7. 內節點在計算時刻節點溫度內節點在計算時刻節點溫度的的差分方程差分方程121222iikkkii,ktttt()xx1iikki,kttt() 22tt