Effect of Microstructural Flexibility on Methane Flow in Kerogen Matrix

Subscriber access provided by UNIV AUTONOMA DE COAHUILA UADEC

C: Surfaces, Interfaces, Porous Materials, and Catalysis

Effect of Microstructural Flexibility on Methane Flow in Kerogen Matrix by Molecular Dynamics Simulations Tianhao Wu, and Abbas Firoozabadi J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.8b12328 • Publication Date (Web): 05 Apr 2019 Downloaded from http://pubs.acs.org on April 5, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 22 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

Effect of Microstructural Flexibility on Methane Flow in Kerogen Matrix by Molecular Dynamics Simulations

Tianhao Wu† and Abbas Firoozabadi*,† † Reservoir

*

Engineering Research Institute, 595 Lytton Avenue Suite B, Palo Alto, CA 94301, USA.

Corresponding author; Tel: 650-326-9172; Fax: 650-472-9285; E-mail: [email protected]

1 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

Page 2 of 22

ABSTRACT We investigate the flow of methane in a kerogen matrix with microstructural flexibilities at various pressures. The matrix is constructed by compressing a collection of 60 type II kerogen macromolecules. In the past, simulations of methane flow in kerogen matrices have been performed assuming rigid molecular structures. We extend the simulations from rigid molecules to flexible molecules. The gas flow simulations are performed based on the boundary-driven method. We introduce a limited number of virtual “nails” to keep the flexible kerogen matrix in place. It is demonstrated that the flexibility of the kerogen microstructure has a significant effect on gas diffusion which is the primary transport mechanism in kerogen that contains micropores. The adsorption in flexible kerogen is higher than in rigid kerogen. Adsorption in flexible kerogen matrix may narrow the main flow path and reduce gas flow in the confined environment under subsurface conditions. On the other hand, the occasional opening of pore throats in the flexible kerogen matrix can provide additional flux. We find that the transport of methane in flexible kerogen is less than that in rigid kerogen due to changes in pore shape, despite the additional transport from pore opening in the flexible kerogen matrix.


Page 3 of 22 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 Shale gas is one of the most important unconventional resources. The matrix of the shale, where the gas resides, is composed of two distinct types of media, organic and inorganic matter, both of which contain nanoscale pores.1 Generally, kerogen is the predominant component of organic matter in most shale gas/oil formations. The complex chemical structure and the broad distribution of nanopores in kerogen affect flow due to pronounced solid-fluid interactions. The composition and the sedimentary process affect the properties of kerogen. Kerogen is usually classified based on maturity indicators, including the vitrinite reflectance, the aromatic-to-aliphatic ratio, and the hydrogen-to-carbon (H/C) and oxygen-to-carbon (O/C) atomic ratios, as described by the van Krevelen diagram. In general, the maturity of kerogen can qualitatively indicate the level of thermal degradation and the state in terms of hydrocarbon generation. Various representative structures of mature and immature kerogen have been constructed based on properties determined from experiments.2-3 The studies based on these kerogen macromolecules include mechanisms of adsorption,4-9 swelling,10-11 and self-diffusion.12-13 The flow of methane in organic-rich shale is a key element of unconventional resource development. The gas transport mechanism in the shale matrix has received significant attention.13-18 However, the investigation of gas transport in a realistic kerogen structure at a molecular scale is in a preliminary stage. Gas flow has been investigated based on graphite slit nanopores by using molecular simulations, which shows high gas flux and separation effect for gas mixtures.19-21 Kerogen matrix is much more complicated than graphite slit nanopores due to the complexities from molecular structure. Collell et al.15 investigated hydrocarbon flow in type II kerogen matrices with rigid microstructure. They simulated methane flow by using the boundary-driven non-equilibrium molecular dynamics (BD-NEMD) and calculated selfdiffusion of the hydrocarbon mixtures by using the equilibrium molecular dynamics (EMD) simulations. Their results show that the permeation mechanism is purely diffusive. The transport coefficients can be 3 ACS Paragon Plus Environment


Page 4 of 22

