Preprint
Article

Study on the Effects of Wettability and Pressure in Shale Matrix Nanopores Imbibition during Shut-in Process by Molecular Dynamics Simulations

Altmetrics

Downloads

131

Views

44

Comments

0

A peer-reviewed article of this preprint also exists.

Submitted:

08 January 2024

Posted:

08 January 2024

You are already at the latest version

Alerts
Abstract
Shut-in after fracturing is generally adopted for wells in shale oil reservoirs, and imbibition occurring in matrix nanopores has been proven as an effective way to improve recovery. In this research, the molecular dynamics (MD) simulation was used to investigate the effects of wettability and pressure on nanopores imbibition during shut-in for a typical shale reservoir Jimsar. The results indicate that the microscopic advancement mechanism of the imbibition front is the competitive adsorption between “interfacial water molecules” at the imbibition front and “adsorbed oil molecules” on the pore wall. The essence of spontaneous imbibition involves the adsorption and aggregation of water molecules onto the hydroxyl groups on the pore wall. The flow characteristics of shale oil suggest that the overall push of the injected water to the oil phase is the main reason for the displacement of adsorbed oil molecules. Thus, shale oil especially the heavy hydrocarbon component in the adsorbed layer tends to slip on the walls. However, the weak slip ability of heavy components on the wall surface is an important reason that restricts the displacement efficiency of shale oil during spontaneous imbibition. The effectiveness of spontaneous imbibition is strongly dependent on the hydrophilicity of the matrix pore's wall. The better hydrophilicity of the matrix pore's wall facilitates higher levels of adsorption and accumulation of water molecules on the pore wall and requires less time for “interfacial water molecules” to compete with adsorbed oil molecules. During the forced imbibition process, the pressure difference acts on both the bulk oil and the boundary adsorption oil, but mainly on the bulk oil, which leads to the occurrence of wetting hysteresis. Meanwhile, shale oil still existing in the pore always maintains a good stratified adsorption structure. Because of the wetting hysteresis phenomenon, as the pressure difference increases, the imbibition effect gradually increases, but the actual capillary pressure gradually decreases and there is a loss in the imbibition velocity relative to the theoretical value. Simultaneously, the decline in hydrophilicity further weakens the synergistic effect of pressure difference because of the more pronounced wetting hysteresis. Thus, selecting an appropriate well pressure enables cost savings and maximizes the utilization of the formation's natural power for EOR.
Keywords: 
Subject: Physical Sciences  -   Atomic and Molecular Physics

1. Introduction

Shale reservoirs exhibit high tightness, low permeability, heterogeneity, and extensive nanopore development [1]. Conventional development techniques typically suffer from inefficiency. However, the utilization of multi-stage horizontal well volume fracturing technology for production purposes represents an effective approach [2,3]. The creation of stimulated fracture networks through hydraulic fracturing enhances the oil and gas seepage capacity of shale reservoirs. The oil present in the fractures formed during hydraulic fracturing can be extracted by utilizing the original formation pressure. Additionally, water imbibition into the matrix pores can displace shale oil present in the nanopores, enabling its extraction to the surface and facilitating EOR [4]. After the fracturing, it is common practice to shut the well for a certain period to allow for the full manifestation of the imbibition effect.
Spontaneous imbibition refers to the process in which porous media imbibe a wetting fluid and displace the original non-wetting fluid [5,6]. In low-permeability tight oil reservoirs, capillary imbibition serves as a significant mechanism for EOR, with capillary force acting as the primary driving force for spontaneous imbibition oil recovery [7]. Numerous scholars have conducted extensive experimental studies on the mechanism of imbibition oil recovery. The findings reveal that spontaneous imbibition within the shale matrix is influenced by various factors, including reservoir fluid, reservoir rock, and physical parameters such as pore structure, pore size distribution, oil-water interfacial tension, solid-liquid interface interaction, emulsification, wettability, temperature, shale oil viscosity, initial water saturation, and permeability, among others [8,9,10,11,12,13]. However, the wettability of the rock wall plays a crucial role in determining the effectiveness of spontaneous imbibition oil recovery. A more hydrophilic core leads to a faster spontaneous imbibition speed and higher recovery. Therefore, surfactants are frequently employed to improve the wettability of oil-wet reservoirs with limited imbibition oil recovery effects [14,15,16]. Secondly, due to the high rock fracture pressure associated with hydraulic fracturing, the pressure within fractures after fracturing is typically higher than the pressure within the matrix pores. As a result, both spontaneous imbibition within shale matrix pores and forced imbibition under the influence of high shut-in pressure take place. Thus, the shut-in pressure also serves as a critical factor that impacts the imbibition for EOR [17]. While traditional experimental methods can assess the imbibition oil recovery effect under various influencing factors, they have limitations in evaluating the oil recovery effects solely by calculating parameter changes from a macro perspective. Experimental research methods targeting various influencing factors include traditional experimental techniques such as weighing and volume observation methods, as well as advanced technologies, such as Nuclear Magnetic Resonance (NMR) and Computed Tomography (CT) technology [18,19,20,21]. Such experimental methods fail to directly depict the evolving characteristics of the oil and water phases within the pores and fall short of elucidating the microscopic mechanism underlying the imbibition oil displacement process.
With the rapid advancement of computer technology, MD simulation methods have gained widespread utilization in the field of oil and gas seepage [22,23]. MD simulation starts from the atomic and molecular scale and employs classical Newtonian mechanics to describe the movement of atoms, thereby capturing the dynamic characteristics of the modeled system [24]. The MD method has reached a relatively mature stage, characterized by minimal artificial influences during the simulation process and high credibility of the obtained results. Consequently, it has emerged as a crucial tool for unraveling the intricate mechanisms governing nanoscale fluid migration. Wang et al. [25] employed MD simulations to investigate the impact of nanoparticles on spontaneous water uptake in highly confined channels. They proposed a competitive mechanism for the spontaneous imbibition of nanofluids in capillaries by integrating the dynamic process of spontaneous imbibition, the contact angle of water in the capillary tube, and the relationship between displacement and time. Similarly, Yang et al. [26] conducted MD simulations to study the imbibition of octane and water into graphite and quartz nanopores. Their findings demonstrated that oil can be imbibed into dense organic nanopores of rocks at a faster rate than predicted by the classical imbibition model (Handy model). In another study, Nabin Kumar Karna et al. [27] performed large-scale atomistic simulations to investigate the capillary water imbibition phenomenon in pore silica nanochannels with varying heights (4-18 nm). They proposed an expression for the evolution of the dynamic contact angle in nanochannels over time, which, when combined with the Bosanquite equation, effectively explained the initial capillary rise. Additionally, Sang et al. [28] examined the imbibition of hydrocarbons in nano-kerogen pores. However, the aforementioned MD simulation studies on imbibition mainly focused on the imbibition of the wetting phase. They did not consider the displacement of the non-wetting phase from the pores. The two-phase imbibition mechanism, involving both the wetting and non-wetting phases, is more complex compared to the sole uptake of the wetting phase. Consequently, subsequent researchers conducted MD studies on the two-phase imbibition of oil and water, exploring various aspects such as the impact of polar components of shale oil on spontaneous imbibition [29] and the microscopic mechanism of water injection huff and puff for EOR [30]. Although some progress has been made, there is still limited research available on imbibition oil displacement, and the microscopic mechanism of the imbibition process involving both oil and water phases remains an area of active exploration.
In this study, we employed MD simulation to investigate the microscopic mechanism of imbibition oil displacement in matrix shale pores with different hydrophilicity and under different shut-in pressures. Initially, we visualized and quantitatively analyzed the spontaneous imbibition process to elucidate the microscopic advancement mechanism of the imbibition front and transportation characteristics of shale oil in different hydrophilic systems. Next, we analyzed the fluid migration characteristics during the forced imbibition process in different hydrophilic systems. Subsequently, based on the simulation results and certain theoretical derivations, we evaluated the synergistic effect of pressure difference and wettability on imbibition oil displacement.

2. Results

2.1. Shale oil migration in different hydrophilic nanopores during imbibition.

