海洋工程装备与技术  2015 , 2 (1): 23-27

海洋油气勘探开发技术与装备

双流体模型中人工扩散的影响数值分析

韩朋飞1, 郭烈锦1*, 李清平2, 姚海元2, 程兵2

1. 西安交通大学动力工程多相流国家重点实验室,陕西 西安 710049
2. 中海油研究总院,北京 100028

Numerical Analysis of the Effect of Artificial Diffusion Terms on Two-Fluid Model

HAN Peng-fei1, GUO Lie-jin1, LI Qing-ping2, YAO Hai-yuan2, CHENG Bing2

1. State Key Laboratory of Multiphase Flow in Power Engineering, Xi'an Jiaotong University, Xi'an, Shaanxi 710049, China
2. CNOOC Research Institute, Beijing 100028, China

中图分类号:  O359

文献标识码:  A

文章编号:  2095-7297(2015)01-0023-05

收稿日期: 2014-12-29

网络出版日期:  2015-02-20

版权声明:  2015 海洋工程装备与技术编辑部 版权所有

基金资助:  国家科技重大专项(2011ZX05026-004-02),国家自然科学基金重点项目(51236007)

作者简介:

作者简介:韩朋飞(1983--),男,博士研究生,主要从事气液两相流方面的研究.

展开

摘要

研究了双流体模型中加入人工扩散项后不同扩散项系数对模拟的影响.对于竖直管的Water Faucet问题,逐步增加人工扩散系数后非物理振荡完全消除.对于水平管流动病态区域的工况,不采用扩散系数在中等精细度网格上计算的结果同实验比较吻合,但无法得到网格无关解.采用Water Faucet问题确定的人工扩散系数进行计算,水动力学不稳定性被抑制,液塞频率被低估.因此,人工扩散的量随不同流动条件发生变化,需要具体确定.

关键词: 双流体模型 ; 人工扩散 ; 段塞流 ; 数值模拟

Abstract

When artificial diffusion terms are introduced in two-fluid model, the effects of different artificial diffusion coefficients are studied. For vertical pipe Water Faucet problem, when the artificial diffusion coefficients are gradually increased, nonphysical oscillations can be eliminated. For horizontal pipe flow in ill-posed region, the calculated results without diffusion terms agree with experimental data but the mesh independent solution cannot be obtained. When the coefficients calibrated in Water Faucet problem are adopted, hydrodynamic instability is suppressed to some extent and slug frequency is underestimated. Therefore, the amount of artificial diffusion varies with different flow conditions and should be determined specially.

Keywords: two-fluid model ; artificial diffusion ; slug flow ; numerical simulation

0

PDF (1452KB) 元数据 多维度评价 相关文章 收藏文章

本文引用格式 导出 EndNote Ris Bibtex

韩朋飞, 郭烈锦, 李清平, 姚海元, 程兵. 双流体模型中人工扩散的影响数值分析[J]. , 2015, 2(1): 23-27 https://doi.org/

HAN Peng-fei, GUO Lie-jin, LI Qing-ping, YAO Hai-yuan, CHENG Bing. Numerical Analysis of the Effect of Artificial Diffusion Terms on Two-Fluid Model[J]. 海洋工程装备与技术, 2015, 2(1): 23-27 https://doi.org/

0 引 言

管道内多相流流动模拟对石油工业的正常生产及安全保障有非常重要的意义.目前,石油开采正在向海洋尤其是深海进发,因此超长距离的管道多相混输对管内流动的预测有着更高和更迫切的需求.

根据以往国内外的研究,一维的漂移流模型和双流体模型是模拟预测超长距离管道流动的有力工具.漂移流模型本质上属于混合物模型,无法对相界面上波动的发展进行预测,因此精确度有限.Lin等[1]很早就从理论上证明了双流体模型可以预测界面上波动的发展,并且Issa等[2]也使用双流体模型在数值模拟中捕捉到了界面波动的实时发展.因此,双流体模型是对长距离管线流动进行准确预测的最好手段.但是,传统的标准双流体模型在某些区域不能保证双曲性,界面上小的波动会被无限放大.这种特征表现在数值模拟当中,就是当网格步长不断减小时,不断有更小的波动被计算分辨,而计算对这些波动的放大系数不是有界的,在波长极小时趋于无穷大,从而无法得到唯一正确的收敛结果.

