Thermodynamic and Structural Characterization of Bulk Organic

Article pubs.acs.org/EF

Thermodynamic and Structural Characterization of Bulk Organic Matter in Chinese Silurian Shale: Experimental and Molecular Modeling Studies Liang Huang,*,†,‡ Zhengfu Ning,†,‡ Qing Wang,†,‡ Rongrong Qi,†,‡ Jing Li,†,‡ Yan Zeng,†,‡ Hongtao Ye,†,‡ and Huibo Qin§ †

State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum (Beijing), Beijing 102249, P. R. China Department of Petroleum Engineering, China University of Petroleum (Beijing), Beijing 102249, P. R. China § State Key Laboratory of Heavy Oil Processing, China University of Petroleum (Beijing), Beijing 102249, P. R. China ‡

S Supporting Information *

ABSTRACT: Silurian marine shale in Sichuan Basin is the most significant target zone for shale gas resources in China. In this work, a combined experimental and molecular simulation study was conducted to characterize the thermodynamic and structural properties of the organic matter in Silurian shale. Organic geochemistry experiments and Fourier transform infrared (FTIR) spectroscopy were performed to provide structural parameters for the main skeleton of the organic matter. A realistic molecular model of the organic matter under typical reservoir conditions was generated by molecular dynamics simulations based on the experimental results and documented analytic data. The thermodynamic and structural properties of the organic matter model were discussed in detail. Clear correlations are found among geochemistry properties and structural parameters of organic matter, independent of organic matter type and maturity. Aromatic units in the organic matter are highly condensed, and the interunit linkages are mainly short methylene groups. Ether groups are the dominant oxygenated compounds, while aromatic sulfur is the main form of organic sulfur. Reasonable consistencies are found on results of compositions of the organic matter fractions and physical density between simulated and available experimental data. The isothermal compressibility and thermal expansion coefficient correspond to the general range of a liquid. In addition to micropores, the organic matter contains a large amount of ultramicropores, which contribute a lot to the high porosity and specific surface area. The porous network is highly connected with few dead pores. Interestingly, the introduction of bitumen fractions has little effect on the spacing of polyaromatic units, but it aggravates the relative slippage of polyaromatic units. Also, separation of lighter compounds is observed in the structure. The carbon dioxide molecules are closer to the oxygenated groups, while the nitrogen molecules and methane molecules are closer to the sulfur functional groups and nitrogen functional groups. This proposed organic matter model can serve as a starting point for further theoretical investigations on gas adsorption and transport mechanisms, representative of the organic matter in Silurian shale at molecular scale.

1. INTRODUCTION As an important alternative for conventional natural gas resource, shale gas has received increasing attention since the beginning of its industrial exploitation.1,2 With the technologic progress of multistage fracturing in horizontal wells, booming exploitation of shale gas has brought a revolution in the global energy market.1,3 Shale gas in reservoirs is composed of free gas in fractures and pores, adsorbed gas on various pore surfaces, and absorbed gas in liquids.4 The adsorbed gas accounts for a large proportion of the total share, which, for example, can be up to 60−85% in the Lewis shale.4 Previous studies have shown that most of the adsorbed gas is stored in the nodules of organic matter (OM),1 and the gas storage capacity is in direct proportion to the total organic carbon (TOC).5−9 Although these organic nodules take up only a small share of the total volumes in shales,10 they possess strong affinity for the hydrocarbons due to high specific surface areas.3,5,11 Wang et al. reported that almost half of the total hydrocarbons are adsorbed in the OM because of its large numbers of adsorption sites for gas molecules.12 © XXXX American Chemical Society

The OM in shale is one kind of mixture of multicompounds, which mainly consists of kerogen, resin, asphaltene, hydrocarbons, carbon dioxide, nitrogen, and water.13−15 The involved compounds and compositions are largely dependent on the type and maturity of the OM.16 The characterization of OM is not only essential to better understand gas storage and transport mechanisms,3,17 but is also fundamental to investigate the feasibility of carbon dioxide storage and methane replacement application.1 To date, various experimental methods, such as chemical degradation,18,19 thermal analysis,20,21 infrared and Raman spectroscopy,22,23 solid-state nuclear magnetic resonance (NMR) spectroscopy,24,25 and combined pyrolysis−gas chromatography−mass spectroscopy (Py−GC−MS),26 have been utilized to characterize the structural properties of OM fractions. These experimental techniques serve as fundamental approaches to obtain the compositional and structural Received: January 12, 2017 Revised: March 29, 2017

A

DOI: 10.1021/acs.energyfuels.7b00132 Energy Fuels XXXX, XXX, XXX−XXX

Article

Energy & Fuels

area. The OM of the investigated shale belongs to type I−II1, which is mainly derived from lower hydrobiont and algae.38,39 Detailed information on the geological setting, well identification, formation thickness, and shale geochemistry of this formation can be found in previous investigations.38,39,43−46 2.2. Organic Geochemistry Experiments. The TOC content of shale samples was measured by a LECO CS230 carbon/sulfur analyzer. Samples were crushed and ground to powders finer than 100 mesh, then treated with HCl to dissolve carbonates, and pyrolyzed up to 540 °C. Due to the absence of vitrinite in Silurian shale, vitrinite-like macerals were targeted for thermal maturity evaluation. The reflectance measurements of vitrinite-like materials were performed on MPV-SP microscope equipped with oil-immersion objective lens and photometer. The measured reflectance was then converted to equivalent vitrinite reflectance by a conversion formula.47 Pyrolysis analyses of samples were carried out with an OGE Workstation. Samples were first crushed into particles smaller than 5 mm, then 100 mg samples were pyrolyzed from 90 to 600 °C to obtain the mass of soluble hydrocarbons (S1), pyrolyzed hydrocarbons (S2), and the maximum temperature for volatilization of hydrocarbons (Tmax). Also, Tmax was utilized to estimate the vitrinite reflectance according to the method proposed by Jarvie et al.48 The classical RockEval parameters, including hydrogen index (HI = S2 × 100/TOC), production index (PI = S1/(S1 + S2)) and generation potential (GP = S1 + S2), were further determined.49 2.3. FTIR Spectroscopy. Chemical structures of studied OM were evaluated by a MAGNA-IR56 FTIR spectrometer. Shale samples were crushed and ground to fine powders smaller than 200 mesh, then a succession of HCL, HNO3, and HF treatments were adopted to remove the carbonate minerals, pyrite minerals, and silicates of grounded samples according to the process provided by Ballice et al.50 Powder samples of 0.5 mg after acid treatments were mixed with 200 mg dry KBr powders, then pressed under vacuum. FTIR spectroscopy was performed at ambient temperature, and spectra was recorded between 4000 and 400 cm−1 with 50 scans at a resolution of 4 cm−1. Lorentz peak-fitting of spectra curves and determination of integrated peak areas were conducted with custom Origin software for the functional group quantifications. 2.4. Generation of Molecular Models. In this section, development of molecular models of individual fraction that composes the OM is discussed, including kerogen, bitumen fractions, and lighter compounds. These fractions are integrated to build the bulk OM model of Longmaxi shale. 2.4.1. Kerogen. Kerogen is composed of polyaromatic units crosslinked by alkyl, ether, and sulfur bridge.28,51 It constitutes the main skeleton of OM and can represent most of the structural characteristics of OM. Detailed analytical data of elemental compositions and functional groups for overmature type II kerogen have been reported by Kelemen et al.42 using XPS method and solid-state C13 NMR spectroscopy. Although the kerogen in Longmaxi shale belongs to overmature type I−II1, not exactly the same as overmature type II, kerogens of different types converge to a similar composition in the overmature stage. Thus, the experimental results of Kelemen et al.42 can still provide general features of the kerogen in Longmaxi shale. These analytic data, combined with the derived structural parameters of the OM in Longmaxi shale based on experimental results, are utilized to derive the kerogen molecule structure. The building process of the kerogen molecule follows the process proposed by Ungerer et al.:52

parameters of OM fractions, which are essential for the generation of molecular models of OM fractions.27−29 Meanwhile, the development of molecular simulation techniques provides an effective tool to predict the thermodynamic properties and structural characteristics of OM fractions at microscopic level, in complement to the macroscopic experiments.25,28,29 Encouraging results for studies on the structures and gases adsorption characteristics of OM fractions with molecular modeling techniques have been reported.30−35 However, unfortunately, few attempts have been made to investigate the OM as an integrated bulk macromolecule under reservoir conditions with the combination of experimental methods and molecular simulation techniques. Silurian shale in Sichuan Basin is currently the most important target zone for shale gas exploration and exploitation in China.36−38 Extensive research has been performed on the characterization of Silurian shale,36−39 while few studies have been concentrated on the molecular characterization of OM fractions.1 In addition, the OM in Silurian shale has a few unique characteristics. The thermal maturity of this OM is very high with the equivalent vitrinite reflectance ranging between 2.0 to 4.5%,38,40 and the sulfur content in kerogen (average 7.69 wt.%)41 is also much higher than documented kerogen such as Duvernay and Green River kerogen.42 As a result, documented molecular models of OM in other shales cannot be directly applied to the OM in Silurian shale. Collell et al.28 proposed an average molecular model of type II OM in the middle of oil formation window using molecular simulation techniques. Although this model is inappropriate for our Silurian sample in the dry gas window, the assumption of atomic balance adopted by Collell et al. during the maturation process is highly provoking to determine the compositions of lighter compounds. However, it should be noted that further adjustments of the compositions of gases compounds are needed since these lighter compounds are inclined to produce from the main skeleton of Silurian OM in the dry gas window. In this work, we aim to gain insights into the structural and thermodynamic properties of the OM in shale of Lower Silurian Longmaxi Formation in Sichuan Basin of China by a combined experimental and molecular simulation method. The experimental techniques enable us to analyze the compositional and structural features of the OM and derive important structural parameters for the OM macromolecule. Then, molecular simulation techniques are utilized to generate the molecular model of bulk OM under typical reservoir conditions. For OM fractions, we mainly consider kerogen, bitumen compounds, and residual lighter compounds. The assumption of atomic balance, associated with the experimental compositions of bitumen fractions and produced gas fractions, is utilized to determine the compositions of OM fractions in our OM model. This molecular model can serve as a starting point for further theoretical investigations of OM in Silurian shale. In the future work, the presented molecular model could be continually improved to better reproduce the physicochemical properties of the OM in Silurian shale based on more detailed analytical data.