The spontaneous imbibition systems consisting of quartz walls WA and WB, along with oil-water, are named spontaneous imbibition system I (SI I) and spontaneous imbibition system II (SI II) respectively. Based on the static test results, symbols I-II represent the strongly hydrophilic system and weakly hydrophilic system respectively. During the spontaneous imbibition process, the fluid pressure surrounding the shale matrix is in equilibrium with the original matrix pore pressure[31]. Therefore, the pressure exerted on both the left and right He plates is maintained at 40MPa.
Firstly, MD simulations were performed on two spontaneous imbibition systems (SI I-II) for a duration of 1 ns. To intuitively reflect the behavior of oil-water phases in different hydrophilic systems during spontaneous imbibition, we calculated and presented the dynamic contact angle and simulated snapshots of imbibition systems I-II in Figure 1. Given the asynchronous advancement of the imbibition front on the upper and lower pore walls, we determine the dynamic contact angle as the arithmetic mean of the arithmetic means of the contact angles θd on the upper and θu lower pore walls. Meanwhile, H-bond exists if the distance of O−O is less than 3.5 Å and simultaneously the angle of H−O···O is less than 30° [32]. As shown in Figure 1a, as the wall wettability changes from strongly hydrophilic to weakly hydrophilic, the advancement distance of the imbibition front decreases significantly. Initially, water molecules located at the entrance penetrate the pore in the form of a meniscus. The progressive movement of the meniscus compels water molecules to continuously occupy nanopore space, leading to the displacement of the oil phase. As shown in Figure 1b, the dynamic contact angle in strongly hydrophilic system I is smaller than that in weakly hydrophilic system II and the dynamic contact angles in systems I and II are both less than 90°, which is in agreement with the results obtained from the static contact angle tests (Figure S5).
To intuitively reflect the axial migration behavior of the oil phase, we calculated the velocity distribution of the oil phase along the Z-axis and the two-dimensional (2-D) density distribution of the oil phase corresponding to the simulated snapshots of the imbibition systems I and II at 0, 500, and 1000 ps as shown in Figure 1a. The experimental test results show that the viscosity of crude decreases with increasing temperature. According to test results provided by the Xinjiang Research Institute of Petroleum Exploration and Development, Jimsar’s low sweet spot shale oil exhibits a high viscosity of tens to hundreds of mPa.S. However, the corresponding experimental data shows that at the reservoir temperature (315K), the viscosity of Jimsar shale oil has low viscosity close to that of water [33,34]. Therefore it can be assumed that the viscosity of shale oil under the reservoir conditions is equal to that of water. Furthermore, assuming the fluid density and viscosity do not vary spatially, The steady-state velocity profile for an incompressible laminar fluid confined between two quartz walls is parabolic and described by the classical Poiseuille equation.
v = Δ P 2 η z 2 w 2 4  
where ΔP is the pressure gradient along the flow direction (The pressure gradient ΔP during spontaneous imbibition is provided by capillary pressure); w and η are the pore width and fluid (shale oil and water) viscosity respectively. Continuous hydrodynamics always assumes the streaming velocity of fluid vanishes at the interface. However, Recent studies have confirmed that pure hydrocarbons, such as octane and decane, exhibit a slip flow behavior when flowing on hydrophilic walls [35,36]. The slip length Ls is defined as the extrapolation distance to the surface location where the fluid velocity is equal to zero.
L S = ± v z surf   d v d z z surf  
where v z surf   is the velocity of the fluid at the position of the fluid-solid interface zsurf. d v d z z surf   is the velocity gradient of the fluid at the fluid-solid interface.
Therefore, taking into account the effect of interface slip, the modified velocity distribution for the H-P equation is as follows:
v = Δ P 2 η z 2 w 2 4 w L s  
Therefore, the Ls can be obtained by fitting the velocity distribution equation (3). Before extracting the velocity distribution curve, the Gibbs dividing surface (GDS) method was employed to determine the fluid-solid interface. However, considering that the imbibition front advances asynchronously on the upper and lower walls, for the accuracy and validity of the results, the velocity data on both sides of the center of the pore were fitted separately. Then the slip length is equal to the arithmetic mean of the slip lengths on the upper and lower walls of the pore. In the Figure 2a-b, the position where z = 0 represents the center of the pore.
From the velocity distribution diagram in Figure 2a,b, it can be observed that the actual velocity matches well with the fitted velocity and the oil phase velocity near the pore wall is much smaller. The occurrence of the slip phenomenon in the oil phase during the migration process proves that the overall push of the injected water to the oil phase is the main reason for the displacement of adsorbed oil molecules and that the oil phase has a larger slip length in more hydrophilic pores. The reason for the slip is that the strong attraction between hydrocarbons and quartz walls leads to the inability of water molecules to effectively peel off the oil phase. During the spontaneous imbibition process, as shown in Figure 2c, the shale oil still existing in the pores of either a strongly hydrophilic system or a weakly hydrophilic system always maintains good layered adsorption characteristics during spontaneous imbibition. Although the methyl and hydroxyl groups of the pore wall in the weakly hydrophilic system are uniformly arranged, the different forces of the methyl and hydroxyl groups on the hydrocarbons lead to the bending of the hydrocarbons that should be adsorbed parallel to the wall. However, the adsorbed oil in the pore walls of the weakly hydrophilic system has a higher adsorption strength (Figure S1), and the bending of the adsorbed hydrocarbons does not affect the final simulation results. To further elucidate the flow behavior of shale oil during spontaneous imbibition, the mass density distribution of shale oil and the number density distribution of different hydrocarbon components after 3ns of structural optimization, and the molecular number changes of different hydrocarbon components over time in the pores during spontaneous imbibition were calculated in Figure 3. It can be found, in agreement with the 2-D density distributions, that a distinct aggregation of hydrocarbon components on both sides of the pore wall occurs (Figure 3a,d). The oil in the pore contains three distinct adsorption layers, each with a width of 5 Å. Figure 3a,d showed that the average bulk density of shale oil calculated by MD simulation from z = -11 Å to z = 11 Å was 0.76 and 0.758 g/cm3 respectively, while the density measured by the experiment was 0.81 g/cm3 [33]. This was because C20 was applied in the MD model to uniformly represent the C20+ molecules of the actual shale oil components, resulting in a lower density calculated by the model than that measured by the experiment. The oil within 5 Å of the quartz wall surface is defined as the first adsorption layer. The second and third adsorption layers are defined analogy. As shown in Figure 3b,e, the heavier hydrocarbon components generally have higher adsorption density in the first adsorption. There is no obvious peak density of C4 in the first adsorption layer in the weakly hydrophilic system II because the pore wall modified by methyl groups of the weakly hydrophilic system II has a stronger adsorption effect on hydrocarbons. Therefore, the pore walls of the weakly hydrophilic system have a stronger adsorption effect on heavy and medium components, leading to a further reduction in the amount of C4 adsorbed near the pore wall and a free state almost throughout the pore channels. During the imbibition period, the number of different hydrocarbon components in the first adsorption layer gradually decreases, and the number of hydrocarbons in the pores with more hydrophilic walls decreases faster (Figure 3c,f).
However, the reduction in the number of different hydrocarbon components in the first adsorption layer is by slip to the outside of the pore or by desorption into a bulk phase fluid. To evaluate the transport characteristics of different hydrocarbon components of adsorbed oil in pores, migration trajectories of C4, C8, and C20 in the first adsorption layer after structural optimization during spontaneous imbibition are labeled as shown in Figure 4. The migration trajectories indicate that the heavy mass components in the adsorbed layer after structural optimization are always adsorbed parallel to the wall during the imbibition process. However, the desorption of light and medium components on the pore wall occurs as they are transported forward.
However, the phenomenon of a single molecule is not enough to reflect the overall motion characteristics of the adsorption layer. Therefore, the mean-squared displacement (MSD) in the X- and Z-directions of different hydrocarbon components in the first adsorption layer after structural optimization during spontaneous imbibition as shown in Figure 5 were calculated [37].
MSD t = r i t r i 0 2
where r i ( t ) is the position of molecule i at time t and r i ( 0 ) is the initial position of molecule i.
For the hydrophilic system I, the MSD in the X-direction of the hydrocarbons is much larger than that in the Z-direction (Figure 5a,b), indicating that the first adsorbed layer of shale oil prefers to slip on the wall surface compared to desorption. As the hydrophilicity of the pore wall decreases, the MSD in the X-direction of the first adsorption layer of shale oil decreases and is similar to the MSD in the z-direction (Figure 5c,d), which means that the slip effect along the wall surface weakens and adsorbed oil molecules are neither easily desorbed nor diffused on the wall surface. Compared to the strongly hydrophilic system, in the weakly hydrophilic system, more adsorption sites are occupied by medium and heavy components, leading to slightly enhanced diffusion of light components, but there is no significant change in the MSD in the z-direction of the medium and heavy components (Figure 5a,c). Both in the strongly hydrophilic and weakly hydrophilic systems, the heavier hydrocarbon component has a smaller MSD in the Z-direction and a larger MSD in the X-direction, which suggests that the heavier component is more inclined to slip along the wall, but its slip effect is diminished. Therefore, the reason why the slip effect on the upper wall of the imbibition system I is more significant is because the adsorption effect of C20 on the lower wall is more significant than on the upper wall (Figure 3b). The slip and diffusion of adsorbed oil on the wall surface are not conducive to imbibition and oil displacement. Therefore, how to effectively strip the adsorbed oil (mainly heavy components) during imbibition is a key issue to be considered to improve oil recovery.

2.2. Microscopic advancement mechanism of imbibition front

