# Simulation of Wellbore Stability during Underbalanced Drilling Operation

Authors: Reda Abdel Azim
Journal: Journal of Applied Mathematics
 (2017)

---

## Abstract

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.


---

## Body

## 1. 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–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].This paper presents a realistic model to evaluate wellbore stability and predict the optimum ECD window to prevent wellbore instability problems.


## 2. 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.


### 2.1. Momentum Conservation

The linear momentum balance equation in terms of total stresses can be written as follows:(1)∇·σ+ρg=0,where σ is the total stress, g 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:(2)ρ=φρl+1-φρs.Equation (1) can be written in terms of effective stress as follows: (3)∇·σ′-pI+ρg=0,where σ is the effective stress, p is the pore pressure, and I 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: (4)σ′=Cε-αTΔT×I,where C is the fourth-order stiffness tensor of material properties, ε is the total strain, αT is the thermal expansion coefficient, and ΔT is the temperature difference. The isotropic elasticity tensor C is defined as(5)C=λδijδkl+2Gδikδjl,where δ is the Kronecker delta and λ is the Lame constant. G is the shear modulus of elasticity. The constitutive equation for the total strain-displacement relationship is defined as follows:(6)ε=12∇u→+∇u→T,where u→ is the displacement vector and ∇ is the gradient operator.


### 2.2. Mass Conservation

The fluid flow in deformable and saturated porous media can be described by the following equation:(7)Ss∂p∂t+β∇·∂u→∂t+∇·q-αT∂T∂t=Q,where β is the Biots coefficient and assumed to be = 1.0 in this study, p is the pore fluid pressure, T is the temperature, αT is the thermal expansion coefficient, q is the fluid flux, and Q is the sink/source, and Ss is the specific storage which is defined by (8)Ss=1-φKs+φKl,where Ks is the compressibility of solid and Kl is the compressibility of liquid. The fluid flux term (q) 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:(9)q=-kμ∇p-ρg→,where k is the permeability of the domain. The Cubic law is used in determining fracture permeability.


### 2.3. Energy Conservation

The energy balance equation for heat transport through porous media can be described as follows:(10)ρcpeff∂T∂t+∇·qT=QT,where qT is the heat flux, QT is the heat sink/source term, and ρcp is the heat storage and equals(11)ρcpeff=φcpρliquid+1-φcpρsolid.In this study, conduction and convection heat transfers are considered during numerical simulation. The heat flux term in (10) can be written as(12)qT=-λeff∇T+cpρliquidv·T,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(13)λeff=φλliquid+1-φλsolid.


### 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: (14)∫ΩwSs∂p∂tdΩ+∫ΩwTα∇·∂u→∂tdΩ+∫Ωwβ∂T∂t-∫Ω∇wT·qHdΩ+∫ΓHqwqH·ndΓ-∫ΩwQHdΩ=0,(15)∫ΓdwbmSs∂p∂tdΓ+∫Γdwα∂bm∂tdΓ+∫Γdwβ∂T∂tdΓ-∫Γd∇wT·bhqHdΩ+∫ΓHqwbhqH·ndΓ+∫ΓdwqH+dΓ+∫ΓdwqH-dΓ=0,(16)∫Ωwcpρ∂T∂tdΩ+∫ΩwcpρqH·∇TdΩ-∫Ω∇wT·-λ∇TdΩ+∫ΓTqw-λ∇T·ndΓ-∫ΩwTQTdΩ=0,(17)∫Γdwbmcplρl∂T∂tdΓ+∫ΓdwcplρlbhqH·∇TdΓ-∫Γd∇wT·-bmλl∇TdΓ+∫ΓTqw-bmλl∇T·ndΓ+∫ΓdwqT+dΓ+∫ΓdqT-dΓ=0,(18)∫Ω∇swT·σ′-αpIdΩ-∫ΩwT·ρgdΩ-∫ΓtwT·t-→dΓ-∫Γdw+T·t→d+dΓ-∫w-T·t→d-dΓ=0,where w is the test function, Ω is the model domain, Γ is the domain boundary, t is the traction vector, superscripts +/- refer to the value of the corresponding parameters on opposite sides of the fracture surfaces, respectively, Ss is the specific storage, n is the porosity, qH is the volumetric Darcy flux, β is the thermal expansion coefficient, QH is the fluid sink/source term between the fractures, qT is the heat flux, cp is the specific heat capacity, bm and bh are mechanical and hydraulic fracture apertures, QT is the heat sink/source term, α is the thermal expansion coefficient, λ is the thermal conductivity, and d 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 p, temperature T, and displacement vector u. All of these variables are approximated by using the interpolation function in finite element space as follows:(19)u=Nuu¯,p=Npp¯,T=NTT¯,where N is the corresponding shape function and u¯, p¯, and T¯ are the nodal unknowns values.


## 2.1. Momentum Conservation

