Flow Rate Ratio Optimization for Waterflooding in Irregular Well Patterns: Insights from Numerical Simulation and 3D Visualization

1. Introduction

As high-quality and easily exploitable oil resources decline, tight reservoirs have become an important alternative to ensure a continued supply of crude oil [1,2]. Horizontal wells with volume fracturing stimulation are effective technologies for developing tight reservoirs [3]. The use of horizontal wells significantly increases oil production that cannot be achieved through a single vertical well. These wells are intersected by artificial fractures that provide a more efficient flow path for crude oil. However, as development progresses, the oil production of horizontal wells consistently declines with the depletion of formation pressure [4]. To replenish formation energy, petroleum engineers have proposed a series of technologies, including waterflooding, gas flooding, water-alternating-gas (WAG) flooding, water injection huff-n-puff, and gas huff-n-puff [5]. Numerous laboratory experiments and numerical simulations have provided substantial evidence that gas injection, particularly CO₂ injection, is a promising enhanced oil recovery (EOR) method due to its favorable attributes, such as high injection capacity and miscibility with oil [6,7]. However, due to the high cost of CO₂ sources, water injection remains the most commonly used method.

Water flooding is an effective and cost-efficient technique to enhance oil recovery following primary production in conventional reservoirs. However, significant reservoir heterogeneity caused by natural and artificial fractures often results in low sweep efficiency of injected water in tight reservoirs [8]. To recover more oil through water flooding, researchers primarily focus on optimizing well patterns [9,10]. The diamond-shaped inverted nine-spot well pattern has been shown to be superior to other well patterns in developing fractured tight reservoirs. However, well patterns in oil fields are often irregular [11,12,13]. A typical irregular horizontal and vertical combined well group from the Fuxian oilfield (Ordos Basin, China) is shown in Figure 1. The newly drilled horizontal well serves as both a production and an exploratory well for understanding the area's geology. Consequently, the spacing between newly drilled horizontal wells is relatively large. Meanwhile, there are typically only 1 to 3 water injection wells around a horizontal well due to cost constraints. Additionally, it is impractical to modify the well pattern by drilling new wells in the current era of low oil prices. Therefore, determining suitable injection parameters to improve oil recovery under the existing irregular well patterns is crucial.

According to the available literature in this field, flow rates of different vertical injectors are typically kept the same within a well group. However, flow rates of different vertical injectors should not be uniform in irregular well groups due to the varying relative positions of vertical injectors and horizontal producers [14,15]. To the best of our knowledge, limited literature has focused on the optimization of flow rates in irregular well groups. Therefore, this study employs streamline simulation and 3D visualization technology to optimize the injection parameters of irregular well groups.

The findings of this research have significant practical implications for the development of tight reservoirs, particularly in mature fields with irregular well patterns. By optimizing injection rates based on spatial relationships between wells, this approach can enhance sweep efficiency without the need for costly infrastructure changes, making it economically viable for operators in low-price oil environments. Furthermore, the methodology can be extended to other unconventional reservoirs, such as shale or low-permeability formations, where uneven fluid distribution is a persistent challenge. From a broader perspective, this work contributes to sustainable hydrocarbon recovery by maximizing resource utilization and delaying well abandonment, thereby reducing the environmental footprint per unit of oil produced.

In this paper, twenty-two cases with varying relative positions of vertical injectors and horizontal producers were initially studied. Based on the simulation results, the controlling area of the horizontal well was divided into three regions where vertical injectors should perform different flow rates. Subsequently, a typical irregular well group with one horizontal well and three vertical wells was studied to determine the flow rate ratio between the injectors from different regions. Finally, analyses of the streamline simulation were conducted to illustrate the mechanism behind the difference in sweep efficiency caused by varying flow rate ratios between the injectors.

2. Reservoir Modeling

2.1. Modeling of a Combined Well Network of One-Injection-one-Production Water Straight

To optimize the injection parameters for the horizontal-straight well group, a conceptual model consisting of one vertical injection well and one horizontal production well was developed. The governing mathematical model for seepage mechanics is presented as follows:

(1)

Seepage model of fluid in Region I:

Fluid flow in Region I comprises three components: (1) Darcy flow through the artificial main fracture along the x-direction, (2) Darcy flow from the fracture in Region II into the artificial main fracture along the y-direction, and (3) fluid transfer from the matrix in Region II into the artificial main fracture along the y-direction. The oil and water phase continuity equations governing fluid flow from the artificial main fracture to the wellbore are expressed as Equation (1-2):

$$-{\displaystyle \frac{𝜕}{𝜕x}}\left({\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\cdot {\displaystyle \frac{𝜕p_{\mathrm{wf}}}{𝜕x}}\right)+{\displaystyle \frac{𝜕\left({\rho }_{\mathrm{wm}}\cdot {\tau }_{\mathrm{mfw}}\right)}{𝜕y}}+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}{\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕p_{\mathrm{wf}}^{2\to 1}}{𝜕y}}\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wf}}{\phi }_{\mathrm{f}}{S}_{\mathrm{wf}}\right)$$

(1)

$$-{\displaystyle \frac{𝜕}{𝜕x}}\left({\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\cdot {\displaystyle \frac{𝜕p_{\mathrm{of}}}{𝜕x}}\right)+{\displaystyle \frac{𝜕\left({\rho }_{\mathrm{om}}\cdot {\tau }_{\mathrm{mfo}}\right)}{𝜕y}}++{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}{\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕p_{\mathrm{of}}^{2\to 1}}{𝜕y}}\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{of}}{\phi }_{\mathrm{f}}{S}_{\mathrm{of}}\right)$$

(2)

In Equation [1,2], variables $p_{\mathrm{of}}^{2\to 1}$ and $p_{\mathrm{wf}}^{2\to 1}$ denote the oil-phase and water-phase pressures (MPa) in the Region II fracture grid at the Region II-I interface, respectively.

The governing equations for matrix-to-fracture fluid transfer toward the artificial main fractures are:

$$\begin{array}{l}{\tau }_{\mathrm{mfw}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)}{𝜕y}}\\ {\tau }_{\mathrm{mfo}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)}{𝜕y}}\end{array}$$

(3)

Substituting Equation (3) into Equations (1) and (2) yields the continuity equation for fluid flow in Region I:

$$\begin{array}{l}-{\displaystyle \frac{𝜕}{𝜕x}}\left({\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\cdot {\displaystyle \frac{𝜕p_{\mathrm{wf}}}{𝜕x}}\right)+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{wm}}\cdot \sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)}{𝜕y}}\right]\\ +{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}{\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕p_{\mathrm{wf}}^{2\to 1}}{𝜕y}}\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wf}}{\phi }_{\mathrm{f}}{S}_{\mathrm{wf}}\right)\end{array}$$

(4)

$$\begin{array}{l}-{\displaystyle \frac{𝜕}{𝜕x}}\left({\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{of}}}}\cdot {\displaystyle \frac{𝜕p_{\mathrm{of}}}{𝜕x}}\right)+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{om}}\cdot \sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)}{𝜕y}}\right]\\ +{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}{\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_o}}{\displaystyle \frac{𝜕p_{\mathrm{of}}^{2\to 1}}{𝜕y}}\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{of}}{\phi }_{\mathrm{f}}{S}_{\mathrm{of}}\right)\end{array}$$

(5)

(2)

Seepage model of fluid in Region II:

Region II represents a dual-porosity system characterized by multidirectional fluid flow, encompassing matrix flow in both x- and y-directions within the region, along with fluid transfer from the matrix to fractures. At the Region II-I boundary, y-directional flow occurs from both matrix and fractures in Region II to Region I, while the Region II-III interface exhibits x-directional flow from Region III to Region II through matrix and fracture networks. Similarly, y-directional transfer takes place from Region IV to Region II at their shared boundary. The governing differential equations describing this coupled matrix-fracture flow system, including the oil-water phase behavior in fractures, are expressed as Equation (6-7):

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\nabla p_{\mathrm{wf}}\right]+\nabla \left({\rho }_{\mathrm{wm}}\cdot {\tau }_{\mathrm{mfw}}\right)-{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕p_{\mathrm{wf}}^{2\to 1}}{𝜕y}}\right]\\ +{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wm}}^{3\to 2}}{𝜕x}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wf}}^{3\to 2}}{𝜕x}}\right)\right]\\ +{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wm}}^{4\to 2}}{𝜕y}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wf}}^{4\to 2}}{𝜕y}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wf}}{\phi }_{\mathrm{f}}{S}_{\mathrm{wf}}\right)\end{array}$$

