Abstract
Aiming to study the characteristics and mechanism of oil-water two-phase flow, based on the method of dissipative particle dynamics (DPD), meso-scale models of pressure-driven flow in quartz slits with width of 5-15 nm were simulated. Water is the driving phase. The results show that, as the liquid flew, it stratified gradually. The stratification came earlier as the width, and the driving force increased. The flow in slit of 5 nm was slug flow, while velocity profiles in other models showed parabolic characteristics. When the liquid stratified, the velocity gradient in the near-wall area decreased as a result of the low momentum exchange rate between from water to oil. The density peaks of fluid appeared at the fluid-solid interface and the oil-water interface, caused by the fluid-solid interaction force and the oil-water repulsion, respectively. The oil phase showed stronger aggregation due to the attractive force from the wall. When the water saturation was 25%, the water phase was a water droplet surrounded by oil. When the water saturation was 75%, oil and water stratified. As the repulsive force between the two phases decreased, the surface tension decreased, the mixing degree increased, and the parabolic characteristics of the velocity profile were more obvious.
1. Introduction
At present, acid fracturing is a necessary process of shale oil and gas development. Once fractured, the rock turns to be metrics around with fracture networks, which is filled with fracturing fluid, oil, water, and any other formation fluids there. The size of microfractures and pores in shale ranges from several nanometers to several micrometers. The flow characteristic at this scale is still unclear [1–4]. Ji [5] carried out experiments on the flow characteristics of oil in nanochannel arrays and unsteady water-flooding oil flow experiments in nanochannel arrays. It is reported that the flow of single-phase oil in the nanochannel array conformed to the linear characteristics described by the traditional Hagen-Poiseuille equation, but the experimental volume flow was smaller than the theoretical volume flow. The boundary layer decreased nonlinearly with pressure gradient and eventually tended to constant. Under the same pressure gradient, the resistance coefficient ratio decreased as the pipe diameter increased. In channels of same size, the resistance coefficient ratio decreased with pressure gradient and tended to constant. Wang et al. [6] established a microscale flow model of gas-water two-phase fluid, multiple microscale effects of gas-water two-phase fluids considered, including the following: (1) gas slip, (2) water slip, and (3) viscosity changes of boundary layer. The results showed that microscale effect was not negligible as the size of the throat was less than 3 μm. Based on the distribution of particle size of real core particles, Feng et al. [7] built a three-dimensional pore structure model of the unconsolidated sandstone reservoir and established a two-phase flow model. The results showed that due to the complex pore structure of the porous medium, the fluid exhibited different flow characteristics when flowing through different pores, which in turn affected the overall pressure of the oil-water two-phase flow. According to the bimodal pore size distribution of shale reservoirs, Cui [8] constructed a Gaussian mixture model of shale. The results showed that there was no slip effect on the capillary force curve, while the phase permeability curve was significantly affected.
Experimental research is the most direct research method. In experiment study of two-phase flow, the fluid interface was difficult to detect, and the method of calculating flow rate in single phase flow failed, which calculated the average flow rate by cumulative volume [9]. Therefore, experimental research on two-phase flow was rarely seen at present. In addition, the driving force of discrete fluid in flow simulation was huge when converted into pressure gradient, which made it difficult to compare with the experimental data [10].
Molecular dynamics is the most widely used discrete medium simulation method at present, which is generally used to simulate physical and chemical phenomena in scale of several nanometers to dozens of nanometers. But for the flow simulation, the calculation of molecular simulation is too large, and a lot of the calculations are useless [10]. In DPD, the molecules are coarse grained into beads, and the parameters independent with flow characteristics are ignored, which greatly reduces the computational complexity of flow simulation. In DPD, one of the most critical steps is to calculate the conservative force between the beads. The conservative force between fluid beads is solved by solubility, while the force of solid is set manually according to experience, as it is insoluble [11].
In this paper, based on DPD method, models of oil-water two-phase flow in quartz slits with width of 5-15 nm were established pressure-driven flow in confined space that was simulated, and then, the flow characteristics were analysed, and the effect of fluid-solid interaction on meso-scale was explained.
2. Model Setup
2.1. Geometric Model
The work flow of the simulation is shown in Figure 1. Modeling principles were reported in details, and this paper adopted the same modeling method [11, 12]. In order to reduce the amount of calculation and improve the simulation efficiency, all molecules in DPD models were coarse grained into beads (Table 1 and Figure 1). The flow was controlled by external driving force and the interaction force between the beads of fluids. For quartz molecules at the fluid-solid interface, the excess oxygen atoms combined with the hydrogen atoms in fluid by hydrogen bonds, thus applied additional forces on fluid beads. So, quartz molecules were coarse grained into two different kinds of beads, and those on surface of the slits were S_wall and the others were S.
The geometric model consisted of four parts. The specific size was shown in Figure 2. For solid walls of top and bottom, the outer wall was composed of S beads, and the inner side was S_wall beads. On left side of the slit was water, and the right was oil. Water saturation was 50%. As in Table 1, the coarse graining density is 3.
2.2. Parameters
The force field consisted of 7 forces, including three basic forces in DPD method, which are dissipative force, random force, and conservative force. Dissipative force and random force were taken as 4.5 and 3, respectively, according to DPD method [10]. Then, they were the spring force on the bead bond and the bond angle. The spring coefficient of the bond between H beads was 300 kcal/(mol·Å2), and the stiffness coefficient of the angle between H-H bonds was 300 kcal/(mol·rad2), so as to maintain as rigid structure of n-hexane. Force between fluid and quartz on fluid-solid interface, namely, fluid-solid interaction, was characterized by LJ96X potential function in the function library of material studio. The characteristic of this function was that the effective action range of 1 μm, and it was attractive within the scope of action, which conformed to the simulated physical scenario.
The most critical is the conservative force. According to the principle of DPD, the conservative force between the same kind of beads is 78. The conservative force between different beads was calculated by solubility parameters [10]:
In which, is the conservative force between the same kind of beads, is the conservative force between different kind of beads, is Flory-Huggins parameter, is average molar volume of oil and water, cm3/mol, is gas constant, 8.314 J/(mol·K), and solubility parameter, (J/mol)0.5. In this paper, the parameter of the conservative force between fluids was 150 [11–15]. As comparisons, the models with oil-water conservation force of 78, 100, and 120 were also simulated. The conservative force between fluid and solid is an empirical value. The conservative force between solid wall and water phase is 140, and that of oil phase is 90 [11], as shown in Table 2.
The external driving force, abbreviated as EF, was applied to beads of water phase. The time interval was 15 fs, totally 200, 000 intervals, which were 3000 ps. The results, mainly the coordinates of the beads, saved every 4,000 intervals. So the result of one model contained 51 steps, the first of which was the initial state of the model. Specific parameters of each model were shown in Table 3.
2.3. Model Stability and Validation
Figure 3 is the curve of temperature and system energy curve of model 1. Temperature and energy were important indicators to measure the stability of the system, which remained stable during the simulation. The total simulation time of model 1 reached 10 ns. The model ran stably during the simulation. As the simulation cost too much time, other models ran for 3000 ps. In models 8, 9, 11, 12, and 13, the temperature rose as the displacement phase, water, got enough energy. However, the results could still support the conclusions as all models reached stability during the simulation.
When the fluid was static, oil showed stronger wettability than water, and the contact angle was 55° (Figure 1), which was in accord with physical reality and indicated that under the force field of the model, physical properties of digital fluid were similar to those of the real fluid [11].
3. Flow Characteristics
3.1. Stratification Phenomenon
Oil-water stratification appeared in all models except models 1 and 4. As it flew, oil accumulated toward the solid wall gradually due to its better wettability and delaminated finally (Figure 4). The stratification process differed from each other as flow conditions changed. The increase of w and EF was conducive to occurrence of stratification, so it delaminated earlier (Figure 5). Such increase could be considered as improvements of flow conditions, while it also brought negative effects. For example, the oscillation of the flow was intensified. For models which were not stratified due to oscillation, the degree of oil-water separation in the model with good flow conditions was higher than that with poor flow conditions. The thickness of oil layers on upper and lower wall of the slit was the same in models with poor flow conditions. It became different as the flow conditions became better.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
(l)
(m)
(n)
(o)
(p)
(q)
(r)
(s)
(t)
3.2. Velocity Profile along -Axis
In models 1 and 2, it flew slowly as w and EF were small (Figures 6(a) and 6(b)). As a surface force, the effect of fluid-solid interaction depended on the ratio of area of fluid-solid interface to volume of fluid. As the volume in nanoslits was extremely small, the influence of fluid-solid interaction, which was characterized by LJ96X, on the flow was obvious. Thermal motion played the main role when the fluid flew slowly, as it accounted for a high proportion in total motion. Therefore, the velocity profile oscillated severely.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
In laminar flow of Newtonian fluid, the motion of fluid transmits from center to side of the wall by viscous force. The viscous force is reflected in rate of microscopic momentum transfer caused by the collision of fluid beads. When w was small, the time required for the transmission of the collision was small. It could be considered that it reached the wall instantaneously, thus showing the characteristics of slug flow.
On oil-water interface, velocity of only a small number of beads was counted, so it was insufficient to represent the whole. Velocity profile near oil-water interface oscillated violently. There was slope mutation of velocity profile at the oil-water interface (Figures 6(c)–6(i)). Because water obtained momentum by EF, oil obtained momentum through bead collision at oil-water interface and transferred it to the boundary through viscous force. Slope of velocity profile of oil was obviously lower than that of water. In the center of the slit, as there was hardly any oil, velocity profile of oil there was almost a horizontal line of zero (Figures 6(c)–6(f), 6(h), and 6(i)).
Similarly, at the edge of the slit, the area oil accumulates, and the velocity profile of water was approximately zero. When a small amount of beads was mixed with the other phase, the velocity profile curve of it oscillated severely (Figure 6(g)).
3.3. Velocity Profile along -Axis
As there was no external force along -axis, all momentum of the beads came from random collisions. So, theoretically speaking, the momentum in direction should always fluctuate around zero. At the oil-water interface and the fluid-solid interface, the velocity profile of a single phase fluctuated violently due to the occasional intrusion of beads into different phases (Figure 7). But velocity profile of the total fluid remained smoothly around zero throughout the slit (Figure 7).
(a) Model 5
(b) Model 6
This result supported the analysis of velocity profile along -axis in the previous section. That is, the random collision of beads without external force led to the oscillation of the velocity profile, both at the oil-water interface and the fluid-solid interface.
3.4. Number Density Profile
The bead density profiles of the oil and water showed different characteristics in different regions of the slit (Figure 8). As the fluid delaminated (or the tended to delaminate), the water concentrated mainly in the center of the slit, so the density of water was zero near the wall, the area where the oil gathered. Oil was the opposite. The peak of oil density profile appeared at the fluid-solid interface, as LJ96X provided additional attraction to oil beads there. The viscosity increased accordingly as the density increased; thus, the flow resistance increased, and the velocity decreased.
At the oil-water interface, oil density decreased while water density increased (Figure 8). Different from model 7, no obvious oil-water interface appeared in model 1 or 4, as the density profile was calculated by averaging over time. In slits of 15 nm, there was a clearer plane stage in density profile, and the density decreased rapidly at the oil-water interface, as the flow was stable, and the two phases mixed within a narrow range (Figure 8).
3.5. Distribution Characteristics of Fluid Molecules
The radial distribution function (RDF or g(r)) reflected the aggregation characteristics of discrete media. In pure fluid, spacing of oil/water bead was about 5 Å. So RDF of oil/water showed the first peak at 5 Å (Figure 9), and the peaks appeared repeatedly at multiples of 5 Å, with decreasing peak values. The first peak was the most important feature to measure the characteristics of fluid accumulation. The peak value of oil was higher than water, indicating that the oil phase had better agglomeration, due to LJ96X on oil beads (Figure 9). RDF of water in the middle of the slit was lower than that on the edge (Figure 10(a)). RDF of oil of upper half and lower half was approximately equal (Figure 10(b)). At the oil-water interface, oil and water repelled each other. Beads there were tended of to escape the interface, so the aggregation degree of beads on both sides of the interface increased, and the value of RDF increased.
(a)
(b)
(c)
(d)
(a) Water phase
(b) Oil phase
4. Analysis of Influencing Factors
4.1. Water Saturation
The water saturation in models 10 and 11 was 25% and 75%, respectively. In both models, oil and water flew in layers. When the water saturation is 25%, the stratification occurred at step 11, and the water was entrained as a droplet by oil (Figure 11). When the water saturation was 75%, the water was continuous, and the oil-water interface fluctuated as it flew (Figure 12). In model 10, the flow rate was low, as the saturation of water, the displacement phase, was small. There was no obvious difference between the slope of oil and oil velocity profile (Figure 13(a)). In model 11, the slope of the velocity profile of oil was obviously smaller than that of water, as water obtained external momentum from EF (Figure 13(b)).
(a) Model 10
(b) Model 11
4.2. Temperature
As the total momentum was increased by the action of EF, the temperature rose, too. The actual temperature difference between models 8, 12, and 13 was 40 K. The flow velocity increased with temperature. The most prominent feature was that the velocity gradient of oil at the edge increased (Figure 14). The increase in velocity led to increase in the rate of momentum exchange. Thus, the momentum of water was transferred to oil more easily. Therefore, it was necessary to maintain high formation pressure in reservoir development. In water-flooding development, cold water would be injected cautiously, and its cooling effect should be fully estimated.
(a) Model 12
(b) Model 13
4.3. Conservative Force between Oil and Water
The repulsive force between the oil-water beads was inflected as two-phase surface tension macroscopically. In reservoir development, the impact of surface tension on development efficiency was critical. In general, it was the direction of development and enhanced oil recovery to reduce surface tension and make the displaced phase miscible with the displaced phase. When the conservative force between H and W beads was reduced to 78, miscibility between was achieved (Figure 15(a)). Compared with model 8, the surface tension of models 14-16 was lower. In model 14, oil-water was completely miscible (Figure 15(a)). In model 15, oil-water stratified and finally divided into five-layer flow (oil-water-oil-water-oil) (Figures 15(b)–15(e)). In model 16, after the stratification, the thickness of the upper oil layer was about 2 times that of the lower layer. There was trapped water in oil, which gradually entered water phase as it flew (Figures 15(f)–15(h)). In model 14, as was completely miscible, the momentum of water, obtained from driving force, was fully transferred to oil beads, and the velocity of oil, water, and the whole coincided completely (Figure 16(a)). The velocity profiles of models 15 and 16 were divided into several segments according to the stratification (Figures 16(b) and 16(c)).
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(a) Model 14
(b) Model 15
(c) Model 16
5. Summary and Conclusions
(i)The oil-water two-phase flow in nanoscale slits tended to stratify. Increase of slit width and driving force promoted stratification. As the wettability of oil was stronger, it assembled close to the solid wall, while water flew in the middle of the slit. Since water was the driving phase, the only momentum source of oil beads came from the collision with water beads at the oil-water interface. Thus, velocity of oil was reduced, and the velocity gradient was smaller than water. The velocity profile exhibited different characteristics in the middle and on edge of the slit(ii)Due to LJ96X at fluid-solid interface and repulsive force at oil-water interface, fluid density reached a peak there. As the width of slit increased, the oil-water interface became clearer gradually. When water saturation was 25%, water was completely surrounded by oil, forming water droplets. When it was 75%, the liquid stratified. The increase of temperature also increased the momentum of the system; thus, the velocity increased(iii)When the conservative force between oil and water was reduced to 78, complete miscibility was achieved, and the efficiency of momentum exchange was the highest. The velocity profiles of oil and water completely overlapped, showing parabolic characteristics. With the increase of the conservative force, which was also the repulsive force between the two, stratification gradually became obvious and was accelerated. As a result, the efficiency of momentum exchange between the two decreased, and the velocity gradient was clearly distinguished. The velocity gradient of oil phase was obviously lower than that of water(iv)The influence of slit width and driving force on flow characteristics was great. A smaller width and higher driving force promoted the stratification. When the width was large, the fluid fluctuated violently and was not easy to delaminate
Data Availability
The data is available upon request from the authors.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This research was funded by the Research Startup Fund of Civil Aviation University of China (No. 2020KYQD41).