鈾礦通風(fēng)尾氣中氣態(tài)放射性核素氡大氣擴(kuò)散數(shù)值模擬

謝東 ，王漢青, ，劉澤華，葉勇軍，熊軍

(1.南華大學(xué) 城市建設(shè)學(xué)院，湖南 衡陽(yáng)，421001；2.中南大學(xué) 能源科學(xué)與工程學(xué)院，湖南 長(zhǎng)沙，410083；3.湖南工業(yè)大學(xué) 土木工程學(xué)院，湖南 株洲，412007)

鈾礦通風(fēng)尾氣中含有大量的氡及氡子體、鈾礦塵和放射性氣溶膠等放射性有害物質(zhì)，這些放射性有害物在隨風(fēng)輸送和彌散的過(guò)程中，不僅對(duì)公眾的集體劑量貢獻(xiàn)大，而且對(duì)地表土壤、植物也造成一定污染[1]。放射性核素在大氣中的遷移擴(kuò)散是一個(gè)非常復(fù)雜的過(guò)程，核素一方面隨大氣運(yùn)動(dòng)而遷移、擴(kuò)散和沉降，一方面還伴隨有自身衰變的物理過(guò)程，其污染程度與污染物濃度、排放口位置、排放高度、排風(fēng)風(fēng)速、當(dāng)?shù)卮髿饴?lián)合頻率等因素相關(guān)。有關(guān)數(shù)據(jù)表明：鈾礦井通風(fēng)排出的放射性廢氣是鈾礦冶系統(tǒng)對(duì)環(huán)境影響最大的污染源項(xiàng)[2]，其中222Rn是鈾礦排風(fēng)口流出物中的主要核素，其次為鈾礦塵238U和226Ra。放射性核素在大氣中遷移和擴(kuò)散的研究在國(guó)際上已進(jìn)行多年，其中大多數(shù)是針對(duì)核電站或核裝置泄漏對(duì)環(huán)境及人類的影響[3-5]和核事故后果評(píng)價(jià)[6-7]的研究，對(duì)鈾礦排風(fēng)井尾氣中放射性核素對(duì)人和環(huán)境的影響研究較少，且前人針對(duì)鈾礦通風(fēng)尾氣污染的研究工作主要以現(xiàn)場(chǎng)實(shí)測(cè)為主[8-10]，采用數(shù)值模擬研究方法和研究結(jié)論很少[11]?；谝陨戏治?，本文作者建立核素氡在大氣中遷移擴(kuò)散的物理模型和計(jì)算模型，對(duì)放射性核素氡遷移擴(kuò)散規(guī)律開(kāi)展數(shù)值模擬研究。

1 數(shù)值模擬

1.1 物理模型

本文研究的鈾礦山物理模型長(zhǎng)×寬×高為 800 m×700 m×190 m，模型建立方法參考文獻(xiàn)[12]。排風(fēng)井長(zhǎng)×寬×高為2.7 m×2.7 m×2.0 m。物理模型如圖1所示。

圖1 鈾礦井通風(fēng)物理模型示意圖Fig.1 Physical model scheme of uranium mine shaft ventilation

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

假定：(1) 計(jì)算區(qū)域的上端自由面沒(méi)有質(zhì)量和溫度的交換；(2) 環(huán)境氣流為連續(xù)、不可壓縮的流體；(3)鈾礦通風(fēng)尾氣排風(fēng)口與環(huán)境氣流沒(méi)有熱量的交換，三維控制方程分別為求解質(zhì)量輸運(yùn)的連續(xù)性方程，動(dòng)量的Navier-Stokes方程及湍流k-ε兩方程模型。考慮上面物理模型簡(jiǎn)化假設(shè)，針對(duì)本文所建立的數(shù)學(xué)模型，其基本形式與各控制方程具有相同的形式，可以使用如下通用方程表示[13]：

式中：ρ為密度；U為速度；φ為通用變量，代表U，V，W，k，ε和1，當(dāng)φ=1時(shí)即為連續(xù)性方程；Γφ為通用變量φ的有效擴(kuò)散系數(shù)；Sφ為源項(xiàng)。k-ε兩方程模型中的幾個(gè)系數(shù)取值如下

考慮氡在大氣中的遷移擴(kuò)散涉及到的物理化學(xué)因素包括平均風(fēng)對(duì)流，湍流擴(kuò)散，放射性核素衰變，干、濕沉降，重力對(duì)核素的作用等，氡氣在大氣中的遷移擴(kuò)散過(guò)程可以描述為：

