摘要： 針對(duì)介觀尺度下蠟晶顆粒聚團(tuán)在管壁形成宏觀尺度下石蠟層的現(xiàn)象，采用過(guò)濾雙流體模型研究水平管道輸送伴有固態(tài)瀝青質(zhì)的稠油的液固兩相流流動(dòng)和傳熱特性.對(duì)比均勻結(jié)構(gòu)曳力模型和亞格子曳力模型的模擬結(jié)果，獲得了顆粒相時(shí)均體積分?jǐn)?shù)、顆粒相壓力和過(guò)濾傳熱系數(shù)的修正系數(shù)等參數(shù)分布規(guī)律.同時(shí)分析了壁面溫度對(duì)稠油溫度分布和傳熱的影響.結(jié)果表明：與均勻模型Gidaspow model相比，過(guò)濾模型能夠揭示蠟晶顆粒非均勻流動(dòng)結(jié)構(gòu)對(duì)液固兩相流的變化規(guī)律.Filtered model Ⅱ的結(jié)果更加接近于試驗(yàn)結(jié)果，F(xiàn)iltered model Ⅱ的結(jié)果隨著蠟晶顆粒體積分?jǐn)?shù)的增加而減小.當(dāng)蠟晶顆粒相體積分?jǐn)?shù)增加到12%時(shí)，加入壁面修正的亞格子模型的誤差僅為4.1%.而ZAMBRANO等人模擬結(jié)果與試驗(yàn)壓降的誤差為12.9%.亞格子尺度模型能夠詳細(xì)地再現(xiàn)瀝青質(zhì)顆粒聚團(tuán)的介觀尺度流動(dòng)行為，為稠油輸送過(guò)程中的蠟晶聚集以及流動(dòng)過(guò)程的預(yù)測(cè)提供了新的方法.

關(guān)鍵詞： 稠油；過(guò)濾雙流體模型；水平管道；液固兩相流；傳熱

中圖分類(lèi)號(hào)： O359" 文獻(xiàn)標(biāo)志碼： A" 文章編號(hào)： 1674-8530（2024）08-0785-09

DOI：10.3969/j.issn.1674-8530.22.0125

收稿日期： 2022-05-18； 修回日期： 2022-11-13； 網(wǎng)絡(luò)出版時(shí)間： 2024-07-19

網(wǎng)絡(luò)出版地址： https：//link.cnki.net/urlid/32.1814.th.20240716.1710.002

基金項(xiàng)目： 國(guó)家自然科學(xué)基金資助項(xiàng)目（51876032）；黑龍江省自然科學(xué)基金重點(diǎn)項(xiàng)目（ZD2019E002）

第一作者簡(jiǎn)介： 陳羽佳（1994—），黑龍江大慶人，女，博士研究生（chenyujia@stu.nepu.edu.cn），主要從事介觀尺度多相流數(shù)值模擬研究.

通信作者簡(jiǎn)介： 王淑彥（1974—），黑龍江大慶人，女，教授（wangshuyan@nepu.edu.cn），主要從事石油多相流模擬與試驗(yàn)研究.

陳羽佳，王淑彥，邵寶力，等. 過(guò)濾雙流體模型的水平管道內(nèi)稠油輸送數(shù)值模擬［J］. 排灌機(jī)械工程學(xué)報(bào)，2024，42（8）：785-793.

CHEN Yujia， WANG Shuyan， SHAO Baoli，et al. Numerical simulation of heavy oil transportation based on filtered two-fluid model in horizontal pipeline［J］. Journal of drainage and irrigation machinery engineering（JDIME），" 2024， 42（8）： 785-793. （in Chinese）

Numerical simulation of heavy oil transportation based

on filtered two-fluid model in horizontal pipeline

CHEN Yujia， WANG Shuyan*， SHAO Baoli， YUAN Zihan， XIE Lei， MA Yimei

（Key Laboratory of EOR， Ministry of Education， Northeast Petroleum University， Daqing， Heilongjiang 163318， China）

Abstract： Aiming at the phenomenon that wax crystal particles at mesoscopic scale aggregate on the pipe wall to form a macroscopic paraffin layer， the filtered two-fluid model was used to study the liquid-solid two-phase flow and heat transfer characteristics of heavy oil accompanied by solid asphaltene transported in horizontal pipelines. By comparing the simulation results obtained from the homogeneous drag model （Gidaspow model）， the Filtered model without wall correction （Filtered model Ⅰ）， and the Filtered model with wall correction （Filtered model Ⅱ）， the distribution laws of parameters such as the time-average asphaltic residues volume fraction， the solid pressures， and the correction to filtered interphase heat transfer coefficient along the radial direction were obtained. At the same time， the influence of pipe wall temperature on the temperature distribution and heat transfer coefficient of heavy oil was analyzed. The results show that compared with the Gidaspow model， the filtered model can show the changing law of the non-uniform flow structure of wax crystal particles on liquid-solid two-phase flow. It is revealed that the predicted pressure drops based on Filtered model Ⅱ corresponds well with the experimental data， and the Filtered model Ⅱ simulation results decrease as wax particle volume fraction increases. Obviously， when the volume fraction of wax particles increases to 12%， the error of the sub-lattice model with wall correction is only 4.1%， while the error between ZAMBRANO simulation and experimental pressure drop is 12.9%. The mesoscopic-scale flow behavior of particle clusters is detailedly reproduced with a filtered model considering wall correction， providing a new strategy for the prediction of wax aggregation and flow process during heavy oil transportation.

Key words： heavy oil;filtered two-fluid model;horizontal pipe;liquid-solid two-phase flow;heat transfer

