Modeling the Thermodynamics of Fluids Treated ... - ACS Publications

Jan 25, 2017 - CNRS 7274), 1 rue Grandville, 54000 Nancy, France. •S Supporting Information. ABSTRACT: It is currently acknowledged that thermody-...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/IECR

Modeling the Thermodynamics of Fluids Treated by CO2 Capture Processes with Peng−Robinson + Residual Helmholtz Energy-Based Mixing Rules Silvia Lasala,*,†,‡ Paolo Chiesa,† Romain Privat,‡ and Jean-Noel̈ Jaubert‡ †

Politecnico di Milano, Department of Energy, via Lambruschini, 4, 20156 Milan, Italy Université de Lorraine, Ecole Nationale Supérieure des Industries Chimiques, Laboratoire Réactions et Génie des Procédés (UMR CNRS 7274), 1 rue Grandville, 54000 Nancy, France



S Supporting Information *

ABSTRACT: It is currently acknowledged that thermodynamic models routinely used to model the behavior of fluids involved in CO2 capture and storage processes are not sufficiently accurate. Such a deficiency is one of the main sources of uncertainty for estimating the costs associated with these technologies and, as a direct consequence, engineers largely oversize the equipment. In order to reduce the uncertainty related to the calculation of thermodynamic properties, this work aims at providing the optimal parameters of the model Peng−Robinson + equation of state (EoS)/ aE,γ‑Wilson mixing rules that are currently able to accurately res correlate the phase behavior of systems encountered in CCS technologies. In particular, this model is optimized in this work over the experimental vapor−liquid equilibria (VLE) data of 23 binary systems, resulting from the binary combination of CO2, H2, N2, O2, Ar, CO, CH4, H2S, and H2O, for which data have been found in the literature. The paper shows that this equation of state is reliably applicable to calculate VLE properties of these binary mixtures and suggests the application of this model to represent the thermodynamics of multicomponent fluids treated by CO2 capture processes. Another outcome of this work is the conceptual framework outlined to enable the optimization of thermodynamic models based on the Maximum Likelihood Method, accounting for the correlation and joint uncertainty of model parameters. been named the Peng−Robinson + EoS/aE,γ res mixing rules, since it joins the standard Peng−Robinson (PR78) cubic equation of state (EoS)7,8 and an advanced class of mixing rules (EoS/aE,γ res ) that originates by coupling an equation of state (in this case, PR78) to the residual part of an activity-coefficient model (γmodel), that actually represents the residual part of an excess Helmholtz energy (aEres) model. This class of mixing rules is thoroughly described in ref 6. The main purpose of this paper is to provide the results obtained from the optimization of the two parameters of the model PR+EoS/aE,γ res mixing rules, thus enabling its application to mixtures containing CO2, N2, Ar, O2, CO, H2, H2O, H2S, and CH4. In particular, such an optimization has been based on the Maximum Likelihood Method (MLM), to allow the statistically sound estimation of parameters and their uncertainty, with the aim of providing, to future users of the model, the actual

1. INTRODUCTION Two connected fields of engineering studies where thermodynamics plays a dominant role are enhanced oil recovery (EOR) and CO2 capture and sequestration (CCS) processes. Their common ground is represented by the necessity of accurately modeling the thermodynamic properties of mixtures containing CO2, N2, CO, Ar, O2, CH4, H2, H2S, and H2O.1−5 The stillunsolved problem lies in the thermodynamic modeling of CO2based mixtures containing highly volatile impurities, such as N2, CO, Ar, O2, CH4, and H2. However, mixtures involved in CCS and EOR processes might not necessarily be composed of, mainly, CO2. This is the case, for instance, of fluids treated by processes for the purification of CO2 where, by the nature of these technologies, there is the necessity of modeling the thermodynamic properties of streams containing CO2, N2, CO, Ar, O2, CH4, H2, H2S, H2O, and other components, in a wide range of possible compositions. Thermodynamic models should thus be able to accurately represent the properties of all the possible combinations of these species. In a previous paper,6 we presented a thermodynamic model accurate in calculating vapor−liquid equilibrium properties of binary systems containing CO2 and N2, O2, or Ar. This model has © 2017 American Chemical Society

Received: Revised: Accepted: Published: 2259

November 1, 2016 December 17, 2016 January 25, 2017 January 25, 2017 DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research ⎧ ai(T ) = ac , iαi(T ) ⎪ ⎪ α (T ) = [1 + m (1 − T /T )]2 i c ,i ⎪ i ⎪ 2 ⎪ mi = 0.37464 + 1.54226ωi − 0.26992ωi ⎪ ⎛ R2T 2 ⎞ ⎪ c ,i ⎪ ac , i = Ωa⎜⎜ P ⎟⎟ ⎝ c ,i ⎠ ⎪ ⎪ ⎪ ⎛ RTc , i ⎞ ⎨ ⎪bi = Ωb⎜⎜ P ⎟⎟ ⎝ c ,i ⎠ ⎪ ⎪ 1/3 1/3 ⎪ X = −1 + (6 2 + 8) − (6 2 − 8) ⎪ 3 ⎪ ⎪ Ω = 8 + 40X ≈ 0.45724 ⎪ a 49 − 37X ⎪ ⎪ Ω = X ≈ 0.07780 ⎪ b ⎩ X+3

uncertainty related to its application. The applied optimization methodology is also presented in detail in this work. In particular, after recalling in section 2 the theory at the basis of the considered thermodynamic model, section 3 presents the procedure of model optimization, applied over vapor−liquid equilibrium (VLE) data and based on the MLM, referring to a more detailed description in Appendix I in the Supporting Information. To help to analyze optimization results, the same section also describes the link between joint confidence regions of optimal parameters and two quantities that can be more easily determined and synthetically reported in optimization results (with respect to confidence regions): the factor of correlation between them and the individual disjoint uncertainty of each parameter. The formula for the mathematical determination of these properties is reported in Appendix II in the Supporting Information. Furthermore, section 4 classifies the 36 binary

(3)

23 of these mixtures are reviewed and analyzed, together with the

At this step, we note that setting s = 1 and selecting a Van Laar expression for aE,γ res /ΛEoS in eq 2 makes it possible to derive the mixing rules used in the well-established PPR78 model.10−12 The pure components properties (critical temperature (Tc), pressure (Pc), and acentric factor (ω)) considered in this work are reported in Table 1.

effects that such features may have in the optimization process.

Table 1. Pure Components Properties Applied in the Modela

mixtures at the basis of the modeling of multicomponent systems composed of CO2, N2, Ar, O2, CO, H2, H2O, H2S, and CH4. In the same section, the characteristics of the VLE data available for

Finally, section 5 reports the optimization results.

2. THE MODEL Thermodynamic models available in the literature contain parameters, more or less based on physical foundation, that require an adjustment over experimental data. The model considered in this work has been presented by the same authors in a previous paper.6 This model is based on the Peng−Robinson a

EoS, P(T , v , z) =

am RT − 2 v − bm v + 2bmv − bm 2

(1)

Tc [K]

Pc [MPa]

ω

H2O H2S CO2 CH4 Ar O2 CO N2 H2

647.096 373.53 304.21 190.564 150.86 154.58 132.92 126.19 33.19

22.064 8.96291 7.383 4.599 4.868 5.043 3.499 3.3958 1.313

0.344861 0.0985 0.22394 0.011548 0.0000 0.02218 0.048162 0.0372 −0.21599

All properties have been collected from the DIPPR database.

In this study, the Soave α-function was selected because, as shown in refs 13 and 14, it holds the thermodynamically consistent positive−monotonically decreasing−convex temperature dependency, within the temperature domains that characterize the systems of the engineering applications treated in this study. Nonetheless, in general, other types of α-functions could be selected15,16 to improve the accuracy in modeling the thermodynamic properties of pure components. However, the application of α-functions should be always supported by the preliminary evaluation of their thermodynamic consistency, in the entire temperature domain of interest. Other characteristics of the applied mixing rules (eq 2) are specified in the following: • parameter s has been set to 2, as proposed by Voutsas et al. when they developed the UMR-PR (Universal Mixing Rule− Peng−Robinson) model,17 • parameter Λ, present in the mixing rule for am, has been considered to be equal to −0.52398, which is the value determined by Michelsen along the derivation of the MHV118 mixing rule, as justified in ref 6, and

and differs from the original formulation of Peng and Robinson uniquely in the applied mixing rules (MR), which, in this case, are the EoS/aE,γ res MR (see ref 6 for more details about its derivation), with the particular selection of the quadratic mixing rule and Lorentz9 combining rule for the covolume: ⎧ NC NC ⎛ b 1/ s + b 1/ s ⎞s i j ⎪b = ⎟ ∑ ∑ zizj⎜⎜ m ⎟ ⎪ 2 ⎪ ⎝ ⎠ i=1 j=1 ⎨ NC ⎛ ⎞ ⎪a a a E, γ ⎪ m = ∑ zi⎜ i ⎟ + res ⎪ bm ΛEoS b ⎩ i=1 ⎝ i ⎠

component

(2)

where pure component energy and covolume parameters, ai and bi, are calculated as in the original PR model: 2260

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research • the Wilson’s activity coefficient model has been applied to derive the expression for the residual part of the excess Helmholtz energy model, aE,γ res :

According to Gibbs’ phase rule, there are two independent variables in each binary VLE data point, leading to six possible combinations of independent variables: (T,P), (x1,y1), (P,x1), (P,y1), (T,x1), and (T,y1). A wide number of methods has been proposed for the reduction of VLE data. They mainly differ in the way the objective function is defined. Objective functions (OFs) applied to optimize models over P−T−x−y binary data always require the calculation of at least one property, as a function of one of these six couples of independent variables. For example, leastsquared functions usually applied to binary systems are

E,γ‐Wilson E,γ‐Wilson ares (T , z) = a E,γ‐Wilson − acomb NC NC NC ⎛Φ ⎞ = −∑ zi ln(∑ zj Ωji(T )) − ∑ zi ln⎜ i ⎟ ⎝z ⎠ j=1 i = 1 =1 i i    a E,γ‐Wilson

E,γ‐Wilson acomb

(4)

where Φi is the volume fraction of the ith molecule, zv Φi = NCi i ∑ j = 1 zjvj

⎡ nexp OF = ⎢∑ (Piexp − Picalc)2 + ⎢⎣ i = 1