approximated by the Onsager’s autocorrelation coefficients, which is linear with adsorption loading and inversely proportional to the alkane chain length. Obliger et al.16 studied the hydrocarbon mixture transport in a mature kerogen matrix by using the external force method (EFM). They showed that the interactions between mixture components and between molecules of the same species are negligible. The permeance, which is proportional to the collective diffusivity, could be estimated from the self-diffusivity of a single component. The diffusivity is found to scale with the adsorption loading and the reciprocal of the alkane length. Obliger et al.13 extended the simulations to different kerogen structures to compute the self-diffusivity by using the EMD simulations, and expressed the transport diffusivity as a function of porosity and loading.16 Okamoto et al.22 simulated methane flow through parallel-plate channels with rough surfaces composed of type III kerogen. The width of the channel was varied from 2.9 nm to 34 nm. They argued that the slip-velocity models could not describe gas flow in the channels when the width is smaller than 6 nm. The kerogen structures by the authors above are kept fixed as rigid solid. However, kerogen structure may have pronounced interaction with fluid molecules, at which conditions the microstructure may change. The flexibility of the structure in the kerogen matrix has been demonstrated to have a significant effect on gas adsorption and kerogen swelling.10,17,23 The volume expansion appears to quadratically increase with the amount of gas adsorbed. The kerogen volume can expand to 5.4% and 11% upon CH4 and CO2 adsorption at 192 atm, respectively.10 Vasileiadis et al.12 studied the gas selfdiffusion, including pure components and a gas mixture in the type II flexible kerogen. They demonstrated that the diffusion across the kerogen is anisotropic and the maximum component of the diffusion coefficient correlates linearly with the limiting pore diameter. To the best of our knowledge, the simulation of gas flow under pressure gradient in kerogen matrix with the consideration of the flexibility of molecular structure has not been reported yet.


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


In this study, we perform a systematic investigation of kerogen matrix construction and gas flow simulations. The central objective of this work is to investigate the effect of the flexibility of the kerogen molecular structure on methane flow in a kerogen matrix. The remainder of this paper is organized into three sections. First, we present the models and the methods of kerogen matrix construction, equilibrium simulation, and gas flow simulation. Second, the results of transport diffusivity are covered, and the mechanisms are analyzed. Finally, conclusions are drawn.

MODELS AND METHODS Kerogen matrix construction

We apply molecular dynamics (MD) simulations using the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS).24 The simulations are run on workstations equipped with Intel E5-2690 v4 CPUs and NVIDIA Tesla P100 GPUs. The kerogen macromolecule unit is adapted from Ungerer et al.3 which reproduces elemental and functional data of kerogen including the H/C, O/C, and N/C ratios, average aromaticity, and average size of the aromatic unit. A type II immature kerogen macromolecule (C252H294O24N6S3) is used to construct the kerogen matrix. We construct a larger kerogen matrix than those in previous studies with side length dimensions ranging from 2.5 nm to 5.0 nm.3,13,15-16 We place 60 kerogen macromolecules at random positions in a cubic box with a side length of 20.0 nm. The periodic boundary condition is applied to all three directions. To create the pore space, we simultaneously and randomly insert 15 dummy atoms. The initial configuration is presented in Fig. 1. The all-atom model with a pcff+ force field is used for kerogen and methane.25 A cut-off of 14 Å is set for short-range interactions.

The

long-range

electrostatic

interactions

are

computed

using

the

Particle−Particle/Particle−Mesh (PPPM) method with an accuracy of 1×10-5.26-27 The Waldman–Hagler



Page 6 of 22

combining rules are used for the Lennard–Jones (LJ) interactions.28 The dummy atoms have an LJ diameter σ = 2.0 nm, energy parameter ε/kB = 10 K, and mass m = 1.0 g/mol. We perform a series of stepwise cooling and heating simulations to create a stable matrix with a skeletal density in agreement with the literature.3,23 The procedures for the MD simulations are listed in Table 1. In the first step, the conformation is relaxed in the NVT ensemble at 300 K. The pressure for the steps of the NPT ensemble is set as 200 bar to compress the structure. The temperature for the heating and cooling process ranges from 300 K to 1000 K. The temperature is finally increased from 300 K to 333.15 K, which is kept constant to represent the reservoir conditions for the remaining procedures in this study. The thermostat29-30 and barostat31 are achieved by coupling to the atom velocities and domain dimensions, respectively. The thermostat is applied to the translational degrees of freedom for the atoms. The barostat is isotropically coupled to the volume of the overall cubic simulation box without changing the shape. The damping factors are set to 100 fs and 1000 fs for the control of temperature and pressure, respectively. The equations of motion are integrated using the velocity-Verlet scheme.

