液化天然氣排放形成的羽流過(guò)程數(shù)值研究

張小斌，厲勁風(fēng)，邱利民（浙江大學(xué)制冷與低溫研究所，浙江省制冷與低溫技術(shù)重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 310027）

?

液化天然氣排放形成的羽流過(guò)程數(shù)值研究

張小斌，厲勁風(fēng)，邱利民
（浙江大學(xué)制冷與低溫研究所，浙江省制冷與低溫技術(shù)重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 310027）

摘要：大量LNG溢出蒸發(fā)后將形成低溫云團(tuán)，對(duì)風(fēng)向下游地面人員可能造成凍傷、燃燒以及缺氧窒息等危險(xiǎn)。基于計(jì)算流體力學(xué)方法，構(gòu)建了深低溫兩相多組分流動(dòng)的Navier-Stokes方程以及湍流封閉方程，考慮了空氣中水蒸氣由于溫度降低到飽和溫度以下而相變傳質(zhì)過(guò)程。由于氣相中氧氣等非液化氣體的存在，修正了計(jì)算水蒸氣相變率的Hertz-Knudsen方程。詳細(xì)給出了計(jì)算過(guò)程邊界條件設(shè)置，評(píng)估了地球自轉(zhuǎn)引起的科里奧利力影響特性。數(shù)值模擬了美國(guó)LLNL完成的Burro系列LNG排放羽流實(shí)驗(yàn)，發(fā)現(xiàn)考慮水蒸氣相變的計(jì)算結(jié)果，相比于未考慮水蒸氣相變的計(jì)算，更接近實(shí)驗(yàn)結(jié)果。研究結(jié)果對(duì)LNG接收終端等安全環(huán)境評(píng)估及設(shè)計(jì)有指導(dǎo)意義。

關(guān)鍵詞：羽流；天然氣；低溫；凝結(jié)；數(shù)值分析

2015-05-29收到初稿，2015-10-16收到修改稿。

聯(lián)系人：邱利民。第一作者：張小斌(1976—)，男，教授。

Received date: 2015-05-29.

Foundation item: supported by the Outstanding Youth Foundation of Zhejiang Province (R15E060001) and the National Natural Science Foundation of China (51576169).

引　言

液化天然氣（LNG）在儲(chǔ)存和運(yùn)輸中由于泄漏事故引起的安全問(wèn)題受到廣泛關(guān)注[1-2]。LNG接收終端儲(chǔ)罐的儲(chǔ)量規(guī)格普遍為10萬(wàn)多立方米，一旦破裂大量LNG（主要為甲烷）將頃刻溢出，于地面或海面快速相變蒸發(fā)后形成低溫云。低溫云團(tuán)大概以相等或稍小于空氣速度向下游傳播，云團(tuán)周邊氣體逐漸與空氣混合稀釋并升溫，同時(shí)將空氣中的水蒸氣液化成飽和液滴，形成羽流。云團(tuán)核心區(qū)在較長(zhǎng)時(shí)間內(nèi)保持低溫以及大于空氣的密度，因此云團(tuán)向下游運(yùn)動(dòng)時(shí)高度可能下降到地面。這將對(duì)地面人員造成凍傷、燃燒等危險(xiǎn)，也可能導(dǎo)致缺氧災(zāi)難。據(jù)報(bào)道，當(dāng)氧氣含量小于17.7%時(shí)，人將出現(xiàn)呼吸困難、頭昏乏力等現(xiàn)象[3]。云團(tuán)到達(dá)地面時(shí)離氣源的距離、組分含量及覆蓋面積等取決于環(huán)境風(fēng)速、濕度，氣源高度、溫度及流量，以及地形結(jié)構(gòu)等。因此，獲知不同條件下LNG儲(chǔ)罐溢流引起的羽流可能對(duì)地面人員造成危險(xiǎn)的安全距離等參數(shù)是建設(shè)LNG儲(chǔ)罐必須的安全條件。

Roberts等[4-5]首先得到基于簡(jiǎn)化假設(shè)的大氣擴(kuò)散問(wèn)題的理論解，該理論用以分析相關(guān)參數(shù)及源項(xiàng)的敏感性分析。高斯(Gaussian)羽流[6-7]被認(rèn)為是最簡(jiǎn)單的羽流情況，指點(diǎn)源連續(xù)無(wú)方向的向無(wú)限空間環(huán)境排放擴(kuò)散組分。高斯羽流理論解已被廣泛應(yīng)用于大工業(yè)氣體排放，工業(yè)及農(nóng)業(yè)養(yǎng)殖場(chǎng)氣味傳輸，火山爆發(fā)以及核、生化工業(yè)等的污染物排放預(yù)測(cè)，并最終形成污染物排放安全工業(yè)標(biāo)準(zhǔn)。