and Ωji(T) is the local composition parameter, expressed as an exponential function of two adjustable parameters, Aij and Aji, a priori temperature-independent and proportional to molecular interaction energies (Aji = NA(εij − εjj)). Ωji(T ) =

Picalc

nexp



i=1

⎥⎦

∑ (y1,expi − y1,calci )2 ⎥

(5)

calc y1,i

where and are determined from bubble point exp calculations, as a function of Texp i and x1,i ; and ⎡ nexp calc 2 OF = ⎢∑ (x1,exp i − x1, i ) + ⎢⎣ i = 1

⎛ Aij ⎞ exp⎜ − ⎟ vi ⎝ T⎠

vj

nexp



i=1

⎥⎦

∑ (y1,expi − y1,calci )2 ⎥

xcalc i

where and ycalc 1,i are determined exp as a function of Texp i and Pi .

Therefore, the application of these mixing rules requires optimization of the two binary interaction parameters contained in the Wilson γ-model: Aij and Aji. This model thus describes, as does most of the available equations of state, the thermodynamics of multicomponent systems by properly combining the effects of pairwise interactions between bodies of equal or different chemical species. Where, in particular, “pairwise interactions between bodies of equal chemical species” are the typical interactions of a pure-component system. Since these parameters, empirically, reflect binary molecular interaction forces, it is recommended to fit model parameters over experimental data which macroscopically reflect the theoretical significance of these parameters. First, it is necessary to optimize binary interaction parameters over data of binary mixtures. If, differently, binary interaction parameters are directly fitted over data of the studied multicomponent system, instead of binary data of its binary subsystems, predictions of the considered multicomponent system clearly improve. However, this fitting method could introduce, in the regressed binary interaction parameters, the unwanted contribution of other types of interactions. Second, it is necessary to optimize binary interaction parameters over data that reflect the nature of the intermolecular forces between molecules, such as saturation data. In fact, as it is wellknown, intermolecular bonds govern the state of aggregation of matter. Vapor pressure, which is an indicator of the tendency of molecules to escape from the liquid state, is a property that thus “testifies to a cohesive force between the particles of matter”. This property was one of the first properties “to be related convincingly to what now we call intermolecular forces”.19

(6)

from PT-flash calculations,

In other cases, each residual of least-squares objective function might be weighted over the uncertainty of the experimental exp property, uexp f i , or the experimental property itself, f i , 1 OF = nexp

⎛ f exp − f calc ⎞2 ∑ ⎜⎜ i exp i ⎟⎟ uf ⎠ i=1 ⎝ i

1 OF = nexp

⎛ f exp − f calc ⎞2 ∑ ⎜⎜ i exp i ⎟⎟ fi ⎠ i=1 ⎝

nexp

or nexp

These objective functions are usually referenced with the expression “weighted least squares” objective functions. Most of the optimization procedures currently applied rely on least squared minimizations by the use of simple or weighted objective functions. In other cases, authors simply minimize OF =

1 nexp

nexp



f iexp − f icalc

i=1

or OF =

1 nexp

nexp

∑ i=1

f iexp − f icalc f iexp

Results obtained by the implementation of different objective functions differ one from each other. To cite an example, the application of OF in eq 5 leads to the minimization of errors between experimental data and calculations of pressure and vapor-phase composition, regardless of deviations between experimental data and calculations of temperatures and liquidphase composition. When thermodynamic models under optimization systematically deviate from data, even whether optimized over them, no one available optimization model could lead to the statistically optimal estimation of parameters. In these cases, the selection of the proper objective function could be performed as to guarantee the minimization of the errors over the thermodynamic properties of interest.

3. MODEL OPTIMIZATION Generally, the optimization of models consists of the determination of parameters that minimize residuals between experimental data and model calculations. For the reasons detailed in the previous paragraph, the optimization of the couple of parameters (A12, A21) of the considered model is performed in this work over vapor−liquid equilibrium data of binary systems. All the measurable properties of binary systems in each ith VLE state (i = 1, ..., ne) are total pressure Pi, temperature Ti, liquid phase mole fraction x1,i and vapor phase mole fraction y1,i. 2261

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

principle, introduced by Fisher in 1922,21 is a methodology that enables the statistically sound estimation of parameters contained in a model. Maximum Likelihood Method (MLM) is based on different general assumptions, which must be considered, in particular, (1) experimental observations at the basis of the optimization are independent and described by the same, a priori unknown, distribution; (2) experimental errors are larger than those introduced by models that, in other words, means that models under calibration do not introduce systematic deviations. More details about the assumptions underlying this method, contextualized to its application to VLE measurements, are reported in the following part of the section, in Appendix I in the Supporting Information and in refs 20 and 22−24. It is worth anticipating here that, when at least one of the assumptions upon which MLM is based is not valid, its application−in terms of its statistical foundation−is equivalent to the application of a weighted least-squared method. In such a case (not-satisfaction of MLM assumptions), the model optimization may be carried out by conveniently applying the objective function that mostly suits our needs: if the aim is to minimize, at most, the residuals between model calculations and experimental data, without accounting for the data uncertainties, it should be applied the (simple or weighted) least-squared method that mostly minimize such residuals; in contrast, if the objective is to minimize residuals by weighting data over experimental uncertainties, the weighted least-squared method should be more conveniently applied. A general treating of the maximum likelihood method is featured in Appendix I in the Supporting Information, together with its specific application in the optimization of models over VLE data of the binary systems applied in this work. The maximum likelihood function minimized in this study for the optimization of the binary interaction parameters, β = (A12,A21) over P−T−x−y data of the equation of state PR+EoS/aE,γ‑Wilson res of binary systems, is given as follows:

However, if a model is able to represent experimental data without introducing systematic errors, parameters may be optimized with statistical optimization methods, accounting for the uncertainty of measurement data. In particular, when parameters under optimization are physically grounded, as it is for the considered equation of state, their statistically optimal estimate reflects the physical meaning of the parameters themselves, thus helping to develop predictive thermodynamic models. The most powerful methods for the estimation of statistically based parameters rely on the maximum likelihood principle. In fact, differently from methods based on least squared minimization, if the assumptions under which this method is based on are valid, results obtained (i.e., optimal parameters) from maximum likelihood methods do not depend on the selected independent variables used to calculate properties within the objective function or, more in general, on the objective function itself.20 With reference to this observation, it is remarked that, as shown in a previous publication of the same authors,6 the studied thermodynamic model enables the accurate representation of the phase behavior of different high-pressure systems studied in this work, such as the binary mixtures CO2−N2, CO2−O2, and CO2− Ar. In this work, the optimization of the equation of state presented in section 2 has thus been performed over VLE data of binary mixtures obtained from the binary combination of CO2, N2, Ar, O2, CO, H2, H2O, H2S, CH4, by applying the maximum likelihood method, presented in the following subsection 3.1 and, in more detail, in Appendix I in the Supporting Information. Section 3.2 explains the concept of correlation between parameters, their individual uncertainties, and their joint confidence region, leading to Appendix II in the Supporting Information that explains the mathematical description of their formulations. 3.1. The Optimization Methodology Determined by Maximum Likelihood Principle. The maximum likelihood ⎧ ⎡ exp calc 2 n ⎪ ⎪ 1 ⎢ exp ⎛⎜ y1, i − y1, i ⎞⎟ ⎢∑ min ⎨ ⎜ ⎟ + A12 , A 21⎪ nexp ⎢ σy i=1 ⎝ ⎠ 1, i ⎪ ⎣ ⎩

⎛ calc 2 (x1,exp − y1,calc (y1,exp ) cov(x1, i , y1, i ) i − x1, i )σy1, i ⎜ i i − ∑⎜ σ 2σ 2 − cov(x1, i , y1, i )2 σy σx1,i 2σy 2 − cov(x1, i , y1, i )2 i = 1 ⎜ σy ⎝ 1,i x1,i y1,i 1, i 1, i nexp

Variances σx1,i2,σy1,i2, and the covariance cov(x1,i, y1,i) are calculated as detailed in Appendix I in the Supporting Information (see eqs (A7)−(A9)). 3.2. The Correlation and Uncertainty of Optimal Parameters. When more parameters are contemporarily optimized in a model, it is necessary to evaluate their mutual correlation, in order to attest the presence of relationships between them and the dependence of the optimal results on the applied optimization method. A high degree of correlation between two parameters underlies that the use, in the model, of relationships between them could lead to represent experimental data as accurate as the use of two individual parameters. The analysis of this characteristic is fundamental to evaluate the occurrence of overparameterized situations (ill-conditioned optimization cases) and to eventually detect functional relationships between parameters. Also, it is always necessary to evaluate the uncertainty of parameters under regression that, in the presence of correlation between parameters, must be assessed by considering, jointly, the optimal couple of parameters.

⎞ 2 ⎤⎫ ⎪ ⎟ ⎥⎪ ⎟ ⎥⎬ ⎟ ⎥⎪ ⎠ ⎦⎪ ⎭

(7)

Information about correlation and uncertainty of optimal parameters are all contained in the so-called “joint confidence region” of optimal parameters (A12,A21), which corresponds to the region in the space (A12,A21) that most probably contains the true values of parameters. In other words, the complete output of a multiparameter optimization process should always provide the optimal discrete value of each parameter, in this case, the optimal values for A12 and A21, and the joint confidence region of these parameters. The contour of these regions associated with a specific confidence level α, OF(β)α, can be determined as shown in Appendix II in the Supporting Information, as a function of the minimum value of the objective function, OF(βopt). These regions are qualitatively represented in Figure 1. More specifically, if the model is linear in its parameters (that is not our case), the joint-confidence region given by eq (A16) in Appendix II in the Supporting Information, is characterized by an exact elliptical shape (Figure 1a). In contrast, if the model is nonlinear in its parameters (as it is in our case), the contour plot 2262

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

