Molecular structures of kerogen control hydrocarbon production in unconventional reservoirs. Significant progress has been made in developing model representations of various kerogen structures. These models have been widely used for the prediction of gas adsorption and migration in shale matrix. However, using density functional perturbation theory (DFPT) calculations and vibrational spectroscopic measurements, we here show that a large gap may still remain between the existing model representations and actual kerogen structures, therefore calling for new model development. Using DFPT, we calculated Fourier transform infrared (FTIR) spectra for six most widely used kerogen structure models.
Philippe F. Weck1, Eunja Kim2, Yifeng Wang1, Jessica N. Kruichak1, Melissa M. Mills1, Edward N. Matteo1 & Roland J.-M. Pellenq3,4,5
- Received: Accepted:Published online:
- © The Author(s) 2017
Limited agreement between the model predictions and the measurements highlights that the existing kerogen models may still miss some key features in structural representation. A combination of DFPT calculations with spectroscopic measurements may provide a useful diagnostic tool for assessing the adequacy of a proposed structural model as well as for future model development. This approach may eventually help develop comprehensive infrared (IR)-fingerprints for tracing kerogen evolution.
Kerogen is a high-molecular weight, carbonaceous polymer material resulting from the condensation of organic residues in sedimentary rocks; such organic constituent is insoluble either in aqueous solvents or in common organic solvents1. In addition to its carbon skeleton, kerogen is also predominantly made of hydrogen and oxygen, as well as residual nitrogen and sulfur. Although kerogen plays a central role in hydrocarbon (i.e., crude oil and natural gas) production from source rocks in geologic environments2, the crucial interplay between its complex nanoscale structure and its properties has only been revealed in recent years3,4,5,6,7.
A detailed understanding of the various kerogen structures associated with different oil- or gas-prone kerogen types and maturity is particularly important to understanding nanofluidic processes underlying crude oil and natural gas extraction7. During its diagenesis, catagenesis, and metagenesis maturation stages, kerogen undergoes successive transformations concomitant with a gradual increase in sp2/sp3 carbon hybridization ratio3,5 before eventually reaching full condensation of its polynuclear aromatic units to form sp2-hybridized graphitic structures as components of a complex network of carbon nanopores.
Existing model structures of kerogen were generally constrained from elemental analyses and functional group data obtained from X-ray photoelectron spectroscopy (XPS) and 13C nuclear magnetic resonance analyses (NMR)4. Recently, Bousige et al. developed kerogen structure models using a molecular dynamics-hybrid reverse Monte-Carlo (MD-HRMC) method, which minimizes the configurational energy while constraining possible molecular configurations from the pair distribution function obtained from inelastic neutron scattering measurements5.
Direct comparison between vibrational properties measured by neutron, Raman or infrared spectroscopies and spectra simulated from classical or quantum mechanical methods has proven to be a powerful tool for model validation5, 8. In this study, infrared signatures for a variety of generic kerogen models proposed recently4, 5 were calculated within the framework of density functional perturbation theory (DFPT) and compared to Fourier transform infrared (FTIR) spectra collected from kerogen samples with different types and maturity. Computed infrared (IR) signatures were also compared to generalized phonon densities of states (GDOS) obtained from recent inelastic neutron scattering experiments for various kerogen samples5. An attempt is made in this study to critically assess the adequacy of existing model representations of kerogen structures, using characteristic IR signatures as a diagnostic tool. We will show that a combination of DFPT calculations with vibrational spectroscopic measurements, especially FTIR measurements, can be a useful tool for developing more realistic models for kerogen structure representation. This approach may eventually lead to the development of a comprehensive IR-fingerprints for tracing kerogen evolution. Note that IR spectroscopy was successfully used in previous investigations to differentiate crude oil from residual fuel oil samples9, 10.
Results and Discussion
The 3D-periodic kerogen structures initially considered in this study were based on generic molecular fragments proposed by Ungerer et al.4 to represent kerogen structures with lacustrine (type I) and marine (type II) depositional origin.
Type I kerogen oil shales are typically the most promising deposits for conventional oil retorting and type II kerogen is found in numerous oil shale deposits. For that reason, the focus of this study will be mainly on type-I and -II kerogen. The models utilized (cf. Fig. 1) contained 118 atoms (type-I) and 89 atoms (type-II) per simulation cell. In terms of oxygen-to-carbon (O/C) and hydrogen-to-carbon (H/C) atomic ratios, these type I (H/C = 1.29; O/C = 0.06) and type II (H/C = 0.57; O/C = 0.06) models correspond to kerogen in upper- and lower-left quadrants of the van Krevelen diagram, mapping H/C as a function of O/C to represent the chemical evolution of kerogen11, 12.
Figure 1 Structures of 3D-periodic type-I (left) and type-II (right) kerogen models used in DFT/DFPT calculations at the GGA/PBE level of theory. The simulation cells are indicated by solid lines. Color legend: grey, C; white, H; purple, N; red, O; yellow, S.
Kerogen maturation over geologic times is accompanied by a decrease of H/C and O/C atomic ratios. In terms of vitrinite reflectance R0, another common maturity indicator, these type I and II models correspond to ~0.6 and ~1.5%R0, respectively; higher values of %R0 indicate higher maturity of the samples.
Typical bands of interest for spectral signatures of kerogen include mainly the aliphatic band (2800–3000 cm−1), aromatic bands (700–900 cm−1, ~1600 cm−1, and 3000–3100 cm−1) and carbonyl and carboxyl absorption band (~1700 cm−1). Results from DFPT linear response and Born effective charge (BEC) calculations carried out at the generalized gradient approximation (GGA) level of theory with the parameterization of Perdew, Burke, and Ernzerhof (PBE) for those type-I and type-II models are shown in Fig. 2, along with FTIR spectra collected for Woodford, Marcellus, and Mancos kerogen samples.
The Mancos sample is a gas-prone, type-III kerogen (humic origin), with a vitrinite reflectance of 0.65%R0. The Marcellus sample is composed of mainly Type-II kerogen, with an estimated vitrinite reflectance value of 2.51%R0. While the simulated frequencies and intensities vary in general from experimental values, agreement between DFPT results and experimental data is generally slightly better using the type-I model than with its type-II counterpart. In particular, compared to all experimentally characterized samples, the type-II model tends to overestimate IR absorbance in the 1200–1400 cm−1 and 1700–2200 cm−1 ranges, while the type-I model simulations feature intensity ratios more in line with experiments.
Therefore, details of the present vibrational analysis given below will focus on the type-I model, although results for the type-II model will also be discussed.
Figure 2 Infrared spectra simulated from density functional perturbation theory (DFPT) at the GGA/PBE level for the type-I and type-II models shown in Fig. 1 and Fourier transform infrared spectra collected for Mancos, Marcellus and Woodford kerogen samples.
In the 2600–2800 cm−1 range, both type-I and type-II models display a similar double-peak structure characteristic of C–H stretch, with the asymmetric and symmetric C–H stretches corresponding mainly to the highest and lowest-frequency peaks, respectively; C–H stretches from carbon rings are also in this frequency domain. The simulated frequencies are redshifted compared to measured symmetric and asymmetric C–H stretches centered around ~2850 and ~2920 cm−1, respectively, for the Woodford and Mancos samples.
This can be ascribed in part to the inability of standard DFT/DFPT to describe accurately long-range van der Waals interactions resulting from dynamical correlations between fluctuating charge distributions and to the fact that the type-I and type-II models used in this study might be more representative of isolated molecular fragments than confined, entangled molecular structures forming kerogen bulk; in the latter case, a blueshift of the simulated C–H stretch frequencies with respect to the present simulations would be expected. It is also worth noting that the type-II model used here still has significant IR absorbance in the 2600–2800 cm−1 range, suggesting that the number of C–H bonds in the model might not be representative of very mature samples with high level of carbon condensation such as Marcellus.
For the type-I model, the highest-frequency mode for ring C–C bond stretching is at 1743 cm−1 (2174 cm−1 in the type-II model), followed by an intense series of peaks below ~1595 cm−1 originating from combinations of ring C–C bond vibration and C–H twisting. Series of H–C–H scissoring and/or rocking modes are predicted around 1430–1495 cm−1. These ring C–C bond vibration/C–H twisting and H–C–H scissoring and/or rocking modes are most likely responsible for the FTIR peaks centered at ~1600 and 1450 cm−1, respectively. The fact that the peak around ~1600 cm−1 subsists for all levels of maturity of the samples, while it becomes sharper and experiences a slight redshift with kerogen maturity, confirms that it is due to ring C–C bond vibration (and becomes more pronounced and more uniform with carbon condensation into larger, regular 2D carbon structures). Since the peak centered around ~1450 cm−1 for Woodford and Mancos samples tends to fade away in the Marcellus sample, this also indicates that significant H–C–H scissoring and/or rocking modes have been suppressed as a result of dehydrogenation and condensation in the most mature samples. A smaller peak observed around ~1370 cm−1 for Woodford and Mancos samples and absent for the Marcellus sample is assigned, based on our type-I calculations, to combinations of CH3 “umbrella” modes and H–C–H wagging modes. Vibrational modes dominating in the ~960–1250 cm−1 range are complex combinations of H–C–H twisting, wagging and rocking, C–H and CH3 wagging, and ring C–C stretch. In particular, calculations suggest that the double-peak structure around ~1030 and ~1090 cm−1 for the Mancos sample (and to some extent the Woodford sample featuring a broad peak centered around ~1070 cm−1) stems predominantly from C–C stretching in hydrocarbon chains and C–C stretching in carbon rings; this would explain why the former peak fades away as hydrocarbon chains disappear in the mature Marcellus sample, while the latter peak subsists in mature samples, although narrower and slightly redshifted as a result of kerogen condensation.
The sharp peaks in the vicinity of ~780 and ~800 cm−1 are assigned to C–C–C ring bending modes and C–O stretching modes in chains and rings, respectively, and the peak near ~740 cm−1 in the Mancos sample is assigned to ring C–N/C–C vibrations. The peak near ~695 cm−1 in both Woodford and Mancos samples is also assigned to C–O stretching modes in chains; this peak tends to disappear in more mature samples. At lower frequencies, a number of complex combinations of carbon skeletal modes from rings and chains contribute to the IR absorbance.
The peak observed at ~520 cm−1 in both Woodford and Mancos samples (as well as in the Marcellus sample to some extent) originates from C–S stretching. The modes below this frequency are due essentially to complex, collective ring vibrational modes.
In addition to the type-I and -II models proposed by Ungerer et al.4, complex kerogen structures based on the sophisticated EFK, MEK, MarK and PY02 models (cf. Fig. 3) recently proposed by Bousige and co-workers5 were also used to simulate IR spectra. These four molecular models with a density of 1.2 g/cm3 were built to represent an immature type-II kerogen from the carbonate-rich Eagle Ford Play (EFK; H/C = 1.19; O/C = 0.10; 0.65%R0), an immature sulphur-rich type-II kerogen from the Middle East (MEK; H/C = 1.45; O/C = 0.08; 0.55%R0), a mature type II kerogen from the clay-rich Marcellus Play (MarK; H/C = 0.46; O/C = 0.11; 2.20%R0), and a shungite from Russia (PY02; H/C = 0.09; O/C = 0.01) sharing many features with excessively mature kerogen5.
Figure 3 Structures of the EFK, MEK, MarK and PY02 models (Top; cubic box size of 50 × 50 × 50 Å3); representative 3D-periodic portions of EFK, MEK, MarK and PY02 models used in the present DFT/DFPT calculations at the GGA/PBE level (Bottom); simulations cells are indicated by solid lines). Color legend: grey, C; white, H; red, O.
The 3D-periodic portions of EFK, MEK, MarK and PY02 models used in the present DFT/DFPT calculations are also displayed. These smaller models have somewhat different H/C and O/C ratios (see Fig. 4) and lower densities than the samples characterized experimentally (i.e., EFK: H/C = 1.45, O/C = 0.03, ρ = 0.7 g/cm3; MEK: H/C = 1.55, O/C = 0.06, ρ = 0.8 g/cm3; MarK: H/C = 0.54, O/C = 0.04, ρ = 0.6 g/cm3; PY02: H/C = 0.13, O/C = 0.01, ρ = 0.5 g/cm3).
The models used in DFT/DFPT calculations (cf. Fig. 3) contained 174 (EFK), 165 (MEK) and 119 (MarK and PY02) atoms per simulation cell. Experimental characteristics of the EFK, MEK and MarK kerogen samples and of the Mancos, Woodford and Marcellus kerogen samples discussed above are summarized in Table 1.
Table 1: Experimental characteristics of the Mancos, Woodford, Marcellus, EFK, MEK and MarK kerogen samples(1).
Figure 5 shows the infrared spectra simulated using DFPT for the representative 3D-periodic portions of the EFK, MEK, MarK and PY02 models depicted in Fig. 3, along with the FTIR spectra of the Mancos, Marcellus and Woodford kerogen samples.
Figure 4van Krevelen diagram displaying H/C versus O/C ratios for kerogen structures of the sophisticated EFK, MEK, MarK and PY02 models from Bousige et al. (ref. 5) and their simplified variants used in the present DFPT calculations, as well as the type-I and -II structures proposed by Ungerer et al. (ref. 4). Isovalue contour lines of the vitrinite reflectance R0, another common maturity indicator, are also represented, along with typical domains (cyan) for kerogen samples with different depositional origins (types I to IV). Maturation increases with decreasing H/C and O/C ratios.
A detailed analysis of the generalized phonon densities of states was carried out for the EFK, MEK, MarK and PY02 samples by Bousige et al. using force-field molecular dynamics simulations and inelastic neutron scattering measurements5. All vibrational eigenmodes are accessible by neutron spectroscopy, unlike in IR spectroscopy where selections rules apply. Therefore, the complex vibrational characteristics of these models will be only discussed succinctly hereafter. Nevertheless, as shown in ref. 5, the GDOS measured by neutron spectroscopy (cf. Fig. 6) are dominated by vibrational modes of hydrogen atoms forming bonds with carbon by s-sp2 or s-sp3 overlap.
While such neutron spectroscopy measurements probing C–H bonding are useful to estimate maturation in terms of sp2/sp3 or H/C ratios, IR spectroscopy can also provide valuable information in terms of C–O, C– S, C–N, C–C and complex carbon-backbone vibrational modes, as discussed above for the analysis of type-I and -II models of Ungerer et al.4. Such information is crucial for constraining functional group distributions and their chemical bonding environments in kerogen.
Figure 5 Infrared spectra simulated from density functional perturbation theory (DFPT) at the GGA/PBE level for representative portions of the EFK, MEK, MarK and PY02 models displayed in Fig. 3 and Fourier transform infrared spectra collected for Mancos, Marcellus and Woodford kerogen samples.
Note that models used in this study cover a wide range of parameter space in the van Krevelen diagram (Fig. 4). We would expect that the three kerogen samples we studied should be reasonably represented by some of these models. However, the comparison of the calculated with the measured FTIR spectra (Figs 2 and 5) shows only a limited agreement between model predictions and experimental measurements, thus highlighting a possible deficiency of the existing models in the representation of actual kerogen structures, especially for functional group distributions in the material.
Previous work shows that nanometer-scale pore structures and pore surface functional groups are the two most important attributes controlling hydrocarbon sorption and release7, 13 As shown in Fig. 6, existing models appear to have a reasonable representation of hydrogen pairing distribution as well as pore size distribution5 but have difficulty in capturing some key characteristics of the measured FTIR spectra, especially related to C–O, C–S, C–N, C–C vibrational modes. The work on other natural carbon materials indicates a critical role of surface functional groups in chemical sorption and desorption. For example, metal sorption on activated carbon has been found to highly depend on surface functional groups14. Surface functional groups can also modify the surface wettability of carbon materials15.
Figure 6 Generalized phonon densities of states (GDOS) of the EFK, MEK, MarK and PY02 samples from inelastic neutron scattering experiments (ref. 5). Error bars are computed from the square root of the neutron count.
It is anticipated that functional groups on kerogen pore surface may directly control hydrocarbon disposition and release in unconventional reservoirs. A similar mechanism may also play a critical role in regulating metal release from shale matrix in unconventional oil/gas production, an issue of a great environmental concern16. Therefore, development of new structural representation models that can better capture both nanoscale pore structures and surface functional groups of actual kerogen is highly desirable. Our work shows that a combination of DFPT calculations with spectroscopic measurements may provide a useful tool for future model development.
We acknowledge that the comparison made above between model predictions and experimental data is relatively rough. Kerogen is known for its structural heterogeneities1. A possible improvement to the existing comparison could perform DFPT simulations for different possible model structures and then use an averaged FTIR spectrum for comparison. It is expected that a large set of different kerogen model structures would be required to obtain a realistic and statistically-meaningful ensemble of kerogen properties to represent samples characterized experimentally. In this sense, it would be helpful to develop a limited number of end-member model structures such that an actual kerogen structure can represented by a combination of those end-member model structures.
The existing classification of kerogen using the van Krevelen diagram focuses mainly on the origin of the carbon materials and their compositional evolution during thermal cracking, with less attention on an associated structural change (in particular, a change in surface functional group distribution). Given the structural heterogeneities, it would be very difficult, if not impossible, to represent actual kerogen structures with one or a definite set of “representative” model structures. The proposed end-member structure concept will help circumvent this difficulty, in which an actual kerogen structure will be described with a probabilistic distribution of the end members.
This approach may also lead us to develop a comprehensive library of infrared fingerprints for tracing kerogen structural evolution.
First-principles total energy calculations were carried out using DFT, as implemented in the Vienna ab initio simulation package17 (VASP). The exchange-correlation energy was calculated using GGA with the PBE parameterization18 (PBE). Such standard functionals were found in previous studies to correctly describe the geometric parameters and vibrational properties of a variety of carbonaceous structures6, 19,20,21. The interaction between valence electrons and ionic cores was described by the projector augmented wave (PAW) method22, 23. The C(2 s2,2p2), N(2 s2,2p3), O(2 s2,2p4) and S(3 s2,3p4) electrons were treated explicitly as valence electrons in the Kohn-Sham (KS) equations and the remaining core electrons together with the nuclei were represented by PAW pseudopotentials. The KS equation was solved using the blocked Davidson iterative matrix diagonalization scheme24. The plane-wave cutoff energy for the electronic wavefunctions was set to 500 eV, ensuring the total energy of the system to be converged to within 1 meV/atom. Electronic relaxation was performed with the residual minimization method direct inversion in the iterative subspace (RMM-DIIS), preconditioned with residuum-minimization. A periodic unit cell approach was used in the calculations. Kerogen models based on the structures reported recently by Ungerer et al.4 and Bousige et al.5 with different depositional origins and maturity levels were used in the calculations. Integrations in the Brillouin zone were carried out at the Γ point. Using these structures, DFPT linear response calculations were carried out at the GGA/PBE level of theory with VASP to determine the vibrational frequencies and associated intensities. The latter were computed based on the Born effective charges tensor, which corresponds to the change in atoms polarizabilities with respect to an external electric field. Compared to simple finite-displacement approaches typically used in DFT methods to describe vibrational properties of many-particle systems, DFPT-based methods can be regarded as considerably more effective since additional physical properties can be derived from the total energy with respect to perturbations25. We extensively tested the accuracy of DFPT approaches, with and without long-range corrections for van der Waals interaction, in previous studies8, 26,27,28,29,30,31.
Benchmark test calculations were carried out using VASP DFPT for simple molecules in the gas phase. For the linear CO2 molecule (D4h group), the IR-active C–O bending (Eu) and asymmetric stretch modes (A2u) were predicted to be at 631 and 2362 cm−1, respectively, i.e., in fair agreement with the measured values of 670 and 2350 cm−1 32.
The computed IR-active C–H modes for the CH4 molecule (Td group) were 1353 and 2851 cm−1 for the bending (T2) and asymmetric stretch modes (T2), respectively. These results are comparable to the experimental values of 1306 and 3019 cm−1 for gaseous CH432. Some of the discrepancies between DFPT and experimental FTIR results might stem in part from the fact that DFPT/DFT methods do not accurately predict the long-range van der Waals interactions.
Functional groups were identified and analyzed using FTIR spectroscopy on kerogen extracted from shale samples of the Mancos (Utha, USA), Woodford (Oklahoma, USA) and Marcellus (Pennsylvania, USA) formations obtained from TerraTek Inc (Schlumberger). The kerogens were lightly crushed using a mortar and pestle, then analyzed. Powdered samples were placed directly on the attenuated total reflection (ATR) attachment-Smart Orbit Diamond (3000–200 cm−1 band) for analysis. Spectra were collected using a Thermo Nicolet 380 FTIR spectrometer and OMNIC software suite. Scans were done on each sample at a resolution of 4 cm−1, with absorbance spectra ranging from 4000 to 400 cm−1. Blank readings were used between every samples.
Sandia National Laboratories is a multi-program laboratory managed and operated by Sandia Corporation, a wholly owned subsidiary of Lockheed Martin Corporation, for the U.S. Department of Energy’s National Nuclear Security Administration under contract DE-AC04–94AL85000. This work was supported by the National Energy Technology Laboratory (NETL) and by Laboratory Directed Research and Development (LDRD) funding (both to Y. Wang) from Sandia National Laboratories. This work was supported by the X-Shale project enabled through MIT’s Energy Initiative in collaboration with Shell and Schlumberger. Additional support was provided by the ICoME2 Labex (ANR-11-LABX-0053) and the A*MIDEX projects (ANR-11-IDEX-0001-02) co-funded by the French programme Investissements d’Avenir managed by ANR, the French National Research Agency. We thank Tuan Anh Ho and Louise J. Criscenti (Sandia National Laboratories) for useful discussions, and Benoit Coasne (CNRS and University Grenoble Alpes) and Colin Bousige (MIT-CNRS) for providing structural information regarding the EFK, MEK, MarK and PY02 models.
P.F.W. and E.K. carried out the density functional perturbation theory calculations and analyzed the infrared spectra. Y.W., J.N.K., M.M.M. and E.N.M. collected and analyzed the experimental infrared data for the Mancos, Woodford and Marcellus kerogen samples. R.J.-M.P. provided the data and structures of the E.F.K., M.E.K., MarK and PY02 models used as input for density functional perturbation theory calculations. P.F.W. wrote the first draft of the manuscript, with the help of E.K., Y.W. and E.N.M. All the authors discussed the results and commented on the manuscript.
Tissot, B. P. & Welte, D. H. Petroleum Formation and Occurrence, Ch. 4, 131–159 (Springer-Verlag, 1984).
Kerr, R. A. Natural gas from shale bursts onto the scene. Science 328, 1624–1626 (2010).
Vandenbroucke, M. & Largeau, C. Kerogen origin, evolution and structure. Org. Geochem. 38, 719–833 (2007).
Ungerer, P., Collell, J. & Yiannourakou, M. Molecular modelling of the volumetric and thermodynamic properties of kerogen: influence of organic type and maturity. Energy Fuels 29, 91–105 (2015).
Bousige, C. et al . Realistic molecular model of kerogen’s nanostructure. Nature Mater. 15, 576–582 (2016).
Weck, P. F., Kim, E. & Wang, Y. F. Van der Waals forces and confinement in carbon nanopores: interaction between CH4, COOH, NH3, OH, SH and single-walled carbon nanotubes. Chem. Phys. Lett. 652, 22–26 (2016).
Ho, T., Criscenti, L. J. & Wang, Y. F. Nanostructural control of methane release in kerogen and its implications to wellbore production decline. Sci. Rep. 6, 28053 (2016).
Johnson, T. J. et al . Time-resolved infrared reflectance studies of the dehydration-induced transformation of uranyl nitrate hexahydrate to the trihydrate form. J. Phys. Chem. A 119, 9996–10006 (2015).
Mattson, J. S., Mark, H. B. Jr., Kolpack, R. L. & Schutt, C. E. A rapid nondestruction technique for infrared identification of crude oils by internal reflection spectroscopy. Anal. Chem. 42, 234–238 (1970).
Mattson, J. S. Fingerprinting of oil by infrared spectroscopy. Anal. Chem.43, 1872–1873 (1971).
van Krevelen, D. W. Coal: Typology, Chemistry, Physics, Constitution , Ch. 25, 238–262 (Elsevier, 1961).
Seewald, J. S. Organic–inorganic interactions in petroleum-producing sedimentary basins. Nature 426 , 327–333 (2003).
Wang et al. Control of pertechnetate sorption on activated carbon by surface functional groups. J. Colloid. Interface Chem. 305, 209–217 (2007)
Deheryan et al. Direct correlation between the measured electrochemical capacitance, wettability and surface functional groups of CarbonNanosheets. Electrochimica Acta 132, 574–582 (2014).
Engle, M. A. & Rowan, E. L. Geochemical evolution of produced waters from hydraulic fracturing of the Marcellus Shale, northern Appalachian Basin—A multivariate compositional data analysis approach: Inter. J. Coal Geol. 126, 45–56 (2014).
Kresse, G. & Furthmuller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B 54, 11169 (1996).
Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 77, 3865 (1996).
Weck, P. F., Kim, E., Balakrishnan, N., Cheng, H. & Yakobson, B. I. Designing carbon nanoframeworks tailored for hydrogen storage. Chem. Phys. Lett. 439, 354–359 (2007).
Miller, G. et al. Hydrogenation of single-wall carbon nanotubes using polyamine reagents: combined experimental and theoretical study. J. Am. Chem. Soc. 130, 2296–2303 (2008).
Chang, K., Kim, E., Weck, P. F. & Tomanek, D. Nanoconfinement effects on the reversibility of hydrogen storage in ammonia borane: a first-principles study. J. Chem. Phys. 134, 214501 (2011).
Blochl, P. E. Projector augmented-wave method. Phys. Rev. B 50, 17953–17979 (1994).
Kresse, G. & Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B 59, 1758–1775 (1999).
Davidson, E. R. Methods in Computational Molecular Physics Vol. 113 (eds Diercksen, G. et al.) Ch. 4, 95–113 (NATO Advanced Study Institute, Plenum, 1983).
Gonze, X. & Lee, C. Dynamical matrices, Born effective charges, dielectric permittivity tensors, and interatomic force constants from density-functional perturbation theory. Phys. Rev. B 55, 10355 (1997).
Weck, P. F. & Kim, E. Assessing Hubbard-corrected AM05 + U and PBEsol + U density functionals for strongly correlated oxides CeO2 and Ce2O3. Phys. Chem. Chem. Phys. 18, 26816 (2016).
Weck, P. F. & Kim, E. Uncloaking the thermodynamics of the studtite to metastudtite shear-induced transformation. J. Phys. Chem. C 120, 16553 (2016).
Weck, P. F., Kim, E., Tikare, V. & Mitchell, J. A. Mechanical properties of zirconium alloys and zirconium hydrides predicted from density functional perturbation theory. Dalton Trans. 44, 18769 (2015).
Johnson, T. J. et al. Dehydration of Uranyl Nitrate Hexahydrate to the Trihydrate under Ambient Conditions as Observed via Dynamic Infrared Reflectance Spectroscopy. Proc. SPIE 9455, 945504 (2015).
Weck, P. F. & Kim, E. Layered uranium(VI) hydroxides: structural and thermodynamic properties of dehydrated schoepite α-UO2(OH)2. Dalton Trans. 43, 17191 (2014).
Weck, P. F. & Kim, E. Solar Energy Storage in Phase Change Materials: First-Principles Thermodynamic Modeling of Magnesium Chloride Hydrates. J. Phys. Chem. C 118, 4618 (2014).13.
Wang, Y. F. Nanogeochemistry: Nanostructures, emergent properties and their control on geochemical reactions and mass transfers. Chem Geol 378, 1–23 (2014).
Stein, S. E. Carbon dioxide and methane, NIST Standard Reference Database 69: NIST Chemistry WebBook, http://webbook.nist.gov, (Date of access: 22/05/2017) (2009).
Publisher's note: AllAboutShale.com remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
The authors declare that they have no competing interests.
Correspondence to Philippe F. Weck. requests for materials should be addressed to P.F.W. (email:firstname.lastname@example.org)
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.