Molecular Dynamics Simulation on the Mechanism of Shale Oil Displacement by Carbon Dioxide in Inorganic Nanopores

Abstract

Shale oil reservoirs feature a considerable number of nanopores and complex minerals, and the impact of nano-pore confinement and pore types frequently poses challenges to the efficient development of shale oil. For shale oil reservoirs, CO₂ flooding can effectively lower crude oil viscosity, enhance reservoir physical properties, and thereby increase recovery. In this paper, the CO₂ displacement process in the nanoscale pores of shale oil was simulated through the molecular dynamic simulation method. The performance disparity of quartz and calcite slit nanopores was discussed, and the influences of nanoscale pore types and displacement rates on CO₂ displacement behavior were further analyzed. The results demonstrate that the CO₂ displacement processes of different inorganic pores vary. In contrast, the displacement efficiency of light oil components is higher and the transportation distance is longer. Intermolecular interaction has a remarkable effect on the displacement behavior of CO₂ in nanopores. On the other hand, it is discovered that a lower displacement rate is conducive to the miscible process of alkane and CO₂ and the overall displacement process of CO₂. The displacement efficiency drops significantly with the increase in displacement velocity. Nevertheless, once the displacement speed is extremely high, a strong driving force can facilitate the forward movement of alkane, and the displacement efficiency will recover slightly.

1. Introduction

China holds an abundance of shale oil resources, making it a crucial sector for increasing domestic crude oil reserves and production. Shale reservoirs are characterized by the extensive development of nanopores, which exhibit low permeability, low porosity, and high heterogeneity [1,2]. This unique microscopic pore structure causes the micro-scale flow of oil and water and the development mechanisms in shale reservoirs to differ significantly from conventional reservoirs. Field experiments and engineering applications have shown that CO₂ can effectively recover oil and gas resources from shale reservoirs, achieving significant production enhancement [3,4,5,6]. Under reservoir conditions, CO₂ exists in a supercritical state, possessing a density similar to that of a liquid while maintaining a viscosity close to that of a gas. This highly diffusive property allows CO₂ to penetrate pores larger than its molecular diameter (0.33 nm), thereby enabling it to mobilize residual oil trapped in the deep nanopores of shale reservoirs through energy enhancement, viscosity reduction, and miscibility mechanisms. This makes CO₂ injection a practical choice for enhancing oil recovery in shale formations [7,8].

Due to the nanoscale characteristics of shale reservoirs, the interactions between fluids and micropore walls become significantly intensified within confined spaces, resulting in variations in the presence of crude oil on shale reservoir surfaces. This, in turn, affects the microscopic displacement processes of CO₂ [9,10]. Currently, experimental challenges in preparing and characterizing specimens for fluid flow in nanopores—alongside the high precision required in measurement and substantial experimental error—limit the consistency and reproducibility of results [11]. The analysis of experimental findings often relies on hypothetical models, and due to the inherent difficulties encountered in experiments, as well as the significant influence of rock walls on fluid flow within nanopores, it is critical to minimize artificial assumptions regarding solid–liquid interactions in simulations. Molecular dynamics (MD) simulations offer unique advantages in this respect, as it can explore fluid flow behavior on a molecular scale without needing conventional assumptions, relying instead on well-established potential energy interaction models in physics [12,13,14]. With advancements in high-performance computing, the challenges associated with the high computational demand and slower computation speed of MD simulations are being addressed. The adsorption strength of CO₂ and crude oil components on rock surfaces plays a crucial role in determining crude oil mobility and CO₂ displacement efficiency [15]. Previous studies have examined the flow behavior of single-component light alkanes (with carbon numbers less than 10) within nanopores, as well as the effectiveness of CO₂ in enhancing recovery from such systems [16,17,18,19].

In most studies, shale oil is often simplified as a single-component oil model, which significantly differs from the actual composition of shale oil. Moreover, another critical factor influencing shale oil displacement behavior, the type of nanopore, also requires further investigation. To address these issues, this study employs three alkane molecules (C₄, C₁₂, C₂₀) to simulate the composition of shale oil through proportional mixing, aiming to investigate the mechanisms of CO₂ static dissolution and the dynamic displacement of crude oil within nanopores. By using equilibrium molecular dynamics simulations, the swelling characteristics of crude oil in a CO₂ system within narrow slits are examined, along with the effects of varying CO₂ concentrations. Additionally, non-equilibrium molecular dynamics simulations are conducted to study the CO₂ displacement of shale oil, thereby revealing the underlying mechanisms of CO₂-driven oil displacement in shale reservoirs.

2. Model Construction

2.1. Model Parameters