(1) Data of functional groups and elemental compositions are gathered based on XPS method and solid-state C13 NMR spectroscopy. (2) A total number of carbon atoms (150−250) is selected, and the quantity of aromatic carbon is determined according to the obtained aromaticity from experiments. (3) The quantities of H, O, N, and S atoms are determined by matching the functional data and elemental compositions.

2. MATERIALS AND METHODOLOGY 2.1. Materials. Ten marine shale samples were collected from the Lower Silurian Longmaxi Formation in Sichuan Basin of China. Six of the studied samples were obtained from the evaluation well N201 in the Changning shale gas demonstration area, while the other four samples are outcrops sampled from the same formation in Chongqing B


Article

Energy & Fuels

units in the molecular structures of asphaltene and resin according to the conclusion of Yen.61 For aromatic hydrocarbons and saturated hydrocarbons, we limit our consideration to representative compounds. We choose toluene and dimethylnaphtalene as the representative aromatic compounds, which are also adopted in Collell’s work.28 For saturated hydrocarbons, naphthenic molecules are considered in this work, since nearly all long aliphatic compounds disappear in the metagenesis stage, and lighter alkyl compounds will be expelled from the bitumen during the extraction process. In this study, we choose fichtelite, a common tricyclic alkane in bitumen, as the representative saturated compound. The specific molecular size of asphaltene and resin, and the compositions of representative compounds are adjusted to match the composition distribution of bitumen compounds in the OM of Longmaxi shales.60 2.4.3. Lighter Compounds. In this study, we consider methane, ethane, carbon dioxide, nitrogen, and water as lighter compounds. Lighter compounds in OM can be sorted into two types: most of lighter compounds are prone to be expelled out of the carbon skeleton with the production of shale gas, while a small number of lighter molecules are blocked in the tiny pockets of OM structure. The former type can be quantitatively analyzed by gas chromatography of produced gas, while the latter one is hard to be directly quantitated because of the loss during the sampling process. Previous scholars have assumed that the atomic balance is preserved during the maturation evolution, namely the produced compounds during the maturation process are not expelled.28,62 We also adopt this assumption to determine the initial content of lighter compounds in the OM model. Detailed analytical data of the immature type I kerogen in Green River shale have been documented by Kelemen et al.42 It is selected for atomic balance analysis, since kerogen is almost the exclusive compound of immature OM. In this study, S/C atomic ratio is excluded for atomic balance due to the unique high sulfur content in the kerogen of Longmaxi shale.41 Furthermore, the produced gas compounds are deducted from the initial content of lighter compounds to represent the realistic distribution of lighter compounds in the reservoir. OM. The bulk OM model of Longmaxi shale is constructed by integrating the kerogen molecules, bitumen fractions and lighter compounds discussed above. First, six kerogen molecules meeting the minimized size requirements are selected to form the skeleton of the OM model. Then, a suite of bitumen fractions, including asphaltene, resin, toluene, dimethylnaphtalene and naphthenic, are introduced by matching their compositions with experimental results60 summarized in section 3.2.2. Furthermore, lighter compounds, including methane, ethane, carbon dioxide, nitrogen, and water, are introduced to generate the initial version of the OM model according to the atomic balance during the maturation evolution. However, it is worth noting that the integrated OM model in this stage cannot represent the actual occurrence structure of the OM under reservoir conditions, since a large number of lighter compounds will be expelled during the dry gas window. To build a more realistic OM model, the produced lighter compounds, with their compositions summarized in section 3.2.3, are deducted from the present model, which is constrained by the gas chromatography results of produced gas from Longmaxi shale.63 Finally, the realistic bulk OM model, with the breakdown of compositions of OM fracions listed in section 3.2.3, is constructed to characterize the OM of Silurian shale under typical reservoir conditions. 2.5. Simulation Details. 2.5.1. Geometry Optimization of OM fractions. The optimal geometry configurations of OM fractions are obtained by molecular mechanics and molecular dynamics simulations using the Forcite module in Materials Studio. The methodology of geometry optimization refers to the work of Ru et al.26 For the molecular mechanics calculation, we adopt the smart minimization algorithm with a fine convergence criterion. The max iterations are set as 5000. The pcff force field is used to calculate the interatomic interactions, which describes the dispersion-repulsion interactions using Lennard-Jones 6−9, and describes the Coulombic interactions using force field assigned point charges.64 The energy of

(4) The number of aromatic units is determined to approximate the experimental result of aromatic carbons per aromatic unit. (5) The initial version of the kerogen structure is generated by connecting the determined units and adjusting the bond lengths, bond angles, and dihedral angles. (6) The generated kerogen structure is refined to match the experimental data of functional groups and elemental compositions. In this study, only organic sulfur groups including thiol, sulfide, and thiophenic are considered for the sulfur-bearing groups. The mineral sulfur forms such as pyrite and sulfate are excluded for discussion. The ratio of sulfur to carbon number is determined by the sulfur content of the kerogen in Longmaxi shale by Fu et al.41 The fractions of individual organic sulfur group result from the experiment-based correlations discussed in section 3.1.1. For oxygenated groups, carbonyl groups, alcoholic groups and ether groups are considered, while oxygenated compounds in residual minerals and water are excluded from discussion. Fractions of oxygenated units are determined by the analysis of FTIR spectra. In the case of nitrogenbearing groups, we mainly focus on pyrrolic, pyridinic and quaternary units, with the corresponding fractions derived from the C13 NMR results.42 The molecular topology of the kerogen is largely inspired from the kerogen fragments of previous models in the literature. 52,53 Aromaticity and average carbon number per aromatic cluster are derived from the experiment-based correlations in section 3.1.1. To better approach a realistic molecular model, branched alkyl and ether segments, in addition to linear alkyl fragments, are also introduced as cross-links between the aromatic units. Moreover, apart from aromatic rings and unsaturated units, we also design monocyclic and polycyclic saturated rings to match the low H/C ratio of the kerogen. For the initial version of kerogen molecule, we have tried to build a flat-shape structure with minimized extension so as to approach the equilibrium conformation. Then adjustments on the functional groups location and polyaromatic units shape are performed to improve the agreements between the model and experimental results. Each change of the molecule structure is checked and optimized by performing a local energy minimization. Thereafter, the molecule model is relaxed through a succession of molecular dynamic simulations to obtain the most stable conformation. Finally, pcff force field is utilized to assign the atom types and charges to verify the structural consistency. For the structure size, a unit of 100−260 carbon atoms for each kerogen molecule, a total of 1700 atoms and a box size of 25 Å of bulk kerogen model are considered as the minimized requirements to capture the major features of kerogen, as proved by previous research.28,52,54−58 2.4.2. Bitumen. Bitumen refers to the fraction of sedimentary OM soluble in dichloromethane.59 Bitumen can be divided into asphaltene, resin, aromatic hydrocarbons, and saturated hydrocarbons.60 The asphaltene and resin molecules are produced from the kerogen skeleton during the early stage of the maturation evolution. Their structures are similar to that of kerogen, but with a lower molecular weight. Both asphaltene and resin are composed of polyaromatic units with polar heteroatoms (O, N, and S), which are similar to those in crude oils of the same type and maturity.28 Compared with asphaltene, resin has a lower aromaticity with more polar fractions, and it also has a lower molecular weight.57 However, the content of resin is much higher than that of asphaltene in bitumen groups.60 To meet the higher weight faction of resin compared with that of asphaltene, either the molecular number or the molecular weight of resin can be increased. Unfortunately, the increase of the molecular number of resin can lead to an extremely large simulation system due to the small factions of resin and asphaltene in the OM. Thus, for the purpose of atomic balance, the molecular weight of resin is adjusted to be higher than that of asphaltene in this study. Also, more polar fractions (nitrogen and oxygen groups) are assigned for resin to represent its higher molecular polar, while sulfur groups are excluded for discussion due to their low concentrations. The nitrogen groups are in the form of quinoline units, while the oxygen groups are represented by ether C


Article

Energy & Fuels Table 1. Organic Geochemical Characteristics of Studied Silurian Shale Samplesa sample

depth (m)

TOC (wt %)

Reqv (%)

S1 (mg/g)

S2 (mg/g)

Tmax (°C)

Rcal (%)

HI (mg/g)

PI

GP (mg/g)

N201−1 N201−2 N201−3 N201−4 N201−5 N201−6 CQ-1 CQ-2 CQ-3 CQ-4

2152.26 2261.20 2297.00 2319.85 2336.46 2377.84 outcrop outcrop outcrop outcrop

2.22 2.65 3.19 2.82 3.01 2.83 2.53 2.50 2.84 1.39

