Simulation of Wellbore Stability during Underbalanced Drilling Operation

The wellbore stability analysis during underbalance drilling operation leads to avoiding risky problems. These problems include (1) rock failure due to stresses changes (concentration) as a result of losing the original support of removed rocks and (2) wellbore collapse due to lack of support of hydrostatic fluid column. Therefore, this paper presents an approach to simulate the wellbore stability by incorporating finite element modelling and thermoporoelastic environment to predict the instability conditions. Analytical solutions for stress distribution for isotropic and anisotropic rocks are presented to validate the presented model. Moreover, distribution of time dependent shear stresses around the wellbore is presented to be compared with rock shear strength to select appropriate weight of mud for safe underbalance drilling.


Introduction
Very recent studies highlighted that the wellbore instability problems cost the oil and gas industry above 500$-1000$ million each year [1].The instability conditions are related to rocks response to stress concentration around the wellbore during the drilling operation.That means the rock may sustain the induced stresses and the wellbore may remain stable without collapse or failure if rock strength is enormous [2].Factors that lead to formation instability are coming from the temperature effect (thermal) which is thermal diffusivity and the differences in temperature between the drilling mud and formation temperature.This can be described by the fact that if the drilling mud is too cold, this leads to decreasing the hoop stress.These variations in hoop stress have the same effect of tripping while drilling which generates swab and surge and may lead to both tensile and shear failure at the bottom of the well.
The interaction between the drilling fluids with formation fluid will cause pressure variation around the wellbore, which results in time dependent stresses changes locally [3].Therefore, in this paper the interaction between geomechanics and formation fluid [4] is taken into consideration to analyze time dependent rocks deformation around the wellbore.
Another study shows that the two main effects causing collapse failure are as follows: (1) poroelastic influence of equalized pore pressure at the wellbore wall and (2) the thermal diffusion between wellbore fluids and formation fluids [3][4][5].
Numerous scientists presented powerful models to simulate the effect of poroelastic, thermal, and chemical effects by varying values of formation pore pressure, rock failure situation, and critical mud weight [3,6].These models mentioned that controlling the component of the water present in the drilling fluid results in controlling the wellbore stability.More or less, there are many parameters that could be controlled during the drilling operation as unfavorable in situ condition [7,8].In addition, mud weight (MW)/equivalent circulation density (ECD), mud cake (mud filtrate), hole inclination and direction, and drilling/tripping practice are considered the main parameters that affect wellbore mechanical instability [9,10].
The factors that affect the mechanical stability are membrane efficiency, water activity interaction between the drilling fluid and shale formation, the thermal expansion, thermal diffusivity, and the differences in temperature between the drilling mud and formation temperature [11,12].

Derivation of Governing Equation for Thermoporoelastic Model
The equations used to simulate thermoporoelastic coupling process are momentum, mass, and energy conservation.These equations are presented in detail in this section.

Momentum Conservation.
The linear momentum balance equation in terms of total stresses can be written as follows: where  is the total stress,  is the gravity constant, and  is the bulk intensity of the porous media.The intensity should be written for two phases, liquid and solid, as follows: Equation ( 1) can be written in terms of effective stress as follows: where  is the effective stress,  is the pore pressure, and  is the identity matrix.This equation for the stress-strain relationship does not contain thermal effects and, to include the thermoelasticity, the equation can be written as follows: where  is the fourth-order stiffness tensor of material properties,  is the total strain,   is the thermal expansion coefficient, and Δ is the temperature difference.The isotropic elasticity tensor  is defined as where  is the Kronecker delta and  is the Lame constant. is the shear modulus of elasticity.The constitutive equation for the total strain-displacement relationship is defined as follows: where  →  is the displacement vector and ∇ is the gradient operator.

Mass Conservation.
The fluid flow in deformable and saturated porous media can be described by the following equation: where  is the Biots coefficient and assumed to be = 1.0 in this study,  is the pore fluid pressure,  is the temperature,   is the thermal expansion coefficient,  is the fluid flux, and  is the sink/source, and   is the specific storage which is defined by where   is the compressibility of solid and   is the compressibility of liquid.The fluid flux term () in the mass balance in (7) can be described by using Darcy's flow equation because the intensity has been assumed constant in this study: where  is the permeability of the domain.The Cubic law is used in determining fracture permeability.