The primary mineral composition of the Lucaogou shale oil reservoir in the Jimsar Sag, Xinjiang oilfield consists of quartz, feldspar, and carbonate rocks primarily composed of calcite [20]. Based on this, a quartz and calcite fracture model is constructed in this study. Molecular dynamics simulations are conducted using the NVT ensemble (namely canonical ensemble, under which condition the molecule number, volume, and temperature are determined) radial simulation in the LAMMPS (Large-scale Atomic/Molecular Massively Parallel Simulator) software (Version 2020), with the OPLS (Optimized Potentials for Liquid Simulations) force field selected. The monomers of quartz and calcite are constructed using monomers from the MS database. The crystal parameters for α-quartz are a, b, and c at 0.499, 0.499, and 0.546 nm, respectively, with angular parameters α, β, and γ at 90°, 90°, and 99°, and the ClayFF force field is applied. The crystal parameters for calcite are a, b, and c at 0.499, 0.499, and 1.706 nm, respectively, with angular parameters α, β, and γ at 90°, 90°, and 120°, and the OPLS force field is selected. The α-quartz crystal uses the (1 0 0) plane, and the structure is expanded after hydroxylation saturation. The calcite wall plane uses the (1 0 4) plane and undergoes expansion. A fracture model is created by leaving a 40 Å gap between two parallel substrates perpendicular to the Z-axis, with the constructed rock wall plane as the matrix. The fracture space is filled with a certain amount of shale oil molecules, where the shale oil components are represented by a mixture of alkanes, including C₄, C₁₂, and C₂₀, mixed in mass fractions of 5%, 45%, and 50%, respectively [21,22].

For the construction of the swelling model, a specified number of alkane molecules are added as the oil component within the fracture, and the system is filled under conditions of 330 K and 40 MPa temperature and pressure. After the filling process, dynamic equilibration is performed for 2 ns. Following equilibration, different concentrations of CO₂ are introduced around the oil components. This configuration serves as the initial model for the swelling process, with the corresponding models depicted in Figure 1.

For the construction of the displacement model, an alkane mixture was used as the oil phase and filled to the left of the slit centerline under temperature and pressure conditions of 330 K and 40 MPa. On the left side of the slit model, a CO₂ atmosphere was introduced under the same temperature and pressure conditions. A rigid pressure plate was positioned on the left side of the CO₂ system to set different displacement velocities, thereby providing the driving force for the CO₂ system. Finally, a vacuum layer with a certain spatial extent was added to the right side of the overall model to eliminate periodicity along the X-axis. A detailed model diagram is shown in Figure 2.

2.2. Model Methodology and Details

The force field component in molecular dynamics simulations primarily involves van der Waals interactions, Coulombic electrostatic interactions, bond stretching, bond angle variations, dihedral torsion, and coupling effects, all of which contribute to the system’s total potential energy. To ensure both computational efficiency and high accuracy, the total potential energy is primarily represented by the following potential energy terms [18]:

(1)$E_{\mathrm{total}}=E_{\mathrm{bonds}}+E_{\mathrm{angles}}+E_{\mathrm{dihedrals}}+E_{\mathrm{nonbonded}}$

In the equation, E_total represents the total potential energy of the system, measured in kcal/mol; E_bonds denotes bond energy; E_angles refers to angle energy; E_dihedrals represents dihedral energy; and E_nonbonded accounts for van der Waals interactions and other non-bonded interactions. The primary differences across the various molecular force fields lie in the functional forms and parameter expressions within the force field. For the OPLS force field, which is central to this study, each potential energy term can be expressed as follows:

(2)$E_{{\mathrm{bonds}}^{\rho }}={\displaystyle \sum }_{{\mathrm{bonds}}^{\rho }}{K}_r{\left(r-r_0\right)}^2$

(3)$E_{{\mathrm{angles}}^{\rho }}={\displaystyle \sum }_{{\mathrm{angles}}^{\rho }}{K}_{\theta }{\left(\theta -{\theta }_0\right)}^2$

(4)$E_{\mathrm{dihedrals}}=\frac{c_1}{2}[1+\cos (\phi )]+\frac{c_2}{2}[1+\cos (2\phi )]+\frac{c_3}{2}[1+\cos (3\phi )]+\frac{c_4}{2}[1+\cos (4\phi )]$

(5)$E_{\mathrm{nonbonded} }={\displaystyle \sum }_{i>j}\left\{4{\epsilon }_{ij}{f}_{ij}\left[{\left({\displaystyle \frac{{\sigma }_{ij}}{r_{ij}}}\right)}^{12}-{\left({\displaystyle \frac{{\sigma }_{ij}}{r_{ij}}}\right)}^6\right]+{\displaystyle \frac{q_i{q}_j{e}^2}{r_{ij}}}\right\},r<r_{\mathrm{cut}}$

In selecting force fields for this study, the CO₂ fluid phase is described using the transferable potential for phase equilibria (TraPPE) force field, which applies a rigid constraint on CO₂ molecules and excludes intramolecular interactions. Shale oil fluid is described using the optimized potentials for liquid simulations (OPLS) force field, which ensures a high computational accuracy while moderately improving simulation efficiency. The interaction parameters between different atom types are determined using the Lorentz–Berthelot mixing rule, with a cutoff radius of 1.2 nm. Long-range electrostatic forces are calculated using the particle-mesh Ewald (PME) method [19].

