**2. Laboratory experiments for model verification and optimization**

Laboratory experiments for model verification and optimization will be performed on large rock specimen. The large-scale testing facility is currently under construction. Meanwhile preliminary experiments on smaller specimens have been conducted. In our project we focus on fracture creation in basement rock. Therefore mainly rocks like basalt, granite and gneiss are going to be tested for model verification and optimization.

#### **2.1. Large-scale**

Blocks of size 30 x 30 x 45 cm3 will be pre-stressed in a massive steel frame to set up realistic primary stress states before hydraulic stimulation (see Figure 1). Stresses in all three directions can be adjusted independently and will be held constant during injection time. After setting up the primary principal stress state with flat-jacks, the injection interval of the borehole will be pressurized by a syringe pump. Injection pressures up to 65 MPa are possible. In order to allow for the verification of the developed numerical model with the experiments, we are going to monitor the borehole-pressure, the deformation of the rock sample and localize acoustic events occurring during crack formation and propagation by means of ultrasonic transducers. Material parameters will be derived from standard rock mechanical tests.

**Figure 1.** Large- scale tri-axial testing facility under construction, able to apply up to 30 MPa vertically and 15 MPa horizontally

#### **2.2. Small-scale**

power generation. Yet this engineering of the heat exchanger has to be improved such that the

The extension to more complex fracture scenarios as well as the integration with other software for risk assessment simulations requires a computer resource moderate modeling of fracture propagation. The extended finite element method (XFEM) forms a good basis for this. It has been applied to various problems within the area of fracture mechanics. The XFEM allows for the consideration of a priori knowledge about the solution of a hydraulic fracture problem into the approximation space through the addition of enrichment functions [10]. It enables, thereby, the accurate approximation of fields that involve jumps, kinks, singularities, and other nonsmooth features within elements [2, 6, 11]. The developed numeric model will be verified against large-scale laboratory experiments. However, the focus of the present paper will lie on the progress in using XFEM for hydraulic fracture modeling. An XFEM approach in combi‐ nation with an explicit and implicit crack description is applied to a plane-strain hydraulic fracture problem. The implicit description is given by three level-set functions defined in [5] and enables a simple evaluation of the enrichments. In contrast, the explicit crack description is used to perform the crack update. Given the explicit interface, the level-set functions for each

The paper is organized as follows: After a short description of the laboratory part of the joint project in Section 2, the governing equations for a hydraulic fracture problem in its basic form are presented in Section 3. Models are discussed for the solid deformation, fluid flow, and fracture propagation. In Section 4, an XFEM formulation with an explicit-implicit interface description is introduced and the discretization of these governing equations is carried out.

Numerical results are presented in Section 5. Finally, Section 6 concludes this paper.

**2. Laboratory experiments for model verification and optimization**

are going to be tested for model verification and optimization.

**2.1. Large-scale**

Laboratory experiments for model verification and optimization will be performed on large rock specimen. The large-scale testing facility is currently under construction. Meanwhile preliminary experiments on smaller specimens have been conducted. In our project we focus on fracture creation in basement rock. Therefore mainly rocks like basalt, granite and gneiss

Blocks of size 30 x 30 x 45 cm3 will be pre-stressed in a massive steel frame to set up realistic primary stress states before hydraulic stimulation (see Figure 1). Stresses in all three directions can be adjusted independently and will be held constant during injection time. After setting up the primary principal stress state with flat-jacks, the injection interval of the borehole will be pressurized by a syringe pump. Injection pressures up to 65 MPa are possible. In order to allow for the verification of the developed numerical model with the experiments, we are going to monitor the borehole-pressure, the deformation of the rock sample and localize acoustic

outcome can be predicted within specified uncertainties.

712 Effective and Sustainable Hydraulic Fracturing

propagation step can be calculated straightforward.

Preliminary tests were performed on concrete samples of smaller size (15 x 15 x 15 cm3 ) and recently been extended to granite and basalt. Acoustic events were recorded during fracture creation and propagation. The experiments indicated the need to lower the compressive energy induced before breakdown. Further, instead of water and light oil, fluids of higher viscosity will be used from now on to enlarge the regime of fracture propagation (see Figure 2c) and to consider the lag of scaling. Optimization of acoustic emission monitoring is continuously ongoing.

### **3. Governing equations**

Hydraulic fracture propagation is based at least on three physical processes. A fluid flow within the fracture imposes a pressure load on the fracture surfaces. As a result, the rock undergoes a (mechanical) deformation and the fracture starts to propagate when a critical condition is reached [1]. Depending on the different modeling assumptions, this critical condition can be defined by the fracture toughness or another stress-based criterion. The following assumptions are usually made for the hydraulic fracture model [1]: I) the fluid flow is governed by the lubrication theory, II) solid deformation is modeled using the theory of

**3.2. Fluid flow**

leak-off *ql*

The fracturing fluid with the dynamic viscosity *μ* is modeled as laminar flow between two

