Bing Liu*,
Xiaoqi Li,
Chao Qi,
Tingyi Mai,
Kaiyun Zhan*,
Li Zhao and
Yue Shen
School of Science, China University of Petroleum, Qingdao 266580, Shandong, China. E-mail: liubing19720115@gmail.com; zhanky@upc.edu.cn
First published on 4th January 2018
A thorough understanding of supercritical CO2 (scCO2) transport through nanochannels is of prime significance for the effective exploitation of shale resources and the mitigation of greenhouse gas emission. In this work, we employed the non-equilibrium molecular dynamics simulations method to investigate the pressure-driven scCO2 transport behavior through silica nanochannels with different external forces and pore sizes. The simulations reveal that the capability of scCO2 diffusion enhances both in the bulk region and the surface adsorbed layer with the increasing of pressure gradient or nanochannel size, in addition, the slip length increases nonlinearly with the external acceleration or nanochannel width increases and finally reaches a maximum value. The negative slippage occurs at lower pressure gradient or within the narrower nanochannel. Overall, it is the combined effect of strong adsorption, surface diffusion and slippage that causes the nonlinear relation between flow rate and pressure gradient or nanochannel size. The present work would provide theoretical guidance for the scCO2 enhanced shale oil/gas recovery, CO2 storage, and mass transport in nanoporous materials.
Confined fluids in nanopores demonstrate microstructural, dynamical and thermophysical behavior that differ significantly from their bulk counterparts due to the molecular asymmetry between fluid–fluid and solid–fluid interactions.17 This asymmetry leads to inhomogeneous local fluid density distributions, where the superimposition of the attractive potentials from two walls induces confined environments.18 Extensive molecular simulation studies were reported on CO2 sorption in the graphite slit-pores18–24 usually treated as the organic matrix of shale reservoirs. In addition, several molecular simulations were addressed on the adsorption structure of CO2 confined in clay-like25,26 and silica17,27 slit pores. The results revealed that the nanopore size dramatically influences the adsorption structure of CO2, which would has potentially significant effects on the transport properties of CO2 inside the nanopore.28 The diffusion of confined CO2 obviously weakens as compared with that in bulk fluid25,27 and CO2 diffusion coefficient reduces with the decrease in the pore size.25 The diffusion and migration of CO2 in silica nanopore occur predominantly along the surface due to hydrogen bonds formed between CO2 and the surface –OH group.27,29 Besides, the widely developed nanopores and ultra-low permeability in shale reservoirs raise concerns on the applicability of Darcy's law to account for mass transport in shale. However, it breaks down for flow through very narrow pores because the continuum flow hypothesis is invalid.23,28 Efforts have been made to moderate for the breakdown of Darcy's law via gas slippage theory,23 for example, the Klinkenberg effect.30,31 Such empirical corrections can not capture the complex adsorption and transport behavior of the fluid in ultra-confining porous materials.32 Falk et al. demonstrated that the non-Darcy behavior results from strong adsorption in kerogen and the breakdown of hydrodynamics at the nanoscale.32 Monteiro et al. suggested a hydrodynamic model of gas flow in nanoporous media by introducing a pressure gradient-dependent permeability of kerogen, which obeys power function relationship.33 Wang et al. reported the pressure-driven flow behavior of scCO2 confined in shale organic nanopores.24 They observed that the streaming velocity profile tends to a plug flow and the magnitude far exceeds the prediction of the no-slip Poiseuille equation. Firouzi and Wilcox investigated the gas slippage and Klinkenberg effects of CO2 confined in carbon slit pores.23 They reported a plug flow occurs at pore size less than 2 nm, a parabolic velocity profile at pore size larger than 10 nm, and the dramatic underestimation of the CO2 permeability using the bulk phase viscosity. In summary, the transport of scCO2 in inorganic nanochannel remains unclear, and there is a strong need for a reliable theoretical framework of scCO2 transport inside nanoporous matrix due to the constraint of the lowest permeability in the fluid path to the overall permeability of the shale.
Herein, non-equilibrium molecular dynamics (NEMD) simulations were conducted to investigate the pressure-driven scCO2 transport through the silica nanochannel. The effects of flow driving pressure gradient and nanochannel size on scCO2 transport were determined by exploring the adsorption, diffusion and flow of scCO2 molecules. The density profiles and residence autocorrelation functions were employed to describe adsorption. The velocity profiles, the slippage effect, and the flow rate as the function of pressure gradients or nanochannel sizes were applied to characterize the flow behavior.
Fig. 1 The MD simulation model: (a) snapshot of simulation model. (b) Schematic images of different region distribution. |
As applying an acceleration or a force on all atoms directly24,36–38 was improper to simulate desorption from surfaces, the whole nanochannel was divided into three segments along X dimension, named G-region, Buffer and S-region, respectively, as shown in Fig. 1(b). A driving pressure gradient was exerted on the system by applying an acceleration along X direction on scCO2 molecules in G-region in order to generate the Poiseuille flow through the nanochannel. A bigger acceleration corresponds to a greater pressure gradient. The effect of pressure gradient on transport behavior was examined by altering the acceleration ranged from 0.001 to 0.05 kcal (Å−1 g−1) as H is 28.50 Å. The effect of nanochannel sizes on transport behavior was studied by altering H from 7.25 to 50.00 Å as the acceleration is 0.01 kcal (Å−1 g−1). All calculations were executed in S-region. The NEMD simulations were performed to simulate mass transport for 6 ns which is enough to obtain a steady scCO2 flow.
Both the EMD and NEMD simulations were performed using large-scale atomic/molecular massively parallel simulator (LAMMPS)39 software in the NVT ensemble with a time step of 1 fs. Visualization of the dynamics process was carried out using the program VMD.40 The Nosé–Hoover41,42 thermostat was employed to control the temperature. The silica nanochannel was modeled with the CLAYFF43 force field. CO2 molecules were described by the three-site EPM2 (ref. 44) model which was successfully demonstrated to reproduce critical point of CO2 in experiment. All the VDW interactions were described by the 12–6 Lennard-Jones potential. The cutoff of non-bonded interaction was set to be 9.5 Å. Long-range electrostatic interaction was calculated by PPPM45 summation method with the accuracy of 0.0001e. The interaction between different types of atoms was determined by the Lorentz–Berthelot rule: and σij = (σi + σj)/2. All parameters were given in the ESI Table S1.†
In order to quantify the residence time of scCO2 molecules remaining near the nanochannel surface, residence autocorrelation function, CR(t), was calculated by:48–50
(1) |
The CR(t) for scCO2 at the different pressure gradients is shown in Fig. 2(b). CR(t) decreases faster as pressure gradient increases, implying a shorter residence time of scCO2 in the adsorption layer at bigger pressure gradients. The high pressure gradients could facilitate the fluid mobility and enhance the occurrence of molecular collisions, and then driving the adsorbed scCO2 molecules away from the surface. The CR(t) plateaus occurs after about 120 ps, indicating that some of the adsorbed scCO2 molecules are not desorbed within the simulation time because of the electrostatic interaction and the hydrogen bond between scCO2 and the surface.51 These scCO2 molecules that can not be desorbed from the surfaces will reduce the scCO2 flow in shale reservoirs but benefit the CO2 sequestration.
(2) |
Fig. 3 illustrates the diffusion coefficients of scCO2 at different pressure gradients and areas in the nanochannel. As pressure gradient increases, the diffusion coefficient of scCO2 in the adsorption layers increases while it decreases firstly and then rises for bulk scCO2. The variation of diffusion coefficients of scCO2 under lower pressure gradients was mainly affected by the variation in the density distribution. As shown in Fig. 2(a), the density near the surface reduces with the increasing pressure gradient as the acceleration is less than 0.007 kcal (Å−1 g−1). This indicates that the partial pressure in adsorbed layers decreases, leading to the increase of scCO2 diffusion coefficient. The increasing density of bulk scCO2 indicates the enhancement of the partial pressure in the central area of nanochannel, resulting in the decrease of scCO2 diffusion coefficient. In this case, the density distribution plays the dominant role in the diffusibility of scCO2 within the nanochannel. As the acceleration is greater than 0.007 kcal (Å−1 g−1), the largely undermined adsorption strength with increasing pressure gradient indicates the increase of the scCO2 velocity within the nanochannel. The raising velocity may lead to the decrease of scCO2 pressure, which induces the increase of the diffusion coefficient. Herein, the increasing velocity plays the pivotal role in scCO2 diffusibility inside the nanochannel. As indicated from the Stokes–Einstein relation, the increase in the diffusion coefficient implies the decrease in the viscosity. Therefore, we conclude that such changes could affect the flow of scCO2 within the nanochannel.
v(z) = az2 + bz + c | (3) |
Fig. 4 (a) The effects of external acceleration on the velocity profiles of scCO2 molecular in nanochannel. (b) Slip length as a function of the external acceleration for scCO2 molecular. |
In Fig. 4(a), the velocity profiles with the parabolic shape24,36,38 are symmetrical about the central axis of the nanochannel. The larger pressure gradients induce the greater axial velocities and sharper parabolic velocity curves. Moreover, we notice that there are non-zero velocities near the surfaces, which leads to the slippage phenomenon. The slip velocity increases as the pressure gradient increases. Here, we employed the slip length Ls to describe the slippage phenomenon, as defined by52
(4) |
Fig. 5 depicts the scCO2 flow rate as a function of pressure gradient, revealing a nonlinear relationship different from Darcy's law where permeability is defined as a material property of the rock. This indicates that Darcy's law dramatically fails to describe scCO2 transport in nanochannel. Form the analysis of Fig. 2 and 3 mentioned above, we conclude that the permeability increases and the viscosity decreases as pressure gradient increases. Hence, the ratio of the permeability to the viscosity, i.e., mobility, increases with the increasing of pressure gradient, characterizing the seepage capacity of scCO2 in nanochannel. It is observed that the mobility is linearly related to pressure gradient shown by the inset in Fig. 5. This corresponds to the hydrodynamic model of gas flow in nanoporous media where permeability is power function of pressure gradient.33
Fig. 6(b) shows the residence autocorrelation functions, CR(t), as a function of the residence time of scCO2 molecules remaining near the surface. As H increases, CR(t) decays fast, implying the shorter time that scCO2 remains in the absorbed layer. As H increases, the overlapping attractive potentials from two surfaces gradually separate from each other and finally form two single potential systems near the surfaces. This weakens the interaction between scCO2 and the surfaces resulting in the shorter residence time of scCO2 on the surface, and also implies a more frequent exchange of scCO2 appeared between the adsorption layer and bulk phase in the center of nanochannel inside the wider nanochannel. For H = 7.25 Å, CR(t) always equals to one indicating that all scCO2 molecules reside in the adsorption layer. As H further increases from 14.5 Å to 50.0 Å, CR(t) plateaus appears after about 120 ps indicating that a part of scCO2 molecules always reside near the surface within the simulation time. Value of CR(t) plateaus decrease as H increases, also showing the stronger adsorption for narrower nanochannels. As H ranges from 35.5 Å to 50.0 Å, CR(t) curves are basically coincident, suggesting that scCO2 adsorption is independent of H larger than 35.5 Å. This is mainly due to that attractive potentials from two surfaces have been effectively divided into two single potential systems near the channel surface. The results also demonstrate that the narrower nanochannel storages CO2 more securely while the wider nanochannel favors scCO2 transport.
Fig. 7 The diffusion coefficient of scCO2 molecular in different adsorbed layers behind various width of nanochannel. |
In Fig. 8(a), the scCO2 velocity profiles exhibit a parabolic shape and are symmetrical about the central axis of the nanochannel. Moreover, as H increases, the velocity curve becomes sharper and the axial velocity of scCO2 increased. This is caused by the intensive constraint on the scCO2 molecules near the surface in the smaller nanochannel, as shown in Fig. 6. However, the bulk scCO2 in the central area of the wider nanochannel are rare constrained.
In Fig. 8(b), the slip length exhibits non-monotonous increase with a maximum value as H increases. For 7.25 Å < H < 14.5 Å, the slip velocity increases fast while the velocity gradient increment is small as H expands, as indicated by the inset in Fig. 8(b). Then the ratio of the slip velocity to the velocity gradient, i.e., slip length, increases rapidly. For 14.5 Å < H < 28.5 Å, both the slip velocity and the velocity gradient increase rapidly. Ratio of slip velocity to velocity gradient enhances gradually and reaches the maximum value as H = 28.5 Å. Therefore, the slip length increases slowly and its maximum value occurs. For 35.5 Å < H < 50.0 Å, the increase of the slip velocity becomes slower while the velocity gradient increase rapidly. The ratio of them weakens gradually. Consequently, the slip length decreases slowly.
Fig. 8(b) also shows a negative slip length54–57 occurring at H less than about 10.0 Å. In this case, the overlapping between attractive potentials from the surfaces is so strong that a part of scCO2 molecules stick on the surfaces. The external pressure gradient cannot completely overcome the attractive interactions between scCO2 and the surface. This leads to the zero velocity near the surfaces, indicated the presence of scCO2 in boundary layer.
Fig. 9 describes the dependence of flow rate on H. As H is less than 30.0 Å, flow rate increases slowly due to the strong constraint on scCO2 molecules from the superimposition of attractive potentials of two surfaces. The attractive potential from the surface has a significant effect on scCO2 transport through nanochannel. As H is larger than 30.0 Å, flow rate increases rapidly because of the less constraint on scCO2 molecules where attractive potentials from two surfaces separate completely. Here, the attractive potential from the surface has little effect on the transport of scCO2 in nanochannel. Therefore, nanochannels with H less than 30.0 Å contribute to CO2 storage in shale. To quantitatively describe the dependence of the flow rate on H, the relationship between them was obtained by fitting data applying the power function polynomial
Q = AH2 + BH3 + CH4 + I | (5) |
The fitting parameters were shown by the inset in Fig. 9. Compared to the continuum flow with slip effect,52 the fourth power of H appears in the fitted formula. This could be mainly due to the fact that both the viscosity and the slip length of scCO2 in nanochannel relate to H, as indicated in Fig. 7 and 8(b).
Our simulations reveal that the structuring adsorption layers of scCO2 molecules occur within the silica nanochannel yet under the driving of pressure gradient. As either pressure gradient or nanochannel width increases, the peaks of scCO2 density profiles decrease and the density of the bulk phase increases. The higher pressure gradients result in the larger shear force, which leads to desorbing scCO2 molecules near the nanochannel surfaces. The wider nanochannel leads to the weaker overlapping of the attractive potentials from the two surfaces, attenuating the adsorption of scCO2 molecules near the surfaces. The higher pressure gradients and the wider nanochannel can promote scCO2 diffusibility for both bulk region and the adsorption layers near the surfaces. The improvement of scCO2 diffusion ability may be mainly determined by the decrease of partial pressure near the surfaces, the increase of velocity of the central area in the nanochannel, and the surface diffusion of scCO2 hopping from one site to another. The increased diffusion ability also indicates the reduction in scCO2 viscosity. The negative slippage occurs at lower pressure gradients or within the narrower nanochannel is mainly due to the stronger adsorption of scCO2 molecules near the surfaces. The slip length exhibits non-monotonic increase and existence of a maximum value with the increase in the pressure gradient or the nanochannel size. The flow rate of scCO2 through the nanochannel shows a nonlinear characteristic with pressure gradient or nanochannel size. It is also found that the mobility is proportional to the pressure gradient. This implies the breakdown of Darcy's law at nanoscale because of the strong adsorption effect, surface diffusion and slippage effect.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/c7ra11746a |
This journal is © The Royal Society of Chemistry 2018 |