多層平板一維非穩(wěn)態(tài)導(dǎo)熱問題的數(shù)值解法

(吉林化工學(xué)院 石油化工學(xué)院，吉林 吉林 132022)

現(xiàn)實生活和實際生產(chǎn)中存在大量的非穩(wěn)態(tài)導(dǎo)熱問題，如室內(nèi)環(huán)境通過墻壁向室外的熱量傳遞過程、鍋爐內(nèi)部通過壁面向外界的熱量傳遞過程以及埋地原油管道向周圍土壤的散熱過程等。由于非穩(wěn)態(tài)導(dǎo)熱問題中物體內(nèi)部溫度隨時間而變化，因此求解起來比穩(wěn)態(tài)導(dǎo)熱要復(fù)雜得多。即使是最簡單的一維問題，最簡單的初始條件和邊界條件，其解析解的形式也相當(dāng)復(fù)雜，將其應(yīng)用于實際問題的求解中將存在一定的困難。在解決導(dǎo)熱問題時主要有理論分析法、實驗方法及數(shù)值計算方法等，其中數(shù)值解法隨著計算機硬件水平的大幅度提高和計算方法的完善，逐漸受到人們的重視。下面我們將針對多層平板一維非穩(wěn)態(tài)導(dǎo)熱問題的數(shù)值解法進行介紹。

1 提出問題及問題分析

1.1 提出問題

由多層材料組成的無限大平板，兩側(cè)分別置于溫度恒定的高溫和低溫環(huán)境中，兩側(cè)環(huán)境溫度分別為Th和Tl，且Th>Tl，各層材料的厚度及物理性質(zhì)如熱導(dǎo)率、比熱容、密度等已知，如圖1所示。由于高溫側(cè)環(huán)境溫度過高，無法測量高溫側(cè)表面溫度Twh，低溫側(cè)表面溫度Twl可實時監(jiān)測，可以獲得低溫側(cè)各個時刻的表面溫度。需確定各層材料中的溫度分布情況。

圖1 多層平板導(dǎo)熱示意圖

1.2 問題分析

對于以上問題的傳熱過程，熱量是由高溫環(huán)境通過多層無限大平板向低溫環(huán)境傳熱的過程，傳熱環(huán)節(jié)包括：高溫環(huán)境與平板高溫側(cè)的對流傳熱(高溫環(huán)境溫度Th已知，平板高溫側(cè)溫度t(0,τ)未知，對流傳熱系數(shù)hh未知)；平板高溫側(cè)通過各層平板向平板低溫側(cè)的導(dǎo)熱(各層平板的厚度及熱導(dǎo)率、比熱容、密度已知，平板內(nèi)及各平板的接觸面溫度未知)；平板低溫側(cè)與低溫環(huán)境間的對流傳熱(低溫環(huán)境溫度Tl已知，平板低溫側(cè)溫度t(δ,τ)已知，對流傳熱系數(shù)hl未知)?？烧J為垂直于X軸的任一平面內(nèi)溫度均勻，因此可看成是一維導(dǎo)熱問題，且溫度隨時間變化，為非穩(wěn)態(tài)導(dǎo)熱問題。

2 數(shù)學(xué)模型及邊界條件

2.1 數(shù)學(xué)模型

如圖1所示的導(dǎo)熱問題，可通過對直角坐標(biāo)系的導(dǎo)熱微分方程(1)去掉Y、Z兩個方向的溫度二階導(dǎo)數(shù)及內(nèi)熱源項，簡化得到微分方程(2)

(1)

(2)

2.2 邊界條件

上述三個公式中，t(x,τ)表示x位置處τ時刻的溫度，δ為多層平板的總厚度。

3 數(shù)值解法及穩(wěn)定性分析

3.1 數(shù)值解法

將一維非穩(wěn)態(tài)導(dǎo)熱微分方程(2)在空間及時間上進行離散處理，將等號左邊非穩(wěn)態(tài)項取向前差分，將等號右邊擴散項取中心差分，則有[1]。

(3)

上式可改寫為

(4)

兩層材料接觸點處溫度計算采用下式[2]：

(5)

上述為一維非穩(wěn)態(tài)導(dǎo)熱問題的顯式差分格式數(shù)值解法，其優(yōu)點是計算工作量較小，缺點是對時間步長和空間步長存在限制。其求解過程是從已知初始溫度分布出發(fā)，根據(jù)邊界條件依次求得以后各個時層上的溫度分布情況。

3.2 數(shù)值穩(wěn)定性分析

如果把計算區(qū)域用N個節(jié)點將其等分，則節(jié)點0和N為對流傳熱面，與微分形式的數(shù)學(xué)描寫相對應(yīng)的顯式離散形式為

(6)

(7)

3.3 程序設(shè)計

由于邊界節(jié)點溫度是隨時間而變化的，兩側(cè)的對流傳熱系數(shù)也是未知的，這就需要采用逐步逼近法編程求解平板內(nèi)的溫度分布情況。雖然邊界節(jié)點溫度是變化的，但是當(dāng)達到穩(wěn)定狀態(tài)后，平板內(nèi)及邊界節(jié)點溫度均會達到一恒定值。而低溫側(cè)邊界節(jié)點溫度是可以通過實驗測定的，如果知道高溫側(cè)邊界節(jié)點溫度，那么問題就可以進行求解。因此我們按照如下步驟進行編程：

(1)給定平板的初始條件。

(2)假設(shè)達到穩(wěn)態(tài)時的高溫側(cè)邊界節(jié)點溫度Twh(達到穩(wěn)態(tài)時，高溫側(cè)邊界節(jié)點溫度接近高溫環(huán)境溫度)，低溫側(cè)邊界節(jié)點溫度已知，可通過串聯(lián)熱阻疊加原則，計算出兩側(cè)的對流傳熱系數(shù)。

(3)將對流傳熱系數(shù)帶入式(6)、(7)確定邊界節(jié)點溫度。

(4)根據(jù)邊界條件，利用式(4)、(5)分別確定各層平板內(nèi)部以及兩層平板接觸面處的溫度值。

(5)循環(huán)執(zhí)行第(3)、(4)步，直到達到穩(wěn)定狀態(tài)，確定高溫側(cè)邊界節(jié)點溫度Twh1，若|Twh-Twh1|小于允許誤差，則退出并輸出數(shù)據(jù)；否則令Twh=Twh1，繼續(xù)執(zhí)行步驟(2)。

4 結(jié)論

上述數(shù)學(xué)模型及時間-空間的離散顯式格式，主要用來解決多層平板的一維非穩(wěn)態(tài)導(dǎo)熱問題，可以用來確定任意時刻平板內(nèi)的溫度分布情況，也可以用來計算各個時刻通過平板的導(dǎo)熱量。通過公式(5)可以有效的處理多層平板接觸處的節(jié)點溫度計算；公式(6)、(7)用來處理平板兩側(cè)均為對流傳熱的邊界節(jié)點計算問題；對于顯示差分格式，時間和空間步長要滿足一定的限制條件，文中也進行了詳細的分析，并給出了需滿足的計算條件；最后給出了兩側(cè)對流傳熱系數(shù)未知情況下的編程詳細步驟，為實際生產(chǎn)和生活解決此類問題提供了理論依據(jù)和解決方法。