COMSOL 6.4 - Failure of a Multilateral Well

Multilateral wells — those with multiple legs that branch off from a single conduit — can produce oil efficiently because the legs can tap multiple productive zones and navigate around impermeable ones. Unfortunately, drilling engineers must often mechanically stabilize multilateral wells with a liner or casing, which can cost millions of dollars. Leaving the wellbore uncased reduces construction costs, but it runs a relatively high risk of catastrophic failure both during installation and after pumping begins.

Figure 1: Geometry for an analysis of a horizontal open-hole multilateral well.

The poroelastic simulations estimate 3D compaction related to pumping by taking subsurface fluid flow with Darcy’s law and coupling it to structural displacements with a stress-strain analysis. This example focuses on elastic displacements brought on by changing fluid pressures when pumping begins. Related analyses for elasto-plastic deformations are straightforward using material laws automated in the Structural Mechanics Module, Nonlinear Structural Materials Module, and Geomechanics Module.

Model Definition: Flow and Deformation Simulation

The modeled geometry (Figure 1) is the lower half of a branching junction, a segment from a larger well network. The junction lies roughly 25 feet from the start of the well. The entire well network extends much further, perhaps hundreds of feet. The well is 8.5 inches in diameter and sits in a cube 80 inches on each side. Pumps move fluid from the reservoir into the well. Fluid exits the geometry only through the well. The displacement at the reservoir edge is constrained. The walls of the well, however, deform freely. The goal is to solve for the change in fluid pressure, stress, strain, and displacement that the pumping causes rather than their absolute values.

Fluid flow

To describe fluid flow, you insert the Darcy velocity into an equation of continuity

where κ is the permeability, μ is the dynamic viscosity, and p equals the pressure of the oil in the pores.

For the flow boundaries, you already know the change in fluid pressure from the well to the reservoir edge. The planar surface adjacent to the well (between the upper and lower blocks) is a symmetry boundary. Because the well is the only exit for the fluid, there is no flow to or from connecting well segments. In summary,

where n is the normal vector to the boundary.

Solid Deformation

The system of equations that describes the quasi-static deformation is

where σ denotes the stress tensor, and F are external body forces. The stress tensor is augmented by the pressure load due to changes in pore pressure. In this model, the Biot–Willis coefficient is equal to one.

The stress-strain relationship for linear materials relates the stress tensor and the strain tensor ε through the elasticity matrix C, which for isotropic materials is a function of Young’s modulus, E, and Poisson’s ratio, ν.

In this geometrically linear example, the components of the strain tensor depend on the displacement vector u, which has directional components u, v, and w:

The tensors σ and ε are linearly related by Hooke’s law σ = C: ε .

For the boundary conditions, the example constrains movement at all external boundaries. The well opening is free to deform. In summary:

Model Data

This model uses the coefficients and parameters listed in Table 1.

Table 1: Model Data.
Variable	Description	Value
ρf	Fluid density	0.0361 lb/in3
κ	Permeability	1·10-13 in2
μ	Fluid dynamic viscosity	1·10-7 psi·s
E	Young’s modulus	0.43·106 psi
ν	Poisson’s ratio	0.16
ρs	Solids density	0.0861 lb/in3
pr	Pressure in reservoir	122.45 psi
pw	Pressure in well	0 psi

Results: Flow and Deformation Simulation

The following figures show results from simulations for coupled fluid flow and reservoir deformation following a poroelastic approach for the horizontal multilateral well reported in Ref. 1.

The isosurfaces in Figure 2 indicate the fluid pressure throughout the well’s lower half. The streamlines show the fluid paths and velocities. Fluid pressure drops from the reservoir toward the well opening. The velocity typically increases toward the well but remains close to zero near the branching legs of the junction.

Figure 2: COMSOL Multiphysics poroelastic analysis of a multilateral well. Results are fluid pressure (the isosurfaces), pressure gradient (streamlines), and fluid velocities (streamline shading).

Results for elastic deformation appear in Figure 3. The surface shading denotes total displacement. The plot illustrates directional displacements by shifting the shading relative to an outline of the original geometry. For a clear view, the displacements are exaggerated. The uncased surface yields slightly because the deformed shading fills in the hollows of the well. The largest displacements occur just above the split in the well.

Figure 3: COMSOL Multiphysics estimates of displacement. Shading indicates the total displacement, and the geometry appears as lines. Even as the deformed shape shifts, those lines remain steady; the shaded image shows movement relative to the geometry outlines.

Failure Criterion

This model allows the evaluation of failures during postprocessing using results from the fluid-flow and solid-deformation simulations shown in the preceding figures. This discussion follows calculations that Ref. 1 uses to map calculations indicating where pumping could compact the reservoir enough (see Figure 3) that the well fails. Refer to Ref. 1 to estimate the critical rock strength required to successfully emplace the well, and also to learn more about calibration to data.

The 3D Coulomb failure criterion in Ref. 1 relates rock failure, the three principal stresses (σ1, σ2, and σ3), and the fluid pressures as follows:

(1)

Here So is the Coulomb cohesion and