In the swelling model, three different concentrations of a CO₂ atmosphere are introduced around mixed alkane molecules, with CO₂ concentrations set at 1 K, 1.5 K, and 2 K based on the number of added CO₂ molecules, to investigate the swelling effects of varying CO₂ concentrations on the shale oil system. An additional model without a CO₂ atmosphere is included as a control group to evaluate the swelling effect of CO₂ on shale oil. The open-source software LAMMPS is used for simulation, with an initial 1 ns equilibrium simulation of the alkane mixture in the pore space under the NVT ensemble, while keeping the rock wall fixed, to emulate the actual presence of alkanes in reservoir pores. Following the introduction of the CO₂ atmosphere, molecular dynamics simulations are conducted for 2 ns for the entire system, with the temperature held at 330 K, a timestep of 1 fs, and data recorded every 1000 steps for analysis.

In the displacement model, a mixed alkane fluid is added to the left side of the slit pore, a CO₂ atmosphere is introduced on the same side, and a vacuum layer is applied to the right side of the slit to eliminate periodicity in the X-direction. Initially, the slit walls are fixed, and partitions are added to both sides of the mixed alkane fluid and CO₂ system. A 1 ns equilibrium simulation is conducted to simulate the adsorption state of the mixed alkanes in actual reservoir conditions. After reaching dynamic equilibrium, the restriction partition on the alkanes is removed, and a metal pressure plate is placed on the left side of the CO₂ system. Different displacement velocities are set for the pressure plate to drive CO₂ movement along the X-axis into the slit. The pressure plate velocities are set at 1 m/s, 2 m/s, 3 m/s, and 4 m/s, with the system temperature maintained at 330 K, a timestep of 1 fs, and data recorded every 1000 steps for analysis. Detailed component data are provided in Table 1.

3. Results and Discussion

3.1. Mechanism of CO₂-Driven Oil Displacement in Quartz Reservoirs

3.1.1. Static Swelling Behavior

Following molecular dynamics simulations conducted using the LAMMPS software, data visualization and processing were performed in OVITO (version 2020) (Open Visualization Tool, a kind of scientific data visualization and analysis solution for particle-based simulations). Using CO₂ at a concentration of 1 K as an example, snapshot images were analyzed to observe the process. Over the course of the simulation, alkane molecules were progressively swollen by CO₂ molecules, gradually dispersing within the slit and forming a mixed-phase state with CO₂. As shown in Figure 3a, the mixed alkanes initially localized at the center of the slit gradually swelled and mixed with CO₂, eventually occupying the entire slit. Figure 3b reveals that by 30 ps, part of the C₄ components in the mixed alkanes had dissolved into the CO₂ phase. By 200 ps, it became apparent that some alkane molecules adsorbed on the wall surface were displaced by CO₂, which subsequently occupied portions of the rock wall previously adsorbed by alkanes. By the end of the simulation, the mixed alkanes and CO₂ had formed a near-complete mixed-phase state. Low molecular weight alkanes exhibited greater dispersion, whereas high molecular weight alkanes were more concentrated and adhered more closely to the rock wall. Furthermore, the high molecular weight alkanes demonstrated an extended configuration within the CO₂ environment [23,24].

During the simulation process, the interaction energies between the rock wall and the components within the fracture were statistically analyzed, with the results shown in Figure 4. The interaction energy between the rock wall and CO₂ is significantly higher than that of the alkane molecules, indicating that CO₂ molecules dominate in the competitive adsorption process with the rock wall. The overall interaction energy of shale oil is also lower than that of CO₂ molecules. The interaction energies of C₁₂ and C₂₀ with the rock wall are similar, but both are significantly higher than that of the C₄ molecules. Snapshot images further reveal that the distribution of C₄ is more dispersed, while C₁₂ and C₂₀ are located closer to the rock wall.

During the simulation, the density distribution curves of each alkane component were also collected, as shown in Figure 5. From the figure, it can be observed that there are no multiple adsorption peaks for the alkane molecules. Instead, there is only a single adsorption peak at the wall surface, indicating that most of the alkane molecules were displaced from the wall by CO₂, forming a miscible state.

The density distribution curves of each component in atmospheres with different CO₂ concentrations are shown in Figure 6. From Figure 6a, it can be observed that as the CO₂ concentration increases, the peak density of CO₂ increases, although the shape of the curve remains largely unchanged. CO₂ is concentrated near the rock wall surface at all concentrations. For the mixed alkanes, as the CO₂ concentration increases, the concentration of alkane molecules surrounding the rock wall decreases. However, the distribution of alkane molecules shows little difference between concentrations of 1.5 K and 2 K.

3.1.2. Dynamic Displacement Process

Molecular dynamics simulations were conducted using the quartz fracture displacement model constructed earlier. The driving velocity of the metal plate on the left side was set at 1 m/s, 2 m/s, 3 m/s, and 4 m/s, and the molecular movement characteristics during the displacement process were recorded. The data on the density distribution of different molecules along the X-axis, interaction energies, and other parameters were analyzed to summarize the characteristics of CO₂ displacing alkanes in the quartz fracture and the effects of displacement velocity.