The driving force for spontaneous advancement of the meniscus at the molecular level arises from non-bonded interactions between the wall and water molecules, which encompasses both van der Waals and Coulomb interactions. The more hydrophilic walls have higher water-wall interactions and more hydrogen bonds are formed (Figure S2a,b). However, the approximately consistent change characteristics of Figure S2a and S2b indicate that Coulomb interaction predominantly governs this process due to the highly polar hydrophilic hydroxyl groups present on the pore wall. These hydroxyl groups can form hydrogen bonds with polar water molecules, resulting in a strong Coulombic interaction between water molecules and the pore wall [38].
Furthermore, as shown in Figure 2a,b, the oil phase has a larger slip length in more hydrophilic pores, however, the slip of shale oil on the wall surface is significantly reduced in the weakly hydrophilic system. Therefore, the wettability limitation of the slip length indicates that the advancement of the imbibition front is closely related to the fluid-solid interface behavior. Previous studies have provided little insight into the microscopic advancement of the imbibition front during the spontaneous imbibition process. Therefore, based on the relevant research on water imbibition into a capillary [39] and visualization results from the MD simulation of present spontaneous imbibition, we have constructed a schematic diagram depicting the advancement process of the imbibition front during spontaneous imbibition, as shown in Figure 6c. The water in the pore is divided into three parts: “adsorbed water molecules” near the pore wall, “bulk water molecules” in the center of the pore, and “interface water molecules” at the oil-water interface. Shale oil can be divided into “adsorbed oil molecules” near the pore wall and the “bulk phase oil molecules” in the center of the pore.
The “bulk water molecules” are almost not affected by the pore wall, but the “adsorbed water molecules” are firmly adsorbed on the pore wall due to the strong Coulomb interaction. Moreover, due to the strong Coulomb interaction of the solid pore wall, the interfacial water molecules move forward while approaching the pore wall, eventually becoming “adsorbed water molecules”. Therefore, the “interfacial water molecules” are continuously spread on the pore wall to realize the continuous advancement of the imbibition front and finally the displacement of shale oil out of the pore. To verify the imbibition front advancement process in Figure 6c, two “interfacial water molecules” (Miw1 and Miw2) and one adsorbed oil molecule (C20H42) (Mabo) in imbibition system I at the initial moment were marked (Figure 6a), and then their movement trajectory from 0 to 1000 ps was observed. The visualization trajectory shows that Miw1 and Miw2 were adsorbed on the pore wall after advancing a certain distance in the pore, eventually becoming immobile adsorbed water molecules. Secondly, the observation of the trajectory of Mabo reveals that it moves approximately horizontally to the right along the pore wall, which is consistent with the trajectory of C20H42 in Figure 4. The results show that the water molecules of the weakly hydrophilic system take a longer time to adsorb to the same wall position as Miw2, which is the reason why the diffusion ability of shale oil of the weakly hydrophilic system is not significant.
To quantify the migration process of the three labeled molecules in Figure 6a, the X and Z coordinate values of Miw1, Miw2, and Mabo at different moments were counted, as shown in Figure 6d,e. For Miw1, although its X and Z coordinate values fluctuate locally, its X coordinate in general appears first gradually to increase and then stabilize, and its Z coordinate generally appears to gradually decrease and then stabilize, indicating that Miw1 flowed along the imbibition direction (X direction) while approaching and adsorbed on the pore wall. The trajectory of Miw2 is similar to that of Miw1, but the X and Z coordinate values are stabilized more quickly, suggesting that it adsorbs to the pore wall more quickly compared to Miw1. The Z coordinate value of Mabo fluctuates near the pore wall, and the X coordinate value continues to increase overall, which proves that Mabo is continuously displaced along the pore wall to the outside of the pore. The coordinate values of Miw1, Miw2, and Mabo are consistent with the trajectory in Figure 6a. Thus, the Microscopic advancement mechanism of the imbibition front during spontaneous imbibition is the competitive adsorption between “interfacial water molecules” and “adsorbed oil molecules”. To determine the reason for the slower diffusion of shale oil at the pore wall in the weakly hydrophilic system, an “interfacial water molecule” in the weakly hydrophilic system was labeled in Figure 6b, which has the same initial position and the same adsorption site on the wall as Miw2 during the imbibition process. A comparison of Figure 6a and b reveals that it takes longer for “interfacial water molecules” in a weakly hydrophilic system to adsorb to the same position as Miw2, which is the reason why adsorbed oils in the weakly hydrophilic system do not have a significant ability to diffuse on the wall.
Because of the asynchronous advancement of the imbibition front, the spontaneous imbibition process in the nanopore can be divided into two stages based on changes in replacement efficiency[40]:
d i s p l a c e m e n t   e f f i c i e n c y = N O s N N s N O s  
where NOs represent the number of oil molecules originally in the pore; NNs represent the number of oil molecules in the pore at time t.
As shown in Figure 7d, we observe the displacement efficiency increases approximately linearly in the first stage I (At 0-1621 ps). However, the increase in displacement efficiency in stages II (1621-2519 ps) slows down because the imbibition front reaches the end on the upper walls at 1621ps so that the competitive adsorption of "interfacial water molecules" of the imbibition front and "adsorbed oil molecules" on the pore wall only occurs on the lower wall. The displacement efficiency in stage I amounts to 69.23%, with a displacement efficiency per unit time of 0.043%ps-1. Moreover, the displacement efficiency in stage II amounts to 30.77%, with a displacement efficiency per unit time of 0.034%ps-1. That is to say, the displacement velocity in the second stage dropped to 0.77 times that in the first stage. The change in displacement efficiency also proves that spontaneous imbibition is the result of the competitive adsorption between "interfacial water molecules" and "adsorbed oil molecules".
In addition, the velocity distribution of water molecules in the X-direction in both spontaneous imbibition systems I and II at 500ps in Figure 7a,b show significant fluctuations at the oil-water interface. However, the velocity fluctuations of water molecules in the rest of the pore were relatively small. Therefore, the transport of bulk-phase water molecules in the pore is less affected by the wall surface and can be approximated as horizontal advancement. In other words, the essence of spontaneous imbibition oil displacement is also that water molecules outside the pore continuously adsorb and accumulate on the pore wall surface (Figure 7c).
Therefore, the essence of spontaneous oil displacement is the adsorption and accumulation of water molecules outside the matrix pores around the hydroxyl groups on the pore wall, which can be described by the radial distribution function g(r). Figure 8a,b depicts the g(r) between water molecules and hydroxyl groups on the pore wall. The first peaks of the water molecules in both systems I and II were about 0.275 nm away from the hydroxyl groups on the wall, which was close to the distance of the hydrogen bonds [41]. It proves that there are strong hydrogen bonding interactions between the water molecules and the hydroxyl groups on the wall and the water molecules gather and adsorb around hydroxyl groups because of the hydrogen bond interactions. The peak of g(r) and the area under the curve increase with time, suggesting a continuous rise in the coordination number of water molecules surrounding the hydroxyl groups on the pore walls, ie., water molecules spread forward along the pore walls. By comparing the g(r) values in Figure 8a and b, it can be found that the adsorption and aggregation of water molecules on the wall in the weakly hydrophilic system are weakened because of the decrease in the density of hydroxyl groups on the wall. The g(r) between oil and water, displayed in Figure 8c-d, is found to be small and remains stable over time, indicating the oil-water interface is stable and there is no significant miscibility between them. Thus, the displacement of oil occurs in a piston-like way. A clear oil-water interface plays a crucial role in preserving optimal interfacial tension, thereby ensuring an effective capillary pressure that serves as the driving force for spontaneous imbibition.

2.3. Forced imbibition under pressure difference

