Transport of Multicomponent Hydrocarbon Mixtures in Shale Organic

Sep 11, 2015 - Novel molecular simulation process design of adsorption in realistic shale kerogen spherical pores. Bo Zhou , Ruina Xu , Peixue Jiang.M...
0 downloads 0 Views 2MB Size
Subscriber access provided by UNIV OF NEBRASKA - LINCOLN

Article

Transport of Multicomponent Hydrocarbon Mixtures in Shales Organic Matter by Molecular Simulations. Julien Collell, Guillaume Galliero, Romain Vermorel, Philippe Ungerer, Marianna Yiannourakou, Francois Montel, and Magali Pujol J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.5b07242 • Publication Date (Web): 11 Sep 2015 Downloaded from http://pubs.acs.org on September 18, 2015

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 free 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 accessible to all readers and 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.

The Journal of Physical Chemistry C 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 16

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 Journal of Physical Chemistry

Transport of Multicomponent Hydrocarbon Mixtures in Shales Organic Matter by Molecular Simulations. Julien Collell,† Guillaume Galliero,∗,† Romain Vermorel,† Philippe Ungerer,‡ Marianna Yiannourakou,‡ Fran¸cois Montel,¶ and Magali Pujol¶ Laboratoire des Fluides Complexes et leurs R´eservoirs, Materials Design S.A.R.L., and TOTAL S.A. E-mail: [email protected]

Abstract

means of molecular simulations. We performed Molecular Dynamics simulations of hydrocarbons permeating through a molecular model representative of oil-prone type II kerogen. Our results show that the permeation mechanisms through this type of material is purely diffusive. Consequently, we have computed the Onsager’s species-specific transport coefficients of a typical condensate-rich gas mixture within kerogen. Interestingly, we have observed that the transport coefficients matrix can be reasonably approximated by its diagonal terms, the so-called Onsager’s auto-correlation coefficients. Inspired by the classical Rouse model of polymer dynamics and surface diffusion theory, we propose a simple scaling law to predict the transport coefficient of linear alkanes in the mixture. In good agreement with our simulations results, the Onsager’s auto-correlation coefficients scale linearly with the adsorption loading and inversely with the alkane chain length. We believe our results and predictions are applicable to other materials, such as carbon-based synthetic microporous membranes, with structural properties close to that of kerogen.

During the last decade, gas recovered from shale reservoirs has jumped from 2 to 40 percent of U.S. natural gas production. However, in response to the drop of gas prices, the oil and gas industry has set its sights on the oil-prone shale plays, potentially more lucrative. This shift from dry to condensate-rich gas has raised the needs for a better understanding of the transport of hydrocarbons mixtures through organic-rich shale reservoirs. At the micrometer scale, hydrocarbons in shales are mostly located in amorphous microporous nodules of organic matter, the so-called kerogen, dispersed in an heterogeneous mineral matrix. In such multiscale materials, a wide range of physical mechanisms might affect the composition of the recovered hydrocarbon mixtures. More specifically, kerogen nodules are likely to act as selective barriers due to their amorphous microporous structure. In this work, we study the transport of hydrocarbons mixtures through kerogen by ∗

To whom correspondence should be addressed LFC-R, UMR-5150 with CNRS and TOTAL, Universit´e de Pau et des Pays de l’Adour, BP 1155, 64013 Pau, France ‡ Materials Design S.A.R.L., 18 rue de Saisset, 92120 Montrouge, France ¶ TOTAL S.A., CSTJF, Avenue Larribau - 64018 Pau †

ACS Paragon Plus Environment

1

The Journal of Physical Chemistry

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

Introduction

Page 2 of 16

sub-nanometer pore space of the kerogen are expected to greatly impact the species-specific transport properties, thus leading to enhanced separation during recovery. From the experimental point of view, even if some recent studies have been carried out in that direction, 11 it is still extremely difficult to shed light on the dynamics of complex fluid mixtures in the kerogen phase, whose mass fraction only amounts to a few percents of the heterogeneous shale material. When experiments are too difficult to perform, Molecular Dynamics (MD) stands as a well-suited tool to take over the investigations. Recent studies have indeed addressed the rich phenomenology related to pure fluid permeation through molecular models representative of kerogen. Botan et al. performed Non Equilibrium Molecular Dynamics simulations of supercritical methane flowing through disordered microporous carbons membranes. 4 Their work demonstrates the predominance of sieving mechanisms when the fluid molecules diameter is comparable to the average pore size. Following a similar approach, Falk et al. studied the permeation of several pure n-alkane fluids through kerogen-like porous structures. 6 They showed that permeation rates scale like the inverse of the n-alkane chain length, which exemplifies the importance of separation mechanisms involved at the scale of the sub-nanometer porous space. Interestingly, these results are consistent, at least qualitatively, with classical permeation and diffusion experiments. 12–15 MD simulations are therefore capable of mimicking experiments provided the molecular models, including that of the microporous solids, are representative of their real counterparts. In this work, we address the transport of multicomponent mixtures in kerogen by means of MD simulations in the context of hydrocarbon production from shale reservoirs. In all our simulations, we have used a realistic molecular model of oil-prone type II kerogen described in previous studies. 16,17 The diffusional nature of mass transfer in shales organic matter has been clearly identified. Based on the results obtained from MD simulations, we propose simple scaling laws that account for the influence of the

