Thermodynamic Modeling of the Equilibrium Partitioning of

Jan 11, 2019 - Received 29 October 2018. Published online 11 January 2019. Published in print 21 February 2019. +. Altmetric Logo Icon More Article Me...
1 downloads 0 Views 1MB Size
Subscriber access provided by ECU Libraries

Fossil Fuels

Thermodynamic Modeling of the Equilibrium Partitioning of Hydrocarbons in Nanoporous Kerogen Particles Jinlu Liu, and Walter G Chapman Energy Fuels, Just Accepted Manuscript • DOI: 10.1021/acs.energyfuels.8b03771 • Publication Date (Web): 11 Jan 2019 Downloaded from http://pubs.acs.org on January 13, 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 41 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

Thermodynamic Modeling of the Equilibrium Partitioning of Hydrocarbons in Nanoporous Kerogen Particles Jinlu Liu and Walter G. Chapman∗ Department of Chemical and Biomolecular Engineering, Rice University, 6100 Main St, Houston,TX, 77005,USA E-mail: [email protected]

Abstract In unconventional shale, petroleum compounds are generated during kerogen maturation and chemical fractionation occurs between fluids that are expelled out of kerogen and retained within kerogen due to different solubilities of hydrocarbons in kerogen. The compositional variation of petroleum compounds from organic matter source (kerogen) to nanoscale pores in the organic matter and to large fractures makes the estimation of fluid storage and the multicomponent compositional distributions in shale a challenging problem. A thermodynamic model that considers the different interactions between hydrocarbons and kerogen is needed to quantify the partitioning of fluids among different phases. In this work, a new thermodynamic model is proposed that predicts both the adsorption of hydrocarbons in nanoscale pores and the dissolution of fluids in organic matter in equilibrium with the bulk fluid phase. Firstly, by taking into account the chemical and structural similarities of kerogen with asphaltene, we leverage the PC-SAFT model of asphaltene to create a crosslinked kerogen matrix model. The swelling ratios of five kerogens of varying maturities and types in different solvents are

1

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

well quantified by phase equilibrium calculation using the new model. Secondly, the new kerogen matrix model is incorporated in interfacial Statistical Associating Fluid Theory (iSAFT) to form a nanoporous kerogen composite model. The inhomogeneous fluid distribution at the fluid-kerogen interface is then able to be solved and the equilibrium partitioning of mixtures of fluids from bulk phase to nanoporous kerogen phase can be determined. The results show the enrichment of aromatic and cyclic molecules in kerogen, the preferential adsorption/absorption of high molecular weight molecules and carbon dioxide, and the significant amount of hydrocarbon storage both in pore space and kerogen matrix. Examples that compare with fluid at the same condition in a graphite pore indicate that graphite as the post matured kerogen provides the largest fractionation effect for fluid mixtures.

2

ACS Paragon Plus Environment

Page 2 of 41

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

Introduction Production of oil and gas from unconventional shale reservoirs has increased dramatically over the past decade and the rapid development is expected to continue in the foreseeable future. 1 The worldwide abundant distribution of shale formations makes unconventional resources an attractive solution to increasing energy demand. Shale reservoirs are classified as "unconventional" because the oil and gas are trapped in the source rock itself, unlike in "conventional" reservoirs where the hydrocarbons have migrated out of the source rock into more porous rock in a reservoir. Shale reservoirs have extremely low permeability due to the small (nanoscale) pore sizes and complex pore structure, which requires hydraulic fracturing to create flow paths to make producing the oil and gas commercially viable. Even with fracturing, the oil recovery factor from shale is still very low and reported to be in the range of 1% to 10%. 2 Reservoir characterization and production forecasting are challenging as hydrocarbons in the source rock environment have unusual storage mechanisms, phase behavior and flow behavior. The uncertainty results from the unique chemical and petrophysical properties of each reservoir ranging from nanoscale to field scale. A significant amount of storage in shale is due to fluid both adsorbed in the nanopores and absorbed/dissolved in the organic matrix. This fluid storage mechanism has been discussed in a number of studies, 3–7 however, few studies have tackled this problem theoretically from a molecular modeling perspective. In this study, we develop a molecular scale model to predict the hydrocarbon storage capacity and fluid distribution in kerogen matrix and nano-pores present in the kerogen. To estimate the amount of fluid stored in the reservoir, pore structures of shale should be firstly understood as pores provide the storage environment of hydrocarbons. Pores can be directly visualized through different kinds of microscopes such as scanning-electron microscope (SEM) and Helium ion microscope (HIM), 8 or indirectly measured by mercury injection capillary pressure(MICP) and nuclear magnetic resonance (NMR). 9 It has been recognized by many studies that the pore sizes in shale range from a few nanometers to several microns. Pores between grains of organic or inorganic matter are called inter-particle 3

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 41

pores and range from tens of nm to 10 µm. These macropores (> 50 nm) are connected to micropores (< 2nm) and mesopores (2 - 50 nm) which are contained within the organic matter kerogen. (See Figure 1(a) for schematic illustration. This type of microstructure can be clearly seen from HIM images of organic rich shale samples. 8 ) Researchers suggest that these micro- and mesopores form as a result of hydrocarbon generation and phase separation in the kerogen. These pores are also classified as secondary organic porosity 10 and in the analysis by King et al., 8 they exist within a narrow spectrum rather than a broad distribution. nanoporous kerogen particle

H pore

kerogen particle

inter-particle pore (macropore)

intra-particle pore (mesopore-micropore)

kerogen molecule

+ solvent hydrocarbon molecules

kerogen swells in organic solvents

kerogen pore model

(a)

z

(b)

Figure 1: (a) Illustration of the interparticle and intraparticle organic pores in kerogen. Macropores between kerogen grains are called inter-particle pores. Meso/micropores are found inside the kerogen grain particle which are called intra-particle pores. The kerogen matrix does not dissolve in any solvent but becomes swollen in organic solvent. (b) Schematic representation of the modeled nanosized kerogen pores. Pore walls are formed by kerogen molecules and fluids can penetrate the pore wall to dissolve in the kerogen matrix.

The existence of nanosized pores related to organic matter raises questions regarding the hydrocarbon thermodynamic and transport behaviors. Hydrocarbons inside the organic nano pores behave differently from bulk fluid because they are affected by the substantial interactions with the pore surfaces as the pore size is within the range of intermolecular forces. Significant surface adsorption happens on the pore surface because hydrocarbons 4

ACS Paragon Plus Environment

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

tend to be attracted by the organic pore wall. The individual fluid-surface interaction of each molecule depends on the molecule’s chemical composition, size, shape, and also properties of the surface. The adsorption and phase behavior of hydrocarbons in shale have been investigated by many researchers either on synthesized materials with known pore structure or shale samples. 11,12 Numerous molecular simulation and modeling studies have also provided detailed insights into the fluid behavior in a nanopore system. 13–17 It has been commonly acknowledged that hydrocarbons inside the organic nanopores have reduced bubble points / critical points and narrowed phase envelopes. Another important phenomenon is that hydrocarbons dissolve in the organic matter, which forms the third hydrocarbon storage environment in shale apart from the free and adsorbed gas/oil. 5–7,18 The dissolution of oil/gas in the organic matter is expected because oil/gas expulsion from organic matter is a result of limited solubility of the hydrocarbons produced from thermal diagenesis. The dissolved fluid is in thermodynamic equilibrium with the expelled fluid and over a long period of time, the chemical reactions inside the source rock continue with temperature increments and new equilibrium states will establish. With fluid properties strongly dependent on the nanopore properties, studies on the characterization of the solid organic matter in shale also emerged since the 1990s. 19–27 The organic matter is mainly composed of kerogen, which is insoluble in any organic solvent. Kerogen originates from the biological degradation of living matter during which the decomposition of biomolecules and polymerization of the decomposed units occur. Because of the different origins, geological environments and thermal alteration history, kerogens have a wide range of chemical compositions and structures. Nevertheless, kerogen is believed to be macromolecules with a three dimensional cross-linked network. It is of particular interest for people to understand the hydrocarbon generation potential of kerogen samples from a reservoir. Solvent swelling experiments have been used to quantify the fluid retention in kerogen and the kerogen-solvent interactions. If a solvent is more favored by the kerogen or they are more chemically compatible, kerogen would absorb more solvent and expand more

5

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