計(jì)算流體力學(xué)(CFD)方法也已被廣泛應(yīng)用于LNG排放形成羽流擴(kuò)散過(guò)程濃度的時(shí)域及地域預(yù)測(cè)。羽流擴(kuò)散CFD計(jì)算一般基于如Fluent，CFX及Star-CD等商業(yè)軟件[8]。另外也有專門(mén)用于羽流排放過(guò)程而開(kāi)發(fā)的數(shù)值計(jì)算軟件如SLAB、DEGADIS 和FEM3等，但這些程序目前只能計(jì)算穩(wěn)定邊界條件，不能計(jì)算組分濃度非穩(wěn)態(tài)波動(dòng)[9]。其中DEGADIS假設(shè)組分濃度在垂直方向滿足指數(shù)定律分布，水平方向滿足修正的高斯分布，同時(shí)假設(shè)垂直方向風(fēng)速滿足指數(shù)分布，則計(jì)算簡(jiǎn)化為一維模型，只能用于分析重氣體排放到平坦、無(wú)障礙地形的情況。SLAB也是基于簡(jiǎn)化的一維模型，DEGADIS也假設(shè)濃度的橫向和高度方向的分布函數(shù)，不能模擬建筑物等對(duì)組分傳播的影響[10-11]。在羽流數(shù)值分析中廣泛使用的商業(yè)軟件為FEM3，這是個(gè)基于有限元方法的模型包，由美國(guó)Lawrence Livermore National Laboratory (LLNL)開(kāi)發(fā)，最初湍流模型采用混合長(zhǎng)理論（一方程模型），現(xiàn)在已經(jīng)包含了κ-ε二方程模型。另外日本三菱重工也為羽流擴(kuò)散的數(shù)值分析開(kāi)發(fā)了CFD有限元軟件STD。

20世紀(jì)80年代LLNL進(jìn)行了多次不同LNG排放條件的羽流實(shí)驗(yàn)，包括Burro系列及Coyote系列[12]，這些實(shí)驗(yàn)得到了不同位置非穩(wěn)態(tài)變化的濃度和溫度，成為驗(yàn)證排放過(guò)程CFD計(jì)算結(jié)果準(zhǔn)確性的基礎(chǔ)數(shù)據(jù)。如Sklavounos等[10-11]基于CFX軟件計(jì)算，通過(guò)實(shí)驗(yàn)對(duì)比評(píng)估了基于標(biāo)準(zhǔn)κ-ε模型，標(biāo)準(zhǔn)κ-ω和SST湍流模型的計(jì)算結(jié)果，發(fā)現(xiàn)這些湍流模型都能產(chǎn)生相同精度的結(jié)果，但κ-ω低估了最大濃度，而其他模型高估。Gavelli等[13]基于Fluent計(jì)算了有建筑物時(shí)的Falcon系列實(shí)驗(yàn)過(guò)程，對(duì)比了濃度和溫度分布，指出DEGADIS計(jì)算結(jié)果不能滿足工程設(shè)計(jì)要求。而Tauseef等[14]利用Fluent模擬了建筑物對(duì)組分離散的影響，對(duì)比了不同湍流模型結(jié)果，發(fā)現(xiàn)Real κ-ε模型產(chǎn)生了最合適和精確的濃度分布。

上述對(duì)于大流量LNG排放到大氣環(huán)境形成復(fù)雜兩相多組分?jǐn)U散羽流過(guò)程的CFD研究，都沒(méi)有考慮空氣和低溫氣體強(qiáng)烈湍流混合后當(dāng)?shù)販囟燃眲〗档?，?dǎo)致水蒸氣飽和形成液滴的相變物理過(guò)程。然而LLNL的實(shí)驗(yàn)明顯觀察到形成水霧特性，因此忽略水蒸氣相變過(guò)程的羽流CFD計(jì)算的準(zhǔn)確性，需要驗(yàn)證。

本文構(gòu)建了考慮水蒸氣相變的基于混合物模型(mixture model)的兩相多組分?jǐn)?shù)值模型，并用Fluent14.5求解控制方程。數(shù)值模擬了Burro系列實(shí)驗(yàn)中LNG排放過(guò)程。

1　兩相多組分流動(dòng)過(guò)程的數(shù)值模型

1.1基本控制方程

對(duì)于混合物模型，將氣液兩相看成是混合物單相，其連續(xù)性方程為

假設(shè)形成的液滴與氣相運(yùn)動(dòng)速度相同，即兩者沒(méi)有動(dòng)量交換，則動(dòng)量方程為