The transport of fluids within disordered microporous carbons is inherent in numerous applications such as gas separation, 1 energy storage 2 or hydrocarbon production from unconventional petroleum reservoirs. A broad variety of carbon based materials are therefore involved in the aforementioned processes, ranging from synthetic microporous carbons 2–4 to natural materials such as coal 5 or the so-called kerogen contained in organic rich shales. 6 Because of the confinement of the fluid species within such narrow and complex porous structures, the investigation of fluid mixtures transport properties in amorphous nanoporous carbons remains challenging. More specifically, this topic has been of growing interest in the context of hydrocarbon production from shale reservoirs. Nowadays about 40 percent of US natural gas comes from gas shales, whereas it was only 2 percent in 2000. 7 In addition, between 2010 and 2013 the oil production from shales has been multiplied by a factor 3 in the US, which reveals a recent trend in the local petroleum industry now leaning towards the more challenging, yet potentially more lucrative, oil-prone shales plays. Consequently, this move from dry methane to condensate-rich gas shale reservoirs has raised the needs for a proper description of fluid mixtures transport in this type of materials. At the micrometer scale, hydrocarbons in shales are mostly located in amorphous microporous nodules of organic matter, the kerogen, where they were generated over geological times. 8 However, the composition of the steam gas recovered from shale wells usually differs from that of the fluid stored in the organic nodules, with a significant enrichment in lighter components that is detrimental to the profitability of the producing wells. Because of the multi-scale nature of the porous space in shale reservoirs, many possible phenomena are likely to influence the composition of the hydrocarbons mixture on its way from the kerogen nodules to the well head. Among this wide variety of mechanisms, adsorption selectivity 9 and molecular sieving effects 10 occurring in the

ACS Paragon Plus Environment

2

Page 3 of 16

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 Journal of Physical Chemistry

where J~i is the diffusive flux of species i in a mixture of n species, cj is the concentration of ~ j is the chemithe mixture in species j and ∇µ cal potential gradient of species j. Lij = Λij /cj where Λij are the elements of the [Λ] matrix of phenomenological Onsager’s coefficients. 20,21 The elements Lij of the matrix [L] can be determined from EMD simulations using Einstein’s relations 22 :

n-alkane chain length and fluid loading on the species-specific diffusion coefficients. This work is organized in two parts : • The first part is devoted to the study of pure methane mass transfer under non equilibrium conditions. We report results from Boundary Driven Non Equilibrium Molecular Dynamics simulations of gas permeation through a slab of kerogen, in which no assumption is made on the nature of mass transfer mechanisms. Upon comparison with results from classical Equilibrium Molecular Dynamics simulations, we show that equilibrium Onsager’s coefficients are the relevant transport coefficients to describe out of equilibrium mass transfer in kerogen host structures.

1 lim Lij = 6Nj ∆t ∆t→∞ * N ! i X (rl,i (t + ∆t) − rl,i (t))

+  l=1 Nj X (rk,j (t + ∆t) − rk,j (t)) · k=1

where rl,i (t) is the position of the center of mass of molecule l = 1, . . . , Ni of species i at time t, rk,j (t) is the position of the center of mass of molecule k = 1, . . . , Nj of species j at time t and h. . . i is the average over multiple time origins. The [L] matrix is thereafter referred to as the Onsager’s matrix, whose diagonal coefficients are referred to as auto-correlation coefficients and the out of diagonal (cross) terms are referred to as cross-correlation coefficients. For pure compounds, the matrix of Onsager’s coefficients simply reduces to a scalar value Lii = L, given by :

• In the second part, we investigate the diffusion properties of hydrocarbon mixtures within kerogen. For linear alkanes compounds, Onsager’s cross-correlation coefficients are found to be negligible compared to the Onsager’s auto-correlation coefficients of the individual compounds. By extending classical scaling laws to the case of mixtures, we propose a general scaling of Onsager’s auto-correlation coefficients which accounts for the n-alkane chain length and the mixture loading.

1 lim L= 6N ∆t ∆t→∞

Theoretical approaches of multicomponent diffusion In the classical decomposition of mass transfer mechanisms, 18 diffusion processes describe the species specific fluxes that deviate from the average net flux. Diffusion processes can be characterized in the frame of the irreversible process thermodynamics 19 which offers a comprehensive and theoretically grounded description of diffusion. Under isothermal approximation, the diffusive fluxes are given by : J~i = −

n X j=1

Lij

cj ~ ∇µj RT

(2)

*

N X l=1

(rl (t + ∆t) − rl (t))

!2 + (3)

Classically, the Onsager’s auto-correlation coefficients Lii of species i can be decomposed as : Lii = Dsi + Cii∗ (4) with : *N + i X 1 Dsi = lim (rl,i (t + ∆t) − rl,i (t))2 6Ni ∆t ∆t→∞ l=1 (5)