(ii) the minimum value of the objective function, in correspondence with the optimal set of parameters, OF(βopt); and (iii) the curvature of the objective function, whose intersection with OF(β)α defines the contour of the confidence region. In particular, given a certain risk level, α: (i) the lower the number of the experimental data, the greater is the extension of the confidence region; (ii) the lower the minimum value of the objective function OF(βopt), the lower is the resulting maximum objective function OF(β)α, relative to the contour that defines the confidence region; and (iii) (for objective functions that account for the uncertainty of experimental data), the higher are the uncertainties associated with experimental data, the lower are gradients of the objective function with respect to the parameters under optimization, thus resulting in more extended confidence regions. In this work, some joint confidence regions will be presented, with the objective to show the effect that some of the mentioned characteristics of the optimization process (number of experimental data, objective function and thermodynamic system under consideration, etc.) have on the uncertainty and correlation of optimal couples of binary parameters. To lighten the presentation of the numerous results that need to be reported in the following, the representation of the joint confidence region for each couple of optimal parameters will be replaced by two calculated quantities that, when considered together, help to visualize the joint-confidence region of the same parameters: • the correlation factor between two optimal parameters A12 and A21, denoted as c and calculated as shown in Appendix II in the Supporting Information; • the approximated individual confidence intervals for the true value of each parameter, calculated irrespective of the value of the second parameter, as shown in Appendix II in the Supporting Information. To better explain the meaning of these properties, we ask the reader to follow the discussion below by making reference to Figure 2, which reports different types of joint confidence regions

Figure 1. Example of elliptical joint confidence regions: (a) elliptical objective function (OF(β)α) contours for a model linear in its parameters, and (b) irregular OF(β)α contours for a nonlinear model, with more than one possible local minimum.

given by eq (A16) in Appendix II in the Supporting Information will provide a proper nonelliptical confidence region, represented by a banana-shaped or quasi-elliptical form (Figure 1b); in this nonlinear case, however, the unknown statistical distribution associated with its properties entails an uncertainty related to the probability level.24 The visualization of the joint-confidence region of parameters is of primary importance to attest the reliability of the optimization results and to study the functional relationship between parameters and other variables, such as the temperature for the binary interaction parameters considered in this study, as shown in section 4. More explicitly, if confidence regions of two couples of parameters, optimized at two different temperatures, do not overlap (or, with a more specific terminology, are not compatible), then it means that the two couples of parameters identify dif ferent points in the space (A12,A21), which are, therefore, temperature-dependent. It can be observed from eq (A16) in Appendix II in the Supporting Information that, once the α-confidence level is fixed, the extension and shape of the confidence region become dependent on (i) the number of experimental measurements (nexp) used in the optimization process;

Figure 2. Joint confidence regions of parameters characterized by different degrees of correlation. 2263

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research of general parameters β1 and β2. To simplify the analysis, regions shown in Figure 2 result from quasi-linear models, as previously shown in Figure 1a. If the optimization problem is nonlinear, the discussion does not change. The only difference is that parameters will be related, in general, by a nonlinear function, as is typically observed for the A12 and A21 parameters studied in this work. The factor of correlation between two parameters quantifies the effect, on the objective function, of their joint variation (i.e., the covariance), with respect to their individual variation. Because of its formulation (see eq (A13) in Appendix II in the Supporting Information), it may vary between −1 and 1: zero c-values correspond to uncorrelated parameters (Figure 2 (a−c)), in case both of them play a relevant role in the optimization problem (Figure 2 (a)) or, differently, one of them is irrelevant and the other controls the result of the optimization (Figure 2 (b,c)); unitary |c|-values correspond to correlated parameters (Figure 2 (d−g)); intermediate |c|-values (0 < |c| < 1) indicate a different degree of correlation between parameter (e.g., Figure 2 (h,i) are characterized by a correlation factor of about 0.8). Furthermore, it is highlighted that drawing joint confidence regions helps to detect ill-conditioned overparameterized optimization problems that are usually reflected in highly elongated shapes preferentially oriented alongside the axis of one of the two parameters (i.e., regions in Figures 2b and 2c). This mutual correlation between parameters (factor c) is therefore reflected in the shape of their joint confidence region and, in particular, in the extension of the region, around a functional relationship. Moreover, the sign of c indicates whether the major axis of the region has a positive or negative slope (e.g., c < 0 in Figures 2d, 2e, and 2h), c > 0 in Figures 2f, 2g, and 2i); the module of c does not provide any information about the slope of this functional relationship. Consequently, the correlation factor does not provide sufficient information about joint confidence regions of parameters, when those are not drawn. To avoid the representation of each confidence region, it is important to associate to the correlation factor c the calculation of the approximated individual confidence intervals for the parameters, u+βi and u−βi, each one calculated independently of the value of the second parameter. Their determination is described in Appendix II in the Supporting Information (eqs (A14, A15)). The joint consideration of these disjoint uncertainties leads to the representation of the rectangular region shown in Figure 3, which, however, well approximates the exact region only when parameters are uncorrelated (i.e., Figures 2a−c).

Knowledge of the correlation factor and of the disjoint uncertainty of the parameters are useful indicators that help to visualize the direction of the slope of the exact joint-confidence region and its extension. For this reason, the table presented at the end of the paper (Table 4), which contains the numerical optimal value of each couple of parameters A12,A21, will also include the correlation factor (c) and individual disjoint confidence intervals.

4. CHARACTERISTICS OF DATA BASIS ON THE MODEL OPTIMIZATION 4.1. The Classification of VLE Systems under Study. Prior to starting with the presentation of the characteristics of data available in the literature, this section preliminarily classifies the 36 binary systems that result from the binary combination of the 9 considered species, based on their phase behavior. In fact, the phase diagrams of these binary mixtures highly differ from system to system. The classification criterion applied in this section was proposed in 1968, by van Konynenburg and Scott,25,26 based on the peculiarities observed in the phase behavior of binary systems. The preliminary part of this classification distinguishes between systems of class 1 and of class 2, based on the ratio between critical temperatures of the pure components that form the binary mixture; in particular, they discriminated binary mixtures based on the continuity (class 1) or discontinuity (class 2) of the locus of critical points. Considering component 1 as being the most volatile, if Tc2/Tc1 < 2, the system is class 1 and the locus of critical points is continuous; otherwise, the system is class 2 and the locus of the critical points is discontinuous. In particular, all class 2 binary systems treated in this paper are characterized by infinite discontinuity (class 2, type III systems): at a specific composition, critical pressure tends toward infinity. The limit Tc,2/Tc,1 = 2 just gives an idea of a general number around which systems begin to exhibit a completely different phase behavior. For systems where Tc,2/Tc,1 ≈ 2 ± 0.3, this criterion is not sufficiently precise to apply such a classification. However, it should be kept in mind that the complexity in modeling the critical region increases as Tc2/Tc1 increases. As previously anticipated, this classification criterion has been applied in this work to preliminarily comprehend the difficulty in measuring and modeling the critical region of these systems. Table 2 shows the ratio between the critical temperatures of the two components of each of the 36 considered binary mixtures. The analysis of these values shows that there are 12 systems of class 1 (40%), 19 systems of class 2 (53%), and 5 systems that have a Tc2/Tc1 ratio close to 2 and an insufficient number of experimental data to attest what class they belong to. For most of these systems, it has been possible to confirm the behavior predicted from the criterion of van Konyenburg and Scott by examining experimental data found in the literature. Since more than half of the considered number of binary systems is characterized by a discontinuous locus of critical points, the measurement of their phase diagrams requires the achievement of extreme experimental conditions, because of the reached very high-pressure level (with respect to the saturation pressure of the less-volatile component). 4.2. The Analysis of Available Data. An extensive review of available VLE data has been carried out considering references provided by the NIST database, integrated in the process simulator Aspen Plus, and works of revision on high-pressure phase-equilibrium measurements.27−31 All specific references are reported in Appendix IV in the Supporting Information.

Figure 3. Exact joint confidence region of parameters (green) and rectangular region obtained by considering the individual (disjoint) confidence intervals for the parameters (yellow). 2264

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research Table 2. Ratio of Critical Temperatures of the Two Components of the Mixture (Tc2/Tc1)a

a

Data given in bold blue font represent the ratio of critical temperatures of class 1 systems; data given in bold red font represent the ratio of critical temperatures of class 2 systems; data given in italic black font represent the ratio of critical temperatures of systems of undefined class (class 1 or class 2).

Table 3. Sets of VLE Data Points (P, T, x, y), for Each Binary Mixture, Available in the Literaturea

a With reference to the heading of the table, “min” and “max” refer to the minimum and maximum experimental values that are reported within the literature for the specific property. x1 and y1 are the molar fractions of the first component reported in the label of each mixture in the column labelled “binary system”.

whose characteristics are synthetically reported in Table 3, leads to the considerations detailed in the following sections, that will be useful in setting up the optimization process and in analyzing optimization results. 4.2.1. Types of Experimental Phase Diagrams. The ranges of temperature, pressure, and molar fractions covered by data points for the considered systems follow their type of phase behavior, as detailed with the following qualitative classification. In Table 3, systems indicated in blue represent systems for which at least one of the considered temperatures is lower than the critical temperature of both of the components contained in the mixture (Figure 4a). In systems presented in red (both normal and italic font styles), the most volatile component is always supercritical; for these systems, the presence of a critical point prevents vapor and liquid molar fractions of the more-volatile component from achieving the unitary value (Figure 4b). Except

Among the considered 36 binary systems, only for 23 of them data are available in a number sufficiently high to enable the optimization of the two parameters of the model. For each of these 23 systems, Table 3 collects available complete data points (i.e., data points for which measured properties are temperature (T), pressure (P), molar composition of both vapor (y) and liquid (x) phases) and presents some of their characteristics. In particular, for each system, it reports: the total number of measured isothermal phase diagrams (NI), the number of isotherms (NIT) at different temperatures (|ΔT| > 0.5 K), the total number of complete data points (NP), the average number of data points for each measured isotherm (the lower integer part of NP/NI), the average number of data points at similar temperatures (T ± 0.5 K) (the lower integer part of NP/NIT), the ranges of temperature, pressure, and molar fractions of the vapor and of the liquid phases. The analysis of the sets of data, 2265

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

Figure 4. Typical isothermal phase diagrams treated in this study.

Figure 5. Joint-confidence regions OF(β)α, with indefinite (red contour, OF(β)α) and limited (blue contour, OF(β)α) extensions, with the indication of the optimal point (in green) characterized by the optimal value of the objective function OF(βopt) (panel (a)); objective function curve as a function of β1, at constant β2, and objective function curve as a function of β1, at constant β2 (panel (b)), with the indication of the intersection points between the objective functions curves and the levels red contour - OF(β)α and blue contour - OF(β)α fixed by the respective contours in (panel (a)).

for CO2−CH4, all of the other mixtures that belong to this class are characterized by a critical point which may achieve pressures higher than 100 bar and, for some of them (shown in italic font), measured isothermal phase diagrams are characterized by the nonexistence of a critical point (Figure 4c). In this respect, it is noted that, besides the indicated maximum pressure being relatively low, the phase diagram of the H2S−H2O system may potentially achieve very high pressures, since bubble and dew curves exhibit almost-infinite slope. Partitioning the 23 mixtures in two groups based on the maximum measured pressure (Pmax < 100 bar and Pmax > 100 bar), we can observe the experimental evidence of the theoretical