Fig. 1. Kerogen macromolecule and initial configuration of the simulation box. (a) Type II kerogen macromolecule; (b) Initial configuration of 60 type II kerogen molecules and 15 dummy atoms.


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


Table 1. Procedure for creating the kerogen matrix. Ensemble

Temperature (K)

Pressure (bar)

Time (ps)

Time Step (fs)

Notes

NVT

300

-

100

1

-

NPT

300 to 1000

200

100

0.1

-

NPT

1000 to 300

200

100

0.1

-

NVT

300 to 333.15

-

100

1

-

1

After removing the dummy atoms

NVT

333.15

-

100

From the simulations, we construct the kerogen matrix with a size of 7.39 × 7.39 × 7.39 nm3 (see Fig. 2). The dummy atoms are removed subsequently, and the pore space is created. After the dummy atom removal, an extra thermostat simulation of the NVT ensemble is performed for 100 ps to obtain a relaxed kerogen model.

Fig. 2. Constructed kerogen matrix with dummy atoms. Domain: 7.39 × 7.39 × 7.39 nm3. The dummy atoms are represented with a diameter of 1.2 nm to avoid overlapping with the kerogen matrix in the adjacent region.

We perform Monte Carlo (MC) simulations by using a probe molecule to calculate the porosity and accessible surface area. The probe molecule has a diameter of a single spherical methane molecule (σ = 0.373 nm) in the TraPPE force field.32 The porosity is defined as the ratio of the accepted insertions (NA) to the total number of insertions (NT) in MC trials as





NA NT

Page 8 of 22

(1)

where NT is 1×106 in this study. The criteria for the acceptance of an MC trial is based on the distance (D) between the mass center of the probe molecule and the nearest atom in kerogen.33 The trial is accepted if

D

 p k

(2)

2

where σp and σk denote the LJ diameters of the probe molecule and the atoms in the kerogen matrix, respectively. The accessible surface area is estimated based on the method proposed by Düren et al.34 The principle of this method is based on rolling a probe molecule at the surface of every atom in the kerogen matrix, which follows the criterion in Eq. (2). The properties of the created kerogen matrix are listed in Table 2. The skeletal density is computed based on bulk density and porosity. The shape of the pore space is generated based on mapping with the probe molecule (see Fig. 3). The results show that the kerogen matrix has an interconnected pore throughout the x-direction, which has the potential to be the main path of gas flow. There are some isolated pores that may provide an additional path for flow in the flexible kerogen. Table 2. Physical properties of the kerogen matrix. Property

Result

Bulk Density

0.96 g/cm3

Skeletal Density

1.03 g/cm3

Porosity

7.0%

Surface Area

360 m2/g


Page 9 of 22 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


Fig. 3. Pore space in the kerogen matrix. There is an interconnected pore throughout the x-direction, which may serve as the main path of gas flow. (a) Pore space and kerogen matrix; (b) Extracted pore structure. Domain: 7.39 × 7.39 × 7.39 nm3.

Equilibrium simulations

To construct the simulation box with two reservoirs of gas molecules, we remove the periodic boundary condition in the x-direction and extend the simulation box. The smaller part of the kerogen molecules crossing the x-boundary is translocated manually to keep the individual kerogen molecules intact (see Fig. 4a and 4b). The pore structures at the x-direction boundaries in Fig. 4b are slightly changed compared to the structure shown in Fig. 4a and Fig. 3. To reduce the effect of the relatively unconsolidated regions near the x-boundaries of the kerogen, we use a region of 5.17 × 7.39 × 7.39 nm3 in the middle of the kerogen matrix covering approximately 70% of the entire matrix. This is the domain we use for the calculation of properties including pore size distribution, total uptake, and gas flux. We place methane molecules at random positions in the reservoirs on the two sides of the kerogen matrix. The size of the entire system is 40.00 × 7.39 × 7.39 nm3. We insert 1568, 3136, 9996, 16762, and 19440 methane molecules into the reservoirs and obtain five pressure conditions. Two sets of simulations are performed, in which the structures of the kerogen matrix are kept as rigid and flexible molecules, respectively. For