Spontaneous imbibition generally occurs under forced pressure (the difference between shut-in pressure and original pore pressure) during the shut-in period [42]. To simulate the forced imbibition process when the pressure difference is 5, 10, and 20 MPa, we apply a matrix pore pressure of 40 MPa to the He plate 2, and different shut-in pressures of 45, 50, and 60 MPa to the He plate 1. The forced imbibition systems corresponding to systems I and II are named FI I and FI II respectively.
The influence of shut-in pressure on imbibition during the shut-in period is vividly described through a series of snapshots at different moments (0, 250, 500, 750, 1000 ps) under varying pressure differences as illustrated in Figure 9a-b. These snapshots are compared with the corresponding spontaneous imbibition simulation snapshots. It can be intuitively observed that as the shut-in pressure increases, the imbibition effect also increases accordingly. Notably, the 2-D density distribution shows that the adsorbed oil still present in the pores remains consistently better-adsorbed structures during forced imbibition under pressure difference (Figure 9c).
To quantitatively characterize the stability of the adsorption layer during forced imbibition under pressure difference, the MSD in the Z-direction of different hydrocarbon components in the first adsorption layer after structural optimization during forced imbibition was calculated, as shown in Figure 10. In forced imbibition under pressure difference, the diffusion performance of the light component changes relatively more, and the diffusion ability of the medium component in the Z-direction changes relatively less (Figure 10a,b,d,e). The MSD of the heavy components in the Z-direction under the effect of pressure difference is not significant (Figure 10c,f), which means that the pressure difference has no significant desorption effect on the heavy components. Therefore, the axial pressure difference leads to changes irregularly in adsorption diffusion properties but does not effectively strip the adsorbed layer of shale oil (especially the heavy component).
Further observation of the imbibition front of the simulated snapshots of forced imbibition in Figure 9a-b shows that with the increase of pressure difference, the degree of curvature of the meniscus at the imbibition front tends to weaken (ie., The value of θd increases with the increase of shut-in pressure), and it is especially obvious at 20 MPa. In weakly hydrophilic system II, the phenomenon of boundary water lagging behind the water in the bulk phase was observed even under the pressure difference of 20 MPa, which means that the morphology of the meniscus at the imbibition front changes from concave to convex on the aqueous side. This may be due to wetting hysteresis caused by differential pressure [43].
Because the uneven distribution of mixed crude oil results in asynchronous advancement of the imbibition front and difficulty in accurately measuring the dynamic contact angle, other methods are needed to quantitatively characterize the wetting hysteresis phenomenon. The wetting hysteresis phenomenon is closely related to the velocity distribution of the fluid of the nanopore in the Z-direction. Thus, the slip velocity and peak velocity variations can quantitatively characterize the wetting hysteresis phenomenon and can be obtained by fitting equation (3), As shown in Figure 11a. Considering the asynchronous advancement of the imbibition front and the accuracy of the results, both the slip and peak velocities are obtained by arithmetic averaging, as shown in equations (6) and(7).
v s l i p = v s l i p l o w e r + v s l i p u p p e r 2
v p e a k = v p e a k l o w e r + v p e a k u p p e r 2
where v s l i p l o w e r , v s l i p u p p e r , and v s l i p are the slip velocities of the shale oil on the lower and upper walls, and their arithmetic mean respectively; v p e a k l o w e r , v p e a k u p p e r , and v p e a k are the peak velocities of shale oil on the lower and uppper sides of the center of the pore, and their arithmetic mean respectively. The fitting results for the systems I and II are presented in Figure 11b and 11c respectively.
Linear fitting again of the velocity obtained by fitting Equation (3) can obtain the slopes Kpeak and Kslip respectively representing the speed changes of peak velocity and slip velocity with pressure difference. Both peak and slip velocities increase faster for the strongly hydrophilic system I than for the weakly hydrophilic system II, representing a more significant promotion of the imbibition in the strongly hydrophilic system by the pressure difference. In both strongly and weakly hydrophilic systems, Kpeak is greater than Kslip, leading to a decrease in the dynamic contact angle (i.e., wetting hysteresis). Further comparison of the ratios of Kpeak and Kslip (i.e., K p e a k K s l i p ) for the different imbibition systems, respectively, the ratio K p e a k K s l i p of the strongly hydrophilic system I is greater than that of the weakly hydrophilic system II, indicating that the difference between the peak and slip velocity in the strongly hydrophilic system changes less than that of the weak hydrophilic system, which means that the wetting hysteresis phenomenon of the weak hydrophilic system is more significant.
To further clarify the wettability and pressure difference correlation because of the wetting hysteresis, We have carried out certain theoretical derivations and analyses. Based on the reference [44], when forced imbibition occurs after shut-in, the equilibrium of forces on the fluid in the confined capillary can be expressed as:
F c + F p + f + F i + F g = 0
where Fc is the capillary force; ΔF is the displacement force. f is the viscous force generated by the friction between fluid molecules and between fluid molecules and the solid pore wall; Fi is the inertial force; F g is gravity; Because the matrix pore width of tight shale reservoirs is extremely small and much smaller than the imbibition height, the influence of F g can be ignored [45]. Secondly, for pores with extremely small radii, viscous force f dominates and the inertial force F i can be ignored [46]. Then equation (8) can be simplified to:
F c + F p + f = 0
Assuming that the width between two walls, length in the Y-direction, and length in the X-direction of the slitpore are w, l, and L respectively, It is possible to convert F p into wl∆P. According to the Young-Laplace equation [47], the Fc in a single pore is:
F c = 2 σ c o s θ w w l = 2 σ l cos θ
where θ is the contact angle of the oil-water-wall three-phase system, which is assumed to be equal to the static contact angle and is not affected by other factors; σ is the oil-water interfacial tension. The total f of oil and water in a single pore is:
f = 2 x l τ w + 2 l L x τ o  
where x is the position of the oil-water interface in the X-direction; τw and τo are the shear stress of water and oil respectively, which can be expressed as [48]:
τ w = 12 μ w w v τ o = 12 μ o w v
where v is the imbibition speed, μw and μo are the viscosity of water and oil respectively. After bringing equation (10-12) into equation (9) and rewriting v as the derivative dx/dt of x with respect to time t, equation (9) is transformed into:
P w 2 + 2 σw cos θ 24 x μ w + L x μ o d x d t = 0  
After integrating equation (13) under the initial condition of x = 0 and t = 0, we get:
12 μ w μ o x 2 + 24 μ 0 L x Δ P w 2 + 2 σ w c o s θ t = 0
Additionally, it should be noted that the shale oil model used in this study does not account for components such as C20+ and asphaltenes, implying that its viscosity could be lower under reservoir conditions. As previously mentioned, the oil displacement efficiency primarily exhibits a linear change during the imbibition process before the imbibition front reaches the outlet. Considering the low viscosity of the Jimsar shale oil under reservoir conditions, for simplification, the viscosity of shale oil under reservoir conditions can be approximated as μ = μw = μo = 0.36 mPa·s, ie., the viscosity of water under the same conditions obtained from the NIST database. Consequently, equation (14) can be simplified as follows:
24 μ L x Δ p w 2 + 2 σ w c o s θ t = 0
Then the imbibition velocity at pressure difference is
v = x t = Δ P w 2 + 2 σ w c o s θ 24 μ L
The pressure differences are equal to 0, 5, 10, and 20MPa respectively in this study. When the pressure difference is equal to 0, Equation (16) is converted to v = σwcosθ⁄12uL, indicating that the spontaneous imbibition velocity is constant during spontaneous imbibition under reservoir conditions. This may be the reason why the oil displacement efficiency changes approximately linearly under reservoir conditions before the imbibition front reaches the outlet.
The imbibition velocity when the pressure difference is equal to 0 is called spontaneous imbibition velocity (vSI), and the imbibition velocity when the pressure difference is equal to 5, 10, and 20MPa is called forced imbibition velocity. Assuming that when capillary imbibition occurs, the dynamic contact angle does not change with pressure difference, the ratio (R) of forced imbibition to spontaneous imbibition velocity at pressure difference is as follows:
R = v F I v S I = Δ P w 2 2 σ w c o s θ + 1  
From equation (17), assuming that the three-phase contact angle θ is not affected by the pressure difference and is always equal to the static contact angle, the R of the weakly hydrophilic system should be greater than that of the strongly hydrophilic system at different pressure differences. However, the R obtained through MD simulations as shown in Table 1 shows that the R of the strongly hydrophilic system is larger than that of the weakly hydrophilic system, which further proves that the pressure leads to the occurrence of the wetting hysteresis phenomenon.
Some parameters of Equation (16) are: The interfacial tension σ under reservoir conditions obtained by MD simulation is 55.44mN/m (Figure S3) (The results indicate a good match between the simulation and experimental value of 53.2 mN/m); The static contact angles θ for the strongly and weakly hydrophilic walls are 53.9° and 78.97° respectively; The length L and width w of the pore is 15 and 6 nm. The imbibition speed at different pressure differences obtained by substituting the value of these parameters mentioned above into Equation (16) is renamed as the theoretical imbibition speed v t h under different pressure differences. Meanwhile, we rename the imbibition velocity obtained through MD simulation in Table 1 as v M D . The calculation method of the imbibition velocity obtained through MD simulation is: taking the average value of the fitted velocity curve in Figure 11a. Although the theoretically calculated and MD simulated values do not necessarily correspond to the real imbition velocity, their change characteristics can quantitatively reflect a certain imbibition mechanism.
Figure 12a,b illustrates the values of v M D and v t h , as well as the difference (Δv) (∆v= v M D - v t h ) between the them at different pressure differences. The results show basically that v M D and v t h exhibit the same order of magnitude, confirming the reliability of the simulation outcomes. In imbibition systems I and II, as the pressure difference increases, both v M D and v t h gradually increase, while Δv continues to decrease. The results demonstrate that an increase in pressure difference can enhance the imbibition effect, while also leading to a reduction in the effective imbibition velocity relative to the theoretical value because of the wetting hysteresis phenomenon. Moreover, ∆v decreases more rapidly in the weakly hydrophilic system relative to the strongly hydrophilic system, suggesting that a more pronounced wetting hysteresis phenomenon leads to a more pronounced loss of the imbibition velocity relative to theoretical values. The loss of imbibition velocity relative to the theoretical velocity is caused by the reduction of capillary pressure due to wetting hysteresis. Consequently, excessively high imbibition pressure impedes the maximization of economic benefits and reduces the underutilization of natural power.

3. Models and Methodology

3.1. Model system