Figure 7 presents a schematic of the quartz fracture model at a displacement velocity of 1 m/s. During the displacement process, CO₂ molecules continually push the mixed alkane molecules forward, gradually displacing and stripping the molecules adsorbed on the pore wall surface. It can be observed from the figure that alkane molecules exhibit significant displacement in the X-direction over time. At the early stage of displacement, a depression forms in the center of the alkane molecules, and an oil film adsorption phenomenon can be observed on the left side of the fracture model. High molecular weight alkanes (C₁₂ and C₂₀) adsorb on the pore wall surface, obstructing CO₂ movement along the wall, negatively impacting the internal displacement within the fracture. CO₂ continuously advances towards the center of the pore, and as it progresses, the distribution range of alkane molecules increases. From Figure 8, it can be observed that the alkane molecules on the left side are squeezed by CO₂, and the displacement front is continually broken through. Thus, a clear CO₂ breakthrough phenomenon is evident, indicating that the CO₂ in the center is in an ineffective displacement state.

Figure 9 presents the density distribution of each alkane in the fracture when the metal plate is moved closer to the rock wall. As the simulation time progresses, the density peak of the alkane molecules gradually shifts forward, indicating the formation of effective displacement behavior during the simulation. When comparing the displacement distances of different alkane molecules at the same time, C₄ exhibits the largest displacement, followed by C₁₂, and C₂₀ shows the smallest displacement (C₄ > C₁₂ > C₂₀). When the metal plate moves closer to the left side of the fracture wall, the density peaks of the alkane molecules are primarily distributed near the nanopore exit, suggesting that during the simulation, CO₂ molecules have created an effective displacement process. Alkanes with different mass fractions exhibit relatively long displacement distances during the replacement process, resulting in effective displacement. However, compared to C₄, C₁₂ and C₂₀ show more residual presence in the pore. C₄ is almost entirely displaced by CO₂ out of the pore, while C₁₂ and C₂₀ show flatter density distributions within the pore, indicating the residual presence of C₁₂ and C₂₀. Therefore, during the displacement process, CO₂ displaces C₄ more thoroughly, while C₁₂ and C₂₀ have a broader and richer distribution along the pore wall surface [25].

Figure 10 presents the interaction energy between alkanes and the rock wall surface, as well as the interaction energy between CO₂ and alkanes. The former reflects the adsorption strength between shale oil and the rock surface, where a greater interaction force indicates a more stable adsorption state. In cases where the interaction energy between shale oil and the rock wall is high, CO₂ has difficulty displacing and separating the alkanes adsorbed on the pore walls, thereby destabilizing the displacement interface. The latter primarily represents the miscibility between CO₂ and the shale oil components. A greater interaction energy between alkanes and CO₂ facilitates miscibility, enabling a miscible displacement process that stabilizes the displacement front. The adsorption strength of alkane molecules on the pore wall is relatively weak. During the displacement process, the interaction between CO₂ and alkanes plays a dominant role. After CO₂ displaces the alkanes from the slit, the reduced interaction between alkanes and CO₂ leads to the increase observed at the end of the curve in the simulation.

3.1.3. Effect of Drive Velocity

Considering the significant impact of displacement velocity on the displacement effect, a series of molecular dynamics simulations were conducted at different CO₂ displacement velocities. The driving velocity of the metal plate was set to 1 m/s, 2 m/s, 3 m/s, and 4 m/s, respectively.

Figure 11 illustrates the final state of the slit model when the partition plate moves from its initial position to the region near the rock wall under different displacement velocities. At lower velocities, the alkanes are displaced over a longer distance, leading to more effective displacement. When the displacement velocity is 1 m/s, the alkanes in the pore are displaced the furthest, and the CO₂ phase in the pore dissolves more alkanes, indicating good miscibility between alkanes and CO₂, which facilitates CO₂ displacing shale oil. Similarly, it is observed that the distribution of C₁₂ and C₂₀ is more widespread and closer to the wall in the pore, reflecting that high molecular weight alkane molecules have a broader distribution, and the displacement in the X-direction is shorter. When the metal partition plate reaches the rock wall, alkane molecules in the pores with lower displacement velocities exhibit larger X-direction displacements compared to those in the pores with higher displacement velocities. This can be attributed to the fact that at a given displacement velocity, the increase in velocity results in a shorter contact time between CO₂ and alkanes, producing a stronger forward motion phenomenon. As a result, most of the CO₂ remains in an ineffective displacement state, leading to a reduction in the displacement efficiency [22].

Figure 12 presents the density distribution of alkanes in the X-direction at four different displacement velocities. It can be observed that, as the velocity decreases, the peak of the alkane density distribution appears farther along the X-direction. At a velocity of 1 m/s, the alkane molecules are distributed as far as 290 Å, while at 4 m/s, the distribution extends to only 260 Å. The distribution of C12 remains relatively widespread and abundant, showing a more extensive distribution at all four displacement velocities. Generally, lower displacement velocities favor the CO₂ displacement process, with the distribution range decreasing from 180 Å to 150 Å. As the displacement velocity increases, the displacement effect of CO₂ on alkanes weakens, resulting in a reduction in the X-direction displacement distance of the alkanes [22].

3.2. Micro-Scale CO₂ Enhanced Oil Recovery Mechanism in the Calcite Reservoir

