Numerical Study of Interaction Between Hydraulic Fracture and Discrete Fracture Network

This paper discusses the interaction between hydraulic fracturing and the pre-existing discrete fracture network (DFN) in a rock mass subject to in–situ stresses. Two–dimensional computational model studies have been used in an initial attempt towards understanding how reservoir response to fluid injection is affected by some of the DFN characteristics and to operational variables such as injection rate.

1. Introduction

Understanding the mechanics of propagation of hydraulic fracture (HF) in naturally fractured reservoirs is critical to both petroleum and geothermal applications. The objective of fluid injection in these applications varies from creating HF to increasing the permeability of the surrounding rock mass, or “stimulation” of the reservoir. During stimulation, several mechanisms can lead to permeability enhancement, including:

• Opening of pre-existing fractures due to the increase in pressure or the decrease in effective normal stress (This mechanism is reversible; in other words, the fracture will close once pressure dissipates and therefore, fluid injection often needs to be accompanied with injecting proppant into the affected fractures.);

• Opening of pre-existing fractures due to slip-induced dilation, which is referred to as hydro-shearing or shear stimulation, and is permanent;

• Extension of the pre-existing fractures and increase in connectivity of the fracture network; and

• Initiation and propagation of HF, which also results in increase in fracture network connectivity.

The focus of this work is on the numerical modelling of fluid injection into the fractured rock mass and interaction between HF and the discrete fracture network. A series of comparative studies were performed to establish the effect of various in-situ parameters, including geometrical properties of the DFN, such as the level of connectivity and fracture size distribution, and operational parameters such as injection rate. In addition to qualitative evaluation of results, the model responses are compared in terms of a series of indices that were evaluated during injection. These indices include:

• Injection pressure, defined as the pressure at the injection point;

• DFN affected surface area, defined as the surface area of the DFN that has experienced a fluid pressure increase due to injection;

• Fracture surface area, calculated as a total area of fractures in the DFN;

• DFN shear stimulated surface area, defined as the area of fractures that have experienced more than 1 mm of slip;

• Leak-off ratio, defined as the volume of fluid leaked into the DFN divided by the total volume of fluid injected;

• Surface area of the HF;

• Average DFN aperture, defined as the volume of fluid injected into the DFN as a fraction of the affected surface area of the DFN; and

• Average HF aperture, defined as the volume of fluid injected into the HF divided by the surface area of the HF.

It is believed that the characteristics of the DFN have a critical effect on the response of a naturally fractured reservoir to fluid injection. Explicit representation of the DFN with realistic characteristics is thus important in the numerical modelling. Numerous realizations of different DFNs have been generated statistically and represented explicitly in the model.

2. Representation of the discrete fracture network

DFNs often are characterized by statistical parameters associated with one or more identified fracture sets. These statistical parameters typically characterize the fracture size distribution, orientation distribution and density of each fracture set.

It is widely accepted that the fracture length distribution usually follows a power law distribution, which relates the probability of occurrence of a fracture with a length of $$l$$ to the negative exponent of the length, i.e., $$n(l)\propto l^{-\alpha }$$. Value of $$\alpha$$ is site specific, but often varies in the range between 2 and 4. In this two-dimensional study, P21 is used as the measure of the fracture density. (P21 is defined as the sum of fracture or trace lengths divided by the area of the sampling or mapping domain — i.e.,$$P21={\displaystyle \sum l_i/}{L}^2$$ , where L is the linear dimension of the DFN domain.)

Considering the computational requirements of the numerical tool used in this study, it was impractical to represent DFNs with the same level of complexity as that observed in the field. Therefore, the DFN realizations were simplified or filtered. The objective of the filtering process was to reduce the geometrical complexity while preserving the relevant characteristics of the DFN. In the adopted approach, DFN realizations were simplified first by disregarding fractures with a length smaller than a prescribed threshold. The minimum fracture length cut-off is determined based on the length scale of the analysed problem. Also, to conform to what is often observed in the field, closely spaced, sub-parallel fractures (sometimes generated by the Poisson process used for generation of the fracture locations in the synthetic DFN) are disregarded. The latter criterion is based on the field observations and the reasoning that the stress field around a fracture prevents occurrence of sub-parallel fractures in its vicinity. Finally, in this study, the variation of fracture orientation about the mean for each fracture set was disregarded.