(1) ACS Paragon Plus Environment

3

The Journal of Physical Chemistry

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

and * N ! i X 1 Cii∗ = (rl,i (t + ∆t) − rl,i (t)) lim 6Ni ∆t ∆t→∞ l=1   + N i X  · (rk,i (t + ∆t) − rk,i (t)) k=1 k6=l

(6)

where Dsi is the self diffusivity, i.e., the diffusion coefficients of a single tagged molecule of species i within the mixture, and Cii∗ is a cross interaction coefficient, which accounts for correlations between two different molecules k and l of the same species i. In the limit of low concentrations, Cii∗ tends to zero and the self diffusivities equal the Onsager’s auto-correlation coefficients. 4,23

Simulations Details Molecular Dynamics MD simulations have been performed using the MedeA/LAMMPS 24 software, with the all atom pcff+ force field, developed by the team of Materials Design. 25 It is based on the class II Polymer Consistent Force Field (pcff ), 26 with refined non-bonded parameters. This force field explicitly accounts for non-bonded interactions (Lennard Jones 9-6 potential 27 with electrostatic interactions) as well as intramolecular interactions (angles, dihedrals, impropers as well as cross-terms). Interactions between unlike species were described using the W¨aldmann and Hagler combining rules. 28 A cutoff radius of 12 ˚ A was applied to Lennard Jones and electrostatic interactions. Temperature was imposed at 300 K with a Nos´e-Hoover thermostat. 29 The time step used for all MD simulations was set to 1 fs. In this work we have employed two different MD simulations schemes: Equilibrium Molecular Dynamics (EMD) and Boundary Driven Non-Equilibrium Molecular Dynamics (BD-NEMD). EMD simulations were run to compute the Onsager’s coefficients on bulk kerogen. Peri-

odic boundary conditions were applied on the system according to the three dimensions. Diffusion coefficients are classically computed from the study of particle trajectories, and are usually determined in the microcanonical (NVE) ensemble. 30,31 However, due to the cutoff of unbounded interactions, the internal energy was not completely conserved during the simulation runs. This resulted in constant increase of the temperature during NVE simulations. To avoid this issue, the system was thermostated, with a damping time constant τ of 1000 fs, supposed to negligibly impact the dynamics of the system. The position of molecules center of mass, involved in equation 2, were exported every tdump time interval along simulations. The Onsager’s coefficients were then calculated in a post-treatment procedure by performing running averages with several time origins taken each tdump time interval. In order to increase the computational efficiency of our post-treatment, we set tdump to the maximal value required to obtain converged results. EMD simulations were performed for pure compounds over 50 millions of time steps and 150 millions of time steps for mixtures. Details of the EMD simulations are given in table 1. BD-NEMD simulations 32 have been run to study pure methane mass transfer. This algorithm aims to mimic gas permeation experiments through a porous membrane. Similarly to Dual Control Volume - Grand Canonical Molecular Dynamics (DCV-GCMD) simulations, 4 a flow is generated by connecting the porous membrane under study to two reservoirs of fluids at different densities. Yet, BD-NEMD is purely deterministic and does not require any Monte Carlo steps of particle insertion / deletion to generate the density gradient. A snapshot of the system simulated is given in figure 1. The porous structure was immersed in two fluid reservoirs (upstream and downstream reservoirs in figure 1) connected by a perturbation region (thin fluid slab) through periodic boundary conditions. The principle of the BDNEMD algorithm relies on the application of an external field F~ex to fluid particles located in the perturbation region to put the system out of equilibrium. It is worth noticing that the same

ACS Paragon Plus Environment

4

Page 4 of 16

Page 5 of 16

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 Journal of Physical Chemistry

Table 1: Details of simulations performed to determine transport properties in molecular models. T is the temperature, t is the simulated time length, τ is the damping time constant of the Nos´ e-Hoover thermostat and tdump is the time interval at which the position of particles were exported to determine Onsager’s coefficient matrix [L]. Simulation

t (ns)

T [K]

τ (fs)

tdump (fs)

50 10 150

300 300 300

1000 100 1000

5 – 20

EMD pure methane BD-NEMD pure methane EMD mixtures external force F~ex was applied to each atom of methane molecules. At steady state, a density difference is created between the two reservoirs and thus a mass flux is generated within the porous membrane. This method is therefore free of any assumption about the mechanisms of mass transfer through the porous membrane. The heat flux introduced by the external field was removed by two independent thermostats applied on each reservoir, with a damping time constant τ of 100 fs, the temperature being measured on the two dimensions perpendicular to the mass flux. At steady state, we assume that the one dimensional mass flux Q across the porous structure is given by : Q=−

c ∆µ × χef f × RT e

absolute phenomenological coefficient for mass transfer χo was extrapolated at zero chemical potential gradient, i.e., when F~ex → ~0 : χo = lim χef f ∆µ→0

Figure 3(a) illustrates the determination of χo from our simulations data and clearly shows the linear relation between Q and ∆µ, which validates the use of equation 7.