Page 10 of 22

the rigid kerogen, all atoms in the matrix are frozen in the simulation box. We perform the simulations in the NVT ensemble for 10 ns with a time step of 1 fs. The trajectories of the molecules in the system are recorded every 1 ps for the last 5 ns to calculate the properties. The pressures are determined based on the data from the two gas reservoirs. The final pressure in the gas reservoir is determined at equilibrium, which is listed in Table 3. One of the snapshots is presented in Fig. 4c.

Fig. 4. Sketch of the final geometries of kerogen construction and equilibrium simulation. (a) Kerogen matrix with periodic boundary condition in all three directions; (b) Translocation of small parts of kerogen to have full individual molecules at the boundary; (c) Simulation box at equilibrium with methane and kerogen matrix. The snapshot in (c) is taken from the simulation with flexible kerogen at the pressure of 63 bar. The configurations of the simulation box at other pressure conditions are similar to those in (c) but with different gas densities, which are not shown here. Domain in (c): 40.00 × 7.39 × 7.39 nm3, the figure is cropped to fit the scale.

Gas flow simulations

To simulate methane flow in the kerogen matrix under physical pressure differences, we implement the boundary driven non-equilibrium molecular dynamics (BD-NEMD) simulations.15 A region with a size of 5.00 × 7.39 × 7.39 nm3 is set as an accelerative zone (see Fig. 5). An external force (1×10-3 kcal⋅mol1⋅Å-1)

in the x-direction is applied to every molecule in the accelerative zone to generate a pressure


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


difference between the two reservoirs. The gas flow simulations are initiated from the final conditions of the equilibrium simulations. The simulations are performed at a constant temperature by using a thermostat, and the velocity in the x-direction is excluded from the temperature calculation. The simulations at each pressure are kept running until the gas flux through the kerogen matrix has reached a steady state. The time step is set to 1 fs. The total simulation time ranges from 15 ns to 25 ns. In the last 9 ns of each pressure point, the coordinates and velocities of the atoms are recorded every 1 ps for the calculations of different properties. The transport diffusivity is calculated every 3 ns, and the error bars are determined based on the standard deviation of the three intervals within the last 9 ns. Subsequently, the diffusion coefficients are computed based on gas flux:

Dt   J / c

(3)

where J is the flux, Dt is the transport diffusivity (diffusion coefficient), and c is the methane concentration. The concentration gradient between the two gas reservoirs is computed.

Fig. 5. Sketch of boundary-driven non-equilibrium molecular dynamics (BD-NEMD) simulations. Domain: 40.00 × 7.39 × 7.39 nm3.

In the rigid kerogen structure, the matrix is fixed in the middle of the simulation box. To account for the flexibility of the kerogen structure for gas flow simulations, we allow the atoms in the kerogen to move but keep the macroscopic position of the kerogen matrix in the middle. Otherwise, the kerogen matrix would be displaced by the gas pressure difference between the two reservoirs. We freeze the atoms in 8



Page 12 of 22

small regions in the matrix by using virtual “nails” to keep the kerogen in place (see Fig. 6). The size of each region is 0.4 × 0.4 × 0.4 nm3. There are 56 atoms in the nail regions in total, which is 0.16% of the total number of atoms in the kerogen matrix. The macroscopic properties are not affected by these fixed atoms. An animation for the demonstration of the methods is provided in Supporting Information Movie S1.

Fig. 6. The virtual nails to keep the kerogen matrix with a flexible structure in the middle of the simulation box. Domain: 40.00 × 7.39 × 7.39 nm3.