稠油具有密度大，黏度高，含蠟量高等特點(diǎn)［1-2］.在管道運(yùn)輸過(guò)程中，壓力或環(huán)境溫度的變化會(huì)導(dǎo)致瀝青質(zhì)重組分通過(guò)相變以固體蠟晶顆粒形式從稠油中析出并沉積在管壁上［3-4］.此時(shí)，管道內(nèi)輸送稠油成了水平管道內(nèi)液固兩相流的流動(dòng).根據(jù)顆粒相在液相中的空間位置分布和流動(dòng)狀態(tài)，液固兩相流在水平管內(nèi)的均勻流動(dòng)特性的流型表現(xiàn)為均質(zhì)流（homogeneous flow）、非均勻流動(dòng)特性的流型表現(xiàn)為非均質(zhì)流動(dòng)（heterogeneous flow）、移動(dòng)床流動(dòng)（moving-bed flow）和穩(wěn)定床流動(dòng)（stationary-bed flow）.水平管道內(nèi)輸送稠油的蠟晶顆粒結(jié)晶聚集的過(guò)程與管道內(nèi)液固兩相流的流動(dòng)狀態(tài)相互關(guān)聯(lián)，瀝青質(zhì)顆粒是黏附性顆粒［5］，稠油在管道內(nèi)低速運(yùn)移過(guò)程中形成了瞬態(tài)介觀尺度顆粒聚團(tuán)，導(dǎo)致管道內(nèi)呈現(xiàn)出多尺度非均勻流動(dòng)結(jié)構(gòu)，進(jìn)一步增加了液固兩相流動(dòng)的復(fù)雜性.管道壁面上附著的瀝青質(zhì)顆粒聚團(tuán)會(huì)減小管道的輸送面積，降低輸送效率.因此，研究稠油及其瀝青質(zhì)顆粒在管道輸送中介觀尺度下的非均勻液固兩相流動(dòng)特性是十分必要的.

隨著計(jì)算機(jī)技術(shù)的迅速發(fā)展以及多相流理論的完善，計(jì)算流體力學(xué)（computational fluid dynamics， CFD）隨之得到快速發(fā)展，數(shù)值模擬方法已成為分析液固兩相流流動(dòng)中復(fù)雜的非均勻特性的重要方法［6］.主要分為歐拉-歐拉模型（雙流體模型）和歐拉-拉格朗日模型，其中雙流體模型將液體視為連續(xù)相，把顆粒相視為擬流體，并假定顆粒與液體在空間任何位置共存且相互滲透，更適用于模擬稠密顆粒流動(dòng)［7］.在液固兩相流系統(tǒng)中，顆粒相分布均勻是由于連續(xù)相液體對(duì)離散相顆粒的攜帶作用，然而顆粒相的非均勻分布使得液相流動(dòng)穿過(guò)介觀尺度下的顆粒聚團(tuán)，而不是簡(jiǎn)單的液相繞流，因而其液固相間的曳力作用機(jī)制發(fā)生了改變.曳力是液固兩相流動(dòng)中最主要的相間作用力，是用來(lái)表示液固兩相動(dòng)量交換非常關(guān)鍵的因素.目前在液固兩相數(shù)值模擬中，常用的基于顆粒相分布均勻體系的曳力模型有Gidaspow模型［8］、Syamlal-O′Brien模型［9］， Huilin-Gidaspow模型［10］等.高精度網(wǎng)格是雙流體模型準(zhǔn)確計(jì)算流動(dòng)特性的必要條件，瀝青質(zhì)顆粒聚團(tuán)析出形成蠟結(jié)晶，非均勻結(jié)構(gòu)的存在又會(huì)大大增加計(jì)算量［11］.因此，均勻曳力模型用來(lái)模擬顆粒聚團(tuán)的特性限制了雙流體方法數(shù)值模擬的準(zhǔn)確性［12］.為了模擬介觀尺度結(jié)構(gòu)，采用可靠的非均勻曳力模型是實(shí)現(xiàn)準(zhǔn)確數(shù)值模擬的關(guān)鍵.

IGCI等［13］提出了過(guò)濾雙流體模型（Filtered two-fluid model），它是一種適用于粗網(wǎng)格的計(jì)算模型，并且是關(guān)聯(lián)型多尺度數(shù)值模擬方法，其通過(guò)對(duì)給定的亞格子尺度下的雙流體模型的高分辨率模擬得到的計(jì)算結(jié)果采用濾波函數(shù)進(jìn)行過(guò)濾，從而推導(dǎo)出相應(yīng)的過(guò)濾雙流體模型方程的曳力系數(shù)和顆粒相應(yīng)力的本構(gòu)方程.濾波函數(shù)是將大尺度流動(dòng)和小尺度流動(dòng)過(guò)濾分離的過(guò)程，大尺度流動(dòng)通過(guò)普通的數(shù)值模擬方法就可以求解，小尺度流動(dòng)則需要建立與大尺度流動(dòng)的關(guān)系，尋求新的模型求解.因此，通過(guò)不同過(guò)濾長(zhǎng)度對(duì)高精度計(jì)算所得流場(chǎng)進(jìn)行時(shí)間和空間上的平均，對(duì)過(guò)濾區(qū)域內(nèi)的曳力系數(shù)、顆粒相壓力和顆粒相黏度進(jìn)行統(tǒng)計(jì)分析，建立相關(guān)物理量與過(guò)濾長(zhǎng)度的關(guān)系，從而得到過(guò)濾曳力系數(shù)和有效應(yīng)力新的本構(gòu)方程.過(guò)濾雙流體模型與雙流體模型的方程形式上完全相同，區(qū)別在于過(guò)濾后的曳力系數(shù)、顆粒相壓力和顆粒相黏度的本構(gòu)方程不同.過(guò)濾雙流體模型是采用周期性邊界條件，在流動(dòng)區(qū)域內(nèi)采用高精度求解方法而得到的計(jì)算模型，并沒(méi)有考慮對(duì)于壁面邊界流動(dòng)的影響.IGCI等［14］提出了考慮壁面修正的本構(gòu)方程，研究發(fā)現(xiàn)考慮了壁面修正的過(guò)濾雙流體模型預(yù)測(cè)結(jié)果更加接近試驗(yàn)數(shù)據(jù)，并且不同過(guò)濾長(zhǎng)度對(duì)應(yīng)的模擬結(jié)果幾乎相同.AGRAWAL等［15］在亞格子過(guò)濾模型中考慮了過(guò)濾相間傳熱系數(shù)，發(fā)現(xiàn)過(guò)濾后的相間傳熱系數(shù)比微觀的相間傳熱系數(shù)小2個(gè)數(shù)量級(jí).過(guò)濾相間傳熱系數(shù)的修正系數(shù)與過(guò)濾長(zhǎng)度和過(guò)濾顆粒相體積分?jǐn)?shù)相關(guān)，表明了小尺度流動(dòng)結(jié)構(gòu)對(duì)相間的傳熱有顯著的影響［16］.