Energy Conservation.
The energy balance equation for heat transport through porous media can be described as follows: where   is the heat flux,   is the heat sink/source term, and   is the heat storage and equals In this study, conduction and convection heat transfers are considered during numerical simulation.The heat flux term in (10) can be written as where V is the velocity of the fluid.The first term on the right hand side of ( 12) is the conduction term and the second term is the convective heat transfer term and  eff is the effective heat conductivity of the porous medium, which can be defined as 2.4.Discretization of the Equations.First one discretizes the thermoporoelastic governing equations by using Greens' theorem [13] to derive equations weak formulations.The weak form of mass, energy, and momentum balance in (1), (7), and (10) can be written as follows, respectively: where  is the test function, Ω is the model domain, Γ is the domain boundary,  is the traction vector, superscripts +/− refer to the value of the corresponding parameters on opposite sides of the fracture surfaces, respectively,   is the specific storage,  is the porosity,   is the volumetric Darcy flux,  is the thermal expansion coefficient,   is the fluid sink/source term between the fractures,   is the heat flux,   is the specific heat capacity,   and  ℎ are mechanical and hydraulic fracture apertures,   is the heat sink/source term,  is the thermal expansion coefficient,  is the thermal conductivity, and  refers to the fracture plane.
Then the Galerkin method is used to spatially discretize the weak forms of ( 14) to (18).The primary variables of the field problem are pressure , temperature , and displacement vector .All of these variables are approximated by using the interpolation function in finite element space as follows:  where  is the corresponding shape function and , , and  are the nodal unknowns values.

Validation of Poroelastic Numerical Model
The verification of poroelastic numerical model against analytical solutions (see Appendix) is presented in this section.A two-dimensional model of circular shaped reservoir with an intact wellbore of 1000 m drainage radius and 0.1 m wellbore radius is used (see Figure 1).The reservoir input data used are presented in Table 1.The numerical model is initiated with drained condition obtained by using Kirsch's problem [14].These conditions with the analytical solution equations for drained condition for the given pore pressure, displacement, and stresses [15,16] are presented in the Appendix.Flow chart describes the solution process for pressure and displacement for poroelastic model and also for temperature for thermoporoelastic frameworks is presented in Figure 2. The numerical results obtained are plotted against the analytical solutions in Figures 3-6.
For the verification purpose, a number of assumptions are made.
Initial State.In this study, zero time (initial state) is assumed to represent drained situation in which pore pressure is stabilized.
Boundary Conditions.They are boundary conditions for the poroelastic model in this model.
Rock and Fluid Properties.In the numerical model, Young's modulus, Poisson's ratio, porosity, permeability, and total system compressibility as well as viscosity of fluid are assumed   to be independent of time and space in order to be consistent with the analytical solutions.
As can be seen from Figure 3, the numerical results match well with the analytical solutions.Due to discontinuity of initial state and the first time step in the numerical procedure a small mismatch is observed between numerical  and analytical solutions for  = 1 hr.It is evident that, for an initial drained condition and horizontal permeability anisotropy, no directional dependence of the change in pore pressure is observed despite the anisotropic horizontal stress state.In Figure 4 the numerical results for displacement have a small mismatch with the analytical solutions.This is due to the method (Patch Recovery Method) that has been used to distribute initial reservoir displacement and calculating the change in in situ stresses with time.
In Figure 5, the numerical results show a good agreement with the exact solutions for different time and orientations.For all cases, as expected,   approaches the maximum horizontal in situ stress (5800 psi) at far field (away from wellbore).The discontinuity of   at wellbore wall is due to the imposed pressure boundary condition.It is assumed that wellbore pressure is equal to the reservoir pressure at zero time in order to simulate drained initial state.It is also observed that as time progresses, the size of the area, which is affected by the change in   , increases.This is due to change in pore pressure.
In Figure 6 numerical results of   match well with that of the analytical solutions for different time.As expected,   approaches the minimum horizontal in situ stress (5500 psi) at far field (away from wellbore).
The results of pore pressure and effective stress after one hour of production are presented in Figures 7-9 in reservoir entire region.These figures (Figures 7-9) clarify how the pressure and stresses are changing from the wall of the wellbore to the reservoir boundary.

Failure Criteria
Shear failure will occur if where   3 the lowest principle is effective stress and  0 is the rock tensile strength.Using Mohr-Coulomb criteria, shear failure criteria are met when where   1 is the highest principle effective stress,  net is the net shear stress,  is the angle of internal friction, and  0 is the rock shear strength.Once the maximum principle stress surpasses the rock shear strength, rock failure takes place at the wellbore.Therefore, evaluating the highest principle stress

Shear yielding
Figure 10: Shear yielding occurs for underbalanced conditions due to the absence of a support pressure on the borehole wall [16].
is important criterion to predict rock failure for analysis of wellbore stability [16].
Drilling with underbalanced technique where the bottom-hole pressure is lower than the formation pore pressure regularly promotes borehole instability.Thus, it is important to design and determine the ideal range of the bottom-hole pressure during underbalanced drilling operation, to avoid generating hydraulic fractures, differential sticking, or undesirable level of formation damage [17] (see Figure 10).