Three different schemes of simulations are performed for gas flow simulations. In Scheme 1, the kerogen structure is kept rigid from the start of the equilibrium simulations to the end of the gas flow simulations. In Scheme 2, the kerogen structure is flexible from the start of the equilibrium simulations to the end. In some previous studies, the kerogen structure is first simulated with flexibility in the stage of equilibrium simulations or is simulated based on other methods, whereas it is assumed to be rigid in gas flow simulations.15-16 Therefore, we perform Scheme 3, which starts from the final condition of the equilibrium simulations with a flexible kerogen structure and freeze the kerogen structure in the gas flow simulations (see Table 3). Five average pressure conditions with different pressure differences are obtained for each scheme.


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


Table 3. Schemes for kerogen structure in equilibrium and gas flow simulations. Equilibrium Simulation Scheme

Kerogen Flexibility

Gas Flow Simulation

Pressure (bar)

Kerogen Flexibility

Average Pressure (bar)

Pressure Difference (bar)

32

2

32 64 1

Rigid

64

5

210

15

430

430

18

565

565

20

32

32

2

210

Rigid

63 2

Flexible

63

5

207

15

419

419

18

549

549

20

32

32

2

207

Flexible

63 3

Flexible

63

5

207

15

419

419

18

549

549

20

207

Rigid

RESULTS There are various evidences in recent literature which suggests that adsorption and transport in kerogen matrices should be based on flexible kerogen molecules.10,21 In the following, the effect of kerogen molecular flexibility will be investigated on a consistent basis. The transport diffusivity at different pressures in the three schemes is presented in Fig. 7. The results are in the same range as the transport diffusivity in the literature: Collell et al.15 provide the transport diffusivity ranging from 0.84×10-8 to 1.06×10-8 m2/s in the rigid kerogen of a smaller matrix than in our work; Falk et al.18 present the transport diffusivity ranging from 2×10-8 to 5×10-8 m2/s in the rigid porous carbon matrix; Obliger et al.16 simulate the gas mixture transport in the rigid kerogen, at which conditions the diffusivity of methane in the mixture is in the range from 0.5×10-8 to 2×10-8 m2/s. Our transport diffusivity decreases as pressure increases, which is consistent with the trends reported in the literature.15,18 The pores in kerogen created by various authors are primarily micropores smaller than 2 nm. We set 0.8



Page 14 of 22

nm as the upper limit of the representative pore size to estimate diffusivity based on Knudsen diffusion.21,35 The results are 1.76×10-7 and 1.4×10-10 m2/s from Knudsen diffusion for a single tube and a porous medium, respectively. The Knudsen diffusion coefficients are different from the results of the molecular simulations. Knudsen diffusion cannot describe methane flow in a kerogen matrix with a complex pore structure.

Fig. 7. Transport diffusivity vs different average pressures for the three schemes at 333.15 K.

The transport diffusivity in Scheme 1 (rigid kerogen) is generally 2 times higher than that in Scheme 2 (flexible kerogen). The difference is attributed to the change in structure of kerogen from adsorption in Scheme 2. Based on the results of the equilibrium simulations, the adsorption in kerogen is computed. Because of the complex pore structure, the volume of the adsorbed phase cannot be accurately defined. We present the total uptake instead of absolute adsorption. The total uptake describes all gas molecules inside the kerogen matrix including adsorbed gas and free gas. In the constructed kerogen matrix, adsorbed gas is dominant, while the free gas is low. The results of total uptake are presented in Fig. 8. The total uptake in flexible kerogen is approximately 1.5 times higher than in the rigid kerogen. The trend is consistent with the studies based on the hybrid method of grand canonical Monte Carlo and molecular dynamics (GCMC-MD).10,23 We present snapshots of a cross-section with a width of 1 nm at the same position in the kerogen matrix for Schemes 1 and 2 (see Fig. 9). The results show that methane molecules


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


