丁 寧， 林 潔

(1. 中國民航大學(xué) 航空工程學(xué)院， 天津 300300； 2. 中國民航大學(xué) 中歐航空工程師學(xué)院， 天津 300300)

在船舶制造行業(yè)，一線人員經(jīng)常需要在雙層底、雙舷側(cè)以及其他密閉艙室內(nèi)實(shí)施焊接、檢測等作業(yè)。這類作業(yè)環(huán)境內(nèi)的溫度往往較高，傳熱方式以自然對流為主，此外，在建筑行業(yè)等其他領(lǐng)域，一線人員也同樣需要在這類環(huán)境中工作。為保護(hù)人體免遭損害，人員在進(jìn)入這類環(huán)境時(shí)往往需要穿著防護(hù)服，因此防護(hù)服的隔熱性能至關(guān)重要。采用數(shù)值方法研究防護(hù)服的隔熱性能，因其成本和靈活性方面的優(yōu)勢，越來越受到重視。

因隔熱性能與熱源密切相關(guān)，學(xué)者們根據(jù)不同的熱源特點(diǎn)，經(jīng)適當(dāng)假設(shè)后建立了多種數(shù)學(xué)模型，并相應(yīng)采用了多種數(shù)值方法進(jìn)行求解。這些熱源主要是石油化工行業(yè)的強(qiáng)熱流閃火環(huán)境和消防行業(yè)的中等熱流輻射環(huán)境[1]，傳熱方式主要是對流和輻射。其中，自然對流產(chǎn)生的換熱占比較小忽略不計(jì)，Torvi[2]對這類熱源條件進(jìn)行了研究，將輻射和對流換熱合并為單一系數(shù)并視為常量給出了簡潔的表述形式，后續(xù)較多研究也均采納了這一表述[3-5]。為考慮自然對流時(shí)的傳熱情況，Mell等[6]采用Howard模型給出的經(jīng)驗(yàn)公式[7]，在常熱流條件下，利用Nusselt數(shù)計(jì)算了自然對流換熱系數(shù)進(jìn)而獲得了織物溫度分布的數(shù)值解，并與解析解進(jìn)行了對比。然而，Howard模型通常適用于常壁溫或常熱流的對流邊界條件，對于變壁溫、變熱流的對流邊界條件是否適用，并未有明確結(jié)論，也尚未見有其他相關(guān)研究成果，而在高溫環(huán)境中，情況完全不同，其傳熱方式以自然對流為主，在達(dá)到穩(wěn)態(tài)前，自然對流換熱系數(shù)一直在變化，對流邊界上的壁溫和熱流也一直在變化，采用Howard模型難以計(jì)算。此外，Su[8]采用了Howard模型計(jì)算自然對流邊界條件下的冷卻問題，從與實(shí)驗(yàn)結(jié)果的對比上看差異較大。從上述分析可以看出，盡管熱防護(hù)服的邊界條件是影響模型精度的重要因素[2]，但對非穩(wěn)態(tài)條件下自然對流換熱系數(shù)的研究并未引起足夠的重視。如陳揚(yáng)等[9]在自然對流傳熱條件下對織物的非穩(wěn)態(tài)傳熱進(jìn)行了模擬，并與實(shí)驗(yàn)值進(jìn)行了對比，其自然對流換熱系數(shù)則直接指定為10 W/(m2·K)，并未給出相關(guān)確定依據(jù)。Piotr等[10]研究了輻射傳熱占優(yōu)的熱防護(hù)服非穩(wěn)態(tài)傳熱問題，其自然對流換熱系數(shù)也未加說明而直接指定。目前，在防護(hù)服隔熱性能非穩(wěn)態(tài)數(shù)值預(yù)報(bào)方面，強(qiáng)迫對流傳熱和輻射傳熱均有較為成熟的公式可以利用，而針對自然對流傳熱不可忽略的熱源環(huán)境，尚未見有計(jì)算非穩(wěn)態(tài)自然對流傳熱系數(shù)的相關(guān)成果。

為此，本文針對置于高溫艙室中的熱防護(hù)服材料，建立了服裝隔熱性能預(yù)報(bào)的數(shù)學(xué)模型，提出了一種計(jì)算非穩(wěn)態(tài)自然對流換熱系數(shù)的方法，利用有限差分法計(jì)算了模擬皮膚外側(cè)的溫度時(shí)程曲線，并與實(shí)驗(yàn)測量值進(jìn)行了對比，結(jié)果顯示文中自然對流換熱系數(shù)計(jì)算方法有效，可明顯改善防護(hù)服隔熱性能的數(shù)值預(yù)報(bào)精度，在壓縮防護(hù)服研發(fā)周期和實(shí)驗(yàn)成本方面意義明顯。