水平管道內(nèi)輸送稠油是一種典型的非均勻液固兩相流，目前對(duì)于稠油中瀝青質(zhì)顆粒團(tuán)聚現(xiàn)象的研究較少.文中基于SUNDARESAN和IGCI提出通過(guò)過(guò)濾雙流體模型對(duì)水平管道內(nèi)稠油輸送液固兩相流動(dòng)的非均勻特性和傳熱特性進(jìn)行研究，數(shù)值模擬不同蠟晶顆粒體積分?jǐn)?shù)條件下管道內(nèi)蠟晶聚集特性，以及壁面修正對(duì)蠟晶顆粒流動(dòng)和沉積的影響；對(duì)比均勻曳力模型與非均勻過(guò)濾曳力模型對(duì)顆粒體積分?jǐn)?shù)、顆粒相擬溫度、傳熱系數(shù)等參數(shù)分布的影響；同時(shí)分析管壁溫度與稠油溫度差的變化對(duì)稠油析出蠟晶顆粒過(guò)程的溫度分布和傳熱特性的影響.

1" 過(guò)濾雙流體模型

通過(guò)對(duì)亞格子尺度下的雙流體模型方程進(jìn)行空間平均，得到過(guò)濾后的雙流體模型方程.過(guò)濾后的小尺度流動(dòng)結(jié)構(gòu)，其長(zhǎng)度尺度小于過(guò)濾長(zhǎng)度，由高分辨率下的過(guò)濾結(jié)果構(gòu)成的殘差項(xiàng)進(jìn)行計(jì)算［17］.此外，過(guò)濾曳力系數(shù)、顆粒相壓力和顆粒相黏度的相關(guān)性也可從高分辨率下的過(guò)濾雙流體模型的模擬中得到.根據(jù)過(guò)濾長(zhǎng)度，對(duì)微觀雙流體模型方程進(jìn)行網(wǎng)格平均，建立過(guò)濾后的雙流體模型方程［18］.過(guò)濾模型的結(jié)果與在同一系統(tǒng)的細(xì)網(wǎng)格模擬中觀察到的宏觀行為相似.

Eulerian-Eulerian模型包括液相和顆粒相的質(zhì)量和動(dòng)量平衡.過(guò)濾后的液相和顆粒相連續(xù)性平衡方程可描述為

t（φl(shuí)ρl）+

SymbolQC@ ·（φl(shuí)ρlvl）=0，（1）

t（φsρs）+

SymbolQC@ ·（φsρsvs）=0，（2）

式中：vs和vl分別為顆粒相和液相的過(guò)濾速度;φs和 φl(shuí)分別為顆粒相和液相的過(guò)濾體積分?jǐn)?shù);ρs和ρl分別為顆粒相密度和液相密度；

SymbolQC@ 是梯度算子（在各個(gè)方向上都有微分）.

液相和顆粒相動(dòng)量守恒方程為

t（φl(shuí)ρlvl）+

SymbolQC@ ·（φl(shuí)ρlvlvl）=

－φl(shuí)

SymbolQC@ p+φl(shuí)ρlg－βls（vl－vs）+

SymbolQC@ ·τl，（3）

t（φsρsvs）+

SymbolQC@ ·（φsρsvsvs）=－φs

SymbolQC@ p+

SymbolQC@ ·ps+

SymbolQC@ ·τs+φsρsg+βls（vl－vs），（4）

式中：βls為過(guò)濾曳力系數(shù);τs和τl分別為過(guò)濾后的顆粒相和液相的剪切應(yīng)力;g為重力加速度;ρs為壓力；

SymbolQC@ p表征壓力分別在x，y，z方向上的微分.

過(guò)濾后的液相和顆粒相能量平衡方程可描述為［15］

ρlCpl（φl(shuí)Tl）t+

SymbolQC@ ·（φl(shuí)vlTl）=

SymbolQC@ ·（klφl(shuí)

SymbolQC@ Tl）－γ（Ts－Tl），（5）

ρsCps（φsTs）t+

SymbolQC@ ·（φsvsTs）=

SymbolQC@ ·（ksφs

SymbolQC@ Ts）+γ（Ts－Tl），（6）

式中：Cpl和Cps分別為液相和顆粒相的比熱容;Tl和Ts分別為過(guò)濾后的液相和顆粒相溫度，量綱為一;γ為相間傳熱系數(shù);kl和ks分別為液相和顆粒相的傳熱系數(shù)；Tl和Ts分別為液相和顆粒相的溫度.

基于計(jì)算網(wǎng)格尺度，Gidaspow的曳力模型僅適用于均勻液固兩相流動(dòng).然而，在實(shí)際的液固兩相流動(dòng)中，顆粒相并非均勻分布于液相中，顆粒聚集形成顆粒聚團(tuán)形成密相，從而會(huì)形成不均勻的液固兩相結(jié)構(gòu)，這種非均勻結(jié)構(gòu)的存在會(huì)對(duì)流場(chǎng)演變產(chǎn)生重要影響［18］.學(xué)者們?cè)噲D在不改變Eulerian-Eulerian模型和均勻曳力模型的前提下，通過(guò)網(wǎng)格細(xì)化來(lái)解析原來(lái)粗網(wǎng)格內(nèi)無(wú)法考慮的介觀尺度結(jié)構(gòu)，但會(huì)存在結(jié)果誤差大和計(jì)算量大的問(wèn)題，也未從理論上解析細(xì)化網(wǎng)格尺度與介觀尺度結(jié)構(gòu)的關(guān)聯(lián)性［19］.引入湍流效應(yīng)或顆粒之間的黏性作用模型都不能從根本上解決這個(gè)問(wèn)題.因此從介觀尺度結(jié)構(gòu)的角度，研究非均勻曳力模型與非均勻結(jié)構(gòu)參數(shù)的關(guān)聯(lián)性［20］.

過(guò)濾雙流體模型通過(guò)過(guò)濾長(zhǎng)度過(guò)濾得到的新的本構(gòu)方程的曳力系數(shù)，是典型的非均勻曳力模型，也被稱(chēng)為亞格子過(guò)濾曳力模型（Filtered model Ⅰ）.并推導(dǎo)出與顆粒體積分?jǐn)?shù)和過(guò)濾長(zhǎng)度有關(guān)的方程［14］