2.01 2.09 2.13 2.27 2.31 2.37 3.51 3.55 1.68 3.53

0.03 0.02 0.01 0.02 0.02 0.01 0.04 0.07 0.01 0.05

0.05 0.06 0.03 0.04 0.04 0.03 0.03 0.04 0.02 0.03

511 517 512 522 521 531 595 598 490 598

2.04 2.15 2.06 2.25 2.22 2.41 3.55 3.60 1.66 3.60

2.25 2.26 0.94 1.42 1.33 1.06 1.19 1.60 0.70 2.16

0.38 0.25 0.25 0.33 0.33 0.25 0.57 0.64 0.33 0.63

0.08 0.08 0.04 0.06 0.06 0.04 0.07 0.11 0.03 0.08

a The equivalent vitrinite reflectance (Reqv): converted from reflectance of vitrinite-like material. The estimated vitrinite reflectance (Rcal): calculated from Tmax. Hydrogen Index (HI): S2 × 100/TOC. Production Index (PI): S1/(S1 + S2). Generation potential (GP): S1 + S2.

Coulomb and van der Waals is calculated based on atoms, and a fine cutoff distance of 15.5 Å is selected. For the molecular dynamics simulation, we adopt the annealing strategy to obtain the lowest energy structures of all annealing cycles.26 Ten annealing cycles with 5 heating ramps per cycle and 2000 steps per ramp are performed. The simulation is carried out using the canonical ensemble (NVT) with the temperature ranged from 300 to 800 K and controlled by the Nose thermostat.26,28 The interatomic interactions are calculated by the pcff force field, while the energy of Coulomb and van der Waals is calculated based on atoms with a fine cutoff distance of 15.5 Å. The total simulation time of 200 ps is determined to ensure the energy convergence for the lowest energy structures of different annealing cycles. After the annealing treatment, the output configuration of each fraction is further optimized by the molecular mechanics simulation with the above setting parameters so as to ensure the structure reaches the lowest energy state. 2.5.2. Structure Optimization of the OM Model. The generation of the initial structure of the OM model is executed by the Amorphous Cell module. The OM fractions after geometry optimization are initially enclosed in a big simulation box with periodic boundary conditions, with the target density being 0.1 g/cm3.28,52 Then the system is optimized by molecular mechanics calculations to obtain the optimal configurations.26 The temperature is 300 K, and the pcff force field is used.64 Five different configurations are generated in this stage in order to investigate the uncertainties caused by different initial configurations. The initial configurations of the OM models are then relaxed by a succession of molecular dynamics simulations.28 During the relaxation process, the pcff force field is used to calculate the interatomic interactions.28,52 The energy summation method for Coulomb and van der Waals is atom based,26 and a fine cutoff distance of 15.5 Å is selected. The pressure is imposed using an Andersen barostat, while the temperature is imposed with a Nose thermostat.26,28 The relaxation process refers to the procedure adopted by Collel et al.28 and Ungerer et al.52 The OM model is first relaxed with a NVT ensemble at 1000 K. Then succeeding NPT simulations are performed with a stepwise decreasing temperature from 1000 to 300 K. The simulation time is 200 ps for the NVT simulation and 300 ps for the NPT simulations. The simulation time of each run is set to be long enough to ensure the density convergence. 2.5.3. Physical Properties Simulation. The isothermal compressibility and the thermal expansion coefficient of the OM model are determined by the fluctuation method. The relaxed OM model is further simulated by molecular dynamics with a NPT ensemble from 800 to 300 K. In this simulation, we adopt the pcff force field, which has been proved to provide good predictions of thermodynamic properties of organic compounds.28,52,65 The van der Waals interactions are atom-based with a fine cutoff distance of 15.5 Å, while the Coulombic interactions are described by the Ewald method. The pressure is imposed using an Andersen barostat, while the temperature is imposed with a Nose thermostat.26,28 A total of 10 000 frames are obtained from the NPT ensemble over 3 ns

simulation for each run to determine the isothermal compressibility and the thermal expansion coefficient. The isothermal compressibility determined by the volume fluctuations can be given by66 β=

1 (⟨V 2⟩ − ⟨V ⟩2 ) ⟨V ⟩kBT

(1)

where V is the volume of the structure, T is the temperature, kB is the Boltzmann constant, and ⟨···⟩ denotes the ensemble average. Similarly, the thermal expansion coefficient can be calculated by66

α=

1 (⟨VE⟩ − ⟨V ⟩⟨E⟩ + P(⟨V 2⟩ − ⟨V ⟩2 )) ⟨V ⟩T 2kB

(2)

where P is the system pressure and E is the total energy of the system. The characterization of the porous structure and the specific surface area is executed by the atom volumes and surface tool in Materials Studio. The pore volumes are detected by probe molecules of given radius. The probe molecules roll over the van der Waals surface of the atoms to determine the surface of solid skeleton.67 The region left by the surface of solid skeleton is identified as pore volume.

3. RESULTS AND DISCUSSION 3.1. Experimental Characterization of OM. 3.1.1. Organic Geochemistry Characteristics. The organic geochemistry results of the samples are shown in Table 1. The TOC content of all the samples except CQ-4 is above 2.0 wt %., corresponding to the range of good organic richness.68 The generation potential (GP), ranging from 0.03 to 0.11 mg/g, indicates a poor hydrocarbon generation potential. All of the shale samples have low S1 and relatively high TOC values indicating the occurrence of indigenous hydrocarbons.69 Equivalent vitrinite reflectance (Reqv) of the samples matches well with calculated vitrinite reflectance (Rcal). The thermal maturity of Longmaxi shale can be mainly classified as overmature level in the dry gas window, which is further confirmed by the values of Tmax (average 540 °C) and PI (average 0.40).70 A typical reservoir condition (360 K, 23.6 MPa) is determined for the Silurian OM according to the reported temperature and pressure gradients.39 To date, many relationships among geochemistry properties and structural parameters of OM have been reported based on experimental work.42,71−73 However, few strong correlations have been observed due to limited experimental data. We extend previous work in the literature by collecting a large number of documented geochemistry properties and structural parameters on samples of different source types and maturities. It is worth noting that these samples cover a wide range of source type from type I to type III, with the maturity ranging D


Article

Energy & Fuels from immature to overmature.24,42,73−76 Based on these collected experimental data, clear regression equations among geochemistry properties and structural parameters are established, as shown in Figure 1. Then these regression equations, combining with the organic geochemistry results of Longmaxi shale in our work, are utilized to derive the important structural parameters of the OM in Longmaxi Formation. The derived structural parameters, including the percent aromatic carbon (fa), the average aliphatic carbon chain length (Cn), the ratio of aromatic sulfur to organic sulfur (fs) and the average number of aromatic carbons per cluster (Ca), are listed in Table 2. As seen in Table 2, the differences of derived structural parameters among all the samples are ignorable, facilitating the construction of an averaged molecular structure. The averaged H/C atomic ratio (H/C) (49.36%) corresponds to the overmature thermal maturity of the OM. The poor hydrogen content of OM can be attributed to three mechanisms: (1) drastic cracking of aliphatic hydrocarbon; (2) dealkylation of aromatics; (3) aromatization of adjacent aliphatic carbon. These structural changes are further substantiated by the high values of fa and low values of Cn. The averaged fa (85.60%) indicates a high condensation degree of aromatic units, which is also proven by the values of Ca. Ca is another parameter that is closely related to the aromatic cluster size. The characterized Ca (average 20.81) is in good agreement with the experimental results for overmature kerogen.42 The values of Cn approach 1.0 in this work, indicating that nearly all aliphatic carbons are in the form of short methylene groups, which is close to the structure of coals.77 The average fs (81.77%) indicates that aromatic sulfur is the dominant form of organic sulfur in the OM, similar to the overmature kerogen by Kelemen et al.42 3.1.2. Chemical Structure Characterization. In order to characterize the chemical functional groups in the studied OM, sample N201−2 (Ro = 2.09%) and CQ-3 (Ro = 1.68%) after demineralization were analyzed by FTIR spectroscopy. The samples with higher thermal maturity are excluded from our discussion, because damage or collapse of functional group bands can result in tiny peak areas and cause large uncertainties in peak area determination, as indicated by Iglesias et al.78 Standard FTIR spectra of these samples are displayed in Figure 2. The band assignments of investigated functional groups are listed in Table 3. As seen in Figure 2, samples N201−2 and CQ-3 show a small absorption peak of aromatic CHx stretching at 3100−3000 cm−1 and an absence of absorption peak of CO groups at 1800−1650 cm−1. Also, the peak intensity of aliphatic CHx stretching at 3000−2800 cm−1 is remarkably weak, while the peak intensity of the aromatic CHx out-of-plane deformation at 900−730 cm−1, as well as the aromatic CC ring stretching at 1650−1550 cm−1, is significantly strong for N201−2 and CQ-3 samples. These results correspond to the evolutionary paths of functional groups with the thermal maturation of OM.16,24,27 The oxygenated groups detach from the OM skeleton, the aliphatic chains shorten and the degree of aromaticity enhances as the thermal maturity increases.16,24,27 For OM of high maturity (N201−2 and CQ-3), the oxygenated groups and aliphatic side chains tend to disappear, and the condensation of aromatic rings becomes the major evolution path, as shown in Figure 2. Apart from qualitative characterization, many semiquantitative indexes have been derived from FTIR spectra for structural characteristics evaluation.82−84 In this section, we focus on these indexes closely related to the chemical structure of OM:

Figure 1. Correlations among geochemistry properties and structural parameters of OM: (a) elemental hydrogen per 100 carbons versus hydrogen index; (b) aromatic carbon content versus elemental hydrogen per 100 carbons; (c) average aliphatic chain length versus aromatic carbon content; (d) aromatic sulfur content on an organic E


Article

Energy & Fuels

cm−1. The derived results of semiquantitative indexes from FTIR spectra are shown in Table 5. 3.1.2.1. Aliphatic Chain Length. The ratio of integrated peak areas of CH2 asymmetric stretching signals (2923 cm−1) and CH3 asymmetric stretching signals (2956 cm−1) is utilized to represent the CH2/CH3 ratio in the OM. This ratio is proposed to estimate the chain length and degree of branching of aliphatic hydrocarbons.82 The CH2/CH3 ratio presents an increase in the metagenesis stage (Table 5), which was also reported by Chen et al.82 The increase of the CH2/CH3 ratio can be attributed to the detachment of methyl branches during the massive gas formation period.73 Thus, for the OM in Longmaxi Formation with high maturity, nearly all the long aliphatic chains and branching chains are removed, and short methylene groups are the main structure of aliphatic carbons, which is also indicated by the value of Cn derived from the experiment-based correlations. Moreover, this conclusion can be further demonstrated by the neglected signal intensity for the vibration of long-chain alkane ((CH2)n, n ≥ 4) in the range from 730 to 700 cm−1 (Figure 2). 3.1.2.2. Aromaticity. Aromaticity in this section is defined as the fraction of aromatic CH peak areas among the sum of aromatic CH and aliphatic CHx peak areas. The value of defined aromaticity is generally less than 1.0 unless the aliphatic CHx vibration signals can be neglected. Aromaticity increases with the increase of thermal maturity.82 For N201−2, the aromaticity is up to 93.45%, indicating that aromatic carbons dominate the organic carbons in the OM of Longmaxi Shale, which agrees with the result of the experiment-based correlation (Table 2). 3.1.2.3. Condensation of Aromatic Rings. Degree of condensation of aromatic rings is represented by the ratio of aromatic CH out-of-plane vibration and aromatic CC stretching. It increases as the thermal maturity increases (Table 5), which can be attributed to the dealkylation of aromatics and the aromatization of aliphatic carbon structure during the maturation evolution.82 The removal of methyl groups in the aromatic rings contributes to the increase of aromatic CH vibration, while the aromatization results in the increase of both aromatic CH vibration and aromatic CC stretching, but generally with a higher contribution to aromatic CH vibration.73 N201−2 has a high value of condensation degree of aromatic rings (0.95), which corresponds to a large aromatic cluster size, as indicated by the Ca values derived from the experiment-based correlation. 3.1.2.4. Oxygenated Group Compounds and Compositions. In this section, we mainly consider the organic oxygenated group compounds, including carbonyl groups, alcoholic groups, and ether groups.16,18,24,27 Oxygenated compounds in residual minerals and water are excluded from discussion because they are subordinate to the whole skeleton of OM and can be removed after successive acid and dehydration treatments.82 As seen in Table 5, the ether groups dominate the oxygenated groups, while the carbonyl groups and alcoholic groups tend to disappear for OM samples of high maturity. For N201−2, nearly all of the carbonyl groups and alcoholic groups are removed, and ether groups are the exclusive oxygenated compounds (Table 5). 3.2. Molecular Models. 3.2.1. Molecular Structure of Kerogen. A representative molecular structure of kerogen is shown in Figure 4, with its compositions and structural parameters listed in Table 6. The formula of the kerogen unit is C194H96O9N4S7, and the molecular weight is 2848 g/mol. To

Figure 1. continued sulfur basis versus aromatic carbon content. The documented experimental data are collected from previous work in the literature.24,42,73−76

Table 2. Structural Parameters Derived from ExperimentBased Correlationsa sample

H/C (%)

fa (%)

Cn

fs (%)

Ca

N201−1 N201−2 N201−3 N201−4 N201−5 N201−6 CQ-1 CQ-2 CQ-3 CQ-4

49.46 49.47 49.28 49.35 49.33 49.29 49.31 49.37 49.24 49.45

85.53 85.53 85.65 85.61 85.61 85.64 85.63 85.59 85.67 85.54

1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00

81.72 81.72 81.80 81.77 81.78 81.80 81.79 81.76 81.82 81.72

20.79 20.79 20.82 20.81 20.81 20.82 20.82 20.81 20.83 20.80

a

Data of these structural parameters are derived from the regression equations established from experimental results of previous work.24,42,73−76

Figure 2. FTIR spectra of samples of different thermal maturity. Functional group regions: (1) hydroxyl group; (2) aromatic CHx; (3) aliphatic CHx; (4) carbonyl groups; (5) aromatic carbon; (6) aliphatic CHx in-plane bending; (7) alcoholic hydroxyl and ether groups; (8) aromatic CHx out-of-plane deformation. Specific peak assignments are listed in Table 3.

Table 3. Band Assignments of the FTIR Spectra (after Painter et al.;79,80 Wang and Griffith;81 Chen et al.82) wavenumber (cm−1)

assignment

3400 3100−3000 3000−2800 1800−1650 1650−1550 1480−1350 1300−1100 1030 900−730

hydrogen-bonded OH aromatic CHx stretching aliphatic CHx stretching carbonyl groups aromatic CC ring stretching aliphatic CHx in-plane bending alcoholic hydroxyl groups ether groups aromatic CHx out-of-plane deformation

(1) aliphatic chain length; (2) aromaticity; (3) degree of aromatic ring condensation; (4) oxygenated group compounds and composition. Definitions of these indexes are listed in Table 4. Integrated peak areas are utilized to calculate these indexes. For these overlapping peaks, curve-fitting of peaks are performed to determine independent peak areas. Figure 3 demonstrates the fitting peaks of different aliphatic CHx stretching signals for N201−2 in the range of 3000−2800 F


Article

Energy & Fuels

Table 4. Definition of Semiquantitative Indexes from the FTIR Spectra (after Iglesias et al.;78 Chen et al.;82 Ganz and Kalkreuth;85 Guo and Bustin;86 Ruan et al.;87 Sun88) index chain length (CL) aromaticity (AR) degree of condensation (DOC) carbonyl fraction (Oca) alcoholic fraction (Oal) ether fraction (Oet)

band region (cm−1)

index description CH2 asymmetric stretching/CH3 asymmetric stretching aromatic CHx/(aromatic CHx + aliphatic CHx) aromatic CHx out-of-plane deformation/CC stretching carbonyl groups/carbon−oxygen groups alcoholic groups/carbon−oxygen groups ether groups/carbon−oxygen groups

2923/2956 [(3100−3000)+(900−730)]/[(3100−3000)+ (900−730)+(3000−2800)+(1480−1350)] (900−730)/(1650−1550) (1800−1650)/[(1800−1650)+(1300−1100)+1030] (1300−1100)/[(1800−1650)+(1300−1100)+1030] 1030/[(1800−1650)+ (1300−1100)+1030]

Figure 3. FTIR curve-fitting for sample N201−2 in the range from 3000 to 2800 cm−1. Peak assignment is according to Painter et al.79

meet the minimized size requirement for a total of 1700 atoms, six kerogen molecules are included to form the main skeleton of the OM model. As seen in Table 6, the compositions and main structural parameters of the proposed kerogen structure match well with the experiment-based results in this work and the documented data in the literature. The atomic compositions of this structure are in fairly good agreement with the target data. The carbon skeleton is composed of eight polyaromatic units containing 17−24 aromatic carbons (average 19.88), which is in accordance with the analytic data. Aromaticity of this structure is slightly lower than the target data, due to adjustments made to meet the balance with aliphatic carbons. Also, these polyaromatic clusters are linked by short linear, branched chains and saturated cycles, and a good match to the ratio of methylene carbons and terminal aromatic carbons (Cn) is obtained. Furthermore, accordant with the FTIR results, all the oxygen atoms are distributed in ether groups. Moderate deviations are observed for nitrogen and sulfur groups because of their low contents in the structure. The four nitrogen atoms are involved in three pyrrolic groups and one pyridinic group, while the seven sulfur atoms are distributed in six aromatic groups and one aliphatic group. 3.2.2. Molecular Structures of Bitumen Fractions. The representative molecular structures of resin and asphaltene are shown in Figure 5. The proposed structures of resin and asphaltene resemble each other, which are also similar to the fractions of the kerogen molecule. The formula of resin model

Figure 4. Molecular structure of the kerogen in Longmaxi shale. Atom representation: black for carbon atoms, white for hydrogen atoms, red for oxygen atoms, yellow for sulfur atoms, and blue for nitrogen atoms.

is C33H29ON with a molecular weight of 455 g/mol, while the formula of asphaltene model is C22H16O, and the molecular weight equals 296 g/mol. It should be noted that the molecular weight of resin model is bigger than that of asphaltene model in this work, different from the experimental result that asphaltene is heavier than resin.57 As discussed in section 2.4, this is because the content of asphaltene is smaller than that of resin in bitumen groups,60 the molecular weight of resin needs to be increased to meet the higher weight faction of resin compared with that of asphaltene. The resin structure contains one ether group and one quinoline group. It is more polar than the asphaltene structure, which involves only one ether group. Also, the aromaticity of asphaltene (0.91) is higher than that of resin (0.70), while the H/C ratio of asphaltene (0.73) is lower than that of resin (0.88). These differences can be well explained by the extra saturated cycles in the structure of resin. The representative compounds of aromatic hydrocarbons and saturated hydrocarbons, associated with the weight fractions of bitumen compounds, are listed in Table 7. The final molecular structure and molecular number of each

