Gas-Liquid Stratified Flow in Pipeline with Phase Change

Abstract

When the natural gas with vapor is flowing in production pipeline, condensation occurs and leads to serious problems such as condensed liquid accumulation, pressure and flow rate fluctuations, and pipeline blockage. This chapter aims at studying phase change of vapor and liquid-level change during the condensing process of water-bearing natural gas characterized by coupled hydrothermal transition and phase change process. A hydrothermal mass transfer coupling model is established. The bipolar coordinate system is utilized to obtain a rectangular calculation domain. An adaptive meshing method is developed to automatically refine the grid near the gas-liquid interface. During phase change process, the temperature drop along the pipe leads to the reduction of gas mass flow rate and the rise of liquid level, which results in further pressure drop. Latent heat is released during the vapor condensing process which slows down the temperature drop. Larger temperature drop results in bigger liquid holdup while larger pressure drop causes smaller liquid holdup. The value of velocity with phase change is smaller than that without phase change while the temperature with phase change is bigger. The highest temperature locates in gas phase. But near the pipe wall the temperature of liquid region is higher than gas region.

1. Introduction

Condensation occurs when the natural gas with vapor is flowing in production pipeline and leads to serious problems such as condensed liquid accumulation, pressure and flow rate fluctuations, and pipeline blockage or corrosion. The pipe flow together with phase change is commonly encountered in various heat and mass transfer processes over the past four decades, for instance, in petroleum and chemical processing industry, steam-generating equipment, nuclear reactors, geothermal fields, heat exchangers, cooling systems, and solar energy system [1, 2, 3, 4]. In petroleum transportation, two-phase flow characterization is a very common and economic technique, where vapor-liquid two-phase stratified flow is often observed in horizontal or slightly inclined systems [5, 6].

There exist several problems in the pipeline network system that the saturated vapor in gas would condense due to pressure and temperature drop [7, 8]. The condensate would attach to the pipe wall as a form of film or droplet [9, 10]. The condensation will decrease the effective cross-sectional area and cause the increase of pressure drop which may lead to system shutdown [11, 16]. Generally, the condensed water accumulates at the lower parts of the pipeline due to the hilly pipeline route topography, which results in a continuous change of liquid holdup along the pipeline [12, 13, 14]. The changing liquid holdup and flow area are bounded to affect the flow patterns which inevitably influence the operating pressure and temperature inversely. Thus, the flow of condensed water and water-bearing gas in production pipelines is a complex process with coupling of hydraulic, thermal, and phase change phenomena [14, 15, 16, 17]. Researchers have investigated the gas-liquid two-phase pipe flow system by experiments or hydrodynamic and thermodynamic models.

It has been observed experimentally that when phase change occurs during the saturated vapor pipeline transportation process, the thermal gradients are created in the wall of the pipeline that lead to severe liquid condensation and stratified vapor-liquid two-phase flow [17, 18]. The fundamental engineering parameters are the pressure drop, liquid holdup, phase fraction, phase flow rate, temperature, thermo-physical properties of the fluids, and pipeline geometry [19].

Not limited to hydraulic parameters, more recent attention has turned to non-isothermal flow in a pipe or plane channel where some numerical studies are also found [20, 21, 22, 23, 24, 25, 26, 27]. The detailed characteristic of heat transfer is taken into consideration in these mentioned models instead of an average value being represented for the temperature profiles in a pipe [28, 29]. According to their studies, the wall temperature distribution is different from the assumption of fully developed isothermal state [30]. The energy transfer model has been taken into the flow progress for the optimization of transportation, estimation of corrosion, or prediction of wax deposition [31, 32]. Concretely, the two-dimensional (2D) momentum and three-dimensional (3D) energy equations for both phases have been established for dynamic and thermal numerical simulation [8, 27, 30, 33]. The smooth or wavy interface between phases was obtained in a different range of flow rates [17]. For such a two-phase non-isothermal stratified flow, analytical and numerical heat transfer solutions limited to laminar flow and without phase change have been obtained for fully developed stratified flow under different thermal boundary conditions [27]. Then, solutions which are more applicable to fully developed turbulent gas-liquid smooth stratified flow have been obtained through the use of high Reynolds model [30]. Recently, the steady-state axial momentum and energy equations coupled with a low Reynolds model were established and solved [34]. The pressure drop, liquid height, and temperature field which are included in the solutions could match well with the experimental data. However, although equation of state (EOS) was utilized in previous one-dimensional (1D) models to calculate the phase fraction [10, 12, 16, 22, 34, 41, 44], the flow rate, temperature, and pressure were not coupled with the varying liquid level.

According to the description of the physical process of gas-liquid two-phase pipe flow, the model can be divided into isothermal model, gas-liquid two-phase pipe flow model coupled with heat and gas-liquid two-phase pipe flow model coupled with heat and mass transfer. In the gas-liquid two-phase pipe flow model under the isothermal condition, it is assumed that all phases are in thermodynamic equilibrium state without considering the heat transfer process between pipe flow and environment, and the physical parameters of gas-liquid two-phase are just the single value function of pressure. In the gas-liquid two-phase flow model coupled with heat, the heat transfer of gas-liquid and surrounding environment is considered. In the gas-liquid two-phase flow model coupled with heat and mass transfer, the coupling effect of flow, heat transfer, and mass transfer are considered simultaneously.

The two-dimensional (2D) or three-dimensional (3D) stratified gas-liquid two-phase model including mass conversation equation, momentum conversation equation, turbulence model equation, boundary conditions, and related auxiliary equations for model closure were applied to describe the flow in pipelines. Differences among them mainly existed in two aspects. On the one hand, different turbulence models were built including Spalart-Allmaras Model (SAM), $k-\epsilon$ series model, Reynolds Stress Model (RSM), Direct Numerical Simulation (DNS), and so on. On the other hand, different gas-liquid interface configuration models were built including volume of fluid (VOF) method lattice Boltzmann method, level set method, smooth interface model, wavy interface model, and curve interface model. The relationship between pressure drop, liquid holdup, velocity distribution, turbulent viscosity, and shear stress has been determined through these formula or methods. The flow pattern and flow regime and even secondary flow have also been discussed under different boundary conditions. However, the energy conservation equation and phase change of fluid components were not taken into consideration in those studies.

