FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 1 of 35

FRACTURE MODELING
1 Introductory Summary............................................................................................................. 2
2 Concepts................................................................................................................................... 3
2.1 Fundamental Laws............................................................................................................... 4
2.2 Constitutive Laws ................................................................................................................. 4
2.3 Fracture Propagation ........................................................................................................... 6
3 Hydraulic Fracturing Models .................................................................................................. 9
3.1 Two-Dimensional (2D) ....................................................................................................... 11
3.2 Pseudo Three-Dimensional (P-3D) .................................................................................... 15
3.3 Planar Three-Dimensional (PL-3D).................................................................................... 18
3.4 Fully Three-Dimensional (3D) ............................................................................................ 19
4 Examples ................................................................................................................................ 20
4.1 Case History ...................................................................................................................... 20
4.2 Model Comparisons ........................................................................................................... 28
FIGURES
Fig. 1. Modes of loading............................................................................................................... 7
Fig. 2. Fracture divided into elements.......................................................................................... 9
Fig. 3. Representation of a planar fracture. ............................................................................... 10
Fig. 4. KGD geometry. ............................................................................................................... 11
Fig. 5. PKN geometry................................................................................................................. 12
Fig. 6. 2D and radial Sneddon cracks. ....................................................................................... 13
Fig. 7. Elliptical profile (P-3D)..................................................................................................... 17
Fig. 8. Example grid (PL-3D model)........................................................................................... 18
Fig. 9. Fracture profile (PL-3D model). ...................................................................................... 19
Fig. 10. Permeability, thickness and stress profile. .................................................................... 20
Fig. 11. Computed values for Young's modulus and Poisson's ratio. ........................................ 21
Fig. 12. Profile of bottomhole, casing and tubing pressures. ..................................................... 24
Fig. 13. Pressure match for bottomhole and casing pressure. .................................................. 24
Fig. 14. Fracture profile.............................................................................................................. 25
Fig. 15. Fracture width profile. ................................................................................................... 25
Fig. 16. Match of net pressure for calibration fracture and main fracture. ................................. 26
Fig. 17. Fracture profile.............................................................................................................. 26
Fig. 18. Reservoir model for final history match......................................................................... 28
DOWELL CONFIDENTIAL

Section 300
May 1998

FRACTURING ENGINEERING MANUAL

Page 2 of 35
Fig. 19.
Fig. 20.
Fig. 21.
Fig. 22.
Fig. 23.
Fig. 24.
Fig. 25.
Fig. 26.

Fracture Modeling

Schlumberger

Dowell

TRIFRAC length and width profile. ................................................................................32

STIMPLAN length and width. ........................................................................................32
FRACPRO length and width profile. ..............................................................................33
GOHFER length and width profile. ................................................................................33
TERRAFRAC length profile. ..........................................................................................34
STIMPLAN length and width profile...............................................................................34
MEYER length and width profile. ...................................................................................35
Ohio state length profile.................................................................................................35
TABLES

Table 1.
Table 2.
Table 3.
Table 4.
Table 5.

Comparison Of Stress ..................................................................................................22

Permeability and Fluid Loss .........................................................................................22
Design Information .......................................................................................................23
Fracture Model Comparison Runs................................................................................30
Fracture Model Comparison Runs................................................................................31

1 Introductory Summary
The prediction of fracture geometry has been a central issue in engineering design
and evaluation of hydraulic fractures, and many models have been developed over
the years. These models determine fracture geometry by attempting to relate many
variables such as rock properties, fluid properties, fluid volume pumped and stress
data. Some models use a fixed fracture height and others continuously calculate the
height during the simulation. Each change aimed at more closely matching the real
conditions requires more sophistication in modeling the fluid flow in the entire
fracture, effect of proppant and elasticity of the entire system. To be practical,
however, the calculations must be made at reasonable increments along the fracture
and computational time must not be excessive. The degree of sophistication of a
model is therefore somewhat controlled by the practical application. The models are
also data limited.
The comparison of different models can be difficult and confusing because of the
way the various authors handle the variety of conditions, what they feel is important,
what assumptions they make and how portions of the model are coupled. Decisions
on how to handle elasticity, fluid flow, type of grid or cross section, vertical stress
differences and toughness, for example, can have a large impact on the calculated
fracture geometry. There is still much work to be done in obtaining meaningful data
for input into the more sophisticated models.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 3 of 35

Basically, there are four types of fracture models either being used or being
developed in the industry today.
Two-Dimensional (2D)

Pseudo Three-Dimensional (P-3D)

Planar Three-Dimensional (PL-3D)

Fully Three-Dimensional (3D).

The designation of the third listed model as Planar is an area of confusion that
should be clarified. Actually, the term planar means that the fracture occurs in a
plane. This condition is true for all fracture models except the fully 3D. Planar was
simply used to name the model that is more advanced than the pseudo 3D, but not
quite as sophisticated as the fully 3D. A fully 3D model would have the capability of
being nonplanar (fracture could curve or change planes) if the correct stress data
and other information were available for input.
This section on fracture models will be limited to a brief discussion of some of the
concepts that must be considered to build a model, as well as a brief discussion of
each model. It is beyond the scope of this section to cover each model in detail
because of the number of models available, and not having the documentation or
code to examine each model. Also, models change as more data become available
from evaluation, in-situ testing and calibration of logs, and from special industry
projects to calibrate the various models based on the best available information. A
successful model is one that has the ability to match the pressure from the treatment
by using realistic variables based on in-situ data, and to calculate a fracture height
consistent with other methods used on the actual treatment.