Table 5. Semiquantitative Indexes Derived from FTIR Spectra sample

Reqv (%)

CL

AR (%)

DOC

Oca (%)

Oet (%)

Oal (%)

N201−2 CQ-3

2.09 1.68

6.05 4.09

93.45 90.28

0.95 0.58

0.00 0.00

100.00 82.18

0.00 17.82

G


Article

Energy & Fuels Table 6. Compositional and Structural Parameters of the Kerogen Structure

Table 8. Compounds Distribution in the Initial OM Model Based on Atomic Balance

property

structural parameter

target data

this model

representative molecule

atomic compositions

H/C O/C N/C S/C aromaticity ( fa, %) aromatic carbon atoms per cluster methylene carbons/terminal aromatic carbons (Cn) carbonyl (mol % O) ether (mol % O) alcoholic (mol % O) pyrrolic (mol % N) pyridinic (mol % N) quaternary (mol % N) aromatic S (mol % organic S) aliphatic S (mol % organic S)

0.494a 0.047b 0.02b 0.035c 85.60a; 82b 20.81a; 20b 1.0a; 1.1b

0.495 0.046 0.021 0.036 81.96 19.88 1.0

0d 100d 0d 62b 15b 23b 81.77a 18.23a

0 100 0 75 25 0 85.71 14.29

kerogen resin asphaltene toluene 1,3-dimetylnaphtalene tricyclo[8.4.0.02,7] tetradecane methane ethane nitrogen carbon dioxide water

C groups

O groups

N groups

S groups

formula

molecular weight (g/mol)

molecular number

C194H96O9N4S7 C33H29ON C22H16O C7H8 C12H12 C14H24

2848 455 296 92 156 198

6 1 1 1 1 1

CH4 C2H6 N2 CO2 H2O

16 30 28 44 18

500 3 13 12 10

between our initial OM model and an immature type I kerogen of Green River Shale42 are listed in Table 9. As seen in Table 8,

a

Table 9. Atomic Balance between the Initial OM Model and an Immature Kerogen of Green River Shale

Data derived from experiment-based correlations in this work. Reference data for overmature kerogen by Kelemen et al.42 cData derived from sulfur content analysis by Fu et al.41 dData derived from FTIR spectra in this work. b

a

phase

H/C

O/C

N/C

OM model immature kerogena

1.53 1.53

0.051 0.051

0.029 0.029

Analytic data documented by Kelemen et al.42

methane molecules dominate the lighter compounds, consistent with the composition analysis of produced gas from Longmaxi shale.63 Also, the atomic compositions of this initial OM model match fairly well with analytic data of the immature kerogen of Green River Shale,42 as shown in Table 9. However, as discussed in section 2.4, the present OM model derived from atomic balance is not the realistic occurrence state of the OM under reservoir conditions, and these lighter compounds prone to release from the OM should be deducted to obtain the final OM model. Table 10 lists the deducted

Figure 5. Molecular model units of resin and asphaltene in Longmaxi shale. Atom representation refers to the kerogen molecular model: (a) resin molecular model and (b) asphaltene molecular model.

compound are adjusted to match the weight fraction distribution with experimental data in Zhang’s work.60 As seen in Table 7, the weight fraction distribution of bitumen compounds of this work is in agreement with that of experimental result. The molecular numbers of resin and asphaltene with respect to six kerogen molecules correspond to the typical atomic balance of the major compounds in OM.28 3.2.3. Molecular Model of OM. Associated with six kerogen molecules and a suite of bitumen compounds corresponding to the realistic weight fraction distribution, lighter compounds are introduced to build the integrated OM model. The initial molecular numbers of lighter compounds are determined based on the assumption of atomic balance during the maturation evolution. The compound distribution in the OM is shown in Table 8, and the comparison results of atomic compositions

Table 10. Deducted Compounds Distribution from the Initial OM Model representative molecule

molecular number

mole fraction (%)

experimental dataa (%)

methane ethane nitrogen carbon dioxide water

460 2 2 2 2

98.29 0.43 0.43 0.43 0.43

98.54 0.48 0.47 0.5 0.24−0.73

a

Analytic data of gas chromatography on produced gas from Longmaxi shale documented by Gao et al.89

Table 7. Weight Composition of Bitumen Compounds in the OM Model bitumen compounds

representative molecule

molecule numbers

weight fraction (%)

experimental dataa (%)

resin asphaltene aromatic hydrocarbons

resin asphaltene toluene 1,3-dimetylnaphtalene tricyclo[8.4.0.02,7]tetradecane

1 1 1 1 1

37.81 24.60 20.61

2.67−70.85 5.17−31.02 1.51−32.73

16.98

17.50−68.29

saturated hydrocarbons a

Experimental results on 24 samples of Longmaxi shale by Zhang.60 H


Article

Energy & Fuels lighter compounds distribution, which is in accordance with the experimental data on Longmaxi shale gas.89 It is worth noting that the deducted molecular number of methane is determined by adjusting the residual methane numbers to match the bitumen contents with experimental data.60 Hence, the final OM model prepared for property characterization in the following sections is composed of kerogen, bitumen fractions, and residual lighter compounds. The compound distribution and atomic distribution of the final OM molecular model are listed in Tables 11 and 12, respectively. Table 11. Compounds Distribution in the Final OM Model representative molecule

molecular number

weight fraction (%)

kerogen resin asphaltene toluene 1,3-dimetylnaphtalene tricyclo[8.4.0.02,7]tetradecane methane ethane nitrogen carbon dioxide water

6 1 1 1 1 1 40 1 11 10 8

86.12 2.29 1.49 0.46 0.79 1.00 3.23 0.15 1.55 2.22 0.73

Figure 6. Kerogen model densities for different initial configurations at 23.6 MPa (from 1000K to 300 K).

3.3. Thermodynamic Properties of OM. As an important physical property of material, density is considered as one of the criteria for rationality evaluation of the proposed model. The simulated densities of the kerogen model containing ten kerogen molecules for different initial configurations following the relaxation process (from 1000 to 300 K) are shown in Figure 6. The average kerogen density of each configuration increases with decreasing temperature and gradually approaches an equilibrium value at low temperature. As seen in Figure 6, the initial configuration of kerogen structure has a big influence on the kerogen density at high temperature. However, this influence becomes weaker with further decreasing temperature. As discussed in section 3.1.1, the typical reservoir condition (360 K, 23.6 MPa) is determined for the Silurian OM according to the reported temperature and pressure gradients.39 The kerogen density range under typical reservoir condition is rather limited (1.226−1.261 g/cm3), consistent with the density range of overmature type II kerogen by Ungerer et al. (1.21−1.28 g/cm3).28 Figure 7 presents the average density of the OM versus temperature during relaxation. The influence of temperature on the average density of the OM follows the same trend as that for kerogen, but with a larger magnitude. Similar with kerogen, density equilibration is also influenced by the initial configuration. However, with the decrease of temperature, the densities of different configurations finally converge into a small range (1.154−1.201 g/cm3) under typical reservoir condition (360 K and 23.6 MPa). It should be noted that the density reduction of the OM, compared with kerogen, is mainly attributed to the inclusion of bitumen compounds and residual

Figure 7. OM model densities for different initial configurations at 23.6 MPa (from 1000K to 300 K).

lighter compounds under reservoir conditions. As shown in Table 12, the H/C ratio (0.650) of our final OM model corresponds to the early stage of gas generation.73 The density range of our OM model (1.172−1.215 g/cm3) at ambient condition is close to the experimental data on pyrite-free samples at the similar maturity by Okiongbo et al. (1.18−1.25 g/cm3).90 The isothermal compressibility of the kerogen and OM models is shown in Figure 8. The average isothermal compressibility of the kerogen and OM models under reservoir conditions is 3.8 × 10−10 and 9.1 × 10−10 Pa−1, respectively, which is close to the bulk modulus K of OM by Ahmadov et al. (K−1 = 2−7.2 × 10−10 Pa−1).91 As seen in Figure 8, the initial configuration has little influence on the isothermal compressibility of the kerogen models, while it has a big influence on the isothermal compressibility of the OM model at high temperature. However, this influence becomes negligible when the temperature drops below 600 K. Moreover, the isothermal compressibility increases with the temperature for both of the two models. For kerogen system, the isothermal compressibility is less dependent on temperature, with its values varying

Table 12. Atomic Composition of the Final OM Model atomic weight fraction (%) phase

C

H

O

N

S

H/C

O/C

N/C

S/C

final model

78.87

4.27

6.77

3.32

6.77

0.650

0.064

0.036

0.032

I


Article

Energy & Fuels

