基于蒙特卡羅法的輻射熱防護(hù)測(cè)試儀的數(shù)值分析

萬(wàn)賢福， 汪 軍， 黃 偉

(1. 浙江省新型紡織品研發(fā)重點(diǎn)實(shí)驗(yàn)室， 浙江 杭州 310001; 2. 東華大學(xué) 紡織學(xué)院， 上海 201620)

基于蒙特卡羅法的輻射熱防護(hù)測(cè)試儀的數(shù)值分析

萬(wàn)賢福1, 2， 汪 軍2， 黃 偉2

(1. 浙江省新型紡織品研發(fā)重點(diǎn)實(shí)驗(yàn)室， 浙江 杭州 310001; 2. 東華大學(xué) 紡織學(xué)院， 上海 201620)

為了探討輻射熱防護(hù)測(cè)試(RPP)儀在試樣表面產(chǎn)生的輻射熱通量分布，采用蒙特卡羅法對(duì)RPP建立了數(shù)值模型，并進(jìn)行了實(shí)驗(yàn)驗(yàn)證. 計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果十分吻合，證明了模型的合理性與正確性. 計(jì)算結(jié)果顯示，在試樣的入射表面，輻射熱通量呈橢圓形分布，能量集中在內(nèi)部區(qū)域，在中心區(qū)域的輻射熱通量分布相對(duì)均勻，兩側(cè)邊緣的輻射熱通量相對(duì)較低. 研究結(jié)果為進(jìn)一步研究RPP測(cè)試試樣的防護(hù)效果提供了理論依據(jù).

輻射熱傳遞； 輻射熱防護(hù)； 輻射熱通量； 蒙特卡羅法

輻射熱防護(hù)性能廣泛應(yīng)用于評(píng)價(jià)火災(zāi)防護(hù)服的防護(hù)熱輻射傷害的能力[1]. 為了獲取這一指標(biāo)，人們?cè)O(shè)計(jì)出輻射熱防護(hù)測(cè)試(RPP)儀進(jìn)行檢測(cè)[2]. RPP儀產(chǎn)生均勻分布的輻射熱通量，并將試樣暴露在此輻射熱通量之下，檢測(cè)穿過(guò)試樣的輻射熱通量，然后基于Stoll曲線(xiàn)預(yù)測(cè)皮膚發(fā)生二度燒傷的時(shí)間[3-4]. 由于Stoll曲線(xiàn)是基于皮膚直接暴露于均勻分布熱通量的情形下的實(shí)驗(yàn)數(shù)據(jù), 因此，RPP儀準(zhǔn)確預(yù)測(cè)熱輻射皮膚燒傷的關(guān)鍵在于儀器產(chǎn)生均勻分布的輻射熱通量. 盡管文獻(xiàn)[2]給出了RPP儀的設(shè)計(jì)規(guī)范，以保證測(cè)試結(jié)果的準(zhǔn)確與一致. 然而在實(shí)際的測(cè)試中，通過(guò)觀察試樣燒傷程度的分布可以發(fā)現(xiàn)，輻射熱通量并非均勻一致. 由于入射輻射熱通量的分布情況是測(cè)試面料的熱防護(hù)效果的一個(gè)重要決定因素，因此，本文嘗試對(duì)RPP儀產(chǎn)生的輻射熱通量進(jìn)行數(shù)值分析，旨在探討試樣表面的入射輻射熱通量的分布情況，為后續(xù)分析RPP儀測(cè)試面料的熱防護(hù)性能做準(zhǔn)備. 由于蒙特卡羅法在大量復(fù)雜問(wèn)題中所表現(xiàn)出的優(yōu)異性[5-7]，本文采用蒙特卡羅法進(jìn)行建模.

1 模 型

1.1 模型描述

RPP儀的結(jié)構(gòu)如圖1所示，基本由石英燈陣列、燈箱、隔熱板、擋板、試樣夾板以及傳感器等組成. 石英燈管為熱輻射源，發(fā)出熱輻射；石英燈管置于燈箱之內(nèi)，燈箱的一側(cè)裝有隔熱板，用于阻隔熱量，以防止試樣一側(cè)過(guò)熱；擋板用于遮擋熱輻射，當(dāng)測(cè)試開(kāi)啟，擋板迅速被提升，從而將試樣直接暴露在石英燈管產(chǎn)生的熱輻射之下；試樣被夾在試樣夾板1和2之間，傳感器則放置在試樣夾板2的后側(cè).隔熱板和試樣夾板都開(kāi)有一定大小的窗口，測(cè)試時(shí)，熱輻射穿過(guò)這些窗口最終抵達(dá)熱通量傳感器表面.