针对这一问题,很多学者采用各种方法进行了研究.一方面,Zanotti等[3]从纯数学的角度上考虑,通过矩阵预处理实现了无条件的双曲性,尽管对算例的模拟显示瞬态结果也是准确的,但是理论上预处理影响了模型的瞬态项,只适用于稳态问题.另一方面,很多学者尝试引入其他物理机制来保证双曲性,得到了一些有用的结论.Chung等[4]发现虚拟质量力的加入并不是完全有益的,有时它会影响计算的稳定性及精确度.Holmas等[5]则指出表面张力系数作用范围有限,人为增大表面张力系数虽然可以抑制小波长不稳定性的发展,但是这部分非物理的能量会转移到长波里,而通过人工增加二阶扩散项耗散掉不稳定的小波长波动的能量是一种合理的方法.Song等[6]通过考虑截面速度不均匀引入动量通量参数,也在广泛的区域内得到了良态的系统.但是,Issa等[7]发现动量通量参数的预测结果同实验偏离.他们经过同实验的对比发现Holmas的人工扩散方法能够在网格不断加密的情况下得到同实验吻合的预测结果[8].但是,他们只对水平管的一个工况进行了研究,提出的系数对竖直管计算的作用还不清楚,需要进一步研究.

本文针对目前保持双流体模型双曲性的研究,对双流体模型模拟竖直管内的流动进行探索,寻求获得正确预测结果所需人工扩散项系数,并验证此组系数计算水平管段塞流时的精确度.

1 控制方程

加入人工扩散项后的一维双流体模型如下:

αlρlt+αlρlulx-xρlφalαlx=0,(1)

αgρgt+ αgρgugx- xρgφagαgx=0,(2)

αlρlult+ αlρlululx- xαlμalulx

=-αlpxlρlgsinθ-αlρlgcosθ hlx-

τwlSwlA+ τiSiA,(3)

αgρgugt+ αgρgugugx- xαgμagugx=-αgpxgρggsinθ-αgρggcosθ hlx- τwgSwgA- τiSiA,(4)

αlg=1,(5)

式中:下标lg分别表示液相和气相;下标wi分别表示在壁面处和相界面处;α表示相含率;ρ表示密度;u表示速度;p表示压力;g表示重力加速度;hl表示液膜厚度;τ表示切应力;θ表示管道倾角;S表示湿周;A表示管道横截面积.气相为理想可压缩气体,液相为不可压缩流体.

式(1)(4)左边第三项二阶偏导即为加入的人工扩散项,φal和φag分别表示液相和气相连续方程的人工扩散项系数,μal和μag分别表示液相和气相动量方程的人工扩散项系数.剪切应力的表达式如下:

τ= 12urur,(6)

式中:f为Fanning摩擦系数;两相的相对速度ur=ug-ul.

对于水平管和竖直管,两相的湿周计算分布如图1所示.

图1   湿周及两相分布示意图

Fig.1   Schematic representation of wetted perimeters and two-phase distributions

水平管相界面上的摩擦系数采用Taitel-Dukler模型[9]:

fi=16ReiRei<21000.046Rei-0.2Rei2100;(7)

竖直管相界面上的摩擦系数采用Wallis[10]的模型:

fi=0.005 1+300hlD,(8)

D为管道直径.水平管气相与壁面间的摩擦系数同样采用Taitel-Dukler的模型[9]:

fg=16RegReg<21000.046Reg-0.2Reg2100;(9)

水平管液相与壁面间的摩擦系数采用Spedding-Hand模型[11]:

fl=24RelRel<21000.0262αlRel-0.139Rel2100;(10)

竖直管液相与壁面间的摩擦系数采用Wallis模型[10]:

fl=16RelRel<23000.079Rel-0.25Rel2300.(11)

雷诺数的表达式为

Rel=αlρlulDμl,Reg=ρgugDgμg,Rei=ρgul-ugDgμg.(12)

2 数值方法

气液相的连续和动量方程式(1)(4)通过有限体积法在同位网格上进行求解,采用Rhie-Chow插值以消除压力速度的失耦.压力修正方程采用两相混合的质量守恒方程导出,其形式如下:

αgρgDρgDt+αlρlDρlDt+∇(αguglul)=0,(13)

式中:D/Dt为所给量的物质导数.

为了提高计算的稳定性,时间离散格式为一阶欧拉隐式后向差分,空间离散格式采用二阶迎风.为了保证计算收敛性,压力与速度的耦合采用SIMPLE算法,并设定全局最大Courant数为0.1,根据Courant数自动选取时间步长,计算收敛的判据为1×10-6.

3 结果及讨论

3.1 Water Faucet问题

Water Faucet问题是检验数值计算结果精确度的常用标准方法,其过程如图2所示.其结构为一段12 m长,内径为1 m的竖直管道,上部入口含气率为0.2,气相和液相的速度分别为0和10 m/s,出口压强固定为105 Pa.气相的密度为1.16 kg/m3,液相的密度为1000 kg/m3,整个系统的温度固定为50 ℃.此问题的解析解为

α(x,t)=αinVinVin2+2g(x-xin)xxin+Vint+12gt2αin其他;(14)

V(x,t)=Vin2+2g(x-xin)xxin+Vint+12gt2Vin+gt其他. (15)

当采取不同的扩散项系数和网格控制容积(cv)数量时,模拟得到的t=0.5 s和t=2.5 s时的结果和解析解的结果对比分别如图3图4所示.从图3可以看出,当未引入人工扩散时,双流体模型是不稳定的病态问题,网格密度为500 cv时,含气率曲线在梯度剧变处发生了大幅振荡.而如图4显示,未引入人工扩散的计算并不能得到稳定的稳态结果.

图2   Water Faucet问题示意图

Fig.2   Schematic representation of Water Faucet problem

图3   不同扩散项系数及网格密度在t=0.5 s时模拟与解析解比较

Fig.3   Comparison of simulation and analytical results for different diffusion coefficients and mesh densities at t=0.5 s

图4   不同扩散项系数及网格密度在t=2.5 s时模拟与解析解比较

Fig.4   Comparison of simulation and analytical results for different diffusion coefficients and mesh densities at t=2.5 s

当扩散系数采用Issa等[8]的系数时,含气率的非物理振荡已被减弱但尚未消失.若进一步增大动量扩散系数,计算没有明显的改善,因此逐步增大质量扩散系数.质量扩散系数从0.1开始逐步增大时的计算结果与解析解的比较如图5所示.从图中可从看到随着质量扩散系数的逐步增大,非物理波动的幅度逐渐减小;当质量扩散系数达到0.7时,计算结果已经不存在非物理振荡.虽然由于扩散的影响,含气率的曲线变得更加平滑,但是能够得到符合物理过程的无振荡解,损失的精度是可以接受的.当网格密度进一步加密到1500 cv和3000 cv时,仍旧没有非物理振荡出现,可以认为此时能够得到所求解问题的网格无关解.

图5   不同质量扩散项系数在t=0.75 s时模拟与解析解比较

Fig.5   Comparison of simulation and analytical results for different diffusion coefficients at t=0.75 s

3.2 水平管段塞流

Woods等[12]列举的内径为0.095 m水平管内各种表观气速uSG和表观液速下uSL的段塞流液塞频率如图6所示.本文选取了气液相入口表观速度分别为15.9 m/s和1.2 m/s远离稳定边界的病态区域工况,分别使用较粗网格不加入人工扩散和人工扩散系数φ=0.7,μ=1×10-5进行计算.

图6   Woods等[2]列举的液塞频率实验数据

Fig.6   Experimental data for slug frequency presented by Woods et al[12]

在网格密度为1500 cv时,计算得到的无量纲液位高度(hl/D)的瞬态曲线随时间的变化如图7所示,图中相邻曲线间的时间间隔为0.4 s.计算初始状态管道内无量纲液位高度为均匀的0.5,随着气液间相互作用的发展,约2 s时相界面的波动开始出现并在约3 s时发展为液塞.由于初始条件设置,液塞前面的液膜厚度较大,液塞长度能够不断增长,因此需等第一个长液塞流出后才能反映真实的流动.从后续的发展曲线来看,并未出现此种超长的液塞.在第一个超长液塞流出管道后,统计得到的液塞频率为0.88 Hz,与0.79 Hz的实验结果相比,误差为11.4%.虽然计算结果同实验相比很接近,但是当加密网格后,进一步被解析的短波波动使得计算很快发散,无法得到网格独立解.