1 傳熱模型

1.1 模型假設(shè)

圖1示出了傳熱系統(tǒng)的組成，為便于后續(xù)推導(dǎo)，構(gòu)成熱防護(hù)服的3層材料由外到內(nèi)依次編號為1、2、3，空氣層編號為4，模擬皮膚層編號為5。

對于圖1所示傳熱系統(tǒng)，假設(shè)：1)熱量沿垂直于模擬皮膚方向傳遞，為一維傳熱問題；2)能量從外界環(huán)境向熱防護(hù)服傳播時(shí)，考慮對流和輻射；3)各材料的熱屬性穩(wěn)定，與溫度無關(guān)；4)材料4的厚度不大于6.4 mm，不計(jì)熱對流的影響[2]；5)恒溫內(nèi)核溫度為37 ℃。

1.2 傳熱系統(tǒng)微分方程

基于以上假設(shè)，對上述瞬態(tài)傳熱系統(tǒng)的數(shù)學(xué)描述如下。各材料內(nèi)部的熱量傳遞滿足：

(1)

式中：Ti(xi,t)為不同材料內(nèi)部溫度分布，分別與相應(yīng)標(biāo)號的材料對應(yīng)，K；xi為水平坐標(biāo)，采用局部坐標(biāo)系，即各材料的坐標(biāo)原點(diǎn)在該材料的左端壁上，x軸正向向右，m；t為時(shí)間，s；λi為不同材料的熱傳導(dǎo)率，W/(m·K)；ρi為不同材料的密度，kg/m3；ci為不同材料的比熱容，J/(kg·K)。

材料1外側(cè)(左側(cè))熱流邊界條件：

(2)

材料1與材料2及材料2與材料3的熱流邊界條件：

(3)

式中：δi為不同材料厚度，m。

材料3與材料4(空氣層)間的熱流邊界條件：

(4)

式中：ε=1/(1/εin+1/εs-1)，εin為材料3發(fā)射率，εs為模擬皮膚發(fā)射率。

材料4(空氣層)與材料5間的熱流邊界條件：

(5)

材料間的接觸邊界溫度邊界條件：

Ti|xi=δi=Ti+1|xi+1=0(i=1，…，4)

(6)

各材料的初始條件為：

Ti|t=0=C0(i=1，…，5)

(7)

T5|x5=δ5=C0

(8)

式中：C0為恒溫內(nèi)核溫度，K。

1.3 非穩(wěn)態(tài)自然對流換熱系數(shù)

實(shí)際上，外部自然對流換熱系數(shù)與初始和穩(wěn)態(tài)自然對流換熱系數(shù)關(guān)系十分密切，在自然對流狀態(tài)下，邊界層內(nèi)的空氣從初始時(shí)的靜止到穩(wěn)態(tài)時(shí)的運(yùn)動(dòng)，其速度變化較為平緩，相應(yīng)自然對流換熱系數(shù)也應(yīng)是平緩、連續(xù)的變化，因邊界層內(nèi)的空氣流動(dòng)是因空氣溫度降低而至空氣收縮和密度增加所致，而邊界層內(nèi)的空氣溫度取決于環(huán)境溫度和防護(hù)服的表面溫度，在環(huán)境溫度恒定的情況下，防護(hù)服的表面溫度是引起邊界層內(nèi)空氣流動(dòng)的直接因素，因此，有理由認(rèn)為非穩(wěn)態(tài)條件下的自然對流換熱系數(shù)是關(guān)于防護(hù)服表面溫度的連續(xù)函數(shù)，該連續(xù)函數(shù)有2個(gè)值已知，分別對應(yīng)于初始時(shí)和穩(wěn)態(tài)值。

為此，首先確定穩(wěn)態(tài)時(shí)的傳熱系統(tǒng)的自然對流換熱系數(shù)，其數(shù)學(xué)描述如下：

各材料內(nèi)部的熱量傳遞滿足：

(9)

材料1外側(cè)(右側(cè))熱流邊界條件：

hc(T1|x1=0)·(Tair-T1|x1=0)

(10)

式中：hc(T1|x=0)為穩(wěn)態(tài)時(shí)外部環(huán)境自然對流換熱系數(shù)，是T1|x=0的函數(shù)[7]，W/(m2·K)，其他參數(shù)如前所述。

材料1與材料2及材料2與材料3的熱流邊界條件見式(3)；材料3與材料4(空氣層)間的熱流邊界條件見式(4)；材料4(空氣層)與材料5間的熱流邊界條件見式(5)；材料間的接觸邊界溫度邊界條件見式(6)。通過解方程可求出穩(wěn)態(tài)時(shí)熱防護(hù)服表面溫度T1|x=0(令其為C1)以及穩(wěn)態(tài)時(shí)的自然對流換熱系數(shù)hc(C1)。

