Sensitivity Analysis of Computations of the Vapor-Liquid Equilibria of Methane + Methanol or Glycols at Gas Hydrate Formation Conditions

The Soave-Redlich-Kwong (SRK-EOS) and Peng-Robinson (PR-EOS) equations of state are used often to describe the behavior of pure substances and mixtures despite difficulties in handling substances, like water, with high polarity and hydrogen bonding. They were employed in studying the binary vapor-liquid equilibria (VLE) of methane + methanol, monoethylene glycol (MEG), and triethylene glycol (TEG). These liquids are used to inhibit the formation of gas hydrates. The investigation focused on the conditions at which methane-water clathrates can form 283.89 K to 323.56 K and 5.01 MPa to 18.48 MPa. The pressure of methane in methanol is overestimated by a factor of two by either the SRK-EOS or the PR-EOS. In the methane + MEG system, the predicted pressures for both equations of state are generally less than experimental pressure except for the highest concentration of methane in MEG calculated by the SRK-EOS. In the methane + TEG system, the predictions of both models are close and trend similarly. Because of the comparative lack of extensive experimental methane + TEG data, the similarity of the methane + TEG computed results can be used as a basis for further study of this system experimentally.


Introduction
The pipelines used in the offshore production of oil and gas can experience serious safety and flow assurance problems because of plugging by solid deposits of gas hydrates, waxes, asphaltenes, and scale [1].The most common and most serious problems are caused by gas hydrates since they are prone to form quickly C. E. Ozigagu, A. J. Duben DOI: 10.4236/mnsms.2019.91001 2 Modeling and Numerical Simulation of Material Science in the deeper and colder waters below which the hydrocarbons are removed [2].
Gas hydrates are clathrate inclusion compounds in which water molecules form hydrogen bonded cages in a lattice structure stabilized by encapsulating a small guest molecule such as methane or ethane [3].A gas hydrate is an ice-like solid that can exist at temperatures up to 25˚C in systems with moderate pressure at 7 MPa.Dissolved water may condense and alter the physical state to an unwanted two-phase flow [4].These conditions are commonly encountered in oil and gas offshore process facilities [5].
Because of the lack of experimental data, simulation software packages are used to perform complex phase equilibria calculations to model systems in the refining and chemical industries.Cubic equations of state are widely used in these packages to generate vapor-liquid equilibrium (VLE) and thermodynamic data for many process fluids and mixtures [11] [12].Because of their simplicity and accuracy, the Soave-Redlich-Kwong [13] (SRK-EOS) and the Peng-Robinson [14] (PR-EOS) cubic equations of state have been preferred.SRK-EOS and PR-EOS work well with nonpolar molecules, like hydrocarbons, but work less well with highly polar and hydrogen bonded fluids [15].Since the solute is nonpolar methane and the solvents have increasing hydrophobic character from methanol to TEG despite their having hydroxyl groups, the solute-solvent systems will illustrate the strengths and weaknesses of these cubic equations of state.
There are several reports on the VLE of binary systems of methane + thermodynamic inhibitors [16] [17] [18] [19].Jou et al. [16] [17] reported VLE of the binary systems of methane + MEG, diethylene glycol (DEG), and TEG and compared the experimental against predictions using PR-EOS.Wang et al. [18] and Abdi et al. [19] reported the VLE of the binary systems of methane + MEG.
There seems to be no available information on a comparative computational study of the VLE of these binary systems at gas hydrate formation conditions using both the SRK-EOS and PR-EOS, especially for the binary system of methane + TEG, despite the importance of TEG in natural gas dehydration.
The objective of this work is to evaluate the quality of the predictions of VLE of methane + methanol, MEG, and TEG at temperature and pressure conditions suitable for gas hydrate formation using the PR-EOS and SRK-EOS.Specifically, Modeling and Numerical Simulation of Material Science to perform a computational sensitivity analysis on the two equations of state in predicting VLE data in comparison to experimental results for binary mixtures and to ascertain the relative quality of the predictions of the PR-EOS versus the SRK-EOS for these binary systems.
Equations of state are used in commercial modeling packages to predict the behavior of systems for which there may be little or no experimental data available.The two equations of state selected for this study are classical, well regarded equations.Of course, there are many others since research on the forms and parameters used in equations of state is continual.Equations of state are used to make predictions for the design of engineering equipment.Depending on the results generated the equipment may be either over-designed making the systems needlessly expensive to manufacture, install, and operate or else may be under-designed with poor operational performance, risks to safety, and liability for any damages caused by failure of the systems.