式中，mm=alml+agmg為混合物黏度；P為壓力；g為重力加速度矢量。

能量方程為

有效熱導(dǎo)率keff=alkl+agkg，hk為k相焓值。

氣相體積含量ag方程為

考慮水蒸氣相變傳質(zhì)，組分i在氣相混合物中的質(zhì)量含量Ygi，由式(5)計(jì)算

上角標(biāo)或下角標(biāo)g、l分別為氣相和液相。對(duì)于湍流，式(5)中的擴(kuò)散流量Ji由式(6)計(jì)算

Sct是湍流有效Schmidt數(shù)，μt/Sct是湍流引起的有效質(zhì)量擴(kuò)散系數(shù)。因此，在湍流條件下的有效擴(kuò)散包含作為已知條件輸入的層流擴(kuò)散系數(shù)Di,m，加上增加的湍流擴(kuò)散效應(yīng)。湍流Sct表示動(dòng)量和質(zhì)量的相對(duì)擴(kuò)散，數(shù)量級(jí)大概為1，且對(duì)分子流體屬性不敏感，因此采用缺省的常數(shù)Sct=0.7是合理的。另外，湍流擴(kuò)散效應(yīng)遠(yuǎn)遠(yuǎn)大于層流擴(kuò)散效應(yīng)，采用常層流擴(kuò)散系數(shù)，Di,m=2.8×10-5m·s-1也是可以接受的。

為封閉式(1)～式(5)，必須計(jì)算相變引起的水汽化和液化量

1.2相變傳質(zhì)模型

假設(shè)氣液界面溫度平衡，界面為幾個(gè)自由分子厚度形成的Knudsen層，水蒸氣為理想氣體分子，基于分子運(yùn)動(dòng)論可得到計(jì)算氣液相變率的Hertz-Knudsen關(guān)系[15]

式中，h為液化或蒸發(fā)系數(shù)。式(7)由純組分條件推導(dǎo)得到，當(dāng)考慮到氣相中存在非液化的氧氣和氮?dú)鈺r(shí)，Pv指的是水蒸氣的分壓力。其他非液化氣體對(duì)水蒸氣相變率的影響，主要通過(guò)修正通用氣體常數(shù)R來(lái)體現(xiàn)[16]

R=8314.34 J·kmol-1·K-1，Mi為組分i的分子量。多組分氣相中，水蒸氣分壓PH2O由式(9)計(jì)算

式中，yH2O為水蒸氣摩爾含量。在實(shí)際計(jì)算中，通常以溫度作為判斷是否發(fā)生相變的標(biāo)準(zhǔn)。利用Clausius-Clapeyron方程，式(8)可重寫(xiě)為

式中，Tsat為飽和溫度，為水蒸氣分壓PH2O的函數(shù)；rg為飽和水蒸氣密度；L為潛熱。需注意的是傳質(zhì)率m單位為kg·m-2·s-1。在數(shù)值計(jì)算中，需要轉(zhuǎn)變?yōu)閱挝惑w積的傳質(zhì)率。假設(shè)生成的液滴直徑相同，則界面積密度為

Vcell為網(wǎng)格體積，d為液滴直徑。由式(10)和式(11)得到

T*為當(dāng)?shù)貙?shí)際溫度。引入系數(shù)C（L·s-1）

則式(12)最后簡(jiǎn)化為

需要說(shuō)明的是，液化與蒸發(fā)過(guò)程系數(shù)C絕對(duì)值不相等，具體值通過(guò)和實(shí)驗(yàn)對(duì)比得到。由于Tsat為水蒸氣分壓力的函數(shù)，因此式(14)需通過(guò)自定義函數(shù)(UDF)的方法，耦合到Fluent的程序中進(jìn)行迭代計(jì)算。

為簡(jiǎn)化計(jì)算，并沒(méi)有考慮可能發(fā)生的水液固相變?？紤]到液態(tài)水在混合物中體積含量很小，因此雖然液固相變時(shí)密度從1000 kg·m-3減小到冰的900 kg·m-3，對(duì)流場(chǎng)的影響認(rèn)為可忽略不計(jì)。同樣的原因，雖然比熱容從水的4200 J·kg-1·K-1減小到冰的2100 J·kg-1·K-1，對(duì)溫度場(chǎng)的影響也認(rèn)為可忽略不計(jì)。

1.3湍流模型