classification already reported in the previous section 4.1, between class I systems (typically low-pressure systems) and class II-Type III systems (typically high-pressure systems). For systems characterized by a continuous locus of critical points (i.e., class I systems), the problem that may arise is derived from the fact that the model easily reproduces these types of experimental data points. In practice, what may happen in this case is that the minimization of residuals too far below the experimental uncertainty leads to optimal parameters that perfectly reproduce the punctual value of fitted experimental data but completely fail in the calculation of properties of a different nature. The calculation of the meaningful minimal value 2266

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

uncertainty of parameters under optimization, as shown earlier in section 4.2.1. 4.2.3. The Uncertainty of Data. Despite the uncertainties of temperature and pressure measurements that are usually reported, experimentalists that measure VLE data rarely mention the uncertainties of each measured molar fraction. As a consequence, the application of the MLM, as well as the application of any method that requires experimental uncertainties, is usually carried out considering constant molar fraction uncertainties possibly declared by the authors.32 However, molar fraction uncertainties are not constant over an entire set of data but generally vary with the considered compositional measurement.33 Moreover, experimental molar fraction uncertainties are, generally, the ones that mostly contribute to the definition of molar fraction variances.33 For this reason, if the objective is to apply MLM, it is necessary to consider uncertainties specific to each molar fraction, rather than overall or average values. Starting from these considerations, MLM has been applied in this work, considering specific molar fraction uncertainties u(zi), estimated as shown in ref 33 from the knowledge of the experimental uncertainty, u(z•), associated with a specific molar fraction z•:

of the objective function (that is the minimum deviation between experiments and model that does not make sense to be further reduced) is dictated by experimental uncertainties. Moreover, when isothermal VLE data of the type represented in Figure 4ctypically low-temperature systems characterized by a discontinuous locus of critical points (i.e., class 2 systems) are not sufficiently numerous and too dispersed, it may happen that the confidence regions of optimal parameters is indefinitely extended or that, in the worst case, parameters under regression have a tendency to diverge. To better comprehend the meaning of this sentence, Figure 5a qualitatively shows a typical case of the unlimited confidence region (red contour). This outcome results from the combination of two effects: (i) the typical curvature of the OF(β) surface, in the space OF−β1−β2, which presents concave and convex regionssee the variation in convexity of the OF(β1)β2opt curve in Figure 5b; (ii) the high value of OF(β)α, that is related (fixing the number of parameters, the α-value) either to the low number of experimental data considered in the optimization process or to the optimal objective function OF(βopt), as explained in section 3.2 and by eq (A16) in Appendix II in the Supporting Information. More precisely, if the model is able to reproduce experimental data without introducing systematic errors, OF(βopt) is dependent only on the dispersion of data. However, this is a typical problem of phase diagrams where the critical region is not measured (or nor measurable) and for which, thus, only the low-pressure portion of the diagram is available. According to the above observations, to optimize models over these types of experimental phase diagrams, those that typically lead to projections of the OF(β) surface represented in Figure 5b (planes OF-β2 and OF-β1), the key point to obtain limited joint confidence regions of optimal parameters (see the blue contour in Figure 5a) is to reduce OF(β)α (see Figure 5b) by increasing the number of data points and reducing their dispersion. Finally, it is observed that the reduced number of experimental data, on the one side, coupled with their possibly ease of modeling (i.e., many combinations of parameters may lead to extremely low values of the objective function) but, on the other side, entail an extremely low sensibility of the objective function to the variation of parameters that may lead to the extreme divergence and, thus, an indefinite value of optimal parameters. More details and practical examples on the optimization of model parameters over the different types of systems will be given in section 5. 4.2.2. The Numerosity of Data. In systems for which NI ≈ NIT (see Table 3), the evaluation of the agreement between experimental data points is not straightforward since, on average, only one measured isothermal phase diagram for each temperature is available. For this reason, the optimization of models over isothermal VLE data is, generally, not recommended for these systems. Maximum values of NI/NIT ≈ 2−2.5 are associated with the systems N2−CH4, CH4−CO2, and N2−CO2; however, these values remain small. For some systems, available data points (NP) are highly limited. In particular, as observed in Table 3, the average number of data points available for each set of VLE data at similar temperatures (NP/NIT ratios) of each binary, it can be concluded that there are systems for which it is not possible to calibrate models over isothermal data, since those are too scarce (see, for example, mixtures CO2−H2S and H2S−H2O). The detected low number of experimental data, for some of the considered systems, sometimes associated with a challenging thermodynamic modeling, imposes the calculation of the

⎛ z2 − z ⎞ u zi = u z•⎜ •i 2 i • ⎟ ⎝ (z ) − z ⎠

(8)

Appendix III in the Supporting Information provides the u(z•) and z• values associated with each set of data considered in this paper.

5. THE OPTIMIZATION OF THE PR+EoS/aE,γ‑Wilson res MIXING RULES EQUATION OF STATE This section presents the results obtained from the optimization of the PR+EoS/aE,γ‑Wilson MR model to correlate data of systems res composed of CO2, N2, Ar, O2, CO, CH4, H2, H2O, and H2S, by Maximum Likelihood Method (MLM). As detailed previously, this optimization is performed over the VLE data of the binary mixtures obtained from the binary combination of these components. As mentioned in section 2, parameters A12 and A21 of the Wilson model are temperature-independent, from the assumption underlying such an activity coefficient model:34 “the average energy of an i−j interaction is independent of temperature, density, and other species present”. However, since these parameters are representative of intermolecular interaction forces 1−2/2−1, a priori temperature-dependent, it is more meaningful to optimize such parameters over isothermal VLE data and to compare, in a second moment, parameters (A12,A21) and their joint confidence regions obtained at different temperatures, in order to evaluate the coherence between the assumed temperature independence of these parameters and the actual temperature dependence of optimal results. 5.1. Selection of Modeling Options Relative to the Symmetry and Temperature-Dependency of Model Parameters. This section shows some preliminary considerations on the optimization of systems characterized by a different level of nonideality, respectively, low and high-pressure VLE systems. In particular, these considerations concern (i) the temperature dependence of parameters and the effect, on optimization results, of (ii) the number and dispersion (and availability of uncertainties) of experimental data and (iii) the selected optimization method (least-squared or MLM). The 2267

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

• the application of MLM on the entire sample of 32 (P, T, x, y) data points leads to a confidence region of optimal parameters, which is as extended as the one resulting from the use of the least-squares objective function (eq 6); however, the two regions do not overlap. Experimental uncertainties of (P, T, x, y) data points from refs 35 and 36 are reported in Appendix III in the Supporting Information. Despite these differences, it should be noted that, generally, all of the optimization results are compatible, since joint confidence regions of the optimal parameters overlap (also see Figure 7).

following observations are particularly focused on the optimization of parameters A12 and A21 of the considered model (see section 2) over VLE data of different types of more or less nonideal systems. 5.1.1. Modeling Options for Almost-Ideal Systems. This section considers the optimization over quasi-ideal systems (such as, in the case of binary combinations of components of typical CCS mixtures: N2−Ar, N2−O2, Ar−CH4, Ar−O2, CO−CH4, N2−CH4, N2−CO, CO2−H2S, ...) whose experimental VLE data, the applied model is particularly capable to reproduce. For these systems, it has been observed that the optimization of the two parameters A12 and A21 over isothermal VLE data leads to different results (optimal couples of parameters and confidence regions), depending on the characteristics mentioned above. This has been well shown in Figure 6, obtained for the

Figure 7. Vapor−liquid equilibrium (VLE) diagrams of the N2(1)− Ar(2) system at 100 K. Points denoted by asterisks (*) represent experimental points (P, T, x, y) by refs 35 and 36; dotted and continuous lines are calculated with the considered model and different parameter couples (A12, A21) that, with reference to points in Figure 6: continuous red, blue, black, and green lines respectively correspond to the following (A12, A21) points: red cross symbol (= black cross symbol); blue cross symbol; solid black circle; and solid red circle in Figure 6; dotted lines (- -) and ( ) are respectively calculated with parameter couples (A12 = 0, A21 = 0) and (A12 = 10, A21 = 10), which do not belong to any confidence region.

Figure 6. Optimal parameters A12, A21, and their 95% joint confidence regions, for the N2(1)−Ar(2) system at 100 K with different optimization methods. Legend of data points: (solid black circle, ●) optimal couple of parameters obtained with OF (eq 6) and 32 experimental points (P, T, x, y) by refs 35 and 36; (solid red circle, ●) optimal couple of parameters obtained with OF (eq 6) and 7 experimental points (P, T, x, y) by ref 36; (black cross symbol, ×) optimal couple of parameters obtained with OF (eq 6) × 106 and 32 experimental points (P, T, x, y) by refs 35 and 36; (red cross symbol, ×) optimal couple of parameters obtained with OF (eq 6) × 106 and 7 experimental points (P, T, x, y) by ref 36; (blue cross symbol, ×) optimal couple of parameters obtained with OF (eq 7) and 32 experimental points (P, T, x, y) by refs 35 and 36. Legend of lines: (dotted black line) 95% joint confidence region of parameters identified by black circles (●) and crosses (×), (dotted red line) 95% joint confidence region of parameters identified by red circles (●) and crosses (×), (continuous blue line) 95% joint confidence region of parameters identified by a blue cross symbol (×).

Moreover, all these contours denote a high degree of correlation between the two parameters. In fact, correlation factors c between parameters A12 and A21, calculated as explained in Appendix II in the Supporting Information (eq (A13)), are close to unitary values for these types of systems. Such a high degree of correlation also explains the reason why estimated optimal binary parameters differ with the applied optimization method. These features characterize systems where the error introduced by the model (and its parameters) is much lower than the one reflected by the inherent dispersion of experimental data. Similar conclusions can be drawn, in fact, for any similarly quasi-ideal low-pressure system. Another remarkable point is that, generally, all of the 95% joint-confidence regions of isothermal optimal parameters of these systems cross the line defined by the equality A12 = A21, thus highlighting the compatible application of this condition to all these systems, at all temperatures. The optimization of temperature-independent parameters, which the Wilson model is based on, does not degrade the results more than deviations that characterize experimental measurements. For this reason, Table 4 only reports the final, temperature-independent, symmetrical parameters A = A12 = A21 for this system N2−Ar and for the other similarly quasi-ideal mixtures N2−O2, Ar−CH4, Ar−O2, CO−CH4, N2−CH4, and CO2−H2S. This table also presents the resulting uncertainty of