All molecular models were constructed through Materials Studio software developed by Accelrys Company in the United States. The molecular model constructed in this article mainly includes solid quartz walls, nanopores composed of quartz walls, and fluid components (shale oil and water).
Scholars [49,50] have conducted mineral components analysis on shale core samples from the Jimsar reservoir by using a device of X-ray diffraction. The results show that the mineral composition of the Lucaogou Formation shale matrix is complex and diverse, but quartz has the highest relative content and plays an important role. Quartz is often used to simulate inorganic nanopores[40,51,52,53]. Therefore, quartz is used in this work to represent shale reservoir matrix pores. The initial quartz unit cell SiO2 Quartz was taken from the structure database of Materials Studio software and then was cleaved along the ( 0 01 ¯ ) crystal face. After being cleaved, the crystal face was replicated to get the final quartz wall. The complexity and randomness of the mineral distribution in shale reservoirs lead to a heterogeneous distribution of wettability as well [54]. To characterize the uneven wettability of shale reservoirs, silica was modified concerning existing research methods. The hydrophilic modification of the wall surface (named WA) was first realized by adding hydrophilic group hydroxyl (-OH) to all chemically unsaturated Si atoms of the ( 0 01 ¯ ) crystal surface (surface in contact with oil), with a surface hydroxyl density of 9.6 nm−2 which was consistent with the results of the crystal chemistry calculations (5.9–18.8 nm−2 ) [55]. Then hydroxyl hydrogens were substituted with a certain amount of hydrophobic group methyl (-CH3) to obtain the other wall surface (named WB) with different wettability. To make the wall surface as homogeneous as possible, all the hydroxyl and methyl groups were evenly distributed on the wall surface of the quartz. The hydroxyl and methyl density on the surface in contact with the oil of the modified quartz wall was 7.2 nm−2 and 2.4 nm-2. The back surface and left and right end surfaces of the quartz walls were all fully hydroxylated. Considering that spontaneous imbibition in oil-wet pores cannot occur effectively, only hydrophilic matrix pores of the reservoir are considered in this paper. Each wall has a thickness of 21.6 Å and dimensions in the X and Y directions of 150Å and 29.5 Å, respectively. Quartz walls with different hydrophilicity are presented in Figure S4.
Before conducting MD imbibition simulations, the static wettability of different walls was first quantitatively characterized (The corresponding results are presented in Figure S5). Then on this basis, the molecular model of imbibition oil displacement was constructed as shown in Figure 1. The nanopore of the imbibition model consisted of two quartz walls with the same modification as shown in Figure S4. The two walls of the pore were placed parallel to the X-Y plane, and the vertical distance between their wall surfaces in the Z direction was 60 Å. For the shale oil in the pore, considering the complexity of shale oil components, a single oil component cannot truly reflect the occurrence and flow state of shale oil in the pore. According to Table S2, in this study, a ternary system composed of n-Butane (C4H10), n-Octane (n-C8H18), n-Eicosane (n-C20H42) was selected to represent shale oil in the lower “sweet spot” of the Lucaogou Formation in the Jimsar Sag, Xinjiang oilfield. n-Butane (n-C4H10), n-Octane (n-C8H18), and n-Eicosane (n-C20H42) respectively represent the light, medium and heavy components in Jimusar shale oil. The molar fractions and number of molecules of each hydrocarbon component of the shale oil model in the pore are shown in Table S2. Water molecules representing the fracture, i.e., fluid surrounding the matrix, are placed to the left of the pore. Before the imbibition oil displacement simulation, the shale oil was fully adsorbed on the pore wall through 3ns equilibrium molecular dynamics (EMD) at reservoir conditions (353.15K and 40MPa), thereby obtaining a reasonable shale oil density. During the imbibition simulation, He plates were placed on the left and right sides of the oil-water-wall system, and the imbibition oil displacement process was controlled by controlling the pressure attached to the He plate. In this work, it is assumed that the shut-in pressure is equal to the fluid pressure around the matrix. The pressure on He plate 2 always maintains the matrix pore pressure (Pm) of 40MPa, while the pressure exerted on He plate 1 represents the shut-in pressure (Ps). When the pressure on He plate 1 is equal to a matrix pore pressure of 40MPa, water is spontaneously imbibed into the pore. When the pressure applied to He plate 1 exceeds 40MPa, it is the process of forced imbibition oil displacement. The simulation time of imbibition oil displacement is 1ns.

3.2. Simulation Details

Different force fields have different characteristics and applicable conditions. Therefore, choosing an appropriate force field for a specific research system is key to the accuracy and reliability of MD simulation results. In this paper, we use different force field potential energy parameters for different molecular groups. The ClayFF force field was developed by Randall T et al [56] for hydrates and multi-component minerals and their interaction with fluid interfaces and has been successfully applied to hydroxylated silica surface-fluid systems [57,58,59,60]. This paper uses the ClayFF force field to describe the silica crystal and surface hydroxyl groups. Referring to the research of Aleksandr Abramov, the remaining surface methyl groups are described by the DREIDING force field [61]. The OPLS-AA force field [62] is used to describe n-hydrocarbons in shale oil. The OPLS force field can accurately describe hydrocarbon molecules’ thermodynamic properties and structural characteristics [63]. The SPC/E potential energy model [64] is used to describe water molecules. SPC/E can reproduce the dynamic physical behavior of water molecules on the interface and has good compatibility with the CLAYF force field and OPLS-AA force field [65,66,67,68]. The helium plate and quartz base plate are considered rigid bodies. Periodic boundary conditions are applied in all directions. This article uses Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) open-source program software to perform MD simulations [69], with a time step of 1fs. The simulation results were visualized using VMD software [70]. Intermolecular non-bonded interactions (Enonbonded) can be described by the standard 12-6 Lennard-Jones (LJ) potential (ELJ ) [54] and Coulomb potential respectively (Ecoulomb) [58,71].
E n o n b o n d e d = E L J + E c o u l o m b = 4 ε i j σ i j r i j 12 σ i j r i j 6 + q i q j 4 π ε 0 r i j
where rij stands for the distance between two atoms i and j; εij and σij represent the LJ well depth and the zero-potential distance between atoms i and j respectively; ε0 represent the dielectric constant of a vacuum; And qi and qj represent the charge of atoms i and j respectively. The interactions between different atoms are calculated using the Lorentz-Berthelot mixing criterion:
This study uses MD simulation to elucidate the microscopic mechanism of water imbibition and oil displacement in the different wettability pores of the shale matrix under different shut-in pressures. The main conclusions are as follows:
σ i j = 1 2 σ i i + σ j j
ε i j = ε i i ε j j

4. Conclusions

This study uses MD simulation to elucidate the microscopic mechanism of water imbibition and oil displacement in the different hydrophilic pores of the shale matrix under different shut-in pressures. The main conclusions are as follows:
(1) The slip phenomenon suggests that the strong attraction between hydrocarbons and quartz walls leads to the inability of water molecules to effectively peel off the oil phase. The heavy hydrocarbon component in the adsorbed layer tends to slip on the walls and the light component tends to desorb from the pore wall. However, the heavier hydrocarbon components have poorer slip ability on the wall. Weak wall hydrophilicity and slip of shale oil on the wall surface can significantly reduce the efficiency of imbibition displacement. How to effectively strip the adsorbed oil (especially heavy components) during imbibition is a key issue to be considered to improve oil recovery.
(2) The displacement of oil occurs in a piston-like way during spontaneous imbibition. The microscopic advancement mechanism of the imbibition front is the competitive adsorption between “interfacial water molecules” at the imbibition front and “adsorbed oil molecules” on the pore wall. The essence of spontaneous imbibition is the adsorption and aggregation of water molecules onto the hydroxyl groups present on the pore wall. The better hydrophilicity of the matrix pore's wall facilitates higher levels of adsorption and accumulation of water molecules on the pore wall and requires less time for “interfacial water molecules” to compete with “adsorbed oil molecules” because interfacial water molecules have a stronger competitive adsorption effect on the more hydrophilic pore.
(3) During the forced imbibition process, as the pressure difference increases, the imbibition effect gradually increases. The pressure difference acts mainly on the bulk oil thus causing wetting hysteresis to occur. Meanwhile, shale oil still existing in the pore always maintains a good stratified adsorption structure. During forced imbibition, the actual capillary pressure gradually decreases because of the wetting hysteresis phenomenon and there is a loss in the imbibition velocity relative to the theoretical value. Simultaneously, the decline in hydrophilicity further weakens the synergistic effect of pressure difference because of the more pronounced wetting hysteresis. Consequently, excessively high shut-in pressure impedes the maximization of economic benefits and reduces the underutilization of natural power.
The study presented in this article investigates the microscopic mechanism and process of imbibition oil displacement in matrix shale pores under different wettabilities and shut-in well pressures. The findings offer valuable insights into enhancing imbibition oil recovery methods, optimizing well stew time, and improving production in shale reservoirs.

Supplementary Materials

The following supporting information can be downloaded at the website of this paper posted on Preprints.org. Figure S1: The interaction energy of oil and wall during structural optimization of the oil-pore systems of imbibition systems I and II. Figure S2:(a) The interaction energy between water molecules and wall, and (b) Hydrogen bond number as a function of time in spontaneous imbibition systems I-II respectively.Figure S3:Molecular model for calculating the interfacial tension of pure oil and water systems. Figure S4: Quartz wall sizes and different modified wall surfaces with different wettability. Figure S5:(a,b) Initial and equilibrium configuration for static contact angle testing respectively. (c) Schematic diagram of static contact angle calculation method. (d) Static wetting angle of WA and WB. Figure S6:The mass fraction of hydrocarbons with different numbers of carbon atoms of (a) JHW05815 Oil Sample, (b) JHW07121 Oil Sample, and (c) J41 Oil Sample. Table S1: The average mass fraction of light, heavy, and medium hydrocarbon components of the Jimsar shale oil. Table S2: The mole fractions and number of each component of the shale oil model.

Author Contributions

Writing—original draft, writing—review and editing, conceptualization, methodology, W.J. and W.L.; data curation, formal analysis, W.J. and N.J.; supervision, W.J. and N.J.; writing—review and editing, W.J. and X.L.; data curation, W.J. and L.W.; validation, K.W.; software, W.L.; data curation, Y.M.; visualization, W.J.; All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Major Science and Technology project of the CNPC in China (grant No. 2021ZZ05).

Institutional Review Board Statement

The study did not involve humans or animals.

Informed Consent Statement

The study did not involve humans.

Data Availability Statement

