Scholarly article on topic 'Quantifying the potential effects of high-volume water extractions on water resources during natural gas development: Marcellus Shale, NY'

Quantifying the potential effects of high-volume water extractions on water resources during natural gas development: Marcellus Shale, NY Academic research paper on "Earth and related environmental sciences"

Share paper
{"Groundwater flow modeling" / "Marcellus Shale" / "New York State" / "High-volume hydraulic fracturing" / "Water quantity"}

Abstract of research paper on Earth and related environmental sciences, author of scientific article — Laura C. Best, Christopher S. Lowry

Abstract Study region The Marcellus Shale, New York State, USA. Study focus Development of natural gas resources within the Marcellus Shale will require large volumes of water if high-volume hydraulic fracturing expands into New York State. Although this region has ample fresh water resources, it is necessary to explore the response of hydraulically connected groundwater and surface water systems to large withdrawals. Because such effects would not be apparent from a typical water budget approach, this study applied groundwater flow modelling under scenarios of high-volume water withdrawals. Emphasis on water quantity, in contrast with other lines of research concerning water quality, introduced an important perspective to this controversial topic. New hydrological insights for the region The potential effects of the withdrawal scenarios on both the water table and stream discharge were quantified. Based on these impact results, locations in the aquifer and stream networks were identified, which demonstrate particular vulnerability to increased withdrawals and their distribution. These are the locations of importance for planners and regulators who oversee water permitting, to reach a sustainable management of the water resources under changing conditions of energy and corresponding water demand.

Academic research paper on topic "Quantifying the potential effects of high-volume water extractions on water resources during natural gas development: Marcellus Shale, NY"

journal of Hydrology: Regional Studies xxx (2014) xxx-xxx

Contents lists available at ScienceDirect