κ-ε二方程湍流模型已經(jīng)被成功地用于大氣羽流的數(shù)值模擬[17-18]。文獻(xiàn)[8]指出，缺省的κ-ε二方程模型參數(shù)對(duì)計(jì)算結(jié)果的濃度分布誤差很小。Sklavounos等[10]利用軟件CFX建模了Coyote實(shí)驗(yàn)，也發(fā)現(xiàn)使用缺省的湍流參數(shù)與實(shí)驗(yàn)結(jié)果吻合很好。因此，本文計(jì)算采用缺省的Realizable κ-ε二方程湍流模型。

1.4科里奧利力影響

由于地球自轉(zhuǎn)，大氣流動(dòng)坐落于有加速的參考坐標(biāo)中，運(yùn)動(dòng)方程有必要考慮科里奧利力的影響。量綱1參數(shù)Rossby數(shù)(Ro)比較了慣性力與科里奧利力的大小[17]

式中，V、D和f分別為特征速度、特征長(zhǎng)度和科里奧利參數(shù)f=2?sin(φ)，φ為緯度角。由于慣性力在空氣流動(dòng)過(guò)程一般占主導(dǎo)作用，因此，Ro通常用來(lái)判斷是否需要考慮科里奧利力。如果Ro較小，則科里奧利力應(yīng)該考慮。以溢出10萬(wàn)立方米 LNG到水面為例，如果風(fēng)速為2 m·s-1，f約10-4s-1，D取天然氣云最大值約1 km，那么Ro約20，說(shuō)明科里奧利力相比慣性力不重要。本文忽略科里奧利力影響，目前也未見(jiàn)考慮科里奧利力影響的羽流CFD數(shù)值計(jì)算報(bào)道。

2　CFD數(shù)值模型邊界條件

為驗(yàn)證上述數(shù)理方程組及相變模型的準(zhǔn)確性，選擇Burro系列[19]實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比。實(shí)驗(yàn)過(guò)程獲得了不同位置的溫度、風(fēng)速和CH4組分非穩(wěn)態(tài)變化，表1給出了實(shí)驗(yàn)中相關(guān)參數(shù)。

表1　Burro 9系列LNG排放實(shí)驗(yàn)參數(shù)[19]Table 1　Test conditions of LNG dispersion in Burro 9[19]

為將上述物理過(guò)程轉(zhuǎn)變?yōu)镃FD數(shù)值模型，進(jìn)行了與文獻(xiàn)[10]同樣的簡(jiǎn)化，包括：

① 由于Burro 9中主要組分皆為甲烷氣體(CH4)，且所占比重較大（Burro 9為87.4%，Burro 9 為83.1%），忽略少量乙烷及丙烷不會(huì)對(duì)流體物性產(chǎn)生較大影響。因此不考慮LNG組分組成，蒸發(fā)的天然氣近似為100% 甲烷氣體；

② 本文模擬LNG在水池中泄漏并快速相變之后以氣態(tài)NG在環(huán)境中擴(kuò)散的過(guò)程，因此忽略LNG在水池中泄漏并快速相變的過(guò)程，整個(gè)水池液面(d=58 m)為氣體CH4進(jìn)口邊界，溫度T=111.5 K，速度根據(jù)表1給定的LNG排放量計(jì)算，水池面積則CH4流速為

③ 忽略現(xiàn)實(shí)中存在的風(fēng)向隨機(jī)變化，羽流關(guān)于水池中心垂直平面對(duì)稱。

基于上述假設(shè)，構(gòu)建了如圖1所示的平面對(duì)稱CFD數(shù)值模型，該模型和文獻(xiàn)[10]幾乎一致，除了后者沒(méi)有利用平面對(duì)稱條件。

圖1　Burro系列實(shí)驗(yàn)的CFD建模邊界條件Fig.1　CFD geometrical model and boundary conditions of Burro serial LNG spill tests

上述計(jì)算域有7個(gè)邊界，邊界條件設(shè)置分別如下：

① 右邊出口為定壓邊界條件，設(shè)置相對(duì)壓力P=0，平均壓力為大氣壓；

② 計(jì)算域頂面高度足夠大，使得底部CH4的進(jìn)入對(duì)頂部流場(chǎng)影響可忽略不計(jì)。因此，理論上，可設(shè)置為定壓邊界，但頂部定壓邊界將改變左邊空氣入口的速度分布。由此，在忽略風(fēng)向隨機(jī)變化的簡(jiǎn)化條件下，頂面設(shè)置為對(duì)稱邊界條件，滑移速度為0，邊界速度與流場(chǎng)速度相等，也沒(méi)有流體穿過(guò)對(duì)稱面；

③ 由假設(shè)③，圖中計(jì)算域正面為對(duì)稱平面；

④ 計(jì)算域的背面，保證計(jì)算域足夠?qū)?250 m)，則CH4進(jìn)入對(duì)該面的影響可忽略不計(jì)，因此相同于計(jì)算域頂面，背面也設(shè)置為對(duì)稱面；