图7   模拟得到的沿管道无量纲液位高度随时间的变化(φ=μ=0)

Fig.7   Time trace of simulated non-dimensional liquid film height along the pipe (φ=μ=0)

采用前文验证的扩散系数在网格密度为3000 cv的细网格上进行模拟,得到的无量纲液位高度的瞬态曲线随时间的变化如图8所示.从图8可以看出,液面高度的波动被抑制,曲线变得平滑,液塞形成的频率降低,统计得到的液塞频率降低为0.58 Hz,与实验相比减小了26.5%,显然此时人工扩散量是过大的.但是,进一步加密网格未能显著改善预测结果,证明人工扩散量在精细网格上仍旧过大,得到的收敛结果对液塞频率的预测偏低.可以预测当扩散量很大时,双流体模型的性能类似于混合物模型.

图8   模拟得到的沿管道无量纲液位高度随时间的变化(φ=0.7,μ=1×10-5)

Fig.8   Time trace of simulated non-dimensional liquid film height along the pipe(φ=0.7,μ=1×10-5)

4 结 语

通过加入人工扩散可以使双流体模型在原本病态的区域内转变为良态问题.经过一系列的数值实验,得出在人工扩散项系数为φ=0.7,μ=1×10-5时,竖直管Water Faucet问题非物理振荡完全消失.对于水平管内流动病态区域的工况,φ=0.7,μ=1×10-5时水动力学不稳定性被抑制,液塞频率降低.水平管人工扩散项的系数需要根据不同工况决定,有待进一步研究.

The authors have declared that no competing interests exist.


参考文献

[1] Lin P Y, Hanratty T J.

Prediction of the initiation of slugs with linear stability theory

[J]. International Journal of Multiphase Flow, 1986, 12(1): 79.

DOI:10.1016/0301-9322(86)90005-4      URL      [本文引用: 1]      摘要

This paper explores the application of linear stability theory to explain the onset of slugging. It is shown that the inviscid Kelvin-Helmholtz theory correctly predicts stability of a stratified flow only for very large liquid viscosities. In general, however, inviscid theory is in error because it ignores the destabilizing effect of liquid inertia. Good agreement is noted between the linear stability analysis and observations of the initiation of slugs in 2.54- and 9.53-cm horizontal pipes at superficial gas velocities less than 3.3 m/s.
[2] Issa R I, Kempf M H W.

Simulation of slug flow in horizontal and nearly horizontal pipes with the two-fluid model

[J]. International Journal of Multiphase Flow, 2003, 29(1): 69.

DOI:10.1016/S0301-9322(02)00127-1      Magsci      [本文引用: 2]      摘要

<h2 class="secHeading" id="section_abstract">Abstract</h2><p id="">A mechanistic approach to the prediction of hydrodynamic slug initiation, growth and subsequent development into continuous slug flow in pipelines is presented. The approach is based on the numerical solution of the one-dimensional transient two-fluid model equations. The advantage of this approach is that the flow field is allowed to develop naturally from any given initial conditions as part of the transient calculation; the slugs evolve automatically as a product of the computed flow development. The need for the many phenomenological models for flow regime transition, formation of slugs and their dynamics can thus be minimized.</p><p id="">It is shown that when the two-fluid model is invoked within the confines of the conditions under which it is mathematically well-posed, it is capable of capturing the growth of instabilities in stratified flow leading to the generation of slugs. The computed rates of growth of such instabilities compare well with the values obtained from Kelvin&ndash;Helmholtz analyses. Simulations are then carried out for a large number of pipe configurations and flow conditions that lead to slug flow. These include horizontal, inclined and V-section pipes. The results of computations for slug characteristics are compared with data obtained from the literature and it is found that the agreement is remarkable given the simplicity of the one-dimensional model.</p>
[3] Zanotti A L, Méndez C G, Nigro N M, et al.

A preconditioning mass matrix to avoid the ill-posed two-fluid model

