Modeling phase equilibrium with a modified Wong-Sandler mixing rule for natural gas hydrates: Experimental validation

Applied Energy 205 (2017) 749–760

Contents lists available at ScienceDirect

Applied Energy journal homepage: www.elsevier.com/locate/apenergy

Modeling phase equilibrium with a modified Wong-Sandler mixing rule for natural gas hydrates: Experimental validation Niraj Thakre, Amiya K. Jana

MARK

⁎

Energy and Process Engineering Laboratory, Department of Chemical Engineering, Indian Institute of Technology, Kharagpur, 721 302, India

H I G H L I G H T S equilibrium models are developed for natural gas hydrates. • Phase modification is proposed in Wong-Sandler (mWS) mixing rule. • AProposed are validated with experimental data. • Patel-Tejamodels model with mWS shows the best performance. •

A R T I C L E I N F O

A B S T R A C T

Keywords: Phase equilibrium modeling Patel-Teja and Peng-Robinson EoS Modified Wong-Sandler mixing rule Natural gas hydrates Thermodynamic promoter Experimental validation

As natural gas hydrate forms at low temperature and high pressure (i.e., near critical region), vapor phase fugacity calculation becomes crucial in this region due to the intermolecular interactions. This work introduces the modified form of the Wong-Sandler mixing rule in the Peng-Robinson (PR) and Patel-Teja (PT) equation of state (EoS) models to revamp the fugacity calculation of vapor phase. This mixing rule which is based on the excess Gibbs free energy leads to precisely capture the polar and asymmetric properties of mixture components. The modified Wong-Sandler mixing rule is formulated by representing its correction factors for energy and covolume parameters in terms of temperature and gas-phase composition. Here, the hydrate phase fugacity is determined by the Chen-Guo model that is coupled with an EoS model for hydrate equilibrium analysis. On the other hand, the liquid phase nonideality is estimated by using the Wilson activity-coefficient model. The proposed modified Wong-Sandler based phase equilibrium models are shown to be better than the existing phase rule based models with the three example natural gas hydrate systems (CH4 + N2, CH4 + C2H6, CH4 + C2H4) using tetrahydrofuran as a thermodynamic promoter. Validating these models with real time data sets, it is further investigated that the PT based model shows a better performance compared to the PR equation of state. This apart, these two models show their superiority over the existing PT-van der Waals (PT-vdW) model, for which, a new parameter set is proposed for the said three systems for their improved performance.

1. Introduction The world is prone to the shortage of energy urging an enormous amount of fuel in upcoming decades. It is expected that the conventional fuels will continue to supply 80% of the fuel demand that could be liable to another major issue concerning global warming [1,2]. As a bridge fuel between conventional and green energy sources, the natural gas hydrates (NGH) could provide a promising solution for these mutually associated problems of fuel demand and global warming [3,4]. This could be achieved by replacing the methane (CH4) encaged in the hydrates at sea-floor with the greenhouse gases like carbon dioxide (CO2) [5,6]. Besides, the technique of entrapment of huge amount of

⁎

gases into a compact form in the hydrates have a major application in the field of gas storage and transportation [7]; separation [8], desalination of water [9], CO2 sequestration [10], refrigeration [11] and so forth. The gas hydrates are first discovered in 1811 as an inclusion compound in which the water molecules form a known ice-like crystalline structure. This tricky solid solution makes itself thermodynamically stable by encapsulating the small gas molecules into the cavities formed by water molecules [1,6]. The nature of guest-gas molecules decides the structure of gas hydrate that are termed as structure I (sI), structure II (sII) and structure H (sH). The details about structure of gas hydrate and mechanism of hydrate formation are well known [4].

Corresponding author. E-mail address: [email protected] (A.K. Jana).

http://dx.doi.org/10.1016/j.apenergy.2017.08.083 Received 1 February 2017; Received in revised form 7 August 2017; Accepted 11 August 2017 0306-2619/ © 2017 Elsevier Ltd. All rights reserved.

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

wide application in pure and mixture gas hydrate systems [13,25–29]. Modeling phase equilibrium for mixture-gas hydrates requires the appropriate mixing rule to determine the mixture parameters that are dependent on the individual component property. The most common mixing rule which has been used for gas hydrate systems is the van der Waals one-fluid type mixing rule [29]. This is widely applicable for non-polar mixtures but it fails to determine the properties of polar and asymmetric component mixtures [30]. The other category of mixing rules, which are made to capture the polar mixture properties, are based on excess Gibbs free energy (GE) and commonly known as GEtype mixing rules. They include the Huran-Vidal [31], Kurihara [32], modified Huran-Vidal [33], linear combination of Vidal and Michelson (LCVM) [34] and the Wong-Sandler [30] mixing rules. Amongst, the Wong-Sandler (WS) mixing rule has gained much attention for its accuracy in prediction of phase equilibrium for strongly nonideal mixtures. The Helmholtz free energy considered in the WS mixing rule ensures the quadratic dependence of second virial coefficient on the composition of gas components as required by statistical thermodynamics [30]. Also, the Helmholtz free energy at infinite pressure is approximated as Gibbs free energy at zero pressure [30]. It is worth noticing that the GE-type mixing rules are immensely used in VLE calculation [35,36] but not in a common practice for predicting gas hydrate formation conditions [37]. As a brief look into the literature for predictions of the NGH formation conditions, the thermodynamic models have undergone successive modifications over the basic models [38]. Chen and Guo [23] have presented the new approach to model the phase equilibrium of gas hydrates using the concepts reported in their previous work [23]. They have illustrated their model for natural gas hydrate having hydrocarbons of different chain length. The predictions for hydrate formation temperature are found in a good competence with the prediction made by Barken and Shenin [23]. Subsequently, Klauda and Sandler [24] have proposed a fugacity based model using the basic of van der Waals and Platteeuw (vdW-P) model and identified the parameters for natural gas hydrates in the THF solution [39]. They have reported the improved predictions with their model over the classical vdW-P model [39]. The extension of the Chen-Guo model for ternary system of sour natural gas hydrate is presented by Sun et al. [40]. Further as a modification in mixing rule for vapor phase fugacity calculation, Ma et al. [32] have configured the Chen-Guo model with the PT equation of state that was complied with Kurihara mixing rule. The model is illustrated with different inhibitor containing systems showing better prediction of gas solubility over the model using the van der Waals mixing rule [32]. It should be noted that the more accurate GE-type Wong-Sandler (WS) mixing rule is rarely applied to the prediction of gas hydrate formation conditions. Djavidnia et al. [37] have demonstrated the WS mixing rule with the PR-EoS for CH4/CO2/CH3OH system. However, based on our knowledge, there is no work published concerning the use of the WS mixing rule in THF based hydrates. Furthermore, the WS mixing rule with the three-parameter PT-EoS has not been formulated so far for any gas hydrate system.In this context, we introduce a modification in the classical Wong-Sandler mixing rule. The proposed approach concerns the formulation of second virial coefficient as it is vital for calculating the mixture parameters introducing correction factors, namely kij and lij , for the mean-value approximation of the mixture parameters, namely energy (a ) and co-volume (b ) parameter respectively. Those factors are expressed as a linear function of temperature and vapor phase composition. This modified WS mixing rule is then compiled with the PR and PT equation of state models and proposed to use in the Chen-Guo model for phase equilibrium analysis of the gas hydrates. To validate these models, we use experimental data sets of the CH4/N2/ THF (System I) [27], CH4/C2H6/THF (System II) [28], and CH4/C2H4/ THF (System III) [41]. Showing the superiority of the proposed modified (mWS) mixing rule over the existing WS and classical vdW in terms of average absolute relative deviation (AARD), it is investigated that the mWS based PT and PR equation of state models outperform the