Recently, attempts have been made to introduce energy equation into the improved model and the detailed solutions about temperature distribution have also been worked out by considering potential energy, kinetic energy, heat transfer, and Joule-Thomson effect [17, 35, 36, 37, 38]. The phase change was ignored in these models. Although equation of state (EOS) was utilized in previous studies to calculate gas condensation of gas-condensate flow in pipelines [12, 13, 14, 15, 17, 18, 19, 20, 39], the flow rate, temperature, pressure were not coupled with liquid level. Turbulent flow is not considered, which would lead to different numerical results in their one-dimensional (1D) model. Moreover, 1D model could not present the detailed distribution of hydraulic and thermal parameters at pipe cross-section.

This chapter mainly introduces the different turbulence models, the interface shape model, and the phase transition model in the process of gas-liquid two-phase stratified flow in the horizontal pipeline. The turbulence model mainly includes the $k-\epsilon$ model, the large eddy simulation (LES) model, and the turbulence model under the low Reynolds number condition at the near wall surface. The interface shape mainly includes the flat interface model, the wavy interface model, the curve interface model, and the wavy-curve interface model. As a result of the phase transition, a varying wavy curve interface will be produced. If the phase transition is large, the varying wavy curve interface may gradually become varying wavy flat interface, and finally lead to the change of flow pattern, no longer the stratified gas-liquid two-phase flow. The phase transformation model mainly includes the empirical formula of water vapor phase transformation, the Peng-Robinson (PR) equation of hydrocarbon mixture, the SRK equation, and the PT flash model of water hydrocarbon impurities. Also, this chapter compares the different research model, enriches the contents of the book, demonstrates the development process of different mathematical models in different research points, and finally introduces the development direction of each point in the future and the relationship between the complex model and research points.

2. Numerical modeling of stratified gas-liquid wavy pipe flow with phase change

Since pressure and temperature drop along the pipeline, the saturated vapor of hydrocarbons would condenses gradually, as shown in Figure 1. The condensed liquid would attach to the pipe wall as liquid film or accumulate at the lower part of the pipeline, which could result in continuous change of the liquid level [17, 29, 30].

Figure 1.

Vapor-liquid two-phase flow coupled with heat and mass transfer.

Thus, some models of stratified vapor-liquid two-phase flow coupled with phase change were proposed. Assumptions could be made as follows: (1) precipitation of condensed liquid is a flash evaporation equilibrium process which occurs in a short moment; (2) regardless of the attachment to inner wall of the pipe, all the condensed liquid accumulates at the bottom of the pipe; (3) two-phase flow in vapor-liquid pipeline has a stratified flow pattern as well as a stable developing flow area in every calculated pipeline segment; (4) the smooth vapor-liquid interface model is adopted to describe the interface shape; (5) the heavy components of hydrocarbon are simplified as pseudo-component ${\mathrm{C}}_7^{+}$ ; and (6) without considering the effect of gravity on the P–T flash process.

Several kinds of flow patterns are likely to form for vapor/condensation flowing in pipeline: stratified flow, slug flow, annular flow, and stratified-dispersed flow. It is hard to calculate the liquid level, the exact position of liquid film, and the migration patterns of condensed liquid (as shown in Figure 1) via one-dimensional model which could not present the detailed distribution of hydraulic and thermal parameters at pipe cross section. Moreover, the 1D model does not consider the turbulent flow which would lead to different numerical results. Hence, the control volume in three-dimensional coordinates is adopted to discretize the calculating area, as shown in Figure 2.

Figure 2.

Schematic illustration of stratified vapor–liquid flow with phase change.

In this section, the governing equations based on physical conservation are chosen and established in three-dimensional coordinates, which include mass conservation equation, momentum conservation equation, energy conservation equation, turbulent flow model, and phase change model.

2.1. Mass conservation equation

The total mass flow rate $G$ is the sum of vapor and liquid mass flow rates and keeps constant during the flowing process:

$$G_G+G_L\equiv G$$ E1

The mass transfer rate of phase-change within pipe segment has been taken into consideration due to the vapor phase gradually condensing or the liquid evaporating during the flow process. The change of vapor mass flow rate $\Delta G_G$ is equal to opposite number of that of liquid $-\Delta G_L$ :

$$\Delta G_G+\Delta G_L=\mathit{\Delta G}\equiv 0$$ E2

The liquid volume fraction $H_L$ , vapor volume fraction $H_G=1-H_L$ , pressure, vapor, and liquid velocities as well as temperature are the principal unknowns in the cross-sectional area along the pipeline. The derivation of mass conservation equation is presented here based on two-fluid approach. The mass conservation equations for each phase within the control volume are given as:

$$H_L{\rho }_L\frac{\partial w_L}{\partial z}+{\rho }_L{w}_L\frac{\partial H_L}{\partial z}+H_L{w}_L{\left.\frac{\partial {\rho }_L}{\partial P}\right|}_T\frac{\partial P}{\partial z}+H_L{w}_L{\left.\frac{\partial {\rho }_L}{\partial T_L}\right|}_P\frac{\partial T_L}{\partial z}=\frac{\Delta G_L}{A_{\mathit{\text{pipe}}}\mathit{\Delta z}}$$ E3

$$H_G{\rho }_G\frac{\partial w_G}{\partial z}+{\rho }_G{w}_G\frac{\partial H_G}{\partial z}+H_G{w}_G{\left.\frac{\partial {\rho }_G}{\partial P}\right|}_T\frac{\partial P}{\partial z}+H_G{w}_G{\left.\frac{\partial {\rho }_G}{\partial T_G}\right|}_P\frac{\partial T_G}{\partial z}=\frac{\Delta G_G}{A_{\mathit{\text{pipe}}}\mathit{\Delta z}}$$ E4

2.2. Momentum conservation equation

The law of momentum conservation is a universal principle for any flow system that the varying rate of momentum is equal to the sum of the forces imposed on the control volume. Considering the compressibility of vapor and liquid, the equation of momentum is as follows:

$$\frac{\partial \left({\rho }_{G,L}{w}_{G,L}\right)}{\partial t}+\text{div}\left({\rho }_{G,L}{w}_{G,L}\overrightarrow{u_{G,L}}\right)=-\frac{\partial p}{\partial z}+\frac{\partial {\tau }_{\mathit{xz},G,L}}{\partial x}+\frac{\partial {\tau }_{\mathit{yz},G,L}}{\partial y}+\frac{\partial {\tau }_{\mathit{zz},G,L}}{\partial z}+F_{G,L}$$ E5

here $\mathit{dp}/\mathit{dz}$ is the pressure gradient in the axial direction, $\mathrm{Pa}/\mathrm{m}$ . In the condition of steady flow state, the axial pressure gradient $\mathit{dp}/\mathit{dz}$ in liquid is balanced by the shear stress at wall ${\tau }_{\mathit{WL}}$ and interface ${\tau }_{\mathrm{int}}$ , respectively, in vapor by ${\tau }_{\mathit{WG}}$ and ${\tau }_{\mathrm{int}}$ .

