Combined Experimental, Theoretical, and Molecular Simulation

Subscriber access provided by UNIV OF DURHAM

Fossil Fuels

Combined experimental, theoretical and molecular simulation approach for the description of the fluid phase behavior of hydrocarbon mixtures within shale rocks Carmelo Herdes, Camille Petit, Andrés Mejía M., and Erich A Muller Energy Fuels, Just Accepted Manuscript • DOI: 10.1021/acs.energyfuels.8b00200 • Publication Date (Web): 05 Apr 2018 Downloaded from http://pubs.acs.org on April 7, 2018

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 28 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

Combined experimental, theoretical and molecular simulation approach for the description of the fluid phase behavior of hydrocarbon mixtures within shale rocks Carmelo Herdesb, Camille Petita, Andres Mejíac and Erich A. Müllera,* a

Department of Chemical Engineering, Imperial College London, U.K. b Department of Chemical Engineering, Bath University, U.K. c Departamento de Ingeniería Química, Universidad de Concepción, Chile * Author to whom correspondence should be addressed ([email protected])

Abstract An experimental, theoretical and molecular simulation consolidated framework for the efficient characterization of the adsorption and fluid phase behavior of multicomponent hydrocarbon mixtures within tight shale rocks is presented. Fluid molecules are described by means of a top-down coarse-grained model where simple Mie intermolecular potentials are parametrized by means of the statistical associating fluid theory (SAFT). A four component (methane, pentane, decane, naphthalene) mixture is used a surrogate model with a composition representative of commonly encountered shale oils. Shales are modelled as a hierarchical network of nanoporous slits in contact with a mesoporous region. The rock model is informed by the characterization of four distinct and representative shale core samples through nitrogen adsorption, thermogravimetric analysis and contact angle measurements. Experimental results suggest the consideration of two types of pore surfaces; a carbonaceous wall representing the kerogen regions of a shale rock and an oxygenated wall representing the clay-based porosity. Molecular dynamics (MD) simulations are performed at constant overall compositions at a temperature of 398.15 K (257 °F) and explore pressures from 6.9 MPa up to 68.95 MPa (1000 to 10000 psi). Simulations reveal that it is the organic nanopores of 1 nm and 2 nm that preferentially adsorb the heavier components, while the oxygenated counterparts show little selectivity between the adsorbed and free fluid. Upon desorption, this trend is intensified, as the gas phase in equilibrium with a carbon nanopore becomes increasing leaner (richer in light components) and almost completely depleted of the heavy components which remain trapped in the nanopores and surfaces of the mesopores. Oxygenated pores do not contribute to this unusual behavior, even for the very tight pores considered. The results presented elucidate the relative importance of considering both the pore size distribution and the heterogeneous nature of the confining surfaces when theoretically describing adsorption and transport of oil through shale rocks and provide a plausible explanation to the abnormal continuous leaning of shale gases seen during field production.

ACS Paragon Plus Environment

1

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 2 of 28

1. Introduction Oil and gas production from unconventional shale resources, enabled by technologies such as horizontal drilling and hydraulic fracturing, has radically changed the worldwide energy and petrochemical landscape1. In spite of the exponential surge in production, many uncertainties surround the chemical and physical mechanisms occurring in these ultra-tight rocks2. A recent NSF workshop 3 highlighted the main challenges associated with the fundamental description of the processes involved with production from the above-mentioned natural resources. Amongst them were key questions regarding the adsorption of hydrocarbon mixtures in nanoporous shale rocks: What are the interactions of the fluids with the rock? How does the extreme nanoporosity affect both the flow and the ultimate desorption of the adsorbed species? What heavy components, and in what proportion are being left behind after primary production? Can the adsorption behavior be altered to control the amount of the heavier hydrocarbons left behind? How can one rationalize the continuous leaning of production fluids even when the reservoir pressures are above the expected bubble (dew) points of the unconfined fluids? Unfortunately, no suitable theory is currently capable of predicting the behavior of confined fluids with the quantitative accuracy and confidence required by the industry. Methods to understand the fundamental physics involved in the ultimate modelling of these systems require both experiments at challenging conditions of extremely high resolution and the visualization of such processes through molecular simulation. This paper deals with these questions and provides a framework to study the adsorption of hydrocarbon mixtures in confined nanoporous media with an enhanced level of confidence by means of using fine-tuned coarse-grained fluid models. In general, reservoir engineering studies rely on accurate knowledge of thermodynamics and transport properties associated with both the porous media and the fluids. In conventional reservoirs, the pores and pore-throats are significantly larger than the molecular dimensions and in these scenarios the interactions with the pore surface have a relatively small effect on the bulk thermodynamic properties. At this scale, continuum theories hold and are successfully employed in a routine fashion. Conversely, in the case of unconventional shale reservoirs, the pore and pore-throats are much smaller; frequently in the nanometer range 4. It is expected that fluid phases confined in such reduced pore space may exhibit physical and chemical behavior different from bulk phase owing to the interaction between the molecules and the pore surface5 . As a result, the thermodynamic properties of the reservoir fluid could be significantly deviated.6 A present challenge in the reservoir engineering community associated with unconventional shale is that the thermodynamic relations at the macroand nano- scales need to be addressed in a self-consistent way but there are currently no experimental methods to determine the alteration of reservoir fluid thermodynamic properties due to confinement in nanopores. Empirical approaches such as the use of macroscopic equation-of-state (EoS) with modified critical properties6, 7 or ad-hoc adsorption regions and pore size distributions8,9,10 have been proposed, however, they run the risk of being physically unsound and inconsistent. Noticeable exception are


2

Page 3 of 28 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

modern density functional treatments11 which hold promise in describing confined fluid density distributions. There is an opportunity for simulation studies to help to develop the next generation of thermodynamic models and engineering tools to support future unconventional shale developments. While our understanding of the effect of confinement on fluid properties is broad 12 , 13 , shale oils present new challenges due to the multicomponent (and sometimes multiphase) nature of the fluid and the inherent chemical and morphological heterogeneity of the confining media. All these conditions and uncertainties pose challenges to the systematic modelling of these systems. Previously published molecular dynamics simulations studies of shale systems have focused mainly on gas phase fluids consisting of either single or binary component systems8 , 14 , 15 , 16 , 17 , 18,19,20,21 . A notable exception are the simulation studies derived from the work by Sing et al.6 that used molecular simulations to determine how the critical properties of single components of methane, ethane, propane, butane, and octane were shifted within nanopores. Although the confinement will have an effect on the bubble (or dew point) of mixtures in shales22, it is not clear how (and if) the difference in criticality between a confined fluid and a bulk one can (or should) be used directly in EoS macroscopic models and simulators23. Complicating matters further, the existing literature debates about the use of surrogate models for shales, in some instances using clays and in others using graphite structures to model the porous regions of shale rocks. It is becoming clear that both types of regions coexist within an inorganic matrix and the issues with regards to the extent of the porosity and its distribution, the rock morphology and connectivity are still being experimentally resolved24. To the best of our knowledge, no previous molecular dynamics simulation studies have investigated how confinement alters the thermodynamics of multicomponent mixtures and multiple phases at representative reservoir conditions on even the most simple of pore models. This manuscript explores the experimental characterization of four shale rock samples as a backdrop to develop a molecular model that captures in a realistic fashion the geometry and energetical heterogeneity of the rock. The experimental results are employed as a basis for the development of a hierarchical pore model, based on the premise that two types of porous media may be present, one carbonaceous in nature and another clay-like. The model parameters for the fluids are obtained by direct parameterization of a molecular-based EOS, while the solid-fluid interactions are mapped back from experimentally-measured macroscopic contact angles. 2. Rock Characterization Four shale core samples were employed, referenced herein as A, B, C and D, respectively. The shales are borehole samples, taken from four different US reservoirs and chosen to provide the widest possible range of physical characteristics. Samples A and B are consolidated rocks, with a lighter gray color, while the latter two, C and D are more brittle and much darker in appearance. Figure 1 shows examples of scanning electron micrographs (SEM) of the samples. The results seem in line with the expected heterogeneity seen in similar reservoirs25. The SEM images show some commonalities,