3.2.1. Static Swelling Behavior

For CO₂ at a concentration of 1 K, as the simulation progresses, alkane molecules are continuously swollen by CO₂ molecules, gradually dispersing within the slit and forming a mixed-phase state with CO₂, as shown in Figure 13. The mixed alkane molecules and CO₂ largely form a miscible phase, with low molecular weight alkanes more dispersed and high molecular weight alkanes more concentrated, presenting a more extended configuration. This result is consistent with that observed in the quartz slit system. However, in comparison to the quartz system, the distribution of alkane molecules in the calcite system is more concentrated over the same simulation time.

The density distribution curves of each component under different CO₂ concentrations are shown in Figure 14. From the CO₂ curve, it can be observed that as the concentration increases, the peak density of CO₂ rises, but the shape of the curve remains relatively unchanged. For the mixed alkanes, as the CO₂ concentration increases, the concentration of alkane molecules around the rock wall decreases. However, the distribution of alkanes shows little variation between 1.5 K and 2 K, which is consistent with the results obtained from the quartz model.

3.2.2. Dynamic Displacement Process

Figure 15 presents snapshot diagrams of the calcite slit model under a displacement velocity of 1 m/s. These snapshots allow for a visual analysis of the displacement characteristics of each alkane molecule. CO₂ molecules continuously advance, displacing and peeling off alkane molecules adsorbed on the pore wall surface. In the early stages of displacement, depressions form between the alkane molecules, and an oil film adsorption phenomenon can be observed on the left side of the slit model. High molecular weight alkanes (C₁₂ and C₂₀) adsorbed on the pore wall surface hinder the advancement of CO₂ along the wall and negatively affect the displacement of alkanes within the slit. CO₂ continues to advance towards the center of the pore, and as it progresses, the number of alkane molecules adsorbed on the alkane wall increases. From Figure 16, it is evident that the alkane molecules on the left side are compressed by CO₂, and the displacement front is progressively overcome. As a result, a clear forward motion (finger-like displacement) of CO₂ is observed, indicating that the CO₂ at the center of the pore is in an ineffective displacement state. Additionally, the finger-like displacement phenomenon observed in the calcite model is more pronounced than that in the quartz model, suggesting that more CO₂ molecules in the calcite model are in an ineffective displacement state.

Figure 17 illustrates the density distribution of different alkanes in the slit as the metal plate moves towards the rock wall. As the simulation progresses, the peak density of the alkane molecules gradually shifts forward, indicating the formation of an effective alkane displacement phenomenon during the simulation. Specifically, when comparing the displacement of different types of alkane molecules at the same moment, it is found that C₄ has the largest displacement distance, followed by C₁₂, while C₂₀ exhibits the smallest displacement distance (C₄ > C₁₂ > C₂₀). When the metal plate approaches the left wall of the slit, the peak density of all alkane molecules is primarily concentrated near the outlet of the nanopore. This indicates that CO₂ molecules have effectively driven the movement of alkane molecules during the simulation, achieving the significant displacement of alkanes with different mass fractions and resulting in effective replacement. However, compared to C₄, C₁₂ and C₂₀ have more residuals in the pore; C₄ is almost completely displaced by CO₂ from the pore, while C₁₂ and C₂₀ exhibit relatively stable density distributions in the pore, indicating their residual presence. Therefore, it can be observed that CO₂ exhibits better displacement efficiency for C₄, while C₁₂ and C₂₀ are more widely distributed and more abundant on the pore wall, with the results consistent with those observed in the quartz wall model.

3.2.3. Effect of Drive Velocity

Considering the significant impact of displacement velocity on the displacement effect, a series of molecular dynamics simulations were conducted to vary the CO₂ displacement velocity. The displacement velocity of the metal plate was set at 1 m/s, 2 m/s, 3 m/s, and 4 m/s, respectively.

Figure 18 presents the final state of the slit model as the separator moves from its initial position towards the rock wall under different displacement velocities. It can be observed that at lower velocities, alkanes are displaced over a greater distance, showing better displacement performance. Specifically, at a displacement velocity of 1 m/s, the alkane displacement distance in the pore is the greatest. Additionally, CO₂ in the pore dissolves more alkanes, indicating a good miscibility between alkanes and CO₂, which facilitates the replacement of shale oil by CO₂. Furthermore, it was observed that in the pore, the distribution of C₁₂ and C₂₀ is more widespread and closer to the wall, reflecting that high molecular weight alkanes have a broader distribution, with shorter displacements in the X-direction.

As the displacement velocity increases, when the metal separator moves closer to the rock wall, the displacement of alkane molecules in the pore in the X-direction is larger at lower velocities. It is speculated that under a given displacement velocity, as the velocity increases, the contact time between CO₂ and alkanes is reduced, resulting in a stronger finger-like displacement phenomenon, causing most CO₂ to be in an ineffective displacement state, which in turn reduces the displacement effect. Simultaneously, the finger-like displacement phenomenon in the calcite system becomes more pronounced as the displacement velocity increases, leading to a poorer final displacement effect.