⑤ 計(jì)算域底面（模擬地面），由于實(shí)際LNG排放時(shí)間較短（約100 s），與地面換熱只影響到很薄的邊界層，對(duì)CH4宏觀擴(kuò)散影響可忽略不計(jì)，因此設(shè)置為絕熱面(Q=0)；

⑥ CH4進(jìn)口為速度進(jìn)口邊界條件，其值為假設(shè)②計(jì)算。入口湍流強(qiáng)度及黏度比都為4%。入口I值實(shí)際是未知值，因此有必要進(jìn)行I對(duì)計(jì)算結(jié)果影響的敏感性分析。

⑦ 計(jì)算域左邊空氣進(jìn)口條件最為復(fù)雜，因?yàn)榭諝馑俣?、溫度和濕度都是高度的函?shù)。M-O相似律[12]給出了精確計(jì)算方法，但該方法較復(fù)雜。另外一種更簡(jiǎn)單的計(jì)算高度方向速度分布的方法是指數(shù)定律

式中，zr為參考高度，Ur為在高度zr的水平風(fēng)速，z為高度坐標(biāo)，指數(shù)p為按照Pasquill- Gifford分類方法由大氣穩(wěn)定性條件決定的系數(shù)[20]。該方法雖然沒(méi)有M-O相似律精確，但在低高度范圍內(nèi)兩者吻合很好，因此，已廣泛應(yīng)用于該類問(wèn)題的CFD0.1計(jì)5算[13,21]。由表1可知，對(duì)在計(jì)算中，利用UDF方法將上式公式指定為邊界上的速度分布。

入口空氣溫度分布，假設(shè)為表1中所給的常數(shù)，文獻(xiàn)[8,10]也作同樣處理。

入口空氣濕度分布，也假設(shè)為常數(shù)，根據(jù)表1中相對(duì)濕度r計(jì)算得到體積含量

式(1)～式(5)以及式(14)加上Realizable κ-ε湍流模型，構(gòu)成了封閉的基于混合物模型的兩相多組分?jǐn)?shù)值框架，利用CFD商業(yè)軟件Fluent14.5進(jìn)行迭代求解。其中連續(xù)性方程和動(dòng)量方程通過(guò)二階迎風(fēng)方案離散，利用Couple方案計(jì)算速度與壓力修正。計(jì)算收斂標(biāo)準(zhǔn)為能量方程殘差為10-6，其他方程殘差為10-3。計(jì)算過(guò)程，純組分氧氣、氮?dú)?、CH4及水蒸氣的熱導(dǎo)率、比熱容及黏度都為溫度的函數(shù)，液態(tài)水的物性為常數(shù)，物性參數(shù)來(lái)自Refprop 8.0。整個(gè)計(jì)算域超過(guò)17.73萬(wàn)個(gè)網(wǎng)格，將CH4進(jìn)口附近局部加密后網(wǎng)格超過(guò)19.12萬(wàn)個(gè)，參照Burro 9實(shí)驗(yàn)，以相對(duì)排放源中心為原點(diǎn)，坐標(biāo)x=397 m、y=48.75 m處的T3傳感器為比較對(duì)象，網(wǎng)格獨(dú)立性檢驗(yàn)計(jì)算表明基于17.73萬(wàn)個(gè)網(wǎng)格數(shù)的CH4濃度值幾乎和后者一致，見(jiàn)圖2。

圖2　不同網(wǎng)格數(shù)量下計(jì)算的Burro 9 T3傳感器的CH4濃度非穩(wěn)態(tài)變化Fig.2　Simulated CH4concentration of Burro 9 T3 sensor with different mesh number

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

首先對(duì)式（14）中系數(shù)C對(duì)相變率的定量影響進(jìn)行評(píng)估。實(shí)驗(yàn)結(jié)果顯示，液化系數(shù)Cc一般要大于蒸發(fā)系數(shù)Ce[15]，Rubel等[22]的測(cè)量更進(jìn)一步指出Ce≈1.2Cc。采用兩組系數(shù)對(duì)，Ce=1，Cc=1.2及Ce=100，Cc=120，兩者相差100倍，對(duì)Burro 9實(shí)驗(yàn)進(jìn)行了數(shù)值計(jì)算，結(jié)果如圖3所示。由圖可見(jiàn)，兩者計(jì)算結(jié)果吻合較一致。因此，計(jì)算中采用液化和蒸發(fā)系數(shù)對(duì)：Ce=100和Cc=120。