The raw/processed data required to reproduce these findings cannot be shared at this time, as the data also forms part of an ongoing study..

Acknowledgments

The authors are grateful for the financial support of the CNPC in China.

Conflicts of Interest

The authors declare that the publication of this paper has no conflict of interest.

References

  1. Xia, D.; Yang, Z.; Gao, T.; Li, H.; Lin, W. Characteristics of micro-and nano-pores in shale oil reservoirs. Journal of Petroleum Exploration and Production 2021, 11,, 157–169. [Google Scholar] [CrossRef]
  2. Dai, C.; Wang, K.; Liu, Y.; Fang, J.; Zhao, M., Study on the reutilization of clear fracturing flowback fluids in surfactant flooding with additives for Enhanced Oil Recovery (EOR). PloS one 2014, 9, (11), e113723. [CrossRef]
  3. Shaibu, R.; Guo, B. The dilemma of soaking a hydraulically fractured horizontal shale well prior to flowback–A decade literature review. Journal of Natural Gas Science and Engineering 2021, 94, 104084. [Google Scholar] [CrossRef]
  4. Hu, Y.; Ren, F.; Li, J.; Wu, Z.; Peng, H.; Hou, J. Effect of dynamic imbibition on the development of ultralow permeability reservoir. Geofluids 2021, 2021, 1–12. [Google Scholar] [CrossRef]
  5. Sobhani, A.; Ghasemi Dehkordi, M., The effect of nanoparticles on spontaneous imbibition of brine into initially oil-wet sandstones. Energy Sources, Part A: Recovery, Utilization, and Environmental Effects 2019, 41, (22), 2746-2756. [CrossRef]
  6. Morrow, N. R.; Mason, G., Recovery of oil by spontaneous imbibition. Current Opinion in Colloid & Interface Science 2001, 6, (4), 321-337. [CrossRef]
  7. Nguyen, D.; Wang, D.; Oladapo, A.; Zhang, J.; Sickorez, J.; Butler, R.; Mueller, B. In Evaluation of surfactants for oil recovery potential in shale reservoirs, SPE Improved Oil Recovery Symposium, 2014; OnePetro.
  8. Hatiboglu, C. U.; Babadagli, T., Experimental and visual analysis of co-and counter-current spontaneous imbibition for different viscosity ratios, interfacial tensions, and wettabilities. Journal of Petroleum Science and Engineering 2010, 70, (3-4), 214-228.
  9. Sukee, A.; Nunta, T.; Haruna, M. A.; Kalantariasl, A.; Tangparitkul, S. Influence of sequential changes in the crude oil-water interfacial tension on spontaneous imbibition in oil-wet sandstone. Journal of Petroleum Science and Engineering 2022, 210, 110032. [Google Scholar] [CrossRef]
  10. Darvishi, H.; Goudarznia, I.; Esmaeilzadeh, F., Effects of rock permeability on capillary imbibition oil recovery from carbonate cores. 2010.
  11. Strand, S.; Puntervold, T.; Austad, T., Effect of temperature on enhanced oil recovery from mixed-wet chalk cores by spontaneous imbibition and forced displacement using seawater. Energy & Fuels 2008, 22, (5), 3222-3225. [CrossRef]
  12. Zhao, H.; Kang, W.; Yang, H.; Zhang, H.; Zhu, T.; Wang, F.; Li, X.; Zhou, B.; Sarsenbekuly, B.; Aidarova, S., Imbibition enhancing oil recovery mechanism of the two surfactants. Physics of Fluids 2020, 32, (4), 047103. [CrossRef]
  13. Wang, F.; Wang, L.; Jiao, L.; Liu, Z.; Yang, K. Experimental Mechanism for Enhancing Oil Recovery by Spontaneous Imbibition with Surfactants in a Tight Sandstone Oil Reservoir. Energy & Fuels 2023. [Google Scholar] [CrossRef]
  14. Wang, Y.; Xu, H.; Yu, W.; Bai, B.; Song, X.; Zhang, J. Surfactant induced reservoir wettability alteration: Recent theoretical and experimental advances in enhanced oil recovery. Petroleum Science 2011, 8, 463–476. [Google Scholar] [CrossRef]
  15. Pal, N.; Saxena, N.; Laxmi, K. D.; Mandal, A. Interfacial behaviour, wettability alteration and emulsification characteristics of a novel surfactant: Implications for enhanced oil recovery. Chemical Engineering Science 2018, 187, 200–212. [Google Scholar] [CrossRef]
  16. Haghighi, O. M.; Zargar, G.; Khaksar Manshad, A.; Ali, M.; Takassi, M. A.; Ali, J. A.; Keshavarz, A., Effect of environment-friendly non-ionic surfactant on interfacial tension reduction and wettability alteration; implications for enhanced oil recovery. Energies 2020, 13, (15), 3988. [CrossRef]
  17. Meng, Q.; Liu, H.; Wang, J., A critical review on fundamental mechanisms of spontaneous imbibition and the impact of boundary condition, fluid viscosity and wettability. Advances in Geo-Energy Research 2017, 1, (1), 1-17. [CrossRef]
  18. Meng, M.; Ge, H.; Ji, W.; Shen, Y.; Su, S. Monitor the process of shale spontaneous imbibition in co-current and counter-current displacing gas by using low field nuclear magnetic resonance method. Journal of Natural Gas Science and Engineering 2015, 27, 336–345. [Google Scholar] [CrossRef]
  19. Chakraborty, N.; Karpyn, Z.; Liu, S.; Yoon, H. Permeability evolution of shale during spontaneous imbibition. Journal of Natural Gas Science and Engineering 2017, 38, 590–596. [Google Scholar] [CrossRef]
  20. Zhou, Z.; Abass, H.; Li, X.; Bearinger, D.; Frank, W. Mechanisms of imbibition during hydraulic fracturing in shale formations. Journal of Petroleum Science and Engineering 2016, 141, 125–132. [Google Scholar] [CrossRef]
  21. Jing, W.; Huiqing, L.; Genbao, Q.; Yongcan, P.; Yang, G. Investigations on spontaneous imbibition and the influencing factors in tight oil reservoirs. Fuel 2019, 236, 755–768. [Google Scholar] [CrossRef]
  22. Chang, Q.; Huang, L.; Wu, X., A Molecular Dynamics Study on Low-Pressure Carbon Dioxide in the Water/Oil Interface for Enhanced Oil Recovery. SPE Journal 2023, 28, (02), 643-652. [CrossRef]
  23. Goodarzi, F.; Zendehboudi, S., Effects of salt and surfactant on interfacial characteristics of water/oil systems: molecular dynamic simulations and dissipative particle dynamics. Industrial & Engineering Chemistry Research 2019, 58, (20), 8817-8834. [CrossRef]
  24. Boldon, L.; Laliberte, F.; Liu, L., Review of the fundamental theories behind small angle X-ray scattering, molecular dynamics simulations, and relevant integrated application. Nano reviews 2015, 6, (1), 25661. [CrossRef]
  25. Wang, X.; Xiao, S.; Zhang, Z.; He, J., Effect of nanoparticles on spontaneous imbibition of water into ultraconfined reservoir capillary by molecular dynamics simulation. Energies 2017, 10, (4), 506. [CrossRef]
  26. Yang, S.; Dehghanpour, H.; Binazadeh, M.; Dong, P. A molecular dynamics explanation for fast imbibition of oil in organic tight rocks. Fuel 2017, 190, 409–419. [Google Scholar] [CrossRef]
  27. Karna, N. K.; Oyarzua, E.; Walther, J. H.; Zambrano, H. A., Effect of the meniscus contact angle during early regimes of spontaneous imbibition in nanochannels. Physical Chemistry Chemical Physics 2016, 18, (47), 31997-32001. [CrossRef]
  28. Sang, Q.; Zhao, X.-Y.; Liu, H.-M.; Dong, M.-Z., Analysis of imbibition of n-alkanes in kerogen slits by molecular dynamics simulation for characterization of shale oil rocks. Petroleum Science 2022, 19, (3), 1236-1249. [CrossRef]
  29. Wang, S.; Wang, J.; Liu, H.; Liu, F., Impacts of polar molecules of crude oil on spontaneous imbibition in calcite nanoslit: a molecular dynamics simulation study. Energy & Fuels 2021, 35, (17), 13671-13686. [CrossRef]
  30. Shun, W.; Jing, W.; Huiqing, L., Forced Imbibition to Enhanced Tight/Shale Oil Recovery: A Molecular Dynamics Simulation Study. 2021.
  31. Kang, W.-L.; Zhou, B.-B.; Issakhov, M.; Gabdullin, M. , Advances in enhanced oil recovery technologies for low permeability reservoirs. Petroleum Science 2022, 19(4), 1622–1640. [Google Scholar] [CrossRef]
  32. Liu, B.; Liu, W.; Pan, Z.; Yu, L.; Xie, Z.; Lv, G.; Zhao, P.; Chen, D.; Fang, W., Supercritical CO2 breaking through a water bridge and enhancing shale oil recovery: A molecular dynamics simulation study. Energy & Fuels 2022, 36, (14), 7558-7568. [CrossRef]
  33. Wei, B.; Gao, K.; Song, T.; Zhang, X.; Pu, W.; Xu, X.; Wood, C., Nuclear-magnetic-resonance monitoring of mass exchange in a low-permeability matrix/fracture model during CO2 cyclic injection: a mechanistic study. SPE Journal 2020, 25, (01), 440-450. [CrossRef]
  34. Tu, J.; Sheng, J. J. Experimental and numerical study of surfactant solution spontaneous imbibition in shale oil reservoirs. Journal of the Taiwan Institute of Chemical Engineers 2020, 106, 169–182. [Google Scholar] [CrossRef]
  35. Yuan, L.; Zhang, Y.; Liu, S.; Zhang, Y.; Chen, C.; Song, Y. Study on the slip behavior of CO2-crude oil on nanopore surfaces with different wettability. International Journal of Heat and Mass Transfer 2024, 218, 124787. [Google Scholar] [CrossRef]
  36. Wang, S.; Javadpour, F.; Feng, Q. Molecular dynamics simulations of oil transport through inorganic nanopores in shale. Fuel 2016, 171, 74–86. [Google Scholar] [CrossRef]
  37. Zhang, C.; Dai, H.; Lu, P.; Wu, L.; Zhou, B.; Yu, C., Molecular dynamics simulation of distribution and diffusion behaviour of oil–water interfaces. Molecules 2019, 24, (10), 1905. [CrossRef]
  38. Wang, Z.; Yu, C.; Zhao, J.; Guo, P.; Liu, H. Molecular dynamics simulation for quantitative characterization of wettability transition on silica surface. Journal of Materials Research and Technology 2022, 19, 4371–4380. [Google Scholar] [CrossRef]
  39. Li, G.; Su, Y.; Guo, Y.; Hao, Y.; Li, L. Frontier enhanced oil recovery (EOR) research on the application of imbibition techniques in high-pressure forced soaking of hydraulically fractured shale oil reservoirs. Geofluids 2021, 2021, 1–17. [Google Scholar] [CrossRef]
  40. Dong, H.; Zhu, Q.; Wang, L.; Yue, X.; Fang, H.; Wang, Z.; Liu, S.; Wei, S.; Lu, X., Effects of Shale Pore Size and Connectivity on scCO2 Enhanced Oil Recovery: A Molecular Dynamics Simulation Investigation. Langmuir 2023, 39, (17), 6287-6299. [CrossRef]
  41. Shi, P.; Luo, H.; Tan, X.; Lu, Y.; Zhang, H.; Yang, X., Molecular dynamics simulation study of adsorption of anionic–nonionic surfactants at oil/water interfaces. RSC advances 2022, 12, (42), 27330-27343. [CrossRef]
  42. Xu, G.; Shi, Y.; Jiang, Y.; Jia, C.; Gao, Y.; Han, X.; Zeng, X., Characteristics and influencing factors for forced imbibition in tight sandstone based on low-field nuclear magnetic resonance measurements. Energy & Fuels 2018, 32, (8), 8230-8240. [CrossRef]
  43. Albers, B. Modeling the hysteretic behavior of the capillary pressure in partially saturated porous media: a review. Acta Mechanica 2014, 225, 2163–2189. [Google Scholar] [CrossRef]
  44. Wu, P.; Nikolov, A. D.; Wasan, D. T., Capillary rise: validity of the dynamic contact angle models. Langmuir 2017, 33, (32), 7862-7872. [CrossRef]
  45. Zhang, L.; Ping, J.; Tang, B.; Kang, L.; Imani, G.; Yang, Y.; Sun, H.; Zhong, J.; Yao, J.; Fan, D., Mathematical Model of Two-Phase Spontaneous Imbibition with Dynamic Contact Angle. Transport in Porous Media 2023, 148, (1), 157-172. [CrossRef]
  46. Martic, G.; Gentner, F.; Seveno, D.; Coulon, D.; De Coninck, J.; Blake, T., A molecular dynamics simulation of capillary imbibition. Langmuir 2002, 18, (21), 7971-7976. [CrossRef]
  47. Sun, E. W.-H.; Bourg, I. C., Molecular dynamics simulations of mineral surface wettability by water versus CO2: Thin films, contact angles, and capillary pressure in a silica nanopore. The Journal of Physical Chemistry C 2020, 124, (46), 25382-25395. [CrossRef]
  48. Wu, P.; Nikolov, A. D.; Wasan, D. T. Two-phase displacement dynamics in capillaries-nanofluid reduces the frictional coefficient. Journal of colloid and interface science 2018, 532, 153–160. [Google Scholar] [CrossRef]
  49. Pan, L.; Dai, F.; Li, G.; Liu, S. A TGA/DTA-MS investigation to the influence of process conditions on the pyrolysis of Jimsar oil shale. Energy 2015, 86, 749–757. [Google Scholar] [CrossRef]
  50. Dong, X.; Xu, W.; Liu, H.; Chen, Z.; Lu, N.; Wang, W. On the replacement behavior of CO2 in nanopores of shale oil reservoirs: Insights from wettability tests and molecular dynamics simulations. Geoenergy Science and Engineering 2023, 223, 211528. [Google Scholar] [CrossRef]
  51. Huang, P.; Shen, L.; Gan, Y.; Maggi, F.; El-Zein, A.; Pan, Z., Atomistic study of dynamic contact angles in CO2–water–silica system. Langmuir 2019, 35, (15), 5324-5332. [CrossRef]
  52. Le, T. T. B.; Striolo, A.; Gautam, S. S.; Cole, D. R., Propane–water mixtures confined within cylindrical silica nanopores: structural and dynamical properties probed by molecular dynamics. Langmuir 2017, 33, (42), 11310-11320. [CrossRef]
  53. Xu, J.; Zhan, S.; Wang, W.; Su, Y.; Wang, H. Molecular dynamics simulations of two-phase flow of n-alkanes with water in quartz nanopores. Chemical Engineering Journal 2022, 430, 132800. [Google Scholar] [CrossRef]
  54. Liu, L.; Ye, Z.; Lai, N., Characterization Method for Non-uniform Wettability of Shale Oil Reservoir. International Core Journal of Engineering 2022, 8, (4), 179-186. [CrossRef]
  55. Koretsky, C. M.; Sverjensky, D. A.; Sahai, N., A model of surface site types on oxide and silicate minerals based on crystal chemistry; implications for site types and densities, multi-site adsorption, surface infrared spectroscopy, and dissolution kinetics. American Journal of Science 1998, 298, (5), 349-438. [CrossRef]
  56. Cygan, R. T.; Liang, J.-J.; Kalinichev, A. G., Molecular models of hydroxide, oxyhydroxide, and clay phases and the development of a general force field. The Journal of Physical Chemistry B 2004, 108, (4), 1255-1266. [CrossRef]
  57. Nie, X.; Chen, J.; Sheng, N.; Zeng, L.; Yang, H.; Wang, C., Effect of water molecules on nanoscale wetting behaviour of molecular ethanol on hydroxylated SiO. 2017. [CrossRef]
  58. Hong, X.; Yu, H.; Xu, H.; Wang, X.; Jin, X.; Wu, H.; Wang, F. Competitive adsorption of asphaltene and n-heptane on quartz surfaces and its effect on crude oil transport through nanopores. Journal of Molecular Liquids 2022, 359, 119312. [Google Scholar] [CrossRef]
  59. Zhan, S.; Su, Y.; Jin, Z.; Wang, W.; Li, L. Effect of water film on oil flow in quartz nanopores from molecular perspectives. Fuel 2020, 262, 116560. [Google Scholar] [CrossRef]
  60. Deng, Y.; Wu, Q.; Li, Z.; Huang, X.; Rao, S.; Liang, Y.; Lu, H. , Crystal face dependent wettability of α-quartz: Elucidation by time-of-flight secondary ion mass spectrometry techniques combined with molecular dynamics. Journal of Colloid and Interface Science 2022, 607, 1699–1708. [Google Scholar] [CrossRef]
  61. Mayo, S. L.; Olafson, B. D.; Goddard, W. A., DREIDING: a generic force field for molecular simulations. Journal of Physical chemistry 1990, 94, (26), 8897-8909. [CrossRef]
  62. Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J., Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. Journal of the American Chemical Society 1996, 118, (45), 11225-11236. [CrossRef]
  63. Wang, S.; Feng, Q.; Javadpour, F.; Xia, T.; Li, Z. Oil adsorption in shale nanopores and its effect on recoverable oil-in-place. International Journal of Coal Geology 2015, 147, 9–24. [Google Scholar] [CrossRef]
  64. Berendsen, H.; Grigera, J.; Straatsma, T., The missing term in effective pair potentials. Journal of Physical Chemistry 1987, 91, (24), 6269-6271. [CrossRef]
  65. Zhao, J.; Yao, G.; Ramisetti, S. B.; Hammond, R. B.; Wen, D. Molecular dynamics investigation of substrate wettability alteration and oil transport in a calcite nanopore. Fuel 2019, 239, 1149–1161. [Google Scholar] [CrossRef]
  66. Quezada, G. R.; Jeldres, M.; Toro, N.; Robles, P.; Toledo, P. G.; Jeldres, R. I. Understanding the flocculation mechanism of quartz and kaolinite with polyacrylamide in seawater: A molecular dynamics approach. Colloids and Surfaces A: Physicochemical and Engineering Aspects 2021, 608, 125576. [Google Scholar] [CrossRef]
  67. Shen, X.; Bourg, I. C. Molecular dynamics simulations of the colloidal interaction between smectite clay nanoparticles in liquid water. Journal of Colloid and Interface Science 2021, 584, 610–621. [Google Scholar] [CrossRef]
  68. Meng, J.; Li, C.; Yan, S.; Zhang, S.; Zhang, H.; Wang, G.; Yang, X. Atomic-level insights into the mechanism of saline-regulated montmorillonite (001)-salt droplet interface wetting: A molecular dynamics study. Applied Clay Science 2022, 224, 106513. [Google Scholar] [CrossRef]
  69. Plimpton, S., Fast parallel algorithms for short-range molecular dynamics. Journal of computational physics 1995, 117, (1), 1-19. [CrossRef]
  70. Humphrey, W.; Dalke, A.; Schulten, K., VMD: visual molecular dynamics. Journal of molecular graphics 1996, 14, (1), 33-38. [CrossRef]
  71. Guo, L.; Tang, G.; Kumar, S., Droplet morphology and mobility on lubricant-impregnated surfaces: A molecular dynamics study. Langmuir 2019, 35, (49), 16377-16387. [CrossRef]