Figure 19 illustrates the density distribution of alkane molecules in the X-direction at four different displacement velocities. It can be clearly observed that as the displacement velocity decreases, the peak position of the alkane molecule density distribution in the X-direction shifts further, indicating a more widespread and concentrated distribution of alkanes at lower velocities. The distribution of C₁₂ and C₂₀ at all velocities shows a broader spread. In general, lower displacement velocities are more favorable for the CO₂ displacement process. However, as the displacement velocity increases, the ability of CO₂ to displace alkanes weakens, leading to a shorter displacement distance of alkanes in the X-direction, which is consistent with the findings in the quartz model.

4. Conclusions

(1) During the static swelling process, as the simulation time progresses, alkane molecules are continuously swollen by CO₂ molecules, gradually dispersing in the nanospace and forming a miscible phase with CO₂. Low molecular weight alkanes are more dispersed, while high molecular weight alkanes are more concentrated, with their distribution more closely aligned with the rock surface, exhibiting a more extended state in CO₂. The interaction energy between the rock surface and CO₂ is significantly greater than that between the alkane molecules, making CO₂ molecules dominate in the competition for adsorption on the rock surface.

(2) During the dynamic displacement process, alkane molecules exhibit noticeable displacement in the X-direction over time. High molecular weight alkanes (C₁₂ and C₂₀) adsorb on the pore wall surfaces, hindering the forward movement of CO₂ along the pore walls, which negatively impacts the displacement of alkanes within the nanopores. CO₂ continuously advances from the pore center, and the CO₂ finger phenomenon becomes evident, with CO₂ at the pore center in an ineffective displacement state. During displacement, C₄ experiences the greatest displacement in the X-direction, followed by C₁₂, and C₂₀ experiencing less displacement (C₄ > C₁₂ > C₂₀). CO₂ displaces C₄ more effectively, while C₁₂ and C₂₀ are more widely distributed and richer on the pore wall surface.

(3) CO₂ molecules dominate the competitive adsorption with the rock surface, with their interaction energy significantly higher than that of alkane molecules. At lower displacement velocities, alkane molecules can move further in the X-direction. As the displacement velocity increases, the contact time between CO₂ and alkanes decreases, leading to intensified fingering phenomena, with most of the CO₂ in an ineffective displacement state, thereby reducing the displacement efficiency. This effect is more pronounced in the calcite model, indicating that more CO₂ molecules are in an ineffective displacement state. In the CO₂ displacement of shale oil, controlling the appropriate displacement velocity is crucial for optimizing the displacement effectiveness. Low-speed displacement is conducive to enhancing alkane mobility and displacement efficiency, while high-speed displacement may reduce efficiency due to fingering effects and decreased effective contact time.

Author Contributions

Conceptualization, C.L., Z.X. and H.X.; methodology, F.Y. and J.Y.; software, J.Y.; validation, L.R. and J.Y.; formal analysis, T.H., C.J. and C.L.; investigation, H.X., F.Y. and C.L.; data curation, J.Y., H.X. and F.Y.; writing—original draft preparation, J.Y., H.X. and F.Y.; writing—review and editing, J.Y. and Z.W.; supervision, Z.W. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The original contributions presented in the study are included in the article, and further inquiries can be directed to the corresponding author.

Conflicts of Interest

Author Chengshan Li, Hongbo Xue, Liping Rao, Fang Yuan, Zhongyi Xu, Tongtong He and Chengwei Ji were employed by the China National Petroleum Corporation. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Table

Figure 1. Initial state diagrams of the swelling model: (a) quartz slit model, and (b) calcite slit model. (Red represents CO2, green represents C4, blue represents C12, and purple represents C20).

Figure 2. Schematic diagram of the initial construction of the slit displacement model (with red representing CO2, green representing C4, blue representing C12, and purple representing C20).

Figure 3. Swelling process evolution of the quartz model: (a) initial moment; (b) 30 ps; and (c) 200 ps. (Red: CO2, green: C4, blue: C12, and purple: C20).

Figure 4. Interaction forces between molecules in the quartz fracture model.

Figure 5. Density distribution of alkane molecules in the quartz slit model: (a) C4; (b) C12; and (c) C20.

Figure 6. Density distribution of oil components and CO2 in the quartz model under different CO2 concentrations: (a) CO2, and (b) oil.

Figure 7. Schematic of the quartz fracture model at a displacement velocity of 1 m/s (red: CO2, green: C4, blue: C12, and purple: C20): (a) 1000 ps; (b) 3000 ps; (c) 7000 ps; and (d) 8000 ps.

Figure 8. Fingering phenomenon of CO2 during the displacement process in the quartz model (green: C4, blue: C12, and purple: C20): (a) 1000 ps; and (b) 4000 ps.

Figure 9. X-direction density distribution of each alkane component in the quartz model: (a) 100 ps, and (b) final state.

Figure 10. Interaction energy distribution: (a) interaction energy between CO2 and each component; and (b) interaction energy between the quartz wall and each component.

Figure 11. Schematic of final displacement results at different displacement velocities in the quartz model: (a) 1 m/s, (b) 2 m/s, (c) 3 m/s, and (d) 4 m/s (red: CO2, green: C4, blue: C12, and purple: C20).