[J]. Journal of Applied Mechanics, 2007, 74(4): 732.

DOI:10.1115/1.2711224      URL      [本文引用: 1]      摘要

ABSTRACT Two-fluid models are central to the simulation of transport processes in two-phase homogenized systems. Even though this physical model has been widely accepted, an inherently non-hyperbolic and non-conservative ill-posed problem arises from the mathematical point of view. It has been demonstrated that this drawback occurs even for a very simplified model, i.e., an inviscid model with no interfacial terms. Lots of efforts have been made to remedy this anomaly and in the literature two different types of approaches can be found. On one hand, extra terms with physical origin are added to model the interphase interaction, but even though this methodology seems to be realistic, several extra parameters arise from each added term with the associated difficulty in their estimation. On the other hand, mathematical based-work has been done to find the way to remove the complex eigenvalues obtained with two-fluid model equations. Preconditioned systems, characterized as a projection of the complex eigenvalues over the real axis, may be one of the choices. The aim of this paper is to introduce a simple and novel mathematical strategy based on the application of a preconditioning mass matrix that circumvents the drawback caused by the non-hyperbolic behavior of the original model. Although the mass and momentum conservation equations are modified, the target of this methodology is to present another way to reach a steady state solution (using a time marching scheme), greatly valued by researchers in industrial process design. Attaining this goal is possible because only the temporal term is affected by the preconditioner. The obtained matrix has two parameters that correct the non-hyperbolic behavior of the model: the first one modifies the eigenvalues removing their imaginary part and the second one recovers the real part of the original eigenvalues. Besides the theoretical development of the preconditioning matrix, several numerical results are presented to show the validity of the method. [To appear in Journal of Applied Mechanics]
[4] Chung M-S, Lee S-J.

A modified semi-implicit method for a hyperbolic two-fluid model

[J]. Appl Numer Math, 2009, 59(10): 2475.

DOI:10.1016/j.apnum.2009.05.005      URL      [本文引用: 1]      摘要

By introducing the interfacial pressure jump terms based on a surface tension into the momentum equations of a two-phase two-fluid model, the mathematical property of the governing equations is changed to a hyperbolic type. Then the eigenvalues of the equation system always become always real values representing the void wave and the pressure wave propagation speeds as shown in the present author's former article: Numerical Heat Transfer – Part B (40) (2001) 83–97 . To solve the interfacial pressure jump terms with void fraction gradients implicitly, the conventional semi-implicit method should be modified by inserting an intermediate calculation process for a void fraction at a fractional time step. This modified semi-implicit method then becomes stable without conventional additive terms. Consequently, by including the interfacial pressure jump terms with the modified semi-implicit method, the numerical calculations of the void discontinuity propagation and water faucet problems can become more stable and sound than those calculated by using virtual mass terms.
[5] Holmas H, Sira T, Nordsveen M, et al.

Analysis of a 1D incompressible two-fluid model including artificial diffusion

[J]. IMA Journal of Applied Mathematics, 2008, 73(4): 651.

DOI:10.1093/imamat/hxm066      URL      [本文引用: 1]      摘要

This article examines a 1D incompressible two-fluid model including artificial tensor diffusion. The aim is to obtain a formulation that provides convergent numerical solutions for all flow conditions within the stratified and the stratified wavy flow regime. With appropriate simplifications, the two-fluid model reduces to one momentum balance, one mass conservation and two algebraic equations. It has previously been established that a formulation that is well posed in possessing exclusively real characteristics can be obtained by including an axial diffusion term in the momentum balance. In this article, however, we demonstrate that this is not sufficient to obtain a system suitable for numerical simulations. Although the unbounded growth rates of the standard two-fluid model are eliminated, linear stability theory predicts that infinitesimal wavelengths still experience finite growth. This entails that grid refinement always will result in new unstable wavelengths being resolved. On the other hand, if artificial axial diffusion is added to both the mass and the momentum equations as suggested here, a cut-off wavelength is established below which all wavelengths are stable. Thus, a numerically converging model is formed, which retains the long-wavelength properties of the standard two-fluid model. The conclusions of the mathematical analysis are substantiated by numerical simulations of 1D gravity waves.
[6] Song J H, Ishii M.