Molecular models Hereafter we briefly describe the molecular models we used in our simulations for the kerogen and fluid phases. More detailed information on these models is available in previous papers. 16,17

(7)

where c is the methane concentration, χef f is the effective mass transfer coefficient linking the finite difference of chemical potential ∆µ/e across the kerogen membrane to the mass flux Q. The finite difference of chemical potential is defined as : µdown − µup ∆µ = e e

(9)

Kerogen model The organic matter in shales comes from the burial of the sedimentary organic matter. The increase of the temperature and the pressure over geological time led to a cracking of the sediments which generated the fluids (mostly hydrocarbons) and residues, mostly composed of kerogen. 33 At the molecular scale, kerogen is mostly made of polyaromatic clusters of carbon atoms crosslinked to each other by short alkyl groups. 34 As a natural material, it exhibits a non negligible content of heteroatoms (N,S,O) (see table 2). Furthermore, the molecular structure of kerogen depends on the type and the maturity of the sedimentary organic matter. The molecular model of kerogen used in this work aims to mimic a type II shale, typical

(8)

where µup is the chemical potential in the upstream reservoir, µdown is the chemical potential in the downstream reservoir and e is the thickness of the kerogen membrane. The chemical potential Chemical potentials in the two bulk reservoirs were determined from additional Monte Carlo simulations in the NVT ensemble, using the MedeA/GIBBS softwarea . The a Medea/Gibbs license from IFP-EN and Laboratory of Chemical Physics, CNRS-Universit´e Paris Sud

ACS Paragon Plus Environment

5

The Journal of Physical Chemistry

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

perturbation region

upstream reservoir

kerogen membrane

Page 6 of 16

downstream reservoir

External Field

Figure 1: Illustration of the simulation box used to perform BD-NEMD simulations. The external field was applied in the perturbation region, creating a density gradient between the high and low density regions. Periodic boundary conditions were applied. Table 2: Some characteristics of the kerogen molecule and structure used in this work. 16,17

of a marine depositional environment, whose maturity corresponds to the middle of the so called oil formation window, which refers to oilprone shales. Starting from molecular models that are consistent with the chemical structure and composition of kerogen molecules deduced from experiments, 35 we obtained the 3D structures by simulated annealing techniques based on MD under typical reservoir pressure, temperature and hydrocarbons loading conditions as described in our previous papers. 16,36 Upon changing the initial positions and velocities of the molecules in these MD simulations, we generated a set of 20 different 3D kerogen structures. As a preliminary step, we have investigated the self-diffusion of methane in each of these 20 kerogen structures. For the present study, we have selected the one kerogen structure that is the most representative of the average over all the generated structures in terms of diffusion properties. Some of its characteristics are provided in table 2.

Structure : H/C O/C N/C S/C Box length Atoms number Adsorption area Average density Porosity

0.905 0.054 0.021 0.008 28.9 ˚ A 1940 110 m2 /g 0.95 g/cm3 16 %

tane and tetradecane) and two compounds for the aromatics (toluene and dimethylnaphtalene). The composition of this mixture is given in table 3. We emphasize that no gas-liquid like phase transition were observed in our simulations, even for subcritical mixtures. Indeed, it is widely acknowledged that the confinement of the interstitial fluid in micropores results in a significant decrease of the critical temperature, 41–44 which explains that no capillary condensation occurs in microporous materials.