The XFEM With An Explicit-Implicit Crack Description For Hydraulic Fracture Problems

http://dx.doi.org/10.5772/56383

715

**Figure 3.** Sketch of a plane-strain hydraulic fracture with varying aperture *w* and fluid front position *Lf*. A fluid lag is shown at the fracture tip. The fluid is injected at the wellbore and flows into the fracture at a constant rate *Q0*. It can

> 3 . <sup>12</sup> *w* m

( ) <sup>3</sup>

¶ - Ñ× Ñ + =

m

<sup>1</sup> (,) <sup>12</sup> *<sup>l</sup> <sup>w</sup> w q Q xy <sup>t</sup>*

and describes the conservation of the fluid mass for a Newtonian fluid. The fluid is injected into the fracture at a constant rate Q0. For a fracture propagating in an impermeable solid, the

the fluid front *Lf* and the crack tip. However, for reasons of simplicity the lag size is not part of the solution. Taking into account the symmetry of the problem, the boundary conditions

Q0 = at the fracture inlet

0

is negligible and, therefore, set to zero. It is assumed that a fluid lag develops between

d

¶ **<sup>p</sup>** (4)

<sup>2</sup> *<sup>q</sup>* (5)

at the fluid front *<sup>f</sup> q* = *q* (6)

<sup>0</sup> **p** = *p* in the fluid lag (7)

. Fluid pressure field is denoted by *p*.

**q p** =- Ñ (3)

parallel plates with a constant injection rate *Q0*. The fluid flux *q* then reads

leak off into the matrix through the fracture surfaces at a rate *ql*

The Reynolds (lubrication) equation is given by

for the fluid flow problem read as follows:

**Figure 2.** Preliminary small-scale fracturing experiments. a) Fractured concrete sample. b) Located acoustic events, projected onto face D of the sample. Different colors correspond to different time of occurrence. Dashed lines corre‐ spond to minimum, mean and maximum (along the direction of projection) height of the two main fracture surfaces, visible on the photo to the left. c) Fluid pressure (black) and flow rate (blue) curve used to fracture the specimen. Col‐ ored regions show the time and pressure regimes during which the same-colored acoustic events to the left occurred.

linear elasticity, and III) the propagation criterion is given by the conventional energy- releaserate approach of linear elastic fracture mechanics (LEFM) theory. The crack propagates when the mode I stress intensity factor reaches the fracture toughness. Each physical process is modeled separately and coupled iteratively. The governing equations are given as follows:

#### **3.1. Deformation**

A homogenous, isotropic and linear elastic solid is modeled with the concept of equilibrium of forces. Far field stresses and the pressure on the interface are imposed as Neumann boundary conditions, body forces are neglected.

$$\begin{array}{ccccc}\nabla \cdot \sigma + f &= 0 & \text{in} & \Omega\\ \mathbf{u} &= \mathbf{u}\_0 & \text{in} & \Gamma\_d\\ \sigma \cdot \mathbf{n} &= \hat{t} & \text{in} & \Gamma\_n \end{array} \tag{1}$$

*f* denotes the body force vector and *u* the displacement defined on the region Ω. The traction *t* ^ is applied at the outer boundary *Γn* and Dirichlet boundary conditions are defined on *Γ<sup>d</sup>* . Hooke's law of elasticity gives the relation between the stress *σ* and the strain *ε*

$$
\sigma\_{\psi} = \mathbf{C}\_{\psi kl} \mathbf{c}\_{kl} \tag{2}
$$

where *C* is the fourth-order stiffness (elasticity) tensor. Since the fracture aperture *w* is not given directly in this formulation, it has to be determined from the displacement field.

#### **3.2. Fluid flow**

linear elasticity, and III) the propagation criterion is given by the conventional energy- releaserate approach of linear elastic fracture mechanics (LEFM) theory. The crack propagates when the mode I stress intensity factor reaches the fracture toughness. Each physical process is modeled separately and coupled iteratively. The governing equations are given as follows:

**Figure 2.** Preliminary small-scale fracturing experiments. a) Fractured concrete sample. b) Located acoustic events, projected onto face D of the sample. Different colors correspond to different time of occurrence. Dashed lines corre‐ spond to minimum, mean and maximum (along the direction of projection) height of the two main fracture surfaces, visible on the photo to the left. c) Fluid pressure (black) and flow rate (blue) curve used to fracture the specimen. Col‐ ored regions show the time and pressure regimes during which the same-colored acoustic events to the left occurred.

A homogenous, isotropic and linear elastic solid is modeled with the concept of equilibrium of forces. Far field stresses and the pressure on the interface are imposed as Neumann

0 in

ˆ in

*f* denotes the body force vector and *u* the displacement defined on the region Ω. The traction

is applied at the outer boundary *Γn* and Dirichlet boundary conditions are defined on *Γ<sup>d</sup>* .

 e