式中：C為空氣中氡的平均濃度(Bq/m3)；D為氡氣擴(kuò)散系數(shù)(m2/s)；U為氡氣運(yùn)動(dòng)速度(m/s)；wz為重力沉降率(m/s)；Λ為雨水沖刷系數(shù)，取值為10-4s-1[15]；λ為放射性核素的衰變常量，取值為2.1×10-6s-1；Z為重力方向的位移；Q為源項(xiàng)。

1.3 邊界條件

結(jié)合鈾礦區(qū)當(dāng)?shù)氐湫偷臍庀?、地理和地形條件，文中選取D類(中性)大氣穩(wěn)定度，確定4種環(huán)境風(fēng)速和 2種下墊面條件進(jìn)行數(shù)值計(jì)算。排風(fēng)口速度為 3.0 m/s，排風(fēng)口處的氡氣濃度為8 500 Bq/m3。

用具有代表性的地面粗糙類別對(duì)應(yīng)的大氣邊界條件為來(lái)流條件，大氣邊界層風(fēng)速剖面指數(shù)分布為

式中：U2為在H2高度的平均風(fēng)速；U1取當(dāng)?shù)仉x地10 m高度處10 min平均風(fēng)速觀測(cè)值；冪指數(shù)a為地面粗糙度指數(shù)，地形越粗糙，地表對(duì)氣流的阻滯作用越強(qiáng)，a越大。文中粗糙度為0.1 m和1.0 m對(duì)應(yīng)的a分別為0.19和0.25。

本文綜合考慮了山體和森林2種下墊面的影響，依據(jù)GB 50009—2001(《建筑結(jié)構(gòu)荷載規(guī)范》)[16]和實(shí)際的地理地形條件，選擇了0.1 m和1.0 m 2種下墊面粗糙度狀況。在使用k-ε模型對(duì)湍流計(jì)算時(shí)，需要給定進(jìn)口邊界上k和ε的估算值。文中入口處來(lái)流的湍流特性借助以下經(jīng)驗(yàn)公式[17]進(jìn)行計(jì)算：

式中：為進(jìn)口處的平均速度；L為湍流積分尺度；I為湍流通量，%。

這里地面和鈾礦排風(fēng)井壁使用壁面邊界條件，采用非平衡壁面函數(shù)以模擬壁面附近的復(fù)雜流動(dòng)現(xiàn)象，出口邊界為充分發(fā)展假定。

1.4 數(shù)值方法

基于有限體積法，控制方程的擴(kuò)散項(xiàng)和對(duì)流項(xiàng)均采用中心差分格式進(jìn)行離散，源項(xiàng)采用迎風(fēng)差分法，壓力和速度的耦合方式采用SIMPLE法。網(wǎng)格采用非結(jié)構(gòu)化網(wǎng)格，經(jīng)過(guò)多次嘗試劃分和試計(jì)算，最終確定網(wǎng)格間距為4 m，排風(fēng)口周圍網(wǎng)格間距為2 m，網(wǎng)格劃分如圖2所示。求解步驟如下：通過(guò)求解連續(xù)性方程、動(dòng)量方程和k-ε湍流方程得到大氣風(fēng)場(chǎng)結(jié)構(gòu)，在此基礎(chǔ)之上，通過(guò) C++語(yǔ)言編譯 UDF程序?qū)腚睔鈹U(kuò)散方程，求解得到核素氡的濃度分布。這里僅介紹風(fēng)向?yàn)槟巷L(fēng)時(shí)(X=0平面)仿真結(jié)果。

圖2 物理模型網(wǎng)格劃分圖Fig.2 Grid diagram of computation model

2 數(shù)值結(jié)果及分析

2.1 大氣風(fēng)場(chǎng)

圖3和圖4所示分別為X=0斷面上大氣風(fēng)場(chǎng)速度等值線圖。由圖3和4可以看出：z=1.0 m地形的高空平均風(fēng)速要大于z=0.1 m地形的高空平均風(fēng)速；隨著水平風(fēng)速的增大，地面及高空風(fēng)速都隨著增大；山體的形狀對(duì)速度分布有較大影響，下風(fēng)向200～300 m有明顯漩渦存在；隨著速度的增大，該漩渦有向下風(fēng)向移動(dòng)的現(xiàn)象。當(dāng)風(fēng)速大于2.0 m/s時(shí)，該位置漩渦不明顯。由于環(huán)境風(fēng)速入口處較小，導(dǎo)致入口速度梯度較大。不同的粗糙度對(duì)山谷處的風(fēng)速影響比較一致，近地面處風(fēng)速大小接近，但粗糙度增大，下風(fēng)向200～300 m產(chǎn)生的渦旋變大，水平風(fēng)速越大，漩渦風(fēng)速越大。