${\tau }_{\mathit{xz},G,L}={\mu }_{m,G,L}\left(\frac{\partial u_{G,L}}{\partial z}+\frac{\partial w_{G,L}}{\partial x}\right)$ and ${\tau }_{\mathit{yz},G,L}={\mu }_{m,G,L}\left(\frac{\partial v_{G,L}}{\partial z}+\frac{\partial w_{G,L}}{\partial y}\right)$ represents viscous stress introduced by molecular viscosity; ${\mu }_m$ denotes molecular dynamic viscosity; $F_{G,L}={\rho }_{G,L}g\text{sin}\theta$ .

According to the theory of CFD, the N-S equation is applicable to any kind of flow. A direct calculation of N-S equation requires a high computer capacity which is not practical in engineering. Hence, an item of Reynolds is introduced, and the final momentum equation for vapor and liquid stratified flow is as follows:

$$\frac{\partial }{\partial x}\left({\Gamma }_{w,G,L}\left(\frac{\partial w_{G,L}}{\partial x}\right)\right)+\frac{\partial }{\partial y}\left({\Gamma }_{w,G,L}\left(\frac{\partial w_{G,L}}{\partial y}\right)\right)=\frac{\mathit{dp}}{\mathit{dz}}-{\rho }_{G,L}g\text{sin}\theta +\frac{\partial {\rho }_{G,L}{w}_{G,L}{w}_{G,L}}{\partial z}$$ E6

where ${\Gamma }_{w,G,L}$ means effective diffusion coefficient, which is the effective viscosity given as the sum of the molecular and eddy viscosity, ${\Gamma }_{w,G,L}={\mu }_{m,G,L}+{\mu }_{t,G,L}$ .

2.3. Turbulent flow model

Following a similar approach as Xiao et al. [36] and Reboux et al. [37, 38] for two-phase turbulent flow, the eddy viscosity is modeled using the large eddy simulation (LES) turbulence model based on the assumption of non-isotropic turbulence. Meanwhile, changes are made to take account of the progressive attenuation of turbulence close to the wall. LES model is applicable to the flow with different Reynolds number. The subgrid scale viscosity is written as:

$${\mu }_{t,G,L}={\rho }_{G,L}{f}_{\text{int},W}{\left(C_S{D}_S\Delta \right)}^2{\left|S\right|}_{G,L}$$ E7

The Smagorinsky constant is $C_S=0.1-0.2$ , while the resolved rate-of-strain tensor is:

$${\left|S\right|}_{G,L}={\left(2S_{\mathit{ij}}{S}_{\mathit{ij}}\right)}_{G,L}^{1/2}={\left[{\left(\frac{\partial w_{G,L}}{\partial x}\right)}^2+{\left(\frac{\partial w_{G,L}}{\partial y}\right)}^2\right]}^{1/2}$$ E8

The filter width $\Delta$ is defined as $\Delta ={\left(\mathit{Vol}\right)}^{1/3}$ , where $\mathit{Vol}$ is the volume of the computational cell. Damping functions have been introduced into the expressions of turbulent viscosity. Near a wall, a wall-damping function is required to perish the eddy viscosity on the wall.

$$D_S=1-\text{exp}\left(\frac{d^{+}}{25}\right);d^{+}=\frac{{\left({\mathit{\rho \tau }}_{\text{int},W}\right)}^{0.5}d}{{\mu }_m};d=\text{min}\left(\frac{D}{2}-\sqrt{x^2+{\left(y-\frac{D}{2}+h_L\right)}^2},\left|d\right|\right)$$ E9

$$f_{\text{int},W}=1-\text{exp}\left[-1.3\times {10}^{-4}{d}^{+}-3.6\times {10}^{-4}{\left(d^{+}\right)}^2-1.08\times {10}^{-5}{\left(d^{+}\right)}^3\right]$$ E10

And $d^{+}$ is the dimensionless distance to pipe wall or vapor-liquid interface. Fulgosi et al. [39] provided an exponential dependence of $f_{\mathrm{int}}$ on the dimensionless distance to the interface or pipe wall. Parameter $d$ means the distance to pipe wall or vapor-liquid interface.

The choice of turbulence model is crucial for this sort of study, due to the presence of these secondary flows. Nallasamy and Rodi explained that the well-known $k-\epsilon$ model assumes an isotropic eddy viscosity, which makes it unsuitable for problems where the anisotropic nature of the turbulent viscosity is important in the calculation of the flow field. This may become significant in stratified gas-liquid flow as the height of the interface approaches the center of the pipe, as stressed by Newton and Behnia. FLUENT offers different options when using this model. Based on the current flow conditions, the linear pressure-strain model and standard wall functions were selected. Turbulence model used by different researchers is shown in Table 1.

Formula name	expression
Duan et al. [5]	Low Reynolds number model: ${\mu }_t=f_{\mu }{C}_{\mu }{k}^2/\epsilon$; $\frac{\partial }{\partial x}\left(\left({\mu }_m+\frac{{\mu }_t}{{\sigma }_k}\right)\frac{\partial k}{\partial x}\right)+\frac{\partial }{\partial y}\left(\left({\mu }_m+\frac{{\mu }_t}{{\sigma }_k}\right)\frac{\partial k}{\partial y}\right)+G_k-\mathit{\rho \epsilon }=0$; $\frac{\partial }{\partial x}\left(\left({\mu }_m+\frac{{\mu }_t}{{\sigma }_{\epsilon }}\right)\frac{\partial \epsilon }{\partial x}\right)+\frac{\partial }{\partial y}\left(\left({\mu }_m+\frac{{\mu }_t}{{\sigma }_{\epsilon }}\right)\frac{\partial \epsilon }{\partial y}\right)+C_1{f}_1{G}_k-C_2{f}_2\frac{\epsilon }{k}=0$; $C_1=1.92$ $C_2=1.3$ $C_{\mu }=0.09{\sigma }_k=1.0$ ${\sigma }_{\epsilon }=1.3$; $f_{\mu }={\left[1-\text{exp}\left(-0.0165{\text{Re}}_y\right)\right]}^2\times \left(1+\frac{20.5}{\text{Re}{}_t}\right)$; $f_1=1+{\left(\frac{0.05}{f_{\mu }}\right)}^3{f}_2=1-\text{exp}\left(-{\text{Re}}_t^2\right)$
Jiang et al. [40]	LES model: ${\mu }_{t,G,L}={\rho }_{G,L}{f}_{\text{int},W}{\left(C_S{D}_S\Delta \right)}^2{\left|S\right|}_{G,L}$