βls=34CDρl（1－φs）φsvl－vsdp（1－φs）－2.65，（7）

β=βls（1+C），（8）

C=－Fr－1.6filterFr－1.6filter+0.4h3D（φs），（9）

Fr－1filter=gΔfilter/ν2t，（10）

式中：CD為曳力系數(shù)；dp為顆粒直徑；C為過(guò)濾曳力系數(shù)修正因子;νt為顆粒黏度系數(shù);Δfilter為過(guò)濾長(zhǎng)度;Frfilter為過(guò)濾長(zhǎng)度的Froud數(shù)，其中

CD=

24Re（1+0.15Re0.687），Relt;1 000，

0.44，Re≥1 000，

h2D（φs）=

2.7φs0.234，φslt;0.001 2，

－0.019φs－0.455+0.963，0.001 2≤φslt;0.014，

0.868exp（－0.38φs）－

0.176exp（－119.2φs），

0.014≤φslt;0.25，

－4.59×10－5·exp（19.75φs）+

0.852exp（－0.268φs），

0.25≤φslt;0.455，

（φs－0.59）（－1 501φs3+

2 203φs2－1 054φs+162），

0.455≤φs≤0.59，

0，0.59lt;φs≤0.65，（11）

h3D（φs）=h2D（φs）gg（φs），（12）

gg（φs）=φs0.24［1.48+exp（－18φs）］，φslt;0.18，

1，φs≥0.18，（13）

式中：h2D（φs）和h3D（φs）分別為2D和3D模型的過(guò)濾曳力系數(shù)修正因子;gg（φs）為由2D系數(shù)推導(dǎo)為3D模型過(guò)濾曳力系數(shù)的附加因子.

顆粒與壁面相互作用導(dǎo)致壁面區(qū)域的顆粒流動(dòng)的差異，需要通過(guò)加入壁面修正系數(shù)來(lái)考慮壁面效應(yīng)對(duì)曳力系數(shù)、顆粒壓力和顆粒黏度的影響.壁面修正系數(shù)不僅和顆粒體積份額有關(guān)，還與壁面的距離有關(guān).并且，推導(dǎo)出加入壁面修正的亞格子過(guò)濾曳力模型（Filtered model Ⅱ）［14］為

rd=g（D/2－r）/νt2，（14）

βfiltered，scaled=βfiltered（φs，rd）βfiltered，core（φs）=

11+5exp（－0.75rd），0≤rd≤gD2νt2，（15）

βfiltered，core（φs）=βfiltered，periodic（φs），（16）

式中：r為半徑；βflitered，scaled為區(qū)域過(guò)濾曳力系數(shù);βfiltered為過(guò)濾曳力系數(shù)；rd為距離壁面的距離，量綱為一;βfiltered，core為3D模型中心區(qū)域過(guò)濾曳力系數(shù);βfiltered，periodic為周期性邊界條件下的過(guò)濾曳力系數(shù).

通過(guò)能量方程（5）和（6）能夠獲得相間傳熱系數(shù)的空間平均值，進(jìn)而可定義過(guò)濾相間傳熱系數(shù)為

γfilt=γ（Ts－Tl）（Ts－Tl）.（17）

對(duì)于不能被求解的小尺度流動(dòng)中的非均勻結(jié)構(gòu)會(huì)影響過(guò)濾相間傳熱系數(shù)的準(zhǔn)確預(yù)測(cè)，因此過(guò)濾傳熱系數(shù)的修正系數(shù)Q可定義為

Q=1－γfiltγ（φs，vs－vl）.（18）

并且，過(guò)濾相間傳熱系數(shù)的修正系數(shù)Q有關(guān)過(guò)濾顆粒相體積分?jǐn)?shù)和過(guò)濾長(zhǎng)度的擬合公式為

Q=Q1（φs）Q2（Δfilter），（19）

Q1（φs）=

0.995 7φsφs+0.001 4，

φs≤0.24;

0.984 6φs2－0.973 9φs+0.240 8φs2－0.986 8φs+0.243 7，

0.24lt;φs≤0.49;

0，" φsgt;0.49；（20）

Q2（Δfilter）=

Δfilter3－2.068 3Δfilter2+0.464 2ΔfilterΔfilter3－2.043 6Δfilter2+0.451 8Δfilter－0.039，

Δfiltergt;0.257;

0，" Δfilter≤0.257.（21）

2" 計(jì)算參數(shù)及邊界條件

文中計(jì)算模型參照了ZAMBRANO等人實(shí)驗(yàn)室中環(huán)路裝置設(shè)備，試驗(yàn)對(duì)表征水平管道輸送稠油的液固兩相流中的壓降進(jìn)行了試驗(yàn)研究［21］.試驗(yàn)裝置由2根直徑為0.024 3 m、長(zhǎng)度為3.0 m的直管組成，直管左側(cè)由U形管連接，直管右側(cè)連接用于制備稠油漿體的20 L攪拌槽，攪拌槽內(nèi)裝有軸向攪拌器用來(lái)制備均勻的漿體.攪拌罐通過(guò)齒輪泵（GP）將稠油泵入管道回路中來(lái)控制流量.采用壓力傳感器（PT）和溫度傳感器（TT）測(cè)量在線(xiàn)壓力和溫度.壓降是由壓差傳感器（DPT）進(jìn)行測(cè)量.本章采用的試驗(yàn)壓降數(shù)據(jù)是對(duì)2.5 m的直管測(cè)量得到的，試驗(yàn)系統(tǒng)圖如圖1所示.

計(jì)算區(qū)域?yàn)橹睆紻=0.024 3 m，數(shù)值管段長(zhǎng)度L=2.5 m+LE，橫截面為圓形的水平直管，LE為管道入口充分發(fā)展的穩(wěn)態(tài)層流長(zhǎng)度，文中選取LE=0.19 m. ZAMBRANO等試驗(yàn)發(fā)現(xiàn)，在中等流速和射流濃度下，稠油形態(tài)不受瀝青質(zhì)固體顆粒逐漸溶解的影響，成為液相，因此將瀝青質(zhì)作為顆粒相，脫瀝青質(zhì)油作為液相.液相和顆粒相的密度分別為906.59 kg/m3和1 120 kg/m3，顆粒直徑為0.5 mm，以上物性參數(shù)設(shè)置均根據(jù)ZAMBRANO等試驗(yàn)設(shè)置［21］.管道入口選擇速度入口，出口選擇壓力出口，出口設(shè)置為大氣壓力.