(1) 工況一：z=0.1 m；U=0.5，1.0，2.0，4.0(m/s)；X=0平面。

(2) 工況二：z=1.0 m；U=0.5，1.0，2.0，4.0(m/s)；X=0平面。

圖3 粗糙度z=0.1 m時(shí)的大氣風(fēng)場(chǎng)速度等值線圖(m/s)Fig.3 Atmospheric wind velocity contours in different wind speeds at surface roughness of 0.1 m(m/s)

2.2 氡氣濃度分布規(guī)律

圖5所示為在z=0.1 m時(shí)X=0斷面上氡氣濃度等值線分布圖。由圖5可以看出：在水平風(fēng)速作用下，氡氣濃度沿下風(fēng)向呈橢圓形分布，隨著距離的增大，氡氣濃度減?。划?dāng)水平風(fēng)速為0.5 m/s時(shí)，局部氡氣濃度高，局部污染嚴(yán)重，受入口氣流湍流及下游漩渦的影響，上風(fēng)向150 m范圍和下風(fēng)向300 m范圍內(nèi)氡氣濃度均大大超過(guò)了對(duì)照點(diǎn)參考值(12 Bq/m3)；風(fēng)速增大，水平氣流對(duì)尾氣射流的作用力增強(qiáng)，擴(kuò)散高度降低，局部污染減弱；氡氣濃度隨著環(huán)境風(fēng)速的增大逐步減小，當(dāng)風(fēng)速增大至2.0 m/s時(shí)，下風(fēng)向180 m基本達(dá)到了對(duì)照值的污染標(biāo)準(zhǔn)。

圖4 粗糙度z=1.0 m時(shí)大氣風(fēng)場(chǎng)速度等值線圖(m/s)Fig.4 Atmospheric wind velocity contours in different wind speeds at surface roughness of 1.0 m (m/s)

圖5 粗糙度z=0.1 m時(shí)氡氣濃度等值線圖(Bq/m3)Fig.5 Radon concentration contours in different atmospheric wind speeds at surface roughness of 0.1 m (Bq/m3)

z=1.0 m時(shí)X=0斷面上氡氣濃度等值線分布圖如圖6所示，粗糙度對(duì)下風(fēng)向氡氣污染影響明顯，相比z=0.1 m的地形，z=1.0 m的地形氡氣濃度沿下風(fēng)向要小很多，這主要是由于排風(fēng)口上游風(fēng)速切變效應(yīng)和下風(fēng)坡地形的影響，核素氡無(wú)法長(zhǎng)距離擴(kuò)散；當(dāng)風(fēng)速為0.5 m/s時(shí)，由于受入口氣流湍流及下游旋渦的影響明顯，沿上風(fēng)向150 m和下風(fēng)向150 m范圍內(nèi)氡氣濃度均超過(guò)對(duì)照值(12 Bq/m3)；隨著風(fēng)速的增大，水平擴(kuò)散作用增強(qiáng)，氡氣濃度大大減小，當(dāng)風(fēng)速為2.0 m/s時(shí)，氡氣污染距離較小，主要污染區(qū)域出現(xiàn)在下風(fēng)向 100 m范圍內(nèi)。

(1) 工況一：z=0.1 m；U=0.5，1.0，2.0，4.0(m/s)；X=0平面。

(2) 工況二：z=1.0 m；U=0.5，1.0，2.0，4.0(m/s)；X=0平面。

2.3 不同粗糙度氡氣濃度分布比較

2種粗糙度4種風(fēng)速條件下氡氣濃度分布如圖7和圖8所示。在z=0.1 m時(shí)，X=0平面上，隨著風(fēng)速的增大，氡氣濃度衰減很快；風(fēng)速越大，氡氣濃度越小，下風(fēng)向100 m內(nèi)污染嚴(yán)重，此時(shí)4種風(fēng)速工況下氡氣濃度均超過(guò)了對(duì)照點(diǎn)參考值(12 Bq/m3)；隨著距離的增大，氡濃度迅速衰減，在400 m處均達(dá)到了對(duì)照點(diǎn)參考值；當(dāng)粗糙度由z=0.1 m增大到z=1.0 m時(shí)，排風(fēng)口附近氡濃度增大，這與大粗糙度地形的摩擦阻力有關(guān)，但隨著風(fēng)速的增大，受排風(fēng)口上游風(fēng)切變和下風(fēng)坡地形的影響，核素氡擴(kuò)散能力較弱，污染主要集中在400 m范圍內(nèi)。

