Digital-Rock Construction of Shale Oil Reservoir and Microscopic Flow Behavior Characterization

Abstract

In shale oil reservoirs, nano-scale pores and micro-scale fractures serve as the primary fluid storage and migration space, while the associated flow mechanism remains vague and is hard to understand. In this research, a three-dimensional (3D) reconstruction of the shale core and micro-pore structure description technique is established; digital core technology for shale reservoirs was developed using X-ray computed tomography (X-CT), scanning electron microscope (SEM) and a focused ion beam scanning electron microscope (FIB-SEM). Microscopic oil–water two-phase flow is mimicked using the lattice Boltzmann method (LBM), a well-acknowledged approach to exploring nanoconfined fluid dynamics. In addition, coupled with digital cores, the flow characteristics of shale reservoirs are characterized. The total porosities of bedding fractures in shale and lamellar shale are 2.042% and 1.085%, respectively. The single-phase oil flow inside bedding fractures follows Darcy’s linear flow principle. This work can deepen the understanding of the microscopic flow characteristics of continental shale reservoirs and provide a reference for similar problems that may be encountered.

1. Introduction

An apparent difference between shale and conventional reservoirs is that the pores in the shale matrix are mainly nano-pores with very low permeability. There are many types of shale pores with complex structures and obvious multi-scale [1]. Because of the micro-nano scale characteristics of shale pores, using traditional methods [2] to study the rock properties generally lacks quantitative description and cannot shed light on the fluid flow behavior and mechanism in the core microstructure. In the early 1990s, Dunsmuir [3] advanced X-ray computed tomography (X-CT) scanning technology and introduced it to the petroleum industry for the first time. Scholars have investigated the flow mechanism in the microstructure porous media in the petroleum engineering industry using digital core technology [4,5,6]. Based on the computed tomography (CT) scanning images of the realistic core sample collected from the field [7], digital core technology and a 3D model of the physical structure are established using reverse reconstruction software. Rosenberg et al. [8] introduced X-CT scanning imaging technology into the field of the digital core and established the 3D digital core of Fontainebleau sandstone. Gao [9] verified the accuracy and reliability of digital cores constructed by CT and suggested that a stricter constraint shall be applied to the alteration principle of average throat radius. Based on CT images, Liu [10] used the pore-scale reconstruction model to simulate water and CO₂ displacement, revealing the influence of pore characteristics on water displacement efficiency and the two-phase flow mechanisms.

Compared to traditional numerical methods, the lattice Boltzmann method (LBM) [11,12] is an efficient method widely used to simulate fluid percolation characteristics in a digital core. The LBM contains different discrete models for mathematical and physical problems one may face in fluid flow modeling [13,14], and its complete parallelism dramatically improves calculation efficiency. In the 1990s, Martys [15] used the LBM model of two-phase immiscible fluid to estimate the relative permeability of 3D digital cores, and its small model size shall take responsibility for the undesirable reliability of calculation results. Walsh [16] established a new partial boundary regression algorithm, which has mass conservation in heterogeneous media and derived the analytical expression for calculating the permeability of different models. Leclaire [17] performed 3D experiments on water absorption and drainage in digital Berea sandstone, and three primary states of fluid invasion are reproduced [18]. Based on the LBM and the multi-component Shan–Chen method, Lautenschlaeger [19] proposed a homogenization method to simulate the fluid–fluid and solid–fluid interactions in the porous media. To solve the problem of inaccurate throat simulation with the help of machine learning, Rabbani [20] combined the advantages of pore network modeling (PNM) and LBM for the first time, and LBM is utilized to calculate throat permeability. Zhao [21] established an improved PNM to simulate single-phase flow in porous media, using LBM to simulate the fluid flow in realistic throats. Suo et al. [22] combined LBM and PNM, PNM modeled the coarse-scale multi-phase transport in porous media; meanwhile, LBM solved the multi-component fluid dynamics.