in volume. Regular Solution Theory of Hildebrand 28 is applied to understand the variation in swelling ratios of kerogen in different solvents assuming that no specific interactions such as hydrogen bonding and strong dipole-dipole interactions exist., 5–7 The free energy change with respect to solvent dissolution will be balanced by the elastic free energy arisen from network expansion. Flory-Rehner equation 29 and its extension to high crosslinking density by Kovac 30 and by Barr-Howell and Peppas 31 are commonly applied to estimate the degree of cross-linking in kerogen. Regular Solution Theory is limited in assuming no specific interactions which results in insufficiency in describing all the experimental results. Earlier experiments with coal indicate that excessive swelling is observed with hydrogen bonding solvents. A more powerful thermodynamic model is needed to quantify the interactions between kerogen and organic solvents and predict the solvent dissolution. Statistical Associating Fluid Theory (SAFT) 32 and the subsequent variants like Perturbed chain-SAFT (PC-SAFT) 33 are able to model complex fluids with hydrogen bonding, to name a few, like crude oil, 34 homopolymers, 35 copolymers, 36 cross-linked polymers, 37,38 branched polymers, 39 ionic liquids, 40 surfactants, 41 gas hydrate, 42 electrolyte solutions 43 etc. Of greatest significance for the purpose of this work is the sucess of PC-SAFT in modeling the phase behavior of asphaltenes in crude oil. 44–48 In SAFT, the Helmholtz free energy is a summation of different contributions, which are the ideal gas, reference fluid, chain formation, dispersion and specific interaction contributions as mentioned before. Molecules are modeled as chains of spherical segments. Then for molecules without hydrogen bonding or dipole/quadrupole moments, the number of segments (m), the segment diameter (σ) and attraction energy between segments (/k) are the three parameters to describe a molecule. Asphaltene components have been modeled as polyatomic molecules in previous studies and asphaltene precipitation pressure from crude oil has been successfully predicted. 46 Another related example of application is the solvent swelling of crosslinked polymers such as PNIPAAM hydrogel 49 and vapor sorption in cross-linked natural rubber. 50 In modeling the swelling of cross-linked polymers in solvents with PC-SAFT EoS, the elastic

6

ACS Paragon Plus Environment

Page 6 of 41

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

free energy from Flory-Huggins polymer solution theory 51 is incorporated with a constraint applied on the extent of network expansion. 52 As discussed above, interfacial properties and phase behavior of confined fluids is affected by molecular size, shape, and specific interactions. In nanopore systems, where the pore size is of the same size range as molecular interactions, density functional theory has been found to accurately describe fluid structure and properties. Within molecular density functional theory (DFT), a SAFT molecular model can be applied to inhomogeneous fluids (iSAFT) with numerous studies in the literature. 53–56 A good approximation of the inhomogeneous fluid free energy functional is essential and those, like iSAFT, that can reduce to a well established and parameterized equation of state in the bulk is of particular interest for practical applications of interfacial phenomena. These applications include phase behavior of surfactants and micelle formation, 56,57 conformation of surface tethered polymers and dendrimer molecules in different solvents, 58,59 block copolymer/nanoparticle composite, 60 hydrophobic and hydrophilic interactions with a surface, 61 etc. For the hydrocarbon storage and distribution problem here, it is particularly necessary to apply iSAFT to solve for the detailed distribution of molecules within the nanosized pores and the kerogen matrix as has been done for graphite pores. 13 In this work we firstly outline the theoretical framework for (1) bulk kerogen swelling using the PC-SAFT EoS and free energy change due to kerogen network elasticity; 52 (2) inhomogeneous fluid adsorption in kerogen pores and absorption in kerogen matrix using DFT. Then we present how molecular modeling parameters for kerogen matrix are developed based on experimental data and reasonable assumptions. Five kerogen samples from literature are characterized in this step. Subsequently, the established kerogen parameters are used as input to form a permeable slit pore that is the basis for our nanoporous kerogen composite model as shown in Figure 1. By assuming the fluid in the nanoporous kerogen particle is in equilibrium with the fluid in macropores or bulk fluid phase, the hydrocarbon mixture partitioning behavior among the bulk phase, intraparticle nanopore space and

7

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 41

kerogen matrix is studied to quantify the fluid in place and the compositional variations of different fluid environments. Specifically we show examples of partitioning behavior for mixtures with the same carbon number and a Eagle Ford shale gas sample. A comparison is also made with the fluid partitioning to an impermeable graphite pore, 13 which has been widely used to model the shale nanopores.

Method In this section, we present the theory and numerical methods for modeling bulk kerogen as a crosslinked polymer network through PC-SAFT EoS and the inhomogneous nanoporous kerogen composite system through DFT.

PC-SAFT model of cross-linked kerogen For polyaromatic molecules such as asphaltenes, PC-SAFT has been demonstrated to be effective in describing the asphaltene phase behavior in crude oils as a function of temperature, pressure and different solvents and gases including CH4 and CO2 . 45,46,62 Asphaltenes are the dispersed and chemically altered fragments of kerogen that can dissolve in organic solvents under certain conditions. Considering kerogen is a macromolecule with polyaromatic islands crosslinked by low molecular weight linkers, as was shown by Ungerer et al.’s atomistic models, 63 we are confident that PC-SAFT is an appropriate choice for modeling the solubility of components in kerogen. We model cross-linking in the kerogen by analogy to a model for cross-linked polymers. With the addition of an elastic free energy term, PC-SAFT has been used to model the swelling behavior and volatile organic solvent sorption of cross-linked natural rubber. 37 Using this same approach, the total free energy in a kerogen network is a sum of different contributions, Atotal = Aid + Aseg + Ach + Adisp + Aelast

8

ACS Paragon Plus Environment

(1)

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

where Aid is the ideal gas contribution, Aseg is the residual free energy of spherical segments, Ach is the change in free energy on bonding segments to form molecules, Adisp is from the long-range attraction between molecules, and Aelast is the free energy change due to crosslinked polymer network swelling when solvent molecules enter the polymer matrix. Details of the segment, chain and dispersion terms are given in early SAFT papers. 32,33 The elastic free energy is taken from Miao et al., 52 where a finite extensibility of the network is considered. We use the same equation by Arndt: 49 "  # Φ−2 3 (V /V0 )2/3 − 1 V Aelastic = nchain × − ln 2/3 kb T Φ 2 1 − (V /Vmax ) V0

(2)

where nchain is the number of polymer chains between two crosslinkers, kb is the Boltzmann constant, Φ is the network functionality which is treated as a fitting parameter, representing the deviation from a perfect polymer network (a perfect network would have a functionality value of 4; a larger functionality indicates more entanglements of the chains). V is the molar volume of each polymer chain after swelling, V0 is the molar volume of dry polymer chain, Vmax is the maximum volume that the swollen polymer network can reach and is calculated by assuming an ideal tetrahedral network, √ Vmax = ncross−linker NA

3

2 sin 2





109.5 π 2 180◦

!3

 mp σ p

(3)

The chain length between cross-links is given by the product of the number of monomer segments in each chain mp and their segment diameter σp . ncross−linker is the number of cross-linker molecules, which can be calculated by the total mass of linkers and the linker molecular weight. ncross−linker =

mcross−linker Mcross−linker

(4)

The number of chains between linkers is determined by the total number of linkers and its

9

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 41

functionality, nchain = ncross−linker ·

Φ 2

(5)

Then the molecular weight of polymer chains between linkers can finally be calculated by

Mc =

mpolymer + mcross−linker nchain

(6)

where mpolymer is the total mass of polymer before adding linker molecules. The elastic contribution to the total pressure of the cross-linked polymer network phase is thus obtained by p

elast

 elast  ∂A =− ∂V T,ni

(7)

Phase equilibrium calculations As kerogen is insoluble in any solvent, the phase equilibrium between the solvent phase (I) and the swollen kerogen phase(II) is obtained by equalizing the total pressure and the chemical potential of each solvent in both phases.

pI = pII

µI(solvent)i = µII (solvent)i

(8)

i = 1, 2, ...N

(9)

where N is the total number of solvent species.

Molecular density functional theory Molecular density functional theory is used to model the kerogen/fluid interface and fluid confined in nanopores. Our permeable pore model is illustrated by Figure 1. The two slabs of kerogen form a slit pore with permeable kerogen matrix walls. The fluid is inhomogeneous in z direction and is present both in the slit pore and the solid matrix. The pore is repetitive in the z direction to mimic a porous structure. The solubility of each fluid molecule is 10

ACS Paragon Plus Environment

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

determined by the bulk phase equilibrium between hydrocarbon fluids and kerogen solid. The kerogen surface is not smooth and we treat the surface roughness of the kerogen by assuming a diffusive kerogen density profile as was suggested by Ravikovitch and Neimark 64 in the quenched-solid density functional theory. 65 The kerogen density profile is given by,    ρbulk 2 (z − z0 ) s ρs (z) = 1 − tanh 2 δ

(10)

is the equilibrium bulk density of kerogen matrix, z0 = 0 is the solid boundary where ρbulk s and δ characterizes the solid roughness – a larger δ value indicates a rougher surface. In this paper, we used δ = 2.1Å throughout the work. The grand potential of the kerogen-fluid interfacial system is a functional of the hydrocarbon and kerogen density distributions

Ωsf [ρs (r); ρi (r)] = Atotal [ρs (r); ρi (r)] −

N Z X

ρi (r)µbulk dr + i

Z

! ρs (r)µbulk dr s

(11)

i=1

where Atotal [ρs (r); ρ(r)] is the total intrinsic Helmholtz free energy functional, ρi (r) is the is the bulk chemical potential of compoent i. density distribution of component i and µbulk i The free energy functional is consistent with the PC-SAFT equation of state; it contains a sum of different contributions, i.e. ideal gas, short range hard sphere repulsion, chain contribution and long range attraction in perturbation theory. The details of each term can be found in Appendix A.