基于Eulerian-Eulerian方法進(jìn)行瞬態(tài)計(jì)算，液相湍流模型采用k-ε模型.顆粒相碰撞恢復(fù)系數(shù)e表征顆粒之間非彈性碰撞的程度，e值在0～1變化，當(dāng)顆粒之間為完全彈性碰撞時(shí)，e=1；當(dāng)顆粒之間為非完全彈性碰撞時(shí)，e=0.e小于1意味著顆粒因非彈性碰撞而損失能量，文中選用0.9.在Eulerian-Eulerian模型中，液相和顆粒相均被視為相互貫穿的連續(xù)相，顆粒相被處理為擬流體，采用無(wú)滑移邊界條件對(duì)顆粒相與壁面間的相對(duì)運(yùn)動(dòng)進(jìn)行模擬計(jì)算，沒(méi)有單獨(dú)考慮顆粒與管壁間的摩擦.在顆粒相的模擬計(jì)算中，采用Schaeffer模型描述了顆粒間的摩擦應(yīng)力.采用高精度的二階迎風(fēng)差分格式對(duì)動(dòng)量方程中的對(duì)流項(xiàng)進(jìn)行離散.采用的壓力-速度耦合的半隱式方法是SIMPLE算法.時(shí)間步長(zhǎng)為1.0×10-4 s，計(jì)算物理時(shí)間為80 s，取70～80 s的數(shù)據(jù)作為時(shí)間平均的計(jì)算樣本，表1為材料性能及數(shù)值模擬參數(shù)，表中dp為顆粒直徑，vin為入口流動(dòng)速度，μl為液相黏度，T為溫度.

文中采用計(jì)算流體力學(xué)（CFD）軟件Fluent17.0進(jìn)行求解，基于C語(yǔ)言的自定義函數(shù)（UDF）亞格子過(guò)濾曳力模型導(dǎo)入Fluent17.0進(jìn)行求解.

采用加入壁面修正的過(guò)濾曳力模型（Filtered model Ⅱ）進(jìn)行了網(wǎng)格分辨率在數(shù)值模擬計(jì)算中獨(dú)立性的驗(yàn)證，分別生成3組不同尺寸的均勻網(wǎng)格以檢驗(yàn)網(wǎng)格獨(dú)立性.粗網(wǎng)格、中等網(wǎng)格、細(xì)網(wǎng)格的網(wǎng)格節(jié)點(diǎn)數(shù)分別為93 450，170 702，434 520，對(duì)應(yīng)的網(wǎng)格尺寸分別為2.5，2.0，1.5 mm.圖2為基于加入壁面修正的過(guò)濾曳力模型下3種不同網(wǎng)格尺寸的顆粒體積分?jǐn)?shù)的變化分布.圖中R為徑向距離；φs0為初始顆粒體積分?jǐn)?shù).由圖可見(jiàn)，在管道底部，粗網(wǎng)格的顆粒體積分?jǐn)?shù)略高于其他網(wǎng)格.然而，中等網(wǎng)格和細(xì)網(wǎng)格的模擬結(jié)果分辨率非常相似.在細(xì)網(wǎng)格的網(wǎng)格尺度下，發(fā)現(xiàn)了該模型的網(wǎng)格獨(dú)立性.因此，在考慮計(jì)算成本和滿(mǎn)足計(jì)算條件的基礎(chǔ)上，本章選用中網(wǎng)格模型進(jìn)行數(shù)值模擬.

3" 數(shù)值模擬結(jié)果與討論

3.1" 水平管道輸送稠油的流動(dòng)特性

圖3為不同數(shù)值模擬模型所得壓降p結(jié)果與試驗(yàn)壓降數(shù)據(jù)對(duì)比分布.初始顆粒體積分?jǐn)?shù)分別為6%，8%和12%，文中應(yīng)用的試驗(yàn)壓降數(shù)據(jù)來(lái)自ZAMBRANO等［21］的文章，并且ZAMBRANO等采用代數(shù)滑移模型（algebraic slip mixture，ASM）均勻模型預(yù)測(cè)得到的水平管道內(nèi)輸送稠油的壓降與試驗(yàn)數(shù)據(jù)進(jìn)行了對(duì)比.ASM model是一種簡(jiǎn)化的兩相流模型，與雙流體模型相比變量較少，將固液兩相看作混合相，不考慮相間作用力.文中基于過(guò)濾雙流體模型，采用3種曳力模型：Gidaspow模型、亞格子過(guò)濾曳力模型（Filtered model Ⅰ）和加入壁面修正的亞格子過(guò)濾曳力模型（Filtered model Ⅱ）.由圖可得，壓降都隨著顆粒相體積分?jǐn)?shù)的增加而增大.對(duì)比不同曳力模型的預(yù)測(cè)結(jié)果，Gidaspow模型與Filtered model Ⅱ結(jié)果比較接近，并且與試驗(yàn)數(shù)據(jù)吻合較好.此外，3種曳力模型的模擬壓降結(jié)果隨初始顆粒體積分?jǐn)?shù)的增加而增加.在相同條件下，F(xiàn)iltered model Ⅰ與試驗(yàn)壓降的誤差在3種曳力模型中是最大的.然而，當(dāng)顆粒體積分?jǐn)?shù)分別為6%和8%時(shí)，Gidaspow模型的誤差小于Filtered model Ⅱ.當(dāng)顆粒體積分?jǐn)?shù)增加到12%時(shí)，F(xiàn)iltered model Ⅱ的誤差則是最小的，僅為4.1%，說(shuō)明考慮了壁面修正的Filtered model Ⅱ在稠密顆粒流動(dòng)中表現(xiàn)較好.與此同時(shí)，與ZAMBRANO等的模擬結(jié)果相比，對(duì)于均勻曳力模型，ASM model與 Gidaspow model的模擬結(jié)果相比，Gidaspow model與試驗(yàn)數(shù)據(jù)誤差高于ASM model的誤差，ASM model的結(jié)果與試驗(yàn)數(shù)據(jù)吻合較好.然而，對(duì)于非均勻曳力模型，F(xiàn)iltered model Ⅱ的結(jié)果更加接近于試驗(yàn)結(jié)果，F(xiàn)iltered model Ⅱ的結(jié)果隨著顆粒體積分?jǐn)?shù)的增加而減小，與ZAMBRANO等模擬結(jié)果的趨勢(shì)則相反，其模擬結(jié)果與試驗(yàn)壓降的誤差為12.9%.這充分證明了過(guò)濾模型在預(yù)測(cè)高顆粒體積分?jǐn)?shù)流動(dòng)方面優(yōu)于均勻模型ASM模型.在稠油及其瀝青渣油流動(dòng)特性顯示為典型的多尺度流動(dòng)過(guò)程，應(yīng)推薦采用過(guò)濾模型.