3

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 4 of 28

such as a lack of macroporosity even at the highest levels of magnification. Less magnified images are presented in the Supporting Information along with chemical analysis from X-ray diffraction (XRD). Rock B seems to have the least mesoporosity, with large areas presumed to be of quartz-like and/or calcium carbonate materials. Sample A is similar, but shows some regions of heterogeneous patches, possibly housing kerogen and clay patches. Sample C shows some clear layering, consistent with larger regions of kerogen material. Sample D is the least consolidated material.

Fig. 1 Scanning Electron Microscopy (SEM) images of the four samples at 1200x magnification. Further images are provided in the Supporting Information. Table 1 summarizes the results obtained from an analysis of nitrogen adsorption on the samples. Nitrogen adsorption curves and experimental details are included in the Supporting Information. The pore size distribution is obtained from these isotherms using standard density functional theory-based packages. The pore size distribution of the samples is a roughly unimodal broad porosity from 3 nm onwards until around 2040 nm, with a rather flat distribution as exemplified by the plots in figure 2. This is consistent with results reported21,24, 26 in the literature for shales. Notwithstanding, recent experiments27 have convincingly reported that the porosity results obtained from pristine samples are radically different than those obtained with crushed samples, suggesting that crushing the samples (as is standard practice in adsorption studies


4

Page 5 of 28 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

including this one) may expose different pore structures to the adsorbent. Results in Table 1 should thus be taken only as a reference.

Fig. 2. Pore size distribution, expressed as the cumulative volume V for each individual pore width, W, for the four samples as back-calculated from nitrogen adsorption isotherms. Table 1. Textural parameters of the samples determined by N2 sorption at 77 K as well as compositional parameters derived from thermogravimetric analyses. Sample

Surface area (m2/g)

Total volume of pores (cm3/g)

Kerogen mass %

Inorganics mass %

Residue mass %

A B C D

10.1 +/- 0.1 48.7 +/- 2.5 14.3 +/- 0.4 3.9 +/- 0.1

0.015 0.055 0.026 0.007

0.05 1.4 8.3 1.9

39.5 40.1 24.2 35.7

41.3 41.7 34.3 38.1


5

Energy & Fuels

100

90

Mass (%)

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 6 of 28

80 Sample A

70

Sample B Sample C Sample D

60

50 0

200

400

600

800

Temperature (°C)

Fig. 3. Thermogravimetric analysis (TGA) of shale samples in air. Figure 3 includes the thermogravimetric curves for the relevant temperature range. Sample C showed a step decrease in the region of 450-550 °C most likely indicative of the presence of a significant amount of kerogen in the sample and a second decrease in mass above 600 °C, most likely suggesting the decomposition of calcite or dolomite or magnesite28. The other three samples showed similar behavior amongst themselves, with only a very a small decrease in mass in the region < 550 °C and a significant drop above 600 °C. All samples stabilized at circa 60 % mass at high temperatures. Taking the mass loss between 200-550 °C as corresponding to organics and kerogen and taking the mass loss between 650-850 °C to be mineral, one can estimate the percentage (mass %) composition of the materials, as presented in Table 1. The residue corresponds to the material left after heating to 900 °C. 3. Contact angle measurements A pendant drop tensiometer is used to determine the fluid-rock interfacial tensions at high pressures and temperatures. Details of the equipment and experimental procedures are given in the Supporting Information. In the case of porous materials, as happens with the shale samples, it is not possible to subtend a drop of fluid on the surfaces due to the fluid flow through the porous material. This flow creates a permanent wetting condition where the contact angle is equal to zero. In order to measure the contact angle in such porous material, the drag force can be balanced by using a needle, as is shown in Figure 4. This experimental technique has been used previously29 for clays and shales, and provides a route to obtain true contact angles.


6

Page 7 of 28 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

Figure 4. (Left) Image of the view cell of the pendant drop tensiometer loaded with a porous rock sample. A drop of the fluid is placed on the surface by the top needle. However, care is taken to “drag” the drop in the upward vertical direction, creating a fluid meniscus between the surface and the needle. The subtended angle can be used as a measure of the solidfluid contact angle. (Right) Actual snapshot from the cell with an overlay of the contact angle. We consider three fluids, which give us a benchmark for understanding the fluid-solid interactions: n-decane, water and CO2; three different temperatures: ambient and two representative of a reservoir (30°C, 80°C and 120°C); and four different pressures: one ambient and three at typical reservoir conditions, 0.1, 13.79 , 20.68 and 27.58 MPa (14.7, 2000, 3000 and 4000 psi). The exceptions were the experiments with CO2, which are performed only at 20°C and 4.83 MPa (700 psi), so to be removed from the fluid critical conditions. Additionally, decane at 120°C was too close to its boiling point (174.1°C), so this particular state point was not considered. Not all permutations were made, as they were deemed too repetitive, so samples A & B were subject to all tests, while C & D only to a subset of the ones deemed most interesting. Water contact angles measurements were performed on samples that were later discarded for fear that water uptake could swell the existing clay regions and/or otherwise alter the original rock. Furthermore, a few benchmark cases with model standard surfaces are produced, namely with graphite and mica. Surfaces were prepared in accordance to the protocol laid out by Zelenev30 as applied to shale solids. In short, this accounts for dry polishing to a shiny surface with a sequential use of decreasingly fine grit sandpaper (mechanical polishing). The preparation of the surface is determinant in the outcome of contact angle measurements, however we include a few results with unpolished surfaces for comparison. The results are qualitatively similar for all the rock samples studied. Decane and CO2 both exhibit complete wetting of the surface and contact angles must be measured as described in the previous section. Water exhibits only partial wetting and the contact


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 28

angle measurement is done without the need of the injection needle. Average values for the contact angle of water at 30 °C (44.2°, 43.8° , 32.3° , 52.4° for samples A-D respectively) agree with the expected values reported31 for other shale rocks. Figure 5 shows a set of characteristic results of the contact angles of decane as a function of pressure at different temperatures for the different samples. The raw data is tabulated in the Supporting Information.

Fig. 5. Comparison of experimental contact angles of decane over shale rock samples as a function of temperature and pressure. (Left 30ºC, Right 80ºC) 4. Molecular model for shale reservoirs While a fully atomistic model of both the fluids and the solid framework could in principle be devised, the level of detail that this implies is at odds with the level of characterization that is known for shale oil systems. The result of this disparity is the development of models where the uncertainty of the actual chemistry hampers the applicability of the results. Furthermore, all-atom descriptions of both solids and fluids result in extremely high computational demands which in turn limit the complexity of the fluid mixture and the solid morphology, both of which turn out to be crucial in understanding the physics of shale oil adsorption and desorption. Notwithstanding, “heroic” simulations of mixtures of gases in atomistically-detailed kerogen models32,33 have been reported16,34,35 and the results are of interest. 4.1 SAFT force field The quality of the intermolecular potentials determines the accuracy of the property prediction of any molecular simulation. In this manuscript, we pay particular attention to the description of fluids with force fields that are traceable to experimental data; i.e. the parameters are obtained directly by fitting the vapor-liquid equilibria curves corresponding to the vapor pressure and saturated liquid density employing a molecular-based equation of state, namely the so-called SAFT-VR-Mie 36 equation of state. The resulting force field is based on the use of the Mie intermolecular potential uMie, which can be represented by:


8

Page 9 of 28 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)