The flow characterises of the DFN are determined by identifying the clusters and evaluating overall DFN connectivity. A cluster is a group of fractures that are connected to each other; no fracture inside a cluster intersects a fracture belonging to a different cluster. A fully connected DFN is defined as the DFN with one cluster extending to the boundaries of the domain. The partially and sparsely connected DFNs were created by decreasing the fracture density and visually inspecting the size of formed clusters relative to the size of the model.

3. Numerical approach

The numerical analyses of this study are carried out using a distinct-element modelling approach. Simulations were completed using distinct element code UDEC [1]. In this approach, the fractured formation is represented by an assembly of intact rock blocks separated by a pre-existing discrete fracture network. The numerical simulations are performed using a fully coupled hydromechanical model. Fluid flow can only occur within the fractures, separating the impermeable blocks. Initially, the formation is dry. The fluid is injected in the centre of the model at constant rate.

The rock blocks are modelled as elastic and impermeable. The pre-existing fractures are represented explicitly. They are discontinuities which deform elastically, but also can open and slip (as governed by the Coulomb slip law) as a function of pressure and total stress.

UDEC can simulate fracture propagation along the predefined planes only. In order to simulate propagation of an HF, the trajectory of the fracture should be defined explicitly in the model prior to simulations. In this model, the HF is assumed to be planar, aligned with the direction of the major principal stress. The two “incipient surfaces” of the plane of the HF initially are bonded with a strength that is equivalent to specified fracture toughness. Propagation of the HF corresponded to breaking of these bonds. Clearly, the assumption of propagation of the HF in a single planar surface is a simplification. In practice, the massive hydraulic fracturing, used for example in shale stimulation, results in a large number of fractures propagating simultaneously or sequentially. Under certain conditions, the mechanical interaction between these fractures can lead to non-planar and complex trajectories as demonstrated by the results of numerical modelling [2, 3] and experimental observations [4]. Also, non-planar fracture geometry may develop as a result of the interaction with in-situ pre-existing fractures and frictional interfaces [5, 6, 7].

Figure 1 shows the geometry and the set-up of the UDEC model. The model represents a 2D horizontal section of a reservoir with a thickness of 350 m. It is assumed that the injection is through a vertical well located at the centre of the model. The core part of the model containing the DFN is embedded into a larger domain with a regular network of pipes with equivalent permeability to that of the core region. The linear dimensions of the full model are twice as large as those of the core part. In this study, the model core has the dimensions of 1000 m × 1000 m. The state of stress in the plane of the model is assumed to be anisotropic, with maximum principal horizontal stress equal to the vertical stress and the minimum principal horizontal stress equal to half of the vertical stress.

Figure 1.

Geometry and model setup.

The applied injection rate is 0.07 m³/s or 70 kg/s. This rate is approximately equal to 26.4 bpm. Considering the assumed thickness of the formation, an injection rate of $$2\times {10}^{-4}$$m³/s/m is applied in the two-dimensional model. Some sensitivity analysis with respect to injection rate is performed. It is assumed that the pre-existing fractures are already open and conductive, with a uniform aperture for each fracture set. The initial apertures of each fracture set are calculated based on their orientation relative to the in-situ principal stresses. The primary fracture set is assigned an initial aperture of $$3\times {10}^{-5}$$m, while the secondary fracture set is assigned an aperture of $$1.1\times {10}^{-5}$$m. The failure criterion of fractures is defined by the Coulomb slip law. The pre-existing fractures are assumed to have zero cohesion and the friction angle of 30°. The dilation angle is assumed to be 7.5°.

4. Results

4.1. Effect of DFN connectivity

The objective of this study was to evaluate how DFN connectivity can affect the way that injection affects propagation of the HF, and the way that the HF interacts with the pre-existing discrete fracture network. In these models, the fracture network is assumed to be static, that is, the propagation of pre-existing fractures is not allowed. However, the HF can propagate once pressure reaches the critical value, or approximately the magnitude of the minimum principal stress.