is the Coulomb friction angle. When properly calibrated, fail = 0 indicates the onset of rock failure; fail < 0 denotes catastrophic failure; and fail > 0 predicts stability. Because this model solves for the pressure change brought on by pumping as well as the stresses, strains, and displacements that the pressure change triggers, it calculates the expression just given using the change in pressure than its absolute value.

Results: Failure Criterion

The values for the fail function appear in Figure 4. When the fail values become increasingly negative, the potential for failure is higher. As expected, the fail function estimates show the greatest potential for failure just above the split in the well.

Figure 4: Values of the fail function calculated with results from a poroelastic model for the branching junction in an open-hole multilateral well. A negative value for the fail function denotes greater potential for failure.

Conclusions

This example couples fluid flow and solid deformation for a poroelastic analysis. The analysis provides estimates of the pressure change brought on by pumping as well as the stresses, strains, and displacements that the pressure drop triggers. Combining the simulation results with a 3D Coulomb failure expression, maps vulnerability to mechanical failure from the pumping. The data and geometry for this model come from petroleum industry analyses by TerraTek (Ref. 1), which in turn use failure criteria to map the potential for failure during emplacement of the well in addition to the potential for failure when the well is pumped.

Reference

1. R. Suarez-Rivera, B.J. Begnaud, and W.J. Martin, “Numerical Analysis of Open-hole Multilateral Completions Minimizes the Risk of Costly Junction Failures,” Rio Oil & Gas Expo and Conference (IBP096_04), 2004.

Subsurface_Flow_Module/Poromechanics/multilateral_well

Modeling Instructions

From the menu, choose .

New

In the window, click

Model Wizard

In the window, click

In the tree, select > > .

Click .

Click

In the tree, select > .

Click

Geometry 1

In the window, under click .

In the window for , locate the section.

From the list, choose .

Block 1 (blk1)

In the toolbar, click

In the window for , locate the section.

In the text field, type 80.

In the text field, type 40.

Locate the section. In the text field, type -40.

In the text field, type -40.

Cylinder 1 (cyl1)

In the toolbar, click

In the window for , locate the section.

In the text field, type 4.25.

In the text field, type 80.

Locate the section. In the text field, type -6.

In the text field, type -40.

Locate the section. From the list, choose .

Cylinder 2 (cyl2)

Right-click and choose .

In the window for , locate the section.

In the text field, type -12.5.

Rotate 1 (rot1)

In the toolbar, click

and choose .

Select the object only.

In the window for , locate the section.

In the text field, type -20.5.

Locate the section. In the text field, type -6.

In the text field, type -12.5.

Difference 1 (dif1)

In the toolbar, click

and choose .

Select the object only.

In the window for , locate the section.

Click to select the

toggle button for .

Select the objects and only.

Form Union (fin)

In the toolbar, click

Global Definitions

Load some parameters from a text file. These will be used to set up the physics.

Parameters 1

In the window, under click .

In the window for , locate the section.

Click

Browse to the model’s Application Libraries folder and double-click the file multilateral_well_parameters.txt.

Geometry 1

Cleanup Log

In the toolbar, click

In the dialog that opens, click to automatically clean up the geometry.

Definitions

Variables 1

In the window, under > click .

In the window for , locate the section.

In the table, enter the following settings:


Name	Expression	Unit	Description
N	2Socos(phi)/(1-sin(phi))	Pa	Fail parameter 1
Q	(1+sin(phi))/(1-sin(phi))		Fail parameter 2
fail	(((solid.sp3+C1(p_r-p))-Q(solid.sp1+C1(p_r-p))+N(1+(solid.sp2-solid.sp1)/(solid.sp3-solid.sp1)))/C2)[1/psi]		Fail expression

Here, solid.sp1, solid.sp2, and solid.sp3 are the principal stresses. Note that the expression fail defined in Equation 1 has the dimension of stress. By appending the operator [1/psi] to the dimensionful expression enclosed within parentheses you extract its value in psi, which is what Ref. 1 uses. This way, you can analyze the risk for rock failure using the same criterion independently of your choice of base unit system.

Materials

Porous Material 1 (pmat1)

In the window, under right-click and choose > .

In the window for , locate the section.

In the εp text field, type 0.3.

Locate the section. Click

Fluid 1 (pmat1.fluid1)

In the window, click .

In the window for , click to expand the section.

In the tree, select > .

Click

Locate the section. In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Density	rho	0.0361[lb/in^3]	kg/m³	Basic
Compressibility of fluid	chif	4e-10[1/Pa]	1/Pa	Basic
Dynamic viscosity	mu	1e-7[psi*s]	Pa·s	Basic

Porous Material 1 (pmat1)

In the window, click .

In the window for , locate the section.

Click

Global Definitions

Porous Matrix

In the window, under > click .

In the window for , type Porous Matrix in the text field.

Locate the section. In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Young’s modulus	E	4.3e5[psi]	Pa	Young’s modulus and Poisson’s ratio
Poisson’s ratio	nu	0.16	1	Young’s modulus and Poisson’s ratio
Density	rho	0.0861[lb/in^3]	kg/m³	Basic
Permeability	kappa_iso ; kappaii = kappa_iso, kappaij = 0	1e-13[in^2]	m²	Basic
Biot-Willis coefficient	alphaB	1	1	Poroelastic material