Table 1.

Different turbulent models.

2.4. Heat transfer model

It is of great significance to study turbulent flow and heat-transfer mechanism due to the frequent occurrence in many industrial applications, such as heat exchangers, vapor turbines, cooling systems, and nuclear reactors [27]. In this study, a quasi-steady state of temperature profile is calculated where axial thermal conduction is neglected. Since the operating temperature is influenced by ambient temperature, fluid properties, and hydraulic parameters such as flow rate and liquid level; a heat-transfer model for vapor-liquid flow is established.

Energy equation means the increase of energy, which is equal to the result of heat flux entering the representative-element volume deducting the work from internal force. Ignoring axial heat conduction and viscous dispersive item, the equation of energy conservation is as follows:

$$\frac{\partial }{\partial x}\left({\Gamma }_{T,G,L}\left(\frac{\partial T_{G,L}}{\partial x}\right)\right)+\frac{\partial }{\partial y}\left({\Gamma }_{T,G,L}\left(\frac{\partial T_{G,L}}{\partial y}\right)\right)=\frac{\partial {\rho }_{G,L}{c}_{P,G,L}{w}_{G,L}{T}_{G,L}}{\partial z}$$ E11

${\Gamma }_{T,G,L}$ means effective diffusion coefficient which is the effective heat transfer coefficient and defined as, ${\Gamma }_{\mathrm{T},G,L}=\frac{{\mu }_{m,G,L}{c}_{P,G,L}}{{\text{Pr}}_{m,G,L}}+\frac{{\mu }_{t,G,L}{c}_{P,G,L}}{{\text{Pr}}_{t,G,L}}$ . $\frac{{\mu }_{m,G,L}{c}_{P,G,L}}{{\text{Pr}}_{m,G,L}}$ and $\frac{{\mu }_{t,G,L}{c}_{P,G,L}}{{\text{Pr}}_{t,G,L}}$ are respectively molecular thermal diffusion and eddy diffusion of heat transfer. ${\text{Pr}}_m=c_P\frac{{\mu }_m}{\lambda }$ is molecular Prandtl number. ${\mathit{Pr}}_t$ is turbulent Prandtl number, which is assumed by Jones and Launder as 0.9 [41]. However, Hishida et al. obtained the conclusion that the value of ${\mathit{Pr}}_t$ increases with the distance increase to pipe wall [42]. Only if $d^{+}>30$ , ${\mathit{Pr}}_t$ is 0.9. When $d^{+}=10$ , the value of ${\mathit{Pr}}_t$ rises up to 1.6 according to the measurement of circular tube of turbulent heat transfer experiments [33]. Another calculation of turbulent Prandtl number is given by Kays and Mansoori based on the direct numerical simulation [27]:

$${\mathit{Pr}}_{t,G,L}={\left[0.5882+0.228\left(\frac{{\mu }_{t,G,L}}{{\mu }_{m,G,L}}\right)-0.0441{\left(\frac{{\mu }_{t,G,L}}{{\mu }_{m,G,L}}\right)}^2\left(1-\text{exp}\left(-5.165\right)\left(\frac{{\mu }_{m,G,L}}{{\mu }_{t,G,L}}\right)\right)\right]}^{-1}$$ E12

2.5. Phase change model

The pressure and temperature drop along pipelines, which incurs the change of ratio of gas mass volume and liquid mass volume. It is prone to evoke P–T flash evaporation in each component. The mass percentage would change accordingly, influencing the parameters like molar mass, density, viscosity, heat capacity and thermal conductivity. The vapor enthalpy would change during the process of phase change and dissipate into the vapor or liquid in the form of phase change heat which results in variation of fluid temperature [22, 32, 45].

Cong Guo et al. considered a modification of the original model, in which the rate of conductive heat through the tube wall due to temperature difference can be calculated [43]. Peneloux et al. proposed the concept of volume translation. They argue that the volume obtained using SRKEOS is a “pseudo-volume” and proposed a method of calculating the “pseudo volume” [44]. Sadegh et al. proposed an equilibrium criterion for the Peng-Robinson equation of state (PREOS) based on the volume translated Peng-Robinson equation of state (VTPREOS) and a translated functional relationship is used based on the theory of Peneloux et al. to discuss the volume translation technique [45]. The study of Li Zhang et al. shows that the condensation heat transfer coefficient reduces with the increase of wall subcooling from around 2 to 14°C. With the rise in the wall subcooling, the heat flux increases, resulting in an increasing rate of steam condensation, which brings forth a thicker condensate film on the tube surface. The thicker condensate film around the tube offers a higher thermal resistance to steam condensation and in turn reduces the condensation heat transfer coefficient [22]. Considering the volume addition of LSI phase due to the coalescence, Kai Yan et al. used the additional velocity and considered all the conditions when some portion of SSI phase can come into LSI phase [46]. Bonizzi proposed a model for calculating the atomization flux and the bulk concentration based on the recommendation by Williams et al. and Pan and Hanratty [47, 48, 49, 50, 51]. Vinesh et al. showed the physical model, considered for phase change which corresponds to hydrodynamically as well as thermally developing vapor-liquid stratified flow in a plane channel, with heating from the top and cooling from the bottom wall [52] (Table 2).