圖3　不同系數(shù)C時(shí)計(jì)算的Burro 9在不同位置的CH4濃度非穩(wěn)態(tài)變化Fig.3　Simulated CH4concentration of Burro 9 with different evaporation and condensation coefficient pairs

圖4和圖5分別給出了考慮相變及未考慮相變時(shí)，計(jì)算的Burro 9溫度和CH4濃度分布與實(shí)驗(yàn)測(cè)量的對(duì)比[19，23]。根據(jù)實(shí)驗(yàn)布置，圖中傳感器位置：以圖1中CH4進(jìn)口中心為坐標(biāo)原點(diǎn): G6: x=137.9 m，y=24.3 m；G15: x=385.25 m，y=107.55 m；T3: x=397 m，y=48.75 m；T4: x=132.32 m，y=46.47 m。這些坐標(biāo)值根據(jù)實(shí)際傳感器位置以及風(fēng)向平均偏轉(zhuǎn)角修正計(jì)算得到。由圖可見(jiàn)，無(wú)論溫度還是濃度，實(shí)驗(yàn)結(jié)果呈無(wú)規(guī)律波動(dòng)，這是由實(shí)驗(yàn)過(guò)程上游風(fēng)速和風(fēng)向隨機(jī)變化導(dǎo)致。而數(shù)值計(jì)算固定風(fēng)向和風(fēng)速，因此結(jié)果曲線光滑?？傮w上，從實(shí)驗(yàn)數(shù)據(jù)依然能捕捉到濃度和溫度的宏觀變化，因此可用于檢驗(yàn)數(shù)值結(jié)果。在圖4中，CFD計(jì)算的不同位置的溫度，無(wú)論是考慮相變或不考慮相變，和實(shí)驗(yàn)值趨勢(shì)一致，但最低溫度都明顯低于實(shí)驗(yàn)值，這和數(shù)值模型與實(shí)驗(yàn)邊界條件不完全一致有關(guān)。首先，實(shí)際風(fēng)向及風(fēng)速隨機(jī)高頻率變化導(dǎo)致實(shí)驗(yàn)時(shí)云團(tuán)中心平面與計(jì)算模型的云團(tuán)中心平面不一致，因此實(shí)驗(yàn)中傳感器位置與對(duì)應(yīng)的計(jì)算模型中的坐標(biāo)，相對(duì)各自云團(tuán)中心平面位置卻不一致（因?yàn)閷?shí)驗(yàn)時(shí)云團(tuán)中心平面位置隨風(fēng)向在變化）。顯然云團(tuán)中心平面上溫度最低，導(dǎo)致測(cè)量溫度比計(jì)算值高。另外，傳感器離地面1米，云團(tuán)和地面的換熱以及地面的輻射對(duì)會(huì)導(dǎo)致云團(tuán)溫度上升?？傮w上，CFD計(jì)算結(jié)果能夠得到可接受的溫度非穩(wěn)態(tài)變化。另外從圖中可發(fā)現(xiàn)，考慮相變的CFD模型，得到的溫度總是大于未考慮相變的CFD模

圖4　計(jì)算的Burro 9 溫度分布與實(shí)驗(yàn)測(cè)量對(duì)比(h=3 m)Fig.4　Comparison of simulated temporal temperature curves with experimental records for Burro 9 (h=3 m)

圖5　計(jì)算的Burro 9 CH4濃度分布與實(shí)驗(yàn)對(duì)比(h=1 m)Fig.5　Comparison of simulated temporal CH4concentration curves with experimental records for Burro 9 (h=1 m)

對(duì)CH4濃度變化，由圖5可見(jiàn)，未考慮相變的CFD模型，幾乎在受云團(tuán)影響的時(shí)間內(nèi)，在所有位置比考慮相變模型的CFD計(jì)算結(jié)果高。這是由于后者溫度普遍比前者高（圖4），使得氣體黏度降低，Re增加，云團(tuán)中湍流強(qiáng)度更加強(qiáng)烈，由此導(dǎo)致氣體分子運(yùn)動(dòng)更活躍，CH4組分?jǐn)U散及流動(dòng)傳質(zhì)加劇，因此CH4更容易被空氣稀釋，濃度相對(duì)降低。總體上，在所有位置計(jì)算的CH4濃度值和實(shí)驗(yàn)測(cè)量的平均值趨勢(shì)和大小基本一致，但G15位置計(jì)算值始終要小于實(shí)驗(yàn)值。由于G15距離CH4源最遠(yuǎn)，風(fēng)向偏移造成的影響最大。