To sum up, the digital core technology based on CT can genuinely reproduce the internal microstructure inside the core, and its visual accuracy advantage is irreplaceable by other methods. LBM has developed into a widely used and efficient numerical simulation method, which will play an increasingly important role in stimulating the flow characteristics of rocks. In this paper, the digital core constructed by the CT scanning method realizes the realistic microscopic pore structure of the shale core, presented in Section 1. After that, in Section 2, scanning electron microscope (SEM) images are processed to obtain the characteristic parameters of the pore structure in the shale matrix. In Section 3, the PNM of different matrix cores was obtained. Subsequently, LBM is used to quantitatively describe fluid flow characteristics in the microscale pore network in Section 4. At last, conclusions are drawn in Section 5.

2. 3D Reconstruction of Digital Shale Core by CT Scanning

The digital core technology has been widely used to advance oilfield development [23,24], including micro-flow mechanism research, oil displacement mechanism research, and reservoir production performance simulation and prediction. Meanwhile, it is realized that CT scanning is a favorable method for reconstructing the digital core, as it can reproduce the micro-pore structure of the actual core with high accuracy.

2.1. Digital Core Reconstruction Method

The core target selected in this paper is the typical shale reservoir. There are six cores divided into pure shale and Lamellar shale, and the basic parameters of the above six cores are summarized in

Table 1. Basic parameter of collected shale cores.

Shale-Type	Well Name	Core No.	Depth/m	Porosity/%	Permeability/mD
Pure shale	Y-47	CZ-Y47-1	2380.25	7.35	0.012
	G-851	CZ-G851-1	2479.91	6.85	0.009
	G-33	CZ-G33-1	2466	7.66	0.019
Lamellar shale	Y-47	ST-Y47-2	2379.95	7.52	0.011
	G-851	ST-G851-2	2478.21	4.72	0.008
	G-33	ST-G33-4	2398.4	7.29	0.005

. CT detection and digital core reconstruction using the above core mainly include the following six steps: a. Drill 25 mm diameter core on the full-diameter core with wire cutters. b. Use a CT scanner to conduct 3D scanning of a 25 mm diameter core with a scanning accuracy of 13.6 μm. c. Processing the measurement results using Avizo software, identifying identifiable bedding cracks by threshold segmentation, and calculating the width, spacing, and porosity of bedding cracks. d. Drill a core with a diameter of 2 mm at the center of the core with a diameter of 25 mm. e. Use a CT scanner to scan the core with a diameter of 2 mm in three dimensions and a scanning accuracy of 1 μm. f. Processing and calculating the measurement results by using Avizo software.

2.2. Digital Core Reconstruction Method

Figure 1 shows the 3D reconstruction results of 25 mm and 2 mm diameter shale cores; the dark flake and strip areas are shale bedding fractures, and the dark spots are matrix pores; and

Table 2. Characteristic parameters of shale bedding fractures.

Sample No.		CZ-Y47-1	CZ-G851-1	CZ-G33-1	ST-Y47-2	ST-G851-2	ST-G33-4
Core diameter: 25 mm	Average fracture width/μm	39.4	45.5	54.2	54.6	36.5	31.2
	Average fracture spacing/mm	4.5	3.3	7.6	10.71	5.8	15.17
	Porosity/%	1.115	1.756	0.908	0.649	0.802	0.262
Core diameter: 2 mm	Average fracture width/μm	3.9	2.4	4.5	2.8	2.4	4.2
	Average fracture spacing/mm	0.6	0.3	0.8	0.5	0.8	1.2
	Porosity/%	0.828	1.019	0.717	0.713	0.382	0.446
Total porosity/%		1.943	2.776	1.625	1.363	1.184	0.708

shows the statistical table of bedding fracture characteristic parameters of shale. The porosity of fractures in shale cores with a width greater than 13.6μm is about 1.153%, and the average fracture spacing is about 5.15 mm. The porosity of fractures in lamellar shale cores with a width greater than 13.6 μm is about 0.571%.

3. Microporous Structure and Digital Core Construction

3.1. Scanning Images of SEM

