姚淑鵬，李玉星，王武昌，宋光春，姜 凱，施政灼

（中國石油大學（華東）山東省油氣儲運安全省級重點實驗室，山東 青島 266580）

作為一種清潔能源，天然氣水合物（簡稱水合物）礦藏受到了越來越多國家的關注，水合物礦藏的開發(fā)和開采成為人們研究的熱點問題［1］。水合物礦藏主要分布在極地和深水區(qū)域，其中，深水區(qū)域約占總儲量的95%［2-3］。目前，世界上試開采的多為成巖的水合物礦藏，而對于總量大卻埋深淺、膠結性差的非成巖水合物礦藏，我國周守為院士提出了固態(tài)流化綠色開采方法［4-5］。

研究固態(tài)流開采過程中管內(nèi)的多相流動，對了解豎直舉升管內(nèi)水合物流動特性有重要意義，但實驗難度大、成本高，不宜推廣。近年來，隨著計算流體力學（CFD）的發(fā)展，數(shù)值模擬成為替代實驗進行相關研究的有效工具。Fatnes［6］采用Ansys CFX軟件對水平管內(nèi)的水合物流動特性進行了模擬。Balakin等［7-8］通過CFD模擬，對CCl3F（R11）水合物在紊流流態(tài)純水體系中的流動沉積現(xiàn)象進行了研究，并通過求解該模型模擬了管內(nèi)水合物顆粒的粒徑分布。江國業(yè)等［9］采用正交實驗，發(fā)現(xiàn)速度對壓降的影響最大。宋光春等［10-11］通過構建群體相平衡模型對R11水合物的流動特性進行了模擬研究，發(fā)現(xiàn)流速、水合物體積分數(shù)及水合物顆粒在管道內(nèi)分布情況都會對水合物漿液的流動壓降產(chǎn)生影響。周守為等［12］通過有限差分迭代法對豎直井筒內(nèi)的水合物的分解臨界點等流動特性進行了研究，對深水淺層水合物礦藏的開采提出了指導性意見。劉艷軍等［13］對水合物漿體垂直管輸?shù)膸嘧兞鲃犹匦赃M行研究，發(fā)現(xiàn)水合物分解對漿體管道運輸具有減阻作用。上述研究沒考慮水合物顆粒在豎直管內(nèi)流動輸送過程中的聚集與破碎，與實際過程存在一定差距，因此值得進一步完善。

本工作引入基于水合物顆粒聚集動力學的群體平衡模型對水合物漿液在豎直管內(nèi)的流動特性進行了模擬，通過群體平衡模型，模擬了水合物顆粒在管道內(nèi)的聚集與破碎，較好地模擬了水合物漿液在豎直管內(nèi)的流動狀態(tài)，進而對不同平均流速、體積分數(shù)下的水合物漿液流動特性進行研究，為深水淺層水合物礦藏開采的流動安全保障和水合物漿液輸送技術提供理論支持。

1 模型建立

1.1 幾何模型

幾何模型的構建基于Balakin等開展的R11水合物漿液流動特性實驗所使用的管道，對它進行了三維模型的構建，管道長2.0 m，直徑4.52 cm。對構建好的幾何模型進行網(wǎng)格劃分，考慮到流體在邊界區(qū)域的物理參數(shù)變化相對比較劇烈，因此，對進口壁面進行了邊界層效應的處理，共設6層邊界層。對其他網(wǎng)格均以1 mm劃分。本次建模共有170 180個六面體網(wǎng)格，經(jīng)過網(wǎng)格獨立性檢驗，網(wǎng)格質(zhì)量為0.917。模型進口截面示意圖見圖1。

圖1 模型進口截面示意圖Fig.1 Schematic diagram of the inlet section of the model.

1.2 模型基本假設

建模過程中的基本假設有:水合物漿液在豎直管內(nèi)流動過程等溫；忽略水合物漿液內(nèi)相間質(zhì)量傳遞，即不考慮水合物的生成及分解；流體是不可壓縮介質(zhì)；水合物顆粒均為連續(xù)性介質(zhì)；不考慮水合物在管壁上的黏附。

1.3 多相流模型及湍流模型

多相流模型選用歐拉-歐拉雙流體模型，由連續(xù)性方程（見式（1））、動量方程（見式（2））以及用來封閉方程組的本構方程組成。