圖6 粗糙度z=1.0 m時(shí)氡氣濃度等值線圖(Bq/m3)Fig.6 Radon concentration contours in different atmospheric wind speeds at surface roughness of 1.0 m (Bq/m3)

圖7 粗糙度z=0.1 m時(shí)氡氣濃度對(duì)比圖Fig.7 Radon concentration contrast graph at surface roughness of 0.1 m

圖8 粗糙度z=1.0 m時(shí)氡氣濃度對(duì)比圖Fig.8 Radon concentration contrast graph at surface roughness of 1.0 m

3 結(jié)論

(1) 建立核素氡擴(kuò)散數(shù)學(xué)模型，利用CFD方法耦合求解得到不同風(fēng)速和下墊面粗糙度下的大氣風(fēng)場(chǎng)結(jié)構(gòu)及核素氡的濃度分布狀況。環(huán)境風(fēng)速和下墊面粗糙度對(duì)核素氡的遷移擴(kuò)散具有重要的影響。當(dāng)大氣風(fēng)速小于0.5 m/s時(shí)，此時(shí)地面粗糙度對(duì)核素氡遷移擴(kuò)散起主導(dǎo)作用。地面粗糙度越大，排風(fēng)口周圍近距離的氡氣濃度越高，局部污染越嚴(yán)重，主要污染范圍在下風(fēng)向100 m以內(nèi)；當(dāng)大氣風(fēng)速大于2.0 m/s時(shí)，大氣風(fēng)速對(duì)核素遷移擴(kuò)散起主導(dǎo)作用。大氣風(fēng)速越大，氡氣遷移擴(kuò)散能力越強(qiáng)，污染范圍越大，局部污染較輕，主要污染范圍在400 m以內(nèi)。

(2) 大氣中鈾礦通風(fēng)井排出的核素氡遷移擴(kuò)散過(guò)程非常復(fù)雜，影響的因素眾多，將在下一步的研究工作中考慮諸如排放口位置、排放高度、排風(fēng)風(fēng)速等參數(shù)的影響。同時(shí)，現(xiàn)場(chǎng)測(cè)試數(shù)據(jù)將被用于進(jìn)一步驗(yàn)證和完善數(shù)值計(jì)算模型。本文研究方法和研究結(jié)論對(duì)于鈾礦區(qū)的環(huán)境評(píng)價(jià)和治理有重要意義，可為鈾礦區(qū)輻射防護(hù)及新建鈾礦山選址提供參考。

[1]周星火, 李向陽(yáng), 王漢青.鈾礦排風(fēng)口放射性流出物對(duì)地表生態(tài)環(huán)境影響[C]//第二次全國(guó)天然輻射照射與控制研討會(huì)論文匯編.北京, 2005∶ 224-228.ZHOU Xinhuo, LI Xiangyiang, WANG Han-qing.The Impact of Uranic Well Ventilation ejection on the Ground Ecosystem Environment[C]//The Second National Natural Radiation and Controllability Proseminar Thesis Compilation.Beijing, 2005∶224-228.

[2]潘自強(qiáng), 王志波, 陳竹舟, 等.中國(guó)核工業(yè)三十年輻射環(huán)境質(zhì)量評(píng)價(jià)[M].北京∶ 原子能出版社, 1990∶ 36-43.PAN Zhiqiang, WANG Zhibo, CHEN Zhuzhou, et al.Radiation environmental impact assessment of the nuclear industry in china in the three decades[M].Beijing∶ Atomic Energy Press, 1990∶36-43.

[3]Ginzburg H M, Reis E.Consequences of the nuclear power plant accident at Chernobyl[J].Public Health Reports, 1991, 101(1)∶32-40.

[4]Terada H, Chino M.Development of an atmospheric dispersion model for accidental discharge of radionuclides with the function of simultaneous prediction for multiple domains and its evaluation by application to the chernobyl nuclear accident[J].Journal of Nuclear Science and Technology, 2008, 45(9)∶920-931.