Taking core fragments, the shale samples polished by argon ion were scanned by scanning electron microscope in the direction parallel to the bedding plane and perpendicular to the bedding plane of the core, and five viewing zones were taken in each direction for scanning. Figure 2 shows the approximate position of the viewing areas in different directions. The characteristic parameters of matrix pore structure can be obtained through image processing and statistics. In this chapter, six cores in the first Section are scanned by SEM. Figure 3, Figure 4 and Figure 5 shows the scanning results of CZ-Y 47-1. Although the resolution of Figure 3d is 20 μm, its HFW value increases accordingly, so these images can be reasonably used in the process of digital core reconstruction.

3.2. Pore Structure Analysis

After processing the SEM images, the parameters such as specific surface area, pore–throat radius, tortuosity, shape factor and pore-–throat coordination number were obtained. In this Section, CZ-Y47-1, one of the six shale samples is taken as an example, after binarizing the SEM images; as shown in Figure 6, the black pixels are matrix pores, and the matrix surface porosity results in different directions are shown in the title. The plane porosity is the ratio of the pore area to the total image area. Because of the anisotropy of shale, the core porosity is about the average porosity of each plane porosity. Firstly, the processed images are statistically analyzed. According to Equation (1), the plane porosity of three different observation directions can be calculated: 1.8%, 3.6%, and 3.4%, as shown in Figure 6:

$${\varphi }_f=\frac{A_{\mathrm{p}}}{A_{\mathrm{all}}}\hspace{1em}\varphi =({\varphi }_{fv1}+{\varphi }_{fv2}+{\varphi }_{fp})/3$$

(1)

where Φ_f is the plane porosity; A_p is the total pore area, μm²; A_all is the total area of view, μm²; Φ is core porosity; Φ_fv1 is the plane porosity in #1 vertical bedding direction; Φ_fv2 is the plane porosity in #2 vertical bedding direction; Φ_fp is the plane porosity in parallel bedding direction.

Then, the average core porosity of the CZ-Y47-1 core is calculated as 2.9% according to Equation (1). According to the above statistical results and Equation (2), the pore diameter distribution is the ratio of the number of pores with different diameters to the total number of pores; it can be estimated that the average pore diameter of the shale matrix is close to 0.11 μm:

$$D_{\mathrm{p}i}=2\sqrt{A_{\mathrm{p}i}/\pi }\hspace{1em}{\overline{D}}_{\mathrm{p}}={\displaystyle \sum _{i=1}^N{D}_{\mathrm{p}i}}/N$$

(2)

where D_pi is the diameter of a single pore, μm; A_pi is the area of a single pore, μm²; ${\overline{D}}_{\mathrm{p}}$ is the average pore diameter, μm; N is the total number of pores.

The throat diameter distribution refers to the ratio of the number of throats with different diameters to the total number:

$$D_{\mathrm{t}i}=2\sqrt{A_{\mathrm{t}i}/\pi }\hspace{1em}{\overline{D}}_{\mathrm{t}}={\displaystyle \sum _{i=1}^M{D}_{\mathrm{t}i}}/M$$

(3)

where D_ti is the diameter of a single throat, μm; A_ti is the area of a single throat, μm²; ${\overline{D}}_{\mathrm{t}}$ is the average throat diameter, μm; M is the total throat number. The tortuosity is the ratio of the actual length of the flow channel to the apparent length through the flow medium, which indicates the complexity of the fluid flow channel:

$${\tau }_{\mathrm{t}}{}_i=\frac{S_i}{L_i}\hspace{1em}{\overline{\tau }}_{\mathrm{t}}={\displaystyle \sum _{i=1}^M{\tau }_{\mathrm{t}i}}/M$$

(4)

where τ_ti is the tortuosity of a single throat; S_i is the actual length of a single throat between two pores, μm; L_i is the linear length of a single throat between two pores, μm; ${\overline{\tau }}_{\mathrm{t}}$ is the average tortuosity. The coordination number refers to the number of connected throats of a single pore, indicating the pore network’s complexity:

$${\overline{CN}}_{\mathrm{t}}={\displaystyle \sum _{i=1}^{N}C{N}_{\mathrm{p}i}}/N$$

(5)

where CN_pi is the coordination number of a single pore. The shape factor is the ratio of the cross-sectional area of the pore and throats to the square of the perimeter, which can quantitatively show the shape characteristics of the pore-throat unit:

$$F_{\mathrm{p}i}=\frac{A_{\mathrm{p}i}}{S_{perim,i}^2}\hspace{1em}{\overline{F}}_{\mathrm{p}}={\displaystyle \sum _{i=1}^N{F}_{\mathrm{p}i}}/N$$

(6)

where ${\overline{F}}_{\mathrm{p}}$ is the average shape factor.

By calculation according to Equations (3)–(6), the average throat diameter of the shale matrix approaches 0.034 μm; the above results are shown in Figure 7 and Figure 8. The average tortuosity of the shale matrix is 6.53, and the average coordination number of the shale matrix is 1.77, and using the above formulas and image statistics, the parameters are summarized in

Table 3. Summary of pore and throat structure parameters.

Shale Type	No.	Sample No.	Porosity/%	Average Pore Diameter/nm	Average Throat Diameter/nm	Tortuosity	Coordination Number	Shape Factor
Pure Shale	1	CZ-Y47-1	2.9	113	34	6.53	1.77	0.037
	2	CZ-G851-1	2.8	128	32	4.78	1.85	0.034
	3	CZ-G33-1	3.3	99	32	5.65	1.89	0.033
	average value		3	113	32.6	5.65	1.837	0.035
Lamellar shale	4	ST-Y47-2	2.6	150	35	5.1	1.89	0.049
	5	ST-G851-2	1.4	163	34	4.98	1.93	0.054
	6	ST-G33-4	2.7	149	31	5.11	1.88	0.058
	average value		2.23	154	33.3	5.06	1.9	0.054

.

3.3. Digital Core Construction Based on SEM Images

Figure 9 and Figure 10 present digital cores of shale matrix with average pore–throat diameters of 130 and 30 nm, respectively. They are reconstructed by the random growth method, which provides a model for subsequent micro-flow simulation with less than 5% error.

4. Construction of Pore Network Model Based on FIB-SEM Tests

PNM, pore network modeling, can reflect the topological structure of the core pore space and has geometric parameters equivalent to the real core pore space [25,26]. Shale cores were scanned by focused ion beam scanning electron microscope (FIB-SEM) technology, and the slice spacing was close to 10 nm. Then, the PNM of different matrix cores is obtained through image processing and reconstruction.

4.1. Construction of PNM

The pixel size is 5.42 nm, and the effective measurement range is 9100 × 6476 × 3990 nm. The digital core of CZ-Y47-1 was obtained by cutting the image using Avizo software. Image segmentation of organic matter, pore, and matrix is performed to extract the structure of different components, as shown in Figure 11.

Reconstructing the binarized images to obtain the segmented 3D photos, as shown in Figure 12, and pore network modeling is generated, as presented in Figure 13.

4.2. Structural Characteristics Analysis of PNM

The classification of reservoir space type depends on the structural characteristics of the core pore network, which can be displayed intuitively and clearly by extracting the digital core pore network model. The microstructure of pores can be evaluated in detail by statistical analysis of structural parameters [27].

According to the 3D digital core of CZ-Y47-1, volume percentage statistics of the pore, organic matter, and core matrix are 1.85%, 4.12%, and 95.03%, respectively. The fractal dimension of the generated PNM is 1.95, and the tortuosity is 8.72. The statistical distribution of core structure characteristics is shown in Figure 14.

The PNM parameters are shown in

Table 4. Pore network parameters (unit, volume: ml; radius: nm).