N2−Ar system at 100 K by applying different optimization methods, as detailed in the caption of the figure. In more detail, the observation of Figure 6 leads to the following considerations: • the increased number of experimental data used to perform the optimization reduces the uncertainty region of the optimal parameters (see the red contour versus the black one), as anticipated in section 3.2; • the use of a more restrictive objective function, with respect to the least-squares objective function (eq 6) obtained by multiplying it by 106, increases the sensitivity of the optimal solution to the order of magnitude of the applied objective function (see the red/black circle versus the red/black cross); 2268

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

Table 4. Optimization Results Model Systems with the PR+EoS/aE,γ‑Wilson MR Equation of State, with a High Accuracy Level res AAD component 2

T [K]

NP

A12 [K]

xCO2

yCO2

ExCO

Ey

Ar

CH4

112.35−150.72

124

13.5

0.4

= A12

1.44

0.97

1.48

1.20

CO

CH4

123.40−178.00

38

14.6

0.3

= A12

0.66

0.53

1.24

1.40

N2

CH4

122.00−177.60

408

36

6

= A12

1.18

1.20

0.87

1.13

Ar

O2

90.00−110.00

45

7.4

0.2

= A12

2.97

2.76

1.50

5.16

N2

O2

100.11−122.70

42

4.7

0.2

= A12

0.52

0.50

0.08

0.36

N2

Ar

95.02−122.89

113

4.5

0.1

= A12

0.84

0.80

0.79

1.31

CO2 CO2 CO2 CO2 CO2 CO2 CO2

CH4 CH4 CH4 CH4 CH4 CH4 CH4

219.27, 219.90 230.00 240.00 250.00 253.00 270.00 288.60

15 57 12 28 5 58 8

112 130.5 100 80 74 49 −127

16 0.1 16 8 37 13 82

211 177.1 196 206 213 212 450

1.43 0.66 0.25 0.44 0.21 0.27 0.35

0.93 0.86 0.50 1.04 0.30 0.63 0.45

1.34 0.79 0.31 0.74 1.57 0.32 0.10

0.91 0.86 0.37 1.19 1.99 0.61 0.16

CO2

H2S

244.3−348.00

94

138

5

= A12

1.00

0.83

0.88

1.22

CO2 CO2 CO2 CO2

CO CO CO CO

223.0, 223.1 243.1 263.1 283.00, 283.10

11 10 9 6

−24 −70 −97 −164

14 15 19 58

239 319 388 535

16 20 38 134

−0.9778 −0.9798 −0.9478 −0.9636

0.55 0.38 0.15 0.11

0.90 1.40 0.65 0.65

0.65 0.51 0.38 0.24

0.89 0.89 0.88 0.61

CO2 CO2 CO2 CO2 CO2 CO2

Ar Ar Ar Ar Ar Ar

223.07 233.15, 233.32 238.08 253.04, 253.15, 253.28 273.12, 273.15, 273.26 288.23

13 17 17 35 32 8

0* 0* −46 −71 −183 −217

9 12 9 12 16 77

282 251 330 319 479 527

17 18 9 18 36 144

−0.8802 −0.9429 −0.9999 −0.9409 −0.9358 −0.9837

0.50 0.26 0.85 0.54 1.00 0.16

0.88 0.61 0.63 0.90 1.21 0.94

0.43 1.06 0.34 0.70 0.63 0.47

0.43 1.52 0.51 1.21 1.08 0.93

CO2 CO2 CO2 CO2 CO2 CO2 CO2 CO2 CO2

O2 O2 O2 O2 O2 O2 O2 O2 O2

218.04 223.15, 223.17 233.15 238.11 243.15 253.12, 253.15 263.15 273.00, 273.09, 273.15 283.15

8 24 12 13 11 27 10 36 6

43 6 22 −44 10 −93 −153 −193 −247

32 5 65 20 44 9 27 4 126

215 262 211 305 215 347 479 528 629

44 8 87 31 52 17 100 6 511

−0.9319 −0.9251 −0.9212 −0.9448 −0.9362 −0.9139 −0.8470 −0.9218 −0.8900

0.57 0.61 0.42 0.40 0.35 0.61 1.06 0.81 0.96

0.98 0.71 1.33 0.69 4.22 1.42 0.71 1.01 1.40

1.23 3.75 6.62 0.27 6.85 3.00 6.94 3.50 7.30

1.25 1.33 2.47 0.52 3.17 2.02 4.79 2.78 6.25

H2S H2S H2S H2S

CH4 CH4 CH4 CH4

222.2, 239.9, 252.00 277.6 310.9 344.3

14 37 37 17

186 163 49 −25

182 23 18 120

266 214 345 418

145 19 15 130

−0.9721 −0.9277 −0.9993 −0.9999

0.68 0.54 0.49 0.27

2.55 0.51 0.82 1.77

0.61 0.83 0.78 0.50

1.16 1.70 1.75 1.18

CO2 CO2 CO2 CO2 CO2 CO2 CO2 CO2 CO2

N2 N2 N2 N2 N2 N2 N2 N2 N2

218.04 220 223.09 238.06 240 250 253.05 270 288.09, 288.30

8 16 9 12 17 11 14 50 13

108 68 40 −13 −23 −2 −21 −110 −196

102 16 35 13 12 19 32 9 35

163 189 229 267 270 226 256 368 497

72 16 49 23 27 20 38 13 88

−0.9630 −0.9541 −0.9512 −0.9379 −0.9285 −0.9789 −0.9690 −0.9989 −0.9793

0.64 0.15 0.40 0.22 0.91 0.10 0.32 0.52 0.57

1.00 0.57 1.04 0.81 1.01 0.69 1.40 0.44 0.53

0.81 0.17 0.28 0.28 0.22 0.17 0.53 0.15 0.29

1.33 0.28 0.41 0.52 0.36 0.39 0.54 0.34 0.55

component 1

u A12 [K]

2269

A21 [K]

u A21 [K]

15 0.7 15 6 57 11 150

c

−0.9390 −0.9810 −0.9940 −0.9980 −0.9694 −0.9960 −0.9939

2

CO2

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research Table 4. continued AAD c

xCO2

yCO2

ExCO

Ey

0.1

−0.9389

0.08

0.47

0.07

0.39

92

8

−0.9821

0.62

1.24

0.44

0.87

12

89

10

−0.9711

0.50

1.08

0.91

1.69

115

9

74

3

−0.9443

1.34

1.45

0.95

1.37

380 335 358 412 547

36 57 21 34 112

586 488 373 195 −4

279 300 66 55 78

−0.8147 −0.9230 −0.9599 −0.9733 −0.9718

0.14 0.24 0.08 0.10 0.25

0.27 0.53 0.23 0.38 0.91

0.06 0.09 0.11 0.13 0.18

0.30 0.43 0.54 0.56 0.50

component 1

component 2

T [K]

NP

A12 [K]

H2O

H2S

298.15−377.40

65

1270.0

N2

H2

100.00−110.30

32

12

11

CO

H2

105, 115, 125

46

11

CH4

H2

100.00−183.10

318

N2 N2 N2 N2 N2

H2S H2S H2S H2S H2S

256.00 277.00 300.00 321.80 344.18

12 12 11 11 8

u A12 [K] 0.6

A21 [K] 580.0

u A21 [K]

2

CO2

Figure 8. Optimization results by MLM for the CO2(1)−CO(2) system at 223.1 K (red),37,38 243.1 K (blue), 263.1 K (green), and 283.1 K (black) (all data taken from ref 38): (a) joint confidence regions of optimal parameters (indicated by a solid circle, ●), and (b) isothermal phase diagrams calculated MR model and optimal parameters reported in panel (a) and experimental data (indicated by an asterisk, *). Information with the PR+EoS/aE,γ‑Wilson res about the considered experimental uncertainties are reported in Appendix III in the Supporting Information. In particular, specific (rather than constant) molar fraction uncertainties have been applied.

(c). Section 5.2 and Table 5 will provide more details on these relationships. This system is the only one, among the eight treated class I systems (see Table 3), that starts to show evidence of some of the features that characterize highly nonideal systems: asymmetrical and temperature-dependent features of binary interaction parameters. 5.1.2. Modeling Options for Highly Nonideal Systems. The evaluation of confidence regions is still more important whenever our objective is to compare parameters optimized over different thermodynamic conditions (for example, parameters optimized over isothermal VLE data of the same mixture but at different temperatures). In fact, a statistically meaningful comparison would require the comparison between confidence regions of parameters, rather than between optimal punctual values. To better appreciate such an observation, in Figure 8, we consider the CO2−CO system and the optimization of parameters A12 and A21 over isothermal VLE data at four different temperatures. The nonintersection of the well-defined 95%-elliptical regions at 223.1 and 263.1 K reveals that, at the considered 95% confidence level, optimal parameters are temperature-dependent. Moreover, note that the cause of the high extension and lesselliptical shape of the contour, relative to a temperature of 283.1 K, mainly lies in the limited number of experimental data around the critical region. In fact, as shown by Figure 9, this is the

these parameters. Since, for these systems, parameters are symmetrical, their calculated uncertainties are exact individual 95% confidence interval (see Appendix II in the Supporting Information). As a further consideration, the comparison between the confidence regions defined by A12 and A21 in Figure 6 and the confidence interval in Table 4 of parameter A = A12 = A21 leads to the observation that the reduction to two symmetrical parameters (instead of two asymmetrical ones), together with the use of a greater number of nonisothermal experimental data, results in an extremely lower uncertainty of parameters. Finally, it is observed that when the temperatures of the experimental data points increase so much that a critical point appears on each isothermal phase diagram, the assumptions of symmetry and temperature independence of parameters A12 and A21 may be incompatible with optimization results obtained considering asymmetrical and temperature-dependent parameters, such as for the CO2−CH4 system. Temperature-specific joint-confidence regions of parameters are not reported, for reasons of brevity; however the optimal parameters and their individual disjoint uncertainties, reported in Table 4, indicate that parameters A12 and A21 are asymmetrical and temperaturedependent. The linear temperature dependence, which characterizes each parameter, implies their strong correlation, which is also visible in Table 4, by means of the reported correlation factor 2270

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