Experimental Background
The experimental results used in the computational analysis of the VLE of the methane + methanol system were reported by Frost et al. [8] and Wang et al. [18].For Frost et al. [8], the experiments were carried out under isothermal conditions at several temperatures at which hydrates can form.The experimental apparatus was capable of handling hydrocarbon -water-hydrate systems at temperatures ranging from 213 K to 353 K and at pressures up to 40 MPa.Wang et al. [18] used a dual cell mercury-free high-pressure PVT system with a 82.7 MPa working pressure and an operating temperature range of 253 -473 K, and they reported gas solubilities in liquids.The composition of the vapor phase was not reported.
The experimental results reported for the methane + methanol binary system are presented in Table 1.
The experimental results used in the computational analysis of the VLE of methane + MEG were reported by Wang et al. [18], Jou et al. [17], and Abdi et al. [19].Wang et al. [18] used the same apparatus for MEG as for methanol.The equipment of Abdi et al. [19] was a constant-temperature chamber that was capable of working with a temperature range of −35˚C to 200˚C and at a pressure of 100 MPa.The experimental results reported for the methane + MEG binary systems are presented in Table 2.The data from all three sources are solubilities of methane in MEG liquid.The vapor phase compositions were not reported.
The experimental results used in the computational analysis of the VLE of methane + TEG were reported by Jou et al. [16] and are reproduced in Table 3.

Computational Analysis
SRK-EOS and the PR-EOS were used in the computational analysis.Since the systems are non-ideal, fugacities are used instead of pressures.The vapor-liquid equilibria using SRK-EOS and PR-EOS were calculated as described by Gmehling et al. [20] and Poling et al. [21] Calculations using the two equations of state were carried out in similar manner after appropriate changes were made to follow the functional forms and parameter values for each equation.Equation (1) is the SRK-EOS [13].Similarly, Equation ( 2) shows the PR-EOS [14].
( ) ( ) Both equations of state are variants of the classic van der Waals equation of state.The form of the PR-EOS differs from the SRK-EOS in the additional correction in the denominator of the second term that accounts for the attractive forces between molecules when volumes are small.These two equations are for pure substances.Parameters a and b are functions of the critical temperature and pressure of the substance.Parameter a also in-Modeling and Numerical Simulation of Material Science cludes additional correction factors specific to the substance.
For mixtures, a and b are functions of the composition of the system.Van der Waals mixing rules were used.
Off-diagonal a ij values are geometric averages of diagonal values including a further binary interaction coefficient K ij while off-diagonal b ij quantities are arithmetic averages.

(
)( ) ( ) Using an arithmetic average for b ij reduces the value of b for the mixture to a weighted average of b parameters for each component by the mole fraction of each.
The same mixing rule was used for both equations of state.VLE is said to exist when the fugacities of each of the components in the liquid state equal those of the components in the vapor state.Each equation of state leads to a computation of fugacity coefficients from which fugacities can be calculated.From each component the fugacity coefficient ( j ϕ ) can be determined from each equation of state: For the SRK-EOS, ( ) ( ) For the PR-EOS, ( ) ( ) In these two equations, Z is the compressibility and A and B are dimensionless coefficients containing a and b for the mixtures, 2 2

;
; In order to apply the equations of state to these binary systems, several pieces of information for each substance are required-critical temperatures and pressures and acentric factors for the corrections to the a parameter.The critical properties and acentric factors for methane, methanol, and MEG as pure components were taken from Reid et al. [22], and Galvão and Francesconi [23].The critical properties and acentric factor for pure component of TEG were taken from Gironi et al. [24].The critical values and acentric factors are tabulated in The binary interaction parameters which appear in the mixing rule of the equation of state given in Equation ( 4) were chosen from the literature.From Frost et al. [8], K ij = 0.01 was used for methane + methanol.Jou et al. [16] [17] provided K ij for methane + MEG and methane + TEG.In these papers, the binary interaction parameter is given as a linear function of temperature, K ij = a 0 + a 1 T.
The parameters are given in Table 5.
Two binary interaction parameters were used for studying methane + MEG using the numbers in Table 5 and the linear equation in temperature because of the differences in temperatures.The experiments in [18] were performed at 293.2 K yielding K ij = −0.02360.The experiments of both [19] and [17] were done at 298.15 K yielding K ij = −0.0179.For methane + TEG [22], K ij = 0.0095 at 298.15 K.