Case Study
This test case has been taken from a field located in southern part of Iran.The operator is considered a well drilled at an approximate depth of 4000 ft. the recorded pore pressure gradient from the DST test analysis is 7.7 lb/gal.The rock mechanical data used for wellbore stability analysis are determined from triaxial test on core samples and presented in Table 2.The wellbore stability analysis has been executed using underbalance technique.Therefore, the reduction in pore pressure during the drilling process will directly affect the horizontal and shear stresses.To avoid either loss of circulation problems or borehole failure, the mud pressure should be less than the formation fracture pressure and greater than its collapse pressure.Therefore, it is mandatory to predict the changes of stresses values as reservoir pressure drops.In this case study, the mud weight recommended to be used is 5 lb/gal.
To do this analysis, a finite element mesh has been generated as the one presented in Figure 1 and it was refined around the wellbore.The boundary conditions have been assigned to the model and Mohr-Coulomb criteria [18] are used for the simulation to predict the stresses and pore pressure distribution with time around the wellbore.

Results and Discussion
Breakout shear failure occurred during underbalance drilling operation; therefore, it is very important to predict the failure at the wellbore wall using failure criteria.Because of pore pressure dissipation, the failure becomes time dependent as the net shear stress increases with time at the borehole wall.
It can be seen from Figure 11 that the net shear stress is the lowest for the time before starting of the underbalance drilling operation.Then, at the wellbore wall, it can be noticed that the net shear stress drops suddenly after 4 s of the drilling operation.In addition, at this time, the net shear stress value is higher than net shear stress for long time of the drilling operation.This effect of short time of drilling operation on the net shear stress value is uncertain, as this time may be too short to allow failure to be devolved.But, in this case study, by comparing the net shear stress value with the rock shear strength, we found its value lower than the rock shear strength (14 MPa).This proves that the failure will not occur using mud weight of 5 lb/gal.Figures 12,13  it can be seen that the net shear stress accumulated at the wellbore wall and increases the probability of failure of the well.If the breakout occurs, it will initiate near the wellbore not at the wall bore wall (see Figure 15).If there is a breakout, the shear forces will cause rock to fall into the wellbore and in this case the well status becomes unstable (wellbore instability).But, in this case study, the net shear stress is too low to cause failure and this well will not suffer from instability problems even for long drilling period.Figures 16,17,and 18 show -stresses at the wellbore wall with the effect of mud cooling.
The cooling process near the wellbore can alter the stresses significantly and leads to increasing the total stresses and the pore pressure drop inside the formation; those increasing in total stresses and pore pressure cause increasing in the effective stresses near the wellbore (see Figures 15,17,and 18).As time increases, the mud temperature will equilibrate with its surroundings so that the formations higher in the section being drilled are subjected to the increased temperature of the mud.Heating process leads to reducing the pore pressure and net shear stresses near the wellbore (see Figures 11 and 13).
The formation cooling increases the pore pressure (see Figure 19) near the wellbore wall at the beginning of the drilling operation (for 4 s, 7.5 s, and 30 s).This is due to thermal osmosis process that results in fluid movement out Figure 20 shows a relationship between the mud weight and accumulated shear stress around the wellbore.It can be seen from the figure that, with using mud weight of 7.5 ppg, the net shear stress (16 MPa) becomes greater than the rock strength (14 MPa).Therefore, to avoid wellbore breakouts, Mohr-Coulomb failure criterion indicates that the safe mud weight used in this case study is between 5.5 and 7.5 ppg.

Conclusion
An integrated thermoporoelastic numerical model has been presented in this paper to predict the stresses distribution and the instability problem around the wall of the wellbore.The model has been validated against the analytical model.Behaviour of the stresses around the wellbore in underbalance drilling operation is very sensitive to the mud weight and mechanical properties of the rock as well.The pore pressure and stresses around the wellbore are significantly affected by the thermal effects.Thus, when the mud temperature is lower than the formation temperature, the pore pressure changes, and the net shear stresses values are increased around the wellbore which increase the probability of occurrence of the instability problem, if its values become greater than the rock shear strength.
, and 14 show stresses at the wellbore wall.Distribution of the net shear stress along -axis for the cooling (  <   ) effect of mud during the underbalance drilling operation is presented in Figure15.From this figure,

Figure 21 :
Figure 21: Schematic of the problem.

Table 1 :
Reservoir inputs used for validation of poroelastic numerical model using circular homogenous reservoir.
Flow chart describes how the nodal unknowns are solved using iterations process.

Table 2 :
Case study input data.

Table 3 :
General description of the problem.