The DFN realizations used in this study are those shown in Figure 2(a). In these figures, different colours represent different clusters. The DFNs have the maximum fracture length of 1000 m and the minimum spacing of 15 m. The angles of the primary and secondary fracture sets relative to the x-axis are 160° and 60°. Figures 2(b) to (d) show the results of injection into the DFNs with various levels of connectivity. These figures show that as the connectivity of the DFN increases, pore pressure propagates to a much larger portion of the DFN. However, as the connectivity decreases, the HF tends to propagate faster. In this case, those fractures that are connected to the HF will experience increase in pressure and aperture.

Figure 3 shows quantitative comparison of the behaviour of the models based on a series of developed indices (defined in Introduction). The history of injection pressure at the injection well (Figure 3(a)) shows that as the connectivity decreases, the injection pressure increases, until it reaches the value of the hydraulic fracturing pressure. Figure 3(b) shows the history of the DFN affected surface area. It suggests that in the fully connected model, this index increases through time at a much higher rate compared to those of the partially and sparsely connected DFNs. This observation is consistent with the contours shown in Figure 2, and is the direct consequence of the presence of a much larger fracture area connected to the injection point. As a result, a greater connected permeability is available for leak-off into the DFN.

Figure 2.

Effect of DFN connectivity.

Figure 2(c), which shows the DFN shear stimulated surface area, indicates that the shear stimulated area in the partially connected DFN is greater than those of the fully and sparsely connected DFNs. This observation can be better interpreted by evaluating the pressure contours shown in Figure 1. These pressure contours indicate that the injection in the fully connected model has resulted in much lower pressures compared to the pressures in the partially and sparsely connected DFN. This is the result of greater conductivity of the fully connected DFN. In the partially and sparsely connected DFNs, injection into disconnected clusters has resulted in greater pressure increase and eventual propagation of the HF. The fractures that are connected to the injection well and the HF experience much greater pressures, which clearly result in fracture slip.

Figure 3.

Effect of DFN connectivity, history of quantitative indices.

However, in the sparsely connected DFN, the shear stimulated area is smaller than that of partially and fully connected DFN. This is due to the fact that the area of the DFN connected to the HF is much smaller than for the partially connected DFN. Therefore, even when all of those fractures were stimulated, the total stimulated area remains lower than that of the partially connected model. The non-monotonic trend of these observations is dominated by the geometry and size of clusters connected to the HF, and may vary for different DFNs.

Figure 2(d) shows the surface area of the HF. Clearly, in the fully connected DFN, no HF is formed, while in the partially and sparsely connected DFNs the length of the HF has increased with time. This graph indicates that as the degree of connectivity decreases, the HF propagates faster.

Figure 2(e) shows the leak-off ratio for fully, partially and sparsely-connected DFN realizations. For the realizations used in this study, the leak-off ratio for the fully-connected DFN is close to one from the onset of injection. The time history of the leak-off ratio in the partially and sparsely connected DFN shows that the leak-off ratio decreases as the connectivity of the DFN decreases. Also, the leak-off ratio for the partially and sparsely-connected DFNs starts to decline noticeably after certain injection time. The oscillation in the leak-off ratio is due to propagation of the HF. For example, the high points may be corresponding to times when a new cluster gets connected to the HF during the HF propagation process.

Finally, Figure 2(f) shows the average DFN aperture. It is difficult to recognize a trend form these graphs, as these observations are greatly dominated by the geometry and size of clusters connected to the injection well and the HF. However, it seems that as the degree of connectivity decreases, the average DFN aperture increases. The increase in the average DFN aperture is the direct consequence of smaller total area of the fractures which have experienced a pressure increase and greater pressures in those fractures. The average HF aperture is shown in Figure 2(g). The same trend is observed as for the average DFN aperture: the average HF aperture increases as the degree of connectivity of the DFN decreases.

4.2. Effect of fracture size distribution