Figure 12. X-direction density distribution of alkane molecules at different displacement velocities in the quartz model ((a): t = 100 ps, and (b): t = final time).

Figure 13. Swelling process evolution in the calcite model: (a) 30 ps; and (b) 200 ps. (Red: CO2, green: C4, blue: C12, and purple: C20).

Figure 14. Density distribution of oil components and CO2 in the calcite model under different CO2 concentrations: (a) CO2, and (b) oil.

Figure 15. Snapshot of the calcite slit model at a displacement velocity of 1 m/s (red: CO2, green: C4, blue: C12, and purple: C20): (a) 1000 ps; (b) 3000 ps; (c) 7000 ps; and (d) 8000 ps.

Figure 16. Finger-like displacement phenomenon of CO2 during the displacement process in the calcite model (green: C4, blue: C12, and purple: C20). (a) 1000 ps; and (b) 4000 ps.

Figure 17. X-direction density distribution of alkane components in the calcite model: (a) 100 ps, and (b) 8000 ps (final).

Figure 18. Schematic of final displacement results at different displacement velocities in the calcite model: (a) 1 m/s, (b) 2 m/s, (c) 3 m/s, (d) 4 m/s (red: CO2, green: C4, blue: C12, and purple: C20).

Figure 19. X-direction density distribution of alkane molecules at different displacement velocities in the calcite model. ((a): t = 100 ps, (b): t = final time).

Table 1

Statistics of component parameters in the displacement model.

Component	Temperature,	Pressure,	X,	Y,	Z,	Density,	Number of Molecules
	K	MPa	Å	Å	Å	g/cm3
Vacuum layer	33040		40	40	60	-	-
CO2			60		60	0.65	1500
C4			-	-	-	0.66	140
C12			-	-	-	0.73	100
C20			-	-	-	0.75	80
α-SiO2			80	40	10	-	-
CaCO3			80	40	10	-	-

References

1. Wu, B.; Wu, C.; Tan, Q.; Chu, Y.; Liang, C.; Li, W.; Zhang, J.; Chen, Y.; Xu, T.; Wang, L. Accumulation conditions and key technologies for exploration and development of Changji shale oil in Jimusar sag of Junggar Basin. Acta Pet. Sin.; 2024; 45, pp. 437-460.