(6)

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{of}}}}\nabla p_{\mathrm{of}}\right]+\nabla \left({\rho }_{\mathrm{om}}\cdot {\tau }_{\mathrm{mfo}}\right)-{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}{\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕p_{\mathrm{of}}^{2\to 1}}{𝜕x}}\right]\\ +{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{3\to 2}}{𝜕x}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{3\to 2}}{𝜕x}}\right)\right]\\ +{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{4\to 2}}{𝜕y}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{4\to 2}}{𝜕y}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{of}}{\phi }_{\mathrm{f}}{S}_{\mathrm{of}}\right)\end{array}$$

(7)

In the Equation, the phase pressures are defined as follows:

$p_{\mathrm{of}}^{3\to 2}$ and $p_{\mathrm{wf}}^{3\to 2}$: oil-phase and water-phase pressures (MPa) in the Region III fracture grid at the Region III-II interface

$p_{\mathrm{om}}^{3\to 2}$ and $p_{\mathrm{wm}}^{3\to 2}$: corresponding phase pressures (MPa) in the Region III matrix grid

$p_{\mathrm{of}}^{4\to 2}$ and $p_{\mathrm{wf}}^{4\to 2}$: oil-water phase pressures (MPa) in the Region IV fracture grid at the Region IV-II boundary

$p_{\mathrm{om}}^{4\to 2}$ and $p_{\mathrm{wm}}^{4\to 2}$: phase pressures (MPa) in the Region IV matrix grid at the same boundary

The governing differential equations for oil-water two-phase flow in the matrix system are:

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\nabla p_{\mathrm{wm}}\right]-\nabla \left({\rho }_{\mathrm{wm}}\cdot {\tau }_{\mathrm{mfw}}\right)={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wm}}{\phi }_{\mathrm{m}}{S}_{\mathrm{wm}}\right)\\ -\nabla \left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\nabla p_{\mathrm{om}}\right]-\nabla \left({\rho }_{\mathrm{om}}\cdot {\tau }_{\mathrm{mfo}}\right)={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{om}}{\phi }_{\mathrm{m}}{S}_{\mathrm{om}}\right)\end{array}$$

(8)

The mass transfer equation governing matrix-fracture crossflow is expressed as:

$$\begin{array}{l}{\tau }_{\mathrm{mfw}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\nabla \left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)\\ {\tau }_{\mathrm{mfo}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\nabla \left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)\end{array}$$

(9)

By substituting Equation (9) into the governing equations (6) -(8), we derive the coupled system of oil-water phase flow equations for the fracture network:

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\nabla p_{\mathrm{wf}}\right]+\nabla \left[{\rho }_{\mathrm{wm}}\cdot \sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\nabla \left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)\right]-{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕p_{\mathrm{wf}}^{2\to 1}}{𝜕y}}\right]\\ +{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wm}}^{3\to 2}}{𝜕x}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wf}}^{3\to 2}}{𝜕x}}\right)\right]\\ +{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wm}}^{4\to 2}}{𝜕y}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wf}}^{4\to 2}}{𝜕y}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wf}}{\phi }_{\mathrm{f}}{S}_{\mathrm{wf}}\right)\end{array}$$

(10)

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{of}}}}\nabla p_{\mathrm{of}}\right]+\nabla \left[{\rho }_{\mathrm{om}}\cdot \sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\nabla \left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)\right]-{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}{\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕p_{\mathrm{of}}^{2\to 1}}{𝜕x}}\right]\\ +{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{3\to 2}}{𝜕x}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{3\to 2}}{𝜕x}}\right)\right]\\ +{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{4\to 2}}{𝜕y}}-G\right)\right]+{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{4\to 2}}{𝜕y}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{of}}{\phi }_{\mathrm{f}}{S}_{\mathrm{of}}\right)\end{array}$$

(11)