Atotal [ρs (r); ρ(r)] = Aid [ρs (r); ρ(r)]+Ahs [ρs (r); ρ(r)]+Ach [ρs (r); ρ(r)]+Adisp [ρs (r); ρ(r)] (12)

At equilibrium, the grand potential as in eq 11 is a minimum, or in other words, fluid adsorbed in the pore and absorbed/dissolved in the kerogen have the same chemical potential

11

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 41

with that of the bulk phase fluid. δ Ωsf [ρs (r); ρ(r)] = 0 for i = 1, 2, ..., N δ ρi (r) equilibrium

kb T ln [ρi (r)] +

δ Ahs [ρs (r); ρ(r)] δ Ach [ρs (r); ρ(r)] δ Adisp [ρs (r); ρ(r)] + + = µbulk i δ ρi (r) δ ρi (r) δ ρi (r)

(13)

(14)

Note that the variation on density is only performed on fluid densities, since the diffusive kerogen density profile is set by eq 10. This means when solving the equations iteratively, only the fluid densities are iterated. For each functional derivative in eq 14, we refer the readers to Sauer and Gross 66 for hard sphere and dispersion part and Liu et al. 13 for chain part.

Results and Discussion PC-SAFT modeling of kerogen swelling in solvents As kerogen is a type of molecule far from being well-defined in chemical composition and structure, information from physical and chemical analysis such as density, composition in element and functional groups only represents the averaged properties of all the kerogen molecules in a sample. In a Van Krevelen diagram, 67 it can be seen that kerogen composition in C, H and O varies depending on their geological environments (i.e. lacustrine, marine, terrestrial) and the maturation stages (temperature and pressure history). Type I kerogen has the highest H/C ratio followed by Type II and III, and Type I is the most oil-prone kerogen. As kerogens mature, both H/C and O/H ratios decrease due to expulsion of hydrocarbons and water, thus the fraction of aromatic carbons increase. The three types of kerogens evolve to similar composition as they become over mature. By combining the solvent swelling experiment with the Flory-Rehner equation, it was demonstrated that the swelling ratios decrease and the cross-linking density of kerogen increases as kerogen becomes

12

ACS Paragon Plus Environment

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

more mature. 23 The decrease in solvent uptake as kerogen matures implies a decreasing potential of kerogen to generate hydrocarbons through pyrolysis. As a direct measurement of kerogen cross-linking is not available, information related to cross-linking can only be revealed by analyzing the experimental measurements with theoretical models. In this study, we model the kerogen swelling ratios in different solvents by calculating the phase equilibrium between solvent and swollen kerogen phase. Firstly the bulk fluid temperature, pressure and composition are set. Then the bulk fluid densities are calculated using PC-SAFT EoS and chemical potentials of each component are determined. After that, the equilibrium density of each component in the kerogen phase, ρkerogen,b , is determined by i eq 8 and eq 9. By solving the phase equilibrium, the swollen kerogen density, ρswollen,bulk s is also obtained. The swelling ratio, Qv , is the ratio of the swollen kerogen volume to the dry kerogen volume, and the swollen kerogen volume includes both the kerogen and solvent volumes. Qv =

ρdry s ρswollen mixture wtkerogen

(15)

In eq 15, wtkerogen is the weight fraction of kerogen in the kerogen/solvent mixture, ρdry s and ρswollen mixture are the dry kerogen mass density and kerogen/solvent mixture mass density, respectively. In a series of experiments by Larsen et al., 20–23 kerogen samples were firstly extracted from shale by acid and CrCl2 demineralization and exposed to excess amount of solvent to swell. The volumetric swelling ratio was measured until a constant value is reached after enough time. In these experiments, although the initial packing and internal structure of kerogen particles were not specified, the kerogen powder swells without external constraint, thus it is assumed the inter- and intra-granular pores expand in the same ratio as the solid matrix. Therefore, the swelling experiments are modeled by bulk phase equilibrium without considering the details of porosity change. As kerogen is believed to be composed mainly of polyaromatics in a 3D cross-linked network, it is modeled as a network of asphaltene-like molecules with low molecular weight linkers added. While asphaltene is soluble in crude oil under certain conditions, the crosslink13

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 41

