Cascading Failure Phenomenon in the Multi-Stage Hydraulically Fractured Wells

Abstract The phenomenon of cascading fracture failure during flowback and initial production from a horizontal multistage hydraulically fractured well is introduced, described, and investigated. First, a simplified analytical model of production from such well is built. This model allows evaluating a range of systems parameters through which the cascading failure evolves and performing a sensitivity study of this effect. Next, while keeping the physical model of the system relatively simple, the critical flow rates causing the motion of proppant pack in fractures are treated as random variables. This assumption brings the next level of sophistication to the model and allows demonstrating nonobvious effects. In particular, the well production rate losses due to cascading fracture failure are estimated. Finally, the proposed hypothesis is validated by conducting numerical simulations of flowback in the estimated conditions of cascading failure. As a practical outcome of this study, recommendations on the mitigation of the well productivity failure caused by cascading failure are formulated and discussed.


I. INTRODUCTION
Over the last decades, horizontal well drilling and multistage hydraulic fracturing have become increasingly popular in the oilfield industry. Typically, the combination of these technologies is required to produce oil and gas from low-permeability reservoirs at economic rates [4,5,15,16]. Created fractures hydraulically connect the well to a reservoir and substantially increase well productivity roughly proportional to the total number of created fractures. In multistage fractured (MSF) wells, hydraulic fractures produce simultaneously from a reservoir to the same well and interact with each other via the shared fluid rate and bottom-hole pressure in the lateral. Such an interactive fracture behavior may lead to interesting and, at a first glance, unexpected results, such as emergent behavior and cascading fracture failure. Similarly to [18], we define this phenomena as follows. A cascading failure is a process in a system of interconnected hydraulic fractures in which the failure of one or few fractures triggers the failure of other fractures and so on. We assume that these fractures might fail if one of the parameters (e.g., pressure or rate) exceeds certain value. In this paper, we refer to this parameter as critical (critical pressure or critical rate). For example, the fractures might have different critical flow rate, below which they can produce without fracture damage. After the flow rate exceeds this value, the fracture fails and stops producing fluid to the well. Such situation may occur, for example, as a result of destabilization and washing out of the proppant from the near-wellbore region of this fracture. Under the high confining rock stress, the near-wellbore region unsupported by proppant will completely close and reduce fracture productivity to zero [10]. Consider a completed horizontal well equipped with a choke at surface and having multiple productive fractures created downhole. The mechanics of the cascading fracture failure process can be illustrated using the scenario presented in Fig. 1. The wellhead choke allows controlling the total production rate from the well. The fractures are different, so the critical parameter of failure is described by some distribution so that there is a varied failure resistance across the created fractures. The flow rate of fractures is induced by the pressure drop in the well with respect to the reservoir fluid pressure, also referred to as drawdown.
Let us assume that initially (Fig. 1,a) the surface choke is set in such a way that the well produces at a relatively small production rate Q 1 . This total production rate Q 1 is the sum of contributions from N hydraulic fractures, which produce with the same flow rate Q 1 /N . At this rate, all fractures produce without any damage. Later, the choke is opened more, and the production rate is increased to Q 2 ( Fig. 1,b). After this change, the weakest fracture in the well is damaged and stops producing (Fig. 1,c). Due to the presence of a wellhead choke, the total flow rate decrease is translated into the decrease of the bottomhole pressure and increase of the drawdown. This mechanism increases the flow rate bringing it to the value Q 3 that gets close to Q 2 . Due to the smaller number of producing fractures, it increases the producing rate at each remaining fracture, which becomes Q 3 /(N −1) > Q 2 /N ( Fig. 1,c). As a result of the production rate per active fracture ramping up, the critical flow rate is reached for two other fractures in the well (red arrows in Fig. 1,c). As the choke opening increases, the well production rate continues increase leading in turn to increased drawdown for these fractures (Fig. 1,d ). Following this fracture, other fractures continue to fail one by one (Fig. 1,e, f ). This process of cascading failure may impact either some of the fractures or all of them, so that in extreme cases the entire well productivity could be lost ( Fig. 1,f ).
In this work, we investigate mechanics of interactive fracture behavior in a horizontal well during flowback and initial production. We evaluate conditions for triggering cascading fracture failure. The content of this paper is organized as follows. First, we explain the essence of the cascading fracture failure phenomenon in MSF horizontal wells using simple illustrations. Next, we build a physically justified analytical model of well production enabled by a large set of productive but damageable fractures. This deterministic model allows predicting onset of cascading fracture failure for the given well, wellhead, fracture, and reservoir parameters. Next, we introduce probabilistic description of the MSF well. Using this description, we perform comprehensive study of the full system behavior and discuss implications of cascading failure events for the potential total loss of production. After this, we verify our analytical models by running accurate numerical simulations of flowback with consequent cascading failure using a rigorous numerical solver.