Journal of Hydrology: Regional 1 t[ Studies

ELSEVIER journal homepage

Quantifying the potential effects of high-volume water extractions on water resources during natural gas development: Marcellus Shale, NY

4 qi Laura C. Best, Christopher S. Lowry *

5 University at Buffalo Department of Geology, 411 Cooke Hall, Buffalo, NY 14620, United States

10 11 12

article info

Article history: Received 13 February 2014 Received in revised form 8 April 2014 Accepted 1 May 2014 Available online xxx


Groundwater flow modeling Marcellus Shale New York State

High-volume hydraulic fracturing Water quantity


Study region: The Marcellus Shale, New York State, USA. Study focus: Development of natural gas resources within the Marcellus Shale will require large volumes of water if high-volume hydraulic fracturing expands into New York State. Although this region has ample fresh water resources, it is necessary to explore the response of hydraulically connected groundwater and surface water systems to large withdrawals. Because such effects would not be apparent from a typical water budget approach, this study applied groundwater flow modelling under scenarios of high-volume water withdrawals. Emphasis on water quantity, in contrast with other lines of research concerning water quality, introduced an important perspective to this controversial topic. New hydrological insights for the region: The potential effects of the withdrawal scenarios on both the water table and stream discharge were quantified. Based on these impact results, locations in the aquifer and stream networks were identified, which demonstrate particular vulnerability to increased withdrawals and their distribution. These are the locations of importance for planners and regulators who oversee water permitting, to reach a sustainable management of the water resources under changing conditions of energy and corresponding water demand.

© 2014 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (

* Corresponding author. Tel.: +1 716 645 4266. E-mail addresses: (L.C. Best), (C.S. Lowry).


2214-5818/© 2014 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (

2 L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

21 1. Introduction

22 Unconventional natural gas production from shale formations provides a significant domestic

23 energy source in the United States (U.S. Energy Information Administration, 2011). Natural gas extrac-

24 tion from tight geologic formations has increased due to technological advancements of horizontal

25 drilling, leading to economic viability of previously untapped reserves (U.S. Department of Energy,

26 2009). The potential expansion of high-volume hydraulic fracturing (HVHF) of the Marcellus and Utica

27 Shale into New York State to extract natural gas resources is a controversial issue for policy makers,

28 industry stakeholders, and community members. Issues surrounding this debate range from socioeco-

29 nomic to logistic to environmental, with emphasis on water quality dictating the direction of scientific

30 research and media attention.

31 Recently, other environmental concerns associated with HVHF in New York have come to the

32 forefront of discussion. This includes a water quantity perspective, which is traditionally less crit-

33 ical in regions that have ample freshwater supplies in humid climates and/or large, proximate

34 freshwater bodies (Rahm and Riha, 2012). HVHF requires large volumes of water which will ulti-

35 mately increase water demand from the regions that will experience development. Increased water

36 demand will prompt regulators to determine from where, and at what rate, this water should be

37 extracted to protect sustainable use for drinking water, agriculture, and other industry demands.

38 Altered stream geochemistry and consequences to stream ecosystems, as a result of decreased

39 stream discharge, are factors beyond the anthropogenic freshwater demands mentioned above

40 that may merit consideration. Although water budgets from the New York State Department of

41 Environmental Conservation (NYSDEC) demonstrate that increased water demands from HVHF in

42 New York would make up a minor fraction of total water use (NYSDEC, 2011), it is unclear how

43 hydraulically linked groundwater-surface water systems might respond to such a development.

44 Water budgets alone may not be sufficient in predicting the spatially variable response of these

45 systems, particularly in identifying areas which present heightened sensitivity to withdrawals.

46 For example, the response of aquifers and streams to increased withdrawals of water might vary

47 as a function of valley width, thickness and depth of aquifers within the valley fill. Addition-

48 ally, smaller streams might be vulnerable to induced changes in groundwater discharge during

49 drought.

50 The projected path of HVHF development of the Marcellus Shale in New York will most likely

51 focus on the Southern Tier of the state, including Broome and Tioga counties (Fig. 1). The major

52 valleys within these counties overlie an unconsolidated glacial valley-fill aquifer network which

53 has been classified as a sole source aquifer since 1985 (U.S. Environmental Protection Agency,

54 2010). Such a designation emphasizes the importance of this groundwater source to the over-

55 lying municipalities, which receive more than half of their drinking water from the aquifer. In

56 this region there is a high degree of hydraulic connectivity between streams and underlying

57 unconsolidated glacial deposits (Randall, 1977; Wolcott and Coon, 2001; Yager, 1993). High-

58 volume withdrawals of water from groundwater may elicit a response from surface water, or

59 vice versa, due to their physical connectivity (Winter et al., 1998). It is therefore necessary to

60 investigate how different development scenarios might affect both the water table and stream

61 flow.

62 This research focuses on the use of groundwater flow modeling to determine if increased

63 water demand associated with HVHF is enough to cause significant change to groundwater lev-

64 els and stream flow within the study area. The objective of this research is to identify scenarios

65 and locations that are particularly vulnerable to high-volume withdrawals of water and may

66 require further evaluation should water permits be requested. A simulated range of develop-

67 ment scenarios demonstrate how varying well pad density, water source, and water volume

68 might affect the groundwater-surface water systems in the Southern Tier of New York. The

69 importance of this research lies in its application to all stakeholders in the HVHF controversy

70 currently underway in New York. Not only will policy makers and regulators benefit from

71 the predictive capacity of computer modeling, but industry, community members and interest

72 groups can better understand how a water quantity perspective is valuable for sustainable energy development.

L.C. Best, C.S. Lowry / Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

Fig. 1. Marcellus Shale extent across New York and Pennsylvania is highlighted in gray (modified from the U.S. Energy Information Administration and the U.S. Census Bureau). Horizontal wells drilled in the Marcellus Shale in Pennsylvania are in purple and demonstrate the projected path of HVHF in New York State (data from, source: Pennsylvania Department q8 of Environmental Protection). Broome and Tioga counties are highlighted in red. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

73 2. Background

74 Hydraulic fracturing is a process that involves the injection of water into the subsurface in order

75 to fracture tight geologic formations. Fracturing creates pathways through which trapped natural

76 gas flows freely into a well and is subsequently harnessed for energy. The combination of hydraulic

77 fracturing and horizontal drilling has led to the growing viability of unconventional shale plays (U.S.

78 Department of Energy, 2009). Horizontal drilling refers to the lateral drilling of a well bore through a

79 target formation. This allows access to a greater volume of gas-bearing rock, making such drilling ven-

80 tures economically feasible (Soeder, 2010). In HVHF, large volumes of water in addition to proppants

81 and other additives serve as the fracturing fluid. The fluid injection and fracturing process progresses

82 in stages along the horizontal extent of the well, with each horizontal well requiring between 1 and

83 5 million gallons of water (Gregory et al., 2011). Only a fraction of injected fluid actually returns to

84 the surface - referred to as flowback - with the unreturned volume remaining in the subsurface. This

85 fraction can vary greatly between wells, company, and target formation with an estimated average of

86 10-40% flowback (Maloney and Yoxtheimer, 2012; NYSDEC, 2011; Rassenfoss, 2011).

87 In arid climates, where freshwater supply is limited, the quantity of water use associated with

88 HVHF is of concern (Nicot and Scanlon, 2012). In humid climates, where freshwater supply is less

89 emphasized in water resource management, increased water demand associated with HVHF is only

90 beginning to receive recognition (Rahm and Riha, 2012). This is in part due to mass balance or water

91 budget approaches in quantifying the impacts of HVHF water demands. Nicot and Scanlon (2012)

92 estimate water use associated with HVHF is less than 1% of water use in Texas, but may account

93 for larger fractions of water use at the county scale. For example, within the Barnett Shale play in

94 Texas, the 2008 fraction of shale gas water use in the counties of Denton, Johnson, Parker, Tarrant,

95 and Wise was 2.8%, 29%, 10%, 1.4%, and 19%, respectively (Nicot and Scanlon, 2012). Such an inves-

96 tigation emphasizes the numerical insignificance of increased water use from the shale gas industry

97 at a statewide scale, although simultaneously recognizing localized (county-scale) significance. The

98 NYSDEC (2011), estimates that HVHF development would increase water demand by 0.24%. While

L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

99 it is important to acknowledge that an increase of less than 1% of increased water demand is small,

100 localized impacts should not be ignored. Groundwater flow modeling offers a different approach to

101 evaluating increased water demand in the Southern Tier of New York State. This approach captures

102 both regional and localized impacts while complying with the dynamic relationship between stream

103 flow and groundwater.

104 The NYSDEC (2011) predicts a peak development of 2462 wells in one year across the state of New 10Q2 York, with four wells most likely developed on one well pad. It is also estimated that about 2.4 to

106 7.8 million gallons (Mgal) will be used for each horizontal well. Accounting for the recycling of flow-

107 back water, approximately 3.6 Mgal of freshwater for each horizontal well will be required, assuming

108 that 15% of the average demand of 4.2 Mgal is recycled flowback water (NYSDEC, 2011). These pro-

109 jections are the basis for setting up the range of development scenarios to simulate in this research.

110 In addition to well density and water volume, water source is also included in the development sce-

111 narios. Although surface water may be the most likely source (NYSDEC, 2011), municipal pumping

112 wells in Pennsylvania do provide some of the water used in HVHF (Rahm and Riha, 2012). Therefore,

113 both groundwater and surface water are accounted for as potential water sources in the development

114 scenarios. Accounting for both groundwater and surface water withdrawals makes this type of inves-

115 tigation applicable to the HVHF development in the short-term as well as future potential long-term

116 changes in water resources, which may involve surface and groundwater.

117 3. Site characterization

118 The aquifer network that underlies Broome and Tioga counties is part of a complex glacial valley-fill

119 system (Fig. 3). The glacial sediments are a legacy of the Late Wisconsin stage of the last Pleistocene

120 glaciation (Aber, 1980; Scully and Arnold, 1981), deposited approximately 16,650 years ago (Cadwell,

121 1973). The aquifer is composed primarily of ice contact deposits overlain by glacial outwash, which

122 was deposited via meltwater streams (Randall, 1978). The unconsolidated glacial deposits, mainly silty

123 sand and gravel, overlie a thin, discontinuous till, which is underlain by fractured, noncalcareous Devo-

124 nian bedrock (Scully and Arnold, 1981). Geographically discontinuous lacustrine silt and clay overlie

125 ice-contact deposits, generating confined aquifers in parts of the network (MacNish and Randall, 1982;

126 Randall, 1978,1986).

127 Previous work within the proposed study area has clearly defined the depositional history, hydro-

128 logic properties, and hydrostratigraphy of the aquifer network (Fleisher, 1986; Kontis et al., 2004;

129 Randall, 1972,1978, 2001; Reynolds and Williams, 1988; Waller and Finch, 1982). Furthermore, por-

130 tions of the aquifer network, particularly sections which underlie of metropolitan area of Binghamton

131 in Broome County, New York, have been previously modeled (Coon etal., 1998; Randall, 1986; Wolcott

132 and Coon, 2001; Yager, 1986, 1993). Considering the extent to which gas ventures will most likely

133 expand, it is desirable to extend the modeled areas to simulate the regional flow paths throughout

134 Broome and Tioga counties.

135 Within these counties there is a high degree of hydraulic connectivity between streams and the

136 underlying aquifer (Randall, 1977; Wolcott and Coon, 2001; Yager, 1993). Additionally, pumping

137 induced recharge from streambed infiltration is significant in the study area (Kontis et al., 2004;

138 Randall, 2001). If municipal pumping rates increase, it becomes important to account for the pos-

139 sibility of added induced recharge. Conversely, groundwater discharge from stratified drift aquifers is

140 the main source of base flow to streams during periods of drought (Randall, 2010). Increased ground-

141 water pumping rates, therefore, would commonly reduce aquifer discharge to streams resulting in

142 reduced stream flow (Randall et al., 1988), although a few broad valleys are drained only by small

143 streams of local origin.

144 4. Methods

145 The most significant groundwater flow occurs within the broad valley drift aquifers, limited to

146 the main glacial valleys. Major streams in this setting are parallel to the axes of the valley walls

147 and would not help to constrain the hydrologic boundaries for groundwater flow. Because there are

148 limited natural hydrologic features for use as boundary conditions, a two-dimensional watershed scale

149 analytic element model (Jankovic and Barnes, 1999) was first constructed in Visual AEM (Craig and

150 Matott, 2009) to develop boundary conditions for the localized area of interest (Hunt et al., 1998).

151 The scope of the first model encompasses the Upper Susquehanna River basin, including the valleys

152 of Broome and Tioga counties. Using constant head boundary conditions from Visual AEM, a three-

153 dimensional finite difference MODFLOW model (Harbaugh, 2005) was built to focus on the valleys of

154 interest (Fig. 2). The extracted constant head boundaries were placed along the perimeter of the model

155 extent and are significant in their simulation of upland recharge to the valley-fill aquifer network.

156 Furthermore, the analytic element model was calibrated to real-time stream discharge measurements

157 in order to approximate net regional groundwater recharge.

158 The finite difference grid was set up in Groundwater Vistas Version 6 (Rumbaugh and Rumbaugh,

159 2011). The grid is comprised of 193 rows and 281 columns of 250 m x 250 m cells, with a total surface

160 area of approximately 3390 km2 (Best, 2013). The model contains five layers in order to represent the

A | Outwash | Ice contact deposits | Lacustrine deposits □ TN, | Bedrock

Fig. 3. Generalized, conceptual cross-sections of a wide (A) and a narrow glacial valley (B) (modified from Randall, 2001).

6 L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

Upper aquifer Lower aquifer [ | Clay layer Uplands


mm 1 L"LJ—

Fig. 4. Cross-section illustrating the model layers and hydraulic conductivity zones derived from the conceptual valley stratigraphy.

161 hydrostratigraphy of the glacial drift while maintaining the reasonable level of simplicity needed to

162 address the research question from a regional perspective. The top elevation of the uppermost layer

163 represents the land surface and was approximated using an imported and resampled digital elevation

164 model. The bottom elevation of the fifth layer represents the bedrock surface, thereby constraining

165 flow within the valley fill thickness ranging between 3 and 120 m thick. The bedrock surface was

166 interpolated from available well logs, both in published literature (Randall, 1972) and public records

167 from NYSDEC. The first upper two layers of the model represent the unconfined aquifer system. The

168 third layer is a clay unit, which serves to confine the lowest two layers. The thickness and elevation

169 of the third layer was also interpolated from well logs (Randall, 1972). Both aquifer systems - upper

170 and lower - were split into two layers apiece, with their interlayer elevation set at half of the aquifer

171 thickness in each cell.

172 There are four hydraulic conductivity units in this model (Fig. 4). The uplands are considered one

173 homogeneous, low-conductivity unit, primarily serving as a transmitting media between the external

174 boundary conditions and the valley walls. Separate hydraulic conductivity units were assigned to

175 the upper and lower aquifer systems. Cells representing the clay confining unit were assigned to the

176 fourth conductivity field. Any cell in the third layer with a thickness greater than 3 m is considered

177 part of the confining unit. The remainder of the third layer, where the confining unit is thin or absent,

178 is part of the upper aquifer hydraulic conductivity unit. Manual calibration indicated that this model

179 was not significantly sensitive to conductivity of the confining unit in layer three at the regional

180 scale. Although there is extensive heterogeneity within the valley drift sequences, it is difficult to

181 capture such a variability at this scale. Therefore, these hydraulic conductivity values better represent

182 regional, effective conductivity. Uniform recharge of 62 mm/year was applied to the top of the model,

183 representing the component of groundwater recharge derived from the infiltration of precipitation

184 falling directly in the valleys. This value was approximated by adding the total volume removed from

185 the system (through municipal pumping) to the net regional recharge estimated from the analytic

186 element model (Best, 2013). Constant head boundaries on the outside of the active model area provide

187 the lateral aquifer recharge derived from overland runoff, tributary infiltration, and interflow. In the

188 baseline model, the constant head contribution to groundwater inflow from the boundary of the model

189 was approximately 42%. Constant head contributions in the withdrawal scenarios were evaluated to

190 ensure that this fraction of groundwater input did not unrealistically increase, results of which will be

191 discussed in the sensitivity analysis. The Streamflow-Routing Package (Prudic et al., 2004) was used

192 to simulate stream flow within the model domain. This package allows for the exchange of water

193 between the stream and the aquifer as well as the passage of water between stream cells.

194 The inverse modeling approach required calibration of hydraulic conductivity for each designated

195 field to 53 head targets. These head targets come from a variety of sources including USGS real-time

196 wells and various one-time head measurements from the NYSDEC water well program, consulting

197 reports, field work, and mine data. Assuming isotropy, hydraulic conductivity was varied according

198 to improvements in root mean squared error (RMSE). The final RMSE was 7.08 m with the range

199 of observed water level variability across the model domain from 215.5 m above sea level to 364.7 m

ID Source mE mN Head (m)

WL1 Field measurements 426,660.63 4,667,182.17 253.6

WL2 Field measurements 426,434.91 4,666,666.31 256.0

WL3 Field measurements 416,934.66 4,660,943.6 255.7

WL4 Field measurements 418,158.8 4,661,904.68 267.8

WL5 Field measurements 411,633.68 4,661,205.93 248.8

WL6 Field measurements 411,227.09 4,661,209.18 248.7

WL7 Field measurements 402,462.43 4,660,976.75 244.2

WL8 Field measurements 394,885.12 4,662,320 241.3

WL9 Field measurements 398,314.48 4,661,751.53 243.0

BM128 USGS real-time well 429,550.6 4,671,664.6 269.1

Ti891 USGS real-time well 374,027 4,673,501.1 307.9

TI1500 NYSDEC Water Wei Program 379,116.6198 4,651,726.53 234.4

TI408 NYSDEC Water Wei Program 377,539.1847 4,651,917.109 240.8

TI1442 NYSDEC Water We l Program 378,714.3878 4,653,498.125 215.5

TI300 NYSDEC Water We l Program 381,308.7472 4,655,198.097 262.1

TI440 NYSDEC Water We l Program 373,698.1993 4,655,230.196 254.5

TI397 NYSDEC Water We l Program 373,853.7178 4,655,696.355 275.5

BM1205 NYSDEC Water We l Program 431,376.6484 4,657,754.885 265.2

TI706 NYSDEC Water We l Program 372,239.5567 4,658,640.912 275.8

TI358 NYSDEC Water We l Program 372,318.7358 4,658,827.746 261.4

TI1026 NYSDEC Water We l Program 372,255.4025 4,658,881.316 260.3

BM1250 NYSDEC Water We l Program 432,044.3492 4,659,796.603 266.4

TI1126 NYSDEC Water We l Program 399,368.0667 4,659,976.618 255.7

TI626 NYSDEC Water We l Program 371,790.8604 4,660,521.916 275.5

TI1206 NYSDEC Water We l Program 371,890.7241 4,661,211.19 237.4

TI583 NYSDEC Water We l Program 397,962.8874 4,673,073.791 287.7

TI515 NYSDEC Water We l Program 390,059.3803 4,673,417.952 285.9

TI401 NYSDEC Water We l Program 390,607.7528 4,673,721.047 266.7

TI508 NYSDEC Water We l Program 377,385.8323 4,673,797.097 298.4

TI1015 NYSDEC Water We l Program 381,182.2544 4,674,355.165 299.6

TI528 NYSDEC Water We l Program 376,725.2601 4,674,718.865 293.8

TI1110 NYSDEC Water We l Program 376,586.0455 4,674,887.912 298.7

BM1262 NYSDEC Water We l Program 432,795.8726 4,674,986.885 291.1

TI431 NYSDEC Water We l Program 376,636.2575 4,675,787.928 305.7

TI553 NYSDEC Water We l Program 376,607.6474 4,676,251.153 304.2

BM1579 NYSDEC Water We l Program 430,996.907 4,676,472.546 276.1

TI506 NYSDEC Water We l Program 386,540.1612 4,677,074.217 288.3

TI1125 NYSDEC Water We l Program 402,380.7665 4,682,666.895 328.3

TI306 NYSDEC Water We l Program 385,560.3204 4,683,630.709 364.7

BM1643 NYSDEC Water We l Program 417,510.1229 4,692,695.742 312.7

This table includes model head targets that are publicly available.

200 above sea level. This error was considered acceptable due to the large model extent, the coarse cell size,

201 and the simplified heterogeneity. The resolution that can be expected for any model must be reflective

202 of that model's scale. Secondly, the largest residuals are generally located near external boundaries. The

203 boundary conditions, therefore, are controlling the sensitivity of those targets to changes in hydraulic

204 conductivity. Lastly, this research is only investigating the differences between the baseline model

205 and scenario simulations. Such a comparison requires less certainty in absolute values of the baseline

206 model because the error is linearly transferred to the applied scenario models (Table 1). q3

207 5. Development scenarios

208 While there are some projections of HVHF development in New York (Davis and Robinson, 2012;

209 NYSDEC, 2011), it is difficult to definitively predict well pad density, the particular source water that

210 will be used, and the volume of water required for each pad. This research required the design and

211 testing of a range of potential development scenarios to produce meaningful simulations. These devel-

212 opment scenarios are not predictive but serve as an objective quantification of possible increased water

8 L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

E . I 3 ;i □ Er-::::::::::::::::::::::::: ------^------------------ mm m Ej ^Si F--......

8Ê: ZSÏ:::: PJ)lyi"H ?. ' ■ Ijpffll jfafa- % D 10% ■ ; d ..TO F^ffiffffl fflH W 11«^ ,T \'u vV. y ■ iffïïïïllll mrW IlnJHH -ftfmfh i ft U N Un ffi tüimiirtttftttiftft 11 ! ! 11 M li—Il 1111 ! 11 E ^ IV V -yy ? > ■ >.>.5 ■ ■ ■ ■ ■ ri-" r, « ® n * « J H J. t* 1 ; ■ . ■ . 20% g"» iiSfiS«; JHfiH WÎiÎI S-S flW' ■ «■« ■ 7* 7* 7* * iSN >•?• il.....miliar

Fig. 5. Illustration of density schemes. The exclusion process (A) considered NYSDEC lands (green), open water (dark blue), wetlands (light blue), and urban areas (shades of red). Of the available units (B), different densities were selected using an even distribution. Results for (C) 5%, (D) 10%, and (E) 20% development density scenarios are presented in this paper. Selected units are in red. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

213 demand. Three variables were included in each scenario: well pad density, source of water for each

214 well pad, and volume of water per well pad. Although the time over which water is extracted is in

215 fact an important variable, this research distributes all water withdrawals over an entire year using

216 a steady state modeling assumption. As a result of the steady state assumption, boundary conditions

217 represent the average annual flow that enters, and exits the model domain. This is to avoid the asso-

218 ciated uncertainty with the time variable and the added modeling complexity in introducing model

219 transience.

220 5.1. Well pad density

221 Well pad density is the percentage of land developed for natural gas extraction. For this research,

222 instead of considering the impact of individual wells, well pads - upon which multiple wells may be

223 drilled - are assumed to be the trending mode of development. This document uses "unit" to describe

224 the surface area encompassing both the well pad and the wells' underground horizontal extent. Each

225 unit can have one well pad, again accommodating multiple wells. Because of this distinction, well

226 spacing requirements are not addressed in this configuration. Land use and land coverage are the

227 limiting factors in delineating available land for development (Fig. 5). Regulations currently proposed

228 (NYSDEC, 2013) would limit the density of well pads to no more than one pad per square mile. At each

229 pad as many as 9 horizontal wells would be allowed. Accordingly, the study area was subdivided into

230 a grid of 1-square-mile (2.6 km2) units (Fig. 5A). Any unit that overlaps NYSDEC land was excluded.

231 Units were then further excluded based on the percentage of land which is considered "unavailable",

232 including wetlands, open water, and developed/urban areas. Any unit with greater than 75% unavail-

233 able land was next excluded. Of the remaining units, some percentage was selected to represent the

234 density of development across the modeled extent for that particular scenario. The range of devel-

235 opment density simulated is between 5% and 20%. Selection from the available units was based on

236 a regular distribution scheme that required numbering of the units. The first unit is located in the

L.C. Best, C.S. Lowry / Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

Fig. 6. Allocation of water source based on a Euclidean function. For each unit, the nearest water source was designated based on proximity to the nearest stream for a stream source (A) and nearest municipal pumping center for a municipal groundwater source (B), indicated by yellow circles. Distributed pumping locations were determined in the same manner using distance to the nearest cell in a valley. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

bottom left of the model extent and the numbering continues from left to right and from bottom to top. A10% development density, for example, would use one out of every 10 units in the grid (Fig. 5D).

5.2. Water source

Both groundwater and surface water were considered potential water sources in this research. Groundwater is pumped from either municipal wells or new, privately operated wells, the latter of which will be referred to as the distributed pumping source hereafter. Surface water withdrawals are taken directly from streams. The location of each source, or the point of withdrawal, was determined using a Euclidean allocation function. This function locates the closest straight-line distance from each well pad to each source type (Fig. 6). Every well pad, therefore, has a closest municipal pumping source, distributed pumping source, and stream source. While the closest stream source was selected based on shortest distance, the point of withdrawal was applied at the end of that stream segment at the point of confluence with the next converging stream. A source combination was also included in the scenario runs; this option allowed each well pad to take half of its required water from its designated municipal source and half from its designated stream source. Although it is unlikely that private groundwater wells will be the primary source of HVHF water, this research attempts to simulate a range of water supply options to not only quantify the potential changes but further understand the sensitivity of this hydrologic system to high-volume withdrawals. If groundwater wells are selected to supply water for HVHF, it should be noted that the locations of those wells would likely be chosen based on their anticipated pumping capacity; because the allocation of water source in this project considered any valley cell as a potential pumping location (weighted equally), the pumping capacity and productivity of wide as opposed to narrow valleys was not included. These development scenarios are not intended to predict the potential locations of future groundwater wells.

5.3. Water volume

The volume of water required for each well pad is the product of the number of wells developed on the site and the volume of water each well requires. Between 4 and 9 wells could be accommodated on each well pad based on New York spacing requirements. Approximately 3-4Mgal of water is required for each well according to predicted averages (NYSDEC, 2011); these volumes account for the fraction of injected water which may be derived from the flowback of previously developed wells. In these simulations, between 12 and 32 Mgal of water represents the range of possible water volumes withdrawn for each well pad. This range allows flexibility in the absolute number of wells or volume of water required per well. For example, if 4 wells are developed on a well pad with each using 8 Mgal of water, the maximum water volume in the scenario range is met. If 8 wells are developed on a well

10 L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

269 pad with each using 4Mgal of water, the maximum water volume in the scenario range is likewise

270 met.

271 5.4. Simulation and comparison

272 There are two modes of comparison between the baseline model and the various withdrawal

273 scenarios. The baseline model simply refers to the calibrated MODFLOW model in which current

274 pre-development pumping conditions are at steady-state, while the various withdrawal scenarios are

275 individual models with different pumping/withdrawal conditions applied to each. Pre-development

276 pumping refers only to current rates of groundwater pumping from municipal water supply wells.

277 Any change in the water table will be evaluated in the form of a head difference map - hydraulic head

278 in every model cell in the scenario simulation is subtracted from its counterpart in the baseline model.

279 Every cell in the model domain is therefore attributed a number, with positive values indicating a rise

280 in the water table across that cell and negative values indicating a decline in the water table across

281 that cell. No change to the water table after pumping/withdrawal simulations is interpreted from any

282 zero-value cell in the model domain. Additionally, any cell with a value within 25 cm of zero change

283 was also considered no change due to model variability. The second mode of comparison between

284 the baseline model and the various scenario simulations is the percent change in stream flow. As a

285 result of uniform groundwater recharge under the steady state modeling assumption any change in

286 stream flow under a given development scenario represents the change in groundwater discharge

287 to streams, or base flow. Although surface water modeling would emphasize change to total stream

288 flow, assessing percent change through this technique does not depend absolutely on the accuracy of

289 stream flow in the baseline model. In this way, changes to both the water table and stream flow as a

290 result of either groundwater pumping or stream withdrawal are quantified.

291 6. Results and discussion

292 6.1. Potential effects on the water table

293 Results of this research indicate that changes to groundwater levels are minimal at low develop-

294 ment densities and with low water volumes extracted for each pad. Simulated development scenarios

295 demonstrate locally increasing drawdown with increasing development density at a set volume of

296 water per pad (12Mgal, Fig. 7). In this case, the water used for HVHF is from a combination of both

297 municipal groundwater and stream water. Other models in this research, which simulate withdrawals

298 from distributed pumping wells and streams, mirror the positive relationship between increased

299 development density and drawdown. Assuming the well pad density is constant, increasing the vol-

300 umes of water extracted for each well pad likewise increases drawdown. One of the main differences

301 between the sources, however, is the spatial distribution of withdrawals and the subsequent con-

302 centration or dispersal of water level change (Fig. 8). It is clear that groundwater levels throughout

303 the model domain experience no detectable change from stream withdrawals. Groundwater with-

304 drawals, however, have spatially discrete effects on the water table, while the rest of the model area

305 remains unchanged. The few areas experiencing drawdown in the municipal pumping and combina-

306 tion source scenarios are directly adjacent to municipal pumping wells. With increasing withdrawal,

307 the cones of depression at municipal wells in narrow glacial valleys are both expanded and deep-

308 ened (Fig. 7, locations I—III). Municipal wells located in the widest glacial valleys and near major

309 rivers, particularly the Susquehanna River, do not experience the same impact (Fig. 7, location IV). The

310 municipal pumping and combination source scenarios produce the same spatial distribution of water

311 table change although there is a difference in the magnitude of change. A distributed pumping source

312 evokes the most widespread drawdown although the extent of drawdown is still limited to narrow

313 valleys (Fig. 8D).

314 Groundwater levels are relatively insensitive to increased water withdrawals although there are

315 two exceptions. First, greater cones of depression are notable around municipal wells when pumping

316 rates increase (Fig. 7). When the burden of water source is instead split between streams and municipal

317 wells, the effect on the water table is lessened. Vulnerable municipal wells appear to be associated

L.C. Best, C.S. Lowry / Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx 11

]0.26—0.50 ] No change

Fig. 7. Simulated drawdowns near municipal wells, assuming half the water use at each HVHF pad is pumped from municipal wells and half from streams. Although water volume per pad is constant at 12Mgal, development density increases A to C. Note locations labeled I—IV, referred to in text.

L.C. Best, C.S. Lowry/ Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

Fig. 8. Head difference maps for 10% development density scenarios with 12 Mgal of water extracted for each well pad. Each frame displays results of using different water sources: (A) combination of municipal pumping and stream withdrawal, (B) stream water only, (C) municipal pumping only, and (D) distributed pumping.

318 with narrow valleys (Fig. 8C). This may be a result of aquifer geometry, area of contributing recharge,

319 and availability of induced recharge from streams. Aquifer geometry refers to both the width and

320 depth of glacial valley fill. The pumping center near Binghamton, NY (Fig. 7, location IV) is an example

321 of a region within the valley aquifer that has municipal wells with the capacity to accommodate

322 the increased pumping rate. These wells are located in a wide valley with thick aquifer deposits.

323 Additionally, proximity to a large stream allows the possibility for greater induced recharge from

324 the stream. The second susceptibility occurs under distributing pumping conditions, during which

325 significant reductions in groundwater elevations are apparent in narrow valleys (Fig. 8D). Again, this

326 is most likely associated with the aquifer geometry and area of contributing recharge.

327 As demonstrated in Fig. 7, increases in both development density and water volume per pad elicit

328 heightened water table responses; this trend was shared by all sources. Although water table change

329 was still undetectable for stream withdrawals at the maximum development tested, heightened res-

330 olution and smaller scale models might allow for better understanding of the connection between

331 streams and groundwater.

332 6.2. Potential effects on stream flow

333 Changes to stream flow in response to high-volume water withdrawals are spatially variable. The

334 most significant reduction to stream flow is concentrated in one region of the model (Fig. 9, cross-

335 sections 7, 8, and 9). Other areas of the model respond relatively uniformly to extraction scenarios,

336 with the percent reduction in stream flow increasing with increasing development density and water

337 volume per pad. Within the minimum development range, extracting water from both municipal

338 pumping wells and streams reduces stream flow by less than 2% throughout most of the stream net-

339 work (Fig. 9A). At the maximum density of development, stream flow is reduced by up to 13% in a

340 localized region (Fig. 9D). Under those same development conditions, however, stream flow reduction

341 still remains under 3% throughout most of the stream network. Although the magnitude of stream flow

L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

Fig. 9. Percentage of stream flow reduction at points along the stream network. The results shown are from scenarios extracting from a combination of municipal pumping wells and streams. The polygon (B) delineates the aquifer network with the numbers indicating the locations at which stream flow was compared (referred to as stream cross-section). Variations in both development density and water volume per pad are given in each bar graph.

342 reduction changes based on water source, the general spatial distribution persists (Fig. 10). Streams

343 throughout the model respond consistently to applied withdrawal scenarios with the exception of

344 stream cross-sections 7, 8, and 9, which exhibit nearly three times the stream flow reduction as com-

345 pared to the rest of the stream segments. The combination source and stream withdrawals produced

346 the greatest response in stream flow whereas distributed pumping scenario results in a less dramatic

347 response (Fig. 10). Extracting from municipal wells causes more spatial variability in stream flow

348 reduction as compared to the combination source (Fig. 10, cross-section 8).

349 There is a positive relationship between stream flow reduction and volume of extracted water

350 which is determined by both well pad density and water volume per pad. Relatively uniform response

351 throughout most of the stream segments emphasizes the markedly greater response at cross-sections

352 7, 8, and 9 (Fig. 9). These locations are in narrow valleys and represent streams with lesser annual

353 discharge. These two factors dictate the capacity of groundwater-surface water exchange when with-

354 drawals from either the aquifer or the streams are applied. Downstream parts of the stream network

355 (Fig. 9, cross-sections 10-13) demonstrate slightly greater sensitivity to combination source with-

356 drawals than the upstream portion of the model area. This potentially demonstrates how effects from

357 a tributary stream might propagate to the main stream with which it converges. The actual magni-

358 tude of stream flow reduction in a high water use scenario may be considered significant only during

359 drought conditions. This underscores the importance of understanding the implications of withdrawal

360 timing and duration on potentially vulnerable valleys. Incorporation of model transience would help

361 address this uncertainty.

362 The spatial distribution of changes to stream flow is consistent between sources, with the exception

363 of the municipal pumping scenarios (Fig. 10, cross-section 8). This location exemplifies an instance

L.C. Best, C.S. Lowry / Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

Fig. 10. Stream flow changes as a result of different source water. Stream cross-section refers to the reference map (Fig. 7B). Density of development is a variable in each graph while water volume per pad is constant at 12 Mgal per well pad.

364 of "shared response" between stream flow and the water table. At this location, the municipal cone

365 of depression is greatest when water is taken only from the municipal well while the stream flow

366 reduction is comparably small. When the burden of water source is shared with withdrawals from the

367 nearby stream, the water table impact is alleviated (Fig. 8A) while the stream flow reduction intensifies

368 (Fig. 10). Intuitively, stream flow is reduced most when water is taken only from the streams. Results

369 demonstrate that the water table is insensitive to stream withdrawals (Fig. 8). It can be inferred that

370 stream-aquifer connectivity distributes the stream withdrawals over a larger area than concentrated

371 pumping schemes, thus resulting in insignificant drawdown. Only when municipal pumping is added

372 (Fig. 10A) water table and stream flow changes simultaneously emerge. Distributed pumping has the

373 least effect on stream flow because of the distribution of water burden. Many low-capacity wells draw

374 uniformly less from overlying streams than fewer high-capacity wells. If stream flow protection is

375 prioritized based on suggested vulnerability, it is important to note that a distributed pumping source

376 causes the least reductions to stream flow.

377 6.3. Sensitivity

378 There are two aspects of this model that are significant in dictating model results: the volume

379 of water input to the system as a result of aquifer recharge and the connectivity of the aquifer and

380 overlying streams as a result of streambed conductance. In order to determine the impacts of these

381 parameters a sensitivity analysis was conducted. The greatest uncertainty in this model is the value

382 estimated for applied recharge, which is associated with infiltration of direct valley precipitation.

383 Recharge is the main parameter that governs how much water is available to the system. Increasing

384 recharge decreases the percent reduction in stream flow, mainly in areas of the stream network that

385 experience the greatest change (Fig. 11A, cross-sections 7-9). As expected, the greatest reduction to

386 streamflow is identified under zero-recharge, or severe drought, conditions.

L.C. Best, C.S. Lowry / Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx 15

4.0 3.0 2.0 1.0 0.0

4.0 3.0 2.0 1.0 0.0

1 2 3 4 5 6 7 8 9 10 11 12 13

Fig. 11. Sensitivity of model results to (A) groundwater recharge and (B) hydraulic conductivity of streambed sediments. The base case parameters are shown in black. Density of development is constant at 10% and water volume per pad is constant at 12 Mgal per well pad. The simulated water source is the municipal pumping and surface water combination scenario.

387 The hydraulic connectivity between surface water and groundwater is primarily controlled by

388 streambed conductance. Hydraulic conductivity of streambed sediments is one of the variables that

389 determines streambed conductance. Similar to recharge sensitivity, increasing the streambed sedi-

390 ment conductivity reduces the changes to stream flow (Fig. 11B). Again, this sensitivity is generally

391 apparent at stream segments which experienced the greatest change. It is crucial for water resource

392 management analyses to consider the range of results possible given the sensitivity of results to a

393 particular model feature.

394 In addition to the evaluation of model sensitivities to the variability in aquifer recharge and

395 streambed conductance, the impact of specified head boundary conditions was evaluated. The model

396 mass balance was analyzed to determine whether constant head contributions to groundwater input

397 would change under withdrawal scenarios. The input volume from the constant head boundary con-

398 ditions increased by less than 1% for each of the source scenarios at maximum development, with the

399 exception of the distributed pumping case. Distributed pumping induced a 9% increase in the constant

400 head input volume. This volume is less than the applied recharge, which supports the use of constant

401 head boundary conditions at the edge of the model domain. Mass balance results demonstrate that

402 these boundary conditions do not supply unrealistic volumes of water to the aquifer under increased

403 pumping conditions.

404 7. Conclusions

405 Although regions that are water-rich encounter fewer water quantity issues as compared to arid

406 regions, possible implications of energy development and subsequent water demands must be con-

407 sidered. This is particularly applicable in areas that have barriers - legal, physical, or economic - to

408 alternate sources of drinking water so both the quality and sustainable supply of existing sources

409 must be safeguarded. Simulating water table and stream flow response to high-volume water with-

410 drawal scenarios is effective in quantifying the potential impacts of increased water demand associated

411 with HVHF expansion into New York State. This research emphasized a regional perspective to first

412 determine whether changes to the water table and/or stream flow could be detected under potential

413 development scenarios. Identification of high-impact scenarios and susceptible model areas demon-

414 strates the utility of regional groundwater flow modeling in assessing a water quantity concern.

415 The range of development scenarios modeled depict impacts to water resources that are most

416 pronounced at municipal pumping centers and along narrow tributary valleys. Cones of depression

417 would deepen around municipal pumping wells, if postulated HVHF water needs were withdrawn

418 partially or entirely from those wells. Additional drawdown around municipal wells in wide valleys

L.C. Best, C.S. Lowry/Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx

419 would be negligible. Significant drawdown is simulated in narrow tributary valleys under pumping

420 scenarios that call for HVHF withdrawals from new private wells at valley sites closest to postulated

421 gas wells. Results demonstrate the capacity for increased pumping is constrained by the contribut-

422 ing recharge area and aquifer geometry. Furthermore, there is diminished opportunity for induced

423 recharge in streams within these narrow valleys. At these locations, distributed pumping wells would

424 draw more water from the aquifer than could be replenished by groundwater recharge.

425 It is important to recognize that both groundwater pumping and stream withdrawals have an

426 impact on stream discharge. The greatest stream flow reductions were geographically limited to a

427 particular section of the stream network (Fig. 9, cross-sections 7-9). Valley width appears to be the

428 limiting factor in determining the magnitude of stream flow reduction. Some reductions were detected

429 on larger streams at locations downstream from those particular cross-sections. As a result of the high

430 hydraulic connectivity between the streams and underlying aquifer, water resource management

431 decisions pertaining to HVHF water demands should fully represent the freshwater system as a single

432 resource.

433 To best understand changes to cones of depression around municipal pumping centers or nearby

434 stream discharge changes, localized fine-scale models are optimal. Furthermore, transient models

435 would allow quantification of variable withdrawal timing and duration. This research presents a

436 necessary foundation for analyzing water resources at a regional scale with the understanding that

437 individual applications would require further high-resolution analysis. Planning and regulation of

438 HVHF will ultimately encounter water permitting decisions. These decisions should conservatively

439 consider the hydraulically connected groundwater-surface water systems which exhibit spatially

440 distributed sensitivities to high-volume withdrawals.

441 Acknowledgments

44Q4 Funding for this project was supported by the Mark Diamond Research Foundation and the Depart-

44Q5 ment of Geology Champion Fund, University at Buffalo. Special thanks to Gary Priscott and Lucas

444 Mahoney from the NYSDEC as well as both Broome and Tioga counties' Department of Health for

445 access to municipal pumping records.

446 References

447 Aber, J.S., 1980. Nature and origin of exotic-rich drift in the Appalachian Plateau. Am. J. Sci. 280,363-384.

448 Best, L.C., (M.S. thesis) 2013. Modeling the Regional Groundwater Flow in Broome and Tioga Counties, NY under High-Volume

449 Water Extraction Scenarios Driven by Potential Natural Gas Development of the Marcellus Shale. State University of New

450 York at Buffalo, Buffalo, 158 p.

45q6 Cadwell, D.H., 1973. Glacial geology of the Chenango River valley near Binghamton, New York. Publ. Geomorphol. Contrib. (3),

452 31-39.

453 Coon, W.F., Yager, R.M., Surface, J.M., Randall, A.D., Eckhardt, D.A., 1998. Hydrogeology and Water Quality ofthe Clinton Street-

454 Ballpark Aquifer near Johnson City, New York. US Geological Survey Open-File Report 97-102., pp. 1-63.

455 Craig, J.R., Matott, L.S., 2009. Visual AEM, (accessed March 2012). 45q7 Davis, J.B., Robinson, G.R., 2012. Ageographic model to assess and limit cumulative ecological degradation from Marcellus Shale

457 exploitation in New York, USA. Ecol. Soc. 17(2).

458 Fleisher, P.J., 1986. Glacial geology and late Wisconsinan stratigraphy, upper Susquehanna drainage basin, New York. Bull. -

459 New York State Museum (1976) 455,121-142.

460 Gregory, K.B., Vidic, R.D., Dzombak, D.A., 2011. Water management challenges associated with the production of shale gas by

461 hydraulic fracturing. Elements 7 (3), 181-186.

462 Harbaugh, A.W., 2005. M0DFL0W-2005. The U.S. Geological Survey Modular Ground-Water Model - The Ground-Water Flow

463 Process. U.S. Geological Survey Techniques and Methods 6-A16.

464 Hunt, R.J., Anderson, M.P., Kelson, V.A., 1998. Improving a complex finite-difference ground water flow model through the use

465 of an analytic element screening model. Ground Water 36 (6), 1011-1017.

466 Jankovic, I., Barnes, R., 1999. High-order line elements in modeling two-dimensional groundwater flow. J. Hydrol. 226 (3-4),

467 211-223.

468 Kontis, A.L., Randall, A.D., Mazzaferro, D.L., 2004. Regional Hydrology and Simulation of Flow in Stratified-Drift Aquifers in the

469 Glaciated Northeastern United States. US Geological Survey Professional Paper (1415 C)., pp. C1-C156.

470 MacNish, R.D., Randall, A.D., 1982. Stratified-Drift Aquifers in the Susquehanna River Basin. In: New York State Department of

471 Environmental Conservation Bulletin 75, New York, 68 p.

472 Maloney, K.O., Yoxtheimer, D.A., 2012. Research article: production and disposal of waste materials from gas and oil extraction

473 from the Marcellus Shale Play in Pennsylvania. Environ. Pract. 14 (4), 278-287.

474 Nicot, J.-P., Scanlon, B.R., 2012. Water use for shale-gas production in Texas, U.S. Environ. Sci. Technol. 46 (6), 3580-3586.

L.C. Best, C.S. Lowry / Journal of Hydrology: Regional Studies xxx (2014)xxx-xxx 17

475 NYSDEC, 2011. Supplemental Generic Environmental Impact Statement on the Oil, Gas and solution Mining Regulatory Program:

476 Well Permit Issuance for Horizontal Drilling And High-Volume Hydraulic Fracturing to Develop the Marcellus Shale and

477 Other Low-Permeability Gas Reservoirs. Revised Draft [Online]. URL:

478 NYSDEC, 2013. Statewide Spacing Unit Sizes and Setbacks, (accessed February 2013).

479 Prudic, D.E., Konikow, L.F., Banta, E.R., 2004. A New Streamflow-Routing(SFR1) Package to Simulate Stream-Aquifer Interaction

480 with M0DFL0W-2000. U.S. Geological Survey Open-File Report 2004-1042,95 p.

481 Rahm, B.G., Riha, S.J., 2012. Toward strategic management of shale gas development: regional, collective impacts on water

482 resources. Environ. Sci. Policy 17,12-23.

483 Randall, A.D., 1972. Records ofWells and Test Borings in the Susquehanna River Basin. New York State Department of Environ-

484 mental Conservation Bulletin 69, New York, 92 p.

485 Randall, A.D., 1977. The Clinton Street-Ballpark Aquifer in Binghamton and Johnson City. New York State Department of

486 Environmental Conservation Bulletin 73, New York, 71 p.

487 Randall, A.D., 1978. A contribution to the late Pleistocene stratigraphy of the Susquehanna River valley of New York. Empire

488 State Geogram 14 (2), 2-15.

489 Randall, A.D., 1986. Aquifer Model of the Susquehanna River Valley in Southwestern Broome County, New York. U.S. Geological

490 Survey Water-Resources Investigations Report 85-4099,38 p.

491 Randall, A.D., 2001. Hydrogeologic Framework of Stratified-Drift Aquifers in the Glaciated Northeastern United States. U.S.

492 Geological Survey Professional Paper (1415-B)., pp. B1-B179.

493 Randall, A.D., 2010. Low Flow of Streams in the Susquehanna River Basin of New York. U.S. Geological Survey Scientific Inves-

494 tigations Report 2010-5063,57 p. At:

495 Randall, A.D., Snavely, D.S., Holecek, T.P., Waller, R.M., 1988. Alternative Sources of Large Seasonal Ground-Water Supplies in

496 the Headwaters of the Susquehanna River Basin, New York. U.S. Geological Survey Water-Resources Investigations, 121 p.

497 Rassenfoss, S., 2011. From flowbackto fracturing; water recycling grows in the Marcellus Shale. J. Petrol. Technol. 63 (7), 48-51.

498 Reynolds, R.J., Williams, J.H., 1988. Continuous seismic-reflection profiling ofglacial drift along the Susquehanna, Chemung, and

499 Chenango Rivers, South-Central New York and North-Central Pennsylvania. Regional aquifer systems of the United States:

500 the northeast glacial aquifers. Am. Water Resour. Assoc. Monogr. Series 11, 83-103.

501 Rumbaugh, J.O., Rumbaugh, D., 2011. Guide to Using Groundwater Vistas Version 6. Environmental Simulations, Inc.

502 Scully, R.W., Arnold, R.W., 1981. Holocene alluvial stratigraphy in the upper Susquehanna River basin, New York. Quaternary

503 Res. (New York) 15(3), 327-344.

504 Soeder, D.J., 2010. The Marcellus Shale: resources and reservations. Eos 91 (32), 277-288.

505 U.S. Department of Energy, 2009. Modern Shale Gas Development in the United States: A Primer, 116 p.

506 U.S. Energy Information Administration, 2011. Review of Emerging Resources: U.S. Shale Gas and Shale Oil Plays, 105 p.

507 U.S. Environmental Protection Agency, 2010. Clinton Street Ballpark Aquifer System,

508 water/aquifer/clinton/clinton.htm (accessed March 2012).

509 Waller, R.M., Finch, A.J., 1982. An Atlas of 11 Selected Aquifers in New York. U.S. Geological Survey Open-File Report 82-0553,

510 247 p.

511 Winter, T.C., Harvey, J.W., Franke, O.L., Alley, W.M., 1998. Ground Water and Surface Water: A Single Resource. U.S. Geological

512 Survey Circular 1139, 79 p.

513 Wolcott, S.W., Coon, W.F., 2001. Simulation of a Valley-Fill Aquifer System to Delineate Flow Paths, Contributing Areas, and

514 Traveltime to Wellfields in Southwestern Broome County, New York. Water-Resources Investigations - U.S. Geological

515 Survey, Report 01-4171,30 p.

516 Yager, R.M., 1986. Simulation of Ground-Water Flow and Infiltration from the Susquehanna River to a Shallow Aquifer at

517 Kirkwoodand Conklin, Broome County, New York. U.S. Geological Survey-Water-Resources Investigations Report 86-4123.

518 Yager, R.M., 1993. Estimation of Hydraulic Conductivity of a Riverbed and Aquifer System on the Susquehanna River in Broome

519 County, New York. US Geological Survey Water-Supply Paper 2387.