LNG船泄漏事故液池擴展計算及不確定性分析

莊學強，孫 迪，高孝洪，李格升

(1.集美大學 輪機工程學院，福建 廈門 361021；2.福建省船舶與海洋工程重點實驗室，福建 廈門 361021；3.武漢理工大學 能源與動力工程學院，武漢 430070)

LNG船泄漏事故液池擴展計算及不確定性分析

莊學強1，2，孫 迪1，2，高孝洪3，李格升3

(1.集美大學 輪機工程學院，福建 廈門 361021；2.福建省船舶與海洋工程重點實驗室，福建 廈門 361021；3.武漢理工大學 能源與動力工程學院，武漢 430070)

根據小孔射流、慣性-質量平衡擴展等理論搭建包括泄漏源強和液池半徑擴展在內的LNG船泄漏液池擴展模型，分析影響液池擴展變化的重要參數及其不確定性，借助MATLAB編程計算這些不確定對預測結果的影響。研究結果表明：導致LNG船泄漏的事件不同，泄漏口的大小可能會有所不同(0～16 m2不等)，當泄漏口面積<2 m2時，液池擴展對泄漏口大小的變化相當敏感；受泄漏口形狀和流動狀態等因素影響，流量系數取值范圍較廣(0.4～1.0)，預測偏差可達2.5倍；因LNG-水攪合程度及沸騰形式不同，液池蒸發速率為0.042～0.25 kg/(m2·s)，預測的最大液池半徑差別較大。

水路運輸；液池；LNG船；泄漏

據統計,自2010年以來，我國液化天然氣(Liquefied Natural Gas,LNG)進口量的年均增長率為30.5%。相關數據表明，到2020年末，我國很可能成為僅次于日本的全球第二大LNG進口國。[1]LNG船作為天然氣供應鏈中重要的海運裝備，承擔著遠距離、大批量輸送天然氣的重任，直接關系到我國能源規劃的順利實施。然而,LNG船載有大量易燃易爆的液體，一旦因遭受碰撞、擱淺和恐怖襲擊等事故而泄漏，將導致災難性后果。

LNG船舶泄漏事故的危害計算主要包括泄漏源強計算、液池擴展計算、可燃蒸氣云擴散計算、火災爆炸危害計算等。[2]在這些計算中,前一步計算的結果是后一步計算的輸入,后一步計算的準確度在很大程度上依賴于前一步的計算精度。目前研究最多的是LNG船泄漏蒸氣云擴散漂移，得到許多計算模型。但是,對LNG船泄漏液池擴展的研究較少，主要原因在于LNG船泄漏液池擴展是一個復雜的傳質與傳熱耦合過程，影響因素有很多，存在諸多不確定性。

1 LNG船泄漏液池擴展計算模型建立

LNG船發生泄漏事故時，推動液池在水面擴展的主動力是泄漏入水的LNG，泄漏源強是液池擴展計算的必要輸入參數，因此一個完整的LNG船泄漏液池擴展計算模型應包括泄漏源強子模型和液池擴展子模型。下面根據小孔射流、質量平衡及慣性擴展理論建立微分方程。

1.1LNG船泄漏源強計算模型

定義Δt為時間步長，i為時間步數。根據伯努利方程可得LNG船泄漏速率qi為

(1)

(2)

ΔQ=qiΔt

(3)

Qi=Qi-1+ΔQ

(4)

式(1)～式(4)中：Chole為泄漏孔的流量系數；Ahole為泄漏孔的流通面積，m2；PT為LNG貨艙液面上方的壓力，Pa；Patm為大氣壓力，101 325 Pa；ρL為LNG液體密度，kg/m3；hi為第i時間步時LNG液面距離泄漏孔中心線的長度，m；ΔQ為第i時間步內的泄漏量,kg；Qi為總泄漏量,kg；AT，i-1為第i-1時間步對應的液面橫截面積，m2。

1.2LNG船泄漏液池擴展模型

LNG船泄漏液池半徑擴展變化過程的計算式可表示為

Mi=Mi-1+qiΔt-Mevap

(5)

(6)

Vi=Mi/ρl

(7)

(8)

1) 若δi>δmin,則液池半徑擴展速度[3]為

(9)

(10)

2) 若δi≤δmin,則液池半徑為

(11)