In general, the formation pressures of fuel gas hydrates are very high and it is not desirable in concern with the safety and economy of the fuel gas storage and transportation system. Some organic compounds such as tetrahydrofuran (THF), cyclopentane, aldehydes, ketones, and so on, lead to form the gas hydrate at atmospheric pressure and moderate temperature range [12,13]. These compounds dramatically lower down the gas hydrate formation pressure when added to the liquid phase and hence, are termed as thermodynamic promoters. Among the various categories of promoters, the THF solution is reported to have highest stabilization effect [12]. Pahlavanzadeh et al. [13] have demonstrated the effect of different promoters for methane gas hydrates and found the tremendous reduction in hydrate formation pressure with THF over other promoters. On the other hand, there are some organic and inorganic compounds that lead to elevate the gas hydrate formation pressure when present in the system and hence, termed as inhibitors [14]. For instance, the salt content in sea water inhibits the hydrate formation [15]. A large category of promoters and inhibitors are utilized in the actual application of the gas hydrate technique to fairly tune up the gas hydrate formation or dissociation conditions [16]. As the natural gas is relatively clean-burning fuel, the researchers have paid much attention towards developing a model for NGH energy generation unit [2–6,16]. The aim is to develop a model for bulk phase prior to the real system that can further accommodate the inherent properties of saline and sediment environments [17,18]. There are certain challenges associated with the efficient operation and automation of the NGH energy generation system, storage and transportation of fuel gases. In order to overcome these challenges, the time-independent (thermodynamic) hydrate quantification is focused [4]. However, the advancements in thermodynamic models for accurate prediction of the formation condition of gas hydrates are lacking in the literature as compared to the abundant experimental data for phase equilibrium of natural gas hydrates. The basic thermodynamic model is proposed by van der Waals and Platteeuw, which uses statistical thermodynamics to derive the hydrate phase model [19]. However, the thermodynamic predictions are crucial in case of gas mixtures due to the hydrate structure transition [20,21]. The modern spectroscopy enabled more accurate predictions for gas mixtures through the bridge of statistical thermodynamics. Uchida et al. [22] have examined the methane and ethane hydrate properties using Raman spectroscopy, X-ray diffraction and gas chromatography to divulge the mechanism of structure transitions and the cage occupancy through Gibbs free energy calculation. Chen and Guo [23] have proposed a more realistic approach to thermodynamic modeling of gas hydrate. The local stability is considered in their model with a better estimation of coordination number. Meanwhile, Klauda and Sandler [24] have derived a fugacity-based model replacing the reference energy parameter calculation and thus reducing the empirical relations. However, the Chen-Guo model has better fascinated the researchers with the simplified calculation of chemical potential and Langmuir constants [23].The phase equilibrium conditions are obtained by equating the hydrate phase fugacity calculated in the Chen-Guo model to the vapor phase fugacity. Since the hydrate formation occurs above the critical point, the ϕ−ϕ approach (i.e., equation of state (EoS)) is common for predicting the vapor-liquid equilibria (VLE) and frequently used in vapor phase fugacity calculation. Among various equations of state, the Peng-Robinson (PR) [13,25] and the Patel-Teja (PT) [26,27] models are extensively adopted. The PR equation of state assumes fixed value of critical compressibility factor that results in deviation of predicted values of densities of the saturated liquid and critical volumes from their experimental values [26]. In order to address these issues, Patel and Teja [26] have incorporated two additional parameters along with critical temperature and pressure to ensure that the energy parameter approaches a zero value (a → 0 ) at high temperatures. In addition, they gave critical compressibility factor as an empirical dependency on the acentric factor. The PR-EoS and PT-EoS both have 750

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

of the formation condition. It can be noted that even a small concentration of THF dominates in deciding the final structure of gas hydrate by forcing the other guest gas molecules to form the sII hydrate that are supposed to make the sI hydrate at pure state. In general, the compounds that form gas hydrate at relatively normal conditions (high temperature and low pressure) lead to decide the hydrate structure in a mixture gas hydrate. The cavities that fabricate the sII hydrate have two different sized structures: 51264 (formed by six hexagons and twelve pentagons) and 512 (formed by twelve pentagons). The number of water molecules in the sII hydrate is 136 molecules per unit cell with number of small (512) and large (51264) cavities as 16 and 8, respectively. The large and small cavities are treated as the basic and linked cavities, respectively. A two-step model was proposed for the formation of gas hydrate by Chen and Guo [23]. The first step involves the quasi-chemical reaction, in which, the basic hydrate (51264) is formed. This quasi-chemical reaction for sII hydrate formation can be represented as:

existing PT-vdW model, for which, parameter sets are either reproduced (System I and II) or newly proposed (System III) when the same is not available in the literature. 2. Thermodynamic phase equilibrium modeling In this work, we have adopted the Chen-Guo model [23] to predict the phase equilibrium of hydrate phase in terms of fugacity, i.e. fih , with the two other phases, namely liquid ( fil ) and vapor ( f iv ). The presence of thermodynamic promoter (THF) in liquid phase imparts nonideality that is evaluated by the Wilson model [29]. In the vapor phase, the intermolecular interaction is liable to the nonideal behavior of natural gas [42]. The unsaturated hydrocarbon in the gas mixture results in the deviations from ideal behavior owing to their polarity and asymmetric property. These sort of vapor phase nonidealities are registered in terms of f iv , that is calculated with two equation of state models, namely PR and PT. For the evaluation of mixture parameters in the said equation of state models, the classical and modified Wong Sandler mixing rules (WS and mWS) are proposed to hybridize with both the PR and PT yielding the PT-WS, PR-WS, PT-mWS and PR-mWS models. Moreover, the van der Waals mixing rule is also coupled with the PT model to evaluate the proposed PR-mWS and PT-mWS with reference to the PT-vdW. The vapor-liquid-hydrate equilibrium (VLHE) can be defined by

f iv = fil = fih

∑

8x i∗ Ci + 136H2 O→ 8(∑ x i∗ Ci )·136H2 O

i

i

(2)

x i∗

is the mole fraction of the component Ci in the basic hydrate. where The basic cavities are fully occupied by gas molecules in the first step. On the other hand, the linked cavities (512) are generated simultaneously in the same step. Subsequently, the second step proceeds with the adsorption of guest molecules into linked cavities that make the hydrate a stable compound as well. It is interesting to mention here that the linked cavities are not fully filled with gas molecules and thus, making the gas hydrates non-stoichiometric [29]. The occupancy of guest molecules depends on the hydrate formation conditions and can be described by Langmuir adsorption theory.The fugacity of i th component in hydrate phase can be calculated by the Chen-Guo model as

(1)

A complete theoretical background of the model equations used in this work can be found in the concerned literature. The Chen-Guo model for hydrate phase fugacity calculation is derived from statistical mechanical model of van der Waals and Platteeuw [23]. The Patel-Teja equation of state for fugacity calculation for fluid phases is a modified form of Peng-Robinson and Soave-Redlich-Kwong equation of state featuring the substance dependent compressibility factor [26]. The Wong-Sandler mixing rule based PT equation of state model is taken from Yang et al. [35]. The expressions for the differential terms used in this model are obtained by solving the fugacity equation. The thermodynamics of hydrate, vapor and liquid phase is briefly described underneath.

fih = x i∗ fio (1− ∑ θj )α j

(3)

where α is the ratio of linked to basic cavities, which is equal to 1/3 and 2 for sI and sII hydrates respectively. θj is the fractional occupation of linked cavities by j th component and it can be calculated as

θj =

2.1. Hydrate phase

cj f jv 1 + ∑j cj f jv

(4)

in which, cj is the Langmuir constant which can be evaluated from the Lennard-Jones potential function:

To represent the thermodynamic equilibrium of hydrate phase, the Chen-Guo model is very effective. The basic assumptions made in this model for gas hydrate formation mechanism are provided elsewhere [23]. Herein, we are mentioning the key assumptions that can be liberated in the future works:

Yj ⎞ cj = Xj exp ⎛⎜ ⎟ T − Zj ⎠ ⎝

(5)

• The solubility of gas in water is very low and thus, the liquid phase nonideality due to solubilized gas can be neglected. • The sensible heat involved in quasi-chemical reaction to form the

Here, Xj , Yj and Zj are the Antoine constants that are determined against the Langmuir constant (cj ).The fio in Eq. (3) stands for the fugacity in vapor phase, which is in equilibrium with the empty basic hydrate. This term is derived as a function of temperature (T ), pressure ( p ) and activity of water (a w ) and can be represented by the following expression:

•

βp fio = fio,T exp ⎛ ⎞ a w−1/ λ2 ⎝T ⎠

• • •