where *C* is the fourth-order stiffness (elasticity) tensor. Since the fracture aperture *w* is not given directly in this formulation, it has to be determined from the displacement field.

in

*d n*

*ij ijkl kl* = *C* (2)

(1)

G

 G

0

Ñ× + = W

*t*

Hooke's law of elasticity gives the relation between the stress *σ* and the strain *ε*

s

=

*f*

× = **u u n**

s

s

**3.1. Deformation**

714 Effective and Sustainable Hydraulic Fracturing

*t* ^

boundary conditions, body forces are neglected.

The fracturing fluid with the dynamic viscosity *μ* is modeled as laminar flow between two parallel plates with a constant injection rate *Q0*. The fluid flux *q* then reads

**Figure 3.** Sketch of a plane-strain hydraulic fracture with varying aperture *w* and fluid front position *Lf*. A fluid lag is shown at the fracture tip. The fluid is injected at the wellbore and flows into the fracture at a constant rate *Q0*. It can leak off into the matrix through the fracture surfaces at a rate *ql* . Fluid pressure field is denoted by *p*.

$$\mathbf{q} = -\frac{w^3}{12\mu} \nabla \mathbf{p}.\tag{3}$$

The Reynolds (lubrication) equation is given by

$$\frac{\partial w}{\partial t} - \frac{1}{12\,\mu} \nabla \cdot \left( w^3 \nabla \mathbf{p} \right) + q\_l = Q\_0 \delta(\mathbf{x}, \mathbf{y}) \tag{4}$$

and describes the conservation of the fluid mass for a Newtonian fluid. The fluid is injected into the fracture at a constant rate Q0. For a fracture propagating in an impermeable solid, the leak-off *ql* is negligible and, therefore, set to zero. It is assumed that a fluid lag develops between the fluid front *Lf* and the crack tip. However, for reasons of simplicity the lag size is not part of the solution. Taking into account the symmetry of the problem, the boundary conditions for the fluid flow problem read as follows:

$$q = \frac{\mathbb{Q}\_0}{2} \text{ at the fracture inlet} \tag{5}$$

$$q = q\_f \quad \text{at the fluid front} \tag{6}$$

$$\mathbf{p} = p\_0 \text{ in the fluid lag} \tag{7}$$

The flow condition at the fluid front *qf* is determined from the pressure gradient and thus, is part of the solution. The pressure in the lag region is set to a constant value *p0*, that is usually chosen to be zero. Finally, the global volume balance condition

$$Q\_0 t = \int\_{\Gamma\_f} w \, d\Gamma + \int\_0^t \int\_{\Gamma\_f} q\_l \, d\Gamma \, d\tau = V \tag{8}$$

The enrichments, that are realized through the partition of unity (PU) concept, are chosen in such a way that they are able to reproduce the asymptotic behavior near the crack tip. In this work only the toughness-dominated solution is considered. Thus, enrichment functions compatible with the classical square root singularity of LEFM are used to enrich the region

The XFEM With An Explicit-Implicit Crack Description For Hydraulic Fracture Problems

The XFEM formulation with an explicit-implicit crack description used in this work is based on the work done by [5]. The basic idea is recalled in this paper, for further details see the

<sup>4</sup> \* \*

*i i j step j k tip k*

14243 1444444444442444444444443

The first term on the right hand side describes the classical FEM-approximation with contin‐

tinuity in the displacement field across the crack path by incorporating step-functions *Ψstep* with additional nodal unknowns *ai* into the enrichment space. The tip region is enriched with a set

{ } { } <sup>4</sup> <sup>1</sup> (, ) ( ), ( ), ( ) ( ), ( ) ( ) *<sup>m</sup> tip <sup>m</sup> r r cos r sin r sin sin r sin cos*

where *r* and *θ* denote local polar coordinates at the crack tip. When propagating in the toughness-dominated regime the assumption of a square root singularity in LEFM requires to choose *λ* = 1/2. In the viscosity-dominated regime the weaker singularity is taken into account

The crack opening is obtained through interpolation of the displacement field *u*(*x*) at the interface nodes by means of (10). Since the interface represents a discontinuity the interpolation is performed by moving the nodes slightly away in normal direction. The opening is defined as the distance between the positive and negative displacement at the interface (see Figure 5).

 q

<sup>=</sup> Y = (11)

 lq 1

*<sup>m</sup>*(*r*, *θ*)that consider the singularity according to the dominating

 lq  l

*wx u x u x* ( ) ( ) ( ). + - = - (12)

 q

*<sup>m</sup>* are introduced into the approximation locally

 lq

*r* q

http://dx.doi.org/10.5772/56383

. The second term accounts for the discon‐

(10)

717

original work. The enriched approximation of the displacements is stated as follows:

() () () () () (, ) ( ) *step tip h m m*

continuous discontinuous

*(x)* and nodal unknowns *ui*

lll

lq

*i I j I k I m*

Î = Î Î = + ×Y + × Y å å åå