本次模擬的水合物漿液的流動為固液兩相流，因此，需要重點考慮液固的耦合。在Fluent 模擬過程中，液固耦合由相間動量交換實現(xiàn)。本次模擬在計算相間動量交換時，主要考慮相間拖曳力，它由摩擦阻力和形狀阻力組成。相間拖曳力采用Gidaspow 模型［13］計算。

當αs≤20%時，采用Wen-Yu 公式（見式（3））計算。

當αs＞20%時，采用Ergun 公式（見式（4））計算。

則相間拖曳力可以由式（5）計算。

本文模擬中，還需要對水合物顆粒的黏度進行求解，見式（6）［14］。

R11水合物漿液與天然氣水合物漿液的表面黏度分別由Roscoe-Brinkman方程（見式（7））［15-16］和 Thomas方程（見式（8））［17］計算。

根據(jù)式（6）～式（8）可以分別編制水合物顆粒黏度的用戶自定義函數(shù)（UDF），用以計算R11水合物和天然氣水合物顆粒相的黏度。

湍流模型則采用Fluent 15.0 軟件自帶的標準k-ε模型。

1.4 群體平衡模型

最初，群體平衡模型用來對人口數(shù)量的動態(tài)平衡進行描述。后來，Hulburt等［18］將群體平衡模型引入到化工領域中，通過跟蹤氣泡、液滴、顆粒等離散相的數(shù)量密度函數(shù)的變化，將離散相的成核、生長、聚集、破碎等微觀行為與它的體積、面積、粒徑等宏觀屬性聯(lián)系起來，這樣既保留了離散相的微觀屬性，又通過平均化的方式表征了連續(xù)相的運動，減少了計算量。

在群體平衡模型引入中，根據(jù)本文假設，可使用式（9）的群體平衡方程。

碰撞頻率方面，主要考慮由差速沉降與流動剪切造成的碰撞，并取兩者碰撞頻率之和作為水合物顆粒的實際碰撞頻率。其中，差速沉降碰撞頻率采用Camp等［19］提出的公式（見式（10））計算。

流動剪切碰撞頻率方面，當水合物顆粒小于Kolmogorov尺度時，其處于湍流耗散區(qū)，在這一區(qū)域內(nèi)水合物的聚集行為主要受渦內(nèi)局部剪切力的影響。此時，水合物顆粒的碰撞頻率可用Saffman等［20］提出的公式（見式（12））計算。

當水合物顆粒大于Kolmogorov尺度時，其處于湍流慣性區(qū)，在這一區(qū)域內(nèi)水合物被主流場牽引運動。此時，水合物顆粒的碰撞頻率可用Abrahamson［21］提出的公式（見式（13））計算。

采用曲線模型計算水合物顆粒間的聚并效率，由于連續(xù)相為水，水合物顆粒間不存在液橋力，故聚并效率計算時主要考慮范德華力與流動剪切力之比。聚并效率按式（14）［22］計算。

對水合物顆粒破碎效應的模擬選用Fluent 15.0軟件自帶的lehr模型。

根據(jù)式（10）～（14）編制UDF，并以此計算群體平衡模型中的關鍵參數(shù)。

1.5 模型求解

通過Fluent 15.0 軟件，對模型進行求解。進口設定為速度進口，出口為壓力出口，出口壓力為0，壁面無滑移?？紤]到水合物開采時，井筒內(nèi)的實際工況重力方向與水合物漿液的流動方向相反，所以本次模擬中設定重力方向與流動方向相反；本次模擬選取二階迎風格式；當各因子的殘差收斂到10-5時，認為已經(jīng)收斂，模擬結束。部分模擬參數(shù)見表1，模擬工況見表 2［23］。

表1 模型參數(shù)表Table 1 Model parameter

表2 模擬工況參數(shù)表Table 2 Simulation condition parameter

1.6 模型驗證

目前，由于豎直管內(nèi)的水合物漿液的流動實驗較少，相關實驗數(shù)據(jù)尚不完備。為驗證本次模型的準確性，依據(jù)Balakin等［7-8］開展的水平管內(nèi)R11水合物漿流動特性的實驗對模型進行驗證，主要以單位壓降作為驗證參數(shù)。模型驗證的工況見表3，單位壓降的驗證結果見表4。由表4可知，模型能夠較準確地模擬水合物漿液的流動狀態(tài)。

表3 模型驗證的工況參數(shù)Table 3 Condition parameters of model verification

表4 單位壓降的驗證結果Table 4 Verification results of unit pressure drop

2 流動特性研究