Figure 9. Application of different parameter couples to the PR+EoS/aE,γ‑Wilson MR model (denoted by cross symbols in panel (a): purple, A; teal, B; and res brown, C), selected within the joint confidence region of parameters optimized at 283.1 K for the CO2−CO system (panel (a)), to calculate the respective isothermal phase diagrams (A, B, C in panel (b)) in the CO2−CO system at 283.1 K.

Figure 10. Optimization results by MLM for the H2S(1)−N2(2) system at 200.2 (red), 228.0 (blue), 256.0 (green), 277.0 (black), 300.0 (teal), 321.8 (brown), and 344.2 K (purple): (a) joint confidence regions of optimal parameters (indicated with a solid circle, ●); (b) isothermal phase diagrams MR model and parameters obtained from correlations with temperature of A12, A21 optimized at 256.0, 277.0, calculated with the PR+EoS/aE,γ‑Wilson res 300.0, and 321.8 K (optimal couples (A12, A21) reported in panel (a), correlations: A12 = −0.0401T2 + 17.335T − 1226, A21 = 0.0517T2 − 29.326T + 4497, with temperature T given in Kelvin); experimental data, collected from refs 39 and 40 are marked with an asterisk (*). For information about the considered experimental uncertainties, see Appendix III in the Supporting Information.

part of the phase diagram where the application of different couples of parameters, selected among those contained within the confidence region, has the greatest effect in the calculation of the critical region. In particular, the major axis of elliptical confidence regions is always headed in the direction imposed by the highest variability of parameters, which is mainly dictated by the inadequate amount of data available around the critical point. This is still more evident when our intent is to optimize the two asymmetrical parameters A12 and A21 for extremely highpressure systems, such as the type III mixtures (according to the classification of van Konynenburg and Scott; see section 4.1), for which the critical point may be absent or not measurable, since it tends to infinite values. This observation has been previously commented through the general case presented in section 4.2.1; in the following, it is presented in the practical case of model optimization on the type III H2S−N2 system. This binary mixture is characterized by a discontinuous locus of critical points, thus resulting in the impossibility of measuring the extreme highpressure vapor−liquid critical points in most of the temperature range (see experimental VLE data at 200.2−321.8 K in Figure 10b). For this reason, the use of VLE data of these systems to optimize the considered two-parameter model might lead to undefined nonelliptical confidence regions (i.e., ill-conditioned optimization results), such as the one represented in Figure 10a,

relative to the optimization over VLE data at 228 K. Furthermore, Figure 10a also shows joint confidence regions relative to the optimization of temperature-specific parameter couples over other isothermal phase diagrams at 256−344 K. Considering Figure 10b as a complement to Figure 10a, the effect of the different distributions of VLE points within each phase diagram on the uncertainty region of optimal parameters, from the lowest to the highest physically achievable saturation pressure (the critical one, in this case), becomes evident. In fact, as observed for the CO2−CO system at 283.1 K, extended confidence regions may result from the combination of the limited number of experimental data around the critical point and the difficulty in modeling this region (see confidence region of parameters optimized over data at 344 K, in Figure 10a), or from the concomitant absence and dispersion of data and the ease in representing the, only available, low-pressure saturation data (see confidence region of parameters optimized over data at 228 K, in Figure 10a, and the discussion in section 4.2.1). Confidence regions relative to parameters optimized over VLE diagrams at intermediate temperatures (mostly, 300 and 322 K) are well-defined, since, on the one side, the characteristics of available data points are appropriate (in terms of number, dispersion, and distribution) for the determination of finite values of optimal parameters in the considered pressure range and, 2271

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

Figure 11. MLM optimization results of the model PR+EoS/aE,γ‑Wilson MR for the CO2(1)−N2(2) system at 220 (red), 240 (blue), and 270 K (purple): res (a) joint confidence regions of optimal parameters (indicated with a solid circle, ●); (b) isothermal phase diagrams calculated with the PR+EoS/ MR model and parameters obtained from correlations with temperature of A12,A21 optimized at 220, 240, and 270 K (optimal couples (A12,A21) aE,γ‑Wilson res reported in panel (a), correlations: A12 = −3.4768T + 825, A21 = 3.5332T − 585, with temperature T given in Kelvin); experimental data, collected from refs 41−45 are marked with an asterisk (*). For information about the considered experimental uncertainties, see Appendix III in the Supporting Information.

the 18 systems that the model is able to accurately represent. In particular, for these systems, the calculated percent average absolute deviation (AAD) values,

on the other side, the thermodynamic model is able to well reproduce experimental data points. However, it is worth remarking that when parameters are not optimized over almost critical experimental data, their confidence regions, if determined as shown in this work (Appendix II in the Supporting Information), do not incorporate the uncertain modeling of the critical region. For this reason, the extrapolation capability of models optimized over low pressure and limited data should be carefully validated. In the nonideal case-studies CO2−CO and N2−H2S treated above, it has been shown that the optimization of parameters over an increased number of experimental data, for phase diagrams with an equal distribution of saturation points, results in reduced confidence regions (see Figures 8b and 10b). However, confidence regions of parameters optimized over phase equilibrium data more concentrated in the region of higher nonideality (the critical one) are more extended. For this reason, it may happen that confidence regions of parameters optimized over few and low-pressure experimental data present a similar extension as the one of parameters optimized over numerous data that also cover critical conditions (see Figure 11). It is thus evident that the information content of the joined confidence region of parameters optimized at 220 K is lower than that of confidence regions at 240 and 270 K. Therefore, one should always consider the fact that the uncertainty of parameters refers to the thermodynamic conditions under which data have been optimized. In conclusion, within the characteristics of numerosity, dispersion, and distribution of experimental data, the optimization of the considered two-parameter PR+EoS/aE,γ‑Wilson MR res model over few and dispersed data or data that do not represent the actual nonideal behavior of the considered system may lead to either unphysical optimal parameters or detrimental undefined results. A final remark, already highlighted, is devoted to the importance of relating each measurement to its specific uncertainty, when attempting at applying MLM. 5.2. Optimized Parameters for PR+EoS/aE,γ‑Wilson MR res Equation of State. According to the considerations reported in the previous section, the model has been optimized over the 23 binary systems for which experimental data are available. In particular, Table 4 and Appendix IV in the Supporting Information report numerical and qualitative results relative to

AADzCO = 2

100 NP

NP calc exp (P , Ti , A opt) − z CO | ∑ |zCO 2 i 2, i

with zi = xi or zi

i

are lower or only acceptably higher (i.e., considering that standard deviations of experimental data are mostly calculated from assumed experimental uncertainties) than the percentage error calculated from experimental variances, EzCO = 2

100 NP

NP

∑ σz

CO2, i

with zi = xi or zi

i

where the two deviations have been calculated for both the vapor and the liquid phase, considering the experimental NP data exp points, zCO , and the calculations zcalc CO2(Pi,Ti,Aopt) performed with 2,i the model, as a function of thermodynamic property i (Pi,Ti) and optimal parameters, Aopt, zcalc CO2(Pi, Ti, Aopt). Standard deviations σzCO , i are calculated from eqs (A7) and (A8) reported in 2 Appendix II in the Supporting Information. The same table reports the individual uncertainty of parameters, uA12 and uA12, and the correlation factor c, described in section 3.2 and calculated as reported in Appendix II in the Supporting Information. For these systems, the following has been reported: • the temperature independence and symmetry (A12 = A21) of parameters for seven almost-ideal systems that the model accurately reproduces (N2−O2, Ar−CH4, Ar−O2, CO−CH4, N2−CH4, CO2−H2S, N2−Ar) • the temperature dependence and asymmetry of parameters for seven nonideal systems that the experimental VLE data the model accurately reproduces (CO2−CH4, CO2−CO, CO2−Ar, CO2−O2, H2S-CH4, CO2−N2, H2S−N2) • that considering temperature-specific and asymmetric binary parameters for four highly nonideal systems (N2−H2, CO−H2, CH4−H2, H2O−H2S) does not improve modeling results more than considering temperature-independent parameters and, particularly for the H2O−H2S mixture (which isothermal phase diagrams are all of the type (c) reported in Figure 4) the optimization of temperature-specific parameters leads to 2272

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

Table 5. Optimal Correlations Aij = aijT2 + bijT + cij of the Temperature-Specific Parameters Reported in Table 4, for PR+ EoS/aE,γ‑Wilson MR EoS res component 1

component 2

NT

a12 [1/K]

b12

c12 [K]

a21 [1/K]

b21

c21 [K]

CO2 CO2 CO2 CO2 H2S CO2 N2

CH4 CO Ar O2 CH4 N2 H2S

7 4 6 9 4 9 5

0 0 0 0 0 0 0.0517

−2.4551 −2.1977 −3.3775 −3.9845 −1.8857 −3.4881 −29.3260

695 466 758 895 635 833 4497

0 0 0 0 0 0 −0.0401

1.4405 3.9843 3.1465 5.3207 1.5425 3.6760 17.3350

−154 −650 −420 −925 −134 −626 −1226

Table 6. Molecular Parameters of H2, N2, CO, Ar, and CH4 Molecular Size Parameter, σ

Molecular Stickiness, εii

Polarizability, α

Dipole Moment, μ

Quadrupole Moment

species

σ [Å]

ref

εii/kB [K]

ref

α [Å ]

ref

μ [D]

ref

present?

ref

H2 N2 CO Ar CH4

2.9860 3.3130 3.2507 3.4784 3.7039

47 51 51 51 51

19.2775 90.96 92.15 122.23 150.03

47 51 51 51 51

0.787 1.710 1.953 1.664 2.448

48 48 48 48 48

0.000 0.000 0.112 0.000 0.000

49 49 49 49 49

yes yes yes no no

50 50 50 50 50

been found, as follows, comparing this system to other similar systems containing hydrogen for which parameters have been, more easily, optimized (i.e., CO−H2, N2−H2, and CH4−H2). Recalling that binary interaction parameters are directly proportional to the differences between molecular local interaction energies of either different or equal components, AWilson ∝ εji − εii (see section 2), the estimation of AWilson for the ij ij mixture Ar−H2,

undefined and nonelliptical confidence regions that make the results of the optimization unreliable. For these systems, temperature-independent and asymmetric binary parameters have been thus optimized. It is specified that the uncertainties of the binary parameters reported in Table 5 only enable the correct evaluation of the uncertainty of model calculations of VLE properties of binary systems. However, validation of the model calculations over multicomponent systems has been recently performed.46 Based on the obtained temperature-specific parameters and their uncertainty, relationships (i.e., parameters aij, bij, cij) that relate the parameters to the temperature,