圖1 RPP儀的結(jié)構(gòu)示意圖Fig.1 Schematic for the RPP apparatus

由于RPP儀僅檢測(cè)熱輻射的影響，而由熱對(duì)流產(chǎn)生的影響則消除，在燈箱壁面設(shè)置有垂直方向的通風(fēng)對(duì)流口，通過(guò)燈箱內(nèi)垂直方向的機(jī)械對(duì)流以消除熱對(duì)流作用于燈箱外的試樣表面. 石英燈管壁面與燈箱內(nèi)壁對(duì)流的存在，使得對(duì)流與輻射耦合，問(wèn)題變得復(fù)雜，為了簡(jiǎn)化，燈管溫度與燈箱內(nèi)壁的溫度通過(guò)測(cè)試手段獲得，成為已知條件，將問(wèn)題轉(zhuǎn)變成純輻射問(wèn)題. 除此之外，本模型為了進(jìn)一步簡(jiǎn)化問(wèn)題還做出如下假設(shè)：

(1) 壁面均為漫射表面；

(2) 空氣與石英燈管內(nèi)氣體不參與熱輻射，它們對(duì)熱輻射透明；

(3) 石英對(duì)于不同波段的熱輻射或者完全透明，或者完全不透明，石英介質(zhì)不參與熱輻射的吸收與散射，熱輻射的發(fā)射、吸收與反射均發(fā)生在氣、固或固、氣界面上；

(4) 壁面的溫度分布均勻；

(5) 鎢絲表面的發(fā)射為漫反射；

(6) 鎢絲和石英燈管兩端的影響可以忽略不計(jì)；

(7) 壁面上的垂直通風(fēng)孔泄漏的熱輻射可以忽略不計(jì).

基于上述假設(shè)，本文采用蒙特卡羅法建立了RPP儀的熱輻射模型. 雖然RPP儀不僅僅在內(nèi)部產(chǎn)生熱輻射用于檢測(cè)，其外殼也對(duì)外部環(huán)境釋放熱輻射，然而外部熱輻射并非用于測(cè)試，這不是本文關(guān)注的問(wèn)題，所以本文模型并不給予考慮. 此外，本文僅關(guān)注由RPP儀產(chǎn)生的用于測(cè)試的熱輻射，并不關(guān)注熱輻射在試樣內(nèi)部的傳播. 因此，本文模型所模擬的問(wèn)題如下：當(dāng)測(cè)試開(kāi)始時(shí)，擋板去除，熱輻射首先由石英燈發(fā)出，然后在4個(gè)相連的腔體(燈箱、隔熱板窗口、擋板室和試樣夾板窗口)內(nèi)傳播，最后到達(dá)試樣表面. 由于腔體壁面吸收熱輻射后，溫度升高，成為了第二輻射源.

蒙特卡羅法把熱輻射看成大量能量粒子的隨機(jī)傳播行為，其模擬的過(guò)程主要包括：各種表面和形體的幾何描述，發(fā)射、反射、折射和散射的概率模型，能量粒子軌跡與表面、形體的相交判斷與求解等. 標(biāo)準(zhǔn)蒙特卡羅法的這些細(xì)節(jié)工作在文獻(xiàn) [5-7]中都已描述，在此不再贅述. 由于石英燈的模擬在本模型中比較復(fù)雜，因此，本文主要對(duì)其進(jìn)行詳細(xì)闡述.

1.2 石英燈模型

圖2 紅外石英燈模型示意圖Fig.2 Schematic for the model of quartz lamp

石英燈主要包括兩個(gè)部分：燈絲與燈管.由于忽略?xún)啥说挠绊?，在形體上燈管用直管描述，而燈絲簡(jiǎn)化為一個(gè)圓柱體，如圖2所示. 現(xiàn)實(shí)中，燈絲是呈螺旋狀旋繞的細(xì)鎢絲，比較復(fù)雜，但相對(duì)于燈箱內(nèi)部腔體，燈絲形體纖細(xì)，其輻射行為可以與一定直徑的圓柱體等效.