primarily appear in the large interconnected pores in Scheme 1. There are some isolated pores in the rigid kerogen where the methane molecules cannot penetrate. However, the methane molecules can flow into these pores during the evolution of the flexible kerogen structures in Scheme 2. This mechanism may provide more space for adsorption, which is one reason that flexible kerogen can absorb more than the rigid kerogen. Moreover, the total uptake reaches an approximate plateau at high pressure. The smaller change in adsorption at high pressure leads to a less incremental effect of kerogen structure change. Additionally, the adsorbed phase makes a substantial contribution to the total flux.14,21 Less incremental adsorption at higher pressure may result in lower additional flux from the adsorbed phase. These factors are the potential reasons that the transport diffusivity changes less at higher pressures than at lower pressures (Fig. 7). The results also show that the large pores are narrowed in the flexible kerogen due to adsorption-induced kerogen structure change (see Fig. 9).36 To gain insight, we calculate the pore size distribution in the kerogen matrix of Schemes 1 and 2 based on 10 snapshots at each pressure condition (see Fig. 10). The pore size distributions in Scheme 1 at different pressure conditions are identical because of the rigid kerogen matrix. In Scheme 2, with the flexible kerogen, the large pores ranging from 0.5 to 1.1 nm are significantly reduced, while the smaller pores from 0.2 to 0.3 nm are slightly increased. Since the large interconnected pore in the center is the main path for gas flow, the transport diffusivity in Scheme 2 with flexible kerogen is significantly reduced compared to that in Scheme 1.



Page 16 of 22

Fig. 8. Total uptake in the kerogen matrix with flexible and rigid structures at 333.15 K. The error bars are smaller than the symbols.

Fig. 9. Adsorption-induced kerogen structure change and its effect on the main flow path. The snapshots of a cross-section with a size of 1.00 × 7.39 × 7.39 nm3 at the same position in the kerogen matrix for the two schemes. (a) Scheme 1 (rigid kerogen, 64 bar); (b) Scheme 2 (flexible kerogen, 63 bar). The red dashed lines are drawn to illustrate the pore space of the main flow path in Schemes 1 and 2. The main flow path in Scheme 2 is significantly narrowed.


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


Fig. 10. Pore size distributions in Scheme 1 (rigid kerogen) and Scheme 2 (flexible kerogen) at five different pressures. The pore size distributions in Scheme 1 under different pressure conditions are identical due to the rigid kerogen matrix.

The transport diffusivity in Scheme 2 (fully flexible kerogen) is higher than in Scheme 3 (flexible kerogen in equilibrium simulation but rigid kerogen in flow simulation) at low pressure and becomes close at high pressure. Methane molecules can penetrate isolated pores in both schemes during the equilibrium simulation. However, when the kerogen structure is frozen in the stage of gas flow simulation in Scheme 3, some molecules can be trapped in isolated pores. We take a snapshot of the same cross-section in the kerogen matrix in Schemes 2 and 3 and trace a methane molecule that is initially trapped in an isolated pore. The results are illustrated in Fig. 11, and the dynamic process is demonstrated in the Supporting Information Movie S2. We observe that the pore throat is occasionally opened, and the molecule is released to the main path in Scheme 2, whereas the methane molecule in Scheme 3 is permanently trapped in the isolated pore. This mechanism can provide an additional flux for Scheme 2 at low pressure. The effect is not appreciable at high pressure, under which condition the flux is dominated by gas flow in the main path due to high loading.



Page 18 of 22

Fig. 11. Dynamic process of gas release/trapping in the kerogen matrix. A slice with a size of 1.00 × 7.39 × 7.39 nm3 is taken as a cross-section from the middle of the kerogen matrix. A methane molecule is marked to show the process, and the red dashed lines and arrows illustrate the trace of the methane molecule. (a to d) Process of gas release from occasionally opened pores in Scheme 2; (e to h) The marked methane molecule is trapped in Scheme 3. The snapshots are taken from the simulations at 207 bar and 210 bar in Scheme 2 and Scheme 3, respectively.

CONCLUSIONS We have constructed a kerogen matrix to investigate methane flow by using molecular dynamics simulations. The main conclusions drawn from this work are as follows: (1) The flexibility of the kerogen microstructure affects the diffusivity and adsorption. (2) Adsorption in the flexible kerogen matrix is higher than in the rigid kerogen matrix, which supports findings presented in recent literature. (3) Adsorption leads to kerogen structural change in the flexible kerogen, which may narrow the pore space and reduce gas flux.


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


(4) The transport diffusivity is higher in the rigid kerogen matrix than in the flexible kerogen matrix. In the matrix that we constructed, this difference is a factor of two between flexible and rigid kerogen. (5) The occasional opening of pore throats in flexible kerogen can provide additional flux at low pressure. This effect is not appreciable at high pressure.