2 Concepts
Modeling fracturing treatments requires a blending of many different components,
such as rock mechanics, fluid mechanics, rheology and heat transfer. Two sets of
laws are required for this process.
The Fundamental Laws dealing with mass, momentum and energy conservation.
These relate to the physical principles.

The Constitutive Laws include rock elasticity and fluid rheology. These describe
the behavior of a system under a certain number of conditions.

Coupling these two sets with the appropriate boundary conditions produces some
very complicated mathematical formulations. To solve the coupled problem requires
discretization of the system (break into small geometric components such as a grid),
and then writing equations in a form that can be solved with digital computing.

DOWELL CONFIDENTIAL

Section 300

FRACTURING ENGINEERING MANUAL

May 1998

Fracture Modeling

Page 4 of 35

Schlumberger

Dowell

2.1 Fundamental Laws

Conservation of Mass: The mass of a system does not change with time. The
conservation of mass for modeling a fracture treatment is used to give the overall
material balance. This material balance is illustrated by Eq. 1.
Vi = V fp + VLP

(1)

Where:
Vi = volume injected
Vfp = fracture volume

VLp = injected volume lost to the formation.

Conservation of Momentum: Two types of forces can be distinguished; (1) the
body forces such as gravity that act on the whole volume, and (2) the surface
forces, such as pressure forces and fluid friction, that act only along the boundary
of the domain.
The conservation of momentum principle relates the time change in the total
momentum of a body to the applied forces (both on the volume and on the
surface). Many fracture treatments are modeled as quasistatic. This implies that
the rates of the change of velocities are negligible, and therefore the summation
of surface and volume forces is zero.
Conservation of Energy: This pertains primarily to the two fundamental laws of
thermodynamics for a system; (1) the change in total energy of a system is equal
to the work of the forces applied on the system plus the rate of heat transfer, and
(2) the internal energy of a system is a function only of its entropy.

2.2 Constitutive Laws

A mechanical system is completely defined by a certain number of variables that
depend on time and position.
mass

temperature

velocity (three components)

stresses (six components).

These variables make a minimum of 11 unknowns while the number of components

in the three conservation laws is only five. Therefore, constitutive equations to solve
the system need to be defined. Some equations used in the less complicated
models can be simplified to the point where they no longer apply, thereby reducing
the number of computations needed to complete a simulation.
Three types of constitutive relations are considered for the system.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 5 of 35

1. Fluid incompressibility is assumed.

2. A relation exists between stress and strain (or rates of strain).
3. Scalar (magnitude only) quantities are related to fluxes.
Incompressibility: For many fluids, the density does not change with time,
pressure or temperature. This assumption of constant density (incompressible
fluid) causes the equation of mass conservation to simplify and become zero.

Stress/Strain: The section on rock mechanics showed the relationship between

stress () and strain () for an elastic solid, giving the coefficient of proportionality
( = E) where E is Young's modulus. Modulus contrasts act mainly to alter the
shape of the fracture rather than the tip positions at any pressure. Analyse to
determine the effect of modulus contrast are extremely time consuming because
they require large finite-element or similar solutions. Formulas for the average
width are used in many cases to scale the width profile.
For the general case of a fracture growing into the surrounding layers of greater
stress and/or fracture toughness than the pay zone, the pressure, height and
shape of the fracture cross section depend on the stress, toughness, thickness
and modulus of the individual layers. The more simple models would assume a
constant modulus in all layers, plane strain elasticity and no flow-induced
pressure drop in the vertical direction.
When considering fluid flow, the stress tensor includes pressure and a viscous
tensor (). Newton's law for viscous incompressible fluids is simplified to
=

(2)

where is the viscosity of the fluid and is the strain rate tensor. Flow in several
of the fracture models is assumed to be in one direction only (x direction) and the
velocity field is therefore unidirectional. This means that both the shear stress
(xy) and rate of deformation:
xy are related as:
xy =
xy

(3)

Coupling conservation relations with the constitutive equation leads to a very

complicated fluid equation. The more sophisticated PL-3D models may attempt
to handle this equation while the 2D and many P-3D models must make
simplifying assumptions, particularly on velocity, before numerically solving the
system. Handling fluid flow in the fracture to determine the pressure distribution
is very important and necessary to determine the fracture displacement, as well
as using a method to determine how a particular model will handle all of the
variables associated with adding the proppant to the system. An accurate
calculation of the pressure distribution within the fracture is critical for any realistic
simulation.

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 6 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Schlumberger

Dowell