圖4為距離管道入口1 m處管道橫截面的顆粒體積分?jǐn)?shù)分布.其中，圖4a為3種初始顆粒體積分?jǐn)?shù)φs0采用Filtered model Ⅱ模型所得顆粒相體積分?jǐn)?shù)徑向分布.由圖可見(jiàn)，隨著初始顆粒相體積分?jǐn)?shù)的增加，管道底部顆粒相體積分?jǐn)?shù)逐漸增大.

圖4b為不同曳力模型下距離管道入口1 m處時(shí)均顆粒體積分?jǐn)?shù)徑向分布.由圖可見(jiàn)，不同曳力模型下的顆粒體積分?jǐn)?shù)都是由管道頂部到管道底部逐漸增大，說(shuō)明了水平管道中稠油輸送過(guò)程中瀝青質(zhì)顆粒的堆積特性.3種曳力模型中，F(xiàn)iltered model Ⅱ的顆粒體積分?jǐn)?shù)最大，尤其是在靠近管道底部位置， Gidaspow模型的顆粒體積分?jǐn)?shù)最小，這充分說(shuō)明了壁面效應(yīng)導(dǎo)致瀝青質(zhì)顆粒在管道底部的聚集.不論是由于重力作用或是Filtered model Ⅱ模型中的壁面效應(yīng)，直觀的現(xiàn)象都是瀝青質(zhì)顆粒聚集在管道底部.由此推斷，Gidaspow模型適用于低顆粒體積分?jǐn)?shù)的均勻液固兩相流動(dòng).但隨著初始顆粒體積分?jǐn)?shù)的增加，形成了介觀尺度流動(dòng)結(jié)構(gòu)即顆粒聚團(tuán).在這種情況下，過(guò)濾曳力模型將顯示出對(duì)物理本質(zhì)更好的預(yù)測(cè).

圖5為在不同曳力模型下沿管道水平流動(dòng)方向顆粒相體積分?jǐn)?shù)的云圖分布.由圖可知， Gidaspow模型的顆粒相分布比過(guò)濾模型的顆粒相分布更均勻.在管道底部，未加入壁面修正的Filtered model Ⅰ的顆粒沉積量相對(duì)小于Filtered model Ⅱ. Filtered model Ⅱ預(yù)測(cè)的顆粒沉降現(xiàn)象最為明顯.由于壁面效應(yīng)和重力作用，在管道壁面附近的顆粒相體積分?jǐn)?shù)較高，并且有顆粒聚集的趨勢(shì).

分析表明，為了提高曳力模型預(yù)測(cè)的精確度，需要將壁面效應(yīng)引入過(guò)濾模型中，特別是對(duì)于稠油輸送的預(yù)測(cè).基于上述結(jié)果，考慮壁面效應(yīng)會(huì)影響顆粒相的分布，壁面修正對(duì)于捕捉介觀瀝青質(zhì)顆粒團(tuán)簇特征具有重要意義，對(duì)于預(yù)測(cè)管道結(jié)蠟厚度的現(xiàn)象更為直觀.

3.2" 蠟晶顆粒聚團(tuán)的流動(dòng)特性

在顆粒流的動(dòng)力學(xué)理論中，稠密顆粒流的顆粒碰撞產(chǎn)生隨機(jī)脈動(dòng).提出用顆粒擬溫度來(lái)表征隨機(jī)脈動(dòng)的能量和強(qiáng)弱，建立顆粒微觀運(yùn)動(dòng)與宏觀性質(zhì)的關(guān)系.顆粒擬溫度的表達(dá)式為θ=v′2/3，其中v′為顆粒脈動(dòng)速度［22］.圖6為3種曳力模型的顆粒擬溫度隨顆粒體積分?jǐn)?shù)變化的分布.

顆粒擬溫度隨顆粒相體積分?jǐn)?shù)的增加而升高.在3組預(yù)測(cè)值中，F(xiàn)iltered model Ⅱ的顆粒擬溫度的平均值最高，Gidaspow模型的顆粒擬溫度的平均值最低.顆粒擬溫度可以用來(lái)表征顆粒脈動(dòng)強(qiáng)度.此外，在顆粒相體積分?jǐn)?shù)較高時(shí)，主要是由顆粒間的碰撞作用引起顆粒脈動(dòng).隨著顆粒體積分?jǐn)?shù)的增加，顆粒間的碰撞增強(qiáng)，產(chǎn)生強(qiáng)烈的顆粒脈動(dòng)，顆粒擬溫度隨之升高.過(guò)濾曳力模型考慮了介尺度結(jié)構(gòu)對(duì)顆粒相分布不均勻的影響，導(dǎo)致過(guò)濾模型顆粒間碰撞比均勻曳力模型Gidaspow模型的顆粒間碰撞更加強(qiáng)烈.考慮了壁面效應(yīng)的Filtered model Ⅱ與Filtered model Ⅰ相比，由于考慮了顆粒與壁面間的碰撞，導(dǎo)致顆粒間和顆粒聚團(tuán)間的碰撞更加劇烈.