The well-posedness of incompressible one-dimensional two-fluid model

[J]. International Journal of Heat and Mass Transfer, 2000, 43(12): 2221.

DOI:10.1016/S0017-9310(99)00287-2      URL      [本文引用: 1]      摘要

A characteristic analysis on the stability of the governing differential equations for an incompressible one-dimensional two-fluid model is presented. The stability criteria are newly proposed in terms of the momentum flux parameters by incorporating the effect of void fraction and velocity profiles. A simplified two-phase flow configuration constructed by using existing correlation for distribution parameter and experimentally correlated velocity profiles is selected to test the validity of the proposed theory. The curve of calculated momentum flux parameters is compared with the curve of stability criteria. The simplified flow is found to be stable within a wide range of void fraction.
[7] Issa R I, Montini M.

Applicability of the momentum-flux-parameter closure for the two-fluid model to slug flow

[C]. 6th International Symposium on Multiphase Flow, Heat Mass Transfer and Energy Conversion, 2009: 712.

[本文引用: 1]     

[8] Issa R I, Montini M.

The effect of surface tension and diffusion on one-dimensional modelling of slug flow instabilities

[C]. 7th International Conference on Multiphase Flow, 2010.

[本文引用: 2]     

[9] Taitel Y, Dukler A E.

A model for predicting flow regime transitions in horizontal and near horizontal gas-liquid flow

[J]. AIChE Journal, 1976, 22(1): 47.

DOI:10.1002/aic.690220105      URL      [本文引用: 2]      摘要

Abstract Models are presented for determining flow regime transitions in two-phase gas-liquid flow. The mechanisms for transition are based on physical concepts and are fully predictive in that no flow regime transitions are used in their development. A generalized flow regime map based on this theory is presented.
[10] Wallis G B.

One-Dimensional Two-Phase Flow

[M]. New York: McGraw-Hill,1969.

[本文引用: 2]     

[11] Spedding P L, Hand N P.

Prediction in stratified gas-liquid co-current flow in horizontal pipelines

[J]. International Journal of Heat and Mass Transfer, 1997, 40(8): 1923.

DOI:10.1016/S0017-9310(96)00252-9      URL      [本文引用: 1]      摘要

The predictive performance of existing models based on momentum balances has been shown to be generally unsatisfactory. An iterative procedure was developed for the prediction of pressure drop and holdup that incorporated new relationships for the interfacial and liquid friction factors in the solution of the phase momentum balance equations for two phase horizontal co-current flow. The method adequately predicted data for the film plus droplet, annular roll wave and stratified type regimes. Successful performance was achieved regardless of fluid properties or pipe diameter. For gas flow rates, where non-uniform stratified flow with an interfacial level gradient occurred, this method of prediction was inaccurate and open channel flow theory recommended. Thus the model did not apply to the intermittent or non-uniform flow regimes.
[12] Woods B, Fan Z, Hanratty T.

Frequency and development of slugs in a horizontal pipe at large liquid flows

[J]. International Journal of Multiphase Flow, 2006, 32(8): 902.

DOI:10.1016/j.ijmultiphaseflow.2006.02.020      URL      [本文引用: 2]      摘要

The evolution of slugs from a highly disturbed stratified flow at the inlet was studied by measuring the holdup at a number of locations along a 20m length of 0.0763m pipe. At superficial gas velocities less than 3m/s incipient slugs form by waves touching the top of the pipe very close to the entry. At U SG >3鈥4m/s incipient slugs form by wave coalescence farther downstream. These incipient slugs grow as long as the height of the stratified flow in front of them is greater than h 0 . When equal to h 0 , the slug can propagate downstream or decay if it has not reached a stable length. Decaying slugs form large amplitude roll waves which propagate downstream at a lower velocity than slugs. Slugs that overtake these waves increase in size. These results can be used in developing a stochastic model for the evolution.

沪ICP备15013849号-1
网站版权所有 © 《海洋工程装备与技术》编辑部
地址:上海市徐汇区广元西路55号上海交通大学徐汇校区慧谷科技楼4楼 邮编:200030 
电话:021-62800061 E-mail:hgzbyjs@163.com
本系统由北京玛格泰克科技发展有限公司设计开发

/