basic hydrate is assumed to be very low and thus the temperature change in the gas hydrate formation is neglected. Irrespective of the mode of occupancy of the guest molecules, the free energy contribution of the water molecules are considered as the same. The molar volume difference in gas hydrate formation is taken as constant. The cavities formed by water molecules are assumed to be perfectly spherical and thus, the interaction of guest gas with those molecules can be derived as pair potential function. The molecular interactions from basic to linked cavities are considered whereas those interactions in basic hydrate as well as in linked cavities are neglected.

(6)

where β equals 4.242 K/MPa and 10.224 K/MPa for sI and sII hydrates, respectively. λ2 is estimated as the number of basic cavities per water molecule and it is 3/23 and 1/17 for sI and sII hydrates respectively. It should be noted that the activity of water is evaluated here by the Wilson model (Section 2.4). The term fio,T in Eq. (6) is a function of temperature, which has the following form for mixture gas hydrate:

⎛ − ∑j Aij θj ⎞ ⎡ ⎛ B′ ⎞ ⎤ fio,T = exp ⎜ ⎟ × ⎢A′exp T −C′ ⎥ T ⎠⎦ ⎝ ⎝ ⎠ ⎣

In the present study, we consider the systems with THF as a thermodynamic promoter which always forms the sII hydrate irrespective

(7)

This equation resembles an Antoine type equation, in which the constants A′, B′ and C′ are determined by fitting the model to the pure gas 751

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

hydrate experimental data. As most of the compounds in pure form make only a specific structure of gas hydrate, the double gas hydrate data is used here for determining the constants for other structures. For this purpose, the mixture with a specific component is taken, which forms another type of structure. Further, Aij is the binary interaction parameter that specifies the interaction between the guest-gas molecules. Note that Aii = Ajj = 0 and Aij = Aji . In order to predict the gas hydrate formation condition, we here follow a search algorithm (Section 3) using which the pressure is computed at a particular temperature. The hydrate phase concentration of each component can be evaluated from Eq. (3) and subsequently the termination criterion for the pressure search algorithm can be simply taken as

∑

C=

(20)

Eq. (17) resembles for PR by replacing the parameter C with B . The fugacity coefficient for pure component leads to the following expression:

ln(ϕ) = z−1−ln(z−B ) +

a z + M⎞ ln ⎛ 2RTN ⎝ z + Q ⎠ ⎜

p b+c M=⎛ −N ⎞ ⎝ 2 ⎠ RT

(8)

(22)

Ωa α (TR) R2Tc2 Pc

b=

Ωb RTc Pc

c=

Ωc RTc Pc

p b+c Q=⎛ + N⎞ ⎝ 2 ⎠ RT

(24)

(25)

2.2.2. Mixing rules The parameters, a , b and c that are obtained from the PR and PT equation of state models for pure components, are utilized to evaluate the parameters for gas mixtures. In this light, we first briefly highlight the van der Waals and Wong-Sandler mixing rule. Subsequently, we propose a modification in the WS mixing rule for better prediction of the gas hydrate equilibria.

(9)

which reduces to the PR model by simply replacing c with b indicating that the PR model has two parameters (a and b ). The parameters, a , b and c have the following forms:

a=

(23)

f iv = ϕi x i p

2.2.1. Equation of state model In this section, the PR and PT equation of state models are represented in a generalized form. The PT model is governed by

RT a − v−b v (v + b) + c (v−b)

(b + c )2 ⎤ N = ⎡bc + ⎢ ⎥ 2 ⎣ ⎦

Subsequently, the vapor phase fugacity can be evaluated as

As stated, the nonideality in vapor phase is taken into account by using an equation of state model compiled with a mixing rule. In the following, they are briefly presented.

2.2.2.1. The van der Waals mixing rule. The mixture parameters for the three-parameter PT-EoS can be expressed by the van der Waals mixing rule as

(10)

am =

∑∑

x i x j aij (26)

i

j

∑

x i bi

(11)

bm =

(27)

i

(12)

cm =

Here, R is universal gas constant, Tc is critical temperature, Pc is critical pressure and TR (=T / Tc ) is reduced temperature. The parameter α (TR) is defined as

α (TR) = [1 + F (1−TR1/2 )]2

(14)

Ωa = 3ζc2 + 3(1−2ζ c )Ωb + Ωb2 + 1−3ζ c

(15)

(28)

aij = (ai aj )1/2 (1−kij )

(29)

Here, kij is the binary interaction parameter that is estimated against binary gas hydrate formation data. x i and x j are the mole fraction of gas components in vapor phase. Accordingly, the fugacity equation for the PT equation of state using the vdW mixing rule is as follows

∑j x j aij b Q + d⎞ RT ln(ϕi ) = −RT ln(z−B ) + RT ⎛ i ⎞− ln ⎛ d ⎝ v−bm ⎠ ⎝ Q−d ⎠ a (b + c ) a + m 2i 2 i + m3 [ci (3bm + cm) + bi (3cm 2(Q −D ) 8d ⎜

Note that Ωb can be evaluated as the smallest positive real root of the following cubic equation:

Ω3b + (2−3ζ c )Ω2b + 3ζc2 Ωb−ζc3 = 0

x i ci

where aij has the following form

(13)

Ωc = 1−3ζ c

∑ i

The values of Ωa and Ωb for PR are constant and given by 0.45724 and 0.0778, respectively. While, in case of PT, these parameters along with Ωc are defined as function of the parameter ζ c as

(16)

⎜

+

+

(A−2BC −B2−B−C ) z

+ (BC + C −A) B = 0

B=

bp RT

⎟

⎜

⎟

⎟

(30)

in which, the parameters, d and Q stand for

(17)

0.5

in which, A , B and C are is computed as

ap A= 2 2 RT

⎜

⎟

Q + d ⎞ ⎛ 2Qd ⎞ ⎤ + bm)] ⎡ln ⎛ + ⎢ ⎝ Q−d ⎠ ⎝ Q 2−d 2 ⎠ ⎥ ⎦ ⎣

The governing equation of PT (Eq. (9)) can be expressed by an equation that is cubic in z as

(C −1) z 2

(21)

where M , N and Q are defined as

2.2. Vapor phase

z3

⎟

−1/2

x i∗ = 1

i

p=

cp RT

b + cm ⎞2⎤ d = ⎡bm cm + ⎛ m ⎢ 2 ⎝ ⎠⎥ ⎣ ⎦

(31)

bm + cm 2

(32)

(18)

Q=v+

(19)

The van der Waals mixing rule and its family of modified forms have shown their applicability in highly nonideal mixtures. However, in the 752

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

low density limits, the models are found to be inconsistent. Besides, these mixing rules fail to have a quadratic dependence of second virial coefficient as required by the statistical mechanics. In order to address these issues, the Wong-Sandler mixing rule [30] is proposed that is briefly described below.

(

am =

• The molecules in a liquid mixture are so closely packed that leave no free volume. • The excess free volume is considered very small and thus, the excess

bm =

Helmholtz free energy at infinite pressure is approximated as the excess Gibbs free energy at zero pressure.

Q 1−D

bm ⎛ E G + ξm ⎜ ⎝

i

D=

∑

xi

ψi =

ξi =

i

i

ai ⎞ bi ⎟ ⎠

(39)

)

(G

E

a

+ ∑i x i ξi bi i

)

(40) (41)

ci bi

(42)

2 ⎛ 3 + ψi + 1 + 6ψi + ψi ⎞ ln ⎜ ⎟ 1 + 6ψi + ψi2 3 + ψi− 1 + 6ψi + ψi2 ⎠ ⎝

1

(43)

with ξm and ψm defined as: 2 ⎛ 3 + ψm + 1 + 6ψm + ψm ⎞ ln ⎜ ⎟ 1 + 6ψm + ψm2 3 + ψm− 1 + 6ψm + ψm2 ⎠ ⎝

1

ξm =

ψm =

∑

x i ψi,

(44) (45)

i

The fugacity coefficient equation for the PR and PT equation of state using WS mixing rule are represented as follows:

lnϕi = −ln(z−B ) +

+

1 ⎛ ∂nt bm ⎞ (z−1) bm ⎝ ∂ni ⎠T ,nj ⎜

⎟