Wilson Wilson AAr − H 2 ∝ εH 2 − Ar − εAr − Ar and A H 2 − Ar ∝ εAr − H 2 − εH 2 − H 2

is performed over comparisons between local interaction energies εii and εij of the different systems upon which their Aij parameters are dependent:

A12 = a12T 2 + b12T + c12

Wilson AiWilson − H 2 ∝ εH 2 − i − εii and A H 2 − i ∝ εi − H 2 − εH 2 − H 2

A 21 = a 21T 2 + b21T + c 21

with i = N2 , CO, CH4

have been optimized in a second step of the optimization process, by minimizing: NT

∑ i=1

3

First, it is observed that argon is characterized by a volatility that lies between that of H2/N2/CO and CH4, with the order being dictated by the magnitude of the parameter εii (also called “molecular stickiness”), which is representative of the molecular attractive strength34 (see Table 6). Observing the analogy between the parameters Ai−H2 and AH2−i optimized for each of the N2−H2, CO−H2, CH4−H2 systems (reported in Table 4) and the intermolecular forces εii (see Figure 12), parameters AAr−H2 and AH2−Ar have been roughly estimated from the linear correlation Ai−H2(εii) and AH2−i(εii) reported in the caption of Figure 12, obtaining

* AijT,opt

− Aijcalc(T *) * u(AijT,opt )

Parameters aij and bij are reported in Table 5, for each mixture treated in Table 4 and characterized by temperature-specific parameters. With regard to the other five systems not mentioned above (N2−CO, H2S−CO, CO2−H2O, CO2−H2, and Ar−H2): (1) for the N2−CO system, it was considered more reliable to apply a zero value for the parameters than the value optimized from experimental data found in the literature (see Appendix III in the Supporting Information); (2) for the other four systems, the convergence of parameters to isothermal optimal values was particularly critical, since temperature-specific parameters have a tendency to diverge. Thus, for the H2S−CO, CO2−H2O, and CO2−H2 mixtures, after observing the functional form that relates parameters to temperature, a temperature-dependent relationship has been optimized, as reported in Table 7. In contrast, a reasonable (nonoptimal) value for parameters of the mixture Ar−H2 has

⎛ ⎞ Wilson εAr−Ar AAr = 122.23 K⎟ ≈ 66 K − H2⎜ ⎝ kB ⎠

and ⎛ εAr−Ar ⎞ A HWilson = 122.23 K⎟ ≈ 82 K ⎜ 2 − Ar ⎝ kB ⎠

Appendix IV in the Supporting Information shows phase diagrams obtained with the application of the model optimized for the 23 treated systems, the parameters of which are reported 2273

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

symmetricity and temperature dependence. In case of observed temperature dependences, the polynomial relationships between parameters and temperature were then determined. For the last five systems (N2−CO, H2S−CO, CO2−H2O, CO2−H2, and Ar−H2), the optimization process was particularly critical, either because experimental data were not adequate or the system was extremely nonideal. For these systems, the analysis of the confidence regions has enabled to confirm the illconditioning of such optimization problems. However, also for these systems, reliable parameters were provided. In this work, the reader has been urged to focus on the need to account for the uncertainty of optimized parameters. In fact, the extension of uncertainty of optimal parameters may be dependent either on systematic inaccuracies introduced by the thermodynamic model or characteristics of experimental data (numerosity, dispersion, distribution, etc.). In the paper, it is shown that, to reduce these uncertainties, it is necessary, above all, to quantify them and study the causes of their variability.

Figure 12. Parameters Ai−H2 (denoted by a cross symbol, ×) and AH2−i (denoted by an open circle, ○) reported in Table 4, as a function of εii/ kB, for i = N2 (red points), CO (black points), and CH4 (blue points). Dotted lines interpolate both Ai−H2 (- - -) and AH2−i (···), as a function of εii. Linear interpolations are defined as Ai−H2 [K] = 1.770(εii/kB [K]) − 150.5 and AH2−i [K] = −0.283(εii/kB [K]) + 116.4.



in Table 4 (for mixtures N2−O2, Ar−CH4, Ar−O2, CO−CH4, N2−CH4, CO2−H2S, N2−Ar, N2−H2, CO−H2, CH4−H2, H2O−H2S), Table 5 (for CO2−CH4, CO2−CO, CO2−Ar, CO2−O2, H2S−CH4, CO2−N2, N2−H2S), and Table 7 (for H2S−CO, CO2−H2O, CO2−H2, Ar−H2, and N2−CO).

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.iecr.6b04190. Principles and general application of the Maximum Likelihood Method (Appendix I); optimal parameters, including their correlation and uncertainties (Appendix II); source and characteristics of experimental data (Appendix III), including bibliographic references of experimental data (Appendix AIII, section 1) and uncertainties of experimental data (Appendix AIII, section 2); phase diagrams of binary mixtures by optimized model (Appendix IV) (PDF)

6. CONCLUSIONS In response to a previous publication by the same authors,6 which showed the improved modeling capability of a novel advanced thermodynamic model (the Peng−Robinson + EoS/aE,γ res mixing rules) in representing VLE properties of the highly nonideal mixtures CO2−N2, CO2−Ar, and CO2−O2, this paper has been devoted to the optimization of the same model, over VLE experimental data of binary systems resulting from the binary combination of CO2, H2, N2, O2, Ar, CO, CH4, H2S, H2O. These 36 binary systems have been first classified on the basis of the criteria proposed by van Konyenburg and Scott; this classification provided an idea of the difficulty in modeling the phase behavior of these mixtures. Among these 36 systems, VLE experimental data have been found in the literature for the following 23 mixtures: N2−O2, Ar−CH4, Ar−O2, CO−CH4, N2−CH4, CO2−H2S, N2−Ar, CO2−CH4, CO2−CO, CO2−Ar, CO2−O2, H2S−CH4, CO2−N2, H2S−N2, N2−H2, CO−H2, CH4−H2, H2O−H2S, N2−CO, H2S−CO, CO2−H2O, CO2−H2, and Ar−H2. The optimization of the binary parameters of this model (A12, A21), for each of the 23 mixtures, has been performed by applying the Maximum Likelihood Method (MLM), evaluating both the correlation and joint uncertainty of optimal parameter couples. Based on these results, it has been reported that, among these 23 mixtures, the optimized model is able to accurately reproduce experimental data (i.e., within their experimental uncertainty) of the first 18 systems listed above. For each of them, the analysis of the uncertainty of parameters has enabled the evaluation of their



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Romain Privat: 0000-0001-6174-9160 Jean-Noël Jaubert: 0000-0001-7831-5684 Notes

The authors declare no competing financial interest.



NOMENCLATURE AND ACRONYMS

Generic Symbols

a = molar Helmholtz energy or energy parameter in cubic EoS b = co-volume parameter in cubic EoS g = molar Gibbs energy n = number of moles z = generic molar fraction x = liquid-phase molar fraction y = vapor-phase molar fraction P = pressure

Table 7. Estimated Correlations Aij = aijT2 + bijT + cij of Parameters for PR+EoS/aE,γ‑Wilson MR EoS res component 1

component 2

a12 [1/K]

b12

c12 [K]

a21 [1/K]

b21

c21 [K]

H2S CO2 CO2 Ar N2

CO H2O H2 H2 CO

0 −0.01002 0 0 0

−4.6205 9.0015 −4.6793 0 0

1535 −1126 1558 66 0

0.02990 0 0 0 0

−14.4720 3.5113 −0.7579 0 0

2006 −745 385 82 0

2274

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

(6) Lasala, S.; Chiesa, P.; Privat, R.; Jaubert, J.-N. VLE Properties of CO2-Based Binary Systems Containing N2, O2 and Ar: Experimental Measurements and Modelling Results with Advanced Cubic Equations of State. Fluid Phase Equilib. 2016, 428, 18−31. (7) Peng, D.-Y.; Robinson, D. B. A New Two-Constant Equation of State. Ind. Eng. Chem. Fundam. 1976, 15 (1), 59−64. (8) Robinson, D. B.; Peng, D.-Y. The Characterization of the Heptanes and Heavier Fractions for the GPA Peng−Robinson Programs; Gas Processors Association: Tulsa, OK, 1978. (9) Lorentz, H. A. Ueber Die Anwendung Des Satzes Vom Virial in Der Kinetischen Theorie Der Gase. Ann. Phys. 1881, 248 (1), 127−136. (10) Jaubert, J.-N.; Mutelet, F. VLE Predictions with the PengRobinson Equation of State and Temperature Dependent kij Calculated through a Group Contribution Method. Fluid Phase Equilib. 2004, 224 (2), 285−304. (11) Jaubert, J. N.; Privat, R.; Mutelet, F. Predicting the Phase Equilibria of Synthetic Petroleum Fluids with the PPR78 Approach. AIChE J. 2010, 56 (12), 3225−3235. (12) Vitu, S.; Privat, R.; Jaubert, J.-N.; Mutelet, F. Predicting the Phase Equilibria of CO2 + Hydrocarbon Systems with the PPR78 Model (PR EOS and kij Calculated through a Group Contribution Method). J. Supercrit. Fluids 2008, 45 (1), 1−26. (13) Lasala, S. Advanced Cubic Equations of State for Accurate Modelling of Fluid Mixtures. Application to CO2 Capture Systems, Ph.D. Thesis; Politecnico di Milano: Milano, Italy, 2016. (14) Le Guennec, Y.; Lasala, S.; Privat, R.; Jaubert, J.-N. A Consistency Test for α-Functions of Cubic Equations of State. Fluid Phase Equilib. 2016, 427, 513−538. (15) Twu, C. H.; Bluck, D.; Cunningham, J. R.; Coon, J. E. A Cubic Equation of State with a New Alpha Function and a New Mixing Rule. Fluid Phase Equilib. 1991, 69, 33−50. (16) Le Guennec, Y.; Privat, R.; Jaubert, J.-N. Development of the Translated-Consistent Tc-PR and Tc-RK Cubic Equations of State for a Safe and Accurate Prediction of Volumetric, Energetic and Saturation Properties of Pure Compounds in the Sub- and Super-Critical Domains. Fluid Phase Equilib. 2016, 429, 301−312. (17) Voutsas, E.; Magoulas, K.; Tassios, D. Universal Mixing Rule for Cubic Equations of State Applicable to Symmetric and Asymmetric Systems: Results with the Peng−Robinson Equation of State. Ind. Eng. Chem. Res. 2004, 43 (19), 6238−6246. (18) Michelsen, M. L. A Modified Huron-Vidal Mixing Rule for Cubic Equations of State. Fluid Phase Equilib. 1990, 60 (1−2), 213−219. (19) Van der Waals, J. D. Van der Waals and the Physics of Liquids; Rowlinson, J. S., Ed.; On the continuity of the gaseous and liquid states, Part I; Elsevier: Amsterdam, New York, 2004. (20) Neau, E.; Peneloux, A. Estimation of Model Parameters. Comparison of Methods Based on the Maximum Likelihood Principle. Fluid Phase Equilib. 1981, 6, 1−19. (21) Aldrich, J. R. A. Fisher and the Making of Maximum Likelihood 1912−1922. Stat. Sci. 1997, 12 (3), 162−176. (22) Péneloux, A.; Deyrieux, R.; Neau, E. Réduction Des Données Expérimentales et Théorie de L’information. J. Chim. Phys. 1975, 72 (10), 1101−1107. (23) Péneloux, A.; Deyrieux, R.; Canals, E.; Neau, E. The Maximum Likelihood Test and the Estimation of Experimental Inaccuracies. Application to Data Reduction for Liquid−Vapour Equilibrium. J. Chem. Phys. 1976, 73 (7−8), 707−716. (24) Draper, N. R.; Smith, H. Applied Regression Analysis, 3rd Edition; Wiley: New York, 1998 (DOI: 10.1002/9781118625590). (25) Konynenburg, P. H. V.; Scott, R. L. Critical Lines and Phase Equilibria in Binary Van Der Waals Mixtures. Philos. Trans. R. Soc., A 1980, 298 (1442), 495−540. (26) Privat, R.; Jaubert, J. N. Classification of Global Fluid-Phase Equilibrium Behaviors in Binary Systems. Chem. Eng. Res. Des. 2013, 91 (10), 1807−1839. (27) Fornari, R. E.; Alessi, P.; Kikic, I. High Pressure Fluid Phase Equilibria: Experimental Methods and Systems Investigated (1978− 1987). Fluid Phase Equilib. 1990, 57 (1−2), 1−33.

