Subscriber access provided by Washington University | Libraries
Fossil Fuels
Compositional modeling of crude oils using C10C36 properties generated by molecular simulation Philippe Ungerer, Marianna Yiannourakou, Alexander Mavromaras, and Julien Collell Energy Fuels, Just Accepted Manuscript • DOI: 10.1021/acs.energyfuels.8b04403 • Publication Date (Web): 04 Mar 2019 Downloaded from http://pubs.acs.org on March 5, 2019
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 48 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
Compositional modeling of crude oils using C10-C36 properties generated by molecular simulation Philippe Ungerer*&† , Marianna Yiannourakou&, Alexander Mavromaras& , Julien Collell§ &Materials Design, S.A.R.L., 42 Avenue Verdier, 92120 Montrouge, France §Total SA, CSTJF, Avenue Larribau, 64 Pau, France *author to whom correspondence should be addressed.
[email protected] KEYWORDS: reservoir fluid characterization, stock tank oil density, petroleum heavy fractions, n-alkanes, iso-alkanes, naphthenes, aromatics, polyaromatic, naphthenoaromatic compounds, thiophenic compounds, ideal gas heat capacity, ideal gas enthalpy, saturation pressure, saturated liquid density, Monte Carlo method, thermodynamic integration.
ABSTRACT.
Due to the lack of detailed analysis in the C10+ fraction and scarcity of reliable thermodynamic properties on polycyclic compounds, it is usually not feasible to relate crude oil properties with the chemical structure of heavy fractions. Over the last decades the description of C10-C36 fractions has mostly relied on average Cn properties determined from observations. We propose an
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
alternative approach in two major steps. In a first step we use Monte Carlo simulation methods to generate Vapor-Liquid Equilibrium (VLE) data on representative hydrocarbons between C10 and C30, from ambient to near-critical temperature. Based on these results, standard liquid density and saturation pressure are correlated for naphthenic hydrocarbons (mono- and poly-cyclic), aromatic hydrocarbons (mono-cyclic, poly-cyclic, naphthenoaromatic), and thiophenic compounds up to C36. In a second step we apply the predicted properties on C10-C36 families to model nine (9) real crude oils. The Cn fractions (n =10 - 36) are described with an exponential distribution, and the concentrations of n-iso/ naphthenes/ aromatics/ NSO compounds are modelled explicitly. Using crude-specific information (e.g. C1-C10 analysis) and general statistics about reservoir fluids (e.g. target region in ternary diagrams), we obtain an excellent agreement of crude oil density, average molecular weights of Cn+ fractions, and SARA analysis (when available). The predicted standard liquid density of Cn fractions increases with carbon number, as observed on real fluids. The higher Cn density observed in aromatic crudes is also well predicted. These results suggest that additional properties (e.g. VLE of live oils) may be predicted with more insight by applying the proposed simulation-based approach in future studies.
ACS Paragon Plus Environment
Page 2 of 48
Page 3 of 48 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
1. INTRODUCTION The composition of reservoir fluids is used to predict pressure-volume behavior, phase equilibria with high pressure gases, transport properties, and thermal properties. As far as volumetric and phase equilibrium properties are concerned, property prediction is generally based on thermodynamic models such as equations of state, which present the advantage of fast computing time for the numerous phase equilibrium computations involved in three-dimensional reservoir studies [1]. When equations of state are not fast enough, simplifications are introduced such as reducing the number of model components by lumping with minimal loss of accuracy[2,3]. However, the reliability of property predictions is limited by two main factors : the difficulty of analyzing the complex composition of crude oils at the molecular level and the limited knowledge of volumetric and thermodynamic properties for all possible compounds. The molecular level analysis of reservoir fluids by gas chromatography (GC) is most often limited to the C1-C10 fraction [1]. Geochemical studies using mass spectrometry have shown that crude oils contain a wide range of polycyclic saturated and aromatic compounds with 2-6 cycles in the C15+ fraction of crude oils [4,5]. However the standard analysis of crude oils by GC alone can just provide a semi-quantitative analysis between C10 and C15, and qualitative composition of each boiling point fraction between C15 and C20 [6]. High temperature gas chromatography (HTGC) provides the overall distribution up to C50 but without indication of chemical families. Finally, the characterization of heavier fractions – resins and asphaltenes – at a molecular level, still presents many challenges [7-10]. As a consequence, the full composition of the C10+ fractions is approximate in most cases, so that empirical methods are used to model the heavy fraction of reservoir fluids. For instance, the percentage of paraffins, naphthenes and aromatics (PNA) in crude oil fractions may be estimated from refractive index, viscosity and specific gravity [11].
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Exponential or log-normal distribution functions with carbon number are used to match experimental data such as the average molecular weight or liquid density of crude oils [12,13] or the liquid drop-out curve in the case of condensate gases [14]. Whatever the detail of analysis, property prediction accuracy is severely limited by the poor knowledge of volumetric and thermodynamic properties of its polycyclic components which have been less extensively studied than n-alkanes [15,16]. Standard properties (normal boiling temperature and saturated liquid density) have been measured or correlated for heavy cyclic compounds [17]. However, EOS parametrization for such compounds would require more data. Indeed, saturation pressure of pure compounds versus temperature is the main information allowing the accurate parametrization of EOS for phase equilibrium computations [18]. This explains why the common practice in reservoir engineering is still to use average pseudocomponents to approximate the Cn fractions [19,20] although these pseudo-components are considered by experts as insufficiently representative [1]. In the present study, we propose to use molecular simulation as a reliable way of predicting volumetric and thermodynamic properties of representative pure compounds in the C10-C36 range and apply the approach to real crude oils. Until recently, Monte Carlo methods and VLE-targeted forcefields were limited in molecular size to approximately 200 g /mol (i.e. C15) with a poor coverage of multifunctional polycyclic molecules [15]. Recent improvements of algorithms for treating polycyclic molecules and forcefield extensions have allowed to increase the coverage of VLE property prediction from simulation [16]. From this work, the average absolute deviations are less than 1% on liquid density and less than 7 K on boiling temperatures, a good level of performance compared with widely used group contribution methods [21,22]. Another advantage is the good performance in predicting density saturation pressures at typical reservoir conditions
ACS Paragon Plus Environment
Page 4 of 48
Page 5 of 48 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
(298 to 550K). Critical coordinates may be also predicted from simulation, but they are less relevant data for heavy fractions that are thermally unstable at critical temperature, so they will not be prioritized. Monte Carlo simulations being computationally intensive, the predicted properties will be used to produce density correlations or to calibrate EOS to avoid increasing computation time in reservoir simulation. In the present work, we intend to show that the new information on the C10-C36 fractions allows to model consistently nine crude oils of different origins with more chemical insight. This requires a specific optimization method to reconcile the analytical and thermodynamic data available on a given crude oil sample. The present article is organized as follows. In a first section section 2, we expose briefly the theoretical background behind the molecular simulation methods and correlative approaches, aiming mainly at the prediction of VLE properties. We expose the basis of the optimization method to model the composition of reservoir fluids. In section 3, we select a set of representative compounds of heavy fractions in the C10-C30 range in order to fill the most important knowledge gaps with properties obtained by molecular simulation. The VLE properties of these compounds are used to create correlations with carbon number that we can be applied to a wider array of compounds, up to C36, quickly and easily. In section 4 we apply the approach on a set of nine real crude oils with different types of analysis, trying to synthesizes several types of input: stock tank oil (STO) density, average molecular weight, detailed C1-C9 analysis from GC, Saturates/Aromatics/Resins/Asphaltenes analysis of C15+ (SARA), pure component properties from experiments and from simulations, and general compositional trends. In section 5, the significance and applicability of the whole approach is discussed. 2. Methods a. Molecular simulation
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Monte Carlo methods used in this study have been exposed in some detail in recent articles [16,23], therefore, we are only providing a brief summary below. i. Monte Carlo simulations Monte Carlo simulation in the isochoric Gibbs ensemble (GEMC) [24] is used to calculate the VLE properties of pure compounds. These simulations involve two coexisting phases without an explicit interface. Periodic boundary conditions (PBC) and minimum image convention are applied to each phase [25]. The total number of molecules is kept constant, but molecules can transfer between phases. The volume of each phase can change, while the total volume is kept constant. The elementary Monte Carlo moves comprise internal moves and molecular transfers with configurational bias [26] or reservoir bias [27,28] . All simulations are done with MedeA GIBBS [29]. Statistical uncertainties are evaluated systematically on computed properties by block averaging [25] Reservoir bias is used with improvements indicated in [16] for branched and polycyclic molecules. At minimum of four (4) simulations (i.e. four different temperatures) with a temperature spacing of 25-40 K are performed for each compound in the two-phase region between the normal boiling temperature and the critical temperature These simulations provide mainly the vaporization enthalpy, the saturated liquid density, and the saturation pressure. The compressibility coefficient Z=PV/RT of the vapor phase is also computed to apply a quality test [30] i.e. check that Z approaches 1.0 for low temperatures (see Supplementary Information). The critical temperature 𝛵𝑐 and critical density 𝜌𝑐 are obtained by fitting a near-critical scaling expression [26]. MC simulations of the liquid phase in the NPT ensemble are essential to predict properties much below boiling temperature. They allow the computation of saturated liquid density
ACS Paragon Plus Environment
Page 6 of 48
Page 7 of 48 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
and intermolecular energy Uext, which are used in thermodynamic integration (TI) to obtain saturation pressures down to standard temperature (298 K) i.e. reduced temperatures as low as 0.3. Convergence analysis in MedeA® is used to automatically detect the interval in which statistical averages and statistical uncertainties can be safely determined. Plots of vaporization enthalpy (see Supplementary Information) and saturated liquid density with temperature are created from the simulation results and analyzed to help detection of possible convergence problems. ii. Thermodynamic integration Thermodynamic integration (TI) allows computing the saturation pressure at temperatures lower than the normal boiling point [31]. For this purpose, the Clausius Clapeyron equation is written in the following form [32]: ∂(ln (𝑃𝑠𝑎𝑡)) ∂(1/𝑇)
𝑇𝛥𝐻𝑣𝑎𝑝
= ― 𝑃𝑠𝑎𝑡𝛥𝑣𝑣𝑎𝑝 ≈ ―
𝛥𝐻𝑣𝑎𝑝 𝑅
= 𝑓
(𝑇1)
(1)
The numerical integration of equation (1) combines simplicity and accuracy to yield Psat(Ti) down to the lower temperature investigated (e.g. 298 K) with single phase NPT simulations . Iii Forcefields We use United-Atom forcefields (TraPPE-UA and AUA) because they allow to study larger molecules than what would be feasible using All-Atom models. In the original TraPPE-UA forcefield, the bond lengths are rigid [33-36] and multi-functional poly-cyclic molecules are typically not covered. For the polycyclic compounds, either naphthenic, aromatic or thiophenic, we use a recent extension of TraPPE-UA [16]. The extension includes, among others, the explicit account of bond stretching and angle bending with harmonic potentials, making thus TraPPE-UA fully flexible. It involves also five additional united-atom groups.
ACS Paragon Plus Environment
7
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 8 of 48
In addition to the extended TraPPE-UA forcefield used for polycyclic compounds, we have been using the Anisotropic United Atoms (AUA) forcefield for n-alkanes, iso-alkanes, alkylbenzenes and alkyl-cyclohexanes. The AUA description [37] is also based on United-Atoms, but the force centers are not located at the carbon atomic centers. We use here the same parametrization as introduced in earlier publications [32,38,39]. b. Correlation of VLE properties When saturation pressures are available from GEMC simulations and TI, we fit the ln(Psat) values onto a second order polynomial expression of 1/T: ln(Psat) = a/T2 + b/T + c
(2)
The three-parameter equation (2) is used to evaluate the normal boiling temperature (Tb) at a reference pressure P0, the critical pressure, Pc , and the acentric factor, .
𝑇𝑏 =
(
―𝑏 ― 𝑏2 ― 4𝑎(𝑐 ― ln (𝑃0)) 2𝑎
𝜔 = ― 1 ― log10 (
𝑃𝑠𝑎𝑡,0.7 𝑇𝑐 𝑃𝑐
)
―1
(3) (4)
)
The correlation of liquid density over the whole scale from ambient to critical temperature is made with a modified Rackett equation involving the critical coordinates Tc and VC: 𝜌𝑆(𝑇) =
𝜌𝐶 (1 ― 𝑍𝐶
𝑇
𝐴
(5)
𝑇 𝑐)
where 𝑍𝐶 and A are component-specific regressed parameters.
ACS Paragon Plus Environment
8
Page 9 of 48 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
Once the standard liquid density at 15°C is described with equation (5) for the reference compounds of a given family (e.g. polycyclic naphthenes) we regress the coefficients 𝑎, 𝑏 of the empirical equation, in which NC is the number of carbon atoms: 𝜌𝑙𝑖𝑞, 𝑠𝑡𝑑 = 𝑎 + 𝑏/𝑁𝐶 (6) c) Compositional modeling of crude oils When the C1-Cm composition is known from GC analysis, the unknown part of the C10-C36 distribution is modelled with an exponential distribution : 𝑤𝑛 = 𝑤𝑚 𝑒𝑥𝑝( ― 𝑏 ∙ (𝑛 ― 𝑚))
(7)
where m is the carbon number starting the regressed distribution, 𝑤𝑚 is its mass fraction and b is a constant. This exponential distribution is found to be suited for modelling single type II reservoir fluid, unaltered by fractionation or biodegradation. In this work, it is assumed that the same b parameter applies to all the sub-families (n+iso-, naphthenes, aromatics and NSO). Based on the measured C10-Cm-1 distribution and exponential Cm-C36 distribution (7) the mass fraction of each pseudo component is given by: 𝑤𝑛,𝑛 ― 𝑖𝑠𝑜 = 𝑋𝑛 ― 𝑖𝑠𝑜𝑤𝑛 (8a) 𝑤𝑛,𝑐𝑦𝑐 = 𝑋𝑐𝑦𝑐𝑤𝑛
(8b)
𝑤𝑛,𝑎𝑟𝑜 = 𝑋𝑎𝑟𝑜𝑤𝑛 (8c) 𝑤𝑛,𝑁𝑆𝑂 = 𝑋𝑁𝑆𝑂𝑤𝑛 (8d) 36
𝑤36 + = 1 ― ∑𝑖 = 1𝑤𝑖
(8e)
ACS Paragon Plus Environment
9
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 10 of 48
The C36+ fraction is described as a mixture of four pseudo-components of the same families (n-iso, cycloalkanes, aromatics, NSO) with identical weights percentages 𝑋𝑛 ― 𝑖𝑠𝑜, 𝑋𝑐𝑦𝑐, 𝑋𝑎𝑟𝑜 and 𝑋𝑁𝑆𝑂 . The properties of these pseudo compounds will be assigned in agreement with asymptotic behaviour of property trends (see section 4). The optimized oil composition is thus fully described by:
the experimental distribution of the light fractions (C1-C5 and C6-C9)
the experimental overall distribution C10-Cm-1
two parameters 𝑤𝑚 and b for the exponential distribution Cm -C36
the percentages 𝑋𝑛 ― 𝑖𝑠𝑜, 𝑋𝑐𝑦𝑐, 𝑋𝑎𝑟𝑜 and 𝑋𝑁𝑆𝑂 applying beyond C10 in the general case (and to the C6-C9 fraction if data are missing)
the average molecular weight of the C36+ fraction
As the weight fractions must sum to 1.0 , the number of independent parameters to be regressed is 5 in the general case (m