Seepage equations for the oil and water phases in the matrix:

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\nabla p_{\mathrm{wm}}\right]-\nabla \left[{\rho }_{\mathrm{wm}}\cdot \sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\nabla \left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wm}}{\phi }_{\mathrm{m}}{S}_{\mathrm{wm}}\right)\\ -\nabla \left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\nabla p_{\mathrm{om}}\right]-\nabla \left[{\rho }_{\mathrm{om}}\cdot \sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\nabla \left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{om}}{\phi }_{\mathrm{m}}{S}_{\mathrm{om}}\right)\end{array}$$

(12)

(3)

Seepage model of fluid in Region III:

Differential equations of seepage for the oil phase and water phase in the matrix:

The governing differential equations for oil-water two-phase flow in the matrix are given by Equations (13) and (14):

$$\begin{array}{l}-{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕p_{\mathrm{m}}}{𝜕x}}\right]-{\displaystyle \frac{𝜕\left({\rho }_{\mathrm{wm}}\cdot {\tau }_{\mathrm{mfw}}\right)}{𝜕x}}\\ -{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wm}}^{3\to 2}}{𝜕x}}-G\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wm}}{\phi }_{\mathrm{m}}{S}_{\mathrm{wm}}\right)\end{array}$$

(13)

$$\begin{array}{l}-{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕p_{\mathrm{m}}}{𝜕x}}\right]-{\displaystyle \frac{𝜕\left({\rho }_{\mathrm{om}}\cdot {\tau }_{\mathrm{mfo}}\right)}{𝜕x}}\\ -{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{3\to 2}}{𝜕x}}-G\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{om}}{\phi }_{\mathrm{m}}{S}_{\mathrm{om}}\right)\end{array}$$

(14)

The governing differential equations for oil-water two-phase flow in the fracture network are given by Equations (15) and (16):

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\nabla p_{\mathrm{nf}}\right]+\nabla \left({\rho }_{\mathrm{wm}}\cdot {\tau }_{\mathrm{mfw}}\right)\\ -{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{wf}}^{3\to 2}}{𝜕x}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wf}}{\phi }_{\mathrm{nf}}{S}_{\mathrm{wf}}\right)\end{array}$$

(15)

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\nabla p_{\mathrm{f}}\right]+\nabla \left({\rho }_{\mathrm{om}}\cdot {\tau }_{\mathrm{mfo}}\right)\\ -{\displaystyle \frac{𝜕}{𝜕x}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{3\to 2}}{𝜕x}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{of}}{\phi }_{\mathrm{f}}{S}_{\mathrm{of}}\right)\end{array}$$

(16)

Among these equations, the equation governing the cross-flow phenomenon from the matrix to the fractures is as follows:

$$\begin{array}{l}{\tau }_{\mathrm{mfw}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{wm}}{B}_{\mathrm{wm}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)}{𝜕x}}\\ {\tau }_{\mathrm{mfw}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{om}}{B}_{\mathrm{om}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)}{𝜕x}}\end{array}$$

(17)

(4)

Seepage model of fluid in Region IV:

Fluid flow in Region IV consists of y-directional matrix flow and y-directional matrix-fracture transfer to Region II, with the governing differential equations for this coupled matrix-fracture system expressed as follows:

The governing differential equations describing oil-water two-phase flow in the matrix system are given by Equations (18) and (19):

$$\begin{array}{l}-{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{wm}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕p_{\mathrm{m}}}{𝜕x}}\right]-{\displaystyle \frac{𝜕\left({\rho }_{\mathrm{wm}}{\tau }_{\mathrm{mfw}}\right)}{𝜕y}}\\ -{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{4\to 2}}{𝜕y}}-G\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wm}}{\phi }_{\mathrm{m}}{S}_{\mathrm{wm}}\right)\end{array}$$

(18)

$$\begin{array}{l}-{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{om}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕p_{\mathrm{m}}}{𝜕x}}\right]-{\displaystyle \frac{𝜕\left({\rho }_{\mathrm{om}}{\tau }_{\mathrm{mfo}}\right)}{𝜕y}}\\ -{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{om}}^{4\to 2}}{𝜕y}}-G\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{om}}{\phi }_{\mathrm{m}}{S}_{\mathrm{om}}\right)\end{array}$$

(19)