固態(tài)流開采過程中，采掘設備將水合物粉碎成細小顆粒后，將水合物顆粒與海水混合，然后通過密封管道輸送至海上平臺。水合物漿液在管道輸送的過程中，流經(jīng)一段較長的豎直舉升管。因此，主要研究了水合物漿液的不同舉升速度，以及與海水不同比例混合的水合物濃度對豎直管內(nèi)水合物漿液流動壓降、速度、濃度分布以及水合物粒徑分布等流動特性的影響。

2.1 豎直管內(nèi)水合物漿液顆粒粒徑分布

水合物漿液在豎直管內(nèi)的流動過程中，會存在水合物顆粒的聚集和破碎，從而導致水合物漿液在流動中產(chǎn)生的不同粒徑分布。圖2為不同入口流速、不同體積分數(shù)下出口截面x軸方向水合物顆粒的粒徑分布。

由圖2可知，9種工況下水合物顆粒粒徑在整個出口截面上均呈現(xiàn)對稱分布，同時，近壁面處水合物顆粒粒徑較大且粒徑變化梯度較大，管道中間部分水合物顆粒粒徑較小且粒徑較為均勻。根據(jù)管道內(nèi)水合物顆粒的聚集動力學可知，水合物顆粒發(fā)生碰撞聚集主要是流動剪切造成的。水合物漿液在管道近壁面處流動剪切最強，使得水合物顆粒在近壁面處發(fā)生碰撞聚集的概率最高，導致水合物顆粒的粒徑增大。而在管道中間部分，水合物漿液的流動剪切作用較弱且水合物顆粒分布較為均勻，所以水合物顆粒發(fā)生碰撞聚集的概率較小，粒徑也較小。

圖2 不同工況下的出口截面X軸方向水合物顆粒的粒徑分布Fig.2 Particle size distribution of hydrate particles in the outlet section X direction under different working conditions.

對比相同入口流速，不同水合物漿液濃度的水合物顆粒分布可知，在豎直管內(nèi)，水合物漿液濃度越高，近壁面處的顆粒粒徑越大，且粒徑分布梯度也越大，不利于輸送的安全性。這是由于水合物漿液濃度越高，水合物顆粒的數(shù)目越多，發(fā)生碰撞的概率隨之增大，在流動剪切較強的近壁面處，水合物顆粒發(fā)生劇烈的碰撞聚集，使水合物顆粒粒徑不斷增大。而在流動剪切作用較弱的管道中間部分，水合物顆粒粒徑碰撞聚集的概率較小，分布較為均勻。對比相同水合物漿液濃度、不同入口流速的漿液粒徑分布可知，在豎直管內(nèi)，水合物漿液流速越高，水合物顆粒的最大粒徑反而降低，這是因為雖然流速越高，流動剪切作用越強，水合物顆粒碰撞聚集的概率越大，但流動剪切越強，較大的水合物顆粒也越容易破碎，能維持的最大聚集粒徑也就越小，這對輸送的安全性有利。

2.2 豎直管內(nèi)水合物漿液濃度分布

不同入口流速、不同體積分數(shù)下出口截面水合物顆粒濃度分布云圖見圖3。

由圖3可知，9種工況的模擬結果都為均勻懸浮流。這是因為在豎直管內(nèi)，由于水合物漿液流動方向為逆重力方向，兩者可認為在一條直線上，所以沒有出現(xiàn)類似水平管中因重力影響而使?jié)舛仍谡麄€截面上分布不對稱的情況［25］，而是對稱分布的。同時，9種工況的濃度分布都呈現(xiàn)出近壁面處高，管道中間部分低的分布狀況。這是由于近壁面處水合物漿液流速較低，而管道中間部分流速較高所致，流速較高會使水合物顆粒的分散系數(shù)增大，從而導致管道中間部分濃度較低且分布較均勻，而四周因為流速降低而造成不均勻分布，出現(xiàn)水合物顆粒的堆積。

當流速相同時，隨著水合物漿液中水合物顆粒濃度的增加，水合物顆粒濃度分布依然保持管道中間部分低，近壁面處高的分布狀態(tài)，管道中間部分濃度均勻且濃度值在初始設置濃度附近，管內(nèi)濃度整體隨之增大。隨著漿液流速的增大，相同濃度下水合物漿液的濃度分布基本相同。這說明在豎直管道中，水合物漿液的流速對水合物漿液的濃度分布影響較小。

2.3 豎直管內(nèi)水合物漿液流動壓降