⎤ z + c2 Bm A ⎡ 1 ⎛ ∂nt2 am ⎞ 1 ∂n b − ⎛ t m ⎞ ⎥ ln ⎜ ⎟ ⎢ 2 2 B ⎢ am nt ⎝ ∂ni ⎠T ,n bm ⎝ ∂ni ⎠T ,nj ⎥ z + c1 Bm j ⎦ ⎣ ⎜

RT lnϕi = RT ln

⎟

(46)

∂ (nt bm) ⎤ 1 v + RT ⎡ ⎥ v−bm ⎢ v−bm ⎣ ∂ni ⎦ 2

∂ (n d ) ⎡ 1 ⎡ ∂ (nt a) ⎤ am ⎡ ∂ (nt dm) − ∂ (nt em) ⎤ ⎤ am ∂tn m n ∂ni ⎦ ⎥ v−dm i ⎣ ∂ni ⎦ − ⎣ ∂ni − +⎢ ln ⎥ v−em (v−dm)(dm−em) ⎢ dm−em (dm−em)2 ⎥ ⎢ ⎦ ⎣

(33)

+

am

∂ (nt em) ∂ni

(v−em)(dm−em)

−RT lnz

(47)

The partial derivatives in Eqs. (46) and (47) are presented in Appendix A. The parameter Q for the PT-EoS model appearing in the definition of partial derivatives of Eq. (47) is formulated as the same that is done earlier for the PR-EoS model (Eq. (42)), while the parameter D has the following form:

(35)

ai GE + ∗ bi RT c RT

x i ξi

where ψ and ξ are defined as

x i x j Bij

j

∑

(

1 1− ξ RT m

(34)

∑∑

(38)

a ∑i ∑j x i x j b− RT ij

Here Q and D stand for

Q=

ij

⎥ ⎦

cm = ψm bm

Comparing with the vdW based model, the WS mixing rule can be neglected in a small range of temperature and pressure but not near the low and high-density limits of fluids. At the lower limit, the fluid phase acts as a perfect gas and the molecular interactions are negligible, whereas the properties of the same resemble to a liquid at the higher limit. These facts are considered in the WS mixing rule and the model reduces to the virial equation of state for the prior case and to an activity coefficient model for the latter. In contrast, the linear mixing rule fails to capture these phenomena. Another issue is the presence of simple (e.g., hydrocarbons and inorganic gases) and complex (e.g., polar, aromatic and associating species) fluid phases. The unique characteristic of the model is that it can be equally applied in both cases without compromising the accuracy. The restriction in the other models is that they perform well only for which they are designed. For instance, the linear mixing rule (vdW) fails for complex mixtures or at the density limits. Likewise, the density dependent mixing rules proposed for complex mixtures fail when applied to the simple fluids [30]. Moreover, these models are inconsistent with the statistical mechanical requirement of quadratic dependence of the second virial coefficient on composition. The natural gas hydrate system may contain a variety of simple and complex mixtures, thus the WS mixing rule is highly recommended for the same. Here, we briefly discuss the classical form of the WS mixing rule formulated to integrate with the PR and PT equation of state models. The mixture parameters need to be evaluated for the PR model as:

bm =

) ⎤⎥ (1−k )

in which, kij denotes the binary interaction parameter that can be adjusted for the subjected binary system. On the other hand, for the PT-EoS, the mixture parameters are expressed using the WS mixing rule as

2.2.2.2. Wong-Sandler mixing rule: the classical model. In this model, the Helmholtz free energy is considered to develop a density-independent mixing rule. The following key assumptions are made in the WS mixing rule:

am D =Q RT 1−D