Throat Number	88	Total Number of Pores	3866
Maximum throat radius	220.26	Maximum pore diameter	848.43
Average throat radius	44.58	Average pore diameter	60.07
Maximum throat volume	8.37 × 107	Maximum pore volume	3.20 × 108
Average throat volume	1.60 × 107	Average pore volume	1.13 × 106
Maximum throat length	1160.05	Maximum coordination number	4
Average throat length	454.48	Average coordination number	1.41

, and it can be observed that the porosity of micropores is relatively large; the increase of micropores enhances the connectivity between some isolated macropores.

5. Microscopic Flow Characteristics

Based on the reconstructed digital core, LBM can be used to simulate the displacement process and quantitatively characterize the flow characteristics of the fluid in the microporous network [28]. The underlying influence of wettability, interfacial tension, and fluid viscosity on the single-phase oil flow characteristics are investigated, expecting to guide shale oil production prediction and development scheme optimization.

5.1. Establishment of the Flow Model

In 1992, Qian [29] developed the D_dQ_m model that laid a foundation for the extensive application of the LBM. The LBM can also be deduced from the continuous Boltzmann Equation:

$$\frac{\partial f}{\partial t}+e·{\nabla }_{x}f+\mathit{a}·{\nabla }_{\xi }f=\mathsf{\Omega }(f)$$

(7)

where $f$ is the particle number density distribution function, e is the particle velocity, a is the particle motion acceleration, and $\mathsf{\Omega }(f)$ is the collision operator.

The collision operator is a complex nonlinear integral form. A linear collision operator is proposed based on Maxwell equilibrium distribution, and the Boltzmann–BGK Equation is obtained:

$$\frac{\partial f}{\partial t}+e·{\nabla }_{x}f+\mathit{a}·{\nabla }_{\xi }f=v(f^{eq}-f)$$

(8)

BGK is a continuous equation, and the discretization of velocity and space obtains the discretized lattice Boltzmann Equation. From the macroscopic statistical theory, it can be considered that all molecules are moving in specified directions, so the rate of particle movement is simplified to several rates along the specified direction {e₀, e₁,…, e_N}, and N is the direction number of the rate. The distribution function f of particle number density is also discretized into {f₀, f₁, …, f_N}, which defines f_α = f_α (r, e, α, t), α = 0, 1, …, N. The Boltzmann Equation for velocity dispersion can be obtained as follows:

$$\frac{\partial f_{\alpha }}{\partial t}+e_{\alpha }·\nabla f_{\alpha }=\frac{1}{{\tau }_0}(f_{\alpha }^{eq}-f_{\alpha })+{(\mathit{a}·{\nabla }_{\xi }f)}_{\alpha }$$

(9)

The equation is discretized in time and space, and the fully discretized LB Equation is obtained:

$$f_{\alpha }(\mathit{r}+{\mathit{e}}_{\alpha }{\delta }_t,t+{\delta }_t)-f_{\alpha }({\mathit{r}}_t,t)=-\frac{1}{\tau }\left[f_{\alpha }({\mathit{r}}_t,t)-f_{\alpha }^{eq}({\mathit{r}}_t,t)\right]+{\delta }_t{F}_{\alpha }({\mathit{r}}_t,t)$$

(10)

The last term on the right side of the equal sign is the external force term, and the above formula is the Lattice–BGK Equation considering external force. The LBM comprises lattice type, evolution equation of particle number density distribution, equilibrium distribution function, and boundary condition. Using the D₃Q₁₉ lattice, the evolution equation is:

$$f_{\alpha }\left(\mathit{r}+{\mathit{e}}_{\alpha }{\delta }_t,t+{\delta }_t\right)=f_{\alpha }\left(r,t\right)-\frac{1}{\tau }\left(f_{\alpha }\left(r,t\right)-f_{\alpha }{}^{eq}\left(r,t\right)\right)$$

(11)

The equilibrium distribution function is:

$$f_{\alpha }{}^{eq}={\omega }_{\alpha }\rho \left[1+\frac{{\mathit{e}}_{\alpha }\cdot \mathit{u}}{c_s{}^2}+\frac{{\left({\mathit{e}}_{\alpha }\cdot \mathit{u}\right)}^2}{2c_s{}^4}-\frac{u^2}{2c_s{}^2}\right]$$