以Burro 9的G6和T3傳感器所得到的CH4濃度分布實(shí)驗(yàn)數(shù)據(jù)為參考，圖6給出了Fluent模擬結(jié)果與FEM3軟件模擬結(jié)果[8]的對(duì)比。由圖6可見(jiàn)，基于Fluent的模擬結(jié)果要比基于FEM3的模擬結(jié)果更接近實(shí)驗(yàn)數(shù)據(jù)。

圖6　計(jì)算的Burro 9 CH4濃度分布與實(shí)驗(yàn)及FEM3軟件結(jié)果對(duì)比 (h=1 m)Fig.6　Comparison of Fluent simulated temporal CH4concentration curves with experimental records and simulation results of FEM3 for Burro 9 (h=1 m)

4　結(jié)　論

本文構(gòu)建了兩相多組分流動(dòng)的Navier-stokes方程以及湍流封閉方程，給出了基于氣體動(dòng)力學(xué)的羽流過(guò)程水蒸氣相變計(jì)算的修正模型。計(jì)算了Burro系列實(shí)驗(yàn)中LNG排放羽流實(shí)驗(yàn)，得到的溫度及濃度非穩(wěn)態(tài)變化與實(shí)驗(yàn)進(jìn)行了對(duì)比，得出如下結(jié)論。

（1）大流量低溫流體如LNG等排放形成羽流過(guò)程的CFD數(shù)值模型，空氣中水蒸氣飽和形成雨滴的物理過(guò)程不可忽略?？紤]水蒸氣相變過(guò)程的計(jì)算結(jié)果，相比于未考慮水蒸氣相變過(guò)程的計(jì)算，與實(shí)驗(yàn)結(jié)果相比更一致。

（2）基于Fluent的數(shù)值計(jì)算，相對(duì)于FEM3等傳統(tǒng)的擴(kuò)散模擬軟件的模擬計(jì)算，與實(shí)驗(yàn)結(jié)果擬合度更高。

（3）構(gòu)建的基于mixture model兩相多組分?jǐn)?shù)值模型及基于氣體動(dòng)力學(xué)的水蒸氣相變計(jì)算模型，能準(zhǔn)確模擬羽流過(guò)程的水蒸氣相變過(guò)程。若能在模型中增加來(lái)流空氣隨機(jī)性波動(dòng)的考慮，將進(jìn)一步提高數(shù)值計(jì)算結(jié)果的準(zhǔn)確性。

References

[1] 錢(qián)新明，劉牧，劉振翼. 隧道內(nèi)液化天然氣管道泄漏火災(zāi)溫度場(chǎng)的數(shù)值模擬 [J]. 化工學(xué)報(bào)，2009，60(12): 3184-3188. QIAN X M，LIU M，LIU Z Y. Numerical simulation of LNG leakage fire temperature distribution from pipeline in tunnel [J]. CIESC Journal，2009，60(12): 3184-3188.

[2] 唐建峰，蔡娜，郭清，等. LNG垂直噴射源連續(xù)泄漏擴(kuò)散的模擬[J]. 化工學(xué)報(bào)，2013，64(3): 1124-1131. DOI: 10.3969/j.issn. 0438-1157.2013.03.048. TANG K F，CAI N，GUO Q，et al. Simulation of LNG diffusion: a continuous vertical jet release [J]. CIESC Journal，2013，64(3): 1124-1131. DOI: 10.3969/j.issn.0438-1157.2013.03.048.

[3] FLYNN T M. Cryogenic Engineering[M]. New York: Marcel Dekker Inc.，2005.

[4] ROBERTS O F T. The theoretical scattering of smoke in a turbulent atmosphere[C]//Proceedings of The Royal Society of London Series A-Containing Papers of A Mathematical and Physical Character. London，England: Royal Soc London，1923，104: 640-654.

[5] SUTTON O G. A theory of eddy diffusion in the atmosphere [C]//Proceedings of The Royal Society of London Series A-Containing Papers of A Mathematical And Physical Character. London，England: Royal Soc London，1932，135: 143-165.

[6] STOCKIE J M. The mathematics of atmospheric dispersion modeling [J]. SIAM Review，2011，53(2): 349-372.

[7] SEINFELD J H. Atmospheric Chemistry and Physics: from Air Pollution to Climate Change[M]// PANDIS S N. New York: JohnWiley & Sons Inc.，1998.

[8] LUKETA-HANLIN A，KOOPMAN R，ERMAK D L. On the application of computational fluid dynamics codes for liquefied natural gas dispersion [J]. Journal of Hazardous Materials，2007，140: 504-517.

[9] OHBA R，KOUCHI A，HARA T，et al. Validation of heavy and light gas dispersion models for the safety analysis of LNG tank [J]. Journal of Loss Prevention in the Process Industries，2004，17: 325-337.