The objective of this study is to evaluate the effect of the fracture size distribution on DFN response to fluid injection. In these studies, the DFNs have identical connectivity characteristics — that is, the three realizations used in the study are fully connected. However, the realizations belong to three different DFNs with different length exponent, α, and the maximum fracture length, l_max. The realizations used in this study are shown in Figure 4, in which both the primary and secondary fracture sets are orientated favourably for shear failure considering the orientation of the fracture sets relative to the direction of the principal stresses.

Figure 5 shows contours of fracture aperture and slip. In Realizations I and II, the injection rate and connectivity are such that pressures remain smaller than the minimum principal stress, (shown in Figure 6(a); thus, the HF never propagates substantially. In Realization III, the injection pressure is very close to the fracturing pressure. The leak-off ratio for all realizations remains close to one (shown in Figure 6(b)).

These figures suggest that, when compared to Realization II, Realization I, which is characterized by a narrower range of fracture size distribution with a maximum fracture length smaller that the DFN region, experiences a much larger shear stimulated surface area. Realization III seems to have the shear stimulated surface area that is smaller than for Realization I, but greater than for Realization II. This observation is quantitatively supported by the graphs of the DFN shear stimulated surface area shown in Figure 6(c). These results indicate that the shear stimulated surface area can be correlated to the exponet $$\alpha$$ and probability of having fractures large relative to the domain size. Larger $$\alpha$$ results in a narrower range of fracture sizes with a higher frequencty for fractures with the length close to $$l_{\min }$$. For example, very large values of $$\alpha$$ lead to a constant fracture size equal to $$l_{\min }$$ However, it should be noted that Realization I has the smallest $$\alpha$$, but the maximum fracture length is also capped to a value almost equal to ¼ of the DFN region size. Thus, the DFN has a fairly narrow fracture size distribution with the both mean and frequency of large fractures much smaller than those of the same distribution with uncapped $$l_{\max }$$. This condition seems to be optimum for shear stimulation.

Figure 6(d) suggests that the affected DFN surface area follows an opposite trend to the DFN shear stimulated surface area. There are two explanations for this.

• Presence of localized flow channels with large length facilitates flow and propagation of pressure front.

• Fracture dilation associated with slip creates additional volume to accommodate injected fluid, resulting in less pressure increase.

Figure 4.

DFN realizations used in the study for the effect of fracture size distribution.

Figure 5.

The effect of the fracture size distribution on response to injection (after fourteen hours of injection).

5. Effect of injection rate

The injection rate and injection pressure along with viscosity of the injected fluid are the operational parameters that can be used to effectively design hydraulic fracturing and DFN stimulation. In this section, it is evaluated how the injection rate affects propagation of the HF and reservoir stimulation. The initially considered range of injection rates is such that it covers injection pressures smaller than the hydraulic fracturing pressure.

The DFN used in this study is shown in Figure 7. Figure 8 shows contours of apertures for four injection rates of $$2\times {10}^{-5}$$m³/s/m, $$4\times {10}^{-5}$$m³/s/m, $$8\times {10}^{-5}$$m³/s/m and $$2\times {10}^{-4}$$m³/s/m. The results in these figures are compared at the time instances corresponding to identical volume injected into the formation. Figure 9(a) shows the history of injection pressure (versus time and injected volume), which indicates that for the considered injection rates, the pressure remains below the hydraulic fracturing pressure. Figure 9(b) shows that for similar injection time, higher injection rates result in a greater DFN affected area. This observation is expected as the volume injected into the DFN is higher for higher rates.

The plot of the DFN affected area versus injected volume shows a reverse trend to that observed for the DFN affected area versus time. That is, for similar injected volume, greater injection rates result in smaller DFN affected areas. For a smaller injection rate, the time required to inject a similar volume is much longer. Thus, during this longer time, the pressure front can propagate to a larger distance from the injection point.

Figure 9(c) shows the history of the DFN shear stimulated surface area. Again, for similar injection time, greater injection rates result in a greater DFN shear stimulated surface area, due to injection of larger volume of fluid. Comparing the DFN shear stimulated surface areas for similar injected volumes indicates that higher rates result in greater DFN shear stimulated surface area, which is contradictory to the trend observed for the total affected area. This observation is contributed to the fact that higher rates lead to higher pressures. The average DFN apertures as functions of time and injected volumes are shown in Figure 9(d), indicate increasing trend with increasing injection rate.

Figure 10 shows contours of apertures for the injection rates equal to and greater than $$2\times {10}^{-4}$$m³/s/m. The results are shown for similar injected volumes. The histories of injection pressures for these injection rates (models shown in Figure 10) are shown in Figure 11(a). These rates are such that they cause an injection pressure that would result in propagation of the HF. The history of leak-off ratio is shown in Figure 11(b). It is clear from the leak-off ratio that the highest injection rate of $$8\times {10}^{-4}$$m³/s/m results in propagation of the HF.

Figure 12(a) shows a similar trend to that observed in Figure 9(b) for DFN affected surface area. However, Figure 12(b) suggests that the trend in the DFN shear stimulated area is complicated and history of this index versus volume shows a decreasing trend with increasing injection rate. This change in trend is contributed to propagation of the HF for the higher injection rates, shown in Figure 12(c). Once the HF propagates, pressures remain roughly equal to the fracturing pressure. Presence of the HF, with an average aperture greater that the aperture of the DFN fractures, results in preferential fluid flow along the HF. As a result of redistribution of flow, or lower leak-off ratio, potential of shear stimulation decreases.

Figure 6.

Effect of fracture size distribution, history of quantitative indices.

Figure 7.

DFN realization used in the study for the effect of injection rate.

Figure 8.

Aperture contours for different injection rates (injected volume of 4.32 m³/s/m).

Figure 9.

Effect of injection rates, history of quantitative indices.

Figure 10.

Aperture contours for elevated injection rates (injected volume of 14 m³/m).

Figure 11.

Effect of elevated injection rates, history of quantitative indices (pressure and leak-off ratio).

These results suggest different effect of the injection rate depending on the induced injection pressures. The effect of the injection rate is evaluated for the states with the same injected volume. In general, for the range of injection pressures smaller than the hydraulic fracturing pressure:

• Smaller injection rates result in a greater total affected surface area; and

• Greater injection rates result in a greater shear stimulated surface area.

For the range of injection pressures equal to the hydraulic fracturing pressure:

• Smaller injection rates result in a greater total affected surface area; and

• Greater injection rates result in a smaller shear stimulated surface area. The results indicate that once the HF is formed, the increase in the injection rates will not be favourable for shear stimulation.

Figure 12.

Effect of elevated injection rates, history of quantitative indices.

6. Conclusions

The response of pre-existing fracture networks, typically encountered in geothermal reservoirs and shale gas formations, to fluid injection, including potential HF propagation, has been studied numerically. This is the first in a proposed series of studies intended to obtain a better understanding of the complex processes involved in DFN stimulation and hydraulic fracturing by fluid injection.

The sensitivity of the models with respect to various in-situ and operational parameters has been evaluated. The results are summarized as follows:

• DFN properties, i.e., density, length distribution and fracture orientation are critical to the overall response of the formation to injection.

• Injection pressure, and the potential for HF propagation, both increase as the connectivity of the DFN decreases.

• For a fully connected netwrok, the potential for shear stimulation increases as exponet $$\alpha$$ increases.

• Injection rate plays a major role in distributing the fluid between the HF and DFN. It is the combined effect of injection rate and effective permeability (affected greatly by DFN connectivity) that leads to the different mechanisms of behavior.

• For a given injected volume, a lower injection rate increases the proportion of the DFN that is affected (i.e., the surface area of the DFN where an increase in fluid pressure is observed).

• For a given injected volume, higher injection rates lead to a greater DFN shear-stimulated surface area, provided that the pressures remain below the hydraulic fracturing pressure. If the HF is propagated, the reverse is true:, i.e., higher injection rates result in a smaller DFN shear-stimulated surface area.

Acknowledgments

The authors acknowledge the financial support of Sandia National Laboratories for this project.