expansion coefficient for both models is negative over the whole simulated temperature range, with its magnitude being typical of a liquid.66 The characteristic of negative thermal expansion coefficient may be attributed to the anisotropic nature of the synthetic macromolecules of kerogen and OM models.95 As seen in formula 2, the thermal expansion coefficient is dependent on both the system volume and the system energy at microscopic scale. With the increase of temperature, the kerogen molecules expand horizontally along the kerogen slice structure, but contract axially along the direction perpendicular to the kerogen slice structure due to the low frequency vibration.96 The horizontal expansion of kerogen molecules increases the volume of the cubic cell, which is consistent with the results of density vs temperature (Figures 6 and 7), while the axial contraction of kerogen molecules results in the negative thermal expansion coefficient by negative vibration mode.97 3.4. Structural Characteristics of OM. 3.4.1. Properties of Porous Structure. Different from these of macroporous materials, the properties of the porous structure of microporous materials are largely dependent on the diameter of the probe molecules. Carbon dioxide and nitrogen are the most common probe molecules for estimating the porous properties in laboratory experiments, because their molecular diameters are much closer to the methane molecule. In this section, to investigate the size effect of probe molecules on the pore structure, we choose carbon dioxide, nitrogen, and methane as the probe molecules to determine the porosity and specific surface area of the structures. In addition, the innate pore structure parameters without considering the influence of probe molecule are also computed for reference. The average porosity and specific surface area of the kerogen and OM structures are listed in Table 13. Remarkably, the porosity and specific surface area decrease with the increase of the probe molecule diameters. Compared with nitrogen and methane, carbon dioxide has the smallest diameter and linear molecule shape, enabling it to enter smaller microspores, which has been proven by the experiment of low-pressure gas adsorption.98 Moreover, the porosity and specific surface area determined by probe molecules are much smaller than these of innate pore structure. This indicates that a large number of ultramicropores exist in the kerogen and OM systems, which cannot be accessed by the probe molecules. This kind of porous structure composed of micropores and ultramicropores has also been reported in disordered nanoporous carbon54,55 and swollen bituminous coals.99,100 Moreover, the differences of both porosity and specific surface area between the total pore regions and external accessible regions are small, which indicates that the pore structures are highly connected, and there are few dead pores in the structures (Figure 10). Compared with the kerogen structure, the OM structure has similar innate porosity, smaller probed porosity, and much higher specific surface area. This suggests that the OM model has a larger number of ultramicropores in comparison to the kerogen model. These ultramicropores result from a combination of factors, including the generation of oils and gas, the irregular stacking of small fractions, and the random thermal motion of the main skeletons. 3.4.2. Stacking of Polyaromatic Units. Stacking of polyaromatic units refers to the tendency for polyaromatic units to form a parallel or quasiparallel structure under π−π interaction among aromatic rings.93 This tendency is especially remarkable for overmature materials comprising a large number

Figure 8. Isothermal compressibility of kerogen model and OM model at 23.6 MPa (from 800 to 300 K).

between 2.8 and 6.5 × 10−10 Pa−1, which correspond to the general range of a liquid.92 The limited change of isothermal compressibility is closely related to the high condensation degree of aromatic sheets of kerogen structure. For OM structure, the isothermal compressibility increases smoothly with the rise of temperature when the temperature is below 600 K and then increases quickly when temperature is higher than 600 K (Figure 8). As the temperature increases, lighter molecules gradually release from the main skeleton, and the whole structure becomes less consolidated. Once after the temperature reaches the catastrophe point, the OM structure changes dramatically, and its isothermal compressibility is closer to that of a compressed gas. The catastrophe point can be estimated by analyzing the turning point on the curve of isothermal compressibility versus temperature.93 The simulated catastrophe point (Tm = 600 K) of the OM model in this work is fairly close to the reported value (673 K) for Chinese Huadian oil shale kerogen with similar element compositions by Ru et al.26 Figure 9 presents the thermal expansion coefficient of the kerogen and OM models. The thermal expansion coefficient

Figure 9. Thermal expansion coefficient of kerogen model and OM model at 23.6 MPa (from 800 to 300 K).

increases with decreasing temperature, consistent with the experimental trend for polymer blends.94 This is because the structure becomes more consolidated with decreasing temperature, which makes the thermal expansion property of the structure more sensitive to temperature.95 Similarly, the thermal expansion coefficient of the kerogen model is higher than that of the OM model. Interestingly, the thermal J


Article

Energy & Fuels Table 13. Comparison of Structural Properties of Kerogen and OM Models by Different Probe Moleculesa specific surface area (m2/g)

porosity (%) b

c

structure

CO2

N2

CH4

VDW

VDW

CO2

N2

CH4

VDWb

VDWc

kerogen OM

17.25 13.07

15.51 11.02

14.49 9.82

39.48 40.03

39.51 40.06

1389 1512

1216 1303

1114 1181

4256 4926

4263 4936

a

Probe molecules diameters: CO2, 0.33 nm; N2, 0.36 nm; CH4, 0.38 nm. bVDW refers to a probe molecule diameter of zero defined over external accessible regions. cVDW is equivalent to a probe molecule diameter of zero defined over the whole pore regions.

Figure 10. Accessible innate pore volume distribution. Skeleton atoms are wrapped by the blue surface. (a) The kerogen model. (b) The OM model.

Figure 11. View of simulation box of the kerogen and OM models at 360 K, 23.6 MPa by molecular dynamics. (a) The kerogen model. (b) The OM model. Molecule representation: pink for water molecules, blue for nitrogen molecules, and purple for carbon dioxide molecules.

of polyaromatic units. Substantial experimental investigations27,101 and molecular simulation work52,58,64 have reported the observation of stacking of polyaromatic units. In our work, we also observe this stacking phenomenon for both the kerogen structure and the OM structure, as shown in Figure 11. For the kerogen model, the kerogen molecules are highly compacted, and the polyaromatic clusters in the kerogen molecules almost present a parallel slice structure. This can be attributed to the stiff linkages between aromatic units in the kerogen molecules. The linkage units are mainly in the form of short aliphatic side chains and cycloalkanes, which prevent the folding of aromatic clusters. For the OM model, the structure compaction is relieved due to the interlude of resin, asphaltene, and heavier hydrocarbons molecules. Moreover, the structure is less ordered compared with the kerogen structure. However,

most of the polyaromatic units still present a quasiparallel structure. The pair distribution functions between carbon atoms in kerogen molecules for the kerogen structure and the OM structure are illustrated in Figure 12. The peaks at the C−C distance below 3 Å correspond to the C−C interaction within the kerogen molecules, while these above 3 Å correspond to the intermolecule C−C interaction of the kerogen molecules. As reported by Ungerer et al.,52 the first peak at 1.37 Å and the minor peak at 1.55 Å represent the C−C bond lengths in aromatic and saturated structures, respectively, whereas the peaks at 2.41 and 2.81 Å correspond to the carbons separated by two and three bonds in aromatic units. We mainly focus on the peak at 3.68 Å, which reflects the average stack spacing between polyaromatic sheets.52 As can be seen, the average spacing for the kerogen model matches fairly well with that of K


Article

Energy & Fuels

kerogen molecules, aggravating the relative slippage of polyaromatic units and intensifying the disorder of the system. The distribution of nonhydrocarbon lighter molecules is shown in Figure 11b. As can be seen, the water molecules form their own clusters in the lower density regions of the system. This behavior has also been reported by Collell et al.28 and Hu et al.102 The localization of the water clusters can be attributed to the overall hydrophobic nature of the carbon skeleton, while the aggregation of water molecules is caused by the hydrogen bonds formed between the water molecules. As seen in Figure 11b, the nitrogen molecules and carbon dioxide molecules are scattered in the nanometered pockets of the system. Most of these molecules are located close to the polar functional groups in the structure, while a small portion of them are stuck in the cabined throats equivalent to their molecular size. Moreover, it is interesting to note that the nitrogen molecules and carbon dioxide molecules tend to appear in pairs, and the carbon dioxide molecules are relatively closer to the oxygenated functional groups, while the nitrogen molecules are relatively closer to the sulfur functional groups and nitrogen functional groups. These observations reflect the influence of surface chemistry on the competitive adsorption behaviors of lighter molecules in the OM structure. The pair distribution functions between methane and atoms in kerogen (C, H, O, N, and S) in the OM model are displayed in Figure 14. These curves have multiple relatively even peaks

Figure 12. Pair distribution functions between carbon atoms in kerogen molecules for the kerogen and OM model. For clear illustration, the curve of the OM model has been shifted +1 in the Y axis.

the OM model, indicating that the introduction of bitumen fractions and lighter compounds has little effect on the spacing of polyaromatic units. Our spacing of polyaromatic sheets is larger than the 3.35 Å for graphite, which may be attributed to either the repulsive effect between polar functional groups or the steric hindrance effect of side chains. The pair distribution functions between the centroids of kerogen molecules for the kerogen model and the OM model are shown in Figure 13. The first peak sites (8.13 Å for the

Figure 14. Pair distribution functions between methane and atoms in kerogen (C, H, O, N, and S) in the OM model. Figure 13. Pair distribution functions between centroids of kerogen molecules for the kerogen and OM model.

when the interaction distances are larger than 6 Å, which indicates that the methane molecules are evenly distributed in the moderately ordered structure of the OM system. The uneven peaks below 6 Å can reflect the relative interaction intensity between the methane molecules and the atoms in the kerogen molecules. The heights of the peaks between methane and nitrogen/sulfur atoms are much larger than these of other atoms, implying that the sulfur functional groups and nitrogen functional groups in the structure have stronger affinity for methane molecules. Two smaller even peaks below 6 Å are observed between the oxygen atoms and the methane molecules, suggesting that the oxygenated functional groups have weaker affinity for methane molecules. This may be well explained by the sterical hindrance of carbon dioxide molecules located between the methane molecules and the oxygen atoms. As shown in Figure 11b, the carbon dioxide molecules are closer to the oxygen atoms than the methane molecules due to the effect of competitive adsorption, which greatly reduces the

kerogen model, 8.93 Å for the OM model) represent the average centroid distances of adjacent kerogen molecules, whose values are over two times of the spacing value of polyaromatic sheets. This suggests that relative slippage exists in the polyaromatic units for both of the two models. Moreover, the OM model has a larger average centroid distance compared with the kerogen model, indicating that the introduced bitumen fractions aggravate the relative slippage of polyaromatic units. 3.4.3. Distribution of Bitumen Fractions and Lighter Molecules. In this section, the distribution characteristics of bitumen fractions and lighter molecules in the OM model are discussed. As discussed in section 3.4.2, the bitumen fractions including resin, asphaltene and heavier hydrocarbon molecules randomly interpenetrate the main skeleton formed by the L