Flux Laws: These include the relationships of pressure drop in a porous medium
to velocity (Darcy's law), of heat flux to temperature (Fourier's law) or the rate of
reaction to the change in concentration (Fick's law) for acid fracturing.

2.3 Fracture Propagation

Modeling the fracturing process is a very complicated task. Injecting fluid into the
formation will modify the stress distribution and pressure within the formation,
creating a fracture in which the injected fluid flows. This injected fluid exchanges
heat with the formation and leaks off into the formation. The poroelastic effect can
be interpreted as a time-dependent back-stress; whenever a fracture profile is
obtained and a pressure is calculated, an additional term must be added to the
pressure. The fluid is also considered to be a multiphase type because of the
proppant being carried. The proppant also alters the fluid viscosity, which affects the
fluid-flow model. A good model needs rigorous calculation of the fluid fronts as well
as handling the settling of the proppant as in the simpler 2D models up to convection
hindered settling and other more complicated accounting of the proppant in some of
the PL-3D models.
Simulating the propagation of a hydraulic fracture requires consideration of the linear
elastic fracture mechanics (LEFM) formulation, fracture fluid flow and continuity
equation. These three sets need to be coupled to simulate the propagation; the
mathematical problem is complex because of the different types of equations and the
presence of a moving boundary (fracture edge).
The strains produced in the formation caused by the deformation from inducing a
hydraulic fracture are actually quite small. The small value of strain allows the
assumption that the formation deforms in a linear elastic manner. Knowing the state
of stress (or pressure) induced by the fluid in the fracture and the confining stresses
(boundary conditions) allows the calculation of the fracture width (displacement).
One simple method is to consider the crack as a uniformly pressurized elliptical
surface with semiaxes a and b. This condition forces the fracture to assume an
elliptical shape (Sneddon). Different cases may be derived from Sneddon's solution.
These cases can relate to the Griffith crack as well as the radial or penny-shaped
crack. This will be covered in more detail in the discussion of 2D models. The
classical techniques in elasticity do have severe limitations because of the stress
singularity near the crack tip, and special methods have to be used to improve the
accuracy in computing the stress intensity factors. Boundary integral techniques
have become popular because the problem is solved only at the boundary, reducing
the dimension of the problem by one, and provides a simple method to determine
stress intensity factors.
However, implementation can be difficult for
nonhomogeneous media with varying elastic moduli.
Remember from the section on rock mechanics that the stress intensity factor (KI) is
for the opening mode of the fracture and describes the magnitude of the stresses
near the crack tip it depends on the crack geometry, surrounding material and
applied loads. Fig. 1 shows the three fundamental modes of loading for a fracture.
DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 7 of 35

Fig. 1. Modes of loading.

The opening mode corresponds to normal tension, in front of the crack. The sliding
mode is associated with transverse shear while the tearing mode is a longitudinal
shear. Usually, only the opening mode applies for a plane fracture. The other
modes may be applicable in complex fracturing, out of plane fractures or when near
natural fissures.
For a given system of stresses, it is necessary to determine whether a fracture will
initiate and then propagate, and in which direction, and from which point in the
created fracture. Griffith (1921) addressed the problem by defining the surface
energy of a fracture. The Griffith criterion stated that there is an equilibrium if the
specific surface energy is related to the change in internal energy (U) by:
2 =

U
A

(4)

where A is the fracture area. For a uniformly loaded fracture, the critical load (c)
and the surface energy are related, such as

2 E

c =
x ( 1 2 )
f

1/ 2

(5)

where xf is the fracture half-length and E and are the elastic constants.

DOWELL CONFIDENTIAL

Section 300
May 1998

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Page 8 of 35

Schlumberger

Dowell

Another important concept related to the surface energy criterion is that of the critical
stress-intensity factor (KIc), also known as fracture toughness. The condition for a
fracture to be in equilibrium requires that the stress intensity factor (KI) associated
with the load be equal to a critical value (KIc). A simple relation can be derived
between the surface energy and KIc for a uniformly pressurized fracture
=

(1 2 )K 2 Ic
E

(6)

This expression indicates (for linear elastic behavior) that the surface energy
criterion and the critical intensity factor are related and form a unified criterion for
propagation.
The calculation of the pressure distribution in the fracture, due to fluid flow, is
necessary to determine the fracture displacements. As mentioned earlier, the fluidflow problems to include non-Newtonian as well as Newtonian fluids present some
problems and require a large amount of computational time. Because of this
problem, the fracture may be discretized into a series of parallel lines (surfaces) and
the flow considered as quasistatic. The PL-3D models do handle the flow problem a
little more rigorously and, consequently, require more time for computation.
The continuity equation is the last relation and simply describes the conservation of
mass previously discussed. The continuity equation may be written for each fracture
element (in the grid) Flow Rate In = Flow Rate Out + Accumulation.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 9 of 35

3 Hydraulic Fracturing Models

When a fracture is created in a rock mass, the fluid applies pressure to the two faces
of the resulting crack. This pressure causes the crack to open. The width at any
point is therefore dependent on the pressure applied to the faces. The condition of
the plane strain and rectangular elements is used to model the fracture for the 2D
cases. Fig. 2 illustrates a boundary element technique applied to a planar fracture.

Fig. 2. Fracture divided into elements.

The simplest representation is when the width is assumed to be piecewise constant
over each element. For problems solved by this method, the width distribution is
typically obtained for a specified distribution of pressure over the elements.
Alternatively, the pressure distribution can easily be found from a specified width
profile. Some P-3D models may use a little more sophisticated grid system, but fluid
flow is again only 1D and the vertical height growth is handled by shape factors and
other methods although elasticity is considered to be 2D. The PL-3D models use a
more complicated grid system (small squares or triangles). The junction of the grid
elements becomes the node for the solution of the problem. This allows the applied
pressure and displacements everywhere on the fracture faces to be related.
Elasticity is also 3D in these PL-3D models.

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 10 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Schlumberger

Dowell

The planar fracture was briefly introduced in the introductory summary to clarify a
model description. Fractures are planar features as shown in Fig. 3.

Fig. 3. Representation of a planar fracture.

The fracture lies in the xz plane, and the width is represented by the y direction.
Fluid may flow in all directions in this fracture, but 3D fluid flow is very complicated
and various models handle this flow differently. Usually flow along the fracture
length is modeled in detail, while horizontal flow across the width of the fracture and
vertical flow are neglected (2D, P-3D). Horizontal flow may be calculated in an
average sense by assuming that an average velocity adequately represents the flow.
An overall fluid-flow model could therefore be obtained by coupling 1D mass
conservation and momentum conservation equations. Volume conservation would
replace mass conservation if only incompressible fluids are considered.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 11 of 35

3.1 Two-Dimensional (2D)

Most 2D models require that a value for fracture height be input so the length and
width can be calculated from the volume and flow characteristics governed by the
code for the particular model. The models that will be discussed here are the (1)
Khristianovic and Zheltov, with later contributions by Geertsma and de Klerk (KGD
model), (2) Perkins and Kern, and later Nordgren (PKN model), and (3) the radial or
penny-shaped model.
A common simplifying assumption is that the lateral effects of a fracture are small
compared to the vertical effects and can be neglected. This condition is termed
plane strain and implies that each cross section acts independently of any other
section, so that the mechanical analysis need only be performed in two dimensions.
This plane-strain assumption is an integral part of the 2D fracture models. Two
cases of plane strain can be distinguished for the 2D formulations.
Horizontal Plane-Strain Geometry: The fracture zone will deform independently
of the upper and lower layers. This will be possible for free slippage on these
layers and represent a fracture with horizontal penetration that is much smaller
than the vertical, and where the fracture shape does not depend on the vertical
position. This describes the KGD fracture geometry as shown in Fig. 4. This
geometry has a constant height with a cross section that is rectangular in shape.

Fig. 4. KGD geometry.

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 12 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Schlumberger

Dowell

Vertical Plane-Strain Geometry: This condition exists for a large confinement

where the fracture is limited to a given zone. Each vertical cross section deforms
independently of the others. The fracture widths in the vertical direction are
coupled through the continuity and fluid-flow equations; however, because there
is no vertical extension (in each of the vertical sections) during simulation, the
pressure is uniform and the cross-sectional shape of the fracture height is
elliptical. This describes the PKN geometry as illustrated in Fig. 5.

Fig. 5. PKN geometry.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 13 of 35

Sneddon's solution for modeling the behavior of a fracture (for linear elastic
assumption) was mentioned earlier during the discussion of concepts. Fig. 6 shows
the 2D and radial Sneddon cracks. These solutions are for a 2D crack having one
dimension of infinite extent, and the other dimension of finite extent (d in Fig. 6). The
radial or penny-shaped crack is defined by the radius (R). The resulting width is
elliptical in shape for both types of cracks, and is proportional to one of the
characteristic dimensions (either d or R).