In this expression, λ is the repulsion parameters of the intermolecular potential, r is the center-to-center distance of the interacting segments, ε is the energy scale corresponding to the potential well depth, σ is the length scale, corresponding loosely with an effective segment diameter. The value of the attractive exponent is fixed to 6 without any loss in generality37. The Mie potential reverts to the well-known LennardJones model if the repulsive exponents is taken as 12. The key point of the selected force field is the use of the EoS to provide for the parametrization of the potential parameters. A review of the SAFT methodology and its applications is available38 and recent publications highlight the successful application to complex fluids and mixtures39,40,41,42,43. Within this framework, the EoS is fitted to the thermophysical properties of real fluids. By performing the optimization, which requires a negligible computer power, one obtains from the theory a set of adjustable parameters, (l, s, e ). The crucial point is that these parameters correspond to the values of the molecular descriptors of the force field, hence the equation of state is effectively a means for performing a “top-down” parametrization. The three parameters obtained are effective values which represent a sensible compromise capable of describing with accuracy the behavior of pure fluids in a wide range of temperatures and pressures. For elongated molecules, the theory allows for a further parameter, m, the chain length, corresponding to the number of tangent isotropic Mie beads that conform a model of a molecule. Since the parameter fitting is done at the level of the EoS, one is implicitly defining a level of molecular resolution which is essentially using a coarse-grained representation to describe the fluid molecules. Coarse-graining the molecular details in the way done herein allows one to smear out some of the irrelevant degrees of freedom, e.g. bond vibrations, changes in electron density, while focusing on the more relevant traits, such as the average length-to-breadth ratio of an elongated molecule, the corresponding fluid density and the average energetics. In this work, the coarse-graining is done at a level where several heavy atoms are grouped in one bead. 4.2 Shale fluids A detailed description of a crude in terms of individual components is clearly an impossible task. Nevertheless, one can obtain important insights by considering the behavior of simple mixtures which share a proportionality resembling the natural oil. Within this mindset, we consider four representative real surrogate components44 for the description of a synthetic reservoir fluid. A C1 molecule stands for the volatiles. Pure methane fluid properties were used to obtain the parameters for this cut described by a single coarse-grained bead (m=1). A C5 molecule represents the light hydrocarbon molecules. A two-tangent coarse-grained bead (m=2) was selected for this cut and was parameterized with pure pentane fluid properties. The heavier compounds were described with a CX molecule, representing a 10-carbon atom backbone. Available n-


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 28

decane fluid properties were used to develop the force field for this coarse-grained chain of 3 tangent beads (m=3). Finally, the aromatic fraction NP was modelled using a coarse-grained dimer (m=2) with parameters obtained from the fluid properties of naphthalene. A light synthetic crude was considered using the above-mentioned components and compositions given in table 2. Table 2 Molecular composition of the synthetic shale oil. Color code is carried out throughout the manuscript. CG models are a guide to the eye and are not to scale. component C1 (methane) C5 (pentane) CX (n-decane) NP (naphthalene)

% molar 42 25 16 17