[10] SKLAVOUNOS S，RIGAS F. Simulation of Coyote series trials(Ⅰ): CFD estimation of non-isothermal LNG releases and comparison with box-model predictions [J]. Chem. Eng. Sci.，2006，61: 1434- 1443.

[11] SKLAVOUNOS S，RIGAS F. Validation of turbulence models in heavy gas dispersion over obstacles [J]. Journal of Hazardous Materials，2004，A108: 9-20.

[12] ERMAK D L，OHAPMAN R C V，GOLDWIRE H C，et al. Heavy gas dispersion test summary report[R]. 1989，ESL-TR-88-22.

[13] GAVELLI F，BULLISTER E，KYTOMAA H. Application of CFD (Fluent) to LNG spills into geometrically complex environments [J]. Journal of Hazardous Materials，2008，159: 158-168.

[14] TAUSEEF S M，RASHTCHIIAN D，ABBASI S A. CFD-based simulation of dense gas dispersion in presence of obstacles [J]. Journal of Loss Prevention in the Process Industries，2011，24: 371-376.

[15] MAREK R，STRAUB J. Analysis of the evaporation coefficient and the condensation coefficient of water [J]. International Journal of Heat and Mass Transfer，2001，44: 39-53.

[16] GRAYSON G，LOPEZ A，CHANDLER F. CFD modeling of helium pressurant effects on cryogenic tank pressure rise rates in normal gravity[C]// AIAA. Cincinnati，Ohio: 2007，MSFC-448.

[17] ALINOT C，MASSON C. k-ε model for the atmospheric boundary layer under various thermal stratifications [J]. Solar Eng.，2005，127: 438-443.

[18] HUSER A，NILSEN P J，SKATUN H. Applications of k-ε model to the stable ABL: pollution in complex terrain [J]. Wind Eng. Ind. Aero，1997，67/68: 425-436.

[19] KOOPMAN R P，BAKER J，CEDERWALL R T. Burro series data report LLNL/NWC 1980 LNG spill tests[R]. 1982，UCID-19075.

[20] CCPS. Consequence analysis of chemical releases[C]// AIChE. New York: Wiley-Blackwell，1999: 80-83

[21] MCBRIDE M A，REEVES A B，VANDERHEYDEN M D，et al. Use of advanced techniques to model the dispersion of chlorine in complex terrain [J]. Process Safety & Environmental Protection，2001，79(2): 89-102.

[22] RUBEL G O，GENTRY J W. Measurement of the kinetics of solution droplets in the presence of adsorbed mono-layers: determination of water accommodation coefficients [J]. Phys. Chem.，1984，88(14): 3142-3148.

[23] CHAN S T，RODEAN H C，ERMAK D L. Numerical simulations of atmospheric releases of heavy gases over variable terrain[C]//Air Pollution Modeling and Its Applications Ⅲ. New York: Plenum Press，1984，5: 295-328.

Numerical study on plume characteristics of liquefied natural gas spills

ZHANG Xiaobin，LI Jingfeng，QIU Limin
（Key Laboratory of Refrigeration and Cryogenic Technology of Zhejiang Province，Institute of Refrigeration and Cryogenics，Zhejiang University，Hangzhou 310027，Zhejiang，China）

Abstract:When a large number of LNG suddenly leaks and evaporates，it will form a low temperature cloud in the wind downstream，which may cause frostbite，burn and oxygen deficit hazard to the ground personnel. Based on the computational fluid dynamics (CFD) method，the two-phase multicomponent flow mathematical framework together with the turbulence closure are built to model the cryogenic flow. The phase-change mass transfer of water vapor in the air due to the temperature depression is considered. Because of the existence of non-liquefied gases such as oxygen，the Hertz-Knudsen equation for calculating the mass transfer rate of water vapor is modified. Detailed methods for setting the boundary conditions of the computational domain are presented and the influence of Coriolis force caused by the earth rotation is evaluated. The Burro experimental series of LNG released by LLNL are simulated and the results are used to evaluate the numerical models. It is found that the results using models with the phase change of water vapor are closer to the experimental results than that without the phase change. The studies are of directive significance for the safety environment assessment and design of LNG received terminal.

Key words:plume flow; NG; cryogenic; condensation; numerical analysis

DOI：10.11949/j.issn.0438-1157. 20150740

中圖分類號(hào)：TQ 021.4

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

文章編號(hào)：0438—1157（2016）04—1225—08

基金項(xiàng)目：浙江省杰出青年基金項(xiàng)目（R15E060001）；國(guó)家自然科學(xué)基金項(xiàng)目（51576169）。

Corresponding author:Prof. QIU Limin，liminqiu@zju.edu.cn