**u x N xu N x xa N x b**

near the crack tip.

**4.1. Standard formulation**

uous shape functions *Ni*

of enrichment functions *Ψtip*

regime. They can be defined as [10]

q

within the enriched region.

by *λ* = 2/3. Additional degrees of freedom *bk*

equates the fracture volume *V* to the volume of injected fluid and the amount lost to the surrounding rock-mass. The integration is performed over the fluid filled part of the crack *Γ<sup>f</sup>* .

#### **3.3. Propagation condition**

Due to the symmetry in loading and geometry, the hydraulic fracture propagates in pure opening mode, i.e. the tensile stress is acting normal to the plane of the crack. The propagation criterion is formulated in the framework of linear elastic fracture mechanics (LEFM) and accounts for the energy required to break the rock. It is characterized by the stress singularity at the tip and a propagation in mobile equilibrium. The LEFM assumption requires the stress intensity factor *KI* to be equal to the fracture toughness *KIC* of the material [12]

$$
\mathbb{K}\_I = \mathbb{K}\_{IC}.\tag{9}
$$

#### **3.4. Asymptotic behavior**

The hydraulic fracture problem characterized by a strong fluid-solid coupling that is mainly confined to a small region near the crack tip where rapid variation in the fluid pressure occurs. Analyzing the physical process at the tip by comparing the work done by the fluid in extend‐ ing a fracture with the energy required to create new crack surfaces leads to understanding of the propagation regime of a fluid-driven fracture. Two limiting regimes can be detected, a toughness- and a viscosity-dominated regime [3]. In the toughness-dominated regime the inverse square root singularity of LEFM captures the effect of the crack tip process on the total fracture. In contrast, the viscosity-dominated regime is characterized by a singularity that is weaker than the singularity predicted by LEFM. Fracture toughness *KIC* may become irrele‐ vant [9].

#### **4. XFEM approximation**

The extended finite element method (XFEM) allows for the consideration of a priori knowl‐ edge about the solution of a hydraulic fracture problem into the approximation space through the addition of enrichment functions [10]. It enables, thereby, the accurate approximation of fields that involve jumps, kinks, singularities, and other non-smooth features within ele‐ ments [2, 6, 11].

The enrichments, that are realized through the partition of unity (PU) concept, are chosen in such a way that they are able to reproduce the asymptotic behavior near the crack tip. In this work only the toughness-dominated solution is considered. Thus, enrichment functions compatible with the classical square root singularity of LEFM are used to enrich the region near the crack tip.

#### **4.1. Standard formulation**

The flow condition at the fluid front *qf* is determined from the pressure gradient and thus, is part of the solution. The pressure in the lag region is set to a constant value *p0*, that is usually

d dd

equates the fracture volume *V* to the volume of injected fluid and the amount lost to the surrounding rock-mass. The integration is performed over the fluid filled part of the crack *Γ<sup>f</sup>* .

Due to the symmetry in loading and geometry, the hydraulic fracture propagates in pure opening mode, i.e. the tensile stress is acting normal to the plane of the crack. The propagation criterion is formulated in the framework of linear elastic fracture mechanics (LEFM) and accounts for the energy required to break the rock. It is characterized by the stress singularity at the tip and a propagation in mobile equilibrium. The LEFM assumption requires the stress

to be equal to the fracture toughness *KIC* of the material [12]

The hydraulic fracture problem characterized by a strong fluid-solid coupling that is mainly confined to a small region near the crack tip where rapid variation in the fluid pressure occurs. Analyzing the physical process at the tip by comparing the work done by the fluid in extend‐ ing a fracture with the energy required to create new crack surfaces leads to understanding of the propagation regime of a fluid-driven fracture. Two limiting regimes can be detected, a toughness- and a viscosity-dominated regime [3]. In the toughness-dominated regime the inverse square root singularity of LEFM captures the effect of the crack tip process on the total fracture. In contrast, the viscosity-dominated regime is characterized by a singularity that is weaker than the singularity predicted by LEFM. Fracture toughness *KIC* may become irrele‐

The extended finite element method (XFEM) allows for the consideration of a priori knowl‐ edge about the solution of a hydraulic fracture problem into the approximation space through the addition of enrichment functions [10]. It enables, thereby, the accurate approximation of fields that involve jumps, kinks, singularities, and other non-smooth features within ele‐

t

= G+ G = ò òò (8)

. *K K I IC* = (9)

*t Qt w q V <sup>l</sup>*

chosen to be zero. Finally, the global volume balance condition

**3.3. Propagation condition**

716 Effective and Sustainable Hydraulic Fracturing

intensity factor *KI*

vant [9].

ments [2, 6, 11].

**3.4. Asymptotic behavior**

**4. XFEM approximation**

0 0

*f f*

G G

The XFEM formulation with an explicit-implicit crack description used in this work is based on the work done by [5]. The basic idea is recalled in this paper, for further details see the original work. The enriched approximation of the displacements is stated as follows:

$$\mathbf{u}^{h}(\mathbf{x}) = \underbrace{\sum\_{i \in I} \mathbf{N}\_{i}(\mathbf{x}) \mathbf{u}\_{i} + \sum\_{j \in I^{step}} \mathbf{N}\_{j}^{\*}(\mathbf{x}) \cdot \Psi\_{step}(\mathbf{x}) \mathbf{a}\_{j} + \sum\_{k \in I^{\#}} \mathbf{N}\_{k}^{\*}(\mathbf{x}) \cdot \left(\sum\_{m=1}^{4} \Psi\_{tip}^{m}(r, \theta) \mathbf{b}\_{k}^{m}\right)}\_{\text{discontinuous}} \tag{10}$$

The first term on the right hand side describes the classical FEM-approximation with contin‐ uous shape functions *Ni (x)* and nodal unknowns *ui* . The second term accounts for the discon‐ tinuity in the displacement field across the crack path by incorporating step-functions *Ψstep* with additional nodal unknowns *ai* into the enrichment space. The tip region is enriched with a set of enrichment functions *Ψtip <sup>m</sup>*(*r*, *θ*)that consider the singularity according to the dominating regime. They can be defined as [10]

$$\left\{ \left. \Psi\_{t\bar{\eta}}^{\mathfrak{m}} (r, \theta) \right\} \right\}\_{\mathfrak{m} = 1}^{4} = \left\{ r^{\lambda} \cos(\lambda \theta), r^{\lambda} \sin(\lambda \theta), r^{\lambda} \sin(\theta) \sin(\lambda \theta), r^{\lambda} \sin(\theta) \cos(\lambda \theta) \right\} \tag{11}$$

where *r* and *θ* denote local polar coordinates at the crack tip. When propagating in the toughness-dominated regime the assumption of a square root singularity in LEFM requires to choose *λ* = 1/2. In the viscosity-dominated regime the weaker singularity is taken into account by *λ* = 2/3. Additional degrees of freedom *bk <sup>m</sup>* are introduced into the approximation locally within the enriched region.

The crack opening is obtained through interpolation of the displacement field *u*(*x*) at the interface nodes by means of (10). Since the interface represents a discontinuity the interpolation is performed by moving the nodes slightly away in normal direction. The opening is defined as the distance between the positive and negative displacement at the interface (see Figure 5).

$$
\mu w(\mathbf{x}) = \mu^+(\mathbf{x}) - \mu^-(\mathbf{x}).\tag{12}
$$

crack update is applied by simply adding new elements to the crack front according to a given

The XFEM With An Explicit-Implicit Crack Description For Hydraulic Fracture Problems

http://dx.doi.org/10.5772/56383

719

(a) (b) (c)

**Figure 6.** (a) Arbitrary crack surface with normal vectors on each element. (b) Local coordinate system at the crack

The implicit interface description is realized by means of the level-set concept. Three level-set functions are defined according to [5]. They are used to define the region to be enriched and

**•** *Φ*1(*x*) is the (un-signed) distance function to the crack path/surface. That is, the level-set

**•** *Φ*2(*x*) is the (un-signed) distance function to the crack tip(s)/front. That is, the level-set value

**•** *Φ*3(*x*) is a signed distance function to crack path/surface that is extended over the entire domain. The sign is based on the direction of the normal vector of the segment that contains

Since the solid deformation and the fluid flow are coupled iteratively, they are solved inde‐ pendently in each iteration step. Solid deformation is discretized with XFEM as follows:

> ˆ *b c T T d td pd*

where the term on the left **BTCB** denotes the stiffness matrix with the gradient operator **B** [4],

pressure on the interface *Γc*. A classical FEM approach is used to solve the fluid flow equation.

ê ú ë û ò òò **<sup>T</sup> B CB u N N** (13)

the traction on the outer boundary *Γb* and **p** the

W GG

ê ú W × = G+ G

^

value at position x is the shortest distance to the crack path/surface.

at position x is the shortest distance to the crack tip(s)/front.

é ù

front. (c) Crack update according to crack extension vectors at the front.

to evaluate the enrichment functions.

**4.4. Discretization of governing equations**

*N* the shape and enrichment functions, *t*

The pressure is approximated by

the nearest point.

extension vector.

**Figure 4.** The enrichment is acting either along the crack path (dashed field) with the step-enrichment Ψ*step* or within a specified region near the crack tip by defining the enrichments Ψ*tip <sup>m</sup>*(*r*, θ).

**Figure 5.** Interpolation of the crack opening along the interface.

#### **4.2. Numerical integration**

The standard approach in the XFEM for numerical integration is a decomposition of the elements into subelements that align with the discontinuity [6]. A Gauss quadrature is then applied on each of these subelements. For a detailed description of the decomposition method in 2D and 3D the reader is referred to [6].

#### **4.3. Explicit-implicit interface description**

The explicit crack description is given by a mesh that is aligned with the interface. For 2D problems the crack is a line and is represented by one dimensional elements in the 2D space. In three dimensions, the crack is a surface and, thus, described through a two dimensional mesh in the 3D space.

Normal and tangential vectors are computed easily on the interface and can be used to define a local coordinate system at the crack tip/front. On the basis of the explicit interface mesh the crack update is applied by simply adding new elements to the crack front according to a given extension vector.

**Figure 6.** (a) Arbitrary crack surface with normal vectors on each element. (b) Local coordinate system at the crack front. (c) Crack update according to crack extension vectors at the front.

The implicit interface description is realized by means of the level-set concept. Three level-set functions are defined according to [5]. They are used to define the region to be enriched and to evaluate the enrichment functions.


#### **4.4. Discretization of governing equations**

**Figure 5.** Interpolation of the crack opening along the interface.

a specified region near the crack tip by defining the enrichments Ψ*tip*

in 2D and 3D the reader is referred to [6].

**4.3. Explicit-implicit interface description**

The standard approach in the XFEM for numerical integration is a decomposition of the elements into subelements that align with the discontinuity [6]. A Gauss quadrature is then applied on each of these subelements. For a detailed description of the decomposition method

**Figure 4.** The enrichment is acting either along the crack path (dashed field) with the step-enrichment Ψ*step* or within

*<sup>m</sup>*(*r*, θ).

The explicit crack description is given by a mesh that is aligned with the interface. For 2D problems the crack is a line and is represented by one dimensional elements in the 2D space. In three dimensions, the crack is a surface and, thus, described through a two dimensional

Normal and tangential vectors are computed easily on the interface and can be used to define a local coordinate system at the crack tip/front. On the basis of the explicit interface mesh the

**4.2. Numerical integration**

718 Effective and Sustainable Hydraulic Fracturing

mesh in the 3D space.

Since the solid deformation and the fluid flow are coupled iteratively, they are solved inde‐ pendently in each iteration step. Solid deformation is discretized with XFEM as follows:

$$
\left[\int\_{\Omega} \mathbf{B}^T \mathbf{C} \mathbf{B} \, d\Omega \right] \cdot \mathbf{u} = \int\_{\Gamma\_b} \mathbf{N}^T \hat{t} \, d\Gamma + \int\_{\Gamma\_c} \mathbf{N}^T p \, d\Gamma \tag{13}
$$

where the term on the left **BTCB** denotes the stiffness matrix with the gradient operator **B** [4], *N* the shape and enrichment functions, *t* ^ the traction on the outer boundary *Γb* and **p** the pressure on the interface *Γc*. A classical FEM approach is used to solve the fluid flow equation. The pressure is approximated by

$$\mathbf{p}^{h}(\mathbf{x}) = \sum\_{i \neq I} \mathbf{N}\_{i}(\mathbf{x}) \mathbf{p}\_{i}. \tag{14}$$

distribution and the crack opening are calculated until convergence is reached. When the propagation condition is met, the crack is updated for the next time step. Otherwise, the fluid front is moved towards the crack tip with a velocity *v* determined from the fluid flow rate *q*.

The XFEM With An Explicit-Implicit Crack Description For Hydraulic Fracture Problems

http://dx.doi.org/10.5772/56383

721

The numerical results for the crack opening and pressure distribution at the wellbore of a plane-strain hydraulic fracture problem are compared to the similarity solution of a small enough toughness parameter in order to allow a significant fluid lag. The boundary condition of zero displacement at infinity is approximated by a finite body and standard finite elements and a local mesh refinement in the area close to the crack interface. Computational evidence of the validity of approximating the infinite medium with a finite block is provided in [13].

**Figure 7.** The numerical results (blue circles) of dimensionless pressure Π (a) and the crack opening Ω (b) at the well‐ bore as well as the dimensionless crack length γ (c) are compared to the analytical solution (red solid line) for various

The results are scaled to dimensionless quantities in the viscosity scaling. For a detailed description of the scaling for the pressure Π, the opening Ω and the crack length γ see the original publication [13]. The domain and the explicit interface are meshed independently with

Figures 7(a)-(c) show a good agreement of the similarity solution for various values of the fluid fraction *ξ<sup>f</sup>* = *L <sup>f</sup>* / *L* . However, especially for high fluid fraction values where the fluid front is close to the fracture front the results reveal inaccuracies. Special attention has to be paid to the crack tip behavior in the case for a vanishing fluid lag when the pressure becomes singular. Depending on the propagation regime crack propagation is governed either by the classical singularity of linear elastic fracture mechanics or by viscous fluid effects which would lead to

**5.2. Numerical results**

values of the fluid front position.

5000 and 3000 elements, respectively.

a weaker singularity than given by LEFM.

The discretized problem formulation reads

$$\frac{1}{12\mu} \int\_{\Gamma\_f} w^3 \nabla \mathbf{N} \nabla \mathbf{N}^T \, d\Gamma = \left[ \mathbf{N} \frac{w^3}{12\mu} \nabla \mathbf{p} \right]\_{\Gamma\_f} - \int\_{\Gamma\_f} \mathbf{N} \, q\_l \, d\Gamma - \int\_{\Gamma\_f} \mathbf{N} \frac{\Delta w}{\Delta t} d\Gamma \tag{15}$$

This formulation is valid for one half of a symmetric crack where *Γ<sup>f</sup>* denotes the fluid filled region. The flow boundary conditions at the fluid front and the fracture inlet correspond to the first term on the right-hand side. Fluid leak-off and the change of volume over time are taken into account by the second and third term on the right-hand side.

#### **5. Hydraulic fracture propagation**

The problem of a fluid driven fracture in an impermeable elastic solid with a fluid lag is considered here. Simulation results are compared with the asymptotic solutions for zero underpressure/time given in [8]. This solution corresponds to the "beginning" of the fluiddriven fracture evolution and provides initial condition for plane-strain fracture propagations. The propagation regime of a fluid driven fracture is controlled by a parameter representing a dimensionless viscosity *M* (dimensionless toughness *K*) defined as

$$M = \frac{\mu' Q\_0}{E'} \left(\frac{E'}{K'}\right)^4, \qquad K = M^{-1/4}.\tag{16}$$

This formulation uses effective parameters [6]

$$E' = \frac{E}{1 - \nu^2}, \qquad \mu' = 12\,\mu, \qquad K' = 4\left(\frac{2}{pi}\right)^{1/2} K\_{IC'} \tag{17}$$

where μ′ denotes the fluid viscosity, Q0 the constant injection rate, E′ the plane-strain elastic modulus with Poisson's ratio ν and K′ the toughness, respectively. The procedure solving the coupled equations follows that described in [8]. Given a fluid front position *Lf* , a solution is sought for the pressure distribution and crack opening.

#### **5.1. Numerical algorithm**

The simulation process is realized through an iterative coupling of the fluid flow and solid deformation. Starting with an initial solution and a guess for the fluid fraction, the pressure distribution and the crack opening are calculated until convergence is reached. When the propagation condition is met, the crack is updated for the next time step. Otherwise, the fluid front is moved towards the crack tip with a velocity *v* determined from the fluid flow rate *q*.

#### **5.2. Numerical results**

() () .

*w w w d qd d*

This formulation is valid for one half of a symmetric crack where *Γ<sup>f</sup>* denotes the fluid filled region. The flow boundary conditions at the fluid front and the fracture inlet correspond to the first term on the right-hand side. Fluid leak-off and the change of volume over time are

The problem of a fluid driven fracture in an impermeable elastic solid with a fluid lag is considered here. Simulation results are compared with the asymptotic solutions for zero underpressure/time given in [8]. This solution corresponds to the "beginning" of the fluiddriven fracture evolution and provides initial condition for plane-strain fracture propagations. The propagation regime of a fluid driven fracture is controlled by a parameter representing a

> 4 <sup>0</sup> 1/4 ., *<sup>Q</sup> <sup>E</sup> <sup>M</sup> K M*


<sup>2</sup> , 12 , 4 <sup>1</sup> , *IC <sup>E</sup> <sup>E</sup> K K*

where μ′ denotes the fluid viscosity, Q0 the constant injection rate, E′ the plane-strain elastic modulus with Poisson's ratio ν and K′ the toughness, respectively. The procedure solving the

The simulation process is realized through an iterative coupling of the fluid flow and solid deformation. Starting with an initial solution and a guess for the fluid fraction, the pressure

 m

æ ö ¢ ¢¢ = == ç ÷ - è ø

coupled equations follows that described in [8]. Given a fluid front position *Lf*

é ù <sup>D</sup> Ñ Ñ G= Ñ - G- G ê ú

*i I* Î

*f f f f*

 m*<sup>t</sup>* <sup>G</sup> <sup>G</sup> G G

*i i*

**p x N xp** <sup>=</sup> å (14)

(16)

(17)

, a solution is

*l*

1/2

*pi*

ê ú <sup>D</sup> ë û <sup>ò</sup> **NN N p N N** ò ò (15)

*h*

taken into account by the second and third term on the right-hand side.

dimensionless viscosity *M* (dimensionless toughness *K*) defined as

m

2

n

sought for the pressure distribution and crack opening.

**5.1. Numerical algorithm**

*E K*

m

The discretized problem formulation reads

**5. Hydraulic fracture propagation**

This formulation uses effective parameters [6]

m

720 Effective and Sustainable Hydraulic Fracturing

<sup>3</sup> 1 <sup>3</sup> 12 12

*T*

The numerical results for the crack opening and pressure distribution at the wellbore of a plane-strain hydraulic fracture problem are compared to the similarity solution of a small enough toughness parameter in order to allow a significant fluid lag. The boundary condition of zero displacement at infinity is approximated by a finite body and standard finite elements and a local mesh refinement in the area close to the crack interface. Computational evidence of the validity of approximating the infinite medium with a finite block is provided in [13].

**Figure 7.** The numerical results (blue circles) of dimensionless pressure Π (a) and the crack opening Ω (b) at the well‐ bore as well as the dimensionless crack length γ (c) are compared to the analytical solution (red solid line) for various values of the fluid front position.

The results are scaled to dimensionless quantities in the viscosity scaling. For a detailed description of the scaling for the pressure Π, the opening Ω and the crack length γ see the original publication [13]. The domain and the explicit interface are meshed independently with 5000 and 3000 elements, respectively.

Figures 7(a)-(c) show a good agreement of the similarity solution for various values of the fluid fraction *ξ<sup>f</sup>* = *L <sup>f</sup>* / *L* . However, especially for high fluid fraction values where the fluid front is close to the fracture front the results reveal inaccuracies. Special attention has to be paid to the crack tip behavior in the case for a vanishing fluid lag when the pressure becomes singular. Depending on the propagation regime crack propagation is governed either by the classical singularity of linear elastic fracture mechanics or by viscous fluid effects which would lead to a weaker singularity than given by LEFM.

**Author details**

, P. Siebert2

versity, Aachen, Germany

44(5), 739-757.

73(4), 456-481.

(2006).

Sons, Ltd, (2007). , 215-247.

601-620.

, K. Willbrand3

, M. Feinendegen2

1 Chair for Computational Analysis of Technical System, RWTH Aachen University, Aa‐

3 Institute for Applied Geophysics and Geothermal Energy, E.ON ERC, RWTH Aachen Uni‐

[1] Adachi, J, Siebrits, E, Peirce, A, & Desroches, J. Computer simulation of hydraulic fractures. International Journal of Rock Mechanics and Mining Sciences, (2007). ,

[2] Belytschko, T, & Black, T. Elastic crack growth in finite elements with minimal re‐ meshing. International Journal for Numerical Methods in Engineering, (1999). , 45,

[3] Desroches, J, Lenoach, B, Papanastasiou, P, & Thiercelin, M. On the modelling of near tip processes in hydraulic fractures. International Journal of Rock Mechanics

[4] Fish, J, & Belytschko, T. A First Course in Finite Elements, chapter 9, John Wiley &

[5] Fries, T. P, & Baydoun, M. Crack propagation with the extended finite element meth‐ od and a hybrid explicit-implicit crack description. International Journal for Numeri‐

[6] Fries, T. P, & Belytschko, T. The extended/generalized finite element method: An overview of the method and its applications. International Journal for Numerical

[7] Garagash, D. Plane-strain propagation of a fluid-driven fracture during injection and shut-in: Asymptotics of large toughness. Engineering Fracture Mechanics, (2006). ,

[8] Garagash, D. Propagation of a plane-strain hydraulic fracture with a fluid lag: Earlytime solution. International Journal of Solids and Structures, 43(18-19): 5811-5835,

and Mining Sciences & Geomechanics Abstracts, (1993). , 30(7), 1127-1134.

cal Methods in Engineering, (2012). , 89(12), 1527-1558.

Methods in Engineering, (2010). , 84(3), 253-304.

2 Institute of Geotechnical Engineering, RWTH Aachen University, Aachen, Germany

, C. Clauser3

The XFEM With An Explicit-Implicit Crack Description For Hydraulic Fracture Problems

and T. P. Fries1

http://dx.doi.org/10.5772/56383

723

N. Weber1

chen, Germany

**References**

**Figure 8.** The pressure distribution Π(ξ) (a) and the crack opening profile Ω(ξ) (b) for various values of the fluid fraction ξ*f* .

The pressure distribution and the crack opening profile along the dimensionless coordinate *ξ* = *x* / *L* are shown in Figures 8(a) and (b) in the viscosity scaling for fluid fraction values ξ<sup>f</sup> = {0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.97, 0.99}.

### **6. Conclusions**

The XFEM with an explicit-implicit crack description has been applied to a plane-strain hydraulic fracture problem. The crack is described explicitly by a line (2D)/triangular (3D) mesh that is aligned with the interface and implicitly by three level-set functions. The enrich‐ ment functions at the tip can be chosen according to the asymptotic behavior of the hydraulic fracture problem. Depending on the propagation regime the stress singularity can be described either by LEFM or by a singularity, which is weaker than predicted by LEFM. However, in this work a partially filled crack with a significant lag is examined and, therefore, crack propagation is governed by LEFM. The results show a good agreement with the known similarity solutions and can be interpreted as an early-time solution that can be used as a starting point in hydraulic fracture simulations.

### **Acknowledgements**

We acknowledge support for this work from the Federal Ministry for the Environment, Nature Conservation and Nuclear Safety, Germany (FKZ 0325167).