T = temperature V = volume Greek Symbols

α = α-function of the EoS ε = molecular energy interaction parameter σ = standard deviation γ = activity coefficient ω = acentric factor ρ = density Φ = volume fraction Ω = local composition parameter μ = experimental uncertainty β = generic parameter under optimization Bold Characters

bold font = indicates a vector Superscripts

E = excess property γ = activity coefficient model exp = experimental quantity calc = calculated quantity Subscripts

m = mixture property c = property at critical conditions res = residual term of the property comb = combinatorial term of the property x or y = characteristic relative to either the liquid or vapor phase i = characteristics relative to pure single component i ij = characteristics relative to the component couple ij Constants and Parameters

R = universal gas constant ΛEoS = parameter of EoS/aE,γ mixing rules r1, r2 = characteristic parameters of cubic equations of state Aij = binary interaction parameter of the EoS/aE,γ‑Wilson model res Acronyms

CCS = CO2 capture and storage EOR = enhanced oil recovery EoS = equation of state MHV1 = modified Huron−Vidal model MR = mixing rule NC = number of components PR = Peng−Robinson model VdW = van der Waals model VLE = vapor−liquid equilibrium



REFERENCES

(1) Effects of Impurities on Geological Storage of CO2; IEAGHG: Cheltenham, U.K., 2011. (2) Li, H.; Yan, J. Evaluating Cubic Equations of State for Calculation of Vapor−Liquid Equilibrium of CO2 and CO2-Mixtures for CO2 Capture and Storage Processes. Appl. Energy 2009, 86 (6), 826−836. (3) Lasala, S.; Chiesa, P.; Di Bona, D.; Consonni, S. Vapour−Liquid Equilibrium Measurements of CO2 Based Mixtures: Experimental Apparatus and Testing Procedures. Energy Procedia 2014, 45, 1215− 1224. (4) Xu, X.; Lasala, S.; Privat, R.; Jaubert, J.-N. E-PPR78: A Proper Cubic EoS for Modeling Fluids Involved in the Design and Operation of Carbon Dioxide Capture and Storage (CCS) Processes. Int. J. Greenhouse Gas Control 2017, 56, 126−154. (5) Jaubert, J.-N.; Avaullee, L.; Pierre, C. Is It Still Necessary to Measure the Minimum Miscibility Pressure? Ind. Eng. Chem. Res. 2002, 41 (2), 303−310. 2275

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276

Article

Industrial & Engineering Chemistry Research

(50) Computational Chemistry Comparison and Benchmark DataBase Standard Reference Database 101; National Institute of Standards and Technology (NIST): Gaithersburg, MD, 2016. (51) Gross, J.; Sadowski, G. Perturbed-Chain SAFT: An Equation of State Based on a Perturbation Theory for Chain Molecules. Ind. Eng. Chem. Res. 2001, 40 (4), 1244−1260.

(28) Dohrn, R.; Brunner, G. High Pressure Fluid-Phase Equilibria: Experimental Methods and Systems Investigated. Fluid Phase Equilib. 1995, 106, 213−282. (29) Christov, M.; Dohrn, R. High-Pressure Fluid Phase Equilibria: Experimental Methods and Systems Investigated (1994−1999). Fluid Phase Equilib. 2002, 202 (1), 153−218. (30) Dohrn, R.; Peper, S.; Fonseca, J. M. S. High-Pressure Fluid-Phase Equilibria: Experimental Methods and Systems Investigated (2000− 2004). Fluid Phase Equilib. 2010, 288, 1−54. (31) Fonseca, J. M. S.; Dohrn, R.; Peper, S. High-Pressure Fluid-Phase Equilibria: Experimental Methods and Systems Investigated (2005− 2008). Fluid Phase Equilib. 2011, 300, 1−69. (32) Sutton, T. L.; MacGregor, J. F. The Analysis and Design of Binary Vapour−Liquid Equilibrium Experiments. Part I: Parameter estimation and consistency tests. Can. J. Chem. Eng. 1977, 55, 602−608. (33) Lasala, S.; Chiesa, P.; Privat, R.; Jaubert, J.-N. Optimizing Thermodynamic Models: The Relevance of Molar Fraction Uncertainties. J. Chem. Eng. Data 2017, 62, 825−832. (34) Elliott, J. R.; Lira, C. T. Introductory Chemical Engineering Thermodynamics; Prentice Hall: Upper Saddle River, NJ, 2012. (35) Narinskii, G. B. Liquid−Vapour Equilibrium in the Argon− Nitrogen System. I. Experimental Data and Their Verification. Russ. J. Phys. Chem. 1966, 40 (9), 1093−1096. (36) Elshayal, I. M.; Lu, B. C.-Y. Measurement of Total Pressures for Ethylene−Propane Mixtures. Can. J. Chem. Eng. 1975, 53 (1), 83−87. (37) Kaminishi, G.; Toriumi, T. Vapor−Liquid Equilibria in the Systems: CO2−CO, CO2−CO−H2 and CO2−CH4. Rev. Phys. Chem. Jpn. 1968, 38 (1), 79−84. (38) Christiansen, L. J.; Fredenslund, A.; Gardner, N. Gas−Liquid Equilibria of the CO2−CO and CO2−CH4−CO Systems. In Advances in Cryogenic Engineering 1995, 309−319. (39) Kalra, H.; Krishnana, T. R.; Robinson, D. B. Equilibrium-Phase Properties of Carbon Dioxide−Butane and Nitrogen−Hydrogen Sulfide Systems at Subambient Temperatures. J. Chem. Eng. Data 1976, 21 (2), 222−225. (40) Besserer, G. J.; Robinson, D. B. Equilibrium-Phase Properties of Nitrogen−Hydrogen Sulfide System. J. Chem. Eng. Data 1975, 20 (2), 157−161. (41) Al-Sahhaf, T. A.; Kidnay, A. J.; Sloan, E. D. Liquid + Vapor Equilibriums in the N2 + CO2 + CH4 System. Ind. Eng. Chem. Fundam. 1983, 22 (4), 372−380. (42) Brown, T. S.; Sloan, E. D.; Kidnay, A. J. Vapor−Liquid Equilibria in the Nitrogen + Carbon Dioxide + Ethane System. Fluid Phase Equilib. 1989, 51, 299−313. (43) Brown, T. S.; Niesen, V. G.; Sloan, E. D.; Kidnay, A. J. Vapor− Liquid Equilibria for the Binary Systems of Nitrogen, Carbon Dioxide, and N-Butane at Temperatures from 220 to 344 K. Fluid Phase Equilib. 1989, 53 (C), 7−14. (44) Somait, F. A.; Kidnay, A. J. Liquid−Vapor Equilibriums at 270.00 K for Systems Containing Nitrogen, Methane, and Carbon Dioxide. J. Chem. Eng. Data 1978, 23 (4), 301−305. (45) Yucelen, B.; Kidnay, A. J. Vapor−Liquid Equilibria in the Nitrogen + Carbon Dioxide + Propane System from 240 to 330 K at Pressures to 15 MPa. J. Chem. Eng. Data 1999, 44, 926−931. (46) Lasala, S.; Chiesa, P.; Privat, R.; Jaubert, J.-N. Measurement and Prediction of Multi-Property Data of CO2−N2−O2−CH4 Mixtures with the “Peng-Robinson + Residual Helmholtz Energy-Based” model. Fluid Phase Equilib. 2017, 437, 166−180. (47) Ghosh, A.; Chapman, W. G.; French, R. N. Gas Solubility in HydrocarbonsA SAFT-Based Approach. Fluid Phase Equilib. 2003, 209 (2), 229−243. (48) Olney, T. N.; Cann, N. M. M.; Cooper, G.; Brion, C. E. E. Absolute Scale Determination for Photoabsorption Spectra and the Calculation of Molecular Properties Using Dipole Sum-Rules. Chem. Phys. 1997, 223 (1), 59−98. (49) Nelson, R. D., Jr.; Lide, D. R., Jr.; Maryott, A. A. Selected Values of Electric Dipole Moments for Molecules in the Gas Phase; National Standard Reference Data Series, NBS 10; U.S. Department of Commerce: Washington, DC, 1967. 2276

DOI: 10.1021/acs.iecr.6b04190 Ind. Eng. Chem. Res. 2017, 56, 2259−2276