Computational Procedure
The calculation of the vapor-liquid equilibria using SRK-EOS described by Gmehling et al. [20] was the model for the calculations performed in this work.
It was implemented in an Excel spreadsheet along with a cubic equation solver needed to calculate the molar volume from the cubic equations of state.Two spreadsheets were prepared, one for each equation of state, because of differences in calculating parameters and in the equations of state themselves.The spreadsheets can calculate the VLE data of any binary system if the liquid phase mole fraction, temperature, pressure, initial vapor composition, critical conditions and acentric factors are available.To test for the accuracy of these spreadsheets, they were used to reproduce the VLE calculation done on the binary system of nitrogen + methane [20].The computational procedure consisted of the following steps: 1) Calculate the reduced temperatures of each component.
2) Calculate EOS parameters pertaining to pure components and for the mixtures which do not depend on composition.
3) For the liquid phase-Performed once since it is assumed that the liquid composition would remain fixed while the system's vapor composition is calculated to self-consistency.a) Calculate mixture parameters for the liquid state.These will depend on composition.
C. E.  ϕ ) for each component in the liquid phase using either Equation (7) or Equation ( 8) depending on the EOS used.The new pressure will be used as input to the next iteration.
g) Test for self-consistency by determining whether the calculated value of S is within 10 −4 of unity.When this has been achieved, terminate the calculation.
The computational procedure is summarized in the following flowchart (Figure 1).

Test Calculation-Binary System of Nitrogen + Methane
The test calculation using nitrogen + methane converged in eight passes using SRK-EOS and seven passes using PR-EOS.Since the calculations were performed using a spreadsheet, each cycle was tallied by hand.This procedure was used in all of the results reported here.Although each pass required manual input of results of the previous pass, the rapid convergence of the calculation indicated that a more sophisticated program was not needed.Results of the test calculations are reported in Table 6.Both equations of state reproduced the experimental conditions well since both substances consist of simple, nonpolar molecules.
The spreadsheet templates for SRK-EOS and PR-EOS were used to calculate the VLE of methane + methanol, methane + MEG, and methane + TEG.All of the calculations converged.In none of the calculations were more than twenty cycles required to achieve convergence.

Binary System of Methane + Methanol
There are two sets of reference experimental data [8] [18] of the methane + methanol system to provide initial data for the calculations using the equations of state.Figure 2 plots the experimental and computational results by showing the resulting system pressures as a function of the mole fraction of methane in the liquid phase.
Frost et al. [8] reported the experimental composition of the vapor phase along with the composition of the liquid phase.The computations also returned  the compositions of both phases.The compositions corresponding to the pressure values from reference [8] in Figure 2 are shown in Figure 3.
The equations of state overestimate the pressure of the system and underestimate the concentration of methane in the vapor phase.The close agreement of the experimental data from [8] and [18] in the range they share validates experiments.The SRK-EOS consistently yields higher pressures than the PR-EOS.The vapor phase composition is very close to being almost exclusively methane in both experiment and in the models.Despite the slight differences in the calculated mole fractions of methane vs. the experimental value (2% -3% less), the corresponding differences in pressure are quite large-two to three times the experimental value.

Binary System of Methane + MEG
Wang et al. [18], Jou et al. [17], and Abdi et al. [19] have all reported experimental VLE data on the methane + MEG system.Experimental and modeled system pressures based on data in the respective experiments are graphically displayed in Figure 4.
Figure 4 shows the consistency of the three sets of experimental results.All of them fall close to the same line.The corresponding sets of lines for the PR-EOS and SRK-EOS modeled results are similarly clustered falling cleanly on the same trend lines.Contrary to the results of the methane + methanol system, the pressures produced from the models based on the two equations of state are underestimated, falling below the experimental line with the PR-EOS producing results at distances further than the SRK-EOS.