SUPPORTING INFORMATION Supporting information is available online at https://pubs.acs.org/. Movie S1: Demonstration of the methods for rigid and flexible kerogens, and flexible kerogen with virtual nails. Movie S2: Dynamic process of gas release/trap in kerogen matrix.

ACKNOWLEDGMENTS The work was supported by the member companies of the Reservoir Engineering Research Institute (RERI). Their support is appreciated.



Page 20 of 22

REFERENCES (1) Wu, T.; Li, X.; Zhao, J.; Zhang, D. Multiscale Pore Structure and Its Effect on Gas Transport in Organic-Rich Shale. Water Resour. Res. 2017, 53, 5438–5450. (2) Bousige, C., et al. Realistic Molecular Model of Kerogen's Nanostructure. Nature Materials 2016, 15, 576–582. (3) Ungerer, P.; Collell, J.; Yiannourakou, M. Molecular Modeling of the Volumetric and Thermodynamic Properties of Kerogen: Influence of Organic Type and Maturity. Energy & Fuels 2014, 29, 91–105. (4) Wang, T.; Tian, S.; Li, G.; Sheng, M. Selective Adsorption of Supercritical Carbon Dioxide and Methane Binary Mixture in Shale Kerogen Nanopores. Journal of Natural Gas Science and Engineering 2018, 50, 181–188. (5) Huang, L.; Ning, Z.; Wang, Q.; Zhang, W.; Cheng, Z.; Wu, X.; Qin, H. Effect of Organic Type and Moisture on CO2/CH4 Competitive Adsorption in Kerogen with Implications for CO2 Sequestration and Enhanced CH4 Recovery. Applied Energy 2018, 210, 28–43. (6) Huang, L.; Ning, Z.; Wang, Q.; Qi, R.; Zeng, Y.; Qin, H.; Ye, H.; Zhang, W. Molecular Simulation of Adsorption Behaviors of Methane, Carbon Dioxide and Their Mixtures on Kerogen: Effect of Kerogen Maturity and Moisture Content. Fuel 2018, 211, 159–172. (7) Zhao, T.; Li, X.; Zhao, H.; Li, M. Molecular Simulation of Adsorption and Thermodynamic Properties on Type II Kerogen: Influence of Maturity and Moisture Content. Fuel 2017, 190, 198–207. (8) Zhao, T.; Li, X.; Ning, Z.; Zhao, H.; Li, M. Molecular Simulation of Methane Adsorption on Type II Kerogen with the Impact of Water Content. Journal of Petroleum Science and Engineering 2017, 161, 302–310. (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 Mesoporous Mater. 2014, 197, 194– 203. (10) Ho, T. A.; Wang, Y.; Criscenti, L. J. Chemo-Mechanical Coupling in Kerogen Gas Adsorption/Desorption. PCCP 2018, 20, 12390–12395. (11) Pathak, M.; Kweon, H.; Deo, M.; Huang, H. Kerogen Swelling and Confinement: Its Implication on Fluid Thermodynamic Properties in Shales. Scientific Reports 2017, 7, 12530. (12) Vasileiadis, M.; Peristeras, L. D.; Papavasileiou, K. D.; Economou, I. G. Transport Properties of Shale Gas in Relation to Kerogen Porosity. The Journal of Physical Chemistry C 2018, 122, 6166–6177. (13) Obliger, A.; Ulm, F.-J.; Pellenq, R. J.-M. Impact of Nanoporosity on Hydrocarbon Transport in Shales’ Organic Matter. Nano Lett. 2018, 18, 832–837. (14) Wu, T.; Zhang, D. Impact of Adsorption on Gas Transport in Nanopores. Scientific Reports 2016, 6, 23629. (15) Collell, J.; Galliero, G.; Vermorel, R.; Ungerer, P.; Yiannourakou, M.; Montel, F.; Pujol, M. Transport of Multicomponent Hydrocarbon Mixtures in Shale Organic Matter by Molecular Simulations. The Journal of Physical Chemistry C 2015, 119, 22587–22595. (16) Obliger, A.; Pellenq, R.; Ulm, F.-J.; Coasne, B. Free Volume Theory of Hydrocarbon Mixture Transport in Nanoporous Materials. The Journal of Physical Chemistry Letters 2016, 7, 3712–3717. (17) Obliger, A.; Valdenaire, P.-L.; Capit, N.; Ulm, F.-J.; Pellenq, R. J.-M.; Leyssale, J.-M. Poroelasticity of Methane-Loaded Mature and Immature Kerogen from Molecular Simulations. Langmuir 2018. 20 ACS Paragon Plus Environment

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