2. Li, Q.; Li, Q.; Han, Y. A Numerical Investigation on Kick Control with the Displacement Kill Method during a Well Test in a Deep-Water Gas Reservoir: A Case Study. Processes; 2024; 12, 2090. [DOI: https://dx.doi.org/10.3390/pr12102090]

3. Alfarge, D.; Wei, M.; Bai, B. Factors affecting CO₂-EOR in shale-oil reservoirs: Numerical simulation study and pilot tests. Energy Fuels; 2017; 31, pp. 8462-8480. [DOI: https://dx.doi.org/10.1021/acs.energyfuels.7b01623]

4. Jin, L.; Hawthorne, S.; Sorensen, J.; Pekot, L.; Kurz, B.; Smith, S.; Heebink, L.; Herdegen, V.; Bosshart, N.; Torres, J. et al. Advancing CO₂ enhanced oil recovery and storage in unconventional oil play-Experimental studies on Bakken shales. Appl. Energy; 2017; 208, pp. 171-183. [DOI: https://dx.doi.org/10.1016/j.apenergy.2017.10.054]

5. Sheng, J. Critical review of field EOR projects in shale and tight reservoirs. J. Pet. Sci. Eng.; 2017; 159, pp. 654-665. [DOI: https://dx.doi.org/10.1016/j.petrol.2017.09.022]

6. Jia, B.; Tsau, J.S.; Barati, R. A review of the current progress of CO₂ injection EOR and carbon storage in shale oil reservoirs. Fuel; 2019; 236, pp. 404-427. [DOI: https://dx.doi.org/10.1016/j.fuel.2018.08.103]

7. Song, Z.; Song, Y.; Li, Y.; Bai, B.; Song, K.; Hou, J. A critical review of CO₂ enhanced oil recovery in tight oil reservoirs of North America and China. Fuel; 2020; 276, 118006. [DOI: https://dx.doi.org/10.1016/j.fuel.2020.118006]

8. Liu, H.; Tao, J.; Meng, S.; Li, D.; Cao, G.; Gao, Y. Application and prospects of CO₂ enhanced oil recovery technology in shale oil reservoir. China Pet. Explor.; 2022; 27, pp. 127-134.

9. Su, Y.; Wang, H.; Zhan, S.; Wang, W.; Xu, J. Research progress on characterization and simulation of shale oil flow in microscale. J. Shenzhen Univ. Sci. Eng.; 2021; 38, pp. 579-589. [DOI: https://dx.doi.org/10.3724/SP.J.1249.2021.06579]

10. Wang, S.; Feng, Q.; Zha, M.; Lu, S.; Qin, Y.; Xia, T.; Zhang, C. Molecular dynamics simulation of the occurrence state of liquid alkanes in shale organic matter pores and fractures. Pet. Explor. Dev.; 2015; 42, pp. 772-778. [DOI: https://dx.doi.org/10.1016/S1876-3804(15)30081-1]

11. Wang, X.; Chen, J.; Ren, D. Research progress and prospect of pore structure representation and seepage law of continental shale oil reservoir. Pet. Reserv. Eval. Dev.; 2023; 13, pp. 23-30.

12. Hong, X.; Xu, H.; Cui, F.; Yu, H.; Wu, Y.; Wu, H.; Wang, F. Application of molecular simulation in unconventional oil and gas development. Chin. J. Comput. Mech.; 2021; 38, pp. 313-320.

13. Zhang, Y.; Li, D.; Xin, G.; Ren, S. Molecular insights into recovery of shale gas by injecting CO₂-rich industrial waste gas in kerogen slit nanopores. Gas Sci. Eng.; 2023; 115, 205025. [DOI: https://dx.doi.org/10.1016/j.jgsce.2023.205025]

14. Sun, J.; Chen, Z.; Wang, X.; Zhang, Y.; Qin, Y.; Chen, C.; Li, W.; Zhou, W. Displacement characteristics of CO₂ to CH₄ in heterogeneous surface slit pores. Energy Fuels; 2023; 37, pp. 2926-2944. [DOI: https://dx.doi.org/10.1021/acs.energyfuels.2c03610]

15. Perez, F.; Devegowda, D. A molecular dynamics study of soaking during enhanced oil recovery in shale organic pores. SPE J.; 2020; 25, pp. 832-841. [DOI: https://dx.doi.org/10.2118/199879-PA]

16. Santos, M.S.; Franco, L.F.; Castier, M.; Economou, I. Molecular dynamics simulation of n-Alkanes and CO₂ confined by calcite nanopores. Energy Fuels; 2018; 32, pp. 1934-1941. [DOI: https://dx.doi.org/10.1021/acs.energyfuels.7b02451]

17. Wang, H.; Su, Y.; Wang, W.; Jin, Z.; Chen, H. CO₂-oil diffusion, adsorption and miscible flow in nanoporous media from pore-scale perspectives. Chem. Eng. J.; 2022; 450, 137957. [DOI: https://dx.doi.org/10.1016/j.cej.2022.137957]

18. Chen, H.; Wu, Z.; Jiang, S.; Chen, Z. Mechanism of CO₂/N₂ oil exchange in tight reservoirs based on molecular dynamics simulation. Bull. Geol. Sci. Technol.; 2024; pp. 1-14. [DOI: https://dx.doi.org/10.19509/j.cnki.dzkq.tb20230456]

19. Wu, Z.; Sun, Z.; Shu, K.; Jiang, S.; Gou, Q.; Chen, Z. Mechanism of shale oil displacement by CO₂ in nanopores: A molecular dynamics simulation study. Adv. Geo-Energy Res.; 2024; 11, pp. 141-151. [DOI: https://dx.doi.org/10.46690/ager.2024.02.06]

20. Dong, X.; Xu, W.; Liu, H.; Chen, Z.; Lu, N.; Wang, W. On the replacement behavior of CO₂ in nanopores of shale oil reservoirs: Insights from wettability tests and molecular dynamics simulations. Geoenergy Sci. Eng.; 2023; 223, 211528. [DOI: https://dx.doi.org/10.1016/j.geoen.2023.211528]

21. Dong, X.; Xu, W.; Liu, R.; Chen, Z.; Lu, N.; Guo, W. Insights into adsorption and diffusion behavior of shale oil in slit nanopores: A molecular dynamics simulation study. J. Mol. Liq.; 2022; 22, 359. [DOI: https://dx.doi.org/10.1016/j.molliq.2022.119322]

22. Dong, X.; Xu, W.; Liu, H.; Chen, Z.; Lu, N. Molecular insight into the oil displacement mechanism of CO₂ flooding in the nanopores of shale oil reservoir. Pet. Sci.; 2023; 20, pp. 3516-3529. [DOI: https://dx.doi.org/10.1016/j.petsci.2023.08.006]

23. Wang, L.; Zhang, Y.; Liu, Y.; Huang, L.; Zou, R.; Meng, Z.; Lei, H. Molecular dynamics analysis on occurrence characteristics of shale oil and competitive adsorption mechanism of CO₂ and oil. J. China Univ. Pet. Ed. Nat. Sci.; 2023; 47, pp. 128-136.

24. Zhang, Y.; Wang, L.; Zou, R.; Zou, R.; Meng, Z.; Huang, L.; Liu, Y.; Lei, H. Effect of cosolvents on CO₂ displacement of shale oil and carbon storage. Pet. Explor. Dev.; 2023; 50, pp. 1318-1326. [DOI: https://dx.doi.org/10.1016/S1876-3804(24)60484-2]

25. Ren, X.; Wang, J.; Dong, H.; Sheng, P.; Xu, H.; Zhang, H. Molecular simulation of CO₂ huff and puff effects in tight oil considering component differences. Fault-Block Oil Gas Field; 2022; 29, pp. 229-233.