Discussion
There are many observations and comments in the literature stating that the equations of state are most appropriate for nonpolar substances [15] [25].This is why these calculations converged for methanol and the glycols.These solvents are polar, but they do not behave like water.Clathrates formed from water require the presence of two hydrogen bonds to form the cage encapsulating the methane molecule.Replacing a hydrogen atom in water with a methyl group to produce methanol breaks the tendency to form lattices.The glycols are too extended and flexible to form rigid lattices like water.
The calculated pressures of methane as a function of the concentration of methane in the liquid solvents are overestimated for methanol and underestimated for the glycols.Another way of describing this phenomenon is the observation that the equations of state require a higher pressure of methane (since the vapor phase is nearly exclusively methane) to achieve the same concentration of methane in methanol as experiment and that they require a lower pressure of methane to achieve the same concentration of methane in the glycols.The pressures of methane needed to achieve the same concentrations are proxies for the solubility of methane in the liquids.Higher pressure to achieve the same concentration implies a lower solubility of the gas.Lower pressure implies a higher solubility.This is one way of expressing Henry's Law.
Henry's Law states that the concentration of a gas in a liquid solvent is proportional to its partial pressure in the vapor phase at low concentrations of the gas.The constant of proportionality is the Henry's Law constant for the system.
The constant is dependent on temperature and the nature of the substances.
There are many forms of Henry's Law depending on the way the relationship is stated and the units of pressure and concentration selected.In this work, Henry's Law will be stated as P i = H x i where H is the Henry's Law constant.
The graphs relating the pressure of methane to its mole fraction in methanol, MEG, and TEG allow the extraction of the Henry's Law constant at infinite dilution.Trend lines for each curve in the Figure 2, Figure 4, and Figure 5 were determined.Both linear and quadratic functions were used to fit the curves.The trend lines were constrained to go through the origin.The quadratic fits were of better quality, and the linear coefficient of the quadratic function fit was selected as the Henry's Law constant.The values from the trend lines are given in Table 7.
A set of values of the Henry's Law constant for methane in methanol have been reported and graphed by Horsch et al. [26] in a paper describing the molecular modeling of hydrogen bonding fluids.There were ten reports from which Henry's law constants were extracted, and they ranged from 75 MPa to 120 MPa at 298 K.The experimental results of [8] [18] fit within this range at the high end using either the linear or the quadratic fit for H.The quadratic fit produces consistent values of the Henry's Law constant for the experimental data in [8] and the calculations using the equations of state based on the same data.Jou et al. reported Henry's Law constants at 298 K for MEG [17] and TEG [16] of 656.1 MPa and 179.2 MPa, respectively.The values calculated here disagree with Jou's.This may be due to how they were determined.In this report, simple linear and quadratic fits to the solubility data were used.Jou reported the Henry's Law constants as part of a set of parameters in the Krichevsky-Ilinskaya equation.The manner in which they were obtained was not described.Nevertheless, the calculation of the Henry's Law constants from Jou's experimental data used in this report yields values close to the numbers in [16] [17].The Henry's Law constants obtained from the equations of state are less than the experimental values but are in the order of the proximity of the calculated curves to the experimental values.SRK-EOS gives closer agreement than PR-EOS to experiment despite its simpler functional form.This may be the result of greater difference in the sizes of the solute and solvent molecules.The extra correction in the PR-EOS may over-correct and is unnecessary when the solute and solvent molecules are very different in size.This may also be the reason why PR-EOS gave better results than SRK-EOS for methane + methanol since the molecules are much closer in size giving the corrections inserted in the PR-EOS greater relevance.

4 )
For the vapor phase-Performed iteratively until self-consistency.The procedure begins with the experimental pressure and vapor composition as the initial condition.If the model were accurate, the model should return the pressure and vapor composition of the initial condition.a) Calculate mixture parameters for the vapor state.These will depend on composition.b) Calculate the vapor phase molar volume (for SRK-EOS) or molar volume and compressibility (for PR-EOS) by solving the pertinent EOS in the form of a cubic equation in the volume.c)Calculate the fugacity coefficients ( Vi ϕ ) for each component in the vapor phase using either Equation(7) or Equation (8) depending on the EOS used.d) Calculate the vapor phase composition (y i ) from the liquid phase composition (x i ) and liquid and vapor fugacity coefficients.calculated vapor phase mole fractions so they sum to one.a new total pressure by multiplying the input pressure by factor S.

Figure 1 .
Figure 1.Flow chart of the computational procedure.

Figure 2 .
Figure 2. Experimental and Modeled System Pressures for Methane + Methanol.

Figure 5
Figure5displays the experimental[16] and modeled system pressures as a function

Figure 3 .
Figure 3. Experimental and Modeled System Compositions for Methane + Methanol.

Figure 4 .
Figure 4. Experimental and Modeled System Pressures for Methane + MEG.

Figure 5 .
Figure 5. Experimental and Modeled System Pressures for Methane + TEG.

Table 1 .
Experimental Vapor-Liquid Equilibrium Data for the Methane + Methanol System.

Table 2 .
Experimental Vapor-Liquid Equilibrium Data for the Methane + MEG.

Table 3 .
Experimental Vapor-Liquid Equilibrium Data for the Methane + TEG System.

Table 6 .
Test Calculation of the Nitrogen + Methane System.

Table 7 .
Henry's Law Constants from Trend lines for Methane in Methanol, MEG, and TEG.