ing in kerogen is one of the reasons that the kerogen does not dissolve in any solvent. Besides the three parameters in the PC-SAFT EoS for non-associating non-polar molecules, another three parameters regarding the cross-linking property are needed. These are the mass of cross-linkers (mcross−linker , here the unit is g/ (100 g of non-crosslinked polymer), molecular weight of cross-linkers (Mcross−linker ) and network functionality (Φ). Determination of these parameters requires some reasonable assumptions to firstly reduce the number of parameters and then fitting the rest of the parameters to match the dry kerogen and swollen kerogen properties. In previous studies of asphaltenes, a correlation has been developed by Punnapala and Vargas 68 to relate the PC-SAFT parameters to asphaltene molecular weight (Mw) and aromaticity (γ), which reads,

m = (1 − γ)(0.057Mw + 0.8444) + γ(0.0101Mw + 1.7296)     93.98 ln(Mw) + γ 4.6169 − σ = (1 − γ) 4.047 − 4.8013 Mw Mw     9.523 234100 ε/kb = (1 − γ) exp 5.5769 − + γ 508 − Mw Mw1.5

(16)

(17) (18)

By taking the advantage of this correlation, the number of parameters is reduced from 6 to 5. To further reduce the number of parameters to fit, we make assumptions of the molecular weight of different kerogens based on the reported kerogen molecular structure by Ungerer et al.. 63 The linker molecular weight, however, is a bit arbitrary as not much chemical analysis of the exact fraction of linker molecules is available to our knowledge. Pawar et al. 69 simulated the maturation process of kerogen by reactive force field and concluded that the generated cross-linkers are mostly -C-C-, -C-O- and -C-S- covalent bonds. We thus assume the average molecular weight of linkers is 32 g/mol. With these assumptions, three parameters (γ,mcross−linker and Φ) are left to be fit. The numerical procedure we take to fit the parameters is to use grid search over a [γ,mcross−linker ] surface and determine the corre-

14

ACS Paragon Plus Environment

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

sponding Φ to match the dry kerogen density. Then this combination of [γ, mcross−linker , Φ] is used to calculate the swelling ratios (according to Eq 8,9 and 15) and deviation from experimental values of the kerogen in different solvents. The set of parameters that gives the minimum average percent deviation from experiments is chosen as the model parameters. An example of parameter fitting result is shown in Figure 2 with only the parameters that give less than 50% errors plotted. We have characterized 5 kerogen samples of Type I and Type II from experiments by Larsen et al. 21,23 and the resulting parameters are shown in Table 2. Only non-polar nonassociating solvents are included ( Table 1) for fitting because specific interactions between the functional groups in kerogen and solvents are possible, such as B2 sample shows enhanced swelling with hydrogen bond acceptors. The solvent parameters are directly taken from PC-SAFT 33 and the binary interaction parameters (kij ) of each solvent with kerogen follow the settings by Punnapala and Vargas 68 of hydrocarbons with asphaltenes. kij 0 s between n-alkanes and kerogen are set to be 0.01 indicating the incompatibility of alkanes with polyaromatics. Different from Punnapala and Vargas, 68 the smallest kij = −0.004 is given to all aromatic molecules instead of saturates as we consider the highest structural compatibility of aromatics with kerogen. The saturated cyclo-alkanes then have an intermediate kij = 0 with kerogen. The same rule is also applied in the rest of this work. Samples A1 to A3 are Type II kerogen from Paris Basin Toarcian shale in the oil generation window that comprise a very short maturation series. They are well characterized to have H/C ratio of 1.48, 1.24 and 1.04, respectively. The four Type II kerogens (IIA to IID) from Duvernay series during maturation reported by Ungerer et al. 63 have H/C ratios from 1.17 to 0.58 and experimentally measured densities range from 1.18 to 1.4 g/cm3 . Another Type II kerogen from Estonian kukersite shale is measured to have a H/C 1.52 and density of 1.11 g/cm3 . Based on a simple correlation of H/C and density, we estimate the densities of A1, A2 and A3 to be 1.10, 1.15 and 1.20 g/cm3 . The molecular models from Ungerer et al. 63 show molecular weight for Type II kerogen in the middle of oil generation window to be 3500 g/mol, so

15

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the same values are taken for A1 to A3 samples. For B1 sample from Green River oil shale, the same molecular weight and density are taken from IA sample of Ungerer et al. 63 as they have almost the same H/C. For B2 sample, the H/C is 1.61 and close to B1 of 1.53, thus it is reasonable to assume the same molecular weight and density with B1. The molecular weight from different research, however, has a fairly large range. For example, while Ungerer et al.’s atomistic models show kerogen molecular weights in a range of 2500 to 4000 g/mol, higher molecular weight models are proposed, such as 6000 g/mol for a Type II Estonian kukersite kerogen, 70 9000 g/mol for Green River kerogen and even 30000 g/mol for Rundle kerogen by Siskin et al.. 19 The large range of heteroatom functionalities and consideration of other elements such as Cl require higher molecular weight, thus adjustment in molecular weights can be investigated to better represent the real molecules. The fitted swelling ratios for each sample is shown in Figure 3. Better agreement with experiments is expected with kij tunning, however, the experimental measurements for Qv using the same method is reported to have an error of ±0.1, thus it may not be necessary to match experiments exactly. From Figure 3, we can see the modeled swelling ratios are in reasonable agreement with experiments, with the largest deviation seen for B1 sample with CS2 and biphenyl, and B2 sample with cyclohexane. Apart from these, all the other deviations are within 0.1. This shows that the PC-SAFT EoS along with an elastic model describes the physics of kerogen interaction with different organic solvents and the relative swelling behavior of different kinds of solvents are well quantified. Considering a series of assumptions being made, it is encouraging that such reasonable descriptions are achieved. Table 2 shows the characterized cross-linking in kerogen. For samples A1 to A3, the aromaticity increases as kerogen matures and the fraction of cross-linkers also increases. For B1 and B2 which have much higher H/C ratio, the cross-linking densities and aromaticity are much lower. This trend follows the chemical change in kerogen maturation and the inverse correlation of aromaticity with H/C ratio.

16

ACS Paragon Plus Environment

Page 16 of 41

Page 17 of 41

0.8 45 0.75

40 35

aromaticity

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

0.7

30 25

0.65

20 15

0.6

10 5

0.55 2

2.2

2.4

2.6

2.8

3

linker fraction Figure 2: Example of parameter fitting for kerogen A1 showing the errors (%) associated with each set of γ and mcross−linker . Colors correspond to different fitting errors. Although the grid search is done in the whole domain, only the points that result in errors less than 50 % are shown in the figure.

(b)

(a)

Figure 3: Fitting results of PC-SAFT model to kerogen swelling experiments for (a) Paris Basin Toarcian kerogen (b) Rundle kerogen. (Lines are from modeling and bars are from experiments 21,23 )

17

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 41

Table 1: PC-SAFT parameters for solvent molecules 33,68 No. 1 2 3 4 5 6 7 8 9 10 11

Mw [g/mol] n-pentane 72.146 n-heptane 100.203 methylcyclohexane 98.182 cyclohexane 84.147 o-xylene 107.167 benzene 78.114 toluene 92.141 tetralin 132.205 chlorobenzene 112.558 carbon disulfide 76.143 biphenyl 154.211 solvent

m [-] 2.6896 3.4831 2.6637 2.5303 3.1362 2.4653 2.8149 3.3131 2.6485 1.6919 3.8877

σ [Å] 3.7729 3.8049 3.9993 3.8499 3.7600 3.6478 3.7169 3.875 3.7533 3.6172 3.8151

ε/k kij with kerogen [K] 231.20 0.01 238.40 0.01 282.33 0 278.11 0 291.05 -0.004 287.35 -0.004 285.69 -0.004 325.07 -0.004 315.04 -0.004 334.82 -0.004 327.42 -0.004

Table 2: Parameters used and fitted for 5 kerogens kerogen A1 A2 A3 B1 B2

source

Mw density γ mcross−linker [g/mol] [g/cc] [-] [g/100 g] II immature 3500 1.10 0.69 2.300 II partially mature 3500 1.15 0.76 2.775 II most mature 3500 1.20 0.82 2.950 3800 1.00 0.41 0.575 I immature 3800 1.00 0.40 0.425 type

Paris Basin Toarcian

Type Type Type

Green River Rundle

Type

18

ACS Paragon Plus Environment

Φ [-] 12.48 11.26 11.23 24.17 23.53

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

Chemical fractionation of mixed components of equal carbon number The main purpose of the solvent swelling experiments done by many researchers so far is to quantify the solubilities of different solvents in a specific kerogen sample, thus to understand the expulsion behavior of different hydrocarbons from the source rock and provide a basis for expulsion modeling. In general, it has been accepted that the expulsion potential of different types of hydrocarbons should follow: n-alkanes > branched alkanes > cyclic > aromatics > polars & resins > asphaltenes. This causes the petroleum expelled from the source rock to have a different composition from the initially generated oil by kerogen, which is called chemical fractionation. The fractionation also happens between produced oil and reservoir fluid. During shale gas production, a compositional change of the produced gas is commonly reported. Tinni et al. 71 reported an increase in fraction of methane even when the pressure is above the dew point of the reservoir fluid methane/butane mixture. Molecular simulation with a modeled kerogen shows the preferential adsorption of n-butane and the subsequent pore blockage by n-butane reduces the produced butane fraction. In this section, we consider three fluid mixtures, each of which have molecules with the same carbon number of 5, 6 or 7, and one group of benzene and its derivatives as examples to show the chemical fractionation of source rock fluid among bulk phase, adsorbed phase in the nanopore and absorbed phase in kerogen. For each mixture, we have equal mole density of each hydrocarbon in the bulk phase, which corresponds to the produced fluid phase. The PC-SAFT parameters are taken from Gross and Sadowski 33 and the binary interaction parameters kij following Punnapala and Vargas 68 are tabulated in Table 3. We take a kerogen A1(refer to Table 2) matrix surrounding a 3 nm pore as an example to show the partitioning of bulk fluids to pores and kerogen matrix. For the carbon number 6 mixture, we also examine how kerogen maturity affects the fractionation. Density profiles of each hydrocarbon mixture at the interface of kerogen and pore space are shown in Figure 4a,4c and 4d, with only the left half of Figure 1 presented considering the symmetry of the system. Because of the diffusive interface and dissolution of molecules to the kerogen matrix, the fluid 19

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

layering close to the surface and the surface excess adsorption is not as strong as in a graphite pore model, which means the solid confinement effect is restrained. This phenomenon is believed to be closer to the real system, especially for liquid rich shales. To illustrate the difference, the distribution of a mixture with carbon number 6 in a graphite pore is also calculated as shown in Figure 4b. The graphite pore can serve as a model of the extreme case as kerogen pore becomes mature, although in reality we only expect gas to be present in a post mature kerogen pore. Fractions of each component in bulk phase, nanopores(adsorbed) and kerogen matrix (dissolved) are shown in the bar charts where colors correspond to component curves in the figure. The nanoporosity is assumed to be 15%, which falls in the range of various experimental measurements 72,73 and molecular simulation of kerogen. 74,75 The mole fraction of each component in the nanoporous kerogen composite is then calculated and also shown in the chart. The deviation from bulk composition is noticeable both in adsorbed phase and dissolved phase. In Figure 4a, we see that benzene is the largest fraction, saturated rings are the intermediate fraction and alkanes are the least. The same trend is followed also for C7 and C8 mixtures. The fractionation has a larger effect in the dissolved phase, and as the kerogen becomes mature (see Figure 4a), this effect becomes more substantial. This means as petroleum expels from the source rock, the aromatic molecules are retained and alkanes will first migrate from source rock to form the second phase in an oil reservoir. Figure 4a also compares the two different density profiles of nanopores in either an A1 or A2 kerogen matrix. For an A2 sample formed pore, kerogen has a higher density which makes the dissolution of hydrocarbons significantly decrease, however, the surface attraction from the kerogen pore wall becomes stronger, which makes the surface adsorption of benzene a little higher than in the A1 sample. Comparing the fluid compositions in a graphite pore with that in a kerogen pore, we see that the same trend of distribution among the five molecules is followed but the fractionation is larger in the graphite pore. This can be viewed as the limiting case of fluid composition after diagenesis.

20

ACS Paragon Plus Environment

Page 20 of 41

Page 21 of 41

2500

18000 16000 14000

density [mol/m3]

density [mol/m3]

1500

H pore

carbon # = 6

kerogen hexane 2-methylpentane cyclohexane methylcyclopentane benzene

2000

solid line : A1 dashed line: A2

1000 5000 4000

12000 10000 8000

hexane 2-methylpentane cyclohexane methylcyclopentane benzene

6000

3000

500

4000

2000 1000

2000

0

0

0

5

composite dissolved

adsorbed

0 bulk A3 A2 A1

A3 A2 A1 A3 A2 A1

10

z/Å

15

10

20

20

30

25

0

30

0

20.0% 20.0% 20.0% 19.3% 18.9% 20.3% 20.0% 21.5% 19.4% 19.1% 20.3% 20.0% 21.3% 19.4% 19.1% 20.3% 20.0% 21.3% 10.4% 8.4% 20.1% 17.5% 43.6% 11.9% 9.9% 20.8% 18.4% 39.0% 12.8% 10.7% 21.0% 18.8% 36.7% 15.4% 14.3% 20.2% 18.9% 31.3% 15.2% 14.0% 20.5% 19.1% 31.1% 15.2% 13.8% 20.7% 19.2% 31.0% 20.0%

5

10

20.0%

15

z/Å bulk

20.0%

graphite pore

17.4%

20.0%

20.0%

15.4%

(a)

21.8%

20.0% 18.1%

20.0% 27.3%

(b) 2500

2500 kerogen

carbon # = 7

carbon # = 8 kerogen

heptane 2000

2000

2-methylhexane

density [mol/m3]

density [mol/m3]

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

methylcyclohexane ethylcyclopentane

1500

toluene 1000

octane ethylcyclohexane o-xylene

1500

ethylbenzene

1000

500

500

0

0 0

5

10

15

20

25

0

30

5

10

15

bulk

20.0%

20.0%

adsorbed

19.5%

19.2%

dissolved 13.1%10.8% 16.9% composite

20

25

30

z/Å

z/Å 20.0% 19.8%

20.0% 19.9%

17.4%

15.6% 14.0% 18.0%

18.4%

20.0%

bulk

25.0%

21.6%

adsorbed

23.5%

25.0% 23.7%

dissolved 11.6%12.1%

41.7%

composite

34.0%

(c)

16.0% 16.5%

25.0% 26.4%

38.7%

34.1%

25.0% 26.3% 37.6%

33.4%

(d)

Figure 4: Density profiles of hydrocarbon mixtures and composition in each phase. (a) carbon number = 6 mixtures in A1 - A3 kerogen pores. (b) carbon number = 6 mixture in a graphite pore. (c) carbon number = 7 mixture in A1 kerogen pore. (d) carbon number = 8 mixture in A1 kerogen pore. Temperature and pressure are fixed to be T = 328 K, P = 50 MPa. The bar charts show the fraction of each component in the bulk, adsorbed, dissolved and composite based on 15 % nanopore porosity. Colors in the bar charts correspond to component curve colors in the figures. 21

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 41

Table 3: PC-SAFT binary interaction parameters used in the model 68 Components N2 CO2 C1 C2 C3 C4+ alkanes saturated rings aromatics kerogen

N2 0

CO2 0 0

C1 0.03 0.05 0

C2 0.04 0.097 0 0

C3 C4+ alkanes 0.06 0.075 0.1 0.12 0 0.03 0 0.02 0 0.015 0

saturated rings aromatics 0.14 0.158 0.13 0.1 0.03 0.029 0.012 0.025 0.01 0.01 0.005 0.012 0 0.007 0

Fractionation of shale gas Composition of produced shale gas is observed to change as a function of time even when the reservoir pressure is higher than the mixture saturation pressure. It has been recognized that a combination of molecular sieving effect, preferential adsorption/absorption and variation in diffusion coefficients will increase light component fraction in produced fluid over time. 71 This has created additional complexity for estimating the fluid in place and calculating the ultimate recovery factor. In this section, we take a gas sample from Eagle Ford shale play 76 to study the partitioning of different components in nanoporous kerogen. The unprocessed pipeline gas composition and its phase envelope is shown in Figure 5. Temperature is fixed at 328 K and pressure varies. At this temperature the fluid is in either supercritical state or gaseous state. In Figure 6 a-c, the adsorption/absorption isotherms of this shale gas sample in a 3 nm graphite pore and a 3 nm pore in a Sample A3 kerogen matrix are presented. The density change in pore space and in kerogen matrix are shown separately in Figure 6b and 6c, with the density in the kerogen matrix much lower than the density in the pore space. Comparing Figure 6c with 6a, we see that the heavy components such as propane, butane and hexane have substantially lower adsorbed amount in the kerogen pore than in a graphite pore. This is due in part to weaker wall attractions in the kerogen pore. Although the permeable pore is not as selective to fluid molecules as an impermeable pore, the fractionation effect from the kerogen matrix is strong, which makes the heavy components 22

ACS Paragon Plus Environment

kerogen 0.16 0.1 0.07 0.06 0.01 0.01 0 -0.004 0

Page 23 of 41

remain in the nanoporous kerogen. As shown in Figure 6d, the selectivity – which is defined as the ratio between fraction in sorbed (for graphite pore only adsorption exists) phase and the fraction in bulk phase – for heavy components is less in nanoporous kerogen composite than in a graphite pore. In both systems, methane has a selectivity smaller than 1 while butane and hexane have very high selectivity in the sorbed phase, which explains that heavy components are most likely to be left in the gas reservoir after pressure depletion. The fractionation effect is moderated in a kerogen composite only from an equilibrium sorption perspective with the absorbed phase contributing more to the total fractionation effect, while in gas production the much faster diffusion of light components will enhance the selectivity of heavy hydrocarbons to be retained in reservoirs. 120

component mole fraction (%) C1 74.60 C2 13.82 C3 5.43 C4+ 4.46 C6+ 0.48 CO2 1.54 N2 0.16

100

Pressure [bar]

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

80

C

60

40

20

0 0

100

200

300

400

Temperature [K]

Figure 5: Phase diagram of an Eagle Ford shale gas sample 76 calculated from the PC-SAFT equation of state. The mixture is in either gas or supercritical state at T = 328 K.

Prediction of shale gas in place It has been pointed out that the absorbed gas in the kerogen matrix makes a significant portion of the total gas in place. Shale gas is mostly stored in organic pores including interparticle and intraparticle pores. The interparticle pore sizes are on the order of microns and this porosity can be relatively easily measured. The intraparticle porosity is still un23

ACS Paragon Plus Environment

Energy & Fuels

3,000

12,000

graphite pore Density [mol/m3]

10,000

Density [mol/m3]

methane 8,000

6,000

butane

ethane

4,000

kerogen pore: dissolved

2,500

methane

2,000

1,500

1,000

butane ethane

propane

propane

hexane

500

2,000

CO2

hexane

N2

0

10

20

30

40

50

N2

CO2

0

0

0

10

Pressure [MPa]

20

30

40

(b) 10

14,000

kerogen pore: adsorbed

50

C1

40

9 12,000

C2

30

8

methane

C3

20

C4

10

10,000

Selectivity

7

8,000

6,000

4,000

50

Pressure [MPa]

(a)

Density [mol/m3]

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 41

C6

0 0

5

10

15

20

6

CO2 N2

5

solid line: graphite pore dashed: kerogen pore

4 3

butane

ethane

propane

2

2,000

1 0

CO2

hexane 0

10

20

30

N2 40

0 0

50

Pressure [MPa]

(c)

10

20

30

Pressure [MPa]

40

50

(d)

Figure 6: (a)Adsorption isotherms of shale gas in a graphite pore;(b) Dissolution of shale gas in kerogen matrix of Sample A3; (c) Adsorption of shale gas in pore space in Sample A3 formed 3 nm pore; (d) Selectivities of each component in a graphite pore v.s. porous kerogen composite with porosity equal to 15%. T = 328 K.

24

ACS Paragon Plus Environment

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

der investigation through different experimental techniques such as FIB-SEM, HIM. In the current work, we are modeling the intraparticle pores and estimating the relative amount of fluid stored in the intraparticle pores and kerogen matrix. Etminan et al. 77 measured the dissolved gas in kerogen and estimated it accounts for 22% of the total gas in place. Before making predictions of the gas in place, the model is firstly benchmarked with molecular simulation results of gas sorption in kerogen matrix so that this composite model along with parameters in Table 2 is further physically justified. Sun et al. 78 and Huang et al. 75 both simulated the CH4 and CO2 sorption in modeled kerogen matrices with kerogen molecules being fixed during simulation. As the kerogen pore size distribution is probed, it is found that as the kerogen density or maturity increases, the volume of the enterable pore, which is defined as pores larger than 4 Å, also increases. 75 In the work of Sun et al., 78 a mature kerogen of density 1.24 g/cm3 is used. To compare with their simulated CH4 and CO2 absolute loadings at various temperatures, we choose to use our model of the A3 sample, which has a density of 1.20 g/cm3. Without considering the polydispersity in pore size, we use a pore size of 4 Å to represent all the enterable pores. A nanoporosity of 0.08 is assumed by considering the increasing enterable pore volume as kerogen maturity and density increases, in which we refer to the measured enterable pore volume fractions of 0.011, 0.018 and 0.061 for kerogens of density 0.98 ± 0.01 g/cm3 , 1.11 ± 0.01 g/cm3 and 1.14 ± 0.02 g/cm3 by Huang et al.. 75 To make it consistent with the fixed matrix simulation, kerogen bulk density is also fixed to be the dry kerogen density in DFT calculation. It can be seen from Figure 7 that after tuning the porosity, predicted loadings of CH4 and CO2 from theory agree with molecular simulation semiquantitatively at all temperatures studied. The difference in slopes of isotherms should be attributed to the difference in treatment of solid. In molecular simulation, sorption of gas molecules is affected by fixed kerogen molecules while theory accounts for the rearrangement of kerogen as gas enters the kerogen matrix. Then shale gas in place is predicted with the bulk gas composition shown before. In Figure 8a, the mass fractions of dissolved gas is shown considering a nanoporosity of 15%

25

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

in Sample A1, A2 and A3. It can be seen the relative amount of dissolved gas decreases as the kerogen maturity increases. In Sample A3 formed composite system, the mass of hydrocarbons dissolved in the kerogen is approximately the same as the amount of hydrocarbons in nanosized pores. Heavy components (C3+ ) dissolve more in the kerogen matrix and light components are more stored in the nanopores. What is interesting is that we find CO2 is very much preferred by the kerogen matrix and makes around 80% of the total CO2 amount in a kerogen composite. The tremendous dissolution of CO2 in kerogen makes it attractive to displace shale gas with CO2 and simultaneously store CO2 in a depleted gas reservoir. In Figure 8b, the mass density of shale gas in three different phases is shown as a function of pressure. The density in the nanopores is always higher than that in the bulk phase due to surface adsorption and confinement effects and the dissolved gas density is much lower than the bulk density except at low pressures. After applying different porosity values, the total density in the kerogen particle is shown as a comparison to the bulk gas density. The amount of free gas in the interparticle pores depends on the macroporosity and the amount of adsorbed/dissolved gas depends on the volume fraction of kerogen. In a shale gas reservoir, the macropore volume and kerogen volume are of similar values, which implies the amount of gas inside kerogen particles makes around 30% of the total gas in place at high pressure, however, after depressurization, the mass of gas inside the kerogen can be greater than that in the interparticle pores. The crossover of the gas density in the kerogen particle with the bulk density curve is due to the strong adsorption/dissolution of heavy components. As we see the slope of ρ − P curve for the dissolved gas is much less steep than the bulk gas, which shows that the recovery factor of the dissolved gas from depressurization is very low. While in reality, the recovery factor of shale gas also depends largely on the pore connectivity and molecular diffusivity. Also the produced gas composition changes with depressurization, this may serve as a maximum recovery factor estimation that depends solely on ρ − P relations.

26

ACS Paragon Plus Environment

Page 26 of 41

Page 27 of 41

3.5

absolute sorption capacity [mmol/g]

2.5

absolute sorption capacity [mmol/g]

2.0

1.5

298 K 1.0

323 K 343 K 0.5

373 K

3.0 2.5 2.0 1.5

323 K

1.0

343 K 373 K

0.5 0.0

0.0

0

5

10

15

20

0

25

5

Pressure [MPa]

10

15

20

25

Pressure [MPa]

(a)

(b)

Figure 7: Comparison of theory predicted and molecular simulation data 79 of (a) pure CH4 (b) pure CO2 absolute sorption isotherms on kerogen. The theory assumes the parameters for the A3 sample kerogen of density 1.20 g/cm3 .

0%

20%

40%

60%

80%

100%

0.45 54% 50%

C1

0.40

41%

A2

61% 57%

C2

0.35

47%

A3 68%

C3

64% 53%

73%

C4

68% 56% 76%

C6

69% 53% 82% 81% 79%

CO2

0.30 0.25 0.20 0.15 0.10 0.05

45% 42%

N2

bulk fluid T = 328 K nanopore density dissolved density (0% porosity) total density (10% porosity) total density (20% porosity) total density (30% porosity)

A1

Density [g/cm3]

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

33%

0.00 0

10

20

30

40

50

Pressure [MPa]

(a)

(b)

Figure 8: (a) Mass fractions of dissolved gas of each component in different kerogens at T = 328 K, P = 50 MPa; (b) ρ − P relations of Eagle Ford shale gas in different phases.

27

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Conclusion In this work we proposed a new composite kerogen model to predict distribution of hydrocarbons and CO2 in nanoporous kerogen based on the PC-SAFT EoS and Flory-Huggins polymer solution theory. The kerogen matrix model is constructed by adding low molecular weight cross linkers to asphaltene molecules, which are recognized as fragments of kerogen. Parameters for kerogen are determined by applying correlations for asphaltene parameters and by fitting to dry kerogen density and swelling ratios of kerogen in different solvents. This new kerogen model is then incorporated in an inhomogeneous fluid model, which describes the nanoporous structure in kerogen particles, to study the hydrocarbon adsorption in nano-pores and dissolution in kerogen matrix. This work demonstrates a new method for calculating the original gas/oil in place and the reservoir fluid PVT relations, which are both important inputs in unconventional reservoir simulation. The major results and conclusions can be summarized as follows: • The model describes the swelling ratios of kerogen of different thermal maturity in normal, branched, cyclic alkanes and aromatic molecules in quantitative agreement with experimental measurements. A large swelling ratio corresponds to a low crosslinking density and high petroleum expulsion potential in less mature kerogen. • Fractionation of molecules from petroleum source kerogen to nanopores and bulk reservoir fluid is studied. Distributions of molecules with similar molecular weight but different chemical structures at the kerogen-fluid interface are shown. Of the adsorbed and dissolved phases, strong deviation in fluid composition from bulk phase is seen in the dissolved phase due to different solubilities of different components, slight deviation of fluid composition in nanopores from bulk phase is attributed to surface adsorption differences on the kerogen solid. • The total sorption of CH4 and CO2 is calculated using this nanoporous kerogen model and semiquantitative agreement is shown in comparison with molecular simulation by 28

ACS Paragon Plus Environment

Page 28 of 41

Page 29 of 41 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

only tuning the kerogen porosity.Adsorption/absorption isotherms of a shale gas mixture are studied and compared with adsorption in a graphite pore. Heavy components (C3+) are strongly retained by the nanopores and kerogen matrix especially in low pressures. A reduction in gas fractionation is seen for nanoporous kerogen among alkanes, however, CO2 has a noticeable increase in partitioning to organic matrix. • Analysis of the amount of each component in different phases shows a significant amount is stored in kerogen matrix, but recovery of the dissolved gas by pressure depletion is expected to be much more difficult than adsorbed gas and free gas. The whole range of parameters in this composite model include (I) parameters in PCSAFT modeling of kerogen (Table 2 and Eqn 16 - 18); (II) PC-SAFT parameters for solvents and their dispersion interaction parameters with kerogen (e.g. Table 1); (III) binary interaction parameters between solvents (e.g. Table 3). (II) and (III) are relatively available from literature on PC-SAFT EoS, while (I) rely largely on the available experimental data since kerogen does not have a defined chemical formula and structure. Kerogen density and solvent swelling measurements are needed to feed enough information to estimate the kerogen parameters. Once parameters are determined, elemental analysis information such as hydrogen index (HI) and fraction of aromatic carbons help justify the validity of parameters. To avoid an extensive complexity, we have made great simplifications with respect to the real material. These include assuming constant values for the solid roughness and the pore size. No specific interactions (e.g., hydrogen bonding sites) have yet been included in our model. Also a slit-shaped pore is not enough to account for the pore connectivity in a three dimensional structure. However, this model still shows reasonable agreement with molecular simulation of gas sorption in kerogen. Once the physical parameters for a particular sample are established, the model can be utilized to characterize the pore size distribution (PSD) based on measured adsorption isotherms. Then the modeling approach will serve as a validation for multicomponent experiments, in which high measurement accuracy is difficult to achieve. In future work, the kerogen model can be extended to have association sites, given 29

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

that data on the functional group concentration and kerogen swelling in hydrogen bonding solvents is available. The water partitioning to clays and organic matter after hydraulic fracturing can then be predicted with association interactions taken into account.

Supporting Information A. Equations for free energy functionals; B. Procedures for solving the fluid density distribution.

Nomenclature A = Helmholtz free energy (total, hs, ch, disp and elast) per mole of molecules kb = Boltzmann’s constant ≈ 1.381 × 10−23 J/K T = Temperature, K Φ = cross-linked polymer network functionality nchain = number of polymer chains between two cross-linkers V = molar volume of each polymer chain when swollen by solvents Vmax = molar maximum volume of a crosslinked polymer network when swollen by solvents V0 = molar volume of a dry polymer chain when not mixed with solvents ncross−linker = number of cross-linker molecules NA = Avogadro’s number ≈ 6.02 × 1023 molecules/mol mp = number of monomer segments in each polymer chain σp = diameter of segments in a polymer chain, Å mcross−linker = mass of linker molecules added to 1 mole of polymer, g Mcross−linker = molecular weight of linker molecules, g/mol Mc = average molecular weight of polymer chains between two linkers, g/mol mpolymer = mass of the polymer before linker molecules are added, g p = pressure 30

ACS Paragon Plus Environment

Page 30 of 41

Page 31 of 41 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

µ = chemical potential N = total number of solvent species ρs (z) = density distribution of segments in the kerogen matrix in the z direction ρs bulk = equilibrium density of kerogen in kerogen/solvent phase z0 = position of the solid-fluid interface where ρs (z0 ) = ρs bulk /2 δ = solid roughness, Å Ω = grand potential r = spacial coordinate of a segment Subscripts i = for solvent component i s = for solid or kerogen sf = for solid and fluid Superscripts ideal = ideal gas hs = hard sphere ch = chain disp = dispersion elast = elastic I = solvent phase II = kerogen/solvent phase

Acknowledgement The authors thank the Robert Welch Foundation (Grant No. C-1241) and the American Chemical Society Petroleum Research Fund (No. ACS-PRF58859-ND6) for financial support. The authors also thank Prof. George Hirasaki, Dr. Philip Singer, Dr. Dilip Asthagiri,

31

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Zeliang Chen and Arjun Valiya Parambathu for insightful discussions.

32

ACS Paragon Plus Environment

Page 32 of 41

Page 33 of 41 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

References (1) Annual Energy Outlook, www.eia.gov/aeo; 2018. (2) EIA/ARI World Shale Gas and Shale Oil Resource Assessment, Shale gas and shale oil resource assessment methodology; 2013. (3) Ross, D. J.; Bustin, R. M. The importance of shale composition and pore structure upon gas storage potential of shale gas reservoirs. Marine and Petroleum Geology 2009, 26, 916–927. (4) Jin, Z.; Firoozabadi, A. Thermodynamic modeling of phase behavior in shale media. SPE Journal 2016, 21, 190–207. (5) Ertas, D.; Kelemen, S. R.; Halsey, T. C. Petroleum expulsion part 1. Theory of kerogen swelling in multicomponent solvents. Energy & Fuels 2006, 20, 295–300. (6) Kelemen, S. R.; Walters, C. C.; Ertas, D.; Freund, H.; Curry, D. J. Petroleum expulsion part 3. A model of chemically driven fractionation during expulsion of petroleum from kerogen. Energy & Fuels 2006, 20, 309–319. (7) Kelemen, S. R.; Walters, C. C.; Ertas, D.; Kwiatek, L. M.; Curry, D. J. Petroleum expulsion part 2. Organic matter type and maturity effects on kerogen swelling by solvents and thermodynamic parameters for kerogen from regular solution theory. Energy & Fuels 2006, 20, 301–308. (8) King, H. E.; Eberle, A. P.; Walters, C. C.; Kliewer, C. E.; Ertas, D.; Huynh, C. Pore architecture and connectivity in gas shale. Energy & Fuels 2015, 29, 1375–1390. (9) Chen, Z.; Singer, P. M.; Kuang, J.; Vargas, F. M.; Hirasaki, G. J. Effects of bitumen extraction on the 2D NMR response of saturated kerogen isolates. Petrophysics 2017, 58, 470–484.

33

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(10) Curtis, M. E.; Cardott, B. J.; Sondergeld, C. H.; Rai, C. S. Development of organic porosity in the Woodford Shale with increasing thermal maturity. Int. J. Coal Geol. 2012, 103, 26–31. (11) Barsotti, E.; Saraji, S.; Tan, S. P.; Piri, M. Capillary Condensation of Binary and Ternary Mixtures of n-Pentane–Isopentane–CO2 in Nanopores: An Experimental Study on the Effects of Composition and Equilibrium. Langmuir 2018, 34, 1967–1980. (12) Barsotti, E.; Tan, S. P.; Piri, M.; Chen, J.-H. Phenomenological Study of Confined Criticality: Insights from the Capillary Condensation of Propane, n-Butane, and nPentane in Nanopores. Langmuir 2018, 34, 4473–4483. (13) 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, 33, 11189–11202. (14) Didar, B.; Akkutlu, I. Pore-size dependence of fluid phase behavior and the impact on shale gas reserves. Unconventional Resources Technology Conference. 2013; pp 1793– 1814. (15) Liu, X.-Y.; Li, J.-T.; Gu, F.; Wang, H.-J. Phase Equilibria of Hydrogen Bonding Fluid in a Slit Pore with Broken Symmetry. Chinese J. Chem. Phys. 2015, 28 . (16) Pizio, O.; Sokołowski, S.; Sokołowska, Z. Phase behavior of binary symmetric mixtures in pillared slit-like pores: A density functional approach. J. Chem. Phys. 2011, 134, 214702. (17) Dong, X.; Liu, H.; Hou, J.; Wu, K.; Chen, Z. Phase equilibria of confined fluids in nanopores of tight and shale rocks considering the effect of capillary pressure and adsorption film. Ind. Eng. Chem. Res. 2016, 55, 798–811.

34

ACS Paragon Plus Environment

Page 34 of 41

Page 35 of 41 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

(18) Pathak, M.; Kweon, H.; Deo, M.; Huang, H. Kerogen Swelling and Confinement: Its implication on Fluid Thermodynamic Properties in Shales. Sci. Rep. 2017, 7, 1–14. (19) Siskin, M.; Brons, G.; Payack, J. F. Disruption of kerogen-mineral interactions in Rundle Ramsay Crossing oil shale. Energy & Fuels 1989, 3, 108–109. (20) Larsen, J. W.; Parikh, H. M.; Michels, R.; Raoult, N.; Pradier, B. Kerogen macromolecular structure. Prepr Symp Am Chem Soc Div Fuel Chem. 2000; pp 211–5. (21) Larsen, J. W.; Li, S. Solvent Swelling Studies of Green River Kerogen. Energy & Fuels 1994, 8, 932–936. (22) Larsen, J. W.; Li, S. An initial comparison of the interactions of Type I and III kerogens with organic liquids. Org. Geochem. 1997, 26, 305–309. (23) Larsen, J. W.; Parikh, H.; Michels, R. Changes in the cross-link density of Paris Basin Toarcian kerogen during maturation. Org. Geochem. 2002, 33, 1143–1152. (24) Ritter, U. Fractionation of petroleum during expulsion from kerogen. J. Geochemical Explor. 2003, 79, 417–420. (25) Savest, N.; Hruljova, J.; Oja, V. Characterization of thermally pretreated kukersite oil shale using the solvent-swelling technique. Energy & Fuels 2009, 23, 5972–5977. (26) Hruljova, J.; Savest, N.; Oja, V.; Suuberg, E. M. Kukersite oil shale kerogen solvent swelling in binary mixtures. Fuel 2013, 105, 77–82. (27) Kilk, K.; Savest, N.; Hruljova, J.; Tearo, E.; Kamenev, S.; Oja, V. Solvent swelling of Dictyonema oil shale. Oil Shale 2010, 27, 26–36. (28) Hildebrand, J. H.; Prausnitz, J. M.; Scott, R. L. Regular and related solutions: the solubility of gases, liquids, and solids; Van Nostrand Reinhold Co., 1970.

35

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(29) Flory, P. J.; Rehner Jr, J. Statistical mechanics of cross-linked polymer networks I. Rubberlike elasticity. J. Chem. Phys. 1943, 11, 512–520. (30) Kovac, J. Modified Gaussian model for rubber elasticity. Macromolecules 1978, 11, 362–365. (31) Barr-Howell, B. D.; Peppas, N. A. Importance of junction functionality in highly crosslinked polymers. Polymer Bulletin 1985, 13, 91–96. (32) Chapman, W. G.; Gubbins, K. E.; Jackson, G.; Radosz, M. New reference equation of state for associating liquids. Ind. Eng. Chem. Res. 1990, 29, 1709–1721. (33) Gross, J.; Sadowski, G. Perturbed-Chain SAFT: An Equation of State Based on a Perturbation Theory for Chain Molecules. Ind. Eng. Chem. Res. 2001, 40, 1244–1260. (34) Vargas, F. M.; Gonzalez, D. L.; Hirasaki, G. J.; Chapman, W. G. Modeling asphaltene phase behavior in crude oil systems using the perturbed chain form of the statistical associating fluid theory (PC-SAFT) equation of state. Energy & Fuels 2009, 23, 1140– 1146. (35) Tumakaka, F.; Gross, J.; Sadowski, G. Thermodynamic modeling of complex systems using PC-SAFT. Fluid Phase Equilibr. 2005, 228, 89–98. (36) Gross, J.; Spuhl, O.; Tumakaka, F.; Sadowski, G. Modeling copolymer systems using the perturbed-chain SAFT equation of state. Ind. Eng. Chem. Res. 2003, 42, 1266–1274. (37) Gushterov, N.; Doghieri, F.; Quitmann, D.; Niesing, E.; Katzenberg, F.; Tiller, J. C.; Sadowski, G. VOC Sorption in Stretched Cross-Linked Natural Rubber. Ind. Eng. Chem. Res. 2016, 55, 7191–7200. (38) Arndt, M. C.; Sadowski, G. Modeling poly (N-isopropylacrylamide) hydrogels in water/alcohol mixtures with PC-SAFT. Macromolecules 2012, 45, 6686–6696.

36

ACS Paragon Plus Environment

Page 36 of 41

Page 37 of 41 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

(39) Tumakaka, F.; Gross, J.; Sadowski, G. Modeling of polymer phase equilibria using Perturbed-Chain SAFT. Fluid Phase Equilib. 2002, 194-197, 541–551. (40) Chen, Y.; Mutelet, F.; Jaubert, J.-N. Modeling the solubility of carbon dioxide in imidazolium-based ionic liquids with the PC-SAFT equation of state. The Journal of Physical Chemistry B 2012, 116, 14375–14388. (41) Stoychev, I.; Galy, J.; Fournel, B.; Lacroix-Desmazes, P.; Kleiner, M.; Sadowski, G. Modeling the phase behavior of PEO- PPO- PEO surfactants in carbon dioxide using the PC-SAFT equation of state: application to dry decontamination of solid substrates. Journal of Chemical & Engineering Data 2009, 54, 1551–1559. (42) Li, X.-S.; Wu, H.-J.; Englezos, P. Prediction of gas hydrate formation conditions in the presence of methanol, glycerol, ethylene glycol, and triethylene glycol with the statistical associating fluid theory equation of state. Ind. Eng. Chem. Res. 2006, 45, 2131–2137. (43) Cameretti, L. F.; Sadowski, G.; Mollerup, J. M. Modeling of aqueous electrolyte solutions with perturbed-chain statistical associated fluid theory. Ind. Eng. Chem. Res. 2005, 44, 3355–3362. (44) David Ting, P.; Hirasaki, G. J.; Chapman, W. G. Modeling of asphaltene phase behavior with the SAFT equation of state. Petroleum Science and Technology 2003, 21, 647–661. (45) Gonzalez, D. L.; Hirasaki, G. J.; Creek, J.; Chapman, W. G. Modeling of asphaltene precipitation due to changes in composition using the perturbed chain statistical associating fluid theory equation of state. Energy & Fuels 2007, 21, 1231–1242. (46) Gonzalez, D. L.; Vargas, F. M.; Hirasaki, G. J.; Chapman, W. G. Modeling study of CO2-induced asphaltene precipitation. Energy & Fuels 2007, 22, 757–762.

37

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(47) Vargas, F. M.; Gonzalez, D. L.; Hirasaki, G. J.; Chapman, W. G. Modeling asphaltene phase behavior in crude oil systems using the perturbed chain form of the statistical associating fluid theory (PC-SAFT) equation of state. Energy & Fuels 2009, 23, 1140– 1146. (48) Panuganti, S. R.; Vargas, F. M.; Gonzalez, D. L.; Kurup, A. S.; Chapman, W. G. PCSAFT characterization of crude oils and modeling of asphaltene phase behavior. Fuel 2012, 93, 658–669. (49) Arndt, M. C. Thermodynamic modelling of hydrogels; Ph.D. thesis at Laboratory of Thermodynamics, Technische Univerität Dortmund, 2014. (50) Gushterov, N. Thermodynamic and Mechanical Characterization of Shape-Memory Natural Rubber ; Ph.D. thesis at Laboratory of Thermodynamics, Technische Univerität Dortmund, 2017. (51) Flory, P. J. Principles of polymer chemistry; Cornell University Press, 1953. (52) Miao, B.; Vilgis, T. A.; Poggendorf, S.; Sadowski, G. Effect of finite extensibility on the equilibrium chain size. Macromol. Theory Simulations 2010, 19, 414–420. (53) Tripathi, S.; Chapman, W. G. Microstructure of inhomogeneous polyatomic mixtures from a density functional formalism for atomic mixtures. J. Chem. Phys. 2005, 122, 094506. (54) Jain, S.; Jog, P.; Weinhold, J.; Srivastava, R.; Chapman, W. G. Modified interfacial statistical associating fluid theory: Application to tethered polymer chains. J. Chem. Phys. 2008, 128, 154910. (55) Emborsky, C. P.; Feng, Z.; Cox, K. R.; Chapman, W. G. Recent advances in classical density functional theory for associating and polyatomic molecules. Fluid Phase Equilib. 2011, 306, 15–30. 38

ACS Paragon Plus Environment

Page 38 of 41

Page 39 of 41 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

(56) Wang, L.; Haghmoradi, A.; Liu, J.; Xi, S.; Hirasaki, G. J.; Miller, C. A.; Chapman, W. G. Modeling micelle formation and interfacial properties with iSAFT classical density functional theory. J. Chem. Phys. 2017, 146, 124705. (57) Marshall, B. D.; Cox, K. R.; Chapman, W. G. Supramolecular assembly and surfactant behavior of triblock rod-coil amphiphiles at liquid interfaces using classical density functional theory. Soft Matter 2012, 8, 7415–7425. (58) Gong, K.; Marshall, B. D.; Chapman, W. G. Modeling lower critical solution temperature behavior of associating polymer brushes with classical density functional theory. J. Chem. Phys. 2013, 139, 094904. (59) Zhang, Y.; Valiya Parambathu, A.; Chapman, W. G. Density functional study of dendrimer molecules in solvents of varying quality. J. Chem. Phys. 2018, 149, 064904. (60) Feng, Z.; Chapman, W. G. Revisited block copolymer/nanoparticle composites: Extension of interfacial statistical associating fluid theory. Macromolecules 2012, 45, 6658– 6668. (61) Ballal, D.; Chapman, W. G. Hydrophobic and hydrophilic interactions in aqueous mixtures of alcohols at a hydrophobic surface. J. Chem. Phys. 2013, 139, 114706. (62) Gonzalez, D. L.; Ting, P. D.; Hirasaki, G. J.; Chapman, W. G. Prediction of asphaltene instability under gas injection with the PC-SAFT equation of state. Energy & Fuels 2005, 19, 1230–1234. (63) Ungerer, P.; Collell, J.; Yiannourakou, M. Molecular Modeling of the Volumetric and Thermodynamic Properties of Kerogen : In fl uence of Organic Type and Maturity. Energy & Fuels 2015, 29, 91–105. (64) Ravikovitch, P. I.; Neimark, A. V. Density functional theory model of adsorption on amorphous and microporous silica materials. Langmuir 2006, 22, 11171–11179. 39

ACS Paragon Plus Environment

Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(65) Gor, G. Y.; Thommes, M.; Cychosz, K. A.; Neimark, A. V. Quenched solid density functional theory method for characterization of mesoporous carbons by nitrogen adsorption. Carbon N. Y. 2012, 50, 1583–1590. (66) Sauer, E.; Gross, J. Classical density functional theory for liquid–fluid interfaces and confined systems: A functional for the perturbed-chain polar statistical associating fluid theory equation of state. Ind. Eng. Chem. Res. 2017, 56, 4119–4135. (67) Van Krevelen, D. Graphical-statistical method for the study of structure and reaction processes of coal. Fuel 1950, 29, 269–284. (68) Punnapala, S.; Vargas, F. M. Revisiting the PC-SAFT characterization procedure for an improved asphaltene precipitation prediction. Fuel 2013, 108, 417–429. (69) Pawar, G.; Meakin, P.; Huang, H. Reactive Molecular Dynamics Simulation of Kerogen Thermal Maturation and Cross-Linking Pathways. Energy & Fuels 2017, 31, 11601– 11614. (70) Lille, Ü.; Heinmaa, I.; Pehk, T. Molecular model of Estonian kukersite kerogen evaluated by 13C MAS NMR spectraâŸĘ. Fuel 2003, 82, 799–804. (71) Tinni, A.; Perez, F.; Devegowda, D.; Truong, T.; Dang, S.; Sondergeld, C.; Rai, C. InSitu Fractionation in Liquids-Rich Shales and its Implications for EOR: Experimental Verification and Modeling Study. Unconventional Resources Technology Conference. 2018; pp 2451–2461. (72) Milliken, K. L.; Rudnicki, M.; Awwiller, D. N.; Zhang, T. Organic matter–hosted pore system, Marcellus formation (Devonian), Pennsylvania. AAPG bulletin 2013, 97, 177– 200. (73) Zargari, S.; Canter, K. L.; Prasad, M. Porosity evolution in oil-prone source rocks. Fuel 2015, 153, 110–117. 40

ACS Paragon Plus Environment

Page 40 of 41

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

(74) Takbiri-Borujeni, A.; Kazemi, M.; Sun, T.; Mansouri-Boroujeni, M. Effect of Kerogen Type and Maturity on Performance of Carbon Dioxide Storage in Shale. SPE Annual Technical Conference and Exhibition. 2017. (75) Huang, L.; Ning, Z.; Wang, Q.; Zhang, W.; Cheng, Z.; Wu, X.; Qin, H. Effect of organic type and moisture on CO 2/CH 4 competitive adsorption in kerogen with implications for CO 2 sequestration and enhanced CH 4 recovery. Applied Energy 2018, 210, 28–43. (76) George, D. L.; Bowles, E. B. Shale gas measurement and associated issues. Pipeline & Gas Journal 2011, 238 . (77) Etminan, S. R.; Javadpour, F.; Maini, B. B.; Chen, Z. Measurement of gas storage processes in shale and of the molecular diffusion coefficient in kerogen. International Journal of Coal Geology 2014, 123, 10–19. (78) Sun, H.; Zhao, H.; Qi, N.; Li, Y. Molecular Insights into the Enhanced Shale Gas Recovery by Carbon Dioxide in Kerogen Slit Nanopores. The Journal of Physical Chemistry C 2017, 121, 10233–10241. (79) Sun, H.; Zhao, H.; Qi, N.; Li, Y. Molecular Insights into the Enhanced Shale Gas Recovery by Carbon Dioxide in Kerogen Slit Nanopores. J. Phys. Chem. C 2017, 121, 10233–10241.

41

ACS Paragon Plus Environment