(12)

where, $c_s$ is the lattice sound speed, $c_s=\frac{c}{\sqrt{3}}$, and ${\omega }_{\alpha }$ is the weight coefficient:

$${\omega }_{\alpha }=\left\{\begin{array}{cc}4/9,& c_{\alpha }{}^2=0\\ 1/9,& c_{\alpha }{}^2=c^2\\ 1/36,& c_{\alpha }{}^2=2c^2\end{array}\right.$$

(13)

Consider the boundary conditions of flow:

$$\left\{\begin{array}{l}\rho ={\displaystyle \int fde={\displaystyle \sum f_{\alpha }}} \alpha =1,2,\cdots ,19\\ u=\frac{{\displaystyle \int fede}}{\rho }=\frac{{\displaystyle \sum f_{\alpha }{e}_{\alpha }} }{\rho } \alpha =1,2,\cdots ,19\\ P=c_s^2\rho \end{array}\right.$$

(14)

Based on the above equations, a single-phase flow LBM model is established considering the interfacial tension and liquid–solid wettability. The viscosity of the oil phase can influence the relaxation time:

$$\tau =\frac{\upsilon M_w}{RT{\delta }_t}+0.5$$

(15)

The oil phase adsorption layer on the pore wall affects the oil phase flow in the pore. The oil–solid interfacial tension (γ) and contact angle (θ) jointly determine the adhesive force. The contact angle is only used to describe the strength of fluid–solid force and does not indicate the existence of multiphase fluid in porous media; additionally, with the increase of the contact angle, the wettability of the wall gradually tends to water.

$$F_s=r{\gamma }_{so}(1+\cos \theta )$$

(16)

$$F_{\alpha }=\left(1-\frac{1}{2\tau }\right)\frac{\left({\mathit{c}}_{\alpha }-\mathit{u}\right)·{\mathit{F}}_{\mathrm{s}}}{\rho c_s^2}{f}_{\alpha }^{eq}$$

(17)

5.2. Percolation Characteristics of Crude Oil in Matrix Pores

Based on the established shale matrix digital cores with average pore–throat diameters of 130 nm and 30 nm, the permeability of the digital cores is calculated, and the matrix permeability of the shale core is measured by a matrix permeability meter. The simulation parameters are shown in

Table 5. Flow simulation parameters.

Physical Parameters	Value
Pressure gradient (MPa/m)	1, 5, 10, 20, 50, 100, 150, 200, 250, 300
Oil viscosity (mPa·s)	0.7, 0.9, 1.1, 1.3
Oil density (g/cm3)	0.83
Contact angle (degree)	10, 30, 50, 70
Oil solid interfacial tension (mN/m)	0.5, 0.7, 0.9, 1.1
The proportion of water-wetted minerals (%)	0, 30, 70, 100

. Compared with the two models, the error between the model verification measurement results and the simulation results is only 6.2%, clarifying the accuracy and reliability of the proposed model. In this paper, the length scale $L_0$ equals 5 × 10⁻⁸ m, the mass scale M₀ equals 1.51 × 10⁻¹⁹ kg, and the time scale $T_0$ equals 1.51 × 10⁻¹⁰ s. Then, the physical pore length ($L_{physics}$) can be calculated by $L_{physics}=L_0\times L_{lattice}$, where $L_{lattice}$ is the lattice length in LB simulation. The physical viscosity ($v_{physics}$) can be calculated by $v_{physics}=v_{lattice}\frac{L_0^2}{T_0}$, where $v_{lattice}$ is lattice viscosity. The physical density (${\rho }_{physics}$) can be calculated by ${\rho }_{physics}={\rho }_{lattice}\frac{M_0}{L_0^3}$, where ${\rho }_{lattice}$ is the lattice density.

With increased hydrophilic mineral proportion, the starting pressure gradient of single-phase oil flow in the shale matrix decreases, as seen in Figure 15.

As illustrated in Figure 16, with the increase in the proportion of hydrophilic minerals, the adequate area oil flow increases, and the pressure decreases in large-size pores. The adsorption force the water-wet surface imposes on oil molecules is less than that with the oil-wet surface; as a result, oil-flowing resistance in water-wet pores is reduced compared to oil-wet pores, leading to relatively high flow efficiency.

It can be seen from Figure 17 that the starting pressure gradient of single-phase oil flow in the shale matrix decreases with the increase of wetting angle.

As can be seen from Figure 18, with the increase of wetting angle, the area where oil can flow increases, and the pressure in large-size pores decreases.

Figure 19 shows that the starting pressure gradient of single-phase oil flow in the shale matrix increases with the increase of interfacial tension.

It can be seen from Figure 20 that with the increase of interfacial tension, the oil-flow area decreases, and the pressure in large-size pores increases. The adsorption force of oil on the pore wall increases, and the flow characteristics become more complicated.

Figure 21 shows that as the viscosity of the oil phase increases, the starting pressure gradient of single-phase oil flow in the shale matrix increases.

It can be seen from Figure 22 that when the viscosity of the oil phase increases, the oil-flow area decreases, and the pressure in large-size pores decreases. The increasing oil viscosity makes it more difficult to flow in the pores.

According to the single-phase starting pressure gradient curve, the average pore-throat diameters are 30 nm and 130 nm, and the starting pressure gradient of flow in the matrix is 0.103–0.285 and 0.045–0.110 MPa/mm, respectively.

5.3. Flow Characteristics of the Oil Phase in Bedding Fractures

Based on this digital core, the percolation characteristics of single-phase oil in the bedding fractures are simulated. As seen from Figure 23, there is no starting pressure gradient in the bedding fracture flow, while the effect of stress is not considered, and the velocity has a linear relationship with the pressure gradient.

Figure 24 shows the relationship between the single-phase flow velocity of crude oil and pressure gradient under different influential factors. The results show that the single-phase flow follows Darcy’s principle, and there is no threshold pressure gradient. The more significant the proportion of hydrophilic minerals is, the greater the wetting angle of the oil phase is, the smaller the interfacial tension is, and the higher the flow velocity of the oil phase is. In addition, the higher the viscosity of the oil phase, the lower the flow rate of the oil phase.

6. Conclusions

In this paper, the CT scans of cores with diameters of 25 mm and 2 mm were carried out, six 3D digital cores were reconstructed, and the characteristics of shales and lamellar shales were obtained. SEM and FIB-SEM were used to rebuild the digital core, obtain the characteristic parameters of the pore structure of the shale matrix, and reconstruct the 3D pore network model of the shale matrix. Considering the viscosity, interfacial tension, and wettability of fluid, the LBM model of shale matrix/bedding fracture flow is established, and the microscopic flow characteristics of shale oil under different conditions are calculated. The following conclusions were drawn:

(1)

Total porosities of bedding fractures in shale and lamellar shale are 2.042% and 1.085%, respectively; the average fracture spacing is 0.5 and 0.9 mm, and the average fracture widths are 6.8 and 5.5 μm.

(2)

The average porosities of the shale and lamellar shale matrix cores obtained by SEM are 3% and 2.23%, the average pore diameters are 113 and 154 nm, the average throat diameters are 32.6 and 33.3 nm, the average tortuosities are 5.65 and 5.06, and the average shape factors are 0.035 and 0.054, respectively.

(3)

The average pore throat diameter is 30 nm, and the starting pressure gradients in the matrix are 0.103~0.285 and 0.045~0.110 MPa/mm, respectively. The mechanism of single-phase flow in shale bedding fractures is the same as Darcy’s linear flow.

(4)

The accuracy of the digital core model constructed by the three modeling methods needs further improvement. Although it can reflect the structural characteristics of the reservoir, the single structure method and mixed numerical reconstruction method still need to be improved to simulate the pore structure and improve the application results accurately.