式(5)～式(11)中：Ri為迭代步為i時LNG液池的半徑，m；Ks為擴展常數，取1.414；δi為迭代步為i時LNG液池的厚度，m；δmin為LNG液池最小穩定厚度，m；g為重力加速度,9.8 m/s2；Mi為迭代步為i時LNG液池的質量，kg；Vi為迭代步為i時LNG液池的體積，m3；Mevap為Δt時間步長內的蒸發量，kg；mi為Δt時間步長內LNG泄漏的速率,kg/s;mevap為單位面積內的蒸發速率，kg/(m2·s)。

2 LNG船泄漏液池擴展不確定性分析

以上主要根據小孔射流、質量平衡及慣性擴展理論建立LNG船泄漏液池擴展模型，國外許多著名機構(如FAY,QUEST和ABS等)也根據類似理論建立相應的模型。然而,LNG船泄漏及液池在水面擴展是一個非常復雜的過程，許多影響因素(如泄漏口大小、流量系數、液池蒸發速率及儲罐內的壓力等)在突發事故中難以確定，因此采用目前的模型對LNG船泄漏液池進行擴展模擬尚存在一定的不確定性。FAY,QUEST和ABS曾針對同一假定事故分別進行模擬計算，得出的結果相差4倍多。[4]

2.1泄漏口大小的不確定性

對于LNG船而言,無論是碰撞、擱淺等操作性事故導致的泄漏，還是受到恐怖襲擊導致的泄漏，泄漏口大小均難以確定。美國Sandia國家實驗室[5]在建立一系列LNG船的三維有限元結構模型分析各種船舶事故和恐怖襲擊可能造成的損害之后認為：由碰撞、擱淺等事故造成的LNG液艙破損面積≤2 m2；由恐怖襲擊造成的LNG液艙破損面積為2～16 m2。由各行業組成的LNG危害評估小組[6]在研究之后認為：由碰撞、擱淺等事故導致LNG船舶泄漏口的寬度為0～1.5 m，建議在進行風險評估時取平均值0.75 m作為最大可信泄漏孔徑。挪威船級社(Det Norske Veritas,DNV)認為LNG船遭受恐怖襲擊的最大可信泄漏孔徑為1.5 m。由此可見,當前在泄漏口大小的確定方面還存在一定的困難。

為定量分析泄漏口大小的不確定性對LNG船泄漏液池擴展的影響，以某138 000 m3LNG實體船為例，設定不同的泄漏口大小，結合已建立的液池擴展模型加以分析比較。該LNG實體船的貨艙為薄膜型貨艙，總艙容為36 650 m3，其剖面結構尺寸見圖1。假設該船因故發生泄漏事故，泄漏口面積為0.5～16.0 m2不等，泄漏孔中心與液面的初始距離為12.35 m，LNG的密度為442.5 kg/m3，液池蒸發速率為0.2 kg/(m2·s)，最小液池厚度為5.2 mm。應用MATLAB軟件對液池擴展模型進行編程計算，計算結果見圖2。

圖1 某薄膜型LNG船貨艙橫剖面尺寸

圖2 不同泄漏口面積的液池擴展情況比較

由圖2可知，泄漏口的面積不同，LNG液池的擴展情況(主要指液池最大半徑和存在時間)是不同的。液池擴展情況對不同泄漏口面積呈現出的敏感度也是不同的:當泄漏口面積<2 m2時，液池擴展對泄漏口大小的變化相當敏感;而當泄漏口面積>7 m2時，液池擴展隨著泄漏口大小的變化相對不明顯。例如:當泄漏口面積為16 m2時，液池最大半徑和存在時間分別為292.6 m及436 s;當泄漏口面積為7 m2時，液池最大半徑和存在時間分別為284.4 m及463.5 s;當泄漏口面積為0.5 m2時，液池最大半徑和存在時間分別為97.6 m及4 676 s;當泄漏口面積為2 m2時，液池最大半徑和存在時間分別為188.6 m及1 263 s。導致船舶泄漏事故發生的主要原因仍是碰撞和擱淺，這就意味著LNG船發生泄漏事故時，其泄漏口面積通常仍<2 m2，即液池擴展變化對泄漏口面積的變化很敏感。因此,在無法確定泄漏口大小的情況下，LNG船泄漏液池擴展計算有很大的不確定性。

2.2流量系數的不確定性