石英的輻射屬性與輻射的波長(zhǎng)密切相關(guān)，通常對(duì)于0.2～4.0 μm波段的輻射，石英可以看成是透明的，而這個(gè)波段之外的熱輻射則幾乎不透明[8]. 因此，本文模型中的能量粒子按照一定的分布賦予了波長(zhǎng)屬性，根據(jù)能量粒子的波長(zhǎng)，決定了能量粒子在入射石英燈管之后的行為. 如果能量粒子處于不透明的波段，那么它將被吸收或反射；而能量粒子處于透明波段，則可能發(fā)生折射. 這取決于入射角、波長(zhǎng)以及石英的光學(xué)屬性. 對(duì)于石英燈管表面的反射率ρ(θi)，其具體的數(shù)學(xué)描述[9]為：

(1)

其中：

2a2=[(η2-κ2-sin2θi)2+4η2κ2]1/2+

(η2-κ2-sin2θi)

2b2=[(η2-κ2-sin2θi)2+4η2κ2]1/2-

(η2-κ2-sin2θi)

θi為入射角. 對(duì)于氣體-石英界面:

η=ηq， κ=κq

而對(duì)于石英-氣體界面:

其中：ηq和κq分別為石英的復(fù)數(shù)折射率的實(shí)部和虛部.ηq和κq取決于能量粒子的波長(zhǎng)，其數(shù)據(jù)可以查表獲得.

對(duì)于石英燈管表面的折射軌跡，其數(shù)學(xué)描述[10]為

(2)

其中：φ為折射角；θi為入射角.

1.3 統(tǒng)計(jì)與計(jì)算

本文模型中，須跟蹤發(fā)射的能量粒子，直到它被吸收或者到達(dá)試樣表面. 由于本文需分析試樣表面入射的輻射熱通量分布，因此對(duì)試樣表面區(qū)域進(jìn)行網(wǎng)格劃分，統(tǒng)計(jì)網(wǎng)格內(nèi)接收到的能量粒子數(shù)n，可以按式(3)計(jì)算輻射熱通量(q).

(3)

其中：Qt為RPP儀內(nèi)部所有熱輻射源發(fā)出的能量總和；N為熱源發(fā)出的能量粒子數(shù)總和；A為試樣表面的某個(gè)網(wǎng)格面積.

2 驗(yàn)證實(shí)驗(yàn)

驗(yàn)證實(shí)驗(yàn)采用一臺(tái)美國(guó)Richmond 公司生產(chǎn)的RPP儀，這臺(tái)儀器基本根據(jù)ASTM F2702標(biāo)準(zhǔn)[11]設(shè)計(jì)，但少數(shù)幾個(gè)地方略有不同：隔熱板窗口的上邊緣比標(biāo)準(zhǔn)略高8mm，燈管的發(fā)光長(zhǎng)度約為178mm，而不是標(biāo)準(zhǔn)中的127mm，燈管間距約1mm，而不是標(biāo)準(zhǔn)中的0.4mm. 由于本文探討的對(duì)象是RPP儀產(chǎn)生的輻射熱通量(即試樣表面的入射輻射熱通量)，而不是通過(guò)試樣的熱通量. 因此，在進(jìn)行驗(yàn)證實(shí)驗(yàn)時(shí)，直接暴露熱通量傳感器以檢測(cè)所在位置的輻射熱通量，而不加載任何試樣. 通過(guò)移動(dòng)熱通量傳感器的位置，獲得不同位置的輻射熱通量分布. 熱通量傳感器直接采用RPP儀自帶的銅片熱通量傳感器，其探測(cè)區(qū)域是直徑為4cm的圓盤(pán)外形，它所測(cè)得的熱通量值實(shí)際上是圓盤(pán)區(qū)域的平均值. 石英燈管壁和燈箱內(nèi)壁的溫度采用紅外熱像儀(FLIR T440, FLIR System Inc, 美國(guó))進(jìn)行測(cè)量.

數(shù)值模型采用面向?qū)ο蟮腣isual C++進(jìn)行編程實(shí)現(xiàn)，運(yùn)行在一臺(tái)工作站 (Dell Precision T5610, 美國(guó))上，配置有兩個(gè)6核CPU(Intel?Xeon?CPU E5—2620@2.10 GHz).采用20個(gè)線(xiàn)程進(jìn)行并行計(jì)算以提升計(jì)算速度. 驗(yàn)證模型采用的幾何參數(shù)根據(jù)儀器的實(shí)際尺寸進(jìn)行設(shè)置，也就是除了少數(shù)尺寸與標(biāo)準(zhǔn)不同外，其他都可以從ASTM F2702中獲得. 輻射屬性參數(shù)如表1所示.模型中呈螺旋狀的鎢絲采用圓柱體進(jìn)行表示，其直徑需要確定，在此直徑下，其輻射效果與螺旋狀鎢絲相當(dāng). 當(dāng)直徑配置為1.2mm時(shí)，中心區(qū)域的輻射熱通量給出合理的結(jié)果，因此，這個(gè)圓柱體的直徑被設(shè)置為1.2mm.