在初始條件下，上述傳熱系統(tǒng)被置于高溫艙室中，邊界層內(nèi)的空氣尚未開始運(yùn)動(dòng)，其熱傳遞形式主要是熱傳導(dǎo)，將空氣層離散，離散步長為Δxair，并設(shè)空氣的熱導(dǎo)率為λair，根據(jù)傅里葉定律有

(11)

整理后可得初始條件下的自然對流換熱系數(shù)為

(12)

(13)

將(13)代入式(3)，方程(1)～(8)可解。

2 傳熱系統(tǒng)數(shù)值算法

2.1 穩(wěn)態(tài)傳熱系統(tǒng)數(shù)值算法

構(gòu)造函數(shù)：f(T1,1)=q0(T1,1)-q′0(T1,1)，

T1,1∈(Tair,C0)

式中：Ti,j表示第i種材料第j個(gè)節(jié)點(diǎn)的溫度，將公式(9)～(10)、(3)～(6)寫成差分格式整理后可得上式中其他符號的含義，如下：

T4,N4=δ5q0/λ5+C0

穩(wěn)態(tài)時(shí)應(yīng)有f(T1,1)=0，可以證明，q′0為T1,1的單增函數(shù)，q0為T1,1的單減函數(shù)，采用迭代法可求解方程f(T1,1)=0，得到T1,1(即C1)，進(jìn)而得到hc(C1)，由式(17)可得h′c(T1|x=0)。

2.2 瞬態(tài)傳熱系統(tǒng)數(shù)值算法

將式(1)寫成差分格式：

(i=1，…，4)

(14)

式中：Ti,j,k表示材料i的節(jié)點(diǎn)j在tk時(shí)刻的溫度。根據(jù)式(7)給出的初始條件，4種材料在t2時(shí)刻的內(nèi)部節(jié)點(diǎn)溫度均可由式(13)確定，t2時(shí)刻材料邊界處的溫度則需單獨(dú)求解。

因T1,2,2為內(nèi)部節(jié)點(diǎn)，則其值為已知，將式(2)左端寫成向前差分格式，整理后可得：

(15)

式(15)是關(guān)于T1,1,2的4次方程，可采用迭代法求解。

將式(3)寫成差分格式，并結(jié)合式(6)，整理后可得

(16)

(17)

將式(4)、(5)寫成差分格式，整理后可得：

(18)

T4,N4,2=-P1T3,N3,2+P2

(19)

式中：Ck3=λ3/Δx3、Ck41=λ4/Δx4、Ck4=λ4/Δx4、Ck5=λ5/δ5，

P1=(Ck3+Ck41)/(Ck4+Ck5)

式(18)是關(guān)于T3,N3,2四次方程，可采用迭代法求解，由式(19)可確定T4,N3,2。

3 算例及分析

各材料的物理屬性[13]見表1，環(huán)境溫度Tair=75 ℃、恒溫內(nèi)核溫度C0=37 ℃，其他材料初始溫度為37 ℃。

表1 材料的物理屬性Tab.1 Thermophysical/geometrical properties of different material

利用上述方法對圖1所示傳熱系統(tǒng)進(jìn)行分析，服裝內(nèi)表面發(fā)射率εin取0.8，模擬皮膚發(fā)射率εs取0.97，服裝外表面發(fā)射率εout取0.1，λ5/δ5=5.433 W/(m·K)，標(biāo)準(zhǔn)大氣壓下空氣的物理屬性由美國國家標(biāo)準(zhǔn)局氣體熱性能表No. 564根據(jù)空氣溫度插值確定。

圖2示出防護(hù)服外表面在不同邊界條件和不同時(shí)間步長下模擬皮膚外側(cè)溫度的計(jì)算結(jié)果。