Energy & Fuels

■

interaction intensity between the oxygen atoms and the methane molecules.

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected].

4. CONCLUSIONS In this work, we characterized the OM in Chinese Silurian shale using organic geochemistry experiments, FTIR spectroscopy, and molecular dynamics simulations. Based on the experimental results and documented analytic data, we developed a realistic molecular model of the bulk OM, composed of kerogen, bitumen fractions, and residual lighter compounds. Thermodynamic and structural properties of the OM were discussed in detail. Experimental characterization results show the aromatic rings in the structure are highly condensed, and nearly all aliphatic carbons are in the form of short methylene groups. Aromatic sulfur dominates the organic sulfur, while ether groups dominate the oxygenated compounds. The density of the OM model is 1.19 ± 0.02 g/cm3 at ambient conditions, close to the experimental data on pyritefree samples at the similar maturity. The thermal expansion coefficient is negative, with its magnitude being typical of a liquid, which is also confirmed by the isothermal compressibility. The porous network of the OM, composed of micropores and ultramicropores, is highly connected with a few dead pores. These ultramicropores, which cannot be accessed by probe molecules, contribute a lot to the high porosity and specific surface area. Besides, a large number of ultramicropores are generated during the thermal maturation evolution of OM. Separation of bitumen fractions and lighter compounds is observed in the OM structure. The bitumen fractions randomly interpenetrate in the skeleton, aggravating the relative slippage of polyaromatic units, but having little effect on the spacing of polyaromatic units. The water molecules aggregate into small clusters in the lower density regions, confirming the hydrophobic characteristic of the OM. The nitrogen molecules and carbon dioxide molecules tend to appear in pairs. The carbon dioxide molecules are relatively closer to the oxygenated groups, while the nitrogen molecules are closer to the sulfur groups and nitrogen groups. The methane molecules are evenly distributed in the moderately ordered structure, and the sulfur groups and nitrogen groups have stronger affinity for methane molecules. The spatial distribution of lighter molecules reflects the surface chemistry effect and competitive adsorption behavior in the OM. The molecular model of the OM provides a reasonable representation of the OM of Chinese Silurian shale in terms of generic compounds, thermodynamic properties and structural characteristics. Thus, it can serve as a starting point for further theoretical investigations of the OM in Silurian shale. This study shows that molecular simulation techniques can be combined with experimental methods to gain insights into the structural and thermodynamic properties of OM.

■

Article

ORCID

Liang Huang: 0000-0002-2761-3128 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We acknowledge the financial support of the National Natural Science Foundation of China (Grant No. 51504265) and the Science Foundation for the Excellent Youth Scholars of China University of Petroleum (Beijing) (Grant No. 2462015YQ0223). We also acknowledge the Key Laboratory of Petroleum Engineering at China University of Petroleum (Beijing) for providing the experimental facilities in this work. Computer time for this study was provided by the HP High Performance Computing Cluster of the State Key Laboratory of Heavy Oil Processing at China University of Petroleum in Beijing.

■

REFERENCES

(1) Wang, X.; Zhai, Z.; Jin, X.; Wu, S.; Li, J.; Sun, L.; Liu, X. Petrol. Explor. Dev. 2016, 43, 841−848. (2) Marshall, C. P.; Love, G. D.; Snape, C. E.; Hill, A. C.; Allwood, A. C.; Walter, M. R.; Van Kranendonk, M. J.; Bowden, S. A.; Sylva, S. P.; Summons, R. E. Precambrian Res. 2007, 155, 1−23. (3) Collell, J.; Galliero, G.; Vermorel, R.; Ungerer, P.; Yiannourakou, M.; Montel, F.; Pujol, M. J. Phys. Chem. C 2015, 119, 22587−22595. (4) Curtis, J. B. Am. Assoc. Pet. Geol. Bull. 2002, 86, 1921−1938. (5) Ross, D. J. K; Bustin, R. M. Mar. Pet. Geol. 2009, 26, 916−927. (6) Weniger, P.; Kalkreuth, W.; Busch, A.; Krooss, B. M. Int. J. Coal Geol. 2010, 84, 190−205. (7) Ji, W.; Song, Y.; Jiang, Z.; Wang, X.; Bai, Y.; Xing, J. Int. J. Coal Geol. 2014, 134, 61−73. (8) Gasparik, M.; Ghanizadeh, A.; Bertier, P.; Gensterblum, Y.; Bouw, S.; Krooss, B. M. Energy Fuels 2012, 26, 4995−5004. (9) Yang, F.; Ning, Z.; Zhang, R.; Zhao, H.; Krooss, B. M. Int. J. Coal Geol. 2015, 146, 104−117. (10) Curtis, M. E; Ambrose, R. J.; Sondergeld, C. H. Structural characterization of gas shales on the micro-and nano-scales. Proceedings of the Canadian unconventional resources and international petroleum conference; Calgary, Alberta, Canada, Oct 19−21, 2010; SPE Paper 137693. (11) Pribylov, A. A.; Skibitskaya, N. A. Russ. J. Phys. Chem. A 2014, 88, 1028−1036. (12) Wang, F.; Reed, R. M. Pore networks and fluid flow in gas shales. Proceedings of the SPE annual technical conference and exhibition; New Orleans, LA, Oct 4−7, 2009; SPE Paper 124253. (13) Ertas, D.; Kelemen, S. R.; Halsey, T. C. Energy Fuels 2006, 20, 295−300. (14) Eseme, E.; Littke, R.; Krooss, B. M.; Schwarzbauer, J. J. J. Geochem. Explor. 2006, 89, 100−103. (15) Zhang, L.; LeBoeuf, E. J. Org. Geochem. 2009, 40, 1132−1142. (16) Kelemen, S. R.; Walters, C. C.; Ertas, D.; Kwiatek, L. M. Energy Fuels 2006, 20, 301−308. (17) Obliger, A.; Pellenq, R.; Ulm, F. J.; Coasne, B. J. J. Phys. Chem. Lett. 2016, 7, 3712−3717. (18) Budinova, T.; Razvigorova, M.; Tsyntsarski, B.; Petrova, B.; Ekinci, E.; Yardim, M. F. Chem. Erde 2009, 69, 235−245. (19) Bajc, S.; Cvetkovic, O.; Ambles, A.; Vitorovic, D. J. J. Serb. Chem. Soc. 2008, 73, 463−478. (20) Li, S.; Yue, C. Fuel 2003, 82, 337−342. (21) Wang, J.; Du, J.; Chang, L.; Xie, K. Fuel Process. Technol. 2010, 91, 430−433. (22) Jehlicka, J.; Beny, C. J. Mol. Struct. 1999, 480, 541−545.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.energyfuels.7b00132. Molecular model of Silurian kerogen in.xyz format. (ZIP) M


Article

