Sie sind auf Seite 1von 11

Potential Contaminant Pathways from

Hydraulically Fractured Shale to Aquifers

by Tom Myers
various forms of water pollution. Two potential pathwaysadvective transport through bulk media and preferential
ow through fracturescould allow the transport of contaminants from the fractured shale to aquifers. There
is substantial geologic evidence that natural vertical ow drives contaminants, mostly brine, to near the surface
region, could reduce the travel time further. Injection of up to 15,000,000 L of uid into the shale generates
advects through the shale. The advection displaces native uids, mostly brine, and fractures the bulk media
requires from 3 to 6 years to reach a new equilibrium reecting the signicant changes caused by fracking the
wells have a reasonable offset from faults.
The use of natural gas (NG) in the United States has
been increasing, with 53% of new electricity generating
capacity between 2007 and 2030 projected to be with NG-
red plants (EIA 2009). Unconventional sources account
for a signicant proportion of the new NG available to
the plants. A specic unconventional source has been
deep shale-bed NG, including the Marcellus shale primar-
ily in New York, Pennsylvania, Ohio, and West Virginia
(Soeder 2010), which has seen over 4000 wells devel-
oped between 2009 and 2010 in Pennsylvania (Figure 1).
Unconventional shale-bed NG differs from conventional
Hydrologic Consultant, 6320 Walnut Creek Road, Reno,
NV 89523; (775) 530-1483; fax: (775) 530-1483; tom_
Received August 2011, accepted February 2012.
2012, The Author(s)
Ground Water 2012, National Ground Water Association.
doi: 10.1111/j.1745-6584.2012.00933.x
sources in that the host-formation permeability is so low
that gas does not naturally ow in timeframes suitable for
development. Hydraulic fracturing (fracking, the industry
term for the operation; Kramer 2011) loosens the forma-
tion to release the gas and provide pathways for it to move
to a well.
Fracking injects up to 17 million liters of uid
consisting of water and additives, including benzene at
concentrations up to 560 ppm (Jehn 2011), at pressures
up to 69,000 kPa (PADEP 2011) into low permeability
shale to force open and connect the fractures. This is
often done using horizontal drilling through the middle
of the shale with wells more than a kilometer long. The
amount of injected uid that returns to the ground surface
after fracking ranges from 9% to 34% of the injected uid
(Alleman 2011; NYDEC 2009), although some would be
formation water.
Many agency reports and legal citations (DiGiulio
et al. 2011; PADEP 2009; ODNR 2008) and peer-
reviewed articles (Osborn et al. 2011; White and Mathes
872 Vol. 50, No. 6GROUND WATERNovember-December 2012 (pages 872882)
2009 to June 2010 and total Marcellus shale wells in New York and West Virginia. There are 4064 wells shown in Pennsylvania,
48 wells in New York, and 1421 wells in West Virginia. Faulting in the area is documented by PBTGS (2001), Isachsen and
McKendree (1977), and WVGES (2011, 2010a, 2010b).
2006) have found more gas in water wells near areas
being developed for unconventional NG, documenting the
source can be difcult. One reason for the difculty is
the different sources; thermogenic gas is formed by com-
pression and heat at depth and bacteriogenic gas is formed
by bacteria breaking down organic material (Schoell
1980). The source can be distinguished based on both
C and H isotopes and the ratio of methane to higher chain
gases (Osborn and McIntosh 2010; Breen et al. 2007).
Thermogenic gas can reach aquifers only by leaking from
the well bore or by seeping vertically from the source.
In either case, the gas must ow through potentially very
thick sequences of sedimentary rock to reach the aquifers.
Many studies which have found thermogenic gas in water
wells found more gas near fracture zones (DiGiulio et al.
2011; Osborn et al. 2011; Breen et al. 2007), suggesting
that fractures are pathways for gas transport.
A pathway for gas would also be a pathway for u-
ids and contaminants to advect from the fractured shale to
the surface, although the transport time would be longer.
Fracking uid has been found in aquifers (DiGiulio et al.
2011; EPA 1987), although the exact source and pathways
had not been determined. With the increasing development
of unconventional NG sources, the risk to aquifers could
be increasing. With so little data concerning the movement
of contaminants along pathways from depth, either from
wellbores or from deep formations, to aquifers, conceptual
analyses are an alternative means to consider the risks.
The intent of this study is to characterize the risk
factors associated with vertical contaminant transport
from the shale to near-surface aquifers through natural
pathways. I consider rst the potential pathways for
contaminant transport through bedrock and the necessary
conditions for such transport to occur. Second, I estimate
contaminant travel times through the potential pathways,
with a bound on these estimates based on formation
hydrologic parameters, using interpretative MODFLOW-
2000 (Harbaugh et al. 2000) computations. The modeling
does not, and cannot, account for all of the complexities
of the geology, which could either increase or decrease
the travel times compared to those considered herein.
The article also does not include improperly abandoned T. Myers GROUND WATER 50, no. 6: 872882 873
boreholes which could cause rapid transport in addition
to natural pathways.
Method of Analysis
Using the Marcellus shale region of southern New
York (Figure 1), I consider several potential scenarios
of transport from shale, 1500 m below ground surface
(bgs) to the surface, beginning with pre-development
steady state conditions to establish a baseline and
then scenarios considering transport after fracking has
potentially caused contaminants to reach formations
above the shale. To develop the conceptual models and
MODFLOW-2000 simulations, it is necessary rst to
consider the hydrogeology of the shale and the details
of hydraulic fracturing, including details of how fracking
changes the shale hydrogeologic properties.
Hydrogeology of Marcellus Shale
Shale is a mudstone, a sedimentary rock consisting
primarily of clay- and silt-sized particles. It forms
through the deposition of ne particles in a low energy
environment, such as a lake- or seabed. The Marcellus
shale formed in very deep offshore conditions during
Devonian time (Harper 1999) where only the nest
particles had remained suspended. The depth to the
Marcellus shale varies to as much as 3000 m in parts
of Pennsylvania, and averages about 1500 m in southern
New York (Soeder 2010). Between the shale and the
ground surface are layers of sedimentary rock, including
sandstone, siltstone, and shale (NYDEC 2009).
Marcellus shale has very low natural intrinsic perme-
ability, on the order of 10
Darcies (Kwon et al. 2004a,
2004b; Neuzil 1986, 1994). Schulze-Makuch et al. (1999)
described Devonian shale of the Appalachian Basin, of
which the Marcellus is a major part, as containing coaly
organic material and appear either gray or black and
being composed mainly of tiny quartz grains <0.005 mm
diameter with sheets of thin clay akes. Median particle
size is 0.0069 0.00141 mm with a grain size distribu-
tion of <2% sand, 73% silt, and 25% clay. Primary pores
are typically 5 10
mm in diameter, matrix porosity
is typically 1% to 4.5% and fracture porosity is typically
7.8% to 9% (Schulze-Makuch et al. 1999 and references
Porous ow in unfractured shale is negligible due
to the low bulk media permeability, but at larger scales
fractures control and may allow signicant ow. The Mar-
cellus shale is fractured by faulting and contains synclines
and anticlines that cause tension cracks (Engelder et al.
2009; Nickelsen 1986). It is sufciently fractured in some
places to support water wells just 6 to 10 km from where
it is being developed for NG at 2000 m bgs (Loyd and
Carswell 1981). Conductivity scale dependency (Schulze-
Makuch et al. 1999) may be described as follows:
K = Cv
K is hydraulic conductivity (m/s), C is the intercept of a
log-log plot of observed K to scale (the K at a sample
volume of 1 m
), v is sample volume (m
), and m is
a scaling exponent determined with log-log regression;
for Devonian shale, C equals 10
, representing the
intercept, and m equals 1.08 (Schulze-Makuch et al.
1999). The very low intercept value is a statistical but
not geologic outlier because it corresponds with very
low permeability values and demonstrates the importance
of fracture ow in the system (Schulze-Makuch et al.
1999). Most of their 89 samples were small because the
deep shale is not easily tested at a eld-scale and no
groundwater models have been calibrated for ow through
the Marcellus shale. Considering a 1-km square area with
30-m thickness, the Kh would equal 5.96 10
(0.0515 m/d). This effective K is low and the shale would
be an aquitard, but a leaky one.
Contaminant Pathways from Shale to the Surface
Thermogenic NG found in near-surface water wells
(Osborn et al. 2011; Breen et al. 2007) demonstrates the
potential for vertical transport of gas from depth. Osborn
et al. (2011) found systematic circumstantial evidence for
higher methane concentrations in wells within 1 km of
Marcellus shale gas wells. Potential pathways include
advective transport through sedimentary rock, fractures
and faults, and abandoned wells or open boreholes. Gas
movement through fractures depends on fracture width
(Etiope and Martinelli 2002) and is a primary concern for
many projects, including carbon sequestration (Annunzi-
atellis et al. 2008) and NG storage (Breen et al. 2007).
Open boreholes and improperly sealed water and gas
wells can be highly conductive pathways among aquifers
(Lacombe et al. 1995; Silliman and Higgins 1990).
Pathways for gas suggest pathways for uids and
contaminants, if there is a gradient. Vertical hydraulic
gradients of a up to a few percent, or about 30 m over
1500 m, exist throughout the Marcellus shale region as
may be seen in various geothermal developments in
New York (TAL 1981). Brine more than a thousand
meters above their evaporite source (Dresel and Rose
2010) is evidence of upward movement from depth to
the surface. The Marcellus shale, with salinity as high
as 350,000 mg/L (Soeder 2010; NYDEC 2009), may
be a primary brine source. Relatively uniform brine
concentrations over large areas (Williams et al. 1998)
suggest widespread advective transport. The transition
from brine to freshwater suggests a long-term equilibrium
between the upward movement of brine and downward
movement of freshwater. Faults, which occur throughout
the Marcellus shale region (Figure 1) (Gold 1999), could
provide pathways (Konikow 2011; Caine et al. 1996)
for more concentrated advective and dispersive transport.
Brine concentrating in faults or anticline zones reects
potential preferential pathways (Wunsch 2011; Dresel and
Rose 2010; Williams 2010; Williams et al. 1998).
In addition to the natural gradient, buoyancy would
provide an additional initial upward push. At TDS equal
to 350,000 mg/L, the density at 25

C is approximately
874 T. Myers GROUND WATER 50, no. 6: 872882
1290 kg/m
, or more than 29% higher than freshwater.
The upward force would equal the difference in weight
between the injected uid and displaced brine. As an
example, if 10,000,000 L does not return to the surface as
owback (Jehn 2011), the difference in mass between the
volume of fracking uid and displaced brine is approxi-
mately 3,000,000 kg, which would cause an initial upward
force. The density difference would dissipate as the salt
concentration in the fracking uid increases due to diffu-
sion across the boundary between the uid and the brine.
In just Pennsylvania, more than 180,000 wells had
been drilled prior to any requirement for documenting
their location (Davies 2011), therefore the location of
many wells is unknown and some have probably been
improperly abandoned. These pathways connect aquifers
through otherwise continuous aquitards; overpressuriza-
tion of lower aquifers due to injection near the well
pathway could cause rapid transport to higher aquifers
(Lacombe et al. 1995). In the short fracking period, the
region that is overpressurized remains relatively close to
the gas well (see modeling analysis below), therefore it
should be possible for the driller to locate nearby aban-
doned wells that could be affected by fracking. This article
does not consider the potential contamination although
unlocated abandoned wells of all types must be considered
a potential and possibly faster source for contamination
due to fracking.
Effect of Hydraulic Fracturing on Shale
Fracking increases the permeability of the targeted
shale to make extraction of NG economically efcient
(Engelder et al. 2009; Arthur et al. 2008). Fracking
creates fracture pathways with up to 9.2 million square
meters of surface area in the shale accessible to a
horizontal well (King 2010; King et al. 2008) and
connects natural fractures (Engelder et al. 2009; King
et al. 2008). No post-fracking studies that documented
hydrologic properties were found while researching this
article (there is a lack of information about pre- and post-
fracking properties; Schweitzer and Bilgesu 2009), but
it is reasonable to assume the K increases signicantly
because of the newly created and widened fractures.
Fully developed shale typically has wells spaced at
about 300-m intervals (Edwards and Weisset 2011; Soeder
2010). Up to eight wells may be drilled from a single
well pad (NYDEC 2009; Arthur et al. 2008), although
not in a perfect spoke pattern. Reducing by half the
effective spacing did not enhance overall productivity
(Edwards and Weisset 2011) which indicates that 300-m
spacing creates sufcient overlap among fractured zones
to assure adequate gas drainage. The properties controlling
groundwater ow would therefore be affected over a large
area, not just at a single horizontal well or set of wells
emanating from a single well pad.
Fracking is not intended to affect surrounding forma-
tions, but shale properties vary over short ranges (King
2010; Boyer et al. 2006) and out-of-formation fracking is
not uncommon. In the Marcellus shale, out-of-formation
fracks have been documented 500 m above the top of the
shale (Fisher and Warpinski 2011). These fractures could
contact higher conductivity sandstone, natural fractures, or
unplugged abandoned wells above the target shale. Also,
uids could reach surrounding formations just because of
the volume injected into the shale, which must displace
natural uid, such as the existing brine in the shale.
Analysis of Potential Transport along Pathways
Fracking could cause contaminants to reach overlying
formations either by fracking out of formation, connecting
fractures in the shale to overlying bedrock, or by
simple displacement of uids from the shale into the
overburden. Advective transport, considered as simple
particle velocity, will manifest if there is a signicant
vertical component to the regional hydraulic gradient.
Numerical modeling, completed with the MODF
LOW-2000 code (Harbaugh et al. 2000), provides ex-
ibility to consider potential conceptual ow scenarios, but
should be considered interpretative (Hill and Tiedeman
2007). The simulation considers the rate of vertical trans-
port of contaminants to near the surface for the different
conceptual models, based on an expected, simplied, real-
istic range of hydrogeologic aquifer parameters.
MODFLOW-2000 is a versatile numerical modeling
code, but there is insufcient data regarding the geology
and water chemistry between aquifers and the deep shale,
such as salinity proles or data concerning mixing of the
brine with fracking uid, to best use its capabilities. As
more data becomes available, it may be useful to consider
simulating the added upward force caused by the brine by
using the SEAWAT-2000 module (Langevin et al. 2003).
Vertical ow would be perpendicular to the general
tendency for sedimentary layers to have higher horizontal
than vertical conductivity. Fractures and improperly
abandoned wells would provide pathways for much
quicker vertical transport than general advective transport.
This article considers the fractures as vertical columns
with model cells having much higher conductivity than
the surrounding bedrock. The cell discretization is ne, so
the simulated width of the fracture zones is realistic. Dual
porosity modeling (Shoemaker et al. 2008) is not justied
because turbulent vertical ow through the fractures is
unlikely, except possibly during the actual fracking that
causes out-of-formation fractures, a scenario not simulated
here. MODFLOW-2000 has a module, MNW (Halford
and Hanson 2002), that could simulate rapid transport
through open bore holes. MNW should be used in
situations where open boreholes or improperly abandoned
wells are known or postulated to exist.
The thickness of the formations and fault would affect
the simulation, but much less than the several-order-of-
magnitude variation possible in the shale properties. The
overburden and shale thickness were set equal to 1500 and
30 m, respectively, similar to that observed in southern
New York. The estimated travel times are proportional
for thicker or thinner sections. The overburden could
be predominantly sandstone, with sections of shale,
mudstone, and limestone. The vertical fault is assumed T. Myers GROUND WATER 50, no. 6: 872882 875
to be 6-m thick. The fault is an attempt at considering
fracture ow, but the simulation treats the 6-m wide fault
zone as homogeneous, which could underestimate the real
transport rate in fracture-controlled systems which could
be highly affected by dispersion. The simulation also
ignores diffusion between the fracture and the adjacent
shale matrix (Konikow 2011).
There are ve conceptual models of ow and trans-
port of natural and post-fracking transport from the level
of the Marcellus shale to the near-surface to consider
1. The natural upward advective ow due to a head
drop of 30 m from below the Marcellus shale to the
ground surface, considering the variability in both shale
and overburden K. This is a steady state solution for
upward advection through a 30-m thick shale zone
and 1500-m overburden. Table 1 shows the chosen K
values for shale and sandstone.
2. Same as number 1, but with a vertical fracture
connecting the shale with the surface, created using
a high-conductivity zone in a row of cells extending
through all from above the shale to the surface. This
emulates the conceptual model postulated for ow into
the alluvial aquifers near stream channels, the location
of which may be controlled by faults (Williams et al.
1998). The fault K varies from 10 to 1000 times the
surrounding bulk sandstone K (K
3. This scenario tests the effect of extensive fracturing
in the Marcellus shale by increasing the shale K
) from 10 to 1000 times its native value over
an extensive area. This transient solution starts with
initial conditions being a steady state solution from
scenario 1. The K
increases from 10 to 1000 times
at the beginning of the simulation, to represent the
relatively instantaneous change on the regional shale
hydrogeology imposed by the fracking. The simulation
estimates both the changes in ux and the time for the
system to reach equilibrium.
4. As number 3, considering the effect of the same
changes in shale properties but with a fault as in
number 2.
5. This scenario simulates the actual injection of 13 to
17 million liters of uid in 5 d into fractured shale
from a horizontal well with and without a fault.
Model Setup
The model domain was 150 rows and columns spaced
at 3 m to form a 450-m square (Figure 2) with 50 layers
bounded with no ow boundaries. The 30-m thick shale
was divided into 10 equal thickness layers from layer 40
to 49. The overburden layer thickness varied from 3 m
just above the shale to layer 34, 6 m from layer 33 to 29,
9 m from layer 28 to 26, 18 m in layer 25, 30 m from
layer 24 to 17, 60 m from layer 16 to 6, 90 m from layer
5 to 3, and 100 m in layers 2 and 1. A 6-m wide column
from layer 39 to the surface is added for some scenarios
in the center two rows to simulate a higher K fault.
Table 1
Sandstone (ss) and Shale (sh) Conductivity (K)
(m/d) and the Steady State Flux (m
/d) for Model
1 Scenarios
Flux K
1.7 0.1 0.00001
1.8 0.5 0.00001
1.9 1 0.00001
1.9 5 0.00001
2.0 10 0.00001
2.0 50 0.00001
2.0 100 0.00001
1.7 0.1 0.00001
9.5 0.1 0.00005
19.0 0.1 0.0001
81.2 0.1 0.0005
135.9 0.1 0.001
291.5 0.1 0.005
340.9 0.1 0.01
394.3 0.1 0.05
401.8 0.1 0.1
409.2 0.1 0.5
40.7 0.001 0.1
186.0 0.005 0.1
339.1 0.01 0.1
988.3 0.05 0.1
1297.3 0.1 0.1
1748.0 0.5 0.1
1826.1 1 0.1
1902.8 5 0.1
1915.4 10 0.1
338.3 0.1 0.01
984.1 0.5 0.01
1292.5 1 0.01
1731.5 5 0.01
1816.0 10 0.01
17.4 1 0.0001
86.3 1 0.0005
176.7 1 0.001
775.1 1 0.005
1292.5 1 0.01
2746.8 1 0.05
3183.2 1 0.1
3650.5 1 0.5
3719.9 1 1
The model simulated vertical ow between constant
head boundaries in layers 50 and 1, as a source and
sink, so that the overburden and shale properties control
the ow. The head in layers 50 and 1 was 1580 and
1550 m, respectively, to create a gradient of 0.019 over
the prole. Varying the gradient would have much less
effect on transport than changing K over several orders
of magnitude and was therefore not done.
Scenario 5 simulates injection using a WELL bound-
ary in layer 44, essentially the middle of the shale, from
columns 25 to 125 (Figure 2). It injects 15 million liters
over one 5-d stress period, or 3030 m
/d into 101 model
cells at the WELL. The modeled K
was changed to its
876 T. Myers GROUND WATER 50, no. 6: 872882
Figure 2. Model grid through layer 44 showing the horizon-
tal injection WELL (red) and DRAIN cells (yellow) used to
simulate owback. There is only one monitoring well because
the off-center well is not used in layer 44.
assumed fracked value at the beginning of the simula-
tion. Simulating high rate injection generates very high
heads in the model domain, similar to that found sim-
ulating oil discharging from the well in the Deepwater
Horizon crisis (Hsieh 2011) and water quality changes
caused by underground coal gasication (Contractor and
El-Didy 1989). DRAIN boundaries on both sides of the
WELL simulated return ow for 60 d after the completion
of (Figure 2), after which the DRAIN was deactivated.
The 60 d were broken into four stress periods, 1, 3, 6, and
50 d long, to simulate the changing heads and ow rates.
DRAIN conductance was calibrated so that 20% of the
injected volume returned within 60 d to emulate standard
industry practice (Alleman 2011; NYDEC 2009). Recov-
ery, continuing relaxation of the head at the well and the
adjustment of the head distribution around the domain,
occurred during the sixth period which lasted for 36,500 d.
There is no literature guidance to a preferred value
for fractured shale storage coefcient, so I estimated S
with a sensitivity analysis using scenario 3. With fractured
equal to 0.001 m/d, two orders of magnitude higher
than the in situ value, the time to equilibrium resulting
from simulation tests of three fractured shale storage
coefcients, 10
, 10
, and 10
/m, varied twofold
(Figure 3). The slowest time to equilibrium was for S =
/m (Figure 3), which was chosen for the transient
simulations because more water would be stored in the
shale and ow above the shale would change the least.
Scenario 1
Table 1 shows the conductivity and ux values
for various scenarios. The steady state travel time
Figure 3. Sensitivity of the modeled head response to the
storage coefcient used in the fractured shale for model layer
39 just above the shale.
for a particle through 1500 m of sandstone and shale
equilibrates with one of the formations controlling the
advection (Figure 4). For example, when the K
1 10
m/d, transport time does not vary with K
. For
at 0.1 m/d, transport time for varying K
ranges from
40,000 to 160 years. The lower travel time estimate is for
similar to that found by Schulze-Makuch et al. (1999).
The shortest simulated transport time of about 20 years
results from both the sandstone and shale K equaling
1 m/d. Other sensitivity scenarios emphasize the control
exhibited by one of the media (Figure 4). If K
is low,
travel time is very long and not sensitive to K
Scenario 2
The addition of a fault with K one to two orders of
magnitude more conductive than the surrounding sand-
stone increased the particle travel rate by about 10 times
(compare Figure 5 with Figure 4). The fault K controlled
the transport rate for K
less than 0.01 m/d. A highly
Figure 4. Sensitivity of particle transport time over
1500 m for varying shale and sandstone vertical K.
Effective porosity equals 0.1. (1)varying K
, K
m/d; (2)varying K
, K
= 0.1 m/d; (3)varying
, K
= 0.1 m/d; (4)varying K
, K
= 0.01 m/d; and
(5)varying K
, K
= 1.0 m/d. T. Myers GROUND WATER 50, no. 6: 872882 877
Figure 5. Variability of transport through various scenarios
of changing the K for the fault or shale. Effective porosity
equals 0.1. (1)varying K
, K
= 0.01 m/d; (2)varying
, K
= 0.1 m/d; (3)no fault; (4)varying K fault,
= 0.1 m/d, K
= 0.01 m/d. Unless specied, the vertical
fault has K = 1 m/d for variable K
conductive fault could transport uids to the surface in
as little as a year for K
equal to 0.01 m/d (Figure 5).
However, a fault did not signicantly change the overall
model ux, so with fault values are not shown in Table 1.
Scenarios 3 and 4
Scenarios 3 and 4 estimate the time to establish a
new equilibrium once the K
changes, due to fracking,
between values specied in scenarios 1 and 2. Equilibrium
times vary by model layer as the changes propagate
through the domain, and ux rate for the simulated
changes imposed on natural background conditions. The
fracking-induced changes cause a signicant decrease in
the head drop across the shale and the time for adjustment
of the potentiometric surface to a new steady state depends
on the new shale properties.
The time to equilibrium for one scenario 3 simulation,
changing from 10
to 10
m/d with K
to 0.1 m/d, varied from 5.5 to 6.5 years, depending
on model layer (Figure 6). Near the shale (layers 39
and 40), the potentiometric surface increased from 23
to 25 m reecting the decreased head drop across the
shale. One hundred meters higher, in layer 20, the
potentiometric surface increased about 20 m. Simulation
of scenario 4, with a fault with K = 1 m/d, decreased
the time to equilibrium to from 3 to 6 years within the
fault zone, depending on model layer (Figure 6). Highly
fractured sandstone would allow more vertical transport,
but advective ow would also increase so that the base
would control the overall rate.
The ux across the upper boundary changed within
100 years for scenario 3 from 1.7 to 345 m
/d, or
0.000008 to 0.0017 m/d, reecting control by K
. There
is little difference in the equilibrium uxes between
scenario 3 and 4 indicating that the fault primarily affects
the time to equilibrium rather than the long-term ow rate.
Figure 6. Monitoring well water levels for specied model
layers due to fracking of the shale; monitor well in the center
of the domain, including in the fault, K of the shale changes
from 0.00001 to 0.01 m/d at the beginning of the simulation.
Scenario 5: Simulation of Injection
The injection scenarios simulate 15 million liters
entering the domain at the horizontal well and the
subsequent potentiometric surface and ux changes
throughout. The highest potentiometric surface increases
(highest injection pressure) occurred at the end of injec-
tion (Figure 7), with a 2400 m increase at the horizontal
well. The simulated peak pressure both decreased and
occurred longer after the cessation of injection with dis-
tance from the well (Figures 7 and 8). The pressure at
the well returned to within 4 m of pre-injection levels in
about 300 d (Figure 7). After injection ceases, the peak
pressure simulated further from the well occurs longer
from the time of cessation, which indicates there is a pres-
sure divide beyond which uid continues to ow away
from the well bore while within which the uid ows
toward the well bore. The simulated head returned to
near pre-injection levels slower with distance from the
well (Figure 7), with levels at the edge of the shale (layer
40) and in the near-shale sandstone (layer 39) requir-
ing several hundred days to recover. After recovering
from injection, the potentiometric surface above the shale
increased in response to ux through the shale adjusting to
the change in shale properties (Figure 8), as simulated in
scenario 3. The scenario required about 6000 d (16 years)
for the potentiometric surface to stabilize at new, higher,
levels (Figure 8). Removing the fault from the simulation
had little effect on the time to stabilization, and is not
Prior to injection, the steady ux for in situ shale
= 10
m/d) was generally less than 2 m
/d and
varied little with K
(Figure 4). Once the shale was
fractured, the sandstone controlled the ux which ranges
from 38 to 135 m
/d as K
ranges from 0.01 to 0.1 m/d
(Figure 9), resulting in particle travel times of 2390 and
616 years, respectively. More conductive shale would
allow faster transport (Figure 4). Adding a fault to the
scenario with K
equal to 0.01 m/d increased the ux to
approximately 63 m
/d and decreased the particle travel
878 T. Myers GROUND WATER 50, no. 6: 872882
Figure 7. Simulated potentiometric surface changes by layer
for specied injection and media properties. The monitoring
point is in the center of the domain. Fault is included. K
0.01 m/d, K
= 0.001 m/d. S (fractured shale) = 0.001/m,
S (ss) = 0.0001/m.
Figure 8. Simulated potentiometric surface changes for lay-
ers within the shale and sandstone. CW is center moni-
toring well and EW is east monitoring well, about 120 m
from the centerline. Fault is included. The line for layer
2, CW plots beneath the line for layer 2, EW. K
0.01 m/d, K
= 0.001 m/d, S (fractured shale) = 0.001/m,
S (ss) = 0.0001/m.
time to 31 years. Approximately, 36 m
/d owed through
the fault (Figure 9). The fault properties control the
particle travel time, especially if the fault K is two or
more orders of magnitude higher than the sandstone.
Simulated owback varied little with K
because it
had been calibrated to be 20% of the injection volume.
A lower storage coefcient or higher K would allow the
injected uid to move further from the well, which would
lead to less owback.
Vertical ux through the overall section with a fault
varies signicantly with time, due to the adjustments in
potentiometric surface. One day after injection, vertical
ux exceeds signicantly the pre-injection ux about
200 m above the shale (Figure 10). After 600 d, the
vertical ux near the shale is about 68 m
/d and in
Figure 9. Comparison of ux for three scenarios. Flowback
is the same for all scenarios. (1): Kss = 0.01 m/d, Ksh =
0.001 m/d, Fault K = 1 m/d; (2): Kss = 0.01 m/d, Ksh =
0.001 m/d, no fault; (3) Kss = 0.1 m/d, Ksh = 0.001 m/d, no
Figure 10. Upward ux across the domain section as a
function of distance above the top of the shale layer. Cross
section is 202,500 m
layer 2 about 58 m
/d; it approaches steady state through
all sections after 100 years with ux equaling about
62.6 m
/d. The 100-year ux is 61.5 m
/d higher than
the pre-injection ux because of the changed shale
The interpretative modeling completed herein has
revealed several facts about fracking. First, MODFLOW
can be coded to adequately simulate fracking. Simulated
pressures are high, but velocities even near the well do
not violate the assumptions for Darcian ow. Second,
injection for 5 d causes extremely high pressure within
the shale. The pressure decreases with distance from the
well. The time to maximum pressure away from the well
lags the time of maximum pressure at the well. The
pressure drops back to close to its pre-injection level T. Myers GROUND WATER 50, no. 6: 872882 879
at the well within 300 d, indicating the injection affects
the ow for signicantly longer periods than just during
the fracking operation. Although the times may vary
based on media properties, the difference would be at
most a month or so, based on the various combinations
of properties simulated. The system transitions within
6 years due to changes in the shale properties. The
equilibrium transport rate would transition from a system
requiring thousands of years to one requiring less than
100 years within less than 10 years.
Third, most of the injected water in the simulation
ows vertically rather than horizontally through the shale.
This reects the higher K
20 m above the well and the no
ow boundary within 225 m laterally from the well, which
emulates in situ shale properties that would manifest at
some distance in the shale.
Fourth, the interpretative model accurately and real-
istically simulates long-term steady state ow conditions,
with an upward ow that would advect whatever conser-
vative constituents exist at depth. Using low, unfractured
K values, the transport simulation may correspond with
advective transport over geologic time although there are
conditions for which it would occur much more quickly
(Figure 4). If the K
is 0.01 m/d, transport could occur
on the order of a few hundreds of years. Faults through the
overburden could speed the transport time considerably.
Reasonable scenarios presented herein suggest the travel
time could be decreased further by an order of magnitude.
Fifth, fracking increases the K
by several orders
of magnitude. Out-of-formation fracking (Fisher and
Warpinski 2011) would increase the K in the overburden,
thereby changing the regional hydrogeology. Vertical ow
could change over broad areas if the expected density
of wells in the Marcellus shale region (NYDEC 2009)
actually occurs.
Sixth, if newly fractured shale or out-of-formation
fractures come close to contacting fault fracture zones,
contaminants could reach surface areas in tens of years,
or less. Faults can decrease the simulated particle travel
time several orders of magnitude.
Fracking can release uids and contaminants from
the shale either by changing the shale and overburden
hydrogeology or simply by the injected uid forcing other
uids out of the shale. The complexities of contaminant
transport from hydraulically fractured shale to near-
surface aquifers render estimates uncertain, but a range
of interpretative simulations suggest that transport times
could be decreased from geologic time scales to as
few as tens of years. Preferential ow through natural
fractures fracking-induced fractures could further decrease
the travel times to as little as just a few years.
There is no data to verify either the pre- or
post-fracking properties of the shale. The evidence for
potential vertical contaminant ow is strong, but there
are also almost no monitoring systems that would
detect contaminant transport as considered herein. Several
improvements could be made.
Prior to hydraulic fracturing operations, the subsurface
should be mapped for the presence of faults and
measurement of their properties.
A reasonable setback distance from the fracking to
the faults should be established. The setback distance
should be based on a reasonable risk analysis of fracking
increasing the pressures within the fault.
The properties of the shale should be veried, post-
fracking, to assess how the hydrogeology will change.
A system of deep and shallow monitoring wells and
piezometers should be established in areas expect-
ing signicant development, before that development
begins (Williams 2010).
This research was funded by the Park Foundation
and Catskill Mountainkeepers. The author thanks Anthony
Ingraffea, Paul Rubin, Evan Hansen, two anonymous
reviews, and the journal editor for helpful comments on
this article.
Alleman, D. 2011. Water used for hydraulic fracturing: Amounts,
sources, reuse, and disposal. In Hydraulic Fracturing of the
Marcellus Shale, National Ground Water Association short
course. Baltimore, Maryland.
Annunziatellis, A., S.E. Beaubien, S. Bigi, G. Ciotoli, M.
Coltella, and S. Lombardi. 2008. Gas migration along
fault systems and through the vadose zone in the Latera
calder (central Italy): Implications for CO
storage. International Journal of Greenhouse Gas Control
2; 353372. DOI: 10.1016/j.ijggc.2008.02.003.
Arthur, J.D., B. Bohm, and M. Layne. 2008. Hydraulic Fractur-
ing Consideration for Natural Gas Wells of the Marcellus
Shale. Cincinnati, Ohio: Ground Water Protection Council.
Boyer, C., J. Kieschnick, R. Suarez-Rivera, R.E. Lewis, and
G. Waters. 2006. Producing gas from its source. Oilelds
Review 18, no. 3: 3649.
Breen, K.J., K. Revesz, F.J. Baldassare, and S.D. McAuley.
2007. Natural gases in ground water near Tioga Junction,
Tioga County, north-central PennsylvaniaOccurrence
and use of isotopes to determine origins. Scientic Inves-
tigations Report Series 2007-5085. Reston, Virginia: U.S.
Geological Survey.
Caine, J.S., J.P. Evans, and C.B. Forster. 1996. Fault zone
architecture and permeability structure. Geology 24, no. 11:
Contractor, D.N., and S.M.A. El-Didy. 1989. Field application
of a nite-element water-quality model to a coal seam with
UCG burns. Journal of Hydrology 109: 5764.
Davies, R.J. 2011. Methane contamination of drinking water
caused by hydraulic fracturing remains unproven. Proceed-
ings of the National Academy of Sciences USA 108: E871.
DiGiulio, D.C., R.T. Wilkin, C. Miller, and G. Oberly. 2011.
DRAFT: Investigation of Ground Water Contamination near
Pavillion, Wyoming. Ada, Oklahoma: U.S. Environmental
Protection Agency, Ofce of Research and Development.
Dresel, P.E., and A.W. Rose. 2010. Chemistry and origin of oil
and gas well brines in western Pennsylvania. 4th ser., Open-
File Report OFOG 1001.0, 48. Harrisburg: Pennsylvania
Geological Survey.
880 T. Myers GROUND WATER 50, no. 6: 872882
Edward, K.L., and S. Weisset. 2011. Marcellus shale hydraulic
fracturing and optimal well spacing to maximize recov-
ery and control costs. Paper 140463 in SPE Hydraulic
Fracturing Technology Conference, January 2426, 2011,
The Woodlands, Texas.
Energy Information Administration (EIA). 2009. Annual Energy
Outlook with Projections to 2030. Washington, D.C.: U.S.
Department of Energy.
(accessed May 23, 2011).
Engelder, T., G.G. Lash, and R.S. Uzcategui. 2009. Joint sets
that enhance production from Middle and Upper Devonian
gas shales of the Appalachian Basin. AAPG Bulletin 93, no.
7: 857889.
Environmental Protection Agency (EPA). 1987. Report to
congress, management of wastes from the exploration,
development, and production of crude oil, natural gas, and
geothermal energy, volume 1 of 3, oil and gas. Washington,
Etiope, G., and G. Martinelli. 2002. Migration of carrier and
trace gases in the geosphere: an overview. Physics of the
Earth and Planetary Interiors 129: 34.
Fisher, K., and N. Warpinski. 2011. Hydraulic fracture-height
growth: real data. Paper SPE 145949 presented at the
Annual Technical Conference and Exhibition, Denver,
Colorado. DOI: 10.2118/145949-MS.
Gold, D. 1999. Lineaments and their interregional relation-
ships. In The Geology of Pennsylvania, ed. C.H. Schultz,
307313. Harrisburg: Pennsylvania Department of Conser-
vation and Natural Resources.
Halford, K.J., and R.T. Hanson. 2002. User guide for the
drawdown-limited, multi-node well (MNW) package for the
U.S. Geological Surveys modular three-dimensional nite-
difference ground-water ow model, Versions MODFLOW-
96 and MODFLOW-2000, 33. Open-File Report 02-293.
Sacramento, California: U.S. Geological Survey.
Harbaugh, A.W., E.R. Banta, M.C. Hill, and M.G. McDonald.
2000. Modow-2000, the U.S. Geological Survey modular
ground-water model User guide to modularization, con-
cepts and the ground-water ow process. Open-File Report
00-92. Reston, Virginia: U.S. Geological Survey.
Harper, J.A. 1999. Devonian. In The Geology of Pennsylva-
nia, ed. C.H. Schultz, 108127. Harrisburg: Pennsylvania
Department of Conservation and Natural Resources.
Hill, M.C., and C.R. Tiedeman. 2007. Effective Groundwater
Model Calibration: With Analysis of Data, Sensitivities,
Predictions, and Uncertainty. Hoboken, New Jersey: John
Wiley and Sons.
Hsieh, P.A. 2011. Application of MODFLOW for oil reser-
voir simulation during the Deepwater Horizon crisis.
Ground Water 49, no. 3: 319323. DOI: 10.1111/j.1745-
Isachsen, Y.W., and W. McKendree. 1977. Preliminary Brittle
Structure Map of New York, Map and Chart Series No. 31.
Albany, New York: New York State Museum.
Jehn, P. 2011. Well and water testingWhat to look for
and when to look for it. In Groundwater and Hydraulic
Fracturing of the Marcellus Shale, National Gound Water
Association short course. Baltimore, Maryland.
King, G. 2010. Thirty years of gas shale fracturing: What have
we learned? Paper SPE 133456 presented at the SPE Annual
Technical Conference and Exhibition, September 1922,
2010, Florence, Italy.
King, G.E., L. Haile, J. Shuss, and T.A. Dobkins. 2008. Increas-
ing fracture path complexity and controlling downward
fracture growth in the Barnett shale. Paper 119896 pre-
sented at the SPE Shale Gas Production Conference,
November 1618, 2008, Fort Worth, Texas.
Konikow, L.F. 2011. The secret to successful solute-transport
modeling. Ground Water 49, no. 2: 144159. DOI:
Kramer, D. 2011. Shale-gas extraction faces growing public and
regulatory challenges. Physics Today 64, no. 7: 2325.
Kwon, O., A.K. Kronenberg, A.F. Gangi, B. Johnson, and
B.E. Herbert. 2004a. Permeability of illite-bearing shale:
1. Anisotropy and effects of clay content and loading.
Journal of Geophysical Research 109: B10205. DOI:
Kwon, O., B.E. Herbert, and A.K. Kronenberg. 2004b. Perme-
ability of illite-bearing shale: 2. Inuence of uid chemistry
on ow and functionally connected pores. Journal of Geo-
physical Research 109: B10206. DOI:10.1029/2004JB00
Lacombe, S., E.A. Sudickey, S.K. Frape, and A.J.A. Unger.
1995. Inuence of leak boreholes on cross-formational
groundwater ow and contaminant transport. Water
Resources Research 31, no. 8: 18711882.
Langevin, C.D., W.B. Shoemaker, and W. Guo. 2003. MODF
LOW-2000, the U.S. Geological Survey modular ground-
water model Documentation of the SEAWAT-2000 version
with the variable-density ow process (VDF) and the
integrated MT3DMS transport process (IMT), 43 p. Open-
File Report 03-426. Tallahassee, Florida: U.S. Geological
Loyd, O.B., and L.D. Carswell. 1981. Groundwater resources of
the Williamsport region, Lycoming County, Pennsylvania.
Water Resources Report 51. Pennsylvania: Department of
Environmental Resources.
Neuzil, C.E. 1994. How permeable are clays and shales? Water
Resources Research 30, no. 2: 145150.
Neuzil, C.E. 1986. Groundwater ow in low-permeability envi-
ronments. Water Resources Research 22, no. 8: 11631195.
New York State Department of Environmental Conservation
(NYDEC). 2009. Draft Supplemental Generic Environmen-
tal Impact Statement on the Oil, Gas and Solution Min-
ing Regulatory ProgramWell Permit Issuance for Hori-
zontal Drilling and High-Volume Hydraulic Fracturing to
Develop the Marcellus Shale and Other Low-Permeability
Gas Reservoirs. Albany, New York: State Department of
Environmental Conservation.
Nickelsen, R.P. 1986. Cleavage duplexes in the Marcellus Shale
of the Appalachian foreland. Journal of Structural Geology
8, no. 3/4: 361371.
Ohio Department of Natural Resources (ODNR). 2008. Report
on the Investigation of the Natural Gas Invasion of Aquifers
in Bainbridge Township of Geauga County, Ohio. ODNR,
Division of Mineral Resources Management.
2011. Methane contamination of drinking water accom-
panying gas-well drilling and hydraulic fracturing. Pro-
ceedings of the National Academy of Sciences 108, no. 20:
Osborn, S.G., and J.C. McIntosh. 2010. Chemical and isotopic
tracers of the contribution of microbial gas in Devonian
organic-rich shales and reservoir sandstones, northern
Appalachian Basin. Applied Geochemistry 25, no. 3:
Pennsylvania Bureau of Topographic and Geologic Survey
(PBTGS). 2001. Bedrock Geology of Pennsylvania (Digital
Files). Harrisburg, Pennsylvania: PA Department of Con-
servation and Natural Resources.
Pennsylvania Department of Environmental Protection (PADEP).
2011. Marcellus Shale.
(accessed June 1, 2011).
Pennsylvania Department of Environmental Protection (PADEP).
2009. Notice of Violation, Re: Gas Migration Investigation,
Dimock Township, Susquehanna County, Letter from S. C.
Lobins, Regional Manager, Oil and Gas Management, to
Mr. Thomas Liberatore, Cabot Oil and Gas Corporation.
February 27, 2009. 4 pp., PADEP, Meadville, Pennsylvania. T. Myers GROUND WATER 50, no. 6: 872882 881
Schoell, M. 1980. The hydrogen and carbon isotopic compo-
sition of methane from natural gases of various origins.
Geochemica et Cosmochiica Act 44, no. 5: 649661.
Schulze-Makuch, D., D.A. Carlson, D.S. Cherkauer, and
P. Malik. 1999. Scale dependence of hydraulic conductivity
in heterogeneous media. Ground Water 37, no. 6: 904919.
Schweitzer, R., and H.I. Bilgesu. 2009. The role of economics
on well and fracture design completions of Marcellus Shale
wells. Paper 125975 in SPE Eastern Regional Meeting,
September 2325, 2009, Charleston, West Virginia.
Shoemaker, W.B., E.L. Kuniansky, S. Birk, S. Bauer, and
E.D. Swain. 2008. Documentation of a conduit ow process
(CFP) for MODFLOW-2005. U.S. Geological Survey
techniques and methods, Book 6, chapter A24, 50 p.
Reston, Virgina: U.S. Geological Survey.
Silliman, S., and D. Higgins. 1990. An analytical solution for
steady-state ow between aquifers through an open well.
Ground Water 28, no. 2: 184190.
Soeder, D.J. 2010. The Marcellus shale: Resources and reserva-
tions. EOS 91, no. 32: 277278.
T.A.L. Research and Development (TAL). 1981. Geology,
drill holes, and geothermal energy potential of the basal
Cambrian rock units of the Appalachian Basin of New York
State. In Prepared for New York State Energy Research and
Development Authority, 54 p.
West Virginia Geological and Economic Survey (WVGES).
2011. Completed WellsMarcellus Shale, West Virginia.
Morgantown, West Virginia: WVGES.
West Virginia Geological and Economic Survey (WVGES).
2010a. Structural Geologic Map (Faults)Topo of the
Onondaga Limestone or Equivalent, West Virginia. Morgan-
town, West Virginia: WVGES.
West Virginia Geological and Economic Survey (WVGES).
2010b. Structural Geologic Map (Folds)Topo of the
Onondaga Limestone or Equivalent, West Virginia. Morgan-
town, West Virginia: WVGES.
White, J.S., and M.V. Mathes. 2006. Dissolved-gas concentra-
tion in ground water in West Virginia, 8. Data Series 156.
Reston, Virginia: U.S. Geological Survey.
Williams, J.H. 2010. Evaluation of well logs for determining
the presence of freshwater, saltwater, and gas above the
Marcellus shale in Chemung, Tioga, and Broome Counties,
New York. Scientic Investigations Report 20105224, 27.
Reston, Virginia: U.S. Geological Survey.
Williams, J.H., L.E. Taylor, and D.J. Low. 1998. Hydrogeology
and groundwater quality of the Glaciated Valleys of
Bradford, Tioga, and Potter Counties, Pennsylvania. Water
Resource Report 68. Pennsylvania Dept of Conservation
and Natural Resources and U.S. Geological Survey.
Wunsch, D. 2011. Hydrogeology and hydrogeochemistry of
aquifers overlying the Marcellus shale. In Groundwa-
ter and Hydraulic Fracturing of the Marcellus Shale,
National Ground Water Association short course. Balti-
more, Maryland.
882 T. Myers GROUND WATER 50, no. 6: 872882