Figure 1. (a) Snapshots of imbibition systems I and II, and (b) corresponding dynamic contact angles at several different moments.
Figure 1. (a) Snapshots of imbibition systems I and II, and (b) corresponding dynamic contact angles at several different moments.
Preprints 95747 g001
Figure 2. (a,b) Velocity distribution of shale oil in the strongly hydrophilic system I and weakly hydrophilic system II respectively. (b) Two-dimensional density distribution of shale Oil at different moments during spontaneous imbibition in systems I and II.
Figure 2. (a,b) Velocity distribution of shale oil in the strongly hydrophilic system I and weakly hydrophilic system II respectively. (b) Two-dimensional density distribution of shale Oil at different moments during spontaneous imbibition in systems I and II.
Preprints 95747 g002
Figure 3. (a,d) The mass density distribution of shale oil in the pores of spontaneous imbibition systems I and II after 3ns of structural optimization respectively. The fuchsia area represents the first adsorption layer. (b,e) The number density distribution of different hydrocarbon components in the pores of spontaneous imbibition systems I and II after 3ns of structural optimization respectively. (c,f) The number of changes of different hydrocarbon components in the first adsorption layer with time during spontaneous imbibition in spontaneous systems I and II respectively.
Figure 3. (a,d) The mass density distribution of shale oil in the pores of spontaneous imbibition systems I and II after 3ns of structural optimization respectively. The fuchsia area represents the first adsorption layer. (b,e) The number density distribution of different hydrocarbon components in the pores of spontaneous imbibition systems I and II after 3ns of structural optimization respectively. (c,f) The number of changes of different hydrocarbon components in the first adsorption layer with time during spontaneous imbibition in spontaneous systems I and II respectively.
Preprints 95747 g003
Figure 4. Migration trajectories of C4, C8, and C20 in the first adsorption layer after structural optimization during spontaneous imbibition. The red molecules in the blue circle are labeled as different hydrocarbon molecules in the first adsorption layer after structural optimization.
Figure 4. Migration trajectories of C4, C8, and C20 in the first adsorption layer after structural optimization during spontaneous imbibition. The red molecules in the blue circle are labeled as different hydrocarbon molecules in the first adsorption layer after structural optimization.
Preprints 95747 g004
Figure 5. (a,b) MSD in X- and Z-directions of different hydrocarbon components in the first adsorption layer after structural optimization in imbibition system I during spontaneous imbibition respectively. (c,d) MSD in X- and Z-directions of different hydrocarbon components after structural optimization in the first adsorption layer in spontaneous imbibition system II during spontaneous imbibition respectively.
Figure 5. (a,b) MSD in X- and Z-directions of different hydrocarbon components in the first adsorption layer after structural optimization in imbibition system I during spontaneous imbibition respectively. (c,d) MSD in X- and Z-directions of different hydrocarbon components after structural optimization in the first adsorption layer in spontaneous imbibition system II during spontaneous imbibition respectively.
Preprints 95747 g005
Figure 6. (a) Migration trajectories of “interfacial water molecules” and “adsorbed oil molecules” in spontaneous imbibition I; (b) Migration trajectories of “interfacial water molecules” in spontaneous imbibition II; (c) Schematic of the advancing imbibition front. (d,e) Position coordinates of two “interfacial water molecules” and one “adsorbed oil molecule” in the X- and Z-directions in spontaneous imbibition I at several different moments respectively. (f,g) Position coordinates of one “interfacial water molecule” in the X- and Z-directions in spontaneous imbibition II at several different moments respectively.
Figure 6. (a) Migration trajectories of “interfacial water molecules” and “adsorbed oil molecules” in spontaneous imbibition I; (b) Migration trajectories of “interfacial water molecules” in spontaneous imbibition II; (c) Schematic of the advancing imbibition front. (d,e) Position coordinates of two “interfacial water molecules” and one “adsorbed oil molecule” in the X- and Z-directions in spontaneous imbibition I at several different moments respectively. (f,g) Position coordinates of one “interfacial water molecule” in the X- and Z-directions in spontaneous imbibition II at several different moments respectively.
Preprints 95747 g006
Figure 7. (a,b) Velocity distribution of water molecules in the X-direction in spontaneous imbibition systems I and II at 500ps respectively. The fuchsia red bar represents the oil-water interface region. (c) Schematic diagram of the migration trajectories of water molecules outside the pore. (d) Imbibition displacement efficiency of spontaneous imbibition I.
Figure 7. (a,b) Velocity distribution of water molecules in the X-direction in spontaneous imbibition systems I and II at 500ps respectively. The fuchsia red bar represents the oil-water interface region. (c) Schematic diagram of the migration trajectories of water molecules outside the pore. (d) Imbibition displacement efficiency of spontaneous imbibition I.
Preprints 95747 g007
Figure 8. (a,b) Radial distribution function g(r) of O (hydroxyl) and O (H2O) atoms in the spontaneous imbibition systems I and II respectively. (c,d) g(r) of O (H2O) and C (hydrocarbons) atoms in the spontaneous imbibition systems I and II.
Figure 8. (a,b) Radial distribution function g(r) of O (hydroxyl) and O (H2O) atoms in the spontaneous imbibition systems I and II respectively. (c,d) g(r) of O (H2O) and C (hydrocarbons) atoms in the spontaneous imbibition systems I and II.
Preprints 95747 g008
Figure 9. (a,b) Snapshots of forced imbibition simulation I and II under different pressure differences at several moments. (c) Two-dimensional density distribution of shale oil at 1000ps during forced imbibition in systems I and II.
Figure 9. (a,b) Snapshots of forced imbibition simulation I and II under different pressure differences at several moments. (c) Two-dimensional density distribution of shale oil at 1000ps during forced imbibition in systems I and II.
Preprints 95747 g009
Figure 10. (a,b,c) MSD in the Z-direction for C4, C8, and C20 in the first adsorption layer after structural optimization in forced imbibition I during forced imbibition respectively. (d,e,f) MSD in the Z-direction for C4, C8, and C20 in the first adsorption layer after structural optimization in forced imbibition II during forced imbibition respectively.
Figure 10. (a,b,c) MSD in the Z-direction for C4, C8, and C20 in the first adsorption layer after structural optimization in forced imbibition I during forced imbibition respectively. (d,e,f) MSD in the Z-direction for C4, C8, and C20 in the first adsorption layer after structural optimization in forced imbibition II during forced imbibition respectively.
Preprints 95747 g010
Figure 11. (a) Schematic diagram of slip velocity and peak velocity; (b) Slip and peak velocity at different pressure differences in forced imbibition I; (c) Slip and peak velocity at different pressure differences in forced imbibition II.
Figure 11. (a) Schematic diagram of slip velocity and peak velocity; (b) Slip and peak velocity at different pressure differences in forced imbibition I; (c) Slip and peak velocity at different pressure differences in forced imbibition II.
Preprints 95747 g011
Figure 12. v M D , v t h and ∆v in (a) strongly hydrophilic I and (b) weakly hydrophilic systems II respectively.
Figure 12. v M D , v t h and ∆v in (a) strongly hydrophilic I and (b) weakly hydrophilic systems II respectively.
Preprints 95747 g012
Figure 13. Snapshot of the initial configuration of the imbibition model.
Figure 13. Snapshot of the initial configuration of the imbibition model.
Preprints 95747 g013
Table 1. Imbibition velocity and under different pressure differences in imbibition systems I and II by MD simulations.
Table 1. Imbibition velocity and under different pressure differences in imbibition systems I and II by MD simulations.
P (MPa) 0 5 10 20
i m b i b i t i o n   v e l o c i t y (m/s) system I 5.12 6.39 7.28 9.52
system II 3.29 3.90 4.29 5.15
R system I 1 1.24 1.42 1.85
system II 1 1.185 1.3 1.56
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.
Copyright: This open access article is published under a Creative Commons CC BY 4.0 license, which permit the free download, distribution, and reuse, provided that the author and preprint are cited in any reuse.
Prerpints.org logo

Preprints.org is a free preprint server supported by MDPI in Basel, Switzerland.

Subscribe

© 2024 MDPI (Basel, Switzerland) unless otherwise stated