Condensate gas mixture model The fluid phase present in such shales corresponds to a typical condensate-rich gas. 17,37 Molecular models aimed to represent the asphaltene/resin fraction, 33,38,39 the hydrocarbon fraction, carbon dioxide fraction and water fraction were introduced. The hydrocarbon distribution is represented with 8 lumped compounds 37,40 : 6 compounds for the linear alkanes (methane, ethane, propane, butane, oc-

Porous structure generation The kerogen porous structure used in this work was obtained following a simulated annealing procedure. Simulation boxes containing 4 kerogen molecules and the fluid phase were first

ACS Paragon Plus Environment

6

Page 7 of 16

The Journal of Physical Chemistry

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 ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 ACS Paragon Plus Environment

Page 8 of 16

Page 9 of 16

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 Journal of Physical Chemistry

most species, the magnitude of Onsager crosscorrelation coefficients Lij are negligible compared to the Onsager’s auto-correlation coefficients Lii , thus the motion induced by the displacement of species j on species i is small compared to the self-induced motion of any species i. The picture is a little different for methane, which exhibits non vanishing Onsager crosscorrelation coefficients. Intuitively, we can deduce that lightweight methane molecules are therefore more sensible to the motion of heavier fluid species in the mixture. Nonetheless, cross-correlation coefficients involving methane are yet almost one order of magnitude lower than the methane auto-correlation coefficient. Thus, the flux in equation 1 arises mostly from the Onsager’s auto-correlation coefficients, and the full Onsager matrix can be approximated by its diagonal elements. In our EMD simulations we also computed the self diffusion coefficients of the alkane species in the condensate-rich gas mixture. Results are reported in figure 5, where we have plotted the mass transport coefficients multiplied by the alkane chain length for the different alkanes present in the mixture. With the exception of methane, it appears that self diffusivities are similar to Onsager’s auto-correlation coefficients, indicating that cross interactions Cii∗ are of low magnitudes. Regarding the methane, the cross interactions Cii∗ represent roughly the third of the Lii . At first, this observation seems to contradict the results previously obtained in the study of pure methane permeating through kerogen. However, the pure methane maximum loading investigated amounted to 6.41 mmol/g, while the total carbon groups loading amounted to 11.8 mmol/g in the case of the condensaterich gas mixture. This denser environment is thus likely to enhance cross-correlated motion of the most abundant species. Let us now introduce Ds,Nc the self diffusivity of a linear alkane containing Nc carbon atoms. As shown in figure 5, the self diffusivities of linear alkanes scale as the inverse of Nc within error bars. Thus, the self diffusivity of a given linear alkane with Nc CHx groups can roughly

pect that longer linear alkane molecules also permeate through kerogen via molecular diffusion. We indeed do not expect any transition to hydrodynamic regime upon increasing the alkane chain length. In the following of this paper we therefore rely on EMD simulations and assume that alkane mass transfer mechanisms in our kerogen model are diffusive only.

Mixture diffusion Hereafter, we present our results on the diffusion of mixtures within kerogen. Two different mixtures have been studied : the first mixture corresponds to the typical condensate-rich gas mixture described in a previous section. 17 The second mixture under study corresponds to a ternary mixture of methane, ethane and propane with a molecular ratio of 10/1/1 : for 10 molecules of methane inserted within the porous structure, 1 molecule of ethane and 1 molecule of propane are also incorporated. Onsager’s coefficients have been determined from EMD simulations on the condensate-rich gas mixture whose composition is given in table 3. Only the results for linear alkanes and carbon dioxide molecules are reported here. When multicomponent mixtures are involved, cross terms appear in the matrix of Onsager’s coefficients, characterizing interactions between unlike species. A three dimensional map of such coefficients is given in figure 4. For

Figure 4: Three dimensional map of the Onsager’s matrix coefficients Lij for the condensate-rich mixture.

ACS Paragon Plus Environment

9

The Journal of Physical Chemistry

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 ACS Paragon Plus Environment

Page 10 of 16

Page 11 of 16

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 Journal of Physical Chemistry

Table 4: Ncopt found for different linear alkane species in the three types of mixtures investigated : pure fluids, ternary mixture of methane, ethane and propane, the condensate-rich gas mixture. Species ethane (C2 ) propane (C3 ) n-butane (C4 ) n-octane (C8 ) n-tetradecane (C14 )

Pure compounds

Ncopt Ternary mixture

Condensate mixture

2.66 5.34 -

2.73 5.07 -

2.57 3.43 4.73 8.55 15.20

rescaled diffusivities for the ternary mixture of methane, ethane and propane (filled symbols in figure 6) follow the same linear trend with loading, and the data points for the mixture and the pure compounds superimpose within uncertainties. In addition, we report the rescaled diffusivity for the condensate-rich gas mixture (filled orange symbol in figure 6). Since the absolute composition of the condensate-rich gas mixture is fixed, the Ncopt scaling factors of the different compounds were determined at a single loading point. As a result, all the rescaled diffusivities equal that of methane, which is the reference compound, and are superimposed on a single data point. In the case of this complex mixture, the determination of the total loading is somewhat ambiguous as carbon-containing compounds other than alkanes, such as carbon dioxide or aromatics, may contribute significantly to the overall steric hindrance into the kerogen porous space. The wide uncertainties on the x-axis thus correspond to different ways of estimating the total loading:

and aromatic compounds (toluene and dimethylnaphtalene); Under those assumptions the rescaled diffusivity of the condensate-rich mixture aligns with the other data points. Based on these observations, we propose a simple scaling law to describe the diffusivity of alcanes in kerogen as follows: 0 Ds,1 × (1 − θc ) (13) Lii ≈ Nc where Lii is the Onsager’s auto-correlation co0 efficient of the species i, Ds,1 is the zero loading self diffusivity of methane, Nc is the number of carbon groups in the alkane i chain. θc is the fractional occupancy, defined as θc = νc /νc∗ , where νc is the total number of carbon groups adsorbed in the porous structure and νc∗ is the total number of carbon groups for which Lii = 0, obtained by linear extrapolation. The scaling proposed in equation 13 is reported in figure 6 (dashed black line). The zero loading self0 diffusivity of methane Ds,1 = 14.6 · 10−9 m2 /s has been obtained from a linear extrapolation of the methane self diffusivities presented in figure 3 and for each species Nc was set to its optimal value Ncopt . The proposed scaling is observed to be consistent with the general trend followed by the simulations results on pure compounds and mixtures over the entire range of loading investigated.

• the lower limit corresponds to the loading that only accounts for the carbon groups of the linear alkanes (CHx groups) ; • the symbol corresponds to the loading that accounts for the carbon groups of the linear alkanes and carbon dioxide molecules; • the upper limit corresponds to the loading that accounts for the carbon groups of the linear alkanes, carbon dioxide molecules

Conclusions In this work we have investigated the transport of hydrocarbons mixtures in a molecular

ACS Paragon Plus Environment

11

The Journal of Physical Chemistry

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

model representative of type II oil-prone kerogen. Molecular simulations were performed to study the impact of both chemical composition and adsorption loading of different fluid mixtures on the mechanisms of mass transfer within shale’s organic matter. Results from BD-NEMD simulations, which are free of any assumption about the mechanisms of permeation, indicate that the absolute mass transfer coefficient of pure methane is comparable to the Onsager’s coefficient computed under equilibrium conditions. This demonstrates that fluid transport inside kerogen is of diffusive nature and that transport coefficients can therefore be obtained directly by means of EMD simulations. EMD simulations were consequently performed for pure methane, ethane and propane, a ternary mixture of the latter compounds and a complex fluid mixture representative of condensate-rich gas found in oil-prone kerogen, which includes linear alkanes ranging from C1 to C14 . We have shown that for all compounds in the mixtures, the Onsager’s cross-correlation coefficients Lij are small compared to the autocorrelation coefficients Lii . As a result, the full matrix of Onsager’s coefficients can reasonably be approximated by its diagonal terms. For a mixture of n species, the n coupled equations given in equation 1 therefore reduce to a system of n independent equations, which drastically simplifies the treatment of multicomponent mixture diffusion in such materials. Furthermore, results from EMD simulations yield that the Onsager’s auto-correlation coefficient reduces to the self diffusivity for most compounds in the range of loading investigated. Transport coefficients can thus be well approximated by self diffusivities, which are accessible to NMR experiments 11 even though, to the best of our knowledge, no such experiment has yet been performed on isolated kerogen. In agreement with the predictions of surface diffusion theory, 10,52–54 we have observed that the Onsager’s auto-correlation coefficients L of pure methane, ethane and propane scale linearly with loading as L ∼ 1 − θc , where θc is the fractional occupancy of carbon groups inside kerogen. This scaling is also valid in the

Page 12 of 16

case of the ternary mixture of methane, ethane and propane. In addition, we have found that the Onsager’s auto-correlation coefficients Lii of the n-alkanes in the condensate-rich mixture scale like the inverse of the chain length as Lii ∼ 1/Nc , in the same way as the Rouse theory of diluted polymer diffusion. 14,15 This should result in a dynamic separation of alkanes with respect to their chain length, when migrating from the microporous kerogen nodules to larger-scale porous structures in the shale rock. Based on these observations, we have proposed a simple formula for the species specific Onsager’s coefficients Lii such that Lii ≈ 0 0 Ds,1 (1 − θc )/Nc . Where Ds,1 is the selfdiffusivity of methane at infinite dilution in kerogen, which can be easily computed from classical EMD simulations. For all our simulations data, the computed Onsager’s coefficients are consistent with the proposed scaling. We believe our results are applicable to the case of other carbon-based materials, such as manufactured microporous carbon membranes with molecular structures close to that of kerogen. In the context of oil and gas production from shale rocks, other factors are likely to significantly influence the transport and separation of hydrocarbons during recovery. For instance, the organic/inorganic interfaces between the kerogen and the surrounding minerals, 55 as well as the presence of water in the inorganic porous space might play an important role. Such problems will be addressed in future works. Acknowledgement The authors gratefully acknowledge TOTAL S.A. for financial support and for a grant for one of us (Julien Collell). Computer time for this study was provided by the computing facilities of TOTAL S.A., and of the MCIA (M´esocentre de Calcul Intensif Aquitain) of the Universit´e de Bordeaux and of the Universit´e de Pau et des Pays de l’Adour.

References 1. Liang, C.; Sha, G.; Guo, S. Carbon membrane for gas separation derived from coal tar pitch. Carbon 1999, 37, 1391–1397. 2. P´ean, C.; Merlet, C.; Rotenberg, B.; Mad-

ACS Paragon Plus Environment

12

Page 13 of 16

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 Journal of Physical Chemistry

den, P. A.; Taberna, P. L.; Daffos, B.; Salanne, M.; Simon, P. On the dynamics of charging in nanoporous carbon-based supercapacitors. ACS Nano 2014, 8, 1576–1583.

12. Rao, M. B.; Sircar, S. Nanoporous carbon membranes for separation of gas mixtures by selective surface flow. Journal of Membrane Science 1993, 85, 253–264.

3. Karger, J.; Pfeifer, H.; Vartapetian, R. S.; Voloshchuk, a. M. Molecular self-diffusion in active carbons. Pure and Applied Chemistry 1989, 61, 1875–1880.

13. Centeno, T.; Fuertes, A. Supported carbon molecular sieve membranes based on a phenolic resin. Journal of Membrane Science 1999, 160, 201–211.

4. Boan, A.; Vermorel, R.; Ulm, F.-J.; Pellenq, R. J.-M. Molecular simulations of supercritical fluid permeation through disordered microporous carbons. Langmuir : the ACS journal of surfaces and colloids 2013, 29, 9985–90.

14. de Gennes, P. G. Scaling concepts in polymer physics; 1979; p 324 p. 15. Rouse, P. E. A Theory of the Linear Viscoelastic Properties of Dilute Solutions of Coiling Polymers. The Journal of Chemical Physics 1953, 21, 1272.

5. Hu, H.; Li, X.; Fang, Z.; Wei, N.; Li, Q. Smallmolecule gas sorption and diffusion in coal: Molecular simulation. Energy 2010, 35, 2939– 2944.

16. Ungerer, P.; Collell, J.; Yiannourakou, M. Molecular Modeling of the Volumetric and Thermodynamic Properties of Kerogen : In fluence of Organic Type and Maturity. Energy & Fuels 2015,

6. Falk, K.; Coasne, B.; Pellenq, R.; Ulm, F.-j.; Bocquet, L. Subcontinuum mass transport of condensed hydrocar- bons in nanoporous media. Nature Communications 2015, 1–25.

17. Collell, J.; Galliero, G. Determination of the thermodynamic correction factor of fluids confined in nano-metric slit pores from molecular simulation. The Journal of chemical physics 2014, 140, 194702.

7. Cueto-Felgueroso, L.; Juanes, R. Forecasting long-term gas production from shale. Proceedings of the National Academy of Sciences 2013, 110, 19660–19661.

18. Bird, R. B.; Stewart, W. E.; Lightfoot, E. N. Applied Mechanics Reviews; 2007; Vol. 55; p R1.

8. Rexer, T.; Mathia, E.; Aplin, A.; Thomas, K. High-pressure Methane Adsorption and Characterization of Pores in Posidonia Shales and Isolated Kerogens. Energy & Fuels 2014, 28, 2886–2901.

19. Groot, S. D.; Mazur, P. Non-equilibrium thermodynamics; Dover publications, 1962. 20. Onsager, L. Reciprical relations in irreversible processes. I. Physical Review 1931, 37, 405– 426.

9. 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. Microporous and Mesoporous Materials 2014, 197, 194–203.

21. Onsager, L. Reciprical relations in irreversible processes. II. Physical Review 1931, 38, 2265– 2279. 22. Krishna, R.; Van Baten, J. M. Hydrogen bonding effects in adsorption of water-alcohol mixtures in zeolites and the consequences for the characteristics of the Maxwell-Stefan diffusivities. Langmuir 2010, 26, 10854–10867.

10. Krishna, R.; Baten, J. M. V. Unified MaxwellStefan description of binary mixture diffusion in micro- and meso-porous materials. Chemical Engineering Science 2009, 64, 3159–3178. 11. Korb, J.; Nicot, B. Dynamics and Wettability of Oil and Water in Oil Shales. The Journal of Physical Chemistry C 2014, 118, 23212– 23218.

23. Nicholson, D. The transport of adsorbate mixtures in porous materials: Basic equations for pores with simple geometry. Journal of Membrane Science 1997, 129, 209–219.

ACS Paragon Plus Environment

13

The Journal of Physical Chemistry

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 16

36. 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.

24. Plimpton, S. Fast Parallel Algorithms for Short Range Molecular Dynamics. Journal of Computational Physics 1995, 117, 1–19. 25. Yiannourakou, M.; Ungerer, P.; Leblanc, B.; Rozanska, X.; Saxe, P.; Vidal-Gilbert, S.; Gouth, F.; Montel, F. Molecular Simulation of Adsorption in Microporous Materials. Oil & Gas Science and Technology Revue dIFP Energies nouvelles 2013, 18.

37. Lagache, M.; Ungerer, P.; Boutin, a. Prediction of thermodynamic derivative properties of natural condensate gases at high pressure by Monte Carlo simulation. Fluid Phase Equilibria 2004, 220, 211–223.

26. Sun, H.; Mumby, S. J.; Maple, J. R.; Hagler, A. T. An ab Initio CFF93 All-Atom Force Field for Polycarbonates. 1994, 2978–2987.

38. Castex, H. R´esines et asphalt`enes : ´evolution en fonction des types de mati`ere organique et de leur enfouissement. Oil & Gas Science and Technology 1985, 40, 169–189.

27. Jones, J. E. On the Determination of Molecular Fields . II . From the Equation of State of a Gas. Proceedings of the Royal Society of London 1924, 4, 463–477.

39. Ungerer, P.; Rigby, D.; Leblanc, B.; Yiannourakou, M. Sensitivity of the aggregation behaviour of asphaltenes to molecular weight and structure using molecular dynamics. Molecular Simulation 2014, 40, 115–122.

28. Waldman, M.; Hagler, A. T. New combining rules for rare gas van der waals parameters. Journal of Computational Chemistry 1993, 14, 1077–1084.

40. Rokosh, C. D.; Lyster, S.; Anderson, S. D. A.; Beaton, A. P.; Berhane, H.; Brazzoni, T.; Chen, D.; Cheng, Y.; Mack, T.; Pana, C. et al. Open File Report : Summary of Alberta s Shale- and Siltstone-Hosted Hydrocarbon Resource Potential ; 2010; p 339.

29. Hoover, W. G. Canonical dynamics: Equilibrium phase-space distributions. Phys. Rev. A 1985, 31, 1695–1697. 30. Allen, M.; Tildesley, D. Computer simulation of liquids. 1989, 400. 31. Frenkel, D.; Smit, B. Understanding Molecular simulation - From Algorithms to Applications, academic p ed.; 1996; p 658.

41. Evans, R.; Marconi, U.; Tarazona, P. Capillary Condensation and Adsorption in Cylindrical and Slit-like Pores. Journal of the Chemical Society, . . . 1986, 1787, 1763–1787.

32. Frentrup, H.; Avenda˜ no, C.; Horsch, M. Transport diffusivities of fluids in nanopores by non-equilibrium molecular dynamics simulation. Molecular Simulation 2012, 38, 37–41.

42. Singh, S. K.; Sinha, A.; Deo, G.; Singh, J. K. Vapor - Liquid Phase Coexistence, Critical Properties, and Surface Tension of Confined Alkanes. 2009, 7170–7180.

33. Tissot, B. P.; Welte, D. H. In Eos, Transactions American Geophysical Union; Springer Berlin Heidelberg,, Ed.; 1984; p 720.

43. Tarazona, P.; Marconi, U.; Evans, R. Phase equilibria of fluid interfaces and confined fluids: non-local versus local density functionals. Molecular Physics 1987, 60, 573–595.

34. Vandenbroucke, M.; Largeau, C. Kerogen origin, evolution and structure. Organic Geochemistry 2007, 38, 719–833.

44. Didar, R.; Texas, A. Pore-Size Dependence of Fluid Phase Behavior and the Impact on Shale Gas Reserves. 2013, 1–22.

35. Kelemen, S. R.; Afeworki, M.; Gorbaty, M. L.; Sansone, M.; Kwiatek, P. J.; Walters, C. C.; Freund, H.; Siskin, M.; Bence, A. E.; Curry, D. J. et al. Direct Characterization of Kerogen by X-ray and Solid-State 13 C Nuclear Magnetic Resonance Methods. Energy & Fuels 2007, 21, 1548–1561.

45. Herrera, L.; Do, D. D.; Nicholson, D. A Monte Carlo integration method to determine accessible volume, accessible surface area and its fractal dimension. Journal of colloid and interface science 2010, 348, 529–36.

ACS Paragon Plus Environment

14

Page 15 of 16

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 Journal of Physical Chemistry

46. D¨ uren, T.; Millange, F.; F´erey, G. Calculating geometric surface areas as a characterization tool for metal-organic frameworks. The Journal of Physical Chemistry C 2007, 111, 15350– 15356. ´ 47. M¨oller, D.; Oprzynski, J.; M¨ uller, A.; Fischer, J. Prediction of thermodynamic properties of fluid mixtures by molecular dynamics simulations: methane-ethane. Molecular Physics 1992, 75, 363–378. 48. Brochard, L.; Vandamme, M.; Pellenq, R. J. M.; Fen-Chong, T. Adsorption-induced deformation of microporous materials: Coal swelling induced by CO 2-CH 4 competitive adsorption. Langmuir 2012, 28, 2659–2670. 49. Jain, S.; Pellenq, R.; Pikunic, J.; Gubbins, K. Molecular modeling of porous carbons using the hybrid reverse Monte Carlo method. Langmuir 2006, 22, 9942–9948. 50. Maginn, E. J.; Bell, A. T.; Theodorou, D. N. Transport diffusivity of methane in silicalite from equilibrium and nonequilibrium simulations. The Journal of Physical Chemistry 1993, 97, 4173–4181. 51. Chempath, S.; Krishna, R.; Snurr, R. Q. Nonequilibrium Molecular Dynamics Simulations of Diffusion of Binary Mixtures Containing Short n -Alkanes in Faujasite. Journal of Physical Chemistry B 2004, 108, 13481– 13491. 52. K¨arger, J.; Ruthven, D. M.; Theodorou, D. N. Diffusion in Nanoporous Materials, wiley-vch ed.; 2012. 53. Krishna, R. DIFFUSION OF ADSORBED SPECIES : A DESCRIPTION BASED ON THE GENERALIZED EQUATIONS. Chemical Engineering Science 1990, 45, 1779–1791. 54. Reed, D.; Ehrlich, G. Surface diffusion, atomic jump rates and thermodynamics. Surface Science Letters 1981, 102, 588–609. 55. Hantal, G.; Brochard, L.; Dias Soeiro Cordeiro, M. N.; Ulm, F. J.; Pellenq, R. J. M. Surface chemistry and atomic-scale reconstruction of kerogen-silica composites. Journal of Physical Chemistry C 2014, ACS Paragon Plus Environment

15

The Journal of Physical Chemistry

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 ACS Paragon Plus Environment

Page 16 of 16