水合物漿液在豎直管內(nèi)的流動過程中，流動壓降是備受關注的重點問題之一，因此，對豎直管內(nèi)9種 工況組合的單位壓降進行了研究，水合物漿液流動壓降見表5。

由表5可知，當水合物濃度相同時，隨著水合物漿液的流速增加，單位壓降也隨之增大。同時，隨著水合物漿液流速的增大，單位壓降的增大梯度也越來越大。以體積分數(shù)為10%的水合物漿液為例，當入口平均流速從1.0 m/s增加到1.5 m/s（增大50%）和2.0 m/s（增大100%）時，單位壓降分別增加了104.9%，237.2%。當水合物流速相同時，隨著水合物漿液濃度的增加，單位壓降也隨之增大。同時，隨著水合物漿液濃度的增加，單位壓降的增大梯度也越來越大。以入口平均流速1.5 m/s為例，當水合物體積分數(shù)由10%增長到20%（增大100%）和30%（增大200%）時，單位壓降分別增大了13.5%，34.2%。

此外，通過對9組實驗數(shù)據(jù)的觀察可知，單位壓降的增大主要是由流速增大導致水合物漿液在豎直管中流動摩阻增大造成的，較高的流速不利于輸送的經(jīng)濟性。

圖3 不同工況下的出口截面濃度分布云圖Fig.3 Distribution cloud diagram of the outlet section concentrations under different working conditions.

表5 水合物漿液流動壓降Table 5 Flow pressure drop of hydrate slurry

3 結論

1）在豎直管道中，水合物顆粒粒徑在整個截面上是對稱分布的，且為近壁面處顆粒粒徑較大，管道中間部分顆粒粒徑較小。同時，水合物顆粒的最大粒徑隨水合物濃度的增加而增大，隨流速的增加而減小。

2）在豎直管道中，水合物漿液的濃度在整個截面上對稱分布，為均勻懸浮流。水合物濃度為近壁面處高，管道中間部分低。同時，隨著水合物濃度的增加，截面上水合物的濃度整體增加，且管道中間部分濃度均勻且濃度值在初始設置濃度附近，而流速對水合物濃度在橫截面上的分布影響不大。

符 號 說 明

a兩水合物顆粒發(fā)生碰撞后的聚并效率