(18) Falk, K.; Coasne, B.; Pellenq, R.; Ulm, F.-J.; Bocquet, L. Subcontinuum Mass Transport of Condensed Hydrocarbons in Nanoporous Media. Nat. Commun. 2015, 6, 6949. (19) Jin, Z.; Firoozabadi, A. Phase Behavior and Flow in Shale Nanopores from Molecular Simulations. Fluid Phase Equilib. 2016, 430, 156–168. (20) Wu, T.; Firoozabadi, A. Molecular Simulations of Binary Gas Mixture Transport and Separation in Slit Nanopores. The Journal of Physical Chemistry C 2018, 122, 20727–20735. (21) Jin, Z.; Firoozabadi, A. Flow of Methane in Shale Nanopores at Low and High Pressure by Molecular Dynamics Simulations. J. Chem. Phys. 2015, 143, 104315. (22) Okamoto, N.; Kobayashi, K.; Liang, Y.; Murata, S.; Matsuoka, T.; Akai, T.; Takagi, S. Slip Velocity of Methane Flow in Nanopores with Kerogen and Quartz Surfaces. SPE J. 2018, 23, 1–15. (23) Tesson, S.; Firoozabadi, A. Methane Adsorption and Self-Diffusion in Shale Kerogen and Slit Nanopores by Molecular Simulations. The Journal of Physical Chemistry C 2018, 122, 23528–23542. (24) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 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 d’IFP Energies nouvelles 2013, 68, 977–994. (26) Hockney, R.; Eastwood, J. Computer Simulation Using Particles; Adam Hilger: New York, 1989, p 185. (27) Beckers, J.; Lowe, C.; De Leeuw, S. An Iterative PPPM Method for Simulating Coulombic Systems on Distributed Memory Parallel Computers. Molecular Simulation 1998, 20, 369−383. (28) Waldman, M.; Hagler, A. T. New Combining Rules for Rare Gas Van Der Waals Parameters. J. Comput. Chem. 1993, 14, 1077–1084. (29) Nosé, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511–519. (30) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695. (31) Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant Pressure Molecular Dynamics Algorithms. J. Chem. Phys. 1994, 101, 4177–4189. (32) Martin, M. G.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 1. United-Atom Description of N-Alkanes. The Journal of Physical Chemistry B 1998, 102, 2569–2577. (33) Herrera, L.; Do, D.; Nicholson, D. A Monte Carlo Integration Method to Determine Accessible Volume, Accessible Surface Area and Its Fractal Dimension. J. Colloid Interface Sci. 2010, 348, 529– 536. (34) Düren, T.; Millange, F.; Férey, G.; Walton, K. S.; Snurr, R. Q. Calculating Geometric Surface Areas as a Characterization Tool for Metal− Organic Frameworks. The Journal of Physical Chemistry C 2007, 111, 15350–15356. (35) Chen, L.; Zhang, L.; Kang, Q.; Viswanathan, H. S.; Yao, J.; Tao, W. Nanoscale Simulation of Shale Transport Properties Using the Lattice Boltzmann Method: Permeability and Diffusivity. Scientific Reports 2015, 5, 8089. (36) Wu, T.; Zhao, H.; Tesson, S.; Firoozabadi, A. Absolute Adsorption of Light Hydrocarbons and Carbon Dioxide in Shale Rock and Isolated Kerogen. Fuel 2019, 235, 855–867.



TABLE OF CONTENTS IMAGE


Page 22 of 22

Effect of Microstructural Flexibility on Methane Flow in Kerogen Matrix

Recommend Documents