Fig. 6 Coarse grained representation of the molecules used herein. The colored spheres represent the Mie beads and are superimposed on an atomistic representation of the actual molecule as a guide to the eye. Coarse-graining is performed through a top-down approach which does not take into account the underlying atomic structure. C1 (gray is Methane; C5 (blue) is n-pentane, CX (green) is n-decane and NP (red) is naphthalene. Color code is carried out throughout the manuscript. CG models are a guide to the eye and are not to scale. The calculated SAFT-g Mie parameters for this study are given in Table 3. A short-cut calculation method for the parameters45, has been employed. The procedure amounts to recognizing that the SAFT EoS can be described in terms of corresponding states, hence the parameters can be intimately related to critical properties; e.g. the repulsive exponent l is directly related to the acentric factor, the energy scale e is relatable to the critical temperature and the approximate segment diameter s relates closely to the density of the saturated liquid. The quality of the fit, in terms of the fluid phase behavior is shown in Figure 7. A more detailed description of the application and validation of


10

Page 11 of 28 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

the SAFT force field to reservoir gases is given elsewhere46. This method, wrapped into a simple webpage47 facilitates the use of the force fields in modern parallel computers48. Table 3 SAFT-g Mie force field parameters employed (per bead).. kB refers to Boltzmann’s constant. Component methane pentane n-decane naphthalene silica (wall) kerogen (wall) graphite (wall) waterb carbon dioxidec

m 1 2 3 2 1 2

s (nm) 0.37523 0.42477 0.45841 0.46228 0.32700 0.39300 0.34000 0.29148 0.28485

e/kB (K) 170.75 317.50 415.19 557.75 128.15 45.81 28.00a 378.87 194.94

l 16.39 16.06 20.92 19.50 12.00 8.40 14.65

a

For contact angle calculations this value varies from 3.75 to 35 K A more generally applicable SAFT model for water can be found in ref 49 c See also ref 50 for an alternative model. b

The unlike parameters, used to describe the cross interactions, were obtained by applying the combining rules suggested by Lafitte et al. (ref. 36) 𝜎#2 =

.33 4.55 6

; 𝜀#2 =

9 9 .33 .55 9 .35

𝜀## 𝜀#2 ; 𝜆#2 − 3 =

𝜆## − 3 𝜆22 − 3

(2)

where the subscripts ii and jj refer to pure component parameters and ij refers to the cross parameter. The energy cross parameter, eij, appears to be the most critical parameter to describe interactions between different components. Although in principle possible, no adjustment in the cross terms is made, e.g. no binary interaction parameters are employed. All bonded segments within a dimer or chain are bonded rigidly at a distance of s (in practice, a stiff spring with equilibrium constant fixed at s provides equivalent results). In the case of CX ( decane ), the three spherical segments are kept elongated by adding a bond angle bending potential, 𝜙=>?@$ , between the three consecutive beads, 𝜙=>?@$ = 𝑘=>?@$ (𝜃 − 𝜃C )6 , where q is the angle subtended by three consecutively bonded spheres, kangle = 22.18 kJ.mol-1.rad-2 and q0 = 157.6 °. All nonbonded interactions ( i.e. 1-3, 1-4 ) are included.


11

Energy & Fuels

800

600

Temperature (K)

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 12 of 28

400

200

0

0

200

400

600

800

1000

1200

3

Density (kg/m )

Fig. 7 Vapor-liquid equilibria (VLE) for the pure components of the synthetic shale oils parameters. Solid lines are smoothed experimental data retrieved from the NIST webbook51, symbols are molecular simulations for methane (black), pentane (blue), decane (green), naphthalene (red) and carbon dioxide (violet) using parameters from Table 3. 4.3 Fluid-solid interactions The solid-fluid interactions are characterized by relating them to fluid-solid contact angles. The procedure for obtaining the contact angle from simulations is well established52 and we employ that method in our analysis of the results. As an example, in Figure 8 we show two equilibrium snapshots, obtained from coarse-grained molecular dynamics simulations, of the contact angle of decane on a generic substrate with two different values of the solid-fluid energy. The wide variation of behavior observed serves as a basis for suggesting that one could, though simulations, find a given a model surface energy that, on average, reproduces the experimental contact angle behavior of a given fluid on a real rock sample. Since our models of fluids give the correct vapor-liquid interfacial tension53, the contact angle is then only a function of the surface-fluid interactions (neglecting the solid-vapor tension) hence can be tuned to reproduce the experimental data. Coarse graining the surface has the advantage that it “smears” out the structural and energetic heterogeneities, giving in a unique parameter set an effective solid-fluid interaction. We perform simulations for coarse grained models of water and decane at different representative solid-fluid energies to obtain a profile of contact angles as a function of surface energy. For each surface, we run eight simulations at different energetics. To keep the number of variables to a minimum we use a rather smooth surface (i.e. the molecular rugosity is not taken into account) characteristic of a graphite plane. Similarly, the size and bonding distances of the surface beads are kept constant. Temperature is kept fixed at the conditions of the experiment and a relatively small quantity of the model fluid is placed in a simulation box with a surface. Periodic boundary conditions are applied in the direction of the surface, but in not in the direction normal to it, where the opposing wall is reflective (hard wall). Since the vapor pressure


12

Page 13 of 28 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

of the fluid is significantly less than atmospheric, the system is pressurized using an appropriate number of nitrogen molecules. In order to reach the macroscopic limit, we use a relatively large simulation cell with over 5000 fluid molecules.

Fig. 8. Snapshots from simulations of droplets on model surfaces (green). Decane (gray) is placed in the simulation cell and the pressure is maintained by filling the cell with nitrogen (blue). The solid-fluid interaction is varied by changing the potential well of the solid particles. Left figure corresponds to a low energy ( e/kB = 3.5 K ) while the right figure corresponds to a rather attractive wall (e/kB = 28 K ) Temperature is 30 oC and pressure 0.1 MPa (14.7 psi). Not all values of the wall energy will produce a drop for which the contact angle could be measured. For example, simulations for CO2 were not carried out, as the thermodynamic conditions were too close to the critical point. This implies that the vapor phase is almost as dense as the liquid and the formation of a stable drop/film is never achieved. From the simulations, a cloud point data algorithm (ref. 52) is employed to back trace the effective average contact angle, which is plotted as a function of the wall energy in figure 9 for decane. The energy-angle plots show a linear trend which has been reported before54 in similar studies. One can observe that by overlaying the experimental data, that there is a correspondence, in the sense that roughly the same value of the surface energy produces the correct (experimental) contact angle. It also becomes clear that on this plot there is hardly any distinction between the average experimental values themselves and the value corresponding to a graphitic wall (εss/kB = 28 K). A similar analysis can be performed on the plot corresponding to the water droplets (See Supporting Information, Figure S21). Here, the correspondence between the different experimental points and the graphite wall energy is not so clear, but again, it clearly points out the fact that although the rock has a very low percentage kerogen (as seen in the TGA analysis), the use of a graphitic wall as a surrogate for the rock surface is not a poor approximation in terms of the average observable energetics. Several researchers have suggested the use of graphitic pores to model shale rocks55,56 sometimes without any appropriate justification23,23, 57 , while other have used a diametrically different approach, e.g. the use of clay models18 which would be inconsistent with our experimental contact angle data and adsorption studies56. We will consider in our models the two extremes, a carbonaceous surface, fitted to the properties of a graphite surface and an oxygenated surface, fitted to the properties of mica. It is most plausible that actual rock samples will be a combination of both pores. While not


13

Energy & Fuels

the focus here, there are theoretical pathways for obtaining coarse grained surface potentials for heterogeneous surfaces 58 . These methods, of course, would require previous knowledge of the atomistic-detail with respect to the energetic and morphological heterogeneities of the rocks.

30 24.5

25

Solid-fluid (wall) energy (K)

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 14 of 28

Sample D Graphite

Sample C Sample B

Sample A

Sample D Sample C Graphite

24

Sample B

23.5

Sample A

23

20

20

22

24

26

15

10

5

0 0

20

40

60

80

100

120

140

Contact angle θ (o)

Fig. 9. Contact angle θ (o) as a function of the solid-fluid (wall) interaction energy for decane at 30 oC and 0.1 MPa. Fluids is a coarse-grained model with fluid-fluid energies given in Table 3. The solid-fluid energy of the surface (wall energy) is varied. Contact angles are calculated from simulations and shown as open triangles. Red circle corresponds to a graphite surface energy. Solid squares are the experimental results on real surfaces. Solid line is a best fit trend line. Insert shows a zoom of the region of interest. 4.4 Model porous media The modelling of the porosity in shales is a challenge since neither the pore size distribution nor the actual composition of the rock surface is well known nor consistently uniform across different reservoirs. With this in mind, one can only attempt to describe the most general traits expected to be important for the adsorption processes. The pore size distribution suggests nanopores in the range of 2 nm with larger mesoporosity in the range up to 50 nm.26,24 The surface itself is most likely composed of a heterogeneous patch-like mixture composed of regions of carbonaceous deposits akin to kerogen, with significant porosity, regions of clay patches 59,60,61, all inside an inorganic matrix. There is little current understanding with respect to the connectivity amongst these regions or the actual path that shale oil takes within these structures. To complicate matters further, and as expected for natural rock formations, there is significant difference between the characteristics of different reservoirs4. For our simulation studies, we have concentrated on the effect of the nanoporosity in the adsorption of shale fluids. It is in this range of confinement where one expects the most profound effect on the phase behaviour of the fluids. On a similar vein, we have


14

Page 15 of 28 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

considered two extreme surface scenarios, where the surfaces are made of oxygen and others where they are carbonaceous in nature, surrogates of clay and kerogen regions17, 62 . The accuracy of these approximations is still under debate, although our contact angle measurements suggest that a carbon-like pore would provide the more accurate local average fluid-solid interactions (Figure 9) while the TGA results imply that carbonaceous surface accounts for only a small fraction of the rock. The apparent paradox will be discussed and explained by the results. The walls are modelled explicitly with Mie beads akin to the ones used for the fluids. The parameters for the beads are imported from standard representations of oxygen atoms63 (clay) and -CH2united atom 64 (kerogen) potentials. Solid walls are built using simple face-centered cubic unit cells, with an exposed (100) crystallographic plane for the oxygenated surfaces and a (111) face for the kerogen face. Accordingly, the densities the solid densities (number of spheres per unit cell) are matched to silica: 2695 kg.m3 and graphite 2054 kg.m3 respectively. The rugosities of the two walls are different, as the (100) surface is somewhat rough while the (111) is smoother in nature. Slit pores are considered herein which in spite of the complexity of the shale pore network are only approximations for the complex pore networks. In Figure 10 we present a sketch of the simulation cell showing the focus on the nanopore and the relation of this nanopore to the connecting mesopore. This hierarchical pore size distribution is a proxy for a bimodal mesoporous and microporous system. The solid blocks that make up the pores are approximately 4 nm in the x direction, 10 nm in the z direction and 5 nm in the y direction (into the plane of the paper). While many definitions can be put in place, the pore width is measured as the distance between the planes formed by the centers of (opposing) inner-facing wall atoms. The fluid is placed in intimate contact with the porous media and the resulting equilibrium distributions showcase not only the phase equilibria but also the surface wetting effects and the distribution between mesoporous and microporous adsorption.

Fig. 10 Conceptual schematic of the porous media. The system is composed of 2 nm nanopores which “feed” into a mesopore, an order of magnitude larger. Periodic images are applied in all Cartesian directions, so the system replicates infinitely. In the direction in and out of the page


15

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 16 of 28

the rock slabs are replicated, in essence infinitely. The simulation unit cell is drawn in red. 5. Molecular simulation details All simulations are performed using the Gromacs 65 software suite and correspond to classical molecular dynamics simulations. Visualizations are rendered using VMD 66. The system was simulated under the NPzzAT ensemble, at constant number of particles and temperature. In the NPzzAT ensemble the pressure coupling is isotropic in the x and y direction (the transversal area A is kept constant throughout the simulation), but different in the z direction. As the simulation cell is elongated in the z direction, any secondary phase that will appear in the system will do so attempting to form a slab or film spanning the x-y plane in order to minimize the interfacial area. With this cell geometry, the z component of the pressure tensor is the bulk (macroscopic) pressure, regardless of the presence of multiple phases and is a constant in the simulation. The volume changes affect only the mesoporous region (i.e. the solid is maintained isometric). Unless stated otherwise, the reported properties are averaged over at least 3×107 steps (Δt = 0.01 ps) after the equilibration of the simulated system is attained, as determined by monitoring both its total energy and density. The Berendsen thermostat and barostat were selected as the coupling algorithms. Periodic boundary conditions and a cut-off of 2.0 nm were applied in all simulations. An initial cell of appropriate size to fit 2000 molecules (C1=840, C5=500, CX=320, NP=340 molecules) was constructed. We report results for a temperature of 398.15 K (125 °C) and spanning pressures of 6.9 MPa up to 34.47 and 68.95 MPa (1000, 5000 and 10000 psi). The bubble point of this model mixture is obtained by performing simulations46 of the bulk fluid and is found at 566.2 kg/m3 associated to a bubble pressure of 13.6 MPa (1976 psi), for this composition at 125 °C. 6. Results and discussion 6.1 Behavior of the confined fluid At 125 °C, for pressures above 13.7 MPa (2000 psi), the synthetic fluid is in a compressed liquid state, however the presence of a surface affects this behavior. We place the synthetic fluid in contact with a nanopore, while keeping both the temperature constant and the bulk pressure (Pzz) constant at 35 MPa (5000 psi). Figure 11 graphically summarizes the results of five of the most revealing scenarios.


16

Page 17 of 28 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

Fig. 11 Equilibrium simulation snapshots representing five of the most relevant scenarios. In all cases temperature is 125 °C and pressure is 35 MPa (5000 psi). Simulation cells show a periodic image in the vertical direction for clarity. Color code for molecules as in Figure 6. In Figure 11 we show visually the effects due to the confinement in two pores of different nature but of equal pore width, specifically, a silica (clay-like) pore and a wax (kerogen-like) pore both of 2 nm in width. The inorganic pore did not show a preferential adsorption for any of the components of the crude, and a similar composition is found between the homogenous liquid bulk and the adsorbed phase in equilibrium. The presence of the pore did not affect the bulk phase of the oil. On the other hand, the organic pore induces a phase separation, in which an adsorbed phase, a liquid slab and gas phase are found in mutual coexistence. Moreover, this pore imposes a dramatic change in the composition and redistribution of the local densities, in essence, the confinement induces a change in the bubble point of the fluid, as it is now phase separated. Along with the phase separation, this pore forces a dramatic change in the composition of the local phases, presenting in one extreme an adsorbed phase rich in decane, a liquid slab coating the walls of the mesoporous region composed mainly of pentane and naphthalene and the much leaner gas phase rich in methane. To explore the ultra-small confinement effect, a smaller organic pore of 1 nm width is studied. The extra confinement of the system biased even further the composition of the adsorbed phase towards the heavier compounds; the light gas remained almost unchanged in its composition (as compared to the 2 nm pore), while the wetting layer in the mesopore wall increased in width. The 1 nm pore becomes even more selective to decane (methane, pentane and naphthalene were poorly adsorbed), to the extreme that one can talk of a decane pore blockage, in which a marked hysteresis is seen when trying to retrieve (desorb) the decane. This is an expected behavior, and in effect long alkanes, e.g. nonane, are routinely used by experimental groups to “block” nanoporosity in activated carbons67. Consequently, the naphthalene and pentane molecules, that were found inside the pore of 2 nm, are now constituents of the liquid slab outside the pore of 1 nm, alongside the excess of decane molecules that have not found room in the pore. The liquid slab, extends up to 10 nm in width. The gas phase, which is still rich in


17

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 18 of 28

methane, presents now a molar composition of C1=0.509, C5=0.32, CX=0.049, NP=0.122. The last case shown in Figure 11 exemplifies the effects on a 2 nm wax pore of the addition of 300 molecules of CO2 to the system (representing a final 13% of the total composition). The effect of the CO2 injection is subtle and a numerical analysis of its distribution is needed to extract appropriate information. The CO2 was not particularly adsorbed in the pore and it is clear that the solid-fluid interactions are stronger for the longer alkanes than for CO2. On the other hand, the absorption of carbon dioxide on organic phases is important, hence a maximum in the concentration of CO2 molecules is detected in the liquid slab outside the pore, as seen in the density profiles of Figure 12.

Fig. 12 Number densities profile per component along the simulation box for a 2nm kerogen pore in the presence of 13% CO2. The pore entrance is at 9.4 nm. 6.2 Depressurization of a confined fluid It is expected that the presence of surfaces will have an effect on the bubble point of a given global crude composition68. However, the extent and nature of the expected shift in the phase boundary depends crucially on the nature of the solid-fluid interaction and the extent of confinement 69 . We evaluate the effect of the depressurization of the synthetic mixture (Table 2) in the presence of a 2 nm organic pore at 125 °C. Under these conditions, for pressures lower than 85.6 MPa (12414 psi), two phases (vaporliquid) are detected inside the simulation box. As described in detail in the previous section, the organic pore induces a phase separation, in which an adsorbed phase, a liquid slab and a gas phase coexist in equilibrium. Along with the phase separation, this pore imposed a dramatic change in the composition of the local phases. It its worth having first a qualitative look (Figure 13) at the simulation cell before (103.4 MPa; 15000 psi), near (82.7 MPa; 12000 psi) and much below (1.34 MPa; 200 psi) the phase transition point, to gain insights about this system. At the high pressure, the pore walls of both the nanopore and the mesopore strongly adsorb the longer alkanes (green spheres), so the bulk fluid is, even at this stage depleted of the paraffins. Very little of the light gases are present in the pore and the lighter gases are confined to the regions far from the surfaces. The other heavy


18

Page 19 of 28 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

compounds, such as the aromatics are also preferentially adsorbed at the wall. Figure 12 shows this information in a quantitative way, in terms of number density profiles across the z direction. The fluid in the bulk region becomes progressively leaner as the pressure drops. At a pressure just below the phase change point (82.7 MPa; 12000 psi), one can observe an accumulation of methane in the center of the mesopore, the region farthest away from the wall, indicating the onset of the formation of the gas phase. At a point well below the transition, the light gases are depleted form the pores. The nanopores remain completely filled with the heavier compounds and the mesopore walls retain a film of the heavy compounds.

Fig. 13 Snapshots from three representative simulation points corresponding to NPzzT simulation of the mixture in contact with a 2 nm pore organic pore at 125 ºC. Left figure corresponds to 103.4 MPa (15000 psi), where the fluid remains in liquid state. Middle figure corresponds to 82.7 MPa (12000 psi), slightly below the bubble point of the mixture (85.6 MPa; 12413 psi). A rearrangement of methane molecules (white) in the bulk of the fluid range is apparent, an indication of the onset of an additional phase. Right figure shows a detail of the system at 1.34 MPa (200 psi), a pressure much below the bubble point. One can analyze the local molar compositions from a direct inspection of the number densities inside the pore, in the liquid bulk and in the gas bulk. Due to the way the simulations are run, the low-pressure simulations cells are very extended in the z direction; up to the µm scale. The calculated partition coefficients Ki=yi/wi where yi is the molar fraction in the bulk (liquid or vapor) phase and wi is the molar fraction in the pore volume, are reported in Table 4, at the studied pressure conditions. Table 4 Local K-values at different pressure levels, P. yi is the molar fraction in the bulk (liquid or vapor) phase and wi is the molar fraction in the pore volume

Pressure (MPa) C1 C5 CX NP

103.4 3.859 3.399 0.047 1.151

Ki=yi/wi 82.7 7.799 3.245 0.028 0.438


1.34 6.186 3.356 0.003 0.050

19

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 20 of 28

6.3 Energetically heterogeneous pores Shales are composed of both inorganic and organic heterogeneous regions. To understand the effect of heterogeneity of the pores on the adsorption, we describe here a simulation where within the same simulation cell we place two adjacent pores of the same width ( 2 nm ) but of different energetics, one corresponding to a kerogen pore and the other to a clay pore. Simulations are performed as in the cases outlined above; the adsorbent is placed in contact with a constant composition mixture at a given pressure and temperature and the system is left to equilibrate. The adsorption on these systems is seen to be additive, in the sense that the characteristics observed for each of the individual pores is replicated in each of the corresponding parts of the pore, i.e. the wax-like portion of the pore exhibits a preferential adsorption for the heavier alkanes while the silica-like portion of the pore is very non-selective. As a consequence of the above the bubble point pressure is now shifted to value of 41.47 MPa (6000 psi) which is intermediate between that calculated when the mixture is in contact with a purely kerogen pore 85.6 MPa (12414 psi) and the bubble point of the homogeneous system, 13.6 MPa (1976 psi). From a visual analysis of the simulation movies (not included), it is observed qualitatively that the kinetics of the adsorption within each of the sections of the pore is distinctively different. In the wax-like pore, the adsorption and surface diffusivity are relatively slow, while on the inorganic regions of the pore, as a consequence of the weak solid-fluid attraction, the diffusivity is rapid. Other than this, it is seen that there is very little synergy between the adsorbents. This observation is particularly useful for modeling purposes.

Fig. 14 Depressurization of the synthetic crude with composition described in Table 2 at 125 ºC exposed to a composite pore structure, formed by two pores; a kerogen-like 2 nm pore (gray) and a clay-like 2 nm pore (orange). Blue points are individual NPzzAT simulation results; the red symbol corresponds to the extrapolated value for the bubble point,


20

Page 21 of 28 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

calculated as the intercept of two fitted straight lines (not shown). Uncertainties in the simulations are of the size of the symbols. Snaphots of equilibrium configurations before and after the bulk bubble point are shown. It is important to point out that the phase change points described here are not true bubble points of the mixture due to the finite size of the simulation cells and the very sharp asymmetry of the adsorption. It is clearly seen in Figure 14 that the organic pore dominates the adsorption, even though it only accounts for half of the available pore volume. This, in itself is a particularly important result, as it is consistent with experimental observations and provides a clear explanation of why, in spite of shales having a large volume fraction of inorganic materials (clays, quartz, calcite, etc.) it is the organic regions which are responsible for the adsorption of hydrocarbons35,70. The use of clays as models for shales18,Error! Bookmark not defined.,14 will most likely severely underestimate the adsorption. 7. Conclusions and perspectives The present paper considerably improves the current understanding on phase behavior of shale reservoir fluids under confinement. A rigorous model for the molecular modelling of shales is developed, informed by experimentation on four representative samples of the shale cores. We showcase here only one set of realizations, with a fixed shale composition and pore size distribution and kerogen fraction. However, we have performed simulations with different oil compositions, and, as expected, the quantitative results vary ( i.e. the degree in which the light components desorb into the production fluid, etc.) However, the general trends shown in the manuscript are maintained: tight shale pores adsorb strongly the heavy ends in detriment of methane which is desorbed upon a pressure drop. The extent of this process is not only governed by the fluid composition, but also by the pore size distribution and kerogen proportion. An appropriate quantitative model would need to take these variables into consideration.

While more complex models may be devised for both the fluid and particularly for the solids16,32,33 , there is a large uncertainty about the appropriate structure and composition of the solid phase. It then becomes a question of the relaive impact that these details have in the actual description of the behavior of fluids in shale rocks. On top of this, there is the need to perform simulations of multi-component multi-phase systems within a reasonable time frame. The use of accurate coarse grained models for the fluid-fluid behavior coupled with physically sound descriptions of plausible pore geometries, as described herein, seems like a promising course of action. The kerogen pore represents an extreme case of a strongly adsorbent material. It is unlikely to be, on its own, an accurate representation of a shale, but the results have relevance as they correspond to the most extreme case possible. Upon decompression, all of the heavier components remain strongly adsorbed in the nanopore or at the surface corresponding to a thin film within the walls of the larger mesopores. The remaining free gas at pressures below the bubble point is very lean (rich in methane and pentane). This suggests that the kerogen regions in shales may be selectively removing the light


21

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 22 of 28

components, leaving the heavier compounds, if they exist, sequestered in the nanopores. More importantly, this selective process takes place even at pressures above the bubble point (in the one phase region) and is progressive, in the sense that as the pressure is decreased, a higher proportion of light gases is seen in the free fluid. The immediate consequence is the leaning of the production gas (even at pressures above the bulk BP) as is frequently seen in production. The mixed-pore scenario is very enlightening. It combines in the same simulation cell a strongly attractive pore with a very inert pore. We show how this inert pore has very little influence on the adsorption properties of the oil. In this case, by coexisting in the same simulation cell we have been able to observe a more realistic scenario, which might resemble more closely the existing downhole conditions of a reservoir. Two important results are found; firstly, there is no evidence of synergy from the presence of the combined system, the total quantity and composition of the adsorbed fluid is an average between that found for the two extreme cases (purely kerogen or purely clay). The large increase in the bubble pressure for carbon-like pores and the correspondingly rather modest increase for the silica-like wall is consistent with the results obtained from simulations of pure alkanes. Secondly; we did observe a significant difference in the mobility of the species within the pores and its quantification is the subject of a subsequent communication. The carbon-like pores adsorb the larger alkanes very strongly and a very slow diffusivity is seen as compared to that in the silica pores. More detailed studies of the effect of the solid-fluid interactions and the geometry of the pores in the overall mass transport are needed. While we have considered two coexisting homogeneous pores, however several questions arise from these results: How realistic is a representation of shale as an intermediate structure between a carbon and a silica pore? How important will the local heterogeneities be; i.e. will a true representation of the surface consist of side-by-side pores of different nature, a homogeneous pore of intermediate surface attraction or a pore with a heterogeneous composition? In any case, a crucial question that must be answered is the appropriateness of the fluid-solid interactions as the results are crucially dependent on the reliability of our estimates of the solid-fluid interactions. At this point we have assumed plausible extreme scenarios, but an improved description is needed. Possibly the simplest experiments that can provide key information are the measurements of the fluid contact angles on rock samples71. These measurements are not unequivocal, and the techniques used have an influence on the results. In spite of this, being able to have access to experimental results of this type allowed us to fine tune and have confidence in the accuracy and predictive capabilities of our molecular models. Supporting Information Nitrogen isotherms pore size distributions, SEM images and photographs of the Samples are provided. Details of the experimental setup and contact angle determination along with tables of the raw data are presented.


22

Page 23 of 28 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

Acknowledgements This work was generously supported by ConocoPhillips. The authors would like to acknowledge the contribution of Jeff Sheremata and José Antonio Torres to the first draft of this manuscript. Further enlightening discussions with Jamal Sandarusi (ConocoPhillips) are gratefully acknowledged. EAM and CH acknowledge support from the U.K. Engineering and Physical Sciences Research Council (EPSRC) research grants EP/E016340, EP/J014958 and EP/I018212 and from the Thomas Young Centre under grant TYC-101. Simulations described herein were performed using the facilities of the Imperial College High Performance Computing Service.


23

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 24 of 28

References 1

Siirola, J. J. The impact of shale gas in the chemical industry. AIChE J. 2014, 60, 810–819. 2 Striolo, A.; Cole, D.R. Understanding Shale Gas: Recent Progress and Remaining Challenges. Energy Fuels 2017, 31, 10300–10310. 3 Yethiraj, A.; Striolo, A. Fracking: What Can Physical Chemistry Offer? J. Phys. Chem. Lett. 2013, 4, 687–690. 4 Clarkson, C.R.; Solano, N.; Bustin, R.M.; Bustin, A.M.M.; Chalmers, G.R.L.; He, L.; Melnichenko, Y.B.; Radliński, A.P.; Blach, T.P. Pore structure characterization of North American shale gas reservoirs using USANS/SANS, gas adsorption, and mercury intrusion. Fuel, 2013, 103, 606-616. 5 Patwardhan, S. D.; Famoori, F.; Gunaji, R.G.; Govindarajan, S.K. Simulation and Mathematical Modeling of Stimulated Shale Gas Reservoirs. Ind. Eng. Chem. Res. 2014, 53, 19788–805. 6 Singh, S.K.; Sinha, A.; Deo, G.; Singh, J.K. Vapor−Liquid Phase Coexistence, Critical Properties, and Surface Tension of Confined Alkanes. J. Phys. Chem 2009, 113, 7170–7180. 7 Zarragoicoechea, G. J.; Kuz, V.A. Critical shift of a confined fluid in a nanopore. Fluid Phase Equilib. 2004, 220, 7–9. 8 Dhanapal, K.; Devegowda, D.; Zhang, Y.; Contreras-Nino, A. C.; Civan, F.; Sigal, R. Phase Behavior and Storage in Organic Shale Nanopores: Modeling of Multicomponent Hydrocarbons in Connected Pore Systems and Implications for Fluids-in-place Estimates in Shale Oil and Gas Reservoirs SPE-169008-MS . 2014 SPE Unconventional Resources Conference; The Woodlands, Texas USA, 1–3 April 2014. 9 Tan, S.P.; Piri, M. Equation-of-state modeling of confined-fluid phase equilibria in nanopores. Fluid Phase Equilib. 2015, 393, 48–63. 10 Travalloni, L., Castier, M.; Tavares, F. W. Phase equilibrium of fluids confined in porous media from an extended Peng–Robinson equation of state. Fluid Phase Equilib. 2014, 362, 335–341. 11 Liu, J.; Wang, L.; Xi, S.; Asthagiri, D.; Chapman W.G. Adsorption and Phase Behavior of Pure/Mixed Alkanes in Nanoslit Graphite Pores: An iSAFT Application. Langmuir 2017 3311189–202. 12 Evans, R. Fluids adsorbed in narrow pores: phase equilibria and structure. J. Phys. Conden. Matt. 1990, 2, 8989–9007. 13 Gelb, L.; Gubbins, K.; Radhakrishnan, R.; Sliwinska-Bartowiak, M. Phase separation in confined systems. Rep. Prog. Phys. 1999, 62, 1573-1659 14 Sharma, A.; Namsani, S.; Singh, J. K. Molecular simulation of shale gas adsorption and diffusion in inorganic nanopores. Molecular Simulation, 2014, 41, 414–422. 15 Jin, Z.; Firoozabadi, A. Methane and carbon dioxide adsorption in clay-like slit pores by Monte Carlo simulations. Fluid Phase Equilibria 2013, 360, 456-465. 16 Collell, J.; Galliero, G.; Gouth, F.; Montel, F.; Pujol, M.; Ungerer, P.; Yiannourakou, M. Molecular simulation and modelisation of methane/ethane mixtures adsorption onto a microporous molecular model of kerogen under typical reservoir conditions. Micropor. Mesopor. Mat. 2014, 197, 194-203.


24

Page 25 of 28 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

17

Firouzi, M.; Alnoaimi, K.; Kovscek, A.; Wilcox, J. Molecular simulation and experimental characterization of the nanoporous structures of coal and gas shale. Int. J. Coal Geol. 2014, 123, 62–68. 18 Zhai, Z.; Wang, X.; Jin, X.; Sun, L.; Li, J.; Cao, D. Adsorption and Diffusion of Shale Gas Reservoirs in Modeled Clay Minerals at Different Geological Depths. Energy Fuels, 2014, 28, 7467–7473 19 Mosher, K.; He, J.; Liu, Y.; Rupp, E.; Wilcox, J. Molecular simulation of methane adsorption in micro- and mesoporous carbons with applications to coal and gas shale systems. Int. J. Coal Geol. 2013, 109, 36–44. 20 Tong, T.; Cao, D. A mesoscale model for diffusion and permeation of shale gas at geological depth. AIChE J 2018, 64, 1059–1066. 21 Psarras, P.; Holmes, R.; Vishal, V.; Wilcox J. Methane and CO2 Adsorption Capacities of Kerogen in the Eagle Ford Shale from Molecular Simulation Accounts Chem. Res. 2017, 50, 1818–1828. 22 Pathak, M.; Cho, H.; Deo, M. Experimental and Molecular Modeling Study of Bubble Points of Hydrocarbon Mixtures in Nanoporous Media. Energy Fuels 2017, 31, 3427–3435. 23 Jin Z.; Firoozabadi, A. Phase behavior and flow in shale nanopores from molecular simulations. Fluid Phase Equilib. 2016, 430, 156–168. 24 King Jr., H. E. ; Eberle, A. P. R.; Walters, C. C.; Kliewer, C. E.; Ertas, D.; Huynh, C. Pore Architecture and Connectivity in Gas Shale. Energy Fuels, 2015, 29, 1375– 1390. 25 Wüst, R.A.J., Nassichuk, B. R.; Bustin, R.M. “Porosity Characterization of Various Organic-Rich Shales From the Western Canadian Sedimentary Basin, Alberta and British Columbia, Canada.” In Electron Microscopy of shale hydrocarbon reservoirs: AAPG memoir 102, p. 81-100 Ed. W. Camp, E. Díaz, B. Wawak, Chapter 9. 2013. 26 Kuila, U.; Prasad, M. Specific surface area and pore-size distribution in clays and shales. Geophys. Prospecting, 2013, 61, 341–362 27 Xu, M.; Dehghanpour, H. Advances in Understanding Wettability of Gas Shales. Energy Fuels, 2014, 28, 4362-4375. 28 Tiwari, P.; Deo, M. Detailed kinetic analysis of oil shale pyrolysis TGA data. AIChE J., 2011, 58, 505–515, 29 Borysenko, A.; Clennell, B.; Sedev, R.; Burgar, I.; Ralston, J.; Raven, M.; Dewhurst, D.; Liu, K. Experimental investigations of the wettability of clays and shales. J. Geophysical Research, 2009, 114, B07202. 30 Zelenev, A.S., Surface Energy of North American Shales and its Role in Interaction of Shale with Surfactants and Microemulsions. Paper SPE 141459, presented at SPE International Symposium on Oilfield Chemistry, held in The Woodlands, TX, 11–13 April 2011. 31 Baojun, B.; Elgmati, M.; Zhang, H.; Wei, M. Rock characterization of Fayetteville shale gas plays. Fuel 2013, 105, 645–52. 32 Bousige, C.; Ghimbeu, C. M.; Vix-Guterl, C.; Pomerantz, A. E.; Suleimenova, A.; Vaughan, G.; Garbarino, G.; Feygenson, M.; Wildgruber, C.; Ulm, F.-J.; Pellenq, R. J. M.; Coasne, B. Realistic molecular model of kerogen’s nanostructure. Nat. Mater. 2016, 15, 576–582.


25

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 26 of 28

33

Vasileiadis, M.; Peristeras, L. D.; Papavasileiou, K. D.; Economou, I. G. Modeling of Bulk Kerogen Porosity: Methods for Control and Characterization. Energy Fuels 2017, 31, 6004-6018, 34 Michalec, L.; Lísal, M. Molecular simulation of shale gas adsorption onto overmature type II model kerogen with control microporosity. Mol. Phys. 2017, 115, 1–18. 35 Collell, J.; Ungerer, P.; Galliero, G.; Yiannourakou, M.; Montel, F., Pujol, M. Molecular Simulation of Bulk Organic Matter in Type II Shales in the Middle of the Oil Formation Window. Energy Fuels, 2014, 28, 7457–7466. 36 Lafitte, T.; Apostolakou, A.; Avendaño, C.; Galindo, A.; Adjiman, C. S.; Müller, E. A.; Jackson, G. Accurate statistical associating fluid theory for chain molecules formed from Mie segments. J. Chem. Phys. 2013, 139, 154-504. 37 Ramrattan, N. S.; Avendaño, C.; Müller, E. A.; Galindo, A. A corresponding-states framework for the description of the Mie family of intermolecular potentials. Mol. Phys. 2015, 113, 932–947. 38 Müller, E. A.; Jackson, G. Force-Field Parameters from the SAFT-γ Equation of State for Use in Coarse-Grained Molecular Simulations Annu. Rev. Chem. Biomolec. Eng. 2014, 5, 405-27. 39 Theodorakis, P. E.; Müller, E.A.; Craster, R.V.; Matar, O. K. Superspreading: Mechanisms and Molecular Design. Langmuir 2015, 31, 2304-2309. 40 Herdes, C.; Santiso, E. E.; James, C.; Eastoe, J.; Müller, E. A. Modelling the interfacial behaviour of dilute light-switching surfactant solutions. J. Colloid Interf. Sci., 2015, 445, 16–23. 41 Mejía, A.; Cartes, M.; Segura, H.; Müller, E. A. Use of Equations of State and Coarse Grained Simulations to Complement Experiments: Describing the Interfacial Properties of Carbon Dioxide + Decane and Carbon Dioxide + Eicosane Mixtures J. Chem. Eng. Data, 2014, 59, 2928–2941. 42 Müller, E. A.; Mejía, A. Resolving Discrepancies in the Measurements of the Interfacial Tension for the CO2 + H2O Mixture by Computer Simulation J. Phys. Chem. Lett., 2014, 5, 1267–1271. 43 Jiménez-Serratos, G.; Herdes, C.; Haslam, A.J.; Jackson, G.; Müller, E.A. Group Contribution Coarse-Grained Molecular Simulations of Polystyrene Melts and Polystyrene Solutions in Alkanes Using the SAFT-γ Force Field. Macromolecules 2017, 50, 4840–53. 44 Reiter, A. M.; Wallek, T.; Mair-Zelenka, P.; Siebenhofer, M.; Reinberger, P. Characterization of Crude Oil by Real Component Surrogates. Energy Fuels 2014, 28, 5565–5571. 45 Mejía, A.; Herdes, C.; Müller, E. A. Force Fields for Coarse-Grained Molecular Simulations from a Corresponding States Correlation. Ind. Eng. Chem. Res. 2014, 53, 4131–4141. 46 Herdes, C.; Totton, T.; Müller E.A. Coarse grained force field for the molecular simulation of natural gases and condensates. Fluid Phase Equilibria 2015, 406, 91100 47 Ervik, Å.; Mejía, A.; Müller E.A. Bottled SAFT: A Web App Providing SAFT-γ Mie Force Field Parameters for Thousands of Molecular Fluids J. Chem. Inf. Model. 2016, 56, 1609-1614.


26

Page 27 of 28 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

48

Ervik, Å.; Jiménez Serratos, G.; Müller, E.A. raaSAFT: A framework enabling coarse-grained molecular dynamics simulations based on the SAFT-γ Mie force field Comp. Phys. Comm. 2017, 212, 161–79. 49 Lobanova, O.; Avendaño, C.; Lafitte, T.; Müller, E. A.; Jackson, G. SAFT-γ force field for the simulation of molecular fluids: 4. A single-site coarse-grained model of water applicable over a wide temperature range. Molecular Physics, 2015 113, 1228– 1249. 50 Avendaño, C.; Lafitte, T.; Galindo, A.; Adjiman, C. S.; Jackson, G.; Müller, E. A. SAFT-γ Force Field for the Simulation of Molecular Fluids. 1. A Single-Site Coarse Grained Model of Carbon Dioxide. J. Phys. Chem. B 2011 115, 11154–11169. 51 E.W. Lemmon, M.O. McLinden, D.G. Friend, Thermophysical properties of fluid systems, in: P.J. Linstrom, W.G. Mallard (Eds.), NIST Chemistry WebBook, NIST Standard Reference Database Number 69, National Institute of Standards and Technology, Gaithersburg, MD 20899, 2015, ⟨http://webbook.nist.gov⟩ (retrieved August 1, 2015). 52 Santiso, E.; Herdes, C.; Müller, E.A. On the Calculation of Solid-Fluid Contact Angles from Molecular Dynamics. Entropy 2013, 15, 3734–3745. 53 Herdes, C.; Ervik, Å.; Mejía, A.; Müller, E.A. Fluid Phase Equilib. 2017, doi:10.1016/j.fluid.2017.06.016. 54 Theodorakis, P. E.; Müller, E.A.; Craster, R.V.; Matar, O.K. Modelling the superspreading of surfactant-laden droplets with computer simulation. Soft Matter 2015, 11, 9254–61. 55 Falk, K., Pellenq, R.; Ulm, F.-J.; Coasne, B. 2015, Effect of Chain Length and Pore Accessibility on Alkane Adsorption in Kerogen. Energy Fuels, 29, 7889–7896. 56 Lin, K.; Yuan, Q. ; Zhao, Y.-P. Using graphene to simplify the adsorption of methane on shale in MD simulations. Comp. Mater. Sci. 2017, 133, 99–107. 57 Jin, Z. Bubble/Dew Point and Hysteresis of Hydrocarbons in Nanopores From Molecular Perspective Fluid Phase Equilibria 2018, 458 177–85. 58 Forte, E.; Haslam, A. J.; Jackson, G.; E.A. Müller. Effective coarse-grained solidfluid potentials and their application to model adsorption of fluids on heterogeneous surfaces Phys. Chem. Chem. Phys. 2014, 16, 19165–19180. 59 Curtis, M. E.; Ambrose, R. J.; Sondergeld, C. H.; Rai, C. S. Structural Characterization of Gas Shales on the Micro-and Nano-Scales. Paper SPE 137693 Presented at the Canadian Unconventional Resources & International Petroleum Conference; Calgary, Alberta, Canada, 19–21 October 2010 2010. CUSG/SPE 137693 60 Bai, B., Elgmati, M., Zhang, H., Wei, M. Rock characterization of Fayetteville shale gas plays. Fuel, 2013, 105, 645–652. 61 Yuan, W. ; Pan, Z. ; Li, X. ; Yang, Y. ; Zhao, C. ; Connell, L. D. ; et al. Experimental study and modelling of methane adsorption and diffusion in shale. Fuel, 2014, 117, 509–519 62 Ambrose, R. J.; Hartman, R. C.; Diaz-Campos, M.; Akkutlu, I. Y.; Sondergeld, C. H. New Pore-scale Considerations for Shale Gas in Place Calculations. 2010 SPE Unconventional Gas Conference, 23-25 February, Pittsburgh, Pennsylvania, USA paper SPE 131772 .


27

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 28 of 28

63

Vlugt, T.J.H.; Krishna, R.; Smith, B. Molecular Simulations of Adsorption Isotherms for Linear and Branched Alkanes and Their Mixtures in Silicalite. J. Phys. Chem. B 1999, 103, 1102-1118 64 Nath, S.K.; de Pablo, J.J. Simulation of vapour-liquid equilibria for branched alkanes. Molecular Phys. 2000, 98, 231-238. 65 Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J.C.; Hess, B. Erik Lindahl, E. GROMACS: High Performance Molecular Simulations Through MultiLevel Parallelism From Laptops to Supercomputers. SoftwareX 1 2015, 1, 19–25. 66 Humphrey, W.; Dalke, A.; Schulten, K. VMD: visual molecular dynamics. J. Molec. Graphics, 1996, 14, 33-38. 67 Barreda, D.; Pérez-Mas, A. M.; Silvestre-Albero, A.; Casco, M. E.; Rudić, S.; Herdes, C.; Müller, E. A.; Blanco, C.; Santamaria, R.; Silvestre-Albero, J.; RodriguezReinoso, F. Unusual flexibility of mesophase pitch-derived carbon materials: An approach to the synthesis of graphene. Carbon, 2017, 115, 539-545. 68 Luo, S.; Nasrabadi, H.; Lutkenhaus, J. L. Effect of confinement on the bubble points of hydrocarbons in nanoporous media AIChE J. 2016, 62, 1772–1780 69 Gubbins, K. E.; Long, Y.; Śliwinska-Bartkowiak, M. Thermodynamics of confined nano-phases. J. Chem. Thermodyn. 2014, 74, 169–183. 70 Wang, S.; Javadpour, F.; Feng, Q. Fast mass transport of oil and supercritical carbon dioxide through organic nanopores in shale. Fuel 2016, 181, 741–758. 71 Mirchi, V.; Saraji, S.; Goual, L.; Piri, M. Dynamic interfacial tension and wettability of shale in the presence of surfactants at reservoir conditions. Fuel 2015, 148, 127–138.


28

Combined Experimental, Theoretical, and Molecular Simulation

Recommend Documents