流量系數用來表征流體流過孔洞時受到的摩擦損失的大小。在對LNG船泄漏事故風險進行評價時，不同機構采用的流量系數不盡相同(例如：DNV和Sandia取0.6;FAY和ABS取1.0)[7]，因此計算得到的危害后果差別很大。實際上，LNG船泄漏事故中流量系數很難確定,這是因為：

1) 泄漏口的形狀往往難以確定。泄漏口的形狀不同，其流量系數的建議取值是不同的。表1為不同泄漏口形狀液體泄漏流量系數建議值。[8]

表1 不同泄漏口形狀液體泄漏流量系數建議值

2) LNG船的液貨在泄漏到外部環境之前會經過內薄膜、絕熱層、外薄膜和船舶的雙層殼等多層阻礙，流經距離至少2 m，正常達到3～4 m，整個流動過程十分復雜，流量系數難以確定。

3) 根據SOTERIOU等[9]的研究結論，LNG泄漏射流的流量系數主要由泄漏孔內部的流動狀態(單相流動、空穴流動和返流流動)來決定,而泄漏孔內部的流動狀態主要受泄漏孔直徑、長度、上游壓力、下游壓力、液體黏度、密度及飽和蒸汽壓等因素的影響。

圖3為射流流動狀態影響因素。在泄漏過程中，隨著LNG儲罐自由液面高度不斷下降，泄漏孔上游壓力不斷下降，泄漏孔內部的流動狀態將發生變化。此外,流量系數在整個泄漏過程也是不斷變化的。以上述LNG船為例，假設泄漏面積為2 m2，流量系數分別取0.4,0.6和1.0，其他條件不變，計算得到液池擴展情況比較見圖4。

圖3 射流流動狀態影響因素

圖4 不同流量系數的液池擴展情況比較

由圖4可知：流量系數的取值對液池擴展變化的影響較大。在其他條件相同的情況下，當流量系數取0.4時，液池的最大半徑和存在時間分別為123.3 m及2 965 s;而當流量系數取1.0時，液池最大半徑和存在時間分別為188.7 m及1 263 s。由此可見，在LNG船泄漏流量系數難以確定的情況下，液池擴展變化計算也存在著不確定性。

2.3蒸發速率的不確定性

LNG液池蒸發速率的大小主要取決于水與LNG液池間的熱通量。根據液體與下墊面間的溫度差，液體的沸騰形式有核沸騰(Nucleate Boiling)、過渡沸騰(Transition Boiling)和膜沸騰(Film Boiling)，沸騰的狀態不同，熱通量是不同的(見圖5)。研究[3]表明，泄漏入水的LNG在液池內大部分時間都處于膜沸騰狀態，只有當液池蒸發進入后期 (蒸發量達69%～84%)時才轉為過渡沸騰狀態(見圖5中的虛箭頭)。因此，許多蒸發模型經常將膜沸騰熱通量28 kW/m2作為LNG液池與水的熱通量，與之對應的蒸發速率為0.042 kg/(m2·s)。該蒸發速率是假定泄漏入水的LNG對水面沒有擾動的情況下得到的理論值，而實際上泄漏入水的LNG或多或少都會對水面造成擾動，LNG與水的攪合將使熱交換急劇增加，蒸發速率也隨之大大提高。表2[10]為不同泄放(攪合)方式下LNG蒸發速率試驗數據統計,可看出泄放的方式不同，液池蒸發速率有較大差別。然而,表2中的數據僅僅是通過小尺度試驗得到的，對于真正的LNG船泄漏事故，泄漏到水面的LNG速度可高達15 m/s,射流可至水面以下13 m處[11]，LNG與水攪合較為激烈，其蒸發速率更加難以確定。正因為如此，不同機構推薦的LNG水面蒸發率是不同的。例如:DNV使用的PHAST模型假設蒸發率為0.182 kg/(m2·s)；ABS假設蒸發率為0.072 kg/(m2·s)；QUEST軟件假設蒸發率為0.2 kg/(m2·s)[12]。圖6為泄漏面積為2 m2，流量系數為1.0，蒸發速率分別為0.042 kg/(m2·s)(理論值),0.072 kg/(m2·s),0.182 kg/(m2·s)和0.200 kg/(m2·s)時，LNG船泄漏液池擴展半徑計算結果比較。由圖6可知:按理論蒸發速率計算會過高地估算最大液池半徑和存在時間；采用不同機構不同的蒸發率取值，計算結果將有較大差別(例如,按ABS蒸發率計算出液池最大半徑和存在時間為293.7 m及1 532 s，而按 QUEST蒸發率計算出液池最大半徑和存在時間為188.7 m及1 263 s)。蒸發率不確定將直接導致液池擴展計算不確定。