aj

) (

a

⎡ bi− RTi + bj− RT Bij = ⎢ 2 ⎢ ⎣

(36)

D=

The constant c∗ equals −0.623225 for the PR-EoS. The excess Gibbs free energy G E can be calculated using the non-random two-liquid model (NRTL) (Appendix B). The second virial coefficient for pure component is defined as

1 ⎛ a ∑ xi ξi bi −A∞E ⎞⎟ ξm RT ⎜ i i ⎝ ⎠

(48)

A∞E

(37)

is the Helmholtz excess free energy at infinite pressure. It is Here, approximated as the Gibbs excess free energy at zero pressure and can be calculated using a suitable activity coefficient model. In this study, we use the NRTL model (Appendix B).

Further, B for mixture (Bij ) is related to the parameters of pure substance as

2.2.2.3. Modified Wong-Sandler mixing rule: the proposed model. The motivation behind the modification in WS mixing rule is to present the

B = b−

a RT

753

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

second virial coefficient (Bij ) more accurately as it decides the contribution of pure state parameters in the calculation of mixture parameters. The classical WS approach approximates Bij in gas mixtures considering the arithmetic mean of individual components (Eq. (38)). In this work, the formulation of Bij is modified by separately defining the correction factors for energy parameter (a ) and co-volume parameter (b ) as they constitute to the formulation of second virial coefficient in dissimilar fashion.The parameter a in cubic equation of state stands for the energy parameter. Jian and Zongcheng [43] have demonstrated that the parameter a is derived from the UNIFAC model by expressing the Gibbs free energy in two parts: combinational part and residual part.

a = bRT

∑

xi

i

E aii 1 Gres + bi RT A RT

Λ21 Λ12 ⎞ lnγ2 = −ln(X2 + Λ21X1) + X1 ⎛ − ⎝ X2 + Λ21X1 X1 + Λ12X2 ⎠ ⎜

2

−

(49)

bi + bj 2

(1−lij )−

RT

Λ21 =

v1 λ −λ exp ⎛− 21 22 ⎞ v2 RT ⎠ ⎝

(57)

(58)

fi = ai fio

(59)

fio

is the standard phase fugacity of i th component and can be where expressed as

vi (p−pis ) ⎞ fio = pis ϕis exp ⎛⎜ ⎟ ⎝ RT ⎠

(50)

ai aj

(56)

In turn, the fugacity in the liquid phase will be

pis

(60)

ϕis

is the saturation pressure, the standard phase fugacity in which, and vi the volume of the i th component in liquid phase. The saturation pressure ( pis ) can be calculated by the Antoine equation and the standard phase fugacity is estimated by the fugacity equation for pure component given in Eq. (21).

Jian and Zongcheng [43] have introduced the correction factor (lij ) for co-volume parameter (b ), in which, the linear mixing rule is adopted. Likewise, the correction factor (kij ) is defined for energy parameter (a ) which is geometric mean of pure substance parameter [36]. As a combination of these two approaches, the second virial coefficient (Bij ) represented in Eq. (38) can be rewritten in a modified form as

Bij =

v2 λ −λ exp ⎛− 12 11 ⎞ v1 RT ⎠ ⎝

ai = Xi γi

ai aj RT

Λ12 =

Here, v1 and v2 are the molar volume of water and THF, respectively under the specific conditions. The molar energy parameters λ12−λ11 and λ21−λ22 can be determined from the binary experimental data. The activity of i th component can be evaluated as

is the residual excess Gibbs free energy. The combinational Here, part contains the parameter a , which can be expressed as the geometric combination of the energy parameters of pure components [36].The other parameter that constitutes the second virial coefficient is the covolume parameter (b ). This b for gas mixtures is determined with the arithmetic combination of pure component parameters [43]. So, with the new formulation, Bij can be represented as

bi + bj

(55)

where γ1 and γ2 are the fugacity coefficients for water and THF, respectively. X1 and X2 stand for mole fraction of respective components in liquid phase. The parameters Λ12 and Λ21 are expressed as

E Gres

Bij =

⎟

3. Simulation algorithm

(1−kij )

(51)

The simulation proceeds with the systematic calculation of three phase fugacity existed in a gas hydrate system. In this particular study, the procedure follows the sII hydrate formation as the systems contain THF as a thermodynamic promoter. The experimental data for hydrate phase equilibrium is available as pressure-temperature relation. This equilibrium is evaluated here for a particular temperature by adjusting the pressure. Accordingly, the computer assisted pressure search algorithm is developed below.

Panagiotopoulos and Reid [44] have suggested the modification in the correction factor (kij ) to be linearly dependent on composition. In addition, Espanani et al. [36] have considered the correction factor (kij ) as a quadratic function of temperature. However, the coefficient of the second power of temperature is very small and thus the expression for kij can be reduced to a linear function of temperature. Herein, we are proposing the correction factors for energy and co-volume parameter to be linearly dependent on temperature as well as vapor phase composition. In this view, the correction factors lij and kij can be formulated as

lij = lij,0 + lij,1 x i + lij,2 T

(52)

kij = kij,0 + kij,1 x i + kij,2 T

(53)

(i) Assign the hydrate structure (sI or sII) and assume the pressure ( p ) at a particular temperature (T ). (ii) Calculate the vapor phase fugacity coefficient (ϕiv ) for each component using Eqs. (30), (46) and (47) for PT-vdW, PR-WS and PT-WS model, respectively. (iii) Using the values of ϕiv calculated in step (ii), determine f iv using Eq. (25) for respective components. (iv) Compute the Langmuir constant (cj ) for each component from Eq. (5) and hence evaluate the corresponding fractional occupancy (θj ) from Eq. (4). (v) Calculate the temperature dependent term ( fTo ) of standard hydrate phase fugacity from the Antoine-type equation (Eq. (7)). (vi) Calculate the activity of water (a w ) using the Wilson model (Eq. (58)). (vii) Using the calculated values of fTo and a w , estimate the standard hydrate phase fugacity ( fio ) using Eq. (6) for each component. (viii) Apply the equilibrium condition as given in Eq. (1) and compute the hydrate phase concentration ( x i∗) for each component using Eq. (3). (ix) For (∑i x i∗−1) ⩽ ε , where ε is the tolerance limit, save calculated p and terminate. (x) If (∑i x i∗−1) ⩾ ε , adjust p by secant method and go to Step (ii).

Now, for the proposed modified WS mixing rule, the new formulation of Bij defined in Eq. (51) is utilized in place of previously defined Eq. (38) for the classical model. While, the fugacity coefficient equations defined in Eqs. (46) and (47) for the classical WS mixing rule is used as the same. 2.3. Liquid phase In this work, the liquid phase fugacity needs to be evaluated for water and THF. As mentioned in one of the assumptions in the ChenGuo model, the solubility of guest gas is negligible and hence, the activity of water is only affected by the addition of THF. The liquid phase fugacity coefficient can be sufficiently described by the Wilson activity coefficient model as

Λ12 Λ21 ⎞ lnγ1 = −ln(X1 + Λ12X2 ) + X2 ⎛ − ⎝ X1 + Λ12X2 X2 + Λ21X1 ⎠ ⎜

⎟

(54) 754

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

Table 1 Thermodynamic properties used in this study.

Table 3 Antoine constants for temperature dependent term of hydrate phase fugacity ( f io,T ).

Species

ωa

ζc a

Fa

Tc (K)b

Pc (MPa)b

Species

A′ × 10−22 Pa −1

B′

C′

CH4 N2 C2H6 C2H4 THF

0.012 0.039 0.100 0.085 –

0.3240 0.3290 0.3170 0.3180 0.4126

0.455336 0.516798 0.561567 0.554369 0.045241

190.56 126.20 305.32 282.34 540.20

4.599 3.398 4.872 5.041 5.190

CH4 N2 C2H6 C2H4 THF

5.260200a 6.816500a 0.039900a 0.038662 616.2550

−12955.00a −12770.00a −11491.00a −11639.16 −24787.50

4.0800a −1.1000a 30.400a 32.249 −109.30

Ref:

Ref: a b

a

Patel and Teja. [26]. Poling et al.[42].

high-pressure phase equilibrium analysis is performed in a transparent sapphire cell equipped with a voltage regulating magnetic stirrer for maintaining the uniformity in concentration and temperature. For all systems, authors have pre-treated the THF solution and gas mixture in the reactor by repeated pressurization (hydrate formation) and depressurization (hydrate dissociation) to avoid the large induction time in hydrate formation due to memory effect. The solution is then stepwise (0.05 MPa) pressurized until the first trace of hydrate appears. If it can last for at least 4 h, the operating temperature and pressure are recorded for an equilibrium point. The experiments are repeated twice and the uncertainties in the measurements of temperature, pressure, gas composition and aqueous phase concentration of THF are reported within ± 0.1 K, ± 0.025 MPa, ± 0.01%, ± 0.005%, respectively (for all systems). These values ensure the reliability of experimental data sets used in the validation of the models. To illustrate the proposed WS mixing rule hybridized with the PR and PT equation of state models, the experimental data sets for the following THF-hydrate systems are used: CH4/N2/THF (System I), CH4/ C2H6/THF (System II) and CH4/C2H4/THF (System III). Further, we develop the classical WS based equation of state models for the said systems. Finally, a systematic comparison is presented below to show the superiority of the proposed phase equilibrium models over the existing PT-vdW model.

The model is simulated following the above algorithm by coding in MATLAB® environment. The adjustable parameters are determined by minimizing the average absolute relative deviation (AARD) as given below for predicted hydrate formation pressure against the experimental one.

%AARD =

1 n

n

∑ i=1

exp mod ⎛ |pi −pi | ⎞ exp ⎜ ⎟ × 100 pi ⎝ ⎠

(61)

Here n stands for the number of experimental data points and the vertical bars return the absolute value of the quantity in between. The simplex search method (fminsearch: an inbuilt tool in MATLAB®) is used for minimizing the objective function (%AARD). 4. Determining model parameters The three-phase fugacity model comprises of the pure substance and interaction parameters. The former includes the inherent thermodynamic properties (Table 1) that are non-adjustable and thus, taken from the literature [23,29,42]. The reported values of Antoine parameters for determining Langmuir constant and temperature dependent term of fugacity (Tables 2 and 3, respectively) are used as the same for all components except C2H4 and THF as they are reported differently elsewhere [28,29]. A single set of these values are determined by regressing against the experiment data and used as the same for phase equilibrium analysis.On the other hand, the adjustable parameters are figured out separately and reported in this work with their optimized values. These include the interaction parameters that occur in the expression of the vdW and WS mixing rule (kij and lij ) (Table 4), Wilson model (λij ) and the temperature dependent term of fugacity ( Aij ) (Table 5). These adjusted values are reported differently for each system. The NRTL parameters for GE calculation in WS mixing rule is determined as a single set for the WS and the mWS mixing rule and reported in Table 6.

5.1. Prediction using the existing phase equilibrium model Here, we consider a phase equilibrium model that consists of the Chen-Guo model for hydrate phase, Wilson model for liquid phase and an equation of state model (PR/PT) combined with a vdW mixing rule for vapor phase. This model was used earlier by Sun et al. [29] and we compare this model considering two different sets of parameters; one set is already reported in literature [28,29] and second set is proposed here for the three said systems. For these systems, the available experimental data sets are used to reproduce the concerned model parameters. The aim is to show an improved performance with these updated parameters compared to the available parameter values. Fig. 1 depicts the phase equilibrium condition for hydrate formation with a binary gas, namely CH4/N2 (System I), varying the CH4 concentration from 4.90 to 46.28 mol%. Here, a fixed concentration of THF (i.e., 6 mol%) is used and it leads to a dramatic reduction of formation pressure from 1-35 MPa to almost 0.1–3.5 MPa. To predict the phase equilibrium of the same system, Sun et al. [29] have demonstrated the Chen-Guo model with PT equation of state that uses the van der Waals mixing rule for computing the gas mixture properties. They have reported the prediction with their model parameters having an AARD of 6.86%. Here, as stated, the same PT-vdW model is reproduced with a new parameter set (Tables 4 and 5) for System I using which the model resembles more accurately with an AARD reduced to 5.69% from 6.86% [29]. The other two systems considered in this study include Systems II and III having guest gases as saturated (C2H6) and unsaturated (C2H4) hydrocarbon, respectively in mixture with CH4 and THF. For System II, the phase equilibrium data are adopted from literature [28] with a varying CH4 concentration from 16.8 to 93.4 mol% (Fig. 2). In the same

5. Results and discussion Prior to discussing the simulation results, here we highlight the experimental conditions and the procedure of phase equilibrium measurements. All the experiments are performed within the temperature range of 277–287 K following an identical procedure [27,28,41]. The Table 2 Antoine parameters to calculate Langmuir constants. Species

X × 105

Y

Z

CH4 N2 C2H6 C2H4 THF

2.3048a 4.3151a 4.2007a 4.2141 0a

2752.290a 2472.370a −762.9513a −764.5044 0a

23.0100a 0.64000a 29.6727a 29.6491 0a

Ref: a

Sun et al. [29].

Sun et al. [29].

755

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

Table 4 Coefficients to find gas phase interaction coefficient (kij and lij ). System

EoS-MR

kij,0

kij,1

kij,2

lij,0

lij,1

lij,2

I

PT-vdW PR-WS PT-WS PR-mWS PT-mWS

0.1291 0.0105 −0.2419 0.0010 0.2769

−0.0011 −0.0023

−0.0205 33.7688

0.1095 0.4175

−0.0009 0.0025

−0.02070 28.5837

PT-vdW PR-WS PT-WS PR-mWS PT-mWS

0.5467 1.0810 0.9155 0.0609 0.9267

2.3045 −0.4098

−0.0011 0.0091

0.1023 0.0020

10.8491 −0.5178

−0.01980 0.01030

PT-vdW PR-WS PT-WS PR-mWS PT-mWS

0.5003 0.7620 0.8101 0.7131 0.7790

−0.4139 −0.4900

0.0001 −0.0001

0.0002 0.1264

−0.0525 0.0140

0.0010 0.0003

II

III

Note: Subscript 2 stands for N2, C2H6, and C2H4 for Systems I, II, and III respectively; 1 and 3 stand for CH4 and THF respectively for all systems.

5.2. Prediction using the proposed model

Table 5 Hydrate phase interaction coefficient ( Aij ) and Wilson parameters. System

EoS-MR

A12

A23

A31

λ12−λ11

λ21−λ22

I

PT-vdW PR-WS PT-WS PR-mWS PT-mWS

131.4395 233.6686 371.4133 229.1903 333.0456

75.7608 134.1485 134.9798 139.7231 148.4424

623.6073 741.4573 740.1255 791.8032 541.7978

1133.6

1693.8

II

PT-vdW PR-WS PT-WS PR-mWS PT-mWS

146.8703 179.0394 417.5436 177.0732 414.7526

88.0186 67.4025 70.4732 65.8366 70.8059

253.0176 344.1920 359.6322 337.2609 362.7103

1092.5

1661.8

III

PT-vdW PR-WS PT-WS PR-mWS PT-mWS

124.9576 261.2918 267.7592 259.2750 221.7457

65.8115 70.6844 72.9160 71.1481 71.4850

246.9535 338.9092 341.9934 358.9689 361.4776

1088.6

1679.7

Now, we compare this PT-vdW (updated) based model with the four proposed thermodynamic models that comprise either of the PR and PT with a mixing rule among classical and modified WS, yielding PR-WS, PR-mWS, PT-WS and PT-mWS. Finding the optimal parameter sets for the said models, their performance is tested with three different systems. In Fig. 1, we first present a systematic comparison for System I (CH4/N2/THF) using the experimental data. It is investigated that the vdW mixing rule based PT-vdW model shows a worse performance, (5.69% AARD). As indicated, the WS based models are tested on the same system, among which the PT-WS model (AARD = 5.64%) performs a little better than the PR-WS (AARD = 5.68%). Since the improvement is not so significant using the WS over the vdW, the PT-WS and PR-WS models are not shown in Fig. 1; but the corresponding AARD values are listed in Table 7. Now, we would like to compile the proposed mWS mixing rule with an equation of state model (PR/PT) for vapor phase. As considered before, the Chen-Guo model is employed for hydrate phase and the Wilson model for liquid phase. For System I, a significant reduction of AARD is harnessed to have the values of 5.59% and 5.36% for the PRmWS and PT-mWS, respectively. It is evident that the PT performs better than the PR. More importantly, the proposed mWS secures a better response than the other two classical mixing rules, namely vdW and WS. Systems II and III are having comparatively more nonideal vapor phase owing to the interactive species like saturated (C2H6) and unsaturated (C2H4) hydrocarbons, respectively. The phase equilibrium data for these two systems is further used to examine the proposed models. In this light, the phase predictions using the WS mixing rules based models are compared with those obtained using the classical vdW based models (Fig. 2). Although, the classical WS based models are not compared in Fig. 2, their corresponding AARD values are listed in

Note: Subscript 2 stands for N2, C2H6, and C2H4 for Systems I, II, and III respectively; 1 and 3 stand for CH4 and THF respectively for all systems.

way, the equilibrium data for System III are available [41] for a vapor phase composition of CH4 from 11.63 to 93.72 mol% in the CH4/C2H4 gas mixture (Fig. 3). The pressure predictions for System II are available in literature [28] using the PT-vdW based model yielding an average error of 4.4% (i.e., 1.34% AARD). The parameters for this system are updated and the AARD is found to be reduced to 1.12% form 1.34% (Table 7). In the similar fashion, the phase equilibrium of System III is predicted using a new parameter set proposed in this work giving 1.72% AARD. Note that for this System III, there are no results available in the literature.

Table 6 NRTL parameters for WS and mWS mixing rule. System

EoS

a12

a21

b12

b21

c12

c21 × 10 4

d12

d21

e12

f12 × 105

I

PR PT

−164 −10.8

−284 −15.0

−130 −8.66

63.0 0.49

165 5.13

37.0 2.62

118 5.79

350 9.31

0.50 0.52

16.4 −3.42

II

PR PT

−1.23 −1.28

−2.19 2.30

−1.66 −1.63

0.95 1.01

1.69 1.69

46.6 0.48

1.39 1.45

4.71 4.37

0.58 0.61

210 200

III

PR PT

−1.29 −1.33

−2.31 −2.38

−1.64 −1.69

1.01 1.03

1.71 1.69

0.48 0.53

1.46 1.45

4.39 4.18

0.62 0.62

200 200

Note: Subscript 2 stands for N2, C2H6, and C2H4 for Systems I, II, and III respectively; 1 and 3 stand for CH4 and THF respectively for all systems.

756

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

4.0

Fig. 1. Predictions of formation conditions for the CH4/N2/ THF hydrate (System I).

Data (4.90% CH4)

3.5

Data (15.99% CH4) Data (29.21% CH4)

Pressure (MPa)

3.0

Data (46.28% CH4) PT-vdW (Sun et al., 2015) PT-vdW (updated) PR-mWS (proposed) PT-mWS (proposed)

2.5

2.0

1.5

1.0

0.5

0.0 278

280

282

284

286

288

Temperature (K) WS based models (Fig. 3), we achieve the lowest AARD of 0.80% (PRmWS) and 0.77% (PT-mWS). We see the superiority of the WS over the vdW mixing rule that endorses the ability of the prior to capture the asymmetric property of the double-bond containing system (System III). Analyzing Figs. 1, 2 and 3, it can be observed that the simulation using the proposed Wong-Sandler (WS and mWS) mixing rule based model have uncertain positive and negative deviations from other models (vdW) in different cases. This can be explained with the fact that the prior is a free energy based model, which depends on the gas phase composition and temperature of the system at which the different equilibrium data points are recorded. Despite, the latter uses a linear mixing rule having empirical formulation. Consequently, for the WS

Table 7 as 0.96% and 0.90% for PR-WS and PT-WS, respectively. An improved performance of the modified WS based models is reflected through the AARD figures appended in Table 7 (0.93% for PR-mWS and 0.83% for PT-mWS). Likewise, a comparative model prediction with respect to System III is presented in Fig. 3. It is evident that the vdW mixing rule based model leads to the least accuracy having a percent AARD of 1.72%. The improved predictions are obtained using the WS based models having percent AARD of 0.88% (PR-WS) and 0.85% (PT-WS). Although, the improvement is quite significant, only the predictions using mWS based models are shown in Fig. 3 and those of the classical (WS) ones are only represented in Table 7 in terms of their AARD values. For the modified 3.2

2.8

Pressure (MPa)

2.4

2.0

1.6

Fig. 2. Predictions of formation conditions for the CH4/C2H6/ THF hydrate (System II).

Data (93.4% CH4) Data (71.6% CH4) Data (55.4% CH4) Data (38.2% CH4) Data (16.8% CH4) PT-vdW (Sun et al., 2010) PT-vdW (updated) PR-mWS (proposed) PT-mWS (proposed)

1.2

0.8

0.4

0.0 276

278

280

282

284

286

288

Temperature (K) 757

290

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

Fig. 3. Predictions of formation conditions for the CH4/C2H4/ THF hydrate (System III).

3.5

Data (11.63% CH4) Data (23.85% CH4)

3.0

Data (42.88% CH4) Data (70.65% CH4)

Pressure (MPa)

2.5

Data (86.24% CH4) Data (93.72% CH4) PT-vdW (updated) PR-mWS (proposed) PT-mWS (proposed)

2.0

1.5

1.0

0.5

0.0 278

280

282

284

286

288

Temperature (K) 6. Conclusions

Table 7 AARD percentage analysis for all the three systems. System

I II III

Reported (PT-vdW) (% AARD)

6.86a 1.34b NA

This work formulates the thermodynamic phase equilibrium model consisting of the Chen-Guo model (hydrate phase) and the Wilson model (liquid phase) with an equation of state (vapor phase) in conjunction with the Wong-Sandler mixing rule. Here, we use two equation of state models, namely PR and PT, and introduce a modification in the WS mixing rule to develop the hybrid model. These proposed PR-mWS and PT-mWS based models are compared with the classical WS based PR-WS and PT-WS, and the existing PT-vdW and its updated form with reference to the experimental data available for three systems, namely CH4/N2 (System I), CH4/C2H6 (System II), and CH4/C2H4 (System III), in presence of THF as a thermodynamic promoter. The comparative study of different equation of state models and mixing rules envisage that the proposed mWS based model is better against the classical WS and vdW based models to predict the phase equilibrium of the NGH systems. Also, the PT equation of state shows a better prediction over the PR and hence, the PT-mWS based model is best suited for the phase equilibrium predictions of natural gas hydrates. The percent AARD values by using this model are obtained as 5.36% and 0.83% against the reported values (6.86% [29] and 1.36% [28]) for Systems I and II, respectively. The model predictions for System III are not available in open literature and for this a fairly good prediction was obtained with a 0.77% AARD using the PT-mWS based model.

This work (% AARD) PT-vdW

PR-WS

PT-WS

PR-mWS

PT-mWS

5.69 1.12 1.72

5.68 0.96 0.88

5.64 0.90 0.85

5.59 0.93 0.80

5.36 0.83 0.77

Ref: NA = Not available so far. a Sun et al. [29]. b Sun et al. [28].

mixing rule based model, the excess Gibbs free energy contribution of components in the mixture may differ with the cases (i.e., data points) depending on the molecular interactions at particular temperature and pressure, despite the same effect in absent in vdW based model.It is now obvious that for all three example systems, the proposed mWS leads to provide the best performance consistently. However, the improvement using the PT over the PR is not so significant as the value of reduced temperature (TR ) is nearly equal to 1.5 in all these cases of gas hydrate formation conditions. It agrees with the observation reported earlier by Patel and Teja [26] that a high value of reduced temperature (TR > 12) significantly differentiate the PR and PT-EoS in terms of their accuracy of predictions. Appendix A Partial derivatives for Wong-Sandler mixing rule 2 ⎡ ∂ (nt bm) ⎤ = 1 ⎜⎛ 1 ∂ (nt Q) ⎟⎞− Q ⎛1− ∂ (nt D) ⎞ ⎥ ⎢ ∂ni ⎦ ∂ni ⎠ 1−D ⎝ nt ∂ni ⎠ (1−D)2 ⎝ ⎣

(A1)

∂ (nt2 bm) ∂ (nt D) ⎞ 1 ∂ (nt2 am) + bm = RT ⎛⎜D ⎟ ∂ni ∂ni ⎠ nt ∂ni ⎝

(A2)

1 ∂ (nt2 Q) a ⎞ = 2 ∑ x j ⎛b− nt ∂ni ⎝ RT ⎠ij j

(A3)

⎜

⎟

758

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

lnγi∞ ∂ (nt D) ai = + ; for PR−WS ∂ni bi RT c∗

(

(A4) a

a ξi bi −RT lnγi∞ nt RT ∑i x i ξi bii ∂ (nt D) i = − ∂ni ξm RT (ξm RT )2

c

)

∂ξm ∂ni

; for PT −WS

(A5)

c

3 bi −3ψm + bi ψm−ψm2 ⎛ 3 + ψi + 1 + 6ψi + ψi2 ⎞ ∂ξ i ln ⎜ nt m = − i ⎟+ ∂ni (1 + 6ψm + ψm2 )2 3 + ψi− 1 + 6ψi + ψi2 ⎝ ⎠

1 1 + 6ψm + ψm2

⎡ ⎢ (ψi−ψm) + ×⎢ ⎢ 3 + ψm + ⎢ ⎣

c c 2 3 i − 3ψm + i ψm − ψm bi bi 2 1 + 6ψm + ψm

1 + 6ψm + ψm2

−

c c 2 3 i − 3ψm + i ψm − ψm b bi (ψi−ψm)− i 2 1 + 6ψm + ψm

3 + ψm− 1 + 6ψm + ψm2

⎤ ⎥ ⎥ ⎥ ⎥ ⎦ (A6)

∂ (nt cm) ∂ (nt bm) = ∂ni ∂ni

∑ i

xi

ci c c + bm i + bm ∑ x i i bi bi bi i

bm ∂ (nt cm) ⎞ ∂ (nt dm) 1 ⎡ ∂ (nt bm) + + = ⎢−⎛ ∂ni ⎠ ∂ni 2 ⎢ ⎝ ∂ni ⎣ ⎜

+ 3bm

∂ (nt cm) ∂ni

+ 3cm

∂ (nt bm) ∂ni

+ cm

∂ (nt cm) ∂ni

⎟

∂ (nt cm) ⎞ bm ∂ (nt em) 1 ⎡ ∂ (nt bm) + − = ⎢−⎛ ∂ni ⎠ ∂ni 2 ⎢ ⎝ ∂ni ⎣ ⎜

∂ (nt bm) ∂ni

(A7)

∂ (nt bm) ∂ni

bm2 + 3bm

+ 6bm cm +

∂ (nt cm) ∂ni

+ 3cm

cm2

∂ (nt bm) ∂ni

+ cm

∂ (nt cm) ∂ni

⎟

bm2

+ 6bm cm +

cm2

⎤ ⎥ ⎥ ⎦

⎤ ⎥ ⎥ ⎦

(A8)

(A9)

Appendix B NRTL model for G E calculation in Wong-Sandler mixing rule and its modified form

GE τ G τ12 G12 ⎞ = x1 x2 ⎛ 21 21 + RT x2 + x1 G12 ⎠ ⎝ x1 + x2 G21 ⎜

⎟

(A10)

where

G12 = exp(−α12 τ12)

(A11)

G21 = exp(−α21 τ21)

(A12)

The τ12 and τ12 are expressed in terms of reduced temperature

τ12 = a12 + b12 Tr1 + c12 Tr21 + d12/ Tr1

(A13)

τ21 = a21 + b21 Tr 2 + c21 Tr22 + d21/ Tr 2

(A14)

and α is the non-randomness parameter that is also dependent on temperature as

α12 = e12 + f12 T

(A15)

α21 = a12

(A16)

[11] Martinez MC, Dalmazzone D, Fürst W, Delahaye A, Fournaison L. Thermodynamic properties of THF+CO2 hydrates in relation with refrigeration applications. AIChE J 2008;54:1088–95. [12] Seo YT, Kang SP, Lee H. Experimental determination and thermodynamic modeling of methane and nitrogen hydrates in the presence of THF, propylene oxide, 1,4dioxane and acetone. Fluid Phase Equilib 2001;189:99–110. [13] Pahlavanzadeh H, Kamran-Pirzaman A, Mohammadi AH. Thermodynamic modeling of pressure–temperature phase diagrams of binary clathrate hydrates of methane, carbon dioxide or nitrogen+ tetrahydrofuran, 1,4-dioxane or acetone. Fluid Phase Equilib 2012;320:32–7. [14] Sujith KS, Ramachandran CN. Natural gas evolution in a gas hydrate melt: effect of thermodynamic hydrate inhibitors. J Phy Chem B 2017;121:153–63. [15] Yang SHB, Babu P, Chua SFS, Linga P. Carbon dioxide hydrate kinetics in porous media with and without salts. Appl Energy 2016;162:1131–40. [16] Li XS, Xu CG, Zhang Y, Ruan XK, Li G, Wang Y. Investigation into gas production from natural gas hydrate: a review. Appl Energy 2016;172:286–322. [17] Li SL, Ma QL, Sun CY, Chen LT, Liu B, Feng XJ, Wang XQ, Chen GJ. A fractal approach on modeling gas hydrate phase equilibria in porous media. Fluid Phase Equilib 2013;356:277–83. [18] Lee S, Seo Y. Experimental measurement and thermodynamic modeling of hydratephase equilibria for the ternary C2H6+ NaCl+ water and C3H8+ NaCl+ water mixtures in silica gel pores. Energy Fuels 2010;24:6074–80. [19] Tsimpanogiannis IN, Diamantonis NI, Economou IG, Papadimitriou NI, Stubos AK. Influence of combining rules on the cavity occupancy of clathrate hydrates using van der Waals–Platteeuw-theory-based modelling. Chem Eng Res Des 2014;92:2992–3007.

References [1] Sloan ED, Koh CA. Clathrate hydrates of natural gas. 3rd ed. Boca Raton, FL: CRC Press; 2008. [2] Vedachalam N, Ramesh S, Srinivasalu S, Rajendran G, Ramadass GA, Atmanand MA. Assessment of methane gas production from Indian gas hydrate petroleum systems. Appl Energy 2016;168:649–60. [3] Chong ZR, Yang SHB, Babu P, Linga P, Li XS. Review of natural gas hydrates as an energy resource: prospects and challenges. Appl Energy 2016;162:1633–52. [4] Sloan ED. Fundamental principles and applications of natural gas hydrates. Nature 2003;426:353–63. [5] Koh DY, Kang H, Lee JW, Park Y, Kim SJ, Lee J, Lee JY, Lee H. Energy-efficient natural gas hydrate production using gas exchange. Appl Energy 2016;162:114–30. [6] Wang XH, Sun YF, Wang YF, Li N, Sun CY, Chen GJ, Liu B, Yang LY. Gas production from hydrates by CH4-CO2/H2 replacement. Appl Energy 2017;188:305–14. [7] Zhao J, Zhao Y, Liang W. Hydrate-based gas separation for methane recovery from coal mine gas using tetrahydrofuran. Energy Technol 2016;4:864–9. [8] Zhong DL, Ding K, Lu YY, Yan J, Zhao WL. Methane recovery from coal mine gas using hydrate formation in water-in-oil emulsions. Appl Energy 2016;162:1619–26. [9] Yang M, Song Y, Jiang L, Liu W, Dou B, Jing W. Effects of operating mode and pressure on hydrate-based desalination and CO2 capture in porous media. Appl Energy 2014;135:504–11. [10] Lee Y, Kim Y, Lee J, Lee H, Seo Y. CH4 recovery and CO2 sequestration using flue gas in natural gas hydrates as revealed by a micro-differential scanning calorimeter. Appl Energy 2015;150:120–7.

759

Applied Energy 205 (2017) 749–760

N. Thakre, A.K. Jana

containing systems. Fluid Phase Equilib 2003;205:291–302. [33] Chin HY, Hsieh MK, Chen YP, Chen PC, Lin ST, Chen LJ. Prediction of phase equilibrium for gas hydrate in the presence of organic inhibitors and electrolytes by using an explicit pressure-dependent Langmuir adsorption constant in the van der Waals-Platteeuw model. J Chem Thermodyn 2013;66:34–43. [34] Haghtalab A, Moghaddam AK. Prediction of minimum miscibility pressure using the UNIFAC group contribution activity coefficient model and the LCVM mixing rule. Ind Eng Chem Res 2016;55:2840–51. [35] Yang T, Chen GJ, Yan W, Guo TM. Extension of the Wong-Sandler mixing rule to the three-parameter Patel-Teja equation of state: application up to the near-critical region. Chem Eng J 1997;67:27–36. [36] Espanani R, Miller A, Jacoby W. Prediction of vapor-liquid equilibria for mixtures of low boiling point compounds using Wong-Sandler mixing rule and EOS/GE model. Chem Eng Sci 2016;152:343–50. [37] Djavidnia A, Izadpanah AA, Sefti MV, Varaminian F. The equilibrium data and thermodynamic modeling of hydrate formation for CO2 and mixtures of CO2 and CH4 in the presence of methanol. Pet Sci Technol 2013;31:2013–21. [38] Shahnazar S, Hasan N. Gas hydrate formation condition: Review on experimental and modeling approaches. Fluid Phase Equilib 2014;379:72–85. [39] Klauda JB, Sandler SI. Phase behavior of clathrate hydrates: a model for single and multiple gas component hydrates. Chem Eng Sci 2003;58:27–41. [40] Sun CY, Chen GJ, Lin W, Guo TM. Hydrate formation conditions of sour natural gases. J Chem Eng Data 2003;48:600–2. [41] Zhang LW, Huang Q, Sun CY, Ma QL, Chen GJ. Hydrate formation conditions of methane+ ethylene+ tetrahydrofuran+ water systems. J Chem Eng Data 2006;51:419–22. [42] Poling BE, Prausnitz JM, John Paul OC, Reid RC. The properties of gases and liquids. 5th ed. New York, US: McGraw-Hill; 2001. [43] Jian C, Zongcheng L. Study on the GE-type mixing rule for cubic equations of state. Tsinghua Sci Technol 1996;1:410–5. [44] Panagiotopoulos AZ, Reid RC. Multiphase high pressure equilibria in ternary aqueous systems. Fluid Phase Equilib 1986;29:525–34.

[20] Kwon M, Youn Y, Seo Y, Lee JW, Lee J, Lee JY, Kim SJ, Lee H. Improved experimental determinations of phase equilibria and structural transitions of mixed gas hydrates under isothermal conditions. Energy Fuels 2013;27:5144–52. [21] Kwon M, Lee JW, Lee H. Temperature-dependent structural transitions in methane–ethane mixed gas hydrates. J Phys Chem C 2014;118:28906–13. [22] Uchida T, Takeya S, Kamata Y, Ikeda IY, Nagao J, Ebinuma T, Narita H, Zatsepina O, Buffett BA. Spectroscopic observations and thermodynamic calculations on clathrate hydrates of mixed gas containing methane and ethane: Determination of structure, composition and cage occupancy. J Phys Chem B 2002;106:12426–31. [23] Chen GJ, Guo TM. A new approach to gas hydrate modelling. Chem Eng J 1998;71:145–51. [24] Klauda JB, Sandler SI. A fugacity model for gas hydrate phase equilibria. Ind Eng Chem Res 2000;39:3377–86. [25] Sheikhi-Kouhsar M, Bagheri H, Raeissi S. Modeling of ionic liquid+polar solvent mixture molar volumes using a generalized volume translation on the PengRobinson equation of state. Fluid Phase Equilib 2015;395:51–7. [26] Patel NC, Teja AS. A new cubic equation of state for fluids and fluid mixtures. Chem Eng Sci 1982;37:463–73. [27] Sun Q, Guo X, Liu A, Dong J, Liu B, Zhang J, Chen G. Experiment on the separation of air-mixed coal bed methane in THF solution by hydrate formation. Energy Fuels 2012;26:4507–13. [28] Sun CY, Chen GJ, Zhang LW. Hydrate phase equilibrium and structure for (methane + ethane+ tetrahydrofuran+ water) system. J Chem Thermodyn 2010;42:1173–9. [29] Sun Q, Guo X, Chapman WG, Liu A, Yang L, Zhang J. Vapor–hydrate two-phase and vapor–liquid–hydrate three-phase equilibrium calculation of THF/CH4/N2 hydrates. Fluid Phase Equilib 2015;401:70–6. [30] Wong DSH, Sandler SI. A theoretically correct mixing rule for cubic equations of state. AIChE J 1992;38:671–80. [31] Boesen RR, Herslund PJ, Sørensen H. Loss of monoethylene glycol to CO2- and H2Srich fluids: modeled using Soave–Redlich–Kwong with the Huron and Vidal mixing rule and cubic-plus-association equations of state. Energy Fuels 2017. http://dx.doi. org/10.1021/acs.energyfuels.6b02365. [32] Ma QL, Chen GJ, Guo TM. Modelling the gas hydrate formation of inhibitor

760