Fig. 6. 2D and radial Sneddon cracks.

The width is also proportional to the net pressure (pf - min) and inversely proportional
to the plane-strain modulus given as
E = E / (1 2 )

(7)

Young's modulus and Poisson's ratio are at in-situ conditions (E is defined by the
tangent Young's modulus Et).
Sneddon's method has been used in different ways to model 2D fractures. The
characteristic dimension, d, is assumed to be the total tip-to-tip fracture length (2xf)
for the KGD model. Since d is assumed to be the total fracture length, then the
infinite dimension has to correspond to the fracture height. The other assumption is
that the characteristic dimension, d, is the fracture height. This is the condition for
the PKN model; since d is the fracture height, then the infinite dimension has to be
the fracture length.

DOWELL CONFIDENTIAL

Section 300

FRACTURING ENGINEERING MANUAL

May 1998
Page 14 of 35

Fracture Modeling

Schlumberger

Dowell

The KGD and PKN models consider the propagation of a vertical fracture of given
height. However, in some cases, the vertical stress is lower than the horizontal
stress and the fracture will propagate in a horizontal or inclined direction. Both
Perkins and Kern as well as Geertsma and de Klerk gave simplified expressions for
the propagation of a radial horizontal fracture. Similarly, vertical fractures may
propagate radially in thick formations where there are no barriers to height growth.
This situation leads to the same equations when the injection interval is small
compared to the fracture extension. The KGD model is valid when h>>xf. The PKN
model is valid when xf >>h. The radial model is most appropriate when the total
length (2xf or 2R) is approximately equal to the height. Again, the three sets of
equations to be coupled are the elasticity, continuity and fluid flow.
Some of the important characteristics and differences pertaining to the 2D models
are
KGD Fracture Model
A fixed fracture height is assumed, and fluid flow is horizontal only (in the
direction of the propagation).
Crack opening is solved in the horizontal plane. As a result, the fracture
width does not vary with the fracture height, except by the boundary
condition set at the wellbore that specifies a constant total injection rate.
Width is constant in the vertical direction because of the plane-strain
condition and individual horizontal planes.
The model gives wider fracture widths and shorter fracture lengths when
compared to the PKN model.
The flow resistance in the narrow rectangular vertical width is what
determines the fluid pressure gradient in the propagating direction.
The excess pressure (net pressure) decreases with time, and in log-log
coordinates has a slope equal to -1/3.

The model is most appropriate when the fracture length is smaller than the
fracture height.
PKN Fracture Model
A fixed fracture height is assumed and fluid flow is horizontal only (in the
direction of the propagation).
Crack opening is solved in the vertical plane. As a result, the fracturing
fluid pressure is constant in vertical cross sections perpendicular to the
direction of the fracture propagation.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 15 of 35

Each vertical cross section deforms individually and obtains an elliptic

shape with the maximum width in the center. The only coupling between
the different vertical cross sections is due to fluid flow in the fracture. There
is no elastic coupling between the planes.
This model is used to describe the behavior of planar fractures that have a
length-to-height ratio greater than three.
The model gives narrower widths and longer fracture lengths compared to
the KGD model.
The excess pressure increases with times; in log-log coordinates, the slope
is positive.
The model is most appropriate when the fracture length is much larger than
the fracture height.
The 2D fracture models require less input and perform the computations very
quickly. This advantage may be used when the fracture is known to be contained in
height, and the value of the height is known. However, the fracture shape is still
rectangular because no difference is calculated in the height at the tip of the fracture
compared to the height at the wellbore. Because of these conditions as well as the
1D fluid flow and noncoupling of the length and width, the 2D model fracture length is
usually longer, the net pressure is higher and the width is wider than actual or when
compared to the various P-3D or PL-3D models.
3.2 Pseudo Three-Dimensional (P-3D)
Most of the P-3D models incorporate basic assumptions regarding the elastic
properties of the rock layers, fluid flow and fracture initiation. These assumptions
reduce the complexity and number of calculations needed for a fracture simulation.
The elasticity is considered in two dimensions (but not in cross section) and fluid flow
is one-dimensional in the direction of the propagation. Many P-3D simulators
assume all layers have the same elastic, reservoir and fluid-loss properties. Some
P-3D simulators allow the user to input a description for each layer and then use
some method for averaging, while others use calculated shape factors to arrive at
the effect the layers will have on the geometry. For example, a different Young's
modulus could be input for 15 layers. The simulator, depending on the choice of the
author, may then use averaging or other factors to determine the actual impact over
the entire interval. The values of 15 different layers are now reduced to one. This
method can have significant implications in sections with thick layers of coal, shale or
other lithology. Most P-3D simulators have problems running when a large contrast
in elastic properties exists. The assumptions and use of averaging and factors allow
the P-3D simulators to generate a solution in a very short time, even on the small
personal computers.
The P-3D simulators divide the fracture into a fixed number of vertical elements in
the grid. Each vertical element extends from the top of the fracture to the bottom.
DOWELL CONFIDENTIAL

Section 300
May 1998
Page 16 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Schlumberger

Dowell

All of the values describing the fracture geometry, pressure in the fracture and fluid
loss are the same throughout that one single element. This is one reason why most
P-3D simulators estimate more vertical height growth when compared to the PL-3D
simulators. Also, there is no pressure drop within the element from the perforations
to the vertical limits of the fracture. The pressure at the vertical tip of a P-3D fracture
simulation is also higher than that of a PL-3D simulation. The P-3D simulators may
produce solutions with excessive vertical height growth; however, the simpler
discretization scheme does significantly reduce the complexity and computational
time required to run a simulation.
Most P-3D simulators handle the fluid flow only in the horizontal direction within the
fracture. This fluid flow in the horizontal direction (the direction of the propagation) is
because of the single grid element in the vertical direction. This 1D flow also limits
the ability to properly describe the proppant transport within the fracture. To
compensate for the problem of proppant transport, some P-3D models use a
correlation that creates an elliptical proppant front.
P-3D simulations usually produce an elliptical fracture profile (side view) because of
the previously discussed methods for handling vertical growth, and because the P3D models are usually lumped models. Lumped models are those where the grid
data are integrated from the tip to the wellbore and averaged. This method typically
produces the elliptical fracture profile as illustrated in Fig. 7.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Fig. 7. Elliptical profile (P-3D).

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 17 of 35

Section 300
May 1998

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Page 18 of 35

Schlumberger

Dowell

3.3 Planar Three-Dimensional (PL-3D)

PL-3D models are able to provide a detailed description of a hydraulic fracture in
both simple and very heterogeneous lithologic sections. These models require all of
the parameters describing the reservoir, elastic and fluid-loss properties for each
layer. Simulation can take 15 min or several hours depending on the complexity of
the problem and/or the computer being used.
PL-3D models divide the fracture into a mesh of vertical and horizontal grid
elements. The elements may take different shapes (usually small squares or
triangles). Fig. 8 shows an example of a mesh composed of quadrilateral elements
with triangular subelements.

Fig. 8. Example grid (PL-3D model).

A very fine mesh may be used around the edges and at the tip for greater accuracy,
while a coarser grid is used in the center of the fracture. The fracture is assumed to
develop as a plane (elasticity equations relate the pressure in the elements on the
crack faces to the width) and to calculate fracture-height growth. The different
values of pressure, fluid loss and other important parameters are calculated both
vertically and horizontally. The PL-3D model has 3D elasticity and 2D fluid flow.
Fracture (crack) propagation is controlled by the criterion of linear elastic fracture
mechanics. The fracture advances by a method where the stress intensity factor (KI)
is maintained nearly equal to the critical stress-intensity factor (KIc) during the crack
extension at each grid node. The theory of linear elasticity also allows problems to
be solved by the superposition of solutions to subproblems. Superposition is a
method whereby a PL-3D model can use part of Sneddon's solution to calculate the
vertical width profile.
Most PL-3D simulators allow the user to control the refinement of the grid. A coarse
grid with fewer elements will require less computing time, but there may be more
error in the solution. A very fine mesh grid will provide a solution with a very low
material balance error. The more complex grid system also allows the description of
the proppant transport in a more detailed manner.
DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 19 of 35

The 2D grid can produce some problems with complex stress profiles. Some PL-3D
simulators have difficulty handling bounding layers with a value of stress lower than
the pay, which are separated from the pay by a bounding layer with a higher stress.
This situation causes the grid elements to become increasingly skewed as the
fracture grows into the lower stress bounding layer. This continues until the
computational errors are extremely large and the simulator terminates the run.
The advantage of any PL-3D simulator is the ability to model the pressure drop
laterally within the fracture. As the fracture grows with each step of the simulation,
all of the parameters are recalculated in each grid element. Using the grid system,
the pressure at the lateral tip of the fracture can be lower than the pressure at the
wellbore. This allows the simulator to model a fracture with a greater vertical height
at the wellbore than at the fracture tip. The profile of this type of fracture is
determined by the grid and is not a lumped solution. Fig. 9 is a profile from a PL-3D
simulation.

Fig. 9. Fracture profile (PL-3D model).

3.4 Fully Three-Dimensional (3D)
Techniques for general 3D fracture propagation (including out-of-plane) have been
presented by some authors. The numerous problems and the system of equations
needed have not been developed to the point of providing a realistic, operational,
fully 3D model. A model that is truly a fully 3D model is needed to simulate special
conditions that cannot be handled by other models. However, at the present time,
no working model exists in the industry except in research.
DOWELL CONFIDENTIAL

Section 300
May 1998
Page 20 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Schlumberger

Dowell

4 Examples
The information used to compile the examples for the different models comes from a
case study where several calibration tests were performed. Extensive stress testing,
core analysis, reservoir evaluation and fluid analysis were also performed. The
reservoir parameters, stress values and other critical parameters are considered to
be the best available because of the many techniques used to obtain and
corroborate the data. Routine treatment designs usually never have such a large
volume of reliable data available with which to work.
Trying to examine all of the input and output data in detail for the examples, and for
each model, is simply too voluminous to cover in this section. Tables and profiles
will be used to show data and comparisons (case history). The actual model
comparisons will show profiles and an output summary for each model simulation.
4.1 Case History
This case history uses the same basic data that will be used later for the model
comparison. However, this case history uses 15 to 22 layers in the simulation,
whereas only 5 layers were used for the model comparisons given in Subsection 4.2
that follows this case history.
The reservoir stress and elastic properties that are used in the examples are given in
Fig. 10 and Fig. 11.

Fig. 10. Permeability, thickness and stress profile.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 21 of 35

Fig. 11. Computed values for Young's modulus and Poisson's ratio.

DOWELL CONFIDENTIAL

Section 300

FRACTURING ENGINEERING MANUAL

May 1998

Schlumberger

Fracture Modeling

Page 22 of 35

Dowell

Table 1 shows the depth, log stress and modified stress that can be used for layer
input. The modified stresses were obtained by modifying the log stresses to more
closely match the measured net pressure response from the calibration fractures.
These modified stress data were the values used for the model simulations.

Table 1. Comparison of Stress

Top of Zone
(ft)

Log Stress
(psi)

Modified Stress
(psi)

Delta Stress
(psi)

9030
9070
9115
9155
9170
9200
9250
9310
9340
9360
9380
9435
9455
9475
9575

7300
7800
7150
6600
6050
5600
5250
5850
6550
7300
5800
6400
7550
8400
7850

7300
8200
7350
6600
6050
5800
5250
6050
6550
7300
6200
6700
7950
8400
7850

0
400
200
0
0
200
0
200
0
0
400
300
400
0
0

Table 2 gives the permeability and fluid-loss coefficient used in the initial simulations.
Table 3 shows the design information for the treatment and the actual volumes used
during the treatment. The average injection rate was 50 bbl/min and a total of
1,168,910 lbm of sand was placed in the fracture. The average treating pressure
was 3000 psi.

Table 2. Permeability and Fluid Loss

Formation Permeability (md)
Initial Total Leakoff Coefficient ( ft min )

0.0065
0.0010

Reservoir Fluid to Filtrate Permeability Ratio

Reservoir to Filter-Cake Permeability Ratio
Leakoff Interval

10.0
100.0
Entire Fracture

The treatment was pumped down the casing/tubing annulus, and the bottomhole
pressure was measured by a pressure gauge run inside the tubing. Fig. 12 shows
the pressure profile.
DOWELL CONFIDENTIAL

Section 300

FRACTURING ENGINEERING MANUAL

Schlumberger

May 1998

Fracture Modeling

Dowell

Page 23 of 35

Note that the second calibration test performed on the well is also included on this
profile, and ends at a little over 100 min. A shut-in time is obtained and then the
main fracturing treatment starts at approximately 250 min. A PL-3D model
(GOHFER) from Marathon Oil Company was used to match the bottomhole pressure
as well as the surface casing (annulus) pressure. The simulation used stress and
elasticity data from the top of the log to the bottom. These data were used to give 22
layers for input. Fig. 13 shows the pressure match (dotted lines) from the start of the
treatment to the point of shutdown. Fig. 14 and Fig. 15 show the fracture length and
fracture width profiles from the simulation.

Table 3. Design Information

Fluid Volume, BBL
Stage

Design

Actual

Fluid Type

1
2
3
4
5

1500
0
1000
3000
300

1571
0
576
2908
293

40-lbm linear
Shut-in
Slickwater
50-lbm x-link
50-lbm x-link

Proppant
Conc.
(lbm/gal)
0
0
0
0
1

6
7

400
400

477
337

40-lbm x-link
40-lbm x-link

500

500

500

10

Proppant Volume,
lbm
Design
Actual

0
0
0
0
12,600

0
0
0
0
11,110

0
1

0
16,800

0
14,150

40-lbm x-link

42,000

41,120

437

40-lbm x-link

63,000

54,910

600

623

40-lbm x-link

100,800

105,330

11

800

755

40-lbm x-link

168,000

133,630

12

900

807

40-lbm x-link

226,800

227,250

13

900

934

40-lbm x-link

264,600

300,100

14

900

935

40-lbm x-link

302,400

292,420

15

300

240

Slickwater

DOWELL CONFIDENTIAL

Proppant
Type
None
None
None
None
100-Mesh
Sand
None
20/40
Ottawa sand
20/40
Ottawa sand
20/40
Ottawa sand
20/40
Ottawa sand
20/40
Ottawa sand
20/40
Ottawa sand
20/40
Ottawa sand
20/40
Ottawa sand
None

Section 300
May 1998
Page 24 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Fig. 12. Profile of bottomhole, casing and tubing pressures.

Fig. 13. Pressure match for bottomhole and casing pressure.

DOWELL CONFIDENTIAL

Schlumberger

Dowell

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Fig. 14. Fracture profile.

Fig. 15. Fracture width profile.

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 25 of 35

Section 300
May 1998
Page 26 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Schlumberger

Dowell

Net pressure was also matched using the P-3D Cleary model (FRACPRO) as run by
Resource Engineering Systems (RES). The net pressure match is shown in Fig. 16.

Fig. 16. Match of net pressure for calibration fracture and main fracture.
The drop in pressure at about 300 min, also shown on the bottomhole pressure plot,
corresponds to the time that the fracture height during the simulation reached a low
stress zone (see Stress Log Interval 9380 to 9455). Fig. 17 is the fracture profile
data from the simulation for the net pressure match.

Fig. 17. Fracture profile.

DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Calculated Parameters
Fracture Length (ft)

Section 300
May 1998
Page 27 of 35

1819

Propped Length (ft)

1441

Fracture Height (ft)

358

Propped Height (ft)

284

Fracture Width (max.) (in.)

0.72

Proppant Concentration (lbm/ft2)

1.88

These are the dimensions at the end of shut-in. Several methods were used to
attempt to determine the fracture height after the treatment. The method that was
accepted as being the most accurate in this case was the Continuous Microseismic
Radiation (CMR) log. This microseismic height log was run four months after the
treatment, and therefore will more accurately indicate the propped fracture height
rather than the created fracture height.
The CMR log indicated the height was from 9125 ft to 9375 ft, allowing 25 ft at both
top and bottom. This makes the propped height range from 250 ft to 300 ft. The
created fracture-height differences for the various models are difficult to evaluate
because created height can be significantly different from the propped height. This
difference can be caused by the way each model handles the information on the
various layers, and the actual width profile that was calculated.
Another
consideration is how the model treats the proppant movement, settling or the many
other complicated aspects present when the proppant is added to the system.
Several postfracture reservoir evaluation techniques were used to analyze the
results of the treatment. The analysis testing was started after 89 days of production
from the well. A reservoir model was used to obtain the final history match of the
reservoir variables. The reservoir model used data obtained from buildup analysis
as well as the production match for both gas and water (two-phase flow). The result
of this analysis is illustrated in Fig. 18.

DOWELL CONFIDENTIAL

Section 300

FRACTURING ENGINEERING MANUAL

May 1998

Fracture Modeling

Page 28 of 35

Schlumberger

Dowell

Fig. 18. Reservoir model for final history match.

Note that conductivity was calculated to 1515 ft, but possibly only 1100 ft may be
truly effective. Total fracture half-length is shown as 1845 ft.
It must be pointed out that FRACPRO was the initial model used for the project, and
continued to be used through the completion of the project. Comparison with other
models (to be shown later) is more difficult because many times the final geometry
interpretation may be influenced by special data or variables input to obtain the
match. This is not wrong, and is the way that models should be calibrated.
However, other models may not have the advantage of using the information pieces
obtained throughout the time of the project. Also, in the case history just discussed,
the data given to run the simulation with GOHFER differed some from the exact input
used for FRACPRO.
4.2 Model Comparisons
This portion will show comparisons for a similar set of input parameters. The critical
output parameters obtained from simulating a common treatment (using the same
data as in Table 2 and Table 3) will be given, as well as profiles or graphs illustrating
the resulting fracture geometry. These comparisons were made by S.A. Holditch &
Associates, Inc., using original data from a project for the Gas Research Institute.
The reader will see some significant differences in some of the output values, but
these serve to show again the differences between the model code. Different
companies running the same model also obtain different results. This difference may
be caused by how each company handled the stress layers, and the fact that in
some cases the various companies were asked to model the calibration fracture and
use fluid-loss coefficients derived from that analysis. This left room for some to use a
wall-building coefficient initially, while others may have started with a total fluid-loss
coefficient. Hopes are that sometime in the future, one inclusive set of data will be
given to everybody and the models run again for comparison.
DOWELL CONFIDENTIAL

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Section 300
May 1998
Page 29 of 35

The comparisons show results for KDG and PKN modes as well as for the P-3D and
PL-3D runs. Table 4 gives the list of output data from the runs using both GDK
(KGD)) and PKN geometry.
Note that since these models use fixed-height assumptions, most have the same
height listed because it is an input parameter. Comparisons are made for two cases
(1) fluid viscosity set at 200 cp, and (2) viscosity based on n' and k' values. These
values give a viscosity of approximately 450 cp at 37.5 sec -1. Case 3 and Case 4
show the difference between fluid-loss values, which is also evident by the calculated
efficiency. For the sake of space, example profiles for these KGD and PKN
geometries will not be shown.
The list of output data from P-3D and PL-3D simulations is given in Table 5. These
simulations use the same data previously shown in Table 2 and Table 3. The same
information concerning the viscosities applies here as it did from the previous
discussion.
The additional data and comparisons for these P-3D and PL-3D simulations involve
the number of layers used in the particular simulation. The example profiles and
graphs for these simulations will be taken from Case 8 of Table 5 (Variable Viscosity,
5-layer). The detailed differences will not be discussed, but left to the reader to
examine.
Table 5 should be self-explanatory for comparing the information in Case 8 with the
figures. However, to eliminate confusion each simulation will carry the figure number
corresponding to the following list
Fig. 19. SAH (TRIFRAC) length and width profile.