The governing differential equations for oil-water two-phase flow in the fracture system are given by Equations (20) and (21):

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{wf}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rwf}}}{{\mu }_{\mathrm{w}}}}\nabla p_{\mathrm{f}}\right]+\nabla \left({\rho }_{\mathrm{wm}}\cdot {\tau }_{\mathrm{mfw}}\right)\\ -{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{4\to 2}}{𝜕y}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{wf}}{\phi }_{\mathrm{f}}{S}_{\mathrm{wf}}\right)\end{array}$$

(20)

$$\begin{array}{l}-\nabla \left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\nabla p_{\mathrm{f}}\right]+\nabla \left({\rho }_{\mathrm{om}}\cdot {\tau }_{\mathrm{mfo}}\right)\\ -{\displaystyle \frac{𝜕}{𝜕y}}\left[{\rho }_{\mathrm{of}}\cdot {\displaystyle \frac{K_{\mathrm{f}}{K}_{\mathrm{rof}}}{{\mu }_{\mathrm{o}}}}\left({\displaystyle \frac{𝜕p_{\mathrm{of}}^{4\to 2}}{𝜕y}}\right)\right]={\displaystyle \frac{𝜕}{𝜕t}}\left({\rho }_{\mathrm{of}}{\displaystyle \frac{{\phi }_{\mathrm{f}}{S}_{\mathrm{of}}}{B_{\mathrm{of}}}}\right)\end{array}$$

(21)

Among these, the equation that describes the mechanism of cross-flow from the matrix to the fractures is as follows:

$$\begin{array}{l}{\tau }_{\mathrm{mfw}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rwm}}}{{\mu }_{\mathrm{w}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{wm}}-p_{\mathrm{wf}}\right)}{𝜕y}}\\ {\tau }_{\mathrm{mfw}}=\sigma {\displaystyle \frac{K_{\mathrm{m}}{K}_{\mathrm{rom}}}{{\mu }_{\mathrm{o}}}}{\displaystyle \frac{𝜕\left(p_{\mathrm{om}}-p_{\mathrm{of}}\right)}{𝜕y}}\end{array}$$

(22)

(5)

The initial and boundary condition equations:

The initial and boundary condition equations are displayed as follows:

initial conditional equation:

$$\begin{array}{c}{p_{\mathrm{w}\mathrm{f}}|}_{\mathrm{t}=0}={p_{\mathrm{w}\mathrm{m}}|}_{\mathrm{t}=0}=p_{\mathrm{w}\mathrm{i}}(x,y,z);{S_{\mathrm{w}\mathrm{f}}|}_{\mathrm{t}=0}={S_{\mathrm{w}\mathrm{m}}|}_{\mathrm{t}=0}=S_{\mathrm{w}\mathrm{i}}(x,y,z)\\ {S_{\mathrm{o}\mathrm{f}}|}_{\mathrm{t}=0}={S_{\mathrm{o}\mathrm{m}}|}_{\mathrm{t}=0}=S_{\mathrm{o}\mathrm{i}}(x,y,z);{{\phi }_{\mathrm{f}}|}_{\mathrm{t}=0}={\phi }_{\mathrm{f}\mathrm{i}}(x,y,z);\\ {{\phi }_{\mathrm{m}}|}_{\mathrm{t}=0}={\phi }_{\mathrm{m}\mathrm{i}}(x,y,z)\\ {K_{\mathrm{f}}|}_{\mathrm{t}=0}=K_{\mathrm{f}\mathrm{i}}(x,y,z);{K_{\mathrm{m}}|}_{\mathrm{t}=0}=K_{\mathrm{m}\mathrm{i}}(x,y,z)\end{array}$$

(23)

boundary condition equation:

$${p_w|}_{\mathrm{\Gamma }}=C_w$$

$${\frac{𝜕p_{\mathrm{o}\mathrm{m}}}{𝜕n}|}_{\mathrm{\Gamma }}=0,{\frac{𝜕p_{\mathrm{o}\mathrm{f}}}{𝜕n}|}_{\mathrm{\Gamma }}=0,{\frac{𝜕p_{\mathrm{w}\mathrm{m}}}{𝜕n}|}_{\mathrm{\Gamma }}=0,{\frac{𝜕p_{\mathrm{w}\mathrm{f}}}{𝜕n}|}_{\mathrm{\Gamma }}=0$$

(24)