圖7為3種曳力模型的顆粒相壓力隨顆粒相體積分?jǐn)?shù)的變化分布.顆粒相應(yīng)力由顆粒相正應(yīng)力和顆粒相剪切應(yīng)力組成，顆粒相壓力代表瀝青質(zhì)顆粒的正壓力.顆粒相壓力是由顆粒的脈動(dòng)速度與顆粒間的碰撞引起的.在稠密顆粒流中，顆粒碰撞產(chǎn)生的正應(yīng)力傳遞了碰撞沖擊力.顆粒相壓力的表達(dá)式為ps=φsρsθ+2ρs（1+e）φ2sg0θ，其中g(shù)0為徑向分布函數(shù).過(guò)濾曳力模型計(jì)算的顆粒相壓力隨顆粒相體積分?jǐn)?shù)的增加先增大后減小，而Gidaspow模型的計(jì)算結(jié)果僅隨顆粒相體積分?jǐn)?shù)的增加而增大.Gidaspow模型和過(guò)濾曳力模型的顆粒相體積分?jǐn)?shù)分別在5%～15%和5%～60%.結(jié)果表明，過(guò)濾曳力模型能夠模擬出介觀尺度流動(dòng)結(jié)構(gòu)強(qiáng)化了顆粒間的相互作用力.

圖8為3種曳力模型的顆粒相黏度隨顆粒相體積分?jǐn)?shù)的變化分布.顆粒相黏度是顆粒碰撞的切向力引起的剪切應(yīng)力.并且，剪切應(yīng)力是瀝青質(zhì)顆粒的平動(dòng)和碰撞的動(dòng)量交換，其表達(dá)式為μs=45φ2sρsdpg0（1+e）θπ+10ρsdsπθ96（1+e）φsg01+45g0φs（1+e）2.剪切應(yīng)力的變化趨勢(shì)與顆粒相壓力相同.與顆粒相壓力相比，剪切應(yīng)力較小，說(shuō)明了考慮了壁面效應(yīng)，沉積在管壁底部的瀝青質(zhì)顆粒在管道稠油流動(dòng)中平衡了剪切應(yīng)力.

3.3" 水平管道內(nèi)輸送稠油的傳熱特性

圖9為不同管壁溫度下水平管道方向截面稠油溫度T的云圖分布，截面對(duì)應(yīng)的水平位置分別為0.5，1.0，1.5，2.0 m，其中管壁溫度分別為303.15，293.15，283.15 K，稠油初始溫度為328.15 K.

由圖9可知，3種管壁溫度下，在管道入口流動(dòng)過(guò)程中，接近管壁位置的溫度開(kāi)始減小，并且溫度降低的趨勢(shì)逐漸向管道中心區(qū)域延伸.管壁溫度低于稠油溫度，管壁與稠油間的對(duì)流換熱強(qiáng)烈，靠近壁面的稠油溫度受到冷卻，導(dǎo)致壁面向油流中心產(chǎn)生了徑向溫度梯度.隨著稠油在水平管道內(nèi)的流動(dòng)，傳熱方式由徑向?qū)醾鬟f到軸向?qū)α鲹Q熱的熱量交換.壁面溫度與稠油溫度不存在溫度差時(shí)，壁面內(nèi)周向的稠油處于均勻熱流密度的流動(dòng)狀態(tài).當(dāng)壁面溫度低于稠油溫度時(shí)，壁面附近均勻熱流的邊界層開(kāi)始不穩(wěn)定，使得壁面溫度與周向稠油溫度的熱量發(fā)生傳遞，進(jìn)而出現(xiàn)稠油溫度邊緣分布不規(guī)則的震蕩現(xiàn)象.溫度差導(dǎo)致的溫降和流動(dòng)過(guò)程中的熱量損失，使得越靠近管道出口處的管壁周?chē)碛蜏囟仍叫?產(chǎn)生溫降區(qū)域開(kāi)始是分散在管壁上的，進(jìn)而管壁周?chē)鷧^(qū)域均產(chǎn)生溫降變化.管壁到管道中心區(qū)域的傳熱速率較慢，導(dǎo)致管道中區(qū)區(qū)域液相稠油的溫度保持恒定.在同一截面位置和不同管壁溫度的條件下，管壁溫度越低，由管壁向管道中心區(qū)域溫度降低的分布區(qū)域越大，這說(shuō)明了管壁溫度與液相稠油的溫差越大，管道內(nèi)稠油傳熱效果更加明顯.因此，管壁溫度對(duì)于水平管道內(nèi)稠油流動(dòng)的傳熱特性以及蠟晶沉積的過(guò)程具有重要影響.

圖10為基于Filtered model Ⅰ和Filtered model Ⅱ 2種亞格子過(guò)濾曳力模型下過(guò)濾傳熱系數(shù)的修正系數(shù)Q隨著過(guò)濾顆粒相體積分?jǐn)?shù)的變化規(guī)律.當(dāng)Q為0時(shí)，表明介觀尺度下的亞格子結(jié)構(gòu)是均勻分布的；然而Q值越來(lái)越大時(shí)，表明介觀尺度下的亞格子結(jié)構(gòu)的分布越來(lái)越不均勻.由圖可知，F(xiàn)iltered model Ⅰ和Filtered model Ⅱ的過(guò)濾顆粒相體積分?jǐn)?shù)的變化范圍分別是8%～14%和1%～24%；其Q值分別從0.97和0.89開(kāi)始逐漸增大.Filtered model Ⅰ隨過(guò)濾顆粒相體積分?jǐn)?shù)的增加而增長(zhǎng)緩慢，然而Filtered model Ⅱ隨著過(guò)濾顆粒相體積增加的增長(zhǎng)速率先快后趨于穩(wěn)定.并且，F(xiàn)iltered model Ⅱ的Q最大值大于Filtered model Ⅰ的Q最大值.由此得出，加入壁面修正的Filtered model Ⅱ模擬所得的水平管道內(nèi)介觀尺度下蠟晶顆粒聚團(tuán)的分布更加不均勻，液相稠油與顆粒相間的過(guò)濾傳熱系數(shù)更小，意味著耗散的熱量較少.

4" 結(jié)" 論

1） Filtered model對(duì)于顆粒聚團(tuán)的預(yù)測(cè)更詳細(xì)，能夠揭示在稠油輸運(yùn)的水平管道內(nèi)固液兩相所呈現(xiàn)出的流動(dòng)為明顯的不均勻流動(dòng)特性.3種曳力模型中，F(xiàn)iltered model Ⅱ所得壓降結(jié)果與試驗(yàn)壓降值吻合較好，尤其在高顆粒體積分?jǐn)?shù)下，考慮了壁面效應(yīng)的Filtered model Ⅱ模型的優(yōu)越性更為明顯.