虛線是防護(hù)服外表為恒溫75 ℃時(shí)的情況，計(jì)算結(jié)果均高于實(shí)驗(yàn)值[13]，趨于穩(wěn)態(tài)時(shí)最大差值為3.92 ℃，但趨于穩(wěn)態(tài)的時(shí)間與實(shí)驗(yàn)值吻合較好，該邊界條件高估了防護(hù)服外表面的熱傳遞，總體上與實(shí)驗(yàn)結(jié)果偏差較大；點(diǎn)劃線是防護(hù)服外表面自然對流系數(shù)按Howard模型計(jì)算時(shí)的結(jié)果，計(jì)算結(jié)果均低于實(shí)驗(yàn)值，非穩(wěn)態(tài)階段最大差值為6.6 ℃，趨于穩(wěn)態(tài)時(shí)溫度值與實(shí)驗(yàn)結(jié)果吻合較好，但在時(shí)間上有較大的滯后，該邊界條件低估了防護(hù)服外表面的熱傳遞，總體上與實(shí)驗(yàn)結(jié)果偏差更大，盡管如此，利用該邊界條件計(jì)算穩(wěn)態(tài)溫度情況時(shí)卻可以獲得令人滿意的結(jié)果；實(shí)線是上述2種情況計(jì)算結(jié)果的加權(quán)平均，經(jīng)該方式處理后的數(shù)據(jù)在一定程度上得到了改善，但與實(shí)驗(yàn)結(jié)果仍有不小的差距，尤其是在穩(wěn)態(tài)階段。從上述分析可以看出，防護(hù)服外表面的傳熱邊界條件對計(jì)算結(jié)果有著非常重要的影響，是影響防護(hù)服隔熱性能模擬精度的關(guān)鍵因素，從對2種典型邊界條件的計(jì)算可以看出，模擬結(jié)果與實(shí)驗(yàn)結(jié)果差別較大，雖在一定程度上反映出了溫度隨時(shí)間的變化情況，但仍然不夠理想。

圖3示出時(shí)間步長為0.001 s時(shí)計(jì)算結(jié)果相對實(shí)驗(yàn)測量結(jié)果的誤差。可以看出，隨著時(shí)間步長的加密(空間步長相應(yīng)加密)計(jì)算結(jié)果迅速收斂，當(dāng)時(shí)間步長為0.001 s時(shí)，計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果高度一致，從圖3可以看出整體誤差小于0.1 ℃，計(jì)算結(jié)果令人滿意。即便是時(shí)間步長為0.1 s的情況，其計(jì)算結(jié)果也優(yōu)于圖2中利用Howard模型所得的結(jié)果。

圖4示出時(shí)間步長為0.001 s時(shí)熱流密度隨時(shí)間的變化。

從圖4可以看出，熱量主要來源于外界環(huán)境的對流傳熱，輻射傳熱僅占較小的比例，起始時(shí)的總熱流密度最大，接近1.7 kW/m2，幾乎全部來自對流傳熱，隨著時(shí)間的變化，自然對流的熱流密度迅速下降，輻射部分的熱流密度變化不大，相應(yīng)總熱流密度迅速下降，趨于穩(wěn)態(tài)時(shí)，自然對流的熱流密度達(dá)到最小值，約為46 W/m2，輻射熱流密度約為14 W/m2，總熱流密度達(dá)到最小值，約為60 W/m2。熱流密度的變化符合實(shí)際情況，也證實(shí)了高溫環(huán)境條件下，熱量主要來源于自然對流傳熱，在自然對流傳熱不可忽略的熱源環(huán)境下，自然對流傳熱系數(shù)計(jì)算的合理性將直接影響模擬結(jié)果的準(zhǔn)確性。

4 結(jié)束語

自然對流傳熱不可忽略甚至占優(yōu)的高溫環(huán)境是熱防護(hù)服一種較為典型的應(yīng)用場景，文中探討了此時(shí)瞬態(tài)自然對流換熱系數(shù)的一種計(jì)算方法，依此方法預(yù)報(bào)了防護(hù)服隔熱性能并與實(shí)驗(yàn)測量值進(jìn)行了對比，驗(yàn)證了該自然對流換熱系數(shù)計(jì)算方法，解決了自然對流傳熱不可忽略甚至占優(yōu)時(shí)，防護(hù)服非穩(wěn)態(tài)隔熱性能預(yù)報(bào)中的一個(gè)問題。對于作業(yè)于非穩(wěn)態(tài)傳熱階段的防護(hù)服，在進(jìn)行防護(hù)服隔熱性能數(shù)值模擬時(shí)，若高估自然對流的影響，則溫度分布會提前進(jìn)入穩(wěn)態(tài)，且穩(wěn)態(tài)溫度高于實(shí)際溫度，這會導(dǎo)致防護(hù)服設(shè)計(jì)的過于保守；若低估自然對流的影響，則溫度分布會延后進(jìn)入穩(wěn)態(tài)，此時(shí)的設(shè)計(jì)將不足以隔熱，嚴(yán)重時(shí)會導(dǎo)致人體受到損傷。

此外，文中假設(shè)材料的物理屬性不隨溫度變化，但當(dāng)環(huán)境溫度很高時(shí)，外側(cè)材料的物理屬性首先受到溫度的影響，下一步應(yīng)將材料物理屬性隨溫度變化時(shí)的情況納入研究，以適應(yīng)更惡劣的環(huán)境溫度。

致謝感謝天津海太科技有限公司馬曉梅總工在本文修改過程中提出的寶貴意見。