Based on the geological characteristics of the Huangjialing block in Fu County and the average size of the water-straight combined well group, a conceptual flow line model, featuring one straight well for water injection and one horizontal well for oil production, was also established, as shown in Figure 2. The one-injection-one-production water-straight joint well group model with varying relative positions was then simulated to investigate the water injection pattern within the network.

As shown in Figure 2, the conceptual model for a single injection and extraction process displays both centrosymmetric and axisymmetric distributions. Therefore, only a subset of the regions, labeled A, B, C, and D, requires a detailed investigation. To demonstrate this, Region A is selected for consideration. By discretizing the fracturing-reforming and non-fracturing-reforming areas and assigning the injection wells Ii-j to the corresponding regions, 22 conceptual models of the one-injection-one-mining process with various relative positions of horizontal and vertical wells were generated. The matrix of I is presented in Equation (1). All wells labeled I1-1 to I6-41 in Equation (1) correspond to potential injection well locations. However, each simulation case adheres to a one-injection-one-production configuration; only one injection well is active at any given time, while the horizontal production well remains fixed. This controlled approach enables a systematic comparison of water injection efficiency across 22 distinct injection well locations, isolating the influence of spatial configuration on development performance.

$$\mathrm{I}=\left[\begin{array}{ccccc}{\mathrm{I}}_{1\text{-}1}& {\mathrm{I}}_{16\text{-}1}& {\mathrm{I}}_{31\text{-}1}& {\mathrm{I}}_{51\text{-}1}& {\mathrm{I}}_{71\text{-}1}\\ {\mathrm{I}}_{1\text{-}11}& {\mathrm{I}}_{16\text{-}11}& {\mathrm{I}}_{31\text{-}11}& {\mathrm{I}}_{51\text{-}11}& {\mathrm{I}}_{71\text{-}11}\\ {\mathrm{I}}_{1\text{-}21}& {\mathrm{I}}_{16\text{-}21}& {\mathrm{I}}_{31\text{-}21}& {\mathrm{I}}_{51\text{-}21}& {\mathrm{I}}_{71\text{-}21}\\ {\mathrm{I}}_{1\text{-}31}& {\mathrm{I}}_{16\text{-}31}& {\mathrm{I}}_{31\text{-}31}& {\mathrm{I}}_{51\text{-}31}& {\mathrm{I}}_{71\text{-}31}\\ {\mathrm{I}}_{1\text{-}41}& {\mathrm{I}}_{16\text{-}41}& & & \end{array}\right]$$

(1)

Hydraulic fracturing is conducted around the horizontal well. The post-fracturing formation can be characterized by a combination of fracture systems and matrix. The conceptual model is depicted in Figure 3.

As shown in Figure 3, the conceptual model can be divided into four regions (I, II, III, and IV). Due to their centrosymmetric nature, only one region (Region I) was selected to study the injection parameters of the horizontal-straight well group. Region I comprises both the fracture system and the matrix system. The water injection well Wi was positioned in the conceptual model. A total of 22 one-injection-one-production configurations with different relative positions of horizontal and straight wells were obtained.

2.2. Establishment of the One-Injection-One-Production Streamline Model

Two common grid division methods include center-block grids and corner-point grids [16,17]. The former is simpler and more accurate for radial geological simulations, while the latter better simulates complex features like faults and pinch-outs [18]. For the conceptual model, corner-point grids enhance streamline distribution accuracy without requiring complex stratigraphic considerations. Therefore, corner-point grids are the most suitable choice for this model. Statistical analysis of the Huangjialing block's water-straight well group determined the model dimensions as 1400 m × 800 m × 40 m. The plane grid was divided into 140 × 80 cells, with a 10 m × 10 m step size. The tight reservoir grid division is illustrated in Figure 4. The reservoir's top depth was measured at 1425 m.

The typical tight reservoir parameters of the Fuxian Oilfield (Ordos Basin, China) are presented in Table 1. Based on these parameters, the grid properties of the one-injection-one-production conceptual model for this block were determined.

The one-injection-one-production conceptual model for this block was established. Figure 5 illustrates the streamline conceptual model corresponding to I_1-1, in which the water injection well is positioned at coordinates (1, 1).