b(L|L') 粒徑為L'的水合物顆粒破碎后產(chǎn)生粒徑為L水合物顆粒的概率

CD曳力系數(shù)

dp固相顆粒粒徑，m

拖曳力，kg·m/s

G絕對速度梯度，s-1

國民在畜禽產(chǎn)品消費成本方面普遍不高，這得益于國內(nèi)廉價的飼料與勞動力，但在推進畜牧業(yè)現(xiàn)代化發(fā)展中，這種觀念可能成為一種主觀障礙。為有效解決問題，需進一步加快物聯(lián)網(wǎng)技術在畜牧業(yè)發(fā)展中的應用。另外，在將RFID物聯(lián)網(wǎng)技術等應用服務于現(xiàn)代畜牧業(yè)發(fā)展中，不可避免地會出現(xiàn)一些行業(yè)類似的問題，如標準制定問題、成本控制問題以及信息安全問題等。但隨著現(xiàn)代技術的進步、政府的有效引導以及市場的不斷完善，此類問題必將迎刃而解。

g重力加速度，m/s2

H哈梅克常數(shù)

L，L'水合物顆粒粒徑，m

M相間動量交換項，kg/(m·s)2

n(L，t) 粒徑為L的水合物顆粒在t時刻的數(shù)量

密度，m-3

p壓力，Pa

R發(fā)生碰撞兩水合物顆粒的調(diào)和半徑，m

S(L') 粒徑為L'水合物顆粒的破碎頻率，s-1

t時間，s

u 速度矢量，m/s

u平均速度，m/s

V沉降速率，m/s

速度，m/s

相間的相對速度，m/s

z與流體性質(zhì)有關的常數(shù)

▽ 梯度

α體積分數(shù)

β相間曳力常數(shù)

β(L-L'，L') 粒徑分別為L-L'和L'的水合物顆粒的碰撞頻率，m3/s

βDS差速沉降碰撞頻率，s-1

ρ密度，kg/m3

τ應力張量，Pa

φ體積濃度

μ黏度，mPa·s

下角標

i，j水相或水合物顆粒相

l 液相

m 混合相

s 固相

TH 湍流耗散區(qū)

TG 湍流慣性區(qū)

參 考 文 獻

［1］ 宋廣喜，雷懷玉，王柏蒼，等. 國內(nèi)外天然氣水合物發(fā)展現(xiàn)狀與思考［J］. 國際石油經(jīng)濟，2013，21（11）:69-76.

［2］ Boswell R，Collett T S. Current perspectives on gas hydrate resources［J］. Energy Environ Sci，2011，4（4）:1206-1215.

［3］ 許紅，黃君權，夏斌，等. 最新國際天然氣水合物研究現(xiàn)狀及資源潛力評估（下）［J］. 天然氣工業(yè)，2005，25（6）:18-23.

［4］ 周守為，陳偉，李清平. 深水淺層天然氣水合物固態(tài)流化綠色開采技術［J］. 中國海上油氣，2014，26（5）:1-7.

［5］ 周守為，陳偉，李清平，等. 深水淺層非成巖天然氣水合物固態(tài)流化試采技術研究及進展［J］. 中國海上油氣，2017，29（4）:1-8.

［6］ Fatnes E D. Numerical simulations of the flow and plugging behaviour of hydrate particles［D］. Bergen:University of Bergen，2010.

［7］ Balakin B V，Hoffmann A C，Kosinski P. Experimental study and computational fluid dynamics modeling of deposition of hydrate particles in a pipeline with turbulent water flow［J］.Chem Eng Sci，2011，66（4）:755-765.

［8］ Balakin B V，Hoffmann A C，Kosinski P. Population balance model for nucleation，growth，aggregation，and breakage of hydrate particles in turbulent flow［J］. AIChE J，2010，56（8）:2052-2062.

［9］ 江國業(yè)，王曉婭，孫鵬. 基于正交試驗設計的水合物漿液流動特性數(shù)值模擬［J］. 科技導報，2014，32（13）:23-27.

［10］ 宋光春，李玉星，王武昌，等. 基于群體平衡理論的管內(nèi)水合物漿流動特性數(shù)值模擬［J］. 化工進展，2018，37（2）:561-568.

［11］ 宋光春，李玉星，王武昌，等. 基于群體平衡理論的水合物聚集動力學模型［J］. 化工進展，2018，37（1）:80-87.

［12］ 周守為，趙金洲，李清平，等. 全球首次海洋天然氣水合物固態(tài)流化試采工程參數(shù)優(yōu)化設計［J］. 天然氣工業(yè)，2017，37（9） :1-14.

［13］ 劉艷軍，唐孝蓉，胡坤. 天然氣水合物漿體分解對其在垂直管中流動特性影響的研究［J］. 化學通報，2018，81（3）:267-273.

［14］ Ding J，Gidaspow D. A bubbling fluidization model using kinetic theory of granular flow［J］. AIChE J，1990，36（4）:523-538.

［15］ 王繼紅. 冰漿的管道輸送熱流動特性［D］. 大連:大連理工大學，2013.

［16］ Pabst W. Fundamental considerations on suspension rheology［J］. P R Soc A，2004，48（1）:6-13.

［17］ 趙鵬飛，王武昌，李玉星，等. 管道內(nèi)水合物漿流動的數(shù)值模型［J］. 油氣儲運，2016，35（3）:272-277.

［18］ Hulburt H M，Katz S. Some problems in particle technology［J］. Chem Eng Sci，1964，19（8）:555-574.

［19］ Camp T R，Stein P C. Velocity gradients and internal work influid motion［J］. J Bsn Soc Civ Eng，1943，30（4）:219-237.

［20］ Saffman P G，Turner J S. On the collision of drops in turbulent clouds［J］. J Fluid Mech，1956，1（1）:16-30.

［21］ Abrahamson J. Collision rates of small particles in a vigorously turbulent fluid［J］. Chem Eng Sci，1975，30（11）:1371-1379.

［22］ van de Ven T G M，Mason S G. The microrheology of colloidal dispersions Ⅶ. Orthokinetic doublet formation of spheres［J］. Colloid Polym Sci，1977，255（5）:468-479.

［23］ Song Guangchun，Li Yuxing，Wang Wuchang，et al. Investigation of hydrate plugging in natural gas+diesel oil+water systems using a high-pressure flow loop［J］. Chem Eng Sci，2017，158（2）:480-489.

［24］ And X L，Logan B E. Collision frequencies between fractal aggregates and small particles in a turbulently sheared fluid［J］.Environ Sci Technol，1997，31（4）:1237-1242.

［25］ 陳鵬，劉福旺，李玉星，等. 水合物漿液流動特性數(shù)值模擬［J］. 油氣儲運，2014，33（2）:160-164.