2） Filtered model Ⅱ的顆粒擬溫度的平均值最高，Gidaspow模型的顆粒擬溫度的平均值最低.顆粒擬溫度越高，表明顆粒間或顆粒聚團(tuán)間的碰撞增強(qiáng)，產(chǎn)生強(qiáng)烈的顆粒速度脈動(dòng).

3） 管壁溫度是影響水平管道內(nèi)輸送稠油傳熱特性的關(guān)鍵因素.壁面溫度低于液相稠油溫度時(shí)會(huì)產(chǎn)生徑向溫度梯度，利于蠟晶顆粒在管壁處沉積，形成高濃度的石蠟顆粒層.管壁與稠油溫度差越大，水平管道內(nèi)導(dǎo)熱和對(duì)流換熱的傳熱效果更加明顯.并且，考慮壁面效應(yīng)的Filtered model Ⅱ的傳熱效率更高.

參考文獻(xiàn)（References）

［1］" MARTNEZ PALOU R， MOSQUEIRA M de L， ZAPATA RENDN B， et al. Transportation of heavy and extra-heavy crude oil by pipeline： a review ［J］. Journal of petroleum science and engineering， 2011， 75（3/4）： 274-282.

［2］" HASAN S W， GHANNAM M T， ESMAIL N. Heavy crude oil viscosity reduction and rheology for pipeline transportation ［J］.Fuel， 2010，89（5）： 1095-1100.

［3］" ILYIN S O， ARININA M P， POLYAKOVA M Y， et al. Rheological comparison of light and heavy crude oils ［J］. Fuel， 2016， 186： 157-167.

［4］" JIN F Y， JIANG T T， YUAN C D， et al.An improved viscosity prediction model of extra heavy oil for high temperature and high pressure ［J］. Fuel，2022，319：123852.

［5］" AIYEJINA A， CHAKRABARTI D P， PILGRIM A， et al. Wax formation in oil pipelines： a critical review ［J］. International journal of multiphase flow， 2011，37（7）： 671-694.

［6］" EESA M， BARIGOU M. CFD Investigation of the pipe transport of coarse solids in laminar power law fluids ［J］. Chemical engineering science， 2009， 64（2）： 322-333.

［7］" OFEI T N， ISMAIL A Y. Eulerian-Eulerian simulation of particle-liquid slurry flow in horizontal pipe ［J］. Journal of petroleum engineering， 2016，2016： 5743471.

［8］" NERI A， GIDASPOW D. Riser hydrodynamics： simulation using Kinetic theory ［J］. AIChE journal， 2000， 46（1）： 52-67.

［9］" SYAMLAL M， O′BRIEN T J. Computer simulation of bubbles in a fluidized bed ［J］.AIChE symposium series， 1989，85：22-31.

［10］" 張赫銘，李文昊，何新林，等.不同管徑水平管道氣液兩相流動(dòng)數(shù)值模擬［J］.排灌機(jī)械工程學(xué)報(bào)，2021，39（5）：488-494.

ZHANG Heming， LI Wenhao， HE Xinlin， et al. Nume-rical simulation of gas-liquid two-phase flow in horizon-tal pipeline with different diameters［J］. Journal of drainage and irrigation machinery engineering， 2021，39（5）：488-494.（in Chinese）

［11］" KAUSHAL D R， THINGLAS T， TOMITA Y， et al. CFD modeling for pipeline flow of fine particles at high concentration ［J］. International journal of multiphase flow， 2012，43： 85-100.

［12］" SINGH P， VENKATESAN R， FOGLER H S， et al. Formation and aging of incipient thin film wax-oil gels ［J］.AIChE journal， 2000，46（5）： 1059-1074.

［13］" IGCI Y， ANDREWS A T， SUNDARESAN S， et al. Filtered two-fluid models for fluidized gas-particle suspensions ［J］. AIChE journal， 2008，54（6）： 1431-1448.

［14］" IGCI Y， PANNALA S， BENYAHIA S， et al. Validation studies on filtered model equations for gas-particle flows in risers ［J］. Industrial and engineering chemistry research， 2012，51（4）： 2094-2103.

［15］" AGRAWAL K， HOLLOWAY W， MILIOLI C C， et al. Filtered models for scalar transport in gas-particle flows［J］. Chemical engineering science， 2013， 95：291-300.

［16］" LEI H， LIAO J W， ZHU L T， et al. CFD-DEM mode-ling of filtered fluid-particle drag and heat transfer in bidisperse gas-solid flows［J］. Chemical engineering science， 2021， 246： 116896.

［17］" MILIOLI C C， MILIOLI F E. On the subgrid behavior of accelerated riser flows for a high stokes number particulate ［J］. Industrial and engineering chemistry research，2011，50（23）： 13538-13544.

［18］" WANG S Y， YANG Q， SHAO B L， et al. Numerical simulation of horizontal jet penetration using filtered two-fluid model in a gas-solid fluidized bed ［J］.Powder technology， 2015，276： 1-9.

［19］" WANG S Y， SHAO B L， LI X Y， et al. Simulations of vertical jet penetration using a filtered two-fluid model in a gas-solid fluidized bed ［J］. Particuology， 2017，31（2）： 95-104.

［20］" 周嶺，韓晨，白玲，等.基于不同曳力模型的鼓泡流化床CFD-DEM數(shù)值模擬與試驗(yàn)研究［J］.排灌機(jī)械工程學(xué)報(bào)，2022，40（1）：49-54.

ZHOU Ling， HAN Chen， BAI Ling，et al.Numerical simulation and experimental study of CFD-DEM in bubbling fluidized bed based on different drag model［J］. Journal of drainage and irrigation machinery enginee-ring， 2022，40（1）：49-54.（in Chinese）

［21］" ZAMBRANO H， SIGALOTTI L D G， KLAPP J， et al.Heavy oil slurry transportation through horizontal pipelines： experiments and CFD simulations ［J］.Interna-tional journal of multiphase flow， 2017，91： 130-141.

［22］" MESSA G V， MALAVASI S. Improvements in the numerical prediction of fully-suspended slurry flow in horizontal pipes ［J］.Powder technology， 2015， 270（part A）： 358-367.

（責(zé)任編輯" 朱漪云）