The linear momentum balance equation in terms of total stresses can be written as follows:(1)∇·σ+ρg=0,where σ is the total stress, g 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:(2)ρ=φρl+1-φρs.Equation (1) can be written in terms of effective stress as follows: (3)∇·σ′-pI+ρg=0,where σ is the effective stress, p is the pore pressure, and I 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: (4)σ′=Cε-αTΔT×I,where C is the fourth-order stiffness tensor of material properties, ε is the total strain, αT is the thermal expansion coefficient, and ΔT is the temperature difference. The isotropic elasticity tensor C is defined as(5)C=λδijδkl+2Gδikδjl,where δ is the Kronecker delta and λ is the Lame constant. G is the shear modulus of elasticity. The constitutive equation for the total strain-displacement relationship is defined as follows:(6)ε=12∇u→+∇u→T,where u→ is the displacement vector and ∇ is the gradient operator.


## 2.2. Mass Conservation

The fluid flow in deformable and saturated porous media can be described by the following equation:(7)Ss∂p∂t+β∇·∂u→∂t+∇·q-αT∂T∂t=Q,where β is the Biots coefficient and assumed to be = 1.0 in this study, p is the pore fluid pressure, T is the temperature, αT is the thermal expansion coefficient, q is the fluid flux, and Q is the sink/source, and Ss is the specific storage which is defined by (8)Ss=1-φKs+φKl,where Ks is the compressibility of solid and Kl is the compressibility of liquid. The fluid flux term (q) 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:(9)q=-kμ∇p-ρg→,where k is the permeability of the domain. The Cubic law is used in determining fracture permeability.


## 2.3. Energy Conservation

The energy balance equation for heat transport through porous media can be described as follows:(10)ρcpeff∂T∂t+∇·qT=QT,where qT is the heat flux, QT is the heat sink/source term, and ρcp is the heat storage and equals(11)ρcpeff=φcpρliquid+1-φcpρsolid.In this study, conduction and convection heat transfers are considered during numerical simulation. The heat flux term in (10) can be written as(12)qT=-λeff∇T+cpρliquidv·T,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(13)λeff=φλliquid+1-φλsolid.


## 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: (14)∫ΩwSs∂p∂tdΩ+∫ΩwTα∇·∂u→∂tdΩ+∫Ωwβ∂T∂t-∫Ω∇wT·qHdΩ+∫ΓHqwqH·ndΓ-∫ΩwQHdΩ=0,(15)∫ΓdwbmSs∂p∂tdΓ+∫Γdwα∂bm∂tdΓ+∫Γdwβ∂T∂tdΓ-∫Γd∇wT·bhqHdΩ+∫ΓHqwbhqH·ndΓ+∫ΓdwqH+dΓ+∫ΓdwqH-dΓ=0,(16)∫Ωwcpρ∂T∂tdΩ+∫ΩwcpρqH·∇TdΩ-∫Ω∇wT·-λ∇TdΩ+∫ΓTqw-λ∇T·ndΓ-∫ΩwTQTdΩ=0,(17)∫Γdwbmcplρl∂T∂tdΓ+∫ΓdwcplρlbhqH·∇TdΓ-∫Γd∇wT·-bmλl∇TdΓ+∫ΓTqw-bmλl∇T·ndΓ+∫ΓdwqT+dΓ+∫ΓdqT-dΓ=0,(18)∫Ω∇swT·σ′-αpIdΩ-∫ΩwT·ρgdΩ-∫ΓtwT·t-→dΓ-∫Γdw+T·t→d+dΓ-∫w-T·t→d-dΓ=0,where w is the test function, Ω is the model domain, Γ is the domain boundary, t is the traction vector, superscripts +/- refer to the value of the corresponding parameters on opposite sides of the fracture surfaces, respectively, Ss is the specific storage, n is the porosity, qH is the volumetric Darcy flux, β is the thermal expansion coefficient, QH is the fluid sink/source term between the fractures, qT is the heat flux, cp is the specific heat capacity, bm and bh are mechanical and hydraulic fracture apertures, QT is the heat sink/source term, α is the thermal expansion coefficient, λ is the thermal conductivity, and d 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 p, temperature T, and displacement vector u. All of these variables are approximated by using the interpolation function in finite element space as follows:(19)u=Nuu¯,p=Npp¯,T=NTT¯,where N is the corresponding shape function and u¯, p¯, and T¯ are the nodal unknowns values.