驗(yàn)證模型中，最終接收能量粒子的物體不是試樣，而是熱通量傳感器的探測(cè)表面，在統(tǒng)計(jì)時(shí)仍然采用式(3)，但這時(shí)A為熱通量傳感器的探測(cè)區(qū)域的面積，而不是網(wǎng)格面積.

表1 模型中的輻射屬性參數(shù)

3 計(jì)算結(jié)果與分析

為了便于闡述，以試樣夾板1外表面的中心點(diǎn)為原點(diǎn)，建立如圖3所示的直角坐標(biāo)系.

圖3 RPP儀的直角坐標(biāo)系統(tǒng)示意圖Fig.3 Cartesian coordinate system in the RPP apparatus

通過(guò)建立的驗(yàn)證模型，在中心區(qū)域輻射熱通量為21 kW/m2的情形下計(jì)算了在x軸與y軸上的輻射熱通量分布，此時(shí)燈管、燈箱內(nèi)壁、擋板室以及試樣夾板1的溫度測(cè)量結(jié)果分別為894， 556， 370和298K. 燈絲的溫度則依據(jù)文獻(xiàn) [1]的實(shí)驗(yàn)結(jié)果設(shè)為1625K. 由于以傳感器探測(cè)區(qū)域?yàn)槟芰苛Ｗ咏邮彰?，統(tǒng)計(jì)區(qū)域相對(duì)較大，發(fā)射106個(gè)能量粒子，輻射熱通量計(jì)算結(jié)果已經(jīng)穩(wěn)定. 在x軸與y軸上的輻射熱通量分布測(cè)試及計(jì)算結(jié)果如圖4所示.由圖4(a)表明，在x= [-20mm, 20mm] 區(qū)間的輻射熱通量分布較為均勻，而且處于相對(duì)高的熱通量范圍，而在這個(gè)區(qū)間之外，距離原點(diǎn)越遠(yuǎn)，輻射熱通量衰減越快.由圖4(b)表明，在y= [-10mm, 10mm] 區(qū)間的輻射熱通量分布較均勻.由圖4可以看出，測(cè)試結(jié)果與理論計(jì)算結(jié)果比較吻合，證明了模型的合理性與正確性.

(a)

(b)圖4 在x軸和y軸上的輻射熱通量分布Fig.4 Distribution of the heat flux on x and y axis

在模型得以驗(yàn)證之后，為了更加精細(xì)地分析試樣表面入射的輻射分布，對(duì)試樣接收面(與試樣夾板1窗口一樣大小，尺寸為152 mm × 64 mm)劃分152×64個(gè)網(wǎng)格. 統(tǒng)計(jì)每個(gè)網(wǎng)格上接收到的能量粒子，計(jì)算每個(gè)網(wǎng)格上的輻射熱通量值，從而獲得輻射熱通量分布. 由于本文的研究對(duì)象是標(biāo)準(zhǔn)的RPP儀，因此，此時(shí)建立模型的幾何參數(shù)完全按照ASTM F2702進(jìn)行設(shè)置，而不再按照驗(yàn)證儀器的實(shí)際尺寸. 此時(shí)計(jì)算的情形同樣是中心區(qū)域輻射熱通量為21 kW/m2. 通過(guò)發(fā)射1010個(gè)能量粒子，最終得到的結(jié)果如圖5所示. 圖5中的每個(gè)點(diǎn)顏色對(duì)應(yīng)相應(yīng)坐標(biāo)下的熱通量值. 圖5中輻射熱通量值在空間上的波動(dòng)源于蒙特卡羅法的統(tǒng)計(jì)學(xué)特性，但事實(shí)上這些波動(dòng)很小(CV值<0.6%)，并不影響所要求的精度.圖5清晰地呈現(xiàn)了輻射熱通量以坐標(biāo)原點(diǎn)為中心的橢圓形分布，離中心越遠(yuǎn)，橢圓上的輻射熱通量值越小. 除邊緣部分的輻射熱通量相對(duì)較低之外，內(nèi)部的熱通量分布相對(duì)均勻，這與實(shí)驗(yàn)中觀察到試樣燒傷部位集中在中心區(qū)域的現(xiàn)象一致. 注意到在圖5中熱通量值的分布梯度比圖4中大，原因是前者為非常小的網(wǎng)格(1 mm × 1 mm)上的熱通量平均值，而后者為直徑4cm的圓形區(qū)域上的平均值，因此，后者要平緩得多.