Fig. 20. NSI (STIMPLAN) length and width graphs (P-3D).

Fig. 21. RES (FRACPRO) length and width profile (P-3D).

Fig. 22. Marathon (GOHFER) length and width profile (PL-3D).

Fig. 23. ARCO (TERRAFRAC) length profile only (PL-3D).

Fig. 24. ARCO (STIMPLAN) length and width profile (P-3D).

Fig. 25. MEYER length and width profile P-3D).

Fig. 26. Ohio State length profile only (P-3D).

Copies of the simulation plots for length and width were not available for MEYER
(Bells) or TEXACO (FRACPRO).

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 30 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Table 4. Fracture Model Comparison Runs

DOWELL CONFIDENTIAL

Schlumberger

Dowell

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Table 5. Fracture Model Comparison Runs

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 31 of 35

Section 300
May 1998
Page 32 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Fig. 19. TRIFRAC length and width profile.

Fig. 20. STIMPLAN length and width.

DOWELL CONFIDENTIAL

Schlumberger

Dowell

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Fig. 21. FRACPRO length and width profile.

Fig. 22. GOHFER length and width profile.

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 33 of 35

Section 300
May 1998
Page 34 of 35

FRACTURING ENGINEERING MANUAL

Fracture Modeling

Fig. 23. TERRAFRAC length profile.

Fig. 24. STIMPLAN length and width profile.

DOWELL CONFIDENTIAL

Schlumberger

Dowell

FRACTURING ENGINEERING MANUAL

Schlumberger

Dowell

Fracture Modeling

Fig. 25. MEYER length and width profile.

Fig. 26. Ohio state length profile.

DOWELL CONFIDENTIAL

Section 300
May 1998
Page 35 of 35