II. BASIC PHYSICAL MODEL
Consider an L-shaped well with N hydraulic fractures connected to its horizontal section and a choke at the wellhead (Fig.1). Assume that hydraulic fractures produce incompressible single-phase fluid at certain rates, depending on pressure drawdown. Because during initial production primarily water-like fluid is flowing back, the single-phase flow assumption can be justified. Individual fractures inflows q k are commingled into the total production rate at the wellhead and the choke: (1) Under such assumptions, the relation between the total mass flow rate Q and the pressure drop across choke ∆p ch = p wh − p whdc is given by (see, for example, [12]) Here p wh is the wellhead pressure, p whdc is the wellhead downstream choke pressure, d ch is the choke opening, d is the pipe diameter, A = πd 2 /4 is the pipe cross-section area, and ρ is the Pressure drop in the wellbore is assumed to be dominated by the hydrostatic term.
Accordingly, the bottomhole pressure p bh is shared between all fractures and given by the sum of the wellhead pressure p wh and the pressure drop at the vertical section of the well: Here g is the gravity acceleration and h is the depth of the horizontal well section relative to the wellhead.
The inflow from the k-th fracture is where J k is the productivity index of the fracture, which can be estimated, for example, using the Carter model [6].
An aggressive flowback strategy, when drawdown pressures and flow rates are high, can result in fracture damage. As mentioned above, the fracture conductivity can be damaged by different mechanisms that include, for example, washing of proppant out of the fracture, proppant crushing and embedment, fines migration, and tensile failure of rock and fracture faces. In this model, we assume that the k-th fracture maintains productivity as long as its rate is below some critical value, q k c or as long as the bottomhole pressure is above some critical value, p k c . Then, fractures can experience complete and irreversible loss of productivity if the failure criterion is met: or, in terms of pressure