圖5 xOy平面上的輻射熱通量分布

4 結(jié) 語(yǔ)

本文采用蒙特卡羅法對(duì)輻射熱防護(hù)測(cè)試(RPP)儀進(jìn)行數(shù)值分析，計(jì)算結(jié)果表明，在試樣的入射表面，輻射熱通量分布并不均勻，高熱通量集中分布在內(nèi)部的橢圓形區(qū)域. 盡管RPP儀測(cè)試中傳感器所檢測(cè)的中心區(qū)域部分輻射熱通量分布相對(duì)均勻，本文所建立的模型對(duì)于RPP儀的改進(jìn)可以提供重要的理論支持.

[1] KELTNER N R. Evaluating thermal protective performance testing [J]. Journal of the ASTM International , 2005(2):191-204.

[2] NFPA 1977—2005 Protective clothing and equipment for wildland fire fighting [S]. Quincy, MA: National Fire Protection Association, 2005.

[3] STOLL A M, GREENE L C. Relationship between pain and tissue damage due to thermal radiation [J]. Journal of Applied Physiology, 1959,14(3): 373-382.

[4] STOLL A M, CHIANTA M A. Method and rating system for evaluation of thermal protection [J]. Aerospace Medicine, 1969,40(11):1232-1238.

[5] HOWELL J R. Thermal radiation heat transfer [M]. 4th ed. New York: Taylor and Francis, 2010.

[6] MODEST M F. Radiative heat transfer [M]. 2nd ed. Oxford: Academic Press, 2003.

[7] 劉偉.輻射介質(zhì)傳熱[M].北京：中國(guó)電力出版社,2009.

[8] WENTINK T, PIANET W G. Infrared emission spectra of quartz [J]. Journal of the Optical Society of America, 1961,51(6):595-600.

[9] TOULOUKIAN Y S, DEWITT D P. Thermal radiative properties-nonmetallic solids [M]. New York: Plenum Press, 1972.

[10] OZISIK M N. Radiative transfer and interactions with conduction and convection [M]. Hoboken, NJ: John Wiley & Sons, Inc, 1980.

[11] ASTM F2702—2015 Standard test method for radiant heat performance of flame resistant clothing materials with burn injury prediction [S].

[12] TURNER T L, ASH R L. Numerical and experimental analyses by of the radiant heat produce produced by quartz heating systems [R]. NASA Technical Paper, 1994:3387.

Numerical Analysis of a Radiant Protective Performance Apparatus Based on the Monte Carlo Method

WANXian-fu1, 2，WANGJun2，HUANGWei2

(1. Zhejiang Provincial Key Laboratory of Novel Textiles Research and Development, Hangzhou 310001, China;2. College of Textiles, Donghua University, Shanghai 201620, China)

In order to analyze the radiant heat flux distribution at the incident surface of test specimen in radiant protective performance(RPP) apparatus, a numerical model based on the Monte Carlo method was developed. And experiments were conducted to validate the model. It is found that the predictions of the model show a good agreement with the experimental results. Further, at the incident surface of the sample, the calculation results show an oval-shaped distributed heat flux area with a high energy concentrated in the interior region. The radiant heat flux distribution is relatively uniform at the central region, and relatively low on both sides of the edge. The results will provide a theoretical guideline for further research on the protective effect of the specimen in RPP testing.

radiative heat transfer；radiant protective performance；radiant heat flux；Monte Carlo method

1671-0444(2016)01-0063-04

2014-11-14

國(guó)家自然科學(xué)基金資助項(xiàng)目(51106021)；浙江省公益基金資助項(xiàng)目(2013C31006)

萬(wàn)賢福(1979—)，男，江西贛州人，副教授，博士，研究方向?yàn)榉b熱濕舒適性.E-mail: wanxianfu@hotmail.com

TK 121 文獻(xiàn)標(biāo)志碼： A