[5]Leelossy á, Mészáros R, Lagzi I.Short and long term dispersion patterns of radionuclides in the atmosphere around the Fukushima Nuclear Power Plant[J].Journal of Environmental Radioactivity, 2011, 102(12)∶ 1117-1121.

[6]姚仁太, 郝宏偉, 黃杰, 等.田灣核電站場(chǎng)外事故后果評(píng)價(jià)系統(tǒng)—TW-NAOCAS[J].輻射防護(hù)通訊, 2006, 26(2)∶ 1-11.YAO Rentai, HAO Hongwei, HUANG Jie, et al.A real time system TW-NAOCAS accident consequence assessment for Tianwan NPP[J].Radiation Protection Bulletin, 2006, 26(2)∶1-11.

[7]Tseng C C, CHANG Nibin.Environmental exposure assessment for emergency response in a nuclear power plant using an integrated source term and 3D numerical model[J].Environment Modeling & Assessment, 2009, 14(6)∶ 661-675.

[8]王衛(wèi)星, 楊亞新, 王雷明, 等.廣東下莊鈾礦田土壤的天然放射性研究[J].中國(guó)環(huán)境科學(xué), 2005, 25(1)∶ 120-123.WANG Weixing, YANG Yaxin, WANG Leiming, et al.Studies on natural radioactivity of soil in Xiazhuang uranium ore field,Guangdong[J].China Environment Science, 2005, 25(1)∶120-123.

[9]胡恭任, 于瑞蓮.閩北 570鈾礦區(qū)周圍大氣環(huán)境氡濃度水平研究[J].地球與環(huán)境, 2004, 32(1)∶ 60-65.HU Gongren, YU Ruilian.Atmospheric radon in No.570 uranium deposit area, north Fujian province[J].Earth and Environment, 2004, 32(1)∶ 60-65.

[10]李舟, 楊忠, 魏濤, 等.貴州省某磷礦區(qū)氡氣及放射性水平的評(píng)價(jià)[J].中國(guó)輻射衛(wèi)生, 2006, 15(3)∶ 262-263.LI Zhou, YANG Zhong, WEI Tao, et al.Evaluation of radon and radioactivity level in phosphorite area in Guizhou[J].Chinese Journal of Radiological Health, 2006, 15(3)∶ 262-263.

[11]謝東, 王漢青, 李向陽(yáng), 等.我國(guó)鈾礦通風(fēng)尾氣放射性核素遷移的研究現(xiàn)狀及展望[J].鈾礦冶, 2008, 27(2)∶ 81-84.XIE Dong, WANG Hanqing, LI Xiangyang et al.Research situation and progress in diffusion of radionuclide radon ejection from uranium mine well ventilation[J].Uranium Mining and Metallurgy, 2008, 27(2)∶ 81-84.

[12]王漢青, 沈小建, 李向陽(yáng), 等.復(fù)雜山地大氣環(huán)境三維建模與仿真方法[J].建筑熱能通風(fēng)空調(diào), 2006, 25(2)∶ 13-16.WANG Hanqing, SHEN Xiaojian, LI Xiangyang, et al.Method of numerical simulation and three dimensional modeling on complex mountainous atmospheric environment[J].Building Energy & Environment, 2006, 25(2)∶ 13-16.

[13]XIE Dong, WANG Hanqing, DENG Qihong, et al.Investigation on different ventilation patterns in an office room[J].Journal of Central South University of Technology, 2007, 14(S3)∶ 66-70.

[14]Wang J S , Chan T L , Cheung C S, et al.Three-dimensional pollutant concentration dispersion of a vehicular exhaust plume in the real atmosphere[J].Atmospheric Environment, 2006, 40(3)∶484-497.

[15]Lutman E R, Jones S R, Hill R A, et al.Comparison between the predictions of a Gaussian plume model and a Lagrangian particle dispersion model for annual average calculations of long-range dispersion of radionuclides[J].Journal of Environmental Radioactivity, 2004, 75(3)∶ 339-355.

[16]GB 50009—2001, 建筑結(jié)構(gòu)荷載規(guī)范[S].GB 50009—2001, Load code for the design of building structure[S].

[17]王福軍.計(jì)算流體動(dòng)力學(xué)分析—CFD軟件原理與應(yīng)用[M].北京∶ 清華大學(xué)出版社, 2004∶ 114-158.WANG Fujun.Computational fluid dynamics analysis—CFD software principles and applications[M].Beijing∶ Tsinghua University Press, 2004∶ 114-158.