## 3. 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 Figure1). 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.Table 1
Reservoir inputs used for validation of poroelastic numerical model using circular homogenous reservoir.
Parameter
Value
Poisson ratio
0.2
Young’s modulus
40 GPa
Maximum horizontal stress
40 MPa (5800 psi)
Minimum horizontal stress
37.9 MPa (5500 psi)
Wellbore pressure (Pw)
6.89 MPa (1000 psi)
Initial reservoir pressure (Pi)
37.9 MPa (5500 psi)
Fluid bulk module (Kf)
2.5 GPa
Fluid compressibility
1.0 × 10−5 Pa−1
Biot’s coefficient
1.0
Fluid viscosity
3 × 10−4 Pa⋅s
Matrix permeability
9.869 × 10−18 m2 (0.01 md)
Wellbore radius
0.1 m
Reservoir outer radius
1000 mFigure 1
Two-dimensional circular reservoir shape used for validation of poroelastic numerical model withσH = 39.9 MPa and σh = 37.9 MPa, Pr = 37.9 MPa, and Δp = 31 MPa.Figure 2
Flow chart describes how the nodal unknowns are solved using iterations process.Figure 3
Pore pressure as a function of radius and time in poroelastic medium withσH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.Figure 4
X-displacement along x-axis as a function of time in poroelastic medium with σH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.Figure 5
X-component of radial stresses as a function of time in poroelastic medium with σH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.Figure 6
X-component of tangential stresses as a function of time in poroelastic medium with σH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.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 Figure3, 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 t = 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 Figure4 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 Figure5, the numerical results show a good agreement with the exact solutions for different time and orientations. For all cases, as expected, σx approaches the maximum horizontal in situ stress (5800 psi) at far field (away from wellbore). The discontinuity of σx 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 σx, increases. This is due to change in pore pressure.In Figure6 numerical results of σy match well with that of the analytical solutions for different time. As expected, σy 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 Figures7–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.Figure 7
Pore pressure contour map after 1 hr of fluid production (for near-wellbore region) withσH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.Figure 8
X-component of radial stresses contour map after 1 hr of fluid production (for near-wellbore region) with σH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.Figure 9
X-component of tangential stresses contour map after 1 hr of fluid production (for near-wellbore region) with σH = 5800 psi and σh = 5500 psi, Pr = 5500 psi, Pw = 1000 psi, kx = 0.01 md, and ky = 0.01 md.


## 4. Failure Criteria

Shear failure will occur if(20)σ3′<-T0,where σ3′ the lowest principle is effective stress and T0 is the rock tensile strength.Using Mohr-Coulomb criteria, shear failure criteria are met when(21)τnet=12cos⁡φσ1′1-sin⁡φ-σ3′1+sin⁡φ>S0,where σ1′ is the highest principle effective stress, τnet is the net shear stress, φ is the angle of internal friction, and S0 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 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).Figure 10
Shear yielding occurs for underbalanced conditions due to the absence of a support pressure on the borehole wall [16].


## 5. 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 Table2. 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.Table 2
Case study input data.
Mechanical parameters
Poisson ratio (ν)
0.2
Bulk Young’s modulus (E)
4.4 GPa
Maximum horizontal stress (σH)
16.8 MPa (2436 psi)
Minimum horizontal stress(σh)
14 MPa (2030 psi)
Wellbore pressure (Pw)
6.89 MPa (1000 psi)
Hydraulic parameters
Initial reservoir pressure (Pi)
11.1 MPa (1610 psi)
Fluid bulk module (Kf)
0.45 GPa
Fluid compressibility
1.0 × 10−5 Pa−1
Biot’s coefficient
1.0
Physical parameters
Fluid viscosity
3 × 10−4 Pa⋅s
Fluid density
1111 kg/m3
Matrix permeability
9.869 × 10−18 m2 (0.01 md)
Porosity (ϕ)
0.1
Wellbore radius
0.15 m
Reservoir outer radius
1000 m
Formation temperature(Tf)
375Ko
Drilling mud temperature (Tm)
330Ko
Thermal osmosis coefficient(KT)
1 × 10−11 m2/s K
Thermal expansion coefficient of fluid(αf)
3 × 10−4 1/K
Thermal diffusivity(CT)
1.1 × 10−6 m2/s
Thermal expansion coefficient of solid(αs)
1.8 × 10−5 1/KTo do this analysis, a finite element mesh has been generated as the one presented in Figure1 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.


## 6. 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 Figure11 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, and 14 show y-stresses at the wellbore wall.Figure 11
Y-component of net shear stresses as a function of time in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 330 K.Figure 12
Y-component of maximum horizontal stresses as a function of time in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm= 330 K.Figure 13
Y-component of minimum horizontal stresses as a function of time in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 330 K.Figure 14
Y-component of net shear stresses at t = 30 s in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 330 K.Distribution of the net shear stress alongy-axis for the cooling (Tm<Tf) effect of mud during the underbalance drilling operation is presented in Figure 15. From this figure, 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 y-stresses at the wellbore wall with the effect of mud cooling.Figure 15
Y-component of net shear stresses as a function of time in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 300 K.Figure 16
Y-component of minimum horizontal stresses as a function of time in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 300 K.Figure 17
Y-component of maximum horizontal stresses as a function of time in poroelastic medium with σH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 300 K.Figure 18
Pore pressure as a function of radius and time in poroelastic medium withσH = 16.8 MPa and σh = 14 MPa, Pr = 11 MPa, Pw = 9.7 MPa, kx = 0.01 md, ky = 0.01 md, and Tm = 300 K.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 Figures15, 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 Figure19) 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 of the formation. Then, the transient response causes pore pressure on the y-axis to decrease.Figure 19
Y-component of net shear stresses at t = 30 s in poroelastic medium with σH = 2436 psi and σh = 2030 psi, Pr = 1610 psi, Pw = 1421 psi, kx = 0.01 md, ky = 0.01 md, and Tm = 300 K.Figure20 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.Figure 20
Relationship between mud weight and net shear stress (iteration process).


## 7. 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.



---

*Source: 2412397-2017-08-15.xml*