A unified production scheme was applied to simulate the 22 one-injection-one-production models established in Section 2.1, using an injection-production ratio of 1.06 and a simulation duration of 6400 days. The simulation results, evaluated based on the production duration to reach the economic limit daily oil production, cumulative oil production, and streamline distribution in the one-injection-one-production models, are summarized in Table 2.

3. Classification of Water Injection Area

In this study, the economic limit was introduced to facilitate the classification of the water injection area. The economic limit of a horizontal well is defined as the balance between the profits from oil production and the investments in initial costs, personnel management, and equipment maintenance. In the Fuxian oilfield, the economic limit is reached when the oil production rate falls below 3 m³/d. The time at which the production of horizontal wells reaches the economic limit is termed the “limit time.”

Streamline simulations of 22 one-injection-one-production cases were conducted to classify the water injection area in this study. In these cases, the horizontal well was operated for 300 days without energy supplementation, followed by water flooding. The injection-production ratio was set to 1.1:1. To ensure the economic limit could be reached in all models, the simulation period was set to 6100 days. To investigate the effect of different relative locations between vertical and horizontal wells on oil recovery, plane interpolation was performed on the production time, cumulative oil production, and water cut data from 22 streamline models at the economic limit. The results were visualized in 3D using Matplotlib. The simulation results, including limit time, cumulative oil production, and water cut, are presented in Figure 6, Figure 7 and Figure 8.

Statistical analysis of the simulation results revealed that the relative location between the water-injection well and the horizontal well significantly influenced production. As the distance between the vertical well and the fracturing zone decreased, the horizontal well exhibited higher water cut, shorter limit time, and lower cumulative oil production. In the direction perpendicular to the horizontal well, the development effect increased monotonically with the distance between the vertical well and the horizontal wellbore. In the parallel direction, the development effect first decreased and then increased.

The distribution of formation pressure at the end of depletion development significantly influences the water-injection response. As inferred from Figure 9, the injected water can easily sweep the fracturing zone due to the significant reduction in formation pressure within the fracturing zone. Figure 10 illustrates that the streamline is densely distributed perpendicular to the wellbore in the initial period, while the opposite trend is observed at the end of depletion development. However, channeling paths can readily form near the fracturing zone [19]. Therefore, injecting large amounts of water near the fractured zone may result in low oil recovery, while the formation pressure in areas far from the fracturing zone remains high. To enhance the controlled reserves of the horizontal well, it is necessary to improve sweep efficiency through water injection. Additionally, higher water-injection intensity can be applied to enhance oil displacement in areas far from the fracturing zone.

By projecting the economic limit onto a plane (Figure 6) and applying the principle of central symmetry, the classification diagram of water injection regions (Figure 11) was obtained. The region around the horizontal well was divided into three water-injection areas based on the economic limit. Reducing the water injection intensity in region #3 or increasing the volumetric flow rate of injected water in region #1 enhances the effectiveness of water flooding.

4. Water Injection Ratio Coefficient Optimization and Mechanism Analysis

4.1. Water Injection Ratio Coefficient Optimization

As shown in Figure 12, a typical irregular horizontal-vertical well group was established in this section to determine the water injection ratio coefficient for each injection well. To optimize the water injection ratio coefficient, a series of numerical simulations were conducted based on the twelve sets of ratio coefficients (see Table 3). In this numerical simulation, the horizontal well was operated for 300 days without energy supplementation, followed by water flooding for 5800 days at an injection-production ratio of 1.1. Cumulative oil production, water cut, and production rate were used to evaluate the water injection ratio coefficient. The simulation results are presented in Figure 13, Figure 14 and Figure 15.

As shown in Figure 13, the highest cumulative oil production after 6,100 days can be achieved with the water-injection ratio coefficient of combination #4, and the cumulative oil production is 33967.08 m³, which is 2997.6 m³ higher than that the combination #12 with a balanced water-injection ratio. Injection wells in region #1 not only supplemented formation energy but also displaced the remaining oil after the depletion development. So, we can increase the water-injection ratio coefficient of region #1 injection wells to get more oil production. Because the Region #3 is closer to the fracturing area of the horizontal well, the formation energy in this area is severely deficient after depletion development. If excessive water is injected into region #3, water channeling would be occurred along with fractures. It would reduce the sweep efficiency and affect the displacement of the remaining oil. Therefore, the injection intensities in region #3 should be reduced. As shown in Figure 14, the horizontal well reaches the high water cut period (water cut> 60%) till the 6100th days in all the 12 combinations. We can also find from Figure 15 that the production rate in all of the 12 combinations is nearly the same. Therefore, the IRC should be set as 6:3:1 to yield the high oil recovery.