Energy & Fuels (23) Zeng, Y.; Wu, C. Fuel 2007, 86, 1192−1200. (24) Wei, Z.; Gao, X.; Zhang, D.; Da, J. Energy Fuels 2005, 19, 240− 250. (25) Tong, J.; Jiang, X.; Han, X.; Wang, X. Fuel 2016, 181, 330−339. (26) Ru, X.; Cheng, Z.; Song, L.; Wang, H.; Li, J. J. Mol. Struct. 2012, 1030, 10−18. (27) Vandenbroucke, M.; Largeau, C. Org. Geochem. 2007, 38, 719− 833. (28) Collell, J.; Ungerer, P.; Galliero, G.; Yiannourakou, M.; Montel, F.; Pujol, M. Energy Fuels 2014, 28, 7457−7466. (29) Bousige, C.; Ghimbeu, C. M.; Vix-Guterl, C.; Pomerantz, A. E.; Suleimenova, A.; Vaughan, G.; Garbarino, G.; Feygenson, M.; Wildgruber, C.; Ulm, F. J.; Pellenq, R. J. M; Coasne, B. Nat. Mater. 2016, 15, 576−582. (30) Carlson, G. A. Energy Fuels 1992, 6, 771−778. (31) Faulon, J. L.; Hatcher, P. G.; Carlson, G. A.; Wenzel, K. A. Fuel Process. Technol. 1993, 34, 277−293. (32) Domazetis, G.; James, B. D. Org. Geochem. 2006, 37, 244−259. (33) Collell, J.; Galliero, G.; Gouth, F.; Montel, F.; Pujol, M.; Ungerer, P.; Yiannourakou, M. Microporous Mesoporous Mater. 2014, 197, 194−203. (34) Zhou, B.; Xu, R.; Jiang, P. Fuel 2016, 180, 718−726. (35) Zhao, T.; Li, X.; Zhao, H.; Li, M. Fuel 2017, 190, 198−207. (36) Wang, Y.; Dong, D.; Cheng, X.; Huang, J.; Wang, S.; Wang, S. Nat. Gas Ind. 2014, 34, 1−7. (37) Wang, Y.; Dong, D.; Yang, H.; He, L.; Wang, S.; Huang, J.; Pu, B.; Wang, S. Sci. China: Earth Sci. 2014, 57, 313−322. (38) Dai, J.; Zou, C.; Liao, S.; Dong, D.; Ni, Y.; Huang, J.; Wu, W.; Gong, D.; Huang, S.; Hu, G. Org. Geochem. 2014, 74, 3−12. (39) Yang, F.; Ning, Z.; Wang, Q.; Zhang, R.; Krooss, B. M. Int. J. Coal Geol. 2016, 156, 12−24. (40) Zou, C.; Dong, D.; Wang, S.; Li, J.; Li, X.; Wang, Y.; Li, D.; Cheng, K. Pet. Explor. Dev. 2010, 37, 641−653. (41) Fu, X.; Qiu, N.; Qin, J.; Tenger; Liu, W.; Wang, X. Petrol. Geol. Exp. 2013, 35, 545−551. (42) Kelemen, S. R.; Afeworki, M.; Gorbaty, M. L.; Sansone, M.; Kwiatek, P. J.; Walters, C. C.; Freund, H.; Siskin, M. Energy Fuels 2007, 21, 1548−1561. (43) Zhou, Q.; Xiao, X.; Tian, H.; Pan, L. Mar. Pet. Geol. 2014, 56, 87−96. (44) Liang, D.; Guo, T.; Chen, J.; Bian, L.; Zhao, H. Mar. Orig. Pet. Geol. 2008, 13, 1−16. (45) Pu, B.; Jiang, Y.; Wang, Y.; Bao, S.; Liu, X. Acta. Pet. Sin. 2010, 31, 225−230. (46) Wang, Y.; Dong, D.; Li, J.; Wang, S.; Li, X.; Wang, L.; Cheng, K.; Huang, J. Acta. Pet. Sin. 2012, 33, 551−565. (47) Chen, S.; Zhu, Y.; Wang, H.; Liu, H.; Wei, W.; Fang, J. Energy 2011, 36, 6609−6616. (48) Jarvie, D. M.; Claxton, B. L.; Henk, F.; Breyer, J. T. Am. Assoc. Pet. Geol. Annu. Meet. Program. 2001, 10, A100. (49) Behar, F.; Beaumont, V.; Penteado, H. D. B. Oil Gas Sci. Technol. 2001, 56, 111−134. (50) Ballice, L.; Yüksel, M.; Saglam, M.; Schulz, H.; Hanoglu, C. Fuel 1995, 74, 1618−1623. (51) Lille, Ü .; Heinmaa, I.; Pehk, T. Fuel 2003, 82, 799−804. (52) Ungerer, P.; Collell, J.; Yiannourakou, M. Energy Fuels 2015, 29, 91−105. (53) Shinn, J. H. Fuel 1984, 63, 1187−1196. (54) Pikunic, J.; Clinard, C.; Cohaut, N.; Gubbins, K. Langmuir 2003, 19, 8565−8582. (55) Pikunic, J.; Llewellyn, P.; Pellenq, R.; Gubbins, K. Langmuir 2005, 21, 4431−4440. (56) Jain, S.; Pellenq, R.; Pikunic, J.; Gubbins, K. Langmuir 2006, 22, 9942−9948. (57) Mullins, O. C.; Sheu, E. Y.; Hammami, A.; Marshall, A. G. Asphaltenes, Heavy Oils, and Petroleomics; Springer: New York, 2007; p 677. (58) Orendt, A.; Pimienta, I.; Badu, S. Energy Fuels 2013, 27, 702− 710.

(59) Vandenbroucke, M. Oil Gas Sci. Technol. 2003, 58, 243−269. (60) Zhang, D. Characteristics of the Organic-rich Shales of the Wufeng-Longmaxi Group in Wulong, Southeastern Sichuan. Master thesis, Chengdu University of Technology, Chengdu, Sichuan, China, 2015. (61) Yen, T. Div. Pet. Chem. Am. Chem. Soc. 1979, 24, 901−909. (62) Loucks, R. G.; Ruppel, S. C. AAPG Bull. 2007, 91, 579−601. (63) Feng, Z.; Liu, D.; Huang, S.; Wu, W.; Dong, D.; Peng, W.; Han, E. Petrol. Explor. Dev. 2016, 43, 705−713. (64) Yiannourakou, M.; Ungerer, P.; Leblanc, B.; Rozanska, X.; Saxe, P.; Vidal-Gilbert, S.; Gouth, F.; Montel, F. Oil Gas Sci. Technol. 2013, 68, 977−994. (65) Ungerer, P.; Rigby, D.; Leblanc, B.; Yiannourakou, M. Mol. Simul. 2014, 40, 115−122. (66) Lagache, M.; Ungerer, P.; Boutin, A.; Fuchs, A. H. Phys. Chem. Chem. Phys. 2001, 3, 4333−4339. (67) Connolly, M. L. Science 1983, 221, 709−713. (68) Chalmers, G. R. L.; Bustin, R. M. Int. J. Coal Geol. 2007, 70, 223−239. (69) Rabbani, A. R.; Kamali, M. R. J. Pet. Geol. 2005, 28, 413−428. (70) Peters, K. E. Am. Assoc. Pet. Geol. Bull. 1986, 70, 318−329. (71) Dereppe, J. M.; Boudou, J. P.; Moreaux, C.; Durand, B. Fuel 1983, 62, 575−579. (72) Wilson, M. A.; Fargue, E. L.; Gizachew, D. A. P. E. A. J. 1994, 34, 210−215. (73) Tissot, B. P.; Welte, D. H. Petroleum Formation and Occurrence; Springer Verlag: New York, 1984; Vol. 511, pp 174−75. (74) Pernia, D.; Bissada, K. K. A.; Curiale, J. Org. Geochem. 2015, 78, 52−61. (75) Wang, Z.; Liu, Y.; Ruan, Z.; Qin, K. J. China Univ. Pet. (Ed. Nat. Sci.) 1992, 16, 66−71. (76) Wang, Q.; Huang, Z.; Chi, M.; Shi, J.; Wang, Z.; Sui, Y. J. Chem. Ind. Eng. (China) 2015, 66, 1861−1866. (77) Solum, M. S.; Pugmire, R. J.; Grant, D. M. Energy Fuels 1989, 3, 187−193. (78) Iglesias, M. J.; Jiménez, A.; Laggoun-Défarge, F.; Suárez-Ruiz, I. Energy Fuels 1995, 9, 458−466. (79) Painter, P. C.; Snyder, R. W.; Starsinic, M.; Coleman, M. M.; Kuehn, D. W.; Davis, A. Appl. Spectrosc. 1981, 35, 475−485. (80) Painter, P. C.; Starsinic, M.; Coleman, M. M. Determination of functional groups in coal by Fourier transform interferometry. In: Ferraro, J. R., Basile, L. J., Eds.. Fourier Transform Infrared Spectroscopy. Academic Press: New York, 1985; pp 169−240. (81) Wang, S.; Griffiths, P. R. Fuel 1985, 64, 229−236. (82) Chen, Y.; Mastalerz, M.; Schimmelmann, A. Int. J. Coal Geol. 2012, 104, 22−33. (83) Kister, J.; Guiliano, M.; Largeau, C.; Derenne, S.; Casadevall, E. Fuel 1990, 69, 1356−1361. (84) Mastalerz, M.; Bustin, R. M. Fuel 1995, 74, 536−542. (85) Ganz, H.; Kalkreuth, W. Fuel 1987, 66, 708−711. (86) Guo, Y.; Bustin, R. M. Int. J. Coal Geol. 1998, 36, 259−275. (87) Ruan, J. S.; Bai, W. Petrol. Geol. Exp. 1988, 10, 80−85. (88) Sun, X. Spectrochim. Acta, Part A 2005, 62, 557−564. (89) Gao, B. Nat. Gas Geosci. 2015, 26, 1173−1182. (90) Okiongbo, K. S.; Aplin, A. C.; Larter, S. R. Energy Fuels 2005, 19, 2495−2499. (91) Ahmadov, R.; Vanorio, T.; Mavko, G. Leading Edge 2009, 28, 18−23. (92) Motakabbir, K.; Berkowitz, M. J. Phys. Chem. 1990, 94, 8359− 8362. (93) DeLapp, R. C.; LeBoeuf, E. J.; Chen, J.; Gu, B. H. J. Environ. Qual. 2005, 34, 842−853. (94) Weng, G.; Sun, Z.; Yang, J.; An, L. Polym. Mater. Sci. Eng. 2002, 18, 154−158. (95) Takenaka, K. Sci. Technol. Adv. Mater. 2012, 13, 1−11. (96) Sleight, A. W. Endeavour 1995, 19, 64−68. (97) Wang, L.; WANG, F.; Yuan, P.; Sun, Q.; Liang, E.; Jia, Y.; Guo, Z. Mater. Res. Bull. 2013, 48, 2724−2729. N


Article

Energy & Fuels (98) Liu, Y.; Zhu, Y.; Li, W.; Xiang, J.; Wang, Y.; Li, J.; Zeng, F. J. Nat. Gas Sci. Eng. 2016, 30, 119−126. (99) Van Niekerk, D.; Mathews, J. P. Fuel Process. Technol. 2011, 92, 729−734. (100) Mathews, J. P.; Van Duin, A. C.; Chaffee, A. L. Fuel Process. Technol. 2011, 92, 718−28. (101) Zhang, L.; Greenfield, M. L. Energy Fuels 2007, 21, 1102− 1111. (102) Hu, Y.; Devegowda, D.; Striolo, A.; Phan, A.; Ho, T. A.; Civan, F.; Sigal, R. F. SPE J. 2014, 20, 112−124.

O


Thermodynamic and Structural Characterization of Bulk Organic

Recommend Documents