III. PROBABILISTIC DESCRIPTION
The merit of probabilistic description is that it is not necessary to specify the exact properties of individual fractures. A modern multistage hydraulically fractured well can be connected with over a hundred hydraulic fractures, and it can be impractical to characterize parameters of every fracture by either measurements or simulations or using a data analytics approach. Instead, it can be more convenient to specify probabilities of some parameters to stay within certain ranges established statistically. Such probabilistic description allows deriving response of the large system of hydraulically interacting fractures in terms of expected values and variances of parameters.
To formulate the probabilistic description, it is convenient to further simplify physical model described in section II. Productivity indices of all fractures are assumed to be the same, Probabilistic analysis is performed during relatively short time interval, such that productivity index J 0 is assumed to be constant despite the √ t dependency typically observed in the fracture productivity behavior [6]. Assume also that the bottomhole pressure is a non-increasing function of time. Each fracture is considered to operate with the constant productivity index J 0 if its individual rate q k is less than some critical value q k c and totally and irreversibly lose productivity if q k at some moment of time exceeds q k c . Under these assumptions, the inflows are modeled by the following relation: Here θ(·) is the Heaviside step function and p k c = p R − q k c /J 0 is the critical bottomhole pressure for fracture failure.
Using equations (2) and (7), one can obtain the following system of algebraic equations describing steady states of the well: It is important to note that during flowback the well is certainly not in steady state regime so that different transient effects (e.g. pressure and rate jumps) can take place as a result of changes in boundary conditions (e.g. changes of choke opening or pressure at the wellhead).
These effects can also have impact on the evolution of cascading failure process. Nevertheless, in this work we consider these transient effects as the next order approximation to the considered problem. We demonstrate that cascading failure can be described and investigated using this steady-state formulation.
Here p bh0 = p whdc + ρgh is the bottomhole pressure for the fully opened choke and Q = N k=1 q k . Sorting the sequence p k c in increasing order, one can obtain new sequence p (k) c . Then, for each interval of the bottomhole pressure (p (8) is reduced to a quadratic equation with respect to p bh . One of the roots of this equation always exceeds reservoir pressure and has no physical meaning. Another root in the dimensionless form is given by Here reservoir pressure p R and rate of N working fractures at minimum possible bottomhole pressure N J 0 (p R − p bh0 ) are chosen as pressure and flow rate scales, Π = J 2 The candidate root p k bh must lie in the corresponding interval (p actual steady state of the well. Accordingly, we will say that for the given choke opening, the well has steady state if there exists at least one Otherwise, it is said that the well has no steady states.
We also suppose that during startup of the well the choke opening increases and the bottomhole pressure decreases. In this case, the system reaches the steady state with the maximum value of pressure first and stays there until further adjustment of choke opening.
Consequently, we study only the steady state with the maximum value of pressure and correspondingly the maximum number of producing fractures k * .    becomes the sequence of order statistics of the statistical sample p k c . Realization of k-th steady state with k fractures working is now a random event. Below we will calculate the probability P k associated with it.
The joint probability distribution function of order statistics p Let A(k, C) be the event describing that k fractures are producing at choke opening C.
The probability of A(k, C) (i.e., the k-th steady state is reached) is given by Here the integration is carried out over any suitable interval containing the support of f .
The product of the first two theta-functions shows that the k-th candidate pressure (9) is between the k-th and k + 1-th critical pressures. The rest factors show that all candidate pressures larger than k-th lie outside of corresponding intervals. Calculation of P k is outlined For illustration purposes, we consider only the case of the uniform distribution f . The probability density function f (x) and cumulative distribution F (x) in this case are given by However, the procedure presented in Appendix A is valid for an arbitrary probability distribution.
Note that the integral in Eq. (11) depends only on the values of candidate pressures p k bh . This means that the probability calculations and physical model are well separated. Therefore, it is possible to use a more sophisticated description of wellbore flow as compared to Eq. (3) and/or inflow performance relationship different from Eq. (7) to produce the sequence of p k bh . The bottomhole pressure, total flow rate, and number of surviving fractures are now random variables with some probability distributions. In the next section, we will qualitatively analyze these distributions.

IV. QUALITATIVE ANALYSIS
We expect the following schematic picture of different probability distributions under varying choke size. If the maximum candidate pressure p N bh is above the value p + the well  expected values and variances Here the rate Q k corresponds to k producing fractures and is given by for the parameters given in Table I absolute or normalized variance that we could tolerate and optimize the µ Q under constrain of σ Q bounded by this tolerance.
Another observation can be made using Figure 5 (color coding and parameters as in Below we will try to use this information about the past performance of the well to correct our expectations about its future behavior.
The probability of the event that k = k 0 fractures will produce at choke opening C = C 0 given that k 1 fractures produced at choke opening C 1 , k 2 at C 2 , . . . , k M at C M is the conditional probability given by Both numerator and denominator of the fraction in (14) are defined similarly as in Eq. (11): Here m 0 = 0 for the numerator and m 0 = 1 for the denominator in (14). Each factor in the outer product is similar to the product of theta-functions in (11). Note that here we should distinguish the candidate bottomhole pressure corresponding to the different choke sizes. From the assumption that the choke size increases in time, it follows that C m+1 < C m and p i bh (C m+1 ) < p i bh (C m ), ∀i, m. Because the fractures lose productivity irreversibly, the number of surviving fractures in future is always smaller or equal to the number of fractures surviving in the past so that k m ≤ k m+1 , ∀m. The independence of the conditional probability on the previous states, except the last one, means that the process of well evolution is Markovian. This allows us to apply to the problem of proper choke management strategy a formalism of Markov decision processes (MDP) [8,11].
We will refer to the pair (number of intact fractures, current choke parameter) as the state of well s = (k, C). In the subsequent derivations, we consider the choke parameter C We apply the well known value iteration algorithm [3] to the MDP problem described above. After the policy is found, evaluation of the number of producing fractures is required to determine the current state s t and define the next action a t . Throughout the paper, we use the assumption that the productivity index of the individual fractures are equal.
Accordingly, the number of producing fractures may be calculated using the measured flow rate. This assumption is essential for application of MDP to our problem because it reflects fundamental assumption of MDP formalism that the environment is fully observable and the current state completely characterizes the process.   to total number of fractures is the vertical coordinate. Color shows the optimal action for the state, either increase choke opening for some value or keep it constant. As one can note, for the majority of states, the optimal policy recommends keeping choke opening constant. In particular, for the choke sizes larger than 20/64 , there are no observed number of producing fractures such that any action will increase the expected value of the final flow rate. However, for the lower openings and relatively high number of surviving fractures, it is reasonable to increase the choke further. Note that the increase of the choke by more than one step at once is never optimal. This conclusion seems to be natural and may be interpreted as follows. Choke adjustments not only increase flow rate but also allow gaining information about statistics of critical rates. Cautious choke opening in small steps helps collect the information gradually without additional risk of damaging fractures and utilize it at the next steps.

VI. NUMERICAL AND ANALYTICAL SIMULATIONS
For additional illustration of cascading failure in the system of interacting hydraulic fractures and to compare predictions of statistical and deterministic approach, we performed numerical and analytical modeling of the fracture-wellbore system as shown in Fig.1 The models simulated single-phase incompressible flow in the horizontal well connected with 100 hydraulic fractures with parameters similar to Case 4 in Table I in [14], with the numerical models of individual fractures. Each fracture was represented by a grid of conductive cells with fixed width simulating a porous medium. Single-phase incompressible flow in these porous cells was described by a standard Darcy's equation [17] with constant viscosity and permeability of the proppant pack. If the flow rate from the fracture exceeded the critical value, all cells connecting the fracture to the well were updated with zero conductivity. The well and fracture models were coupled using iterative algorithm based on Picard iterations (or fixed point iterations method) [9] that ensured balance of pressures and flow rates between them within certain tolerance.
In the analytical model, we performed fine timestepping such that there could be more     fig. 8 shows that there is good agreement between all of them. We would argue that all of them can be used interchangeably as needed and as allowed by the requirements of the problem that needs to be analyzed. This can be useful for engineering purposes, when, for example, it is necessary to analyze details of the mechanism of fracture conductivity damage inside fractures for some specific realization of fracture properties. The statistical or Monte-Carlo approach can be used first to identify the specific scenario or the trajectory in the total rate versus choke opening space. Then, analytical model can be used to convert rate versus choke curve into the specific realization of critical rates distribution. The latter can be simulated with the numerical model to analyze details of the fracture damage mechanisms and investigate options to re-engineer fracture properties.

VII. CONCLUSION
In this work, the emergent behaviour leading to the formation of the cascading failure mechanism in multistage hydraulically fractured wells during well startup and initial production is described and investigated. This phenomenon can be explained by existence of a feedback loop due to the hydraulic connection between individual producing fractures connected to the same wellbore. In our study, we presented several possible scenarios of evolution of cascading failure. The failure of weaker hydraulic fractures in such systems might increase drawdown on the stronger fractures, which survive. In turn, some of surviving fractures are weaker than others, so further increase of drawdown might lead to the cascading failure of an ever-increasing amount of fractures. We also evaluated specific conditions triggering the cascading failure effect. To the best of our knowledge, the phenomenon of cascading failure of hydraulic fractures has not been studied in the previously published literature.
The hydraulic fracturing process is associated with a lack of robust measurement techniques and hence it is intrinsically uncertain. This provides grounds for treating the hydraulically fractured systems using a probabilistic approach and assuming that the failure criterion is represented as a random variable. Based on the proposed description, the analy-sis of the full system behavior is carried out. As an outcome, several choke opening strategies during the well startup could are formulated. Overall, these strategies deal with the trade-off between the short-term goal of maximizing production and long-term goal of preserving a larger number of hydraulic fractures. It is shown that moderate choke openings favor the well operating in the safe zone with productivity of all fractures preserved at the cost of relatively low production rates. Intermediate values of choke opening might lead to damage of some fractures; however, the production rates will be higher. Finally, the aggressive and large choke opening in extreme cases can lead to damage of all fractures. However, if some fractures survive with this strategy, a high production rate could be achieved. For the particular strategy aiming to maximize production rate, we propose a way to calculate choke management policy taking into account the statistical nature of fractures properties. Given the current number of producing fractures and choke opening, the policy defines the next choke adjustment maximizing the expected value of the final flow rate.
Finally, the evolution of cascading failure is demonstrated using numerical simulations that are based on fewer assumptions than analytical model and probabilistic description. It is shown that numerical results coincide qualitatively and quantitatively with the probabilistic description for the system with large number of hydraulic fractures.
First, let us make the following change of variables in (11): Next, the product of theta-functions in (11) may be transformed to The latter transformation takes into account that π k+1 > π k and consequently θ(π k−1 − ω k )θ(ω k − π k ) ≡ 0.
Here index l = k + 2, N + 1 and the starting member of the sequence is g k,k+2 (ω k+2 ) = θ(ω k+2 − π k+1 )(ω k+2 − π k+1 ). (A5) Using the previous definition, one can rewrite the equation for probability as One can check by the direct substitution into (A4) that the functions g k,l (ω l ) are polynomials multiplied by appropriate theta function: and utilizing the definition given by (A4) one can obtain Here the argument C m 0 in the term g km 0 ,k m 0 +1 (ω k m 0 +1 ; C m 0 ) shows that the values π l (C m 0 ) have been used to construct the polynomial g km 0 ,k m 0 +1 by formulas (A7), (A8).
Because the equations (B1) and (B2) differ only by the change of index m 0 → m 0 + 1 and replacement of the factor it is possible to repeat the calculation exactly in the same way and finally get Here k M +1 = N, π N +1 = 1.
Surprisingly, the conditional probability depends only on the last observed state of the well.
Probably, there exist a more elegant (than direct calculation presented here) way to prove this fact.