Formula name	Expression
Cong Guo et al. [43]	The rate of conductive heat through the tube wall: $Q_2=\mathrm{\Delta }m⌊C_{p,2-2}\left(T_{v,2-2}-T_{f,2-2}\right)+h_{\mathit{fg}}⌋-C_{p,2-2}{m}_{1-1}\left(T_{f,2-2}-T_{f,1-1}\right)$; The temperature of the condensate around circumferential wall: $T_f=T_{w,\mathit{\text{in}}}+0.31\left(T_v-T_{w,\mathit{\text{in}}}\right)$
Peneloux et al. [44]	The “pseudo-volume” obtained using SRKEOS: $v_{\mathit{\text{actual}}}=v_{\mathit{\text{SRKEOS}}}-c$; The definition of “pseudo volume”: $\tilde{V}=V+\sum _{i=1}^n{c}_i{z}_i$
Sadegh et al. [45]	The “pseudo partial volume” can be defined as: ${\tilde{v}}_i={\left(\frac{\partial \tilde{V}}{\partial z_i}\right)}_{T,P,z_j}=v_i+c_{i}j\ne i;i=1,n$; With this definition “pseudo fugacity coefficients” can be defined as: $\text{ln}{\tilde{\varphi }}_i={\int }_0^P\left(\frac{{\tilde{v}}_i}{\mathit{RT}}-\frac{1}{P}\right)\mathit{dP}i=1,n$
Kai Yan et al. [46]	The expression of the additional velocity ${\overrightarrow{\mathbf{U}}}_{\mathit{Add}}$ : ${\overrightarrow{\mathbf{U}}}_{\mathit{Add},s}=\left\{\begin{array}{l}-{\alpha }_3\text{max}\left[{\left({\overrightarrow{\mathbf{U}}}_3-{\overrightarrow{\mathbf{U}}}_m\right)}_s,0\right],\text{if}\frac{\partial {\alpha }_2}{\partial s}>0\\ -{\alpha }_3\text{min}\left[{\left({\overrightarrow{\mathbf{U}}}_3-{\overrightarrow{\mathbf{U}}}_m\right)}_s,0\right],\text{if}\frac{\partial {\alpha }_2}{\partial s}<0\end{array}\right.$
Williams et al. & Pan and Hanratty [47, 48]	The atomization flux: $R_A=\frac{k_A{u}_G^2{\left({\rho }_G{\rho }_L\right)}^{1/2}}{{\sigma }_{\mathit{GL}}}\left({\mathrm{\Gamma }}_Z-{\mathrm{\Gamma }}_{Z,C}\right)$; In which $k_A\approx 2.0\times {10}^{‐6}$ [8, 9].
Bonizzi et al. [51]	The flux of droplet deposition is usually expressed as: $⟨R_D⟩=k_D{C}_B$; $C_B=\frac{{\mathit{EW}}_L}{Q_{G}S}={\rho }_L\frac{{\alpha }_D}{{\alpha }_G}$

Table 2.

Different phase change models.

The equation of state in P–T Flash involves Peng-Robinson (PR) equation, which was proposed by Peng-Robinson in 1976. It can predict the molar volume more accurately than SRK equation and can be applied for polar compounds. Apart from that, the equation is applicable to vapor and liquid at the same time and widely used in the calculation of phase equilibrium.

During the P–T flash calculation process, the required parameters are the component of light hydrocarbon, gasification rate, density, molar mass, and enthalpy. The relation between input and output is explained in Figures 3 and 4.

Figure 3.

Input and output of P–T flash calculation.

Figure 4.

The illustration of the way to obtain the mole fraction by P–T flash in each grid cell.

The fluid in the pipe includes $\mathit{Nc}$ sorts of hydrocarbon, $i=1,2,\dots ,\mathit{Nc}$ . The total energy keeps constant during the phase change process. When vaporization occurs, excess heat of vapor phase is absorbed from the liquid phase. When condensation occurs, excess heat is released from the vapor phase. $L_H$ represents the latent heat of vaporization or condensation, which is calculated by [32]:

$$L_H=\left\{\begin{array}{l}{E}_G-E_L,\text{when}\text{condensation}\text{occurs}\\ E_L-E_G,\text{when}\text{vaporization}\text{occurs}\end{array}\right.$$ E13

Therefore, the latent heat may result in the temperature change of the vapor-liquid system. The enthalpy difference in the process of vaporization or condensation can be obtained by virtue of P–T flash calculation. The latent heat helps to retard the temperature change. The process can be calculated as follows:

$${\int }_{A_{G,L}}{\rho }_{G,L}{w}_{G,L}{L}_H{l}_{\eta }{l}_{\xi }\mathit{\text{dηdξ}}=c_{P,G}{G}_G\Delta T_G+c_{P,L}{G}_L\Delta T_L$$ E14

The parameters of flow state and related dimensionless parameters are given by:

$$G_G={\int }_{A_G}{\rho }_G{w}_G{l}_{\eta }{l}_{\xi }\mathit{\text{dηdξ}},G_L={\int }_{A_L}{\rho }_L{w}_L{l}_{\eta }{l}_{\xi }\mathit{\text{dηdξ}}$$ E15

The properties of fluids are calculated by:

$$\left(e,X,Y,{\rho }_G,{\rho }_L,{\rho }_{\mathit{MIX}},M_G,M_L,M_{\mathit{MIX}},E_G,E_L,E_{\mathit{MIX}}\right)=\mathit{\text{PREOS}}\left(Z,P,T\right)$$ E16

$${\mu }_{m,G}={\mu }_{m,G}\left(T_G,{\rho }_G\right),{\lambda }_G={\lambda }_G\left(T_G\right),c_{P,G}=c_{P,G}\left(P_G,T_G,M_G\right)$$ E17

$${\mu }_{m,L}={\mu }_{m,L}\left(P_L,T_L,M_L,{\rho }_L,X\right),{\lambda }_L={\lambda }_L\left({\rho }_L,T_L\right),c_{P,L}=c_{P,L}\left({\rho }_L,T_L\right)$$ E18

Both gas and liquid phases are not ideal fluids. Vapor phase is a mixture of multi-component light hydrocarbons. Thus, the property of vapor phase is a combination of their quality weighting. The PREOS equation is currently acknowledged as the most accurate formula to calculate the density of mixed vapor. If the pressure, temperature, and relative density of light hydrocarbon are known, the viscosity can be calculated using experimental formulas. Both the thermal conductivity and the heat capacity at constant pressure for vapor mixture are related to temperature and pressure, which also can be calculated using the experimental formulas. For the liquid phase, the density is calculated by PR EOS. The other physical properties including viscosity, thermal conductivity, and specific heat capacity are respectively a function of temperature, density, and other critical properties of components.

The relationship between liquid level $h_L$ and liquid circulation area $A_L$ is as follows:

$$A_{\mathit{\text{pipe}}}{H}_L=A_L=\frac{1}{4}{D}^2\text{arccos}\left(1-\frac{2h_L}{D}\right)-\frac{1}{4}{D}^2\left(1-\frac{2h_L}{D}\right)\sqrt{1-{\left(1-\frac{2h_L}{D}\right)}^2}$$ E19

The value of $H_L$ as well as the calculating method differs in each cross section, as depicted in Figure 6.

As shown in Figure 6, the evaporation fraction, density, and molar weight of vapor and liquid phase can be obtained with the known mole fraction of each component since the temperature and pressure are the same as inlet data in the first cross section. The liquid hold-up can be calculated as follows:

$$H_L=\frac{1}{1+\frac{{\rho }_L{M}_G}{{\rho }_G{M}_L}\frac{e}{\left(1-e\right)}}$$ E20

The temperature, pressure, and mole fraction vary in the second and other subsequent cross section due to the influence of thermal conduction, pressure drop, and phase change but can be derived from the former section. Hence, it requires the calculation of the evaporation fraction, density, and molar weight in each grid cell. The liquid holdup can be calculated as follows:

$$H_L=\frac{\sum _j\frac{{\rho }_{\mathit{MIX},j}{M}_{L,j}}{M_{\mathit{MIX},j}{\rho }_{L,j}}{A}_j\left(1-e_j\right)}{\sum _j\frac{{\rho }_{\mathit{MIX},j}}{M_{\mathit{MIX},j}}\left[\frac{M_{L,j}}{{\rho }_{L,j}}{A}_j\left(1-e_j\right)+\frac{M_{G,j}}{{\rho }_{G,j}}{A}_j{e}_j\right]}$$ E21

Then, the liquid level $h_L$ can be obtained by Eq. (19).

2.6. Wavy gas-liquid interface model

$\mathit{\Delta y}$ is effective level shift of the gas-liquid interface. It increases the length scale of turbulence, which further increases the shear stress near the rough pipe wall. The shift is caused by the increase of flow resistance by the rough interface, which varies the velocity distribution at the interface.

Dimensionless shift at wavy gas-liquid interface is defined as:

$$\Delta y_{\text{int},G,L}^{+}={\rho }_{G,L}\Delta y_{\text{int},G,L}{w}_{\tau ,\text{int},G,L}/{\mu }_{G,L}$$ E22

It can be calculated by (Cebeci and Smith, 1974):

$$\Delta y_{\text{int},G,L}^{+}=0.9\left(\sqrt{{\epsilon }_{\text{int},G,L}^{+}}-{\epsilon }_{\text{int},G,L}^{+}{e}^{-{\epsilon }_{\text{int},G,L}^{+}/6}\right),4.54<{\epsilon }_{\text{int},G,L}^{+}<2000$$ E23

The abovementioned correlation is based on the data fit, the lower limit corresponds to a smooth surface, where $\Delta y_{\text{int},G,L}^{+}\approx 0$ . ${\epsilon }_{\text{int},G,L}^{+}$ is the dimensionless equivalent roughness, expressed by:

$${\epsilon }_{\text{int},G,L}^{+}={\rho }_{G,L}{\epsilon }_{\mathrm{int}}{w}_{\tau ,\text{int},G,L}/{\mu }_{G,L}$$ E24

The ${\epsilon }_{\mathrm{int}}$ and $u_{\tau ,\text{int}}$ are calculated by:

$$-3.6{\text{log}}_{10}\left[\frac{6.9}{{\text{Re}}_G}+{\left(\frac{{\epsilon }_{\mathrm{int}}}{3.7D_G}\right)}^{1.11}\right]=\sqrt{\frac{1}{f_{\mathrm{int}}}}$$ E25

$$w_{\tau ,\text{int},G,L}={\left(\overline{{\tau }_{\mathrm{int}}}/{\rho }_{G,L}\right)}^{1/2}$$ E26

where ${\text{Re}}_G={\rho }_G{u}_G{D}_G/{\mu }_G$ , $D_G=4A_G/\left(S_G+S_{\mathrm{int}}\right)$ . $S_G$ is the gas phase wetted perimeter length; m. $f_{\mathrm{int}}$ and $\overline{{\tau }_{\mathrm{int}}}$ are calculated as:

$$f_{\mathrm{int}}=\frac{2\overline{{\tau }_{\mathrm{int}}}}{{\rho }_G{\left(\overline{w_G}-\overline{w_L}\right)}^2}$$ E27

$$\overline{{\tau }_{\mathrm{int}}}=\overline{{\tau }_{\text{int},G,L}}=\frac{1}{a}{\int }_0^a{\Gamma }_{w,\text{int},G,L}\frac{l_{\eta }}{l_{\xi }}{\left.\frac{\partial w_{G,L}}{\partial \xi }\right|}_{\xi =\pi }\mathit{d\eta }$$ E28

The flow geometry of stratified flow in circular receiver is very complex. In order to simplify the model, many researchers use different models in their researches. The following are some models.

2.7. Boundary conditions

The boundary conditions included in vapor-liquid two-phase stratified pipe flow with heat transfer and phase change involve vapor-liquid interface condition, wall boundary condition, symmetrical boundary condition and inlet boundary condition. As for the vapor-liquid interface, in order to illustrate the mutual influences among flow, heat transfer and phase change, the equal interfacial shear stress has been prescribed as the vapor-liquid interface condition which is related to the fluid properties and velocity distribution and can be calculated as:

$${\left.{\tau }_{\text{int},G}={\left({\mu }_m+{\mu }_t\right)}_G\frac{\partial w_G}{\partial n}\right|}_{\mathrm{int}}={\left.{\left({\mu }_m+{\mu }_t\right)}_L\frac{\partial w_L}{\partial n}\right|}_{\mathrm{int}}={\tau }_{\text{int},L}$$ E29

The temperature and heat flux of two phases are respectively equal at vapor liquid interface (Table 3).

$$q_{\text{int},G}=q_{\text{int},L},T_{\text{int},G}=T_{\text{int},L}$$ E30

Table 3.

Different gas-liquid interface models.

At the pipe wall, the non-slip condition is applied for velocity in both two phases.

$$w_{W,G,L}=0$$ E31

The temperature boundary condition at pipe wall of vapor and liquid phases are convective heat transfer and its coefficient of pipeline outer wall remains constant.

$${\alpha }_{G,L}\left(T_W-T_{G,L}\right)={\lambda }_{G,L}{\left.\frac{{\mathit{dT}}_{G,L}}{\mathit{dn}}\right|}_W$$ E32

The gradients of velocity and temperature are zero at the symmetrical boundary.

$${\left.\frac{\partial w_{G,L}}{\partial n}\right|}_{\mathit{\text{symmetrical}}}=0,{\left.\frac{\partial T_{G,L}}{\partial n}\right|}_{\mathit{\text{symmetrical}}}=0$$ E33

At the inlet of pipeline, the velocity and temperature filed of two phases are respectively equal to the pipe inlet values.

$$w_{G,L}=w_{\mathit{\text{inlet}},G,L},T_{G,L}=T_{\mathit{\text{inlet}}}$$ E34

3. Results and discussion

Based on the theory of flow and heat transfer, turbulent flow and phase equilibrium, the model is solved by multi-physical field coupling numerical simulation. The non-circular liquid and vapor domains in stratified pipe flow can be simply modeled with the bipolar coordinate system, which is helpful in solving the problem caused by the inhomogeneity of boundaries. Bipolar cylindrical coordinate is composed of two orthogonal circles in rectangular coordinate. As the flow field in both phases is bounded by a circular pipe wall and a plane interface, the calculation domain has been converted to rectangle form from the anomalous physical domain by adopting the bipolar coordinate system.

With the increase of axial distance, the liquid level in pipeline changes constantly and leads to the change of flow area in both phases. The grid size changes adaptively along with the flow area, where the flow area is determined by the height of gas-liquid interface. Variable-size grid has advantages in calculating the changing interface. The grid number remains unchangeable. The location of the gas-liquid interface is obtained by the secant method, and the convergence condition is that the conservation of the mass flow rate of the gas-liquid phase and the total mass flow rate equal to the inlet mass flow rate. In this way, the interface is detected.

Vapor-liquid two phase flow and heat transfer coupled with phase change have been simulated in this section. The simulated pipeline is with inner diameter of 100 mm and total length of 6000 m. The superficial velocities of vapor and liquid are respectively $U_{\mathit{SG}}=6.223\mathrm{m}/\mathrm{s}$ and $U_{\mathit{SL}}=0.016\mathrm{m}/\mathrm{s}$ . The pressure gradient is about $\mathit{dp}/\mathit{dz}=-21\mathrm{Pa}/\mathrm{m}$ and the liquid holdup at pipe inlet is $H_L=0.0343$ . The mass flow rate of the fluid at pipe inlet is 5.892 kg/s. The fluid of vapor-liquid mixture at pipe inlet also has a pressure of 10.343 MPa and a temperature of 48°C. The convective heat transfer coefficient is about 10 $\mathrm{W}/\left({\mathrm{m}}^2\cdot \mathrm{K}\right)$ and the ambient temperature is 15°C. The simulated pipeline is divided into 600 segments where the cross-sectional area of each segment has a mesh grid scale of $84\times 122$ for numerical calculation. The fluid is a mixture of multi-component hydrocarbons, which are shown in Table 4.

Compound	CH4	C2H6	C3H8	i-C4H10	n-C4H10	i-C5H12	n-C5H12	n-C6H14	${\mathbf{C}}_{\mathbf{7}}^{+}$
Mole percent (%)	78.03	4.73	5.98	3.05	3.54	2.85	0.54	0.69	0.59

Table 4.

Chemical composition of the light hydrocarbons used in current study.

3.1. Mole fraction, density distribution, and liquid level along the pipeline

The mole fractions of each component, density distribution, temperature distribution, and liquid level along the pipeline are obtained in the condition of condensation production.

Mole fractions of each component at pipe inlet in both phases are shown in Table 5. In vapor phase, the mole fraction of methane is larger than all the other light hydrocarbons ( ${\mathrm{C}}_2^{+}$ ). In liquid phase, the mole fractions of the other light hydrocarbons are bigger than that in vapor phase, but the methane is still the main component.

Compound	CH4	C2H6	C3H8	i-C4H10	n-C4H10	i-C5H12	n-C5H12	n-C6H14	${\mathbf{C}}_{\mathbf{7}}^{+}$	Total
Mole fraction in vapor phase (%)	80.54	4.69	5.61	2.71	3.05	2.25	0.41	0.44	0.30	100
Mole fraction in liquid phase (%)	43.31	5.35	11.05	7.82	10.30	11.21	2.32	4.10	4.54	100

Table 5.

Mole fraction of each component at pipe inlet in both phases.

Mole fractions of each component in both vapor and liquid phase are shown in Figure 5. The content of methane in vapor phase increases when flowing in the pipe while the content of the other light hydrocarbons become less and less in vapor phase. During the condensing process which is dominated by temperature drop, the methane keeps evaporating; while during the evaporating process which is dominated by pressure drop, the other light hydrocarbons keep condensing. Meanwhile, the bigger the molar mass is, the faster the condensing rate is.

Figure 5.

Mole fraction of each component changes along pipeline. (left) vapor phase; (right) liquid phase.

Six pipe cross sections located at every 1000 m along the pipeline are selected to illustrate the change of density and temperature distribution, as shown in Figures 6 and 7. It is exactly because the pressure on the cross section is the same, so the uneven distribution of the temperature leads to the uneven distribution of the fluid density. In the two phases, the density distribution is opposite to the temperature distribution. The high temperature means low density. The temperature value distributed at every single cross section can be ranked in descending order: the interior of liquid phase, most parts in vapor phase, vapor phase near the top wall. However, the descending rank of density is liquid phase, vapor phase near the top wall, and the interior of vapor phase. Along the pipeline, the temperature in the area referenced above decreases gradually while the density increases gradually. Therefore, the temperature of the sequential cross sections tends to be the same and the density distribution within the two phases gradually becomes uniform.

Figure 6.

Density (dimensionless) distribution depicted in 6 pipe cross sections along the pipeline (every 1000 m).

Figure 7.

Temperature (K) depicted in 6 pipe cross sections along the pipeline (every 1000 m).

Figure 8 illustrates the selected 12 pipe cross sections located at every 500 m along the pipeline, where the varying trend of liquid level are presented in three-dimensional coordinate system. The minimum liquid level is 7.82 mm at pipeline inlet. The liquid level at outlet is about 16.18 mm and keeps declining trend, which can also be found in Figure 9(b).

Figure 8.

Liquid level depicted in 12 pipe cross sections along the pipeline (every 500 m).

Figure 9.

Pressure gradient, liquid level, fluid mass flow rate, and temperature along the pipeline. (a) pressure gradient, (b) liquid level, (c) mass flow rate, and (d) bulk temperature.

3.2. Pressure gradient, liquid level, fluid mass flow rate, and temperature along the pipeline

The pressure gradient, liquid level, fluid mass flow rate, and temperature along the pipeline obtained in the condition of condensation production have been compared with that found in literature.

When the phase change behavior is considered along pipe flow, the vapor in gas phase starts to condense to liquid which begins from the pipe inlet due to the significant temperature drop at pipe wall. The pressure drop fit well with the simulated results presented by Sadegh, as shown in Figure 9(a) [16].

During the condensing process, vapor mass flow rate gradually reduces, as shown in Figure 9(c), and the liquid mass flow rate increases due to the constant total mass flow rate. The increase of liquid mass flow rate leads to further rise of liquid level, as shown in Figure 9(b).

The liquid holdup firstly increases until it reaches a maximum value and then gradually decreases. The reason behind this is as follows: The increase of liquid holdup results from the liquid precipitation caused by dominant temperature drop. Due to the large difference between fluid temperature and ambient temperature, the amount of liquid precipitation is greater than liquid evaporation. On the contrary, the decrease of liquid holdup is leaded by liquid evaporation due to dominant pressure drop. Being same to liquid holdup, the liquid mass flow rate maintains the same trend, that is, gradually increasing to reach a maximum value and then gradually reduced, which is depicted in Figure 9(c). As the total mass flow is constant, the mass flow rate of the vapor phase decreases first and then increases. When compared with the process of evaporation, the precipitation process caused by temperature drop is transient and intense, which is related to the temperature difference between the inside and outside of the pipeline and to the convective heat transfer coefficient.

The tendency of temperature drop is similar to that in the existing research [16]. But there also exists difference between vapor bulk temperature, liquid bulk temperature, and the total bulk temperature, which cannot be revealed by one-dimensional model. The liquid bulk temperature is always higher than the vapor bulk temperature while the vapor bulk temperature is almost always equal to the total bulk temperature. Latent heat is revealed during the vapor condensing process which slows down the temperature drop, as shown in Figure 9(d).

3.3. Velocity and temperature distribution at pipe length of 3000 m

Through solving the model, with phase change happening, it can be obtained that the pressure gradient is 21.29 Pa/m and the liquid level is 16.02 mm when axial distance reaches 3000 m.

Figure 10(a) shows that the velocity of vapor phase slows down while approaching either the pipe wall or the vapor-liquid interface because of the hindering effects and fluid viscosity. The velocity of liquid phase keeps increasing from pipe wall to the interface. The maximum velocity at the pipe cross section occurs within the vapor phase.

Figure 10.

Velocity and temperature distribution at pipe length of 3000 m. (a) Velocity profile at pipe cross section; (b) temperature profile at pipe cross section; (c) velocity profile at vertical centerline; and (d) temperature profile at vertical centerline.

In Figure 10(b) shows that the temperatures of both vapor and liquid phase drop while approaching the pipe wall because of the lowest ambient temperature and the convective heat transfer effects. The temperature of liquid phase keeps increasing from pipe wall to the interface. The maximum temperature at the pipe cross section exists within the liquid phase near the interface.

Figure 10(c) shows that the temperature at pipe wall of vapor phase is lower than that of liquid phase when convective heat transfer exists due to the smaller heat carried by the vapor phase than liquid phase. Thus, lower specific heat capacity results in bigger temperature drop at pipe wall of vapor phase. The thermal conductivity of the liquid phase is greater than vapor phase, hence, the temperature gradient in liquid phase is smaller than that in vapor phase, and the bulk average temperature of the liquid phase is higher than the vapor phase. Heat is transferring from liquid phase to vapor phase through the interface, which makes the temperature drop of the liquid phase and reduces the temperature difference between the two phases.

Figure 10(d) reveals the velocity distribution at the centerline of the pipe. By the dragging force of the interface, the velocity of liquid phase reaches the maximum value at the interface while the velocity of vapor phase reaches the maximum value at the location between the interface and its bulk center. The liquid phase slows down while approaching the pipe wall because of the hindering of the pipe wall and its high viscosity.

4. Conclusion

The vapor-liquid two-phase pipe flow and heat transfer are studied by virtue of numerical simulation in light hydrocarbon transportation pipeline coupled with hydraulics, thermodynamics, and phase change. A three-dimensional non-isothermal vapor-liquid stratified flow model including phase change model in bipolar coordinate system has been established, where LES turbulence model is utilized to simulate the turbulence flow and the wall attenuation function is used to describe the inadequacy performance of vapor-liquid interface. The vapor phase and the liquid phase are both considered to be compressible and the PR equation of state is chosen for the vapor-liquid equilibrium calculation where the multi-component hydrocarbon flash calculation is used to evaluate the physical properties, gasification rate, and enthalpy departure of the phases. The P–T flash calculation has been applied to predict the varying liquid level and the multi-component mass fraction in each phase during the process of vapor/liquid stratified pipe flow. The axial pressure gradient, liquid holdup, velocity, and temperature fields have been presented. The fluid mass flow rate, mole fraction, density distribution, and liquid level along the pipeline are also given out.

The simulation results indicate that the influence of pressure and temperature on liquid holdup is different. During the light hydrocarbon transportation process in pipeline, the temperature drop leads to the reduction of vapor mass flow rate and the rise of liquid level as well as mass flow rate. Larger temperature drop results in bigger liquid holdup while larger pressure drop causes smaller liquid holdup due to the change of physical properties and phase equilibrium. After the increase of liquid holdup caused by dominant temperature drop reaching the maximum value, then the decrease of liquid holdup maintains its trend till the pipe outlet, which results from liquid evaporation due to dominant pressure drop.

The highest velocity locates in vapor phase while the highest temperature locates in liquid phase. The liquid bulk temperature is always higher than the vapor bulk temperature. The vapor bulk temperature is almost always equal to the total bulk temperature, which cannot be revealed by one-dimensional model. Latent heat is revealed during the vapor condensing process which slows down the temperature drop. The average velocity of liquid is lower than that of vapor, but the temperature of liquid is higher than vapor.

When the fluid flows in the pipeline, the content of methane in vapor phase increases all the time while the content of the other light hydrocarbons ( ${\mathrm{C}}_2^{+}$ ) become less and less. The condensing rates of the other light hydrocarbons have a positive correlation with their molar mass. The temperature value distributed at every single cross section can be ranked in descending order: the interior of liquid phase, most parts in vapor phase, vapor phase near the top wall. However, the descending rank of density is liquid phase, vapor phase near the top wall, and the interior of vapor phase. Along the pipeline, the temperature in the area referenced above decreases gradually while the density increases gradually. Therefore, the temperature of the sequential cross-sections tends to be the same and the density distribution within the two phases gradually becomes uniform. As for the varying trend of liquid level presented in three-dimensional coordinate system, it has the same trend of the liquid holdup which firstly increases until it reaches a maximum value and then gradually decreases.

Thus, models in this chapter can be utilized to accurately predict pressure gradient, velocity, temperature field, liquid holdup, fluid physical properties, and mole fraction, which are essential to the determination of pipe size, design of downstream equipment, and guarantee of flow assurance.

Acknowledgments

This work was supported by the National Natural Science Foundation of China [Grant number 51474228]; and the Beijing Scientific Research and Graduate Joint Training Program [Grant number ZX20150440].