4.2. Mechanism Analysis of Injection Intensities Classification

To investigate the mechanism of injection intensity classification, streamline distributions under 12 water injection ratio coefficients were simulated using FrontSim (see Figure 16). The distribution range of streamlines indicates the effectiveness of water flooding. The density of streamlines correlates with the sweep efficiency of water flooding. Higher streamline density corresponds to greater sweep efficiency [20].

As shown in Figure 16, the distribution of streamlines remains similar across the 12 different water injection ratio coefficients, suggesting that the classification of water injection intensities has minimal impact on the controlled areas of water flooding. In Figure 16(1) and 14(2), streamlines are densely concentrated around I₁, indicating that the reservoir displacement around I₁ is sufficient under high water injection intensities. However, the streamlines near I₂ and I₃ are sparsely distributed, implying relatively low sweep efficiency in the water flooding process. In Figure 16(3), 14(4), and 14(5), the streamlines near the horizontal well and water injection wells are densely packed, indicating that sweep efficiency is improved and the water flooding effect is more favorable under these three water injection ratio coefficients. Compared to the other combinations, the streamline distribution in combination 4 (IRC = 6:3:1) is denser in the central area of the horizontal well, suggesting that it effectively supplements energy to the central region and enhances the water injection development process. In combination #12, with balanced water injection, the controlled area of I₁ is smaller than in combination #4. Moreover, the water flooding front of I₁ in combination #12 is positioned far from the horizontal well at day 6100, indicating that the balanced water injection method yields less effective results.

Combinations #1 (IRC = 7.5:1.5:1), #4 (IRC = 6:3:1), and #12 (IRC = 1:1:1) were selected to demonstrate the streamline distribution at days 2800 and 4200 (see Figure 17). In combination #1, the water flooding front of I₁ reaches the horizontal well by day 2800. However, by day 4200, water channeling in I₁ results in a rapid increase in water cut. The water injection volumes for each well in combination #12 are nearly identical, but because I₂ and I₃ are located close to the fracturing zone and formation energy in these areas decreases significantly after depletion development, a large pressure gradient forms after water injection, causing the water injection fronts of wells I₂ and I₃ to reach the horizontal well quickly. As a result, the oil cannot be displaced effectively. Additionally, due to insufficient water injection at well I₁, the water flooding front advances slowly, leading to sparse streamlines in the middle of the horizontal well and failing to effectively supplement the formation energy in this area. In combination #4, the water injection ratios for each well are relatively balanced. By day 2800, the water flooding fronts of wells I₂ and I₃ are near the horizontal well. Moreover, as seen in the streamline distribution at day 4200, the oil saturation around the horizontal well has not significantly decreased, indicating that the injected water effectively displaces crude oil in the reservoir towards the production well without causing water channeling.

5. Conclusions

This study develops a novel classification method for water injection regions in irregular horizontal-vertical well groups, demonstrating that optimized, region-specific water injection strategies can significantly enhance oil recovery compared to conventional balanced injection approaches. Through comprehensive numerical and streamline simulations, we establish three distinct injection regions around horizontal wells, with Region #1 requiring the highest injection intensity and Region #3 the lowest, as validated by the superior performance of the 6:3:1 injection ratio coefficient (achieving 2997.6 m³ cumulative production after 6100 days). These findings directly address our research question by proving that classification based on post-depletion pressure and flow field distributions can effectively guide waterflooding strategies to enhance formation energy while preventing premature water breakthrough. The practical implications are substantial, offering field operators a scientifically-grounded framework for designing injection schemes in irregular well patterns, particularly in similar heterogeneous reservoirs. However, we acknowledge that our current model assumes homogeneous reservoir conditions, and future work should incorporate reservoir heterogeneity, real-time dynamic adjustment of injection ratios, and potential integration with machine learning for more complex field applications. This research provides both theoretical advancements and practical methodologies that can be immediately implemented in field operations while paving the way for more sophisticated, data-driven approaches to waterflood optimization in various reservoir types.