圖5 典型的沸騰-熱通量曲線

表2 不同泄放方式下LNG蒸發速率試驗數據統計

圖6 不同蒸發速率的液池擴展情況比較

2.4貨艙壓力不確定分析

目前許多模型都假設在LNG船泄漏過程中作用于貨艙LNG液面上的壓力保持不變。實際上，LNG貨艙是相對封閉的空間，艙內壓力的變化主要取決于LNG泄漏體積率、真空破壞閥補氣率和LNG體積蒸發率三者間的相互作用。當LNG船正常卸貨或泄漏孔洞較小時，卸貨/泄漏體積率等于補氣閥補氣率與LNG體積蒸發率之和，艙內壓力基本上保持不變;但當泄漏孔洞較大時，尤其是在泄漏初期，LNG泄漏體積率明顯大于貨艙的真空破壞閥補氣率與貨艙內LNG的蒸發率之和，貨艙內作用于LNG液面的壓力顯著下降，出現較高的真空度，艙內壓力的下降反過來又會影響LNG泄漏速率和艙內LNG的蒸發率，泄漏率和蒸發率改變又將使艙內壓力發生新的變化，這種相互作用使得貨艙內的壓力變得難以確定。為說明LNG船泄漏過程中貨艙內壓力的變化情況，根據理想氣態方程和Hashemi-Wesson定理[13]建立LNG船泄漏時LNG泄漏體積率、真空破壞閥補氣率及貨艙內LNG體積蒸發率相互作用微分方程組(限于篇幅此處不再給出)，以前述薄膜型貨艙為例，設定泄漏孔徑為4 m，求解貨艙內壓力隨泄漏時間的變化過程，計算結果見圖7。

由圖7可知:在LNG泄漏過程中，艙內壓力并不是穩定不變的，而是會發生復雜變化，在泄漏剛開始時急劇下降，隨后逐漸回升;貨艙在整個泄漏過程中基本上都處于真空狀態。艙內壓力的這種變化將使得平均泄漏速率下降，液池的最大半徑變小，液池的持續時間變長。

圖7 艙內壓力隨泄漏時間變化過程

3 結束語

為對LNG船泄漏液池擴展變化進行計算，根據小孔射流伯努利方程、質量平衡及液體慣性擴展理論建立包括泄漏源強和液池半徑擴展在內的LNG船泄漏液池擴展模型。由于缺乏事故統計數據及事故尺寸的現場試驗數據，模型中有些重要參數尚存在一定的不確定性，影響LNG船泄漏液池擴展計算。主要的不確定參數有泄漏口大小、流量系數、液池蒸發速率和儲罐內壓力，其不確定性表現為：導致LNG船泄漏的事件不同(恐怖襲擊、碰撞或擱淺)，泄漏口大小可能在0～16 m2不等；流量系數受泄漏口形狀、LNG泄漏流動狀態和評估機構主觀因素等影響，取值范圍較廣(0.4～1.0)，計算結果可相差2.5倍；LNG在水面大部分時間內都處于膜沸騰狀態，當液池蒸發快結束時轉為過渡沸騰狀態，蒸發速率在液池存在時間內是變化的，LNG泄放入水的方式不同，與水的攪合激烈程度不同，蒸發速率會有很大不同(常取值0.042～0.25 kg/(m2·s))；LNG船貨艙是相對封閉的空間，當發生事故性泄漏時(尤其是泄漏口徑較大時)，儲罐內會產生一定的真空度，真空度的大小及其對泄漏速率的影響是很難確定的。此外，水面波浪、風向及風速等其他許多因素對LNG液池擴展的影響目前也很難確定。對LNG液池擴展變化的研究有待進一步深入進行。

[1] 趙敏,黃東風,佘孝云.從BP世界能源統計年鑒看中國能源發展[J].能源與環境,2014(6):16-17.

[2] 莊學強,鄭超瑜.液化天然氣船水上泄漏事故源強計算方法[J].大連海事大學學報,2010,36(4):125-128.

[3] HISSONG D W.Keys to Modeling LNG Spills on Water [J].Journal of Hazardous Materials,2007,140(3):465-477.

[4] KYTOMAA H,GAVELLI F.Studies of LNG Spills Over Water Point up Need for Improvement [J].Oil & Gas Journal,2005,103(18):61-65.

[5] HIGHTOWER M M,LUKETA-HANLIN A J,ATTAWAY S W.Breach and Safety Analysis of Spills Over Water from Large Liquefied Natural Gas Carriers [C]∥Listmail.Energy.Ca.Gov,2008.

[6] PITBLADO R M,BAIK J,HUGHES G J,et al.Consequences of LNG Marine Incidents[C].CCPS Conference Orlando,2004:1-19.

[7] BAIK J,RAGHUNATHAN V,WITLOX H.Consequence Modeling of LNG Marine Incidents[C].7th Professional Development Conference & Exhibition,2006:1-7.

[8] 王琳,季潔,賀寶龍,等.LNG地面泄漏蒸發速率的計算[J].天然氣與石油,2014,32(5):29-32.

[9] SOTERIOU C,ANDREWS R,SMITH M.Direct Injection Diesel Sprays and the Effect of Cavitation and Hydraulic Flip on Atomization[J].Simulation & Modeling,1995,120(16):84-93.

[10] WOODWARD J L,PITBLADO R M.LNG Risk Based Safety-Model & Consequence Analysis [M].John Wiley & Sons,Inc.,Hoboken,New Jersey,2010:134-145.

[11] DAHLSVEEN J,KRISTOFFERSON R,SAETRAN L R.Jet Mixing of Cryogen and Water[C].Second Symposium on Turbulence and Shear Flow Phenomena,Stockholm,Sweden,2001:329-334.

[12] 于桂峰,張彬,吳宛青.LNG船舶泄漏水面液池擴展的數值分析研究[J].大連海事大學學報,2008,34(S2):170-172.

[13] KHELIFI-TOUHAMI M S,BENBRIK A,LEMONNIER D,et al.Laminar Natural Convection Flow in a Cylindrical Cavity Application to the Storage of LNG[J].Journal of Petroleum Science and Engineering,2010,71(3-4):126-132.

PoolSpreadCalculationofLNGCarrierSpillsandUncertaintyAnalysis

ZHUANGXueqiang1,2，SUNDi1,2，GAOXiaohong3，LIGesheng3

(1.Marine Engineering Institute,Jimei University,Xiamen 361021,China; 2.Fujian Provincial Key Laboratory of Naval Architecture and Ocean Engineering,Xiamen 361021,China；3.School of Energy and Power Engineering,Wuhan University of Technology,Wuhan 430070,China)

A model for pool spread of LNG carrier spills is proposed based on the theory of orifice flow and inertia-gravity balance.The uncertainties of important parameters are analyzed,and the degree of the impacts of them on predicting results is calculated by means of the software MATLAB.The research reveals the following facts:Different kind of accidents (collision，grounding or terrorist attack)may result different LNG carrier leakage and the possible breach size can vary from 0 to 16 m2; The pool spread is sensitive to the breach size when it is less than 2 m2; The discharge coefficient,which is dependent on the breach shape and the liquid property,varies significantly,and the prediction deviation can be great,possibly up to 250%; The average evaporation rates also show wide variations (0.042～0.25 kg/(m2·s)) depending on the effects of the LNG-water turbulence and the boiling region.

waterway transportation; liquid pool; LNG carrier; release

U698.6;U674.13+3.3

A

2017-01-15

福建省科技重大平臺資助項目(2014H2001)；福建省科技(重點)項目(2017Y0065)；廈門市科技計劃項目(3502Z20173030)；廈門南方海洋研究中心項目(14GZB66NF30)；集美大學博士啟動基金(ZQ2013006)

莊學強(1974—),男,福建泉州人，副教授，碩士生導師，博士，從事水上危險品運輸安全與防污染研究。E-mail:13015917291@163.com

孫 迪(1979—)，女，遼寧沈陽人，講師，博士，主要研究方向為柴油機關鍵摩擦副摩擦監測技術。E-mail:sundi6329@sina.com

1000-4653(2017)02-0073-05