Updates to the Regional Groundwater-Flow Model of the New Jersey Coastal Plain, 1980–2013

Scientific Investigations Report 2023-5066
Prepared in cooperation with the New Jersey Department of Environmental Protection
By:  and 

Links

Acknowledgments

This update of a previously published regional groundwater-flow model of the New Jersey Coastal Plain is derived directly from the published work of Mary Martin Chepiga (Martin, 1998) and Lois Voronin (Voronin, 2004) of the U.S. Geological Survey (USGS); the authors are grateful for their insights and their assistance in developing the updated model. Debra Buxton (USGS, retired) completed most of the work to build the three-dimensional hydrogeologic framework for the groundwater-flow model. Liam Kenefic (USGS) assisted with revisions of the framework. Technical reviews of the report and model by Eric Roman (New Jersey Geological and Water Survey) and Leon Kauffman and Paul Misut (USGS) improved the accuracy of the report. The verification by Wesley Zell (USGS) of the model data release accompanying this report is appreciated. The authors also appreciate the additional technical input from members of the New Jersey Geological and Water Survey, particularly Jeffrey Hoffman (New Jersey State Geologist) and Steven Domber.

Abstract

A 21-layer three-dimensional transient groundwater-flow model of the New Jersey Coastal Plain was developed and calibrated by the U.S. Geological Survey (USGS) in cooperation with the New Jersey Department of Environmental Protection to simulate groundwater-flow conditions during 1980–2013, incorporating average annual groundwater withdrawals and average annual groundwater recharge. This model is the third version of the New Jersey Coastal Plain regional groundwater-flow model that was initially developed as part of the USGS Regional Aquifer System Analysis (RASA) program. The model simulates groundwater flow in 11 aquifers and 10 intervening confining units of the New Jersey Coastal Plain to provide a regional overview of groundwater conditions. Averaged groundwater withdrawal data for 1980 to 2013 were used in the model. The 11 aquifers in New Jersey are, from shallowest to deepest, the Holly Beach water-bearing zone and the confined Cohansey aquifer in Cape May County; the Rio Grande water-bearing zone; the Atlantic City 800-foot sand; the Piney Point, Vincentown, and Wenonah-Mount Laurel aquifers; the Englishtown aquifer system; and the upper, middle, and lower aquifers of the Potomac-Raritan-Magothy (PRM) aquifer system.

The model was developed with the MODFLOW–2005 numerical code and the UCODE parameter estimation technique and calibrated using water-level and base-flow observations. A total of 3,453 water-level observations from 392 wells in New Jersey and 48 wells in Delaware from 1983 to 2013 were used in model calibration, which includes historical water-level trends for 29 wells in New Jersey during 1980–2013 presented in time-series hydrographs. In addition, derived observations also were included by calculating the vertical gradient at 33 pairs of nested observation wells in New Jersey, for a total of 210 observations. Changes in water levels over time were calculated for 134 wells in New Jersey and four wells in Delaware where water levels had varied substantially (approximately 10 ft) over the 30-year span of synoptic water-level measurements, for a total of 767 observations. A total of 1,485 base-flow observations in 47 surface-water basins in New Jersey from 1980 to 2013 were used in model calibration.

Updates to the groundwater-flow model include the conversion to a fully three-dimensional model from the previous quasi-three-dimensional model. The new model will allow for potential future uses such as particle tracking or simulation of variable-density groundwater flow that could not be accomplished with earlier versions of the model. Spatially and temporally variable recharge estimated by using a soil-water balance model resulted in a spatially and temporally finer discretization. The Rio Grande water-bearing zone was added to the model as an aquifer layer to refine estimates of simulated flow in Atlantic and Cape May Counties, New Jersey. Hydrogeologic parameters were updated to include the confining units in New Jersey and corresponding hydrogeologic units in Delaware and eastern Maryland.

The simulated water levels for the New Jersey Coastal Plain aquifers were compared to water-level measurements made during 1980–2013. The average residual for 4,243 water-level observations for New Jersey (simulated water levels minus measured water levels) is 1.5 feet. The simulated water-level contours for the confined aquifers for 2013 were compared to potentiometric surfaces produced from water levels measured during 2013. Simulated water levels generally matched the 2013 potentiometric surfaces of the confined aquifers in the areas of large withdrawals. Hydrographs of wells in the confined Coastal Plain aquifers of New Jersey show that simulated water levels generally match the magnitude and seasonal variation of the observed water levels. Hydrographs of base-flow for the 47 streamgaging stations in New Jersey indicate that most of the simulated and estimated data match reasonably well.

Groundwater withdrawals are an important resource for water supply, agricultural, industrial, and commercial needs in the New Jersey Coastal Plain. Groundwater withdrawals from the New Jersey Coastal Plain aquifers have resulted in persistent, regionally extensive cones of depression in the Englishtown aquifer system and Wenonah-Mount Laurel aquifer in Ocean and Monmouth Counties; Wenonah-Mount Laurel and upper, middle, and lower PRM aquifers in Camden County; and Atlantic City 800-foot sand in Atlantic County. Because hydrologic stresses and water-management needs change with time, periodic updates to the groundwater-flow model are required to provide current information about hydrologic conditions in the New Jersey Coastal Plain and to maintain its usefulness as a tool to manage water resources and develop water-resource strategies. The current updates will support the continued application of this model as a tool for evaluating the regional effects of changes in groundwater withdrawals and of current and potential future water-management strategies on groundwater levels in the New Jersey Coastal Plain.

Introduction

The U.S Geological Survey (USGS) Regional Aquifer System Analysis (RASA) program was begun in 1978 in response to a Congressional mandate to develop quantitative appraisals of the major groundwater systems of the United States (Martin, 1998). The initial New Jersey Coastal Plain groundwater-flow model, constructed by Martin (1998) as part of the New Jersey RASA program, was based on the hydrogeologic framework of the New Jersey Coastal Plain developed by Zapecza (1989). Martin (1998) adapted the hydrogeologic framework presented in Zapecza (1989) for the New Jersey RASA model and represented the Coastal Plain sediments as 10 aquifers and 9 intervening confining units. The RASA model of Martin (1998) was updated by Voronin (2004), who (1) rediscretized the model parameters with a finer cell size, (2) rediscretized stream cells to more accurately represent the streams, (3) used a spatially variable recharge rate based on results of USGS studies conducted in New Jersey since Martin (1998), and (4) updated groundwater-withdrawal data to 1998.

Since its original development, the New Jersey Coastal Plain groundwater-flow model has been used to evaluate the regional effects of groundwater withdrawals on water levels in the confined aquifers in the New Jersey Coastal Plain and as a tool to assess water-management strategies for the region. Studies conducted in New Jersey that have used the model to address water-supply management issues, depletion of streamflow from increased groundwater withdrawals, and continued declining water levels include Navoy (1994), Pope (2006), Watt and Voronin (2006), and Gordon (2007). In some studies (for example, Spitz and others, 2008; Spitz and DePaul, 2008), the regional effect of reducing or increasing groundwater withdrawals on water levels in Water Supply Critical Areas 1 and 2—two areas designated by the New Jersey Department of Environmental Protection (NJDEP) where excessive water use threatens the long-term sustainability of the water supply (fig. 1)—were evaluated. Critical Area 1, designated in 1985, encompasses parts of Middlesex, Monmouth, and Ocean Counties. Regulated aquifers in Critical Area 1 are, in order of increasing depth, the Wenonah-Mount Laurel aquifer, the Englishtown aquifer system, and the upper and middle Potomac-Raritan-Magothy (PRM) aquifers. To improve the management of groundwater resources of the PRM aquifer system in southwestern New Jersey, Critical Area 2 was designated in 1993. This management area encompasses Camden, most of Burlington and Gloucester, and parts of Atlantic, Cumberland, and Salem Counties. Regulated aquifers in Critical Area 2 are, in order of increasing depth, the upper, middle, and lower PRM aquifers.

Map with various color lines showing critical areas and cross sections, large purple
                     box showing model grid outline.
Figure 1.

Map showing A, location of the study area, B, generalized hydrogeologic section A–A’ through the Coastal Plain of southern New Jersey, and C, generalized hydrogeologic section B–B’ through southern Cape May County, New Jersey.

To maintain the usefulness of the New Jersey Coastal Plain groundwater-flow model as a tool for managing water resources and developing water-resource strategies in light of changing hydrologic conditions, the USGS, in cooperation with the NJDEP, revised the model by updating the hydrogeologic framework, hydraulic parameters, and groundwater withdrawals. The updated hydrogeologic framework includes a three-dimensional framework of the New Jersey Coastal Plain that extends into Delaware and parts of eastern Maryland. These updates may result in (1) refined boundary flows for local models, (2) improved estimates of the simulated interaction between water levels in southwestern New Jersey and withdrawals in Delaware, (3) simulation of the unconfined Kirkwood-Cohansey aquifer system and confined Rio Grande water-bearing zone as model layers, and (4) the ability to conduct particle tracking to determine regional sources of flow and times of travel.

Purpose and Scope

This report describes the development and calibration of a numerical groundwater-flow model of the New Jersey Coastal Plain that extends into Delaware and parts of eastern Maryland. The model simulates transient regional groundwater conditions incorporating annual groundwater withdrawals and recharge from 1980 to 2013 and will be used to assess changes in water levels over the simulated period. The fully three-dimensional flow model described in this report is a 21-layer model that represents 11 aquifers and 10 intervening confining units in the New Jersey Coastal Plain. The groundwater-flow model is based on the USGS New Jersey RASA model (Martin, 1998), which was rediscretized in the late 1990s by Voronin (2004). These two previous quasi-three-dimensional models represented 10 aquifers in the New Jersey Coastal Plain but did not explicitly include the intervening confining units or the Rio Grande water-bearing zone. The confining units were previously represented as leakage by using a vertical leakance parameter. The model simulates the groundwater-flow system in the New Jersey Coastal Plain and offshore sediments that contain water with a chloride concentration less than 10,000 milligrams per liter (mg/L). Model design and input data, boundary conditions, and groundwater withdrawals used in the model calibration are described. The hydrogeologic parameters for each aquifer and confining unit, including the correlated units in Delaware and Maryland, also are discussed. A total of 3,453 annual water-level observations from 392 wells in New Jersey and 48 wells in Delaware from 1983 to 2013 were used in model calibration. Hydrographs of simulated and observed water levels in 29 observation wells and simulated and interpreted potentiometric surfaces for 2013 conditions are shown. Hydrographs of simulated and estimated annual base flows in 47 surface-water basins are presented. The calibration of the model, conducted by trail-and-error adjustment and by using UCODE (Poeter and others, 2014), is presented. Model parameter sensitivities and limitations are discussed. Estimates of recharge in New Jersey, applied to the groundwater-flow model by using a soil-water balance (SWB) code by Westenbroek and others (2010), is discussed in the appendix. The model archive for the groundwater-flow model described in this report contains the input and output files and the executable file needed to run the model (Carleton and others, 2023). Digital data describing the extents and thicknesses of the units of the updated hydrogeologic framework can be found in Carleton and others (2023).

Location and Extent of the Model Area

The model area encompasses the Coastal Plain Physiographic Province of New Jersey and a small part of southeastern Pennsylvania, and parts of the Coastal Plain in Delaware and eastern Maryland (fig. 1). The model area covers approximately 16,537 square miles (mi2) and is bounded on the north and northwest by the updip limit of the Coastal Plain sediments (Fall Line), which separates the Coastal Plain sediments from the Piedmont Physiographic Province sediments; on the east and southeast by the Atlantic Ocean; and on the southwest by Delaware Bay, most of Delaware, and part of eastern Maryland. This investigation focused on the Counties of Atlantic, Burlington, Camden, Cape May, Cumberland, Gloucester, Monmouth, Ocean, Salem, and parts of Mercer and Middlesex in New Jersey, and includes New Castle, Kent, and most of Sussex Counties in Delaware and parts of Cecil, Kent, Queen Anne’s, Caroline, Wicomico, and Worcester Counties in Maryland. Topography in New Jersey is relatively flat; altitudes range from sea level along estuaries, bays, and the Atlantic Ocean to nearly 400 feet (ft) in northeastern Monmouth County, New Jersey.

Previous Investigations

Martin (1998), Pope and Gordon (1999), and Voronin (2004) describe simulated groundwater flow in the New Jersey Coastal Plain from a regional perspective. These three studies, along with other studies describing regional simulated flow in the New Jersey Coastal Plain or the North Atlantic Coastal Plain, are listed in table 1. Table 1 also includes regional studies completed in Critical Areas 1 and 2 and in the Atlantic City area of the New Jersey Coastal Plain.

Table 1.    

Previous regional simulation studies conducted for the New Jersey Coastal Plain and Northern Atlantic Coastal Plain.

[Northern Atlantic Coastal Plain includes Delaware and parts of Maryland, New Jersey, New York, North Carolina, and Virginia; Water Supply Critical Area 1 includes parts of Middlesex, Monmouth, and Ocean Counties in New Jersey; Water Supply Critical Area 2 includes Camden and parts of Atlantic, Burlington, Cumberland, Gloucester, Ocean and Salem Counties in New Jersey]

Area of study Simulated aquifers Report name Reference
Northern Atlantic Coastal Plain Ten regional North Atlantic Coastal Plain aquifers Geohydrology and simulation of ground-water flow in the Northern Atlantic Coastal Plain aquifer system Leahy and Martin, 1993
Water Supply Critical Area 1 in New Jersey Upper and middle Potomac-Raritan-Magothy aquifers, Englishtown aquifer system, and Wenonah-Mount Laurel aquifer Hydrogeology, simulation of regional ground-water flow, and saltwater intrusion, Potomac-Raritan-Magothy aquifer system, northern Coastal Plain of New Jersey Pucci and others, 1994
Water Supply Critical Area 2 in New Jersey Upper, middle, and lower Potomac-Raritan-Magothy aquifers Ground-water flow and future conditions in the Potomac-Raritan-Magothy Aquifer system, Camden area, New Jersey Navoy and Carleton, 1995
New Jersey Coastal Plain Ten regional New Jersey Coastal Plain aquifers Ground-water flow in the New Jersey Coastal Plain Martin, 1998
New Jersey Coastal Plain Ten regional New Jersey Coastal Plain aquifers Simulation of groundwater flow and movement of the freshwater-saltwater interface in the New Jersey Coastal Plain Pope and Gordon, 1999
Cape May and parts of Atlantic and Ocean Counties in New Jersey Atlantic City 800-foot sand Ground-water flow and quality in the Atlantic City 800-foot sand, New Jersey McAuley and others, 2001
New Jersey Coastal Plain Ten regional New Jersey Coastal Plain aquifers Documentation of revisions to the Regional Aquifer System Analysis model of the New Jersey Coastal Plain Voronin, 2004
Water Supply Critical Area 2 in New Jersey Upper and middle Potomac-Raritan-Magothy aquifers, Englishtown aquifer system, and Wenonah-Mount Laurel aquifer Recovery of ground-water levels from 1988 to 2003 and analysis of 2003 and full-allocation withdrawals in Critical Area 2, Southern New Jersey Spitz and DePaul, 2008
Water Supply Critical Area 1 in New Jersey Upper and middle Potomac-Raritan-Magothy aquifers, Englishtown aquifer system, and Wenonah-Mount Laurel aquifer Recovery of ground-water levels from 1988 to 2003 and analysis of potential water-supply management options in Critical Area 1, east-central New Jersey Spitz and others, 2008
Salem and Gloucester Counties in New Jersey Upper, middle and lower Potomac-Raritan-Magothy aquifers, and Englishtown aquifer system Simulated effects of allocated and projected 2025 withdrawals from the Potomac-Raritan-Magothy aquifer system, Gloucester and northeastern Salem Counties, New Jersey Charles and others, 2011
Atlantic, Burlington, Gloucester, Camden, Ocean, and Cape May Counties in New Jersey Kirkwood-Cohansey aquifer system, Rio Grande water-bearing zone, and Atlantic City 800-foot sand Simulated effects of alternative withdrawal strategies on groundwater flow in the Kirkwood-Cohansey aquifer system, the Rio Grande water-bearing zone, and the Atlantic City 800-foot sand in the Great Egg Harbor and Mullica River Basins, New Jersey Pope and others, 2012
Ocean and southern Monmouth Counties in New Jersey Kirkwood-Cohansey aquifer system, Rio Grande water- bearing zone, Atlantic City 800-foot sand, Piney Point and Vincentown aquifers Simulated effects of groundwater withdrawals from aquifers in Ocean County and vicinity, New Jersey Cauller and others, 2016
Northern Atlantic Coastal Plain 19 regional North Atlantic Coastal Plain aquifers Documentation of a groundwater flow model developed to assess groundwater availability in the Northern Atlantic Coastal Plain aquifer system from Long Island, New York, to North Carolina Masterson and others, 2016a.
Table 1.    Previous regional simulation studies conducted for the New Jersey Coastal Plain and Northern Atlantic Coastal Plain.

In addition to the previous regional simulation studies listed in table 1, the rediscretized RASA model of the New Jersey Coastal Plain (Voronin, 2004) was used by Pope (2006) to simulate the effects of increased withdrawals from the Atlantic City 800-foot sand on water levels, by Watt and Voronin (2006) to simulate sources of water to wells in the Wenonah-Mount Laurel aquifer, and by Gordon (2007) to produce water budgets for confined aquifers in the New Jersey Coastal Plain under various withdrawal scenarios.

Eight synoptic water-level studies of the confined aquifers of the New Jersey Coastal Plain document water-level data and potentiometric-surface maps for those aquifers. Potentiometric-surface maps in these reports show water levels in the Coastal Plain at 5-year intervals from 1978 through 2013: 1978, Walker (1983); 1983, Eckel and Walker (1986); 1988, Rosman and others (1995); 1993 and 1998, Lacombe and Rosman (1997 and 2001, respectively); 2003, DePaul and others (2009); 2008, DePaul and Rosman (2015); and 2013, Gordon and others (2021).

Conceptual Hydrogeologic Model

The hydrogeologic framework of the New Jersey Coastal Plain is based on the hydrogeologic framework presented in Zapecza (1989). The framework consists of a southeastward-dipping and -thickening wedge of unconsolidated deposits of sand, silt, and clay of Cretaceous to Neogene age underlain by basement rocks and overlain by locally occurring Quaternary sediments (fig. 1). Coastal Plain sediments were deposited in various shelf, marginal marine, nearshore or coastal beach, and deltaic environments, the extent of which fluctuated in response to changes in sea level. Units composed of substantially less permeable sediments (predominantly clays and fine-grained silts) form the confining units, and coarser, more permeable sand and gravel units, which readily produce water, form the aquifers. These deposits are less than 50 ft thick along the updip limit of the Coastal Plain sediments (Fall Line) in New Jersey and thicken to more than 6,500 ft in southern Cape May County. Coastal Plain sediments generally strike northeast-southwest and dip gently from 10 to 60 feet per mile to the southeast (Zapecza, 1989); overlying Quaternary deposits are flat. The sediments composing the New Jersey Coastal Plain aquifers and confining units generally crop out near the Fall Line parallel to strike; the aquifer units generally transition into confined aquifers, except the Piney Point aquifer, which is confined throughout the study area. The aquifers and confining units in the New Jersey Coastal Plain range in age from Cretaceous to Quaternary (fig. 2).

Chart showing various types of hydrogeologic units.
Figure 2.

Hydrogeologic units for the New Jersey Coastal Plain, and the correlated hydrogeologic units in the New Jersey Coastal Plain groundwater-flow model of the Delaware and Maryland Coastal Plain and the regional hydrogeologic units in the North Atlantic Coastal Plain aquifer system. Unit names for New Jersey are from modified from Zapecza (1989). Unit names for Delaware are modified from Vroblesky and Fleck, (1991) and Masterson and others (2016b). Unit names for Maryland are modified from Andreasen and others (2013) and from Masterson and others (2016b). Units in Maryland are not the primary focus of this report. Unit names in the North Atlantic Coastal Plain groundwater-flow model are from Masterson and others (2016b). Letter and number in parentheses (for example, L1) in the New Jersey column is the model-layer designation in New Jersey, Delaware, and Maryland.

Figure 1 is a generalized hydrogeologic cross section through the New Jersey Coastal Plain; aquifers and confining units are shown in figure 2, which also shows the correlation of the corresponding units in Maryland and Delaware with those in New Jersey. Detailed discussions in Zapecza (1989) and Sugarman and others (2005) describe the hydrogeology of New Jersey; Vroblesky and Fleck (1991) describe the hydrogeology of Delaware; and Andreasen and others (2013) describe the hydrogeology of Maryland. In addition, the hydrogeology of Delaware, Maryland, and New Jersey is described in Trapp and Meisler (1992) and Masterson and others (2015).

Martin (1998) and Voronin (2004) adapted the hydrogeologic framework in Zapecza (1989) for use in the New Jersey RASA model and represented the Coastal Plain sediments as 10 major aquifers and 9 intervening confining units. In the RASA groundwater-flow models of Martin (1998) and Voronin (2004), the Atlantic City 800-foot sand and the Rio Grande water-bearing zone were simulated by using one aquifer layer. The downdip area of this model layer was referred to by Martin (1998) as the confined Kirkwood aquifer; the updip area of the model layer represents the lower part of the unconfined Kirkwood-Cohansey aquifer system. The model layer above the lower part of the unconfined Kirkwood-Cohansey aquifer system and the confined Kirkwood aquifer represents the unconfined upper Kirkwood-Cohansey aquifer system. Martin (1998) states that the Kirkwood-Cohansey aquifer system was subdivided into an upper and lower aquifer in updip areas to more accurately represent the vertical head distribution in the unconfined aquifer system and to provide a lateral connection between the confined Kirkwood aquifer and the lower Kirkwood-Cohansey aquifer. This framework also was used in the rediscretized RASA model of Voronin (2004).

The updated model described in this report represents aquifers and intervening confining units with a hydrogeologic framework of 21 layers—11 aquifers and 10 confining units. This representation includes the 10 aquifers simulated by Martin (1998) and Voronin (2004) as well as the Rio Grande water-bearing zone, which is represented as a separate aquifer layer to more accurately estimate flow in Atlantic and Cape May Counties, where this aquifer is used as a source of water. Each aquifer layer is represented in the fully three-dimensional groundwater-flow model by specifying its hydraulic parameters; confining units are simulated as layers by specifying their hydraulic parameters.

Only the freshwater parts of the aquifer (where the chloride concentration is less than 10,000 milligrams per liter (mg/L) are represented in this study. The locations of the 10,000-mg/L isochlors for selected aquifers in the New Jersey were obtained from simulation results from the New Jersey Coastal Plain SHARP (Essaid, 1990) model (Pope and Gordon, 1999). SHARP, a quasi-three-dimensional finite-difference model that can simulate freshwater and saltwater flow, was used by Pope and Gordon (1999) to simulate the groundwater-flow system in the New Jersey Coastal Plain, including the location and movement of the freshwater-saltwater interface. The simulated freshwater-saltwater interface was used to define the downdip limit of those aquifers and adjacent confining units where a downdip limit is not defined and is discussed further in the section on boundary conditions. The locations of the 10,000-mg/L isochlors for aquifers in Maryland and Delaware were determined from data published in Meisler (1989), Vroblesky and Fleck (1991), and Trapp and Meisler (1992).

Aquifer and Confining-Unit Characteristics

Below is a brief description of the hydrogeologic properties of the 21 model layers, including the corresponding units in Delaware and Maryland. Additional data on the hydraulic properties of the aquifers and confining units in the New Jersey Coastal Plain are summarized in Martin (1998). The hydrologic properties of the units in Maryland and Delaware were obtained primarily from the Maryland Geologic Survey (Andreasen and others, 2013) and also Masterson and others (2016a).

Estimates of hydraulic conductivity were also made for 3,621 wells in the New Jersey Coastal Plain from results of specific-capacity tests conducted with pumping rates of 250 gallons per minute (gal/min) for a minimum of 8 hours. These data are stored in the USGS National Water Information System (NWIS) database (U.S. Geological Survey, 2017). The transmissivity of the aquifer in the vicinity of these wells was estimated from specific-capacity data by using the Theis equation as presented in Heath (1983). The hydraulic conductivity for an aquifer at a well location was estimated from the calculated transmissivity value determined by using the length of the screened interval. The range of horizontal hydraulic conductivity estimated for the aquifers is summarized in table 2 and the well locations are shown in figure 3.

Table 2.    

Horizontal hydraulic conductivities estimated from specific-capacity data for wells in selected aquifers in the New Jersey Coastal Plain.

[ft/d, feet per day; PRM, Potomac-Raritan-Magothy]

Specific-capacity tests1
Aquifer Model layer Number of wells Range of estimated horizontal hydraulic conductivity2 (ft/d) Average estimated hydraulic conductivity (ft/d)
Holly Beach water-bearing zone 1 14 21–104 57
Kirkwood-Cohansey aquifer system 3 1,791 1–1,817 110
Rio Grande water-bearing zone 5 3 17–47 32
Atlantic City 800-foot sand 7 94 5–178 50
Piney Point 9 52 1–1,135 43
Vincentown 11 42 3–220 33
Wenonah-Mount Laurel 13 211 1–733 36
Englishtown aquifer system 15 181 2–249 26
upper PRM 17 517 1–776 102
middle PRM 19 514 4–1,926 119
lower PRM 21 202 6–1,114 137
Table 2.    Horizontal hydraulic conductivities estimated from specific-capacity data for wells in selected aquifers in the New Jersey Coastal Plain.
2

Hydraulic conductivity for a well was estimated from specific-capacity data by using the equation in Heath (1983) which was used to calculate the transmissivity of the aquifer in the vicinity of the well. The transmissivity was then divided by the length of the screened interval of the well to yield the hydraulic conductivity. These data are stored in the U.S. Geological Survey National Water Information System (NWIS) database (U.S. Geological Survey, 2017).

Map with various colored dots showing locations of specific-capacity tests.
Figure 3.

Map showing location of sites of specific-capacity tests, New Jersey Coastal Plain.

Horizontal hydraulic conductivity can vary spatially for each aquifer (fig. 3). Table 2 indicates the order-of-magnitude variability observed in the range of horizontal hydraulic conductivity for some aquifers. Ranges of horizontal hydraulic conductivity values for 11 aquifers in the New Jersey Coastal Plain (table 2) indicate that hydraulic conductivities are greatest at wells screened in the Kirkwood-Cohansey aquifer system and the upper, middle and lower Potomac-Raritan-Magothy aquifers. In addition to the specific-capacity data, the horizontal hydraulic conductivities of the aquifers, vertical hydraulic conductivities of the confining units, and aquifer storage coefficients from selected previously published groundwater-flow models of the New Jersey Coastal Plain are summarized in table 3.

Table 3.    

Hydraulic properties of aquifers and confining units from selected groundwater-flow models in the New Jersey Coastal Plain.

[PRM, Potomac-Raritan-Magothy; ft/d, feet per day; ft, feet; >, greater than; <, less than; —, no data]

Aquifer Model layer Horizontal hydraulic conductivity (ft/d) Thickness1 (ft) Storage coefficient Reference
Holly Beach water-bearing zone 1 126–219 15–225 6.0×10−2 Spitz, 1998
Unconfined Kirkwood-Cohansey aquifer system 3 50–60 50–>400 Pope and others, 2012
Rio Grande water-bearing zone 5 50 >0–140 Pope and others, 2012
Atlantic City 800-foot sand 7 50 40–200 1.1×10−5–1.0×10−4 Pope and others, 2012; McAuley and others, 2001
Piney Point 9 23 >0–>200 3.0×10−4 Rush, 1968
Vincentown 11 20 20–>140 Nicholson and Watt, 1997
Wenonah-Mount Laurel 13 13–19 <25–>120 1.5×10−5–3.5×10−4 Nemickas, 1976
Englishtown aquifer system 15 6–70 40–140 Charles and others, 2011
upper PRM 17 18–90 50–>200 1.0×10−4 Navoy and Carleton, 1995; Martin, 1998
middle PRM 19 7–90 <50–>150 1.0×10−4 Navoy and Carleton, 1995; Martin, 1998
lower PRM 21 15–150 >0–250 1.0×10−4 Charles and others, 2011; Martin, 1998
Confining unit Model layer Vertical hydraulic conductivity (ft/d) Thickness (ft) Storage coefficient Reference
Confining unit overlying the estuarine sand and confining unit overlying the confined Cohansey aquifer 2 4.0×10−3 20–150 Spitz, 1998
Confining unit overlying the Rio Grande water-bearing zone 4 1.00×10−5 Pope and others, 2012
Confining unit overlying the Atlantic City 800-foot sand 6 1.0×10−5 <100–>300 Pope and others, 2012
Composite confining unit2 8 1.0×10−5–1.4×10−5 3100 McAuley and others, 2001
Composite confining unit4 10 5.7×10−3 <50–1,190 Rosenau and others, 1969
Composite confining unit5 12 5.6×10−2 60–90 Nichols, 1977
Marshalltown-Wenonah 14 1.5×10−5 20–>180 Nichols, 1977
Merchantville-Woodbury 16 4.3×10−6–3.5×10−3 100–350 Nichols, 1977; Pucci and others, 1994
Confining unit between upper and middle PRM 18 1.8×10−5–8.0×10−1 50–>200 Pucci and others, 1994; Charles and others, 2011
Confining unit between middle and lower PRM 20 2.1×10−5–1.7×10−1 >0–<50 Charles and others, 2011.
Table 3.    Hydraulic properties of aquifers and confining units from selected groundwater-flow models in the New Jersey Coastal Plain.
1

Thickness from Zapecza (1989), except for model layers 1 and 2 (Spitz, 1998), model layers 3 and 5 (Pope and others, 2012), and model layer 6 (McAuley and others, 2001).

2

Corresponds to the basal Kirkwood confining unit.

3

Thickness in vicinity of Atlantic City, New Jersey, from McAuley and others (2001).

4

Corresponds to the Manasquan-Shark River confining unit.

5

Martin (1998) referred to this part of the composite confining unit as the Navesink-Hornerstown confining unit.

Data on the hydraulic properties of confining units within the modeled area in New Jersey are limited, and most estimates of these hydraulic properties are primarily from the groundwater-flow-model studies listed in table 3. Data on the hydraulic properties of hydrogeologic units within the model area in Delaware and counties in eastern Maryland also are limited and are summarized in table 4. Data presented in table 4 are compiled from the literature of previous estimates of the hydraulic properties of aquifers and confining units within the Delaware and Maryland Coastal Plain and from Andreasen and others, 2013.

Table 4.    

Summary of hydraulic properties of aquifers and confining units in Delaware and counties in eastern Maryland.

[ft, feet; ft2/day, feet squared per day; ft/d, feet per day; —, no data]

Aquifer or confining unit Model layer Average or maximum thickness1 (ft) Delaware or county in Maryland2 Mean or range of transmissivity3 (ft2/d) Horizontal hydraulic conductivity3 (ft/d) Vertical hydraulic conductivity3 (ft/d) Storage coefficient3
Upper Chesapeake confining unit 2 2 50 Worcester
50 Wicomico
Manokin aquifer 3 10–20 Wicomico 480–7,440 54
Up to 195 Worcester 4,820–14,800
430–150 Delaware
St. Mary's confining unit 4 Up to 230 Worcester
Delaware 41.0×10-4
Choptank aquifer 5 Wicomico 440–510
Up to 170 Worcester
Lower Chesapeake confining unit 6 Up to 155 Worcester
Calvert aquifer system 7 90 Queen Anne's 30-50
Calvert confining unit 8 Queen Anne's 0.015 0.02
Piney Point aquifer 9 50 Caroline 100–4,670
Queen Anne's 1,203 1.6×10−4
Up to 240 Delaware 5800–5,350 52.8×10−4– 3.0×10−4
Nanjemoy confining unit 10 Queen Anne's 0.011–10.3 2.3×10-4 – 8.33
Aquia (Rancocas) aquifer6 11 Kent 350–8,090 2.0×10−4–4.0×10−4
Up to 289 Queen Anne's 180–5,600 3.0×10−4
Delaware7 530–2,767 0.15–110 1.9×10−4–4.4×10−4
Severn confining unit 12 Up to 78 Kent
Monmouth (Mount Laurel) aquifer6 13 Kent 220–340 1.2×10−3
Delaware7 330–5,600 4.2 1.0×10−4–7.7×10−3
Matawan confining unit 14 54 Delaware8
Matawan aquifer 15 Queen Anne's 931
Magothy confining unit 16 152 Delaware8
Magothy aquifer 17 45 Cecil 3,100–3,440 1.0×10−4
Kent 500 3.0×10−4
Queen Anne's 5,800–10,000 2.2×10−4
Delaware7 410–1,760 3.7×10−5–3.0×10−3
Potomac confining unit 18 119 Delaware8
Potomac-Patapsco aquifer 19 134 Cecil 100–4,150 150 9.6×10−3
318 Kent 50–2,440 40 7.3×10−3
Queen Anne's 1,000–8,800 50 1.0×10−4−2.2×10−4
Delaware7 455–5,350 1.0×10−4–2.2×10−4
Potomac confining unit 20 199 Delaware8
Patuxent aquifer system
 
21
 
Cecil 380–11,250 1.0×10−4–1.1×10−4
Up to 500 Queen Anne's 800
Table 4.    Summary of hydraulic properties of aquifers and confining units in Delaware and counties in eastern Maryland.
2

County in Maryland within modeled area.

3

Data from aquifer tests summarized in Andreasen and others (2013).

6

Aquifers in parentheses are the correlative aquifer in Delaware.

7

Values of transmissivity, horizontal hydraulic conductivity, and storage coefficient in Dugan and others (2008).

8

Average thickness of regional hydrogeologic units in the Northern Atlantic Coastal Plain aquifer system in Masterson and others (2016a).

Simulation of Groundwater Flow

Groundwater flow in the New Jersey Coastal Plain was simulated under transient conditions for the period 1980–2013. In this section, the previous versions of the New Jersey Coastal Plain models are described briefly, and the development of the three-dimensional groundwater-flow model, model calibration, simulation of groundwater flow, the sensitivity of model parameters to water levels and base flows, and model limitations are discussed.

Groundwater-Flow System

Groundwater flow in the New Jersey Coastal Plain is controlled by topography, the hydraulic properties of the sediments, and hydrologic stresses. Flow within the aquifers is predominantly horizontal, although some vertical flow exists. Hydraulic gradients in the confining units generally are vertical because the horizontal hydraulic conductivity is orders of magnitude less than that in the aquifers. The major source of recharge to the aquifers of the New Jersey Coastal Plain is infiltration of precipitation, although leakage from surface-water bodies, lateral flow from adjacent areas, and in confined aquifers, flow from downdip areas also occurs. Most of the groundwater recharge is discharged to nearby surface-water bodies; a smaller amount infiltrates through confining units to recharge the underlying confined aquifers. In addition to the water that discharges from the aquifers to surface-water bodies (streams that flow to the ocean and bays), some groundwater discharges directly to saltwater bodies as submarine groundwater discharge, or is removed by groundwater evapotranspiration, withdrawals from wells, and lateral flow to adjacent areas. Discharge areas in the New Jersey Coastal Plain include the Atlantic Ocean, Raritan Bay, Delaware River, Delaware Bay, and streams. Because groundwater withdrawals in some areas of the New Jersey Coastal Plain are substantial, large cones of depression developed in several aquifers, and recharge and discharge locations shifted. These changes in regional groundwater-flow patterns were investigated in the eight water-level synoptic studies completed for the New Jersey Coastal Plain in cooperation with the NJDEP and are listed in the “Previous Investigations” section of this report.

Design of Previous Models of the New Jersey Coastal Plain

The New Jersey RASA model (Martin, 1998) was developed in the 1980s from a modified version of the Trescott model code (Trescott, 1975). The grid spacing ranges from 6.25 to 9.375 mi2 in onshore areas to 47.5 mi2 in offshore areas. The model used a constant recharge rate of 20 inches per year (in/yr) and withdrawals averaged over time at pumping intervals from 1896–1980. The model was rediscretized in the 1990s (Voronin, 2004) with a smaller grid spacing, recharge based on water-budget equations that used base flow estimated from hydrograph separation of streamflow data for several surface-water basins in of the New Jersey Coastal Plain, and updated withdrawals to 1998, and was run in MODFLOW–96 (Harbaugh and McDonald, 1996). The grid spacing ranged from 0.25 to 0.31 mi2 in onshore areas to 3.16 mi2 in offshore areas. Martin (1998) and Voronin (2004) modeled the aquifers of the New Jersey Coastal Plain using a quasi-three-dimensional representation of the aquifers that represented the confining units by using a vertical leakance parameter (vertical hydraulic conductivity divided by thickness). Voronin (2004) simulated the aquifer system using a 10-layer model and simulated streams using the drain (DRN) and river (RIV) packages of MODFLOW–96 (Harbaugh and McDonald, 1996). Martin (1998) simulated the aquifer system and streams by using an 11-layer model in which the top layer was a constant-head layer that represented the streams. The major differences between the two previous versions of the New Jersey RASA model and the version described in this report are summarized in table 5. The development, design, and calibration of the model described in this report are discussed in the following sections.

Table 5.    

Characteristics of three groundwater-flow models developed for the New Jersey Coastal Plain.

[RASA, Regional Aquifer System Analysis; in/yr, inches per year; SWB, soil-water balance]

Reference Model computer code1 Total model layers2 Withdrawal conditions Most recent water-level data used Minimum cell size (square miles) Boundary flows derived from: Recharge to cells in New Jersey
Martin (1998) Trescott (1975) 11 1896–1980 1978 6.25 Leahy and Martin (1993) Constant of 20 in/yr
Voronin (2004) MODFLOW–96 10 1968–1998 1998 0.25 Leahy and Martin (1993)3 0.01–20 in/yr; varies by cell but not by year
Gordon and Carleton (2023) MODFLOW–2005 21 41978–2013 2013 0.25 Masterson and others (2016a) Used SWB5 model to vary by cell and by year.
Table 5.    Characteristics of three groundwater-flow models developed for the New Jersey Coastal Plain.
1

Reference for MODFLOW–96 is Harbaugh and McDonald (1996); for MODFLOW–2005 is Harbaugh (2005).

2

The groundwater-flow models of Martin (1998) and Voronin (2004) are quasi-three-dimensional models. Confining units are not explicitly modeled but are simulated by a vertical leakance. The vertical leakance is calculated as the hydraulic conductivity divided by thickness.

3

Boundary flows in Delaware were modified to account for the increase in withdrawals.

4

Stress periods 1 through 3 are transitional stress periods in this model that generally correspond to the time periods 1970–76 and 1977–79, and 1980, respectively, to set antecedent conditions in the model. The withdrawal rates for these three time periods are based on a percentage of the average annual 1978–81 withdrawals in Voronin (2004).

5

The SWB model uses a soil-water-balance code (Westenbroek and others, 2010) for estimating an annual recharge for each year simulated in inches per year. The recharge rates for the time periods 1970–76 and 1977–79 (model stress periods 1 and 2), used to set antecedent conditions in the model, are based on a percentage of the 1980 average recharge.

Model Development

The three-dimensional groundwater-flow model developed to simulate the groundwater-flow system in the New Jersey Coastal Plain is based on the USGS computer program MODFLOW–2005 (Harbaugh and others, 2005). Updates to the two previous versions of the RASA model include (1) a 21-layer hydrogeologic framework that extends into Delaware and part of Maryland; (2) spatially variable recharge estimated by using the SWB model (Westenbroek and others, 2010); (3) updated hydraulic parameters to include the confining units in New Jersey and the hydrogeologic units in Delaware and eastern Maryland; (4) updated boundary flows from the groundwater-flow model of the North Atlantic Coastal Plain (NACP) (Masterson and others, 2016a); (5) groundwater-withdrawal data for 1980–2013 for the New Jersey Coastal Plain and for 1980–2010 for the modeled areas in Delaware and eastern Maryland; and (6) the Rio Grande water-bearing zone and the confining units above and below it, added as a model layer to refine the estimate of simulated flow in Atlantic and Cape May Counties, New Jersey. The updip part of model layer 5, where the Rio Grande water-bearing zone is absent (fig. 1), is represented in this model by a thin layer that allows for a hydraulic connection between the overlying and underlying model layers. The groundwater-flow model is used to (1) simulate current (2013) conditions of groundwater flow in the New Jersey Coastal Plain and the confined aquifers of Delaware, (2) improve the representation of the groundwater system of the New Jersey Coastal Plain, and (3) provide the NJDEP and other water managers with the information (Carleton and others, 2023) that can be used to calculate groundwater-budget components and evaluate water-resource management decisions.

Spatial and Temporal Discretization

The finite-difference grid for the numerical model consists of 145 rows, 245 columns, and 21 layers. The model grid size is variable, with a grid spacing of 0.25 mi2 (the smallest) in the northern and southwestern parts of the New Jersey Coastal Plain, 0.31 mi2 in the southeastern part of the Coastal Plain, and as much as 3.16 mi2 in the offshore areas. The model area is approximately 16,597 mi2 and extends from the Fall Line (fig. 1), which is the transition from the consolidated rocks of the upland Piedmont Physiographic Province to the unconsolidated sediments of the Coastal Plain Physiographic Province, to the approximate downdip and (or) freshwater-saltwater interface in each model layer. (The location and extent of the Fall Line in southeastern New York is not shown in figures.) The grid is oriented approximately parallel to the Fall Line and to the strike of the New Jersey Coastal Plain hydrogeologic units (Zapecza, 1989).

A schematic diagram of the model layers used to represent the hydrogeologic units in the New Jersey Coastal Plain and the corresponding hydrogeologic units in Delaware and eastern Maryland is shown in figure 4. Confining units are modeled as individual layers by assigning them appropriate hydrologic properties.

Various colors used in a graph to show confining units and aquifers, chart below showing
                           units simulated.
Figure 4.

Schematic representation of aquifers, confining units and boundary conditions in the New Jersey Coastal Plain groundwater-flow model.

The model area extends vertically from the water table to consolidated bedrock and includes 11 distinct major aquifers and 10 confining units of the New Jersey Coastal Plain. The hydrogeologic framework described by Zapecza (1989) was revised for this study to improve the simulation of the groundwater-flow system. Extents of each model layer in Delaware and eastern Maryland are based on hydrogeologic data from Meisler (1989), Vroblesky and Fleck (1991) and Andreasen and others (2013). The hydrogeologic framework is described in more detail in the data release for this report (Carleton and others, 2023). This data release contains an ArcGIS shapefile containing the top, bottom, thickness, and extent of each layer in the model used to define the aquifers and confining units incorporated into the New Jersey Coastal Plain model (Carleton and others, 2023).

If a hydrogeologic layer does not extend across a given model layer in the MODFLOW-2005 finite-difference solution method but instead the hydrogeologic layer pinches out, but the model layer above and (or) below is an active area in the model, an area was created in the model layer with the limited extent to make that layer continuous with the above and (or) below layer, if needed, and was given a small thickness (1 foot or more). The hydraulic properties of the overlying model layer were assigned so that the layer would be continuous, thereby achieving a hydraulic connection between the layers above and below. For example, model layer 11 which represents the Vincentown aquifer, an aquifer with a small area, but the areas of the overlying and underlying layers extend beyond the area of the Vincentown aquifer, therefore an area is needed to simulate the flow from overlying aquifers to underlying aquifers to more accurately represent the vertical head distribution in the model layers (fig. 1). Additionally, in the updip area of the Piney Point aquifer (model layer 9) where the aquifer pinches out and does not outcrop, an area is needed to simulate the flow from overlying units to underlying units to more accurately represent the vertical head distribution in the model layers (fig. 1). The same is true for the updip area in model layer 5, which in the downdip area (fig. 1) represents the confined Rio Grande water-bearing zone, and the overlying and underlying confining units (model layers 4 and 6, respectively). In the updip parts of those three model layers, an area is needed to simulate the flow from the overlying model layer 3 to the underlying model layer 7 to more accurately represent the vertical connection in the model layers.

A total of 36 stress periods were simulated—1 steady-state period and 35 transient periods. Stress period 1 is steady-state and stress period 2 is transient and 3 years long. Stress periods 1 and 2 are transitional stress periods that generally correspond to 1970–76 and 1977–79, respectively. Stress periods 3 to 36 are 1 year long and generally represent the period 1980–2013.

Lateral and Lower Boundary Conditions

The model is bounded to the north and northwest and below by the contact between the unconsolidated Coastal Plain deposits and consolidated bedrock that underlies the Coastal Plain sediments. This boundary was represented numerically as a no-flow boundary condition. Lateral boundaries for each layer in the groundwater-flow model are shown in a generalized schematic diagram in figure 4. Although water may flow into the model area along the Fall Line between the subsurface fractured rock and the Coastal Plain sediments, this amount is assumed to be negligible and, therefore, is not considered in this report.

In the two previous New Jersey RASA groundwater-flow models (Martin, 1998; Voronin, 2004), the lateral boundaries in the northeast and southwest are specified-flux boundaries that originally were derived from a model of the NACP (Leahy and Martin, 1993). Some modifications were made to the boundary fluxes near Delaware Bay in the rediscretized RASA model (Voronin, 2004). The Leahy and Martin (1993) model area extended from Long Island, New York, to North Carolina and includes the New Jersey Coastal Plain. In the model described in this report, the lateral model boundaries in the northeast and southwest are specified-flux boundaries derived from the NACP aquifer system model (Masterson and others, 2016a), which also extends from Long Island, New York, to North Carolina, by using the MODFLOW Flow and Head Boundary (FHB) package (Leake and Lilly, 1997). Specified fluxes are applied in areas where the aquifers and confining units intersect the flow boundary. Flows from the NACP model (Masterson and others, 2016a) are incorporated into the FHB package for stress periods 1 through 32. Because the end of the NACP model simulation coincides with stress period 32 of the current New Jersey Coastal Plain model, fluxes for stress period 32 are repeated for stress periods 33 through 36. The FHB package is used to specify lateral boundary fluxes in each active cell along column 1 of the model grid, which traverses eastern Maryland and a small part of southwestern Delaware (fig. 1). Some active cells along column 1 where the flux was zero were not included in the FBH model file. Specified fluxes are also input in New Jersey using the FHB package along the northeastern edge of the active model area for some layers in column 245 for the offshore extent of layers 7 and 8 and in layers 16 to 19 for the offshore part of the upper and middle Potomac-Raritan-Magothy aquifers (layers 17 and 19), the Merchantville-Woodbury confining unit (layer 16), and the confining unit between the upper and middle Potomac-Raritan-Magothy aquifers (layer 18). The remaining layers have a no-flow or a constant-head boundary along northeastern limit of the model area (fig. 4).

The southeastern (downdip) boundary of the aquifers and confining-unit model layers in New Jersey, Delaware, and Maryland are simulated as no-flow boundaries except as described in the “General Head Boundary” section of this report. The downdip extents of the aquifers and confining-unit model layers in Delaware and Maryland are from Meisler (1989), Vroblesky and Fleck (1991), and Andreasen and others (2013). The southeastern downdip boundary of the aquifers and confining-unit layers in New Jersey are generally located at the downdip limit of freshwater in the aquifer (10,000 mg/L isochlor) determined by Pope and Gordon (1999). Pope and Gordon (1999) used the SHARP computer model (Essaid, 1990) to simulate the groundwater-flow system, including the location and movement of the freshwater-saltwater interface in nine aquifers and eight intervening confining units in the Coastal Plain. Although the interface position is not static in the model, the actual amount of movement cannot be quantified because of the model’s coarse cell size. Because the SHARP model is quasi-three-dimensional, the downdip limit of the 10,000-mg/L isochlor in those aquifers was assumed to be the same for the overlying confining unit. However, the Englishtown aquifer system and Vincentown aquifer are not continuous throughout the New Jersey Coastal Plain and contain freshwater throughout their confined extent in New Jersey. However, chloride concentrations exceeded 15,000 mg/L at an observation well in the Englishtown aquifer system in Sandy Hook in northeastern Monmouth County (DePaul and Rosman, 2015).

General Head Boundary

The General Head Boundary (GHB) package (McDonald and Harbaugh, 1988) was used to represent flow in the updip direction from limited parts of the downdip lateral boundaries of layers 7, 9, 11, 13, 15, 17, 19, and 21 (fig. 4). However, in model layer 11, the GHB boundary was used in New Jersey only. The GHB head is constant for all stress periods and is set to 0.0 ft (NAVD 88); no attempt was made to estimate the equivalent freshwater head where groundwater density is greater than that of freshwater. A single value of bed conductance of 1 square foot per second (ft2/s) was applied to all GHB cells.

The GHB boundaries were included to prevent the occurrence of unrealistically low heads downdip from withdrawals; because the boundaries are generally far from the withdrawals, they have only a minor effect on the overall flow system. In 2013, flow into the model from GHB cells is about 0.6 percent of the total flow budget.

Upper Boundary Conditions

The upper boundaries in the model represent the streams and recharge in onshore areas and the Delaware and Raritan Bays and Atlantic Ocean (fig. 1) in offshore areas. In offshore areas, the upper boundary is a constant freshwater equivalent water level specified by using the Time-Variant Specified-Head (CHD) package (Harbaugh, 2005) and is applied across the top active layer of model cells used to represent the Atlantic Ocean and the Delaware and Raritan Bays. CHD cells represent the outcrop of layers 2 through 5 in the Delaware Bay, the outcrop in layers 14 through 19 (excluding layer 18) in Raritan Bay, and the outcrop of layers 2, 3, 7, and 10 through 12 in the Atlantic Ocean (fig. 4).

Recharge

The Recharge (RCH) package of the USGS finite-difference groundwater model MODFLOW–2005 was used to apply recharge to the topmost active layer of each onshore model cell. This upper boundary is spatially variable recharge from precipitation estimated by using the Soil-Water Balance (SWB) model (Westenbroek and others, 2010). The SWB model is based on a modified Thornthwaite-Mather method (Thornthwaite and Mather, 1957) that incorporates spatially distributed land cover, soil properties, and daily climate data to produce a detailed spatial and temporal variability of recharge. SWB model-calculated recharge rates from 1980 to 2013 were used as the initial recharge values for New Jersey in the groundwater-flow model. The SWB model is described in appendix 1. Recharge was adjusted as part of the model parameter-estimation process to allow for more variability than that calculated by the SWB model and adjusted for comparison with base-flow estimates at various surface-water drainage basins in the New Jersey Coastal Plain. Average annual precipitation in southern New Jersey (Office of the New Jersey State Climatologist, 2021) during 1970–76 and 1977–79 were 1.23 and 1.31 times the 1980 annual average, respectively. Therefore, the recharge rates applied in model stress periods 1 (associated with 1970–76) and 2 (1977–79) are 1.23 and 1.31 times those applied in stress period 3 (1980). Stress periods 1 and 2 are used in the model to set antecedent conditions for the other years.

Because the surficial aquifer in Delaware and Maryland is not simulated, cell-by-cell vertical flows into or out of active cells in model layer 1 in Delaware and Maryland were calculated from the NACP model (Masterson and others, 2016a) and applied with the recharge.

Drain Package

The MODFLOW Drain package (DRN) is used to represent the rivers and streams in the model where groundwater flows (drains) out of the model area (fig. 5). The DRN package was used for continually gaining stream reaches in the New Jersey Coastal Plain which are not a source of groundwater recharge except on a localized scale (for example, upstream from a dam). Use of the DRN package allows for the representation of headwater reaches, where the streambed elevation may be above the water table in places, preventing simulated small streams from being a source of potentially unlimited volumes of water, as would be the case with the MODFLOW River package (RIV) (discussed below). The drain elevation is the average elevation of a stream or stream tributary in the topmost active layer in each onshore model cell. Streams in New Jersey are those used in the previous model version from Voronin (2004). Model cells containing streams in New Jersey were identified by intersecting the model-grid polygons with streams digitized from USGS 1:24,000-scale topographic maps. Streams in Maryland and Delaware were identified by intersecting the grid with NHDPlus (U.S. Geological Survey, 2005), a comprehensive geospatial dataset developed by the U.S. Environmental Protection Agency and the USGS.

Map showing various colors representing model layers with river and stream reaches.
Figure 5.

Map showing location of stream cells in the model area, New Jersey Coastal Plain groundwater-flow model.

The elevations of drain boundaries are kept constant for all stress periods in the transient model. The drain elevations are long-term averages and do not include annual or seasonal variations. Seasonal variations cannot be considered when annual stress periods are used; neglecting the variation of average elevation from year to year is minor in the context of surface-water elevations averaged over model cells of 0.25 mi2 or larger.

The DRN package also requires specification of streambed conductance, an aggregate parameter that represents the product of the streambed hydraulic conductivity, the length of the stream reach, and the stream width divided by the streambed thickness. The same streambed conductance is applied to all drains and was adjusted during calibration.

River Package

The MODFLOW River package (RIV) package was used to represent the Delaware River (fig. 5). The RIV package (unlike the DRN package) allows the Delaware River to be a source of recharge to the aquifer near pumping centers. This is appropriate because the Delaware River carries substantial flow, is tidal in the simulated area, and can provide unlimited flow to the aquifer. The river stage is constant for all stress periods and increases from 0.0 ft (NAVD 88) at the southwestern extent of model river cells in Delaware Bay to 1.0 ft near Trenton by 0.1 ft in segments of approximately the same length. A single value of riverbed conductance was applied to all river cells, representing the product of the riverbed hydraulic conductivity, length, and width of the river reach in each model cell divided by riverbed thickness. Riverbed conductance was modified during model calibration.

Hydraulic Properties

The hydraulic properties hydraulic conductivity and storage coefficient are described in the “Hydraulic Properties of Aquifers and Confining Units” section near the beginning of this report. As previously mentioned, data for the hydraulic properties of confining units are limited. The hydraulic properties are specified by parameter zones within each layer that are assumed to have similar properties. The initial values are comparable with those from Martin (1998), Voronin (2004) and Masterson and others (2016a). During model calibration, the range of hydraulic properties shown in tables 2 and 3 were used to provide supplemental information for values in areas of New Jersey. The hydraulic conductivities and storage coefficients in table 4 were used to provide supplemental information for values in Delaware and Maryland. Hydraulic-conductivity and storage-coefficient parameters were adjusted manually and by using the parameter-estimation software UCODE-2014 (Poeter and others, 2014) during model calibration as needed to improve the model fit between observed data and simulated results.

Groundwater Withdrawals

Groundwater-withdrawal data for the New Jersey Coastal Plain were tabulated and mapped to assess the volume of water pumped from wells in each of the aquifers. These wells include those used for public supply, large-scale agricultural (irrigation), and commercial or industrial purposes. Withdrawals from small-capacity wells, such as those used for domestic supply, are not incorporated into the model. This is a limitation of the model in the outcrop areas of confined aquifers or the Kirkwood-Cohansey aquifer system in areas where a substantial number of domestic-supply wells are located, but domestic self-supply withdrawals do not make up a substantial part of the groundwater withdrawals. Statewide production from domestic wells accounts for about 3 percent of total withdrawals in New Jersey (New Jersey Department of Environmental Protection, 2017b).

Withdrawal data for 1980–1998 were obtained from Voronin (2004) and were reviewed for updates that may have been made since 2003, such as an update to an aquifer code or surveyed well location. Withdrawal data for 1999–2003 were obtained from data reported to the NJDEP. These data were quality reviewed and incorporated into a water-use database available in the data release accompanying this report (Carleton and others, 2023). This withdrawal information includes permitted data only—that is, only data from wells in which daily withdrawals meet or exceed 100,000 gallons for a period of more than 30 days in a consecutive 365-day period. Withdrawal data for 2004–13 were obtained from the NJWATr database developed by the USGS and maintained by the NJDEP to track water withdrawals, use, treatment, and discharge in New Jersey (New Jersey Department of Environmental Protection, 2017a). Total withdrawals from New Jersey Coastal Plain aquifers input to the model for 2004 to 2013 are shown in figure 6. The largest withdrawals occurred in the unconfined Kirkwood-Cohansey aquifer system, the Atlantic City 800-foot sand, and the upper, middle, and lower Potomac-Raritan-Magothy aquifers (model layers 3, 7, 17, 19, and 21).

Bar graph showing groundwater withdrawals by years.
Figure 6.

Graph showing estimated groundwater withdrawals from New Jersey Coastal Plain aquifers, 2004–13.

Annual withdrawal data for 1980–2013 for Delaware and Maryland were obtained from a dataset of pumping rates for 1980–2010 from the groundwater-flow model developed for a groundwater-availability study of the NACP by Masterson and others (2016a; Jack Monti, Jr., U.S. Geological Survey, written commun., 2014). Agricultural withdrawal data for Delaware and Maryland may be less complete and less well documented (Masterson and others, 2016a) than New Jersey agricultural withdrawals. Masterson and others (2016a) provide a detailed description of the methods used to tabulate and estimate groundwater withdrawals in their study. Withdrawals from the unconfined surficial aquifer in Delaware and Maryland were not included in the New Jersey Coastal Plain groundwater-flow model. The unconfined surficial aquifer in Delaware and eastern Maryland was not included in the model hydrogeologic framework because model layer 1 represents the Holly Beach water-bearing zone in Cape May County, New Jersey, which is not continuous into Delaware.

Stress periods 1 (steady-state) and 2 (transient) are transitional stress periods that generally correspond with 1970–76 and 1977–79, respectively. Stress periods 3 and 4 generally represent 1980 and 1981, respectively. Withdrawals for stress periods 1, 2, and 3 are 73, 86, and 93 percent, respectively, of those in stress period 4, and account for the estimated increase in withdrawals from 1970 to 1980. Withdrawals for stress period 4 generally are an average of 1978–81 withdrawals in Voronin (2004) and represent 1981 withdrawal rates in the groundwater-flow model. Average annual withdrawals were used in the model for stress periods 5 to 36, which are 1 year long, and the average annual withdrawals are equivalent to total annual withdrawals during each of these stress periods. Stress periods 5 through 36 represent each year from 1982 to 2013, respectively.

Model Sensitivity Analysis and Calibration

Model calibration was conducted using manual adjustments of parameters through trial and error and parameter estimation techniques. The UCODE-2014 (Poeter and others, 2014) was used to evaluate the sensitivity of parameters to improve the fit of simulation results to observation data and identify correlated parameters that may not be estimated separately. Parameter estimation is the process of determining which parameters can be statistically estimated and then varying those parameters to minimize the weighted least squares objective function (Hill and Tiedeman, 2007). The objective function measures this fit by quantitatively comparing simulated and observed values. This section of the report describes the observations, residuals, final parameter values, correlations between parameters, and parameter sensitivity achieved by the model calibration and the sensitivity of water levels and base flow to the model parameters.

Effective use of hydrogeologic data and observations to constrain a model is likely to produce a model that improves the accuracy of the representation of groundwater flow and, consequently, simulation results. The match of observed to simulated values is used to evaluate how well a model represents an actual system. The simulated and observed values compared for the New Jersey Coastal Plain groundwater-flow model are water-level altitudes (heads) and groundwater discharge (base flow) to gaged streams which involved the adjustment of model parameters so that residuals (the differences between simulated and observed water levels or base flows) are minimized. During calibration, initial model parameters, such as hydraulic conductivity and storage coefficient, are changed to minimize residuals. All observed and simulated water-level observations and simulated and estimated base-flow observations are available in the data release accompanying this report (Carleton and others, 2023).

Water-Level Observations

The model was calibrated by using water-level observations that include observations derived from selected water-level measurements described below. The water-level observations in New Jersey are a subset of measurements collected during the seven quinquennial New Jersey Coastal Plain water-level synoptic studies (1983, 1988, 1993, 1998, 2003, 2008, and 2013) and were selected to provide a spatial distribution in 10 of the simulated aquifers. The water-level synoptic studies did not include the unconfined Holly Beach water-bearing zone and the unconfined part of the Kirkwood-Cohansey aquifer system. Water levels in correlated aquifers in Delaware were obtained from the Delaware Geological Survey (Delaware Geological Survey, 2017). Most of the water levels were measured during October–December after large summer withdrawals had decreased but before water levels had recovered to those measured during spring, thereby representing an approximation of the annual average water level. Water levels were measured at least once during 1983–2013 in 392 wells in New Jersey and 48 wells in Delaware (fig. 7), yielding a total of 2,533 water-level observations in 440 wells. Hydraulic water levels at each well were calculated by subtracting the water level, in feet below land surface, from the land-surface altitude, in feet above the North American Vertical Datum of 1988 (NAVD 88); the result of this calculation is shown in this report as the water-level altitude.

Map with various colors showing wells.
Figure 7.

Map showing location of wells with measured water levels, New Jersey Coastal Plain.

Average annual water levels were calculated for 29 wells in New Jersey (fig. 7) where continuous water levels were measured: 19 wells with observations for every year during 1980–2013 and 10 wells with data for most of that period (average of 27 years). A total of 920 observations from the 29 wells were included in the model calibration; these data are shown in time-series hydrographs shown in this report in the “Calibration of Water Levels” section.

In addition to the water-level observations used in the model calibration, derived observations also were included by calculating the vertical gradient at nested observation wells in New Jersey. Vertical gradients were calculated for 33 pairs of nested wells, for a total of 210 observations. (Nested well pairs are wells that are near each other but are open to different aquifers.) For each sequential synoptic water-level measurement period, vertical gradients calculated by subtracting water levels in colocated wells screened in different aquifers provided information for calibrating the hydraulic conductivity of the confining units. By attempting to match the gradient across the confining unit, the parameter-estimation function improved the accuracy of the vertical-conductance value, for example, in a situation where simulated water levels in both aquifers were somewhat high but the gradient was exact. Also, the change in water level over time was calculated by subtracting water levels in selected wells in which measurements had been made in sequential water-level synoptic studies; the difference was used in calibrating storage in aquifers and confining units. Changes in water levels over time were calculated for 138 wells, 134 wells in New Jersey and four wells in Delaware, where water levels had varied substantially (approximately 10 ft) over the 30-year span of synoptic water-level measurements, for a total of 767 observations.

Water-level observations and vertical-gradient and change-over-time observations calculated from water-level measurements were all assigned the same weight, except for annual observations from the 29 wells with a continuous water-level recorder. Most water-level hydrographs for selected wells had 34 observations, whereas many synoptically measured wells had at most 7 observations; therefore, synoptic water-level observations and calculated observations were assigned a weight of 5, and hydrograph observations were assigned a weight of 1.1. Hill and Tiedeman (2007) describe methods for assigning weights on the basis of estimated error ranges so that observations with substantial potential errors exert a smaller effect on the parameter estimation function than the more accurate measurements do. This approach was not practical for the model in this study because water-level observations were collected under a wide range of conditions, often unknown. For example, some water-level measurements were made in wells with varying degrees of recent and nearby pumping, whereas others were made in unused wells known to be thousands of feet from current or recent withdrawals; still other water-level measurements were made at production wells that had been pumped as recently as 1 hour earlier. Attempting to quantify potential errors based on proximity in time or space to known and unknown withdrawals for synoptic water-level measurements conducted over a 30-year period was infeasible. Most land-surface-elevation estimates were determined from light detection and ranging (lidar) data considered to be accurate to within 2 ft, a modest error compared to the uncertainties related to proximal withdrawals at the time of measurement.

Base-Flow Observations

The New Jersey Coastal Plain model was calibrated by using a total of 1,485 estimates of base flow from 47 stream basins in New Jersey. The locations of the 47 basins are shown in figure 8. The surface-water basins were selected because of location and availability of data during the period 1980—2013. Annual base flow was calculated by using the hydrograph-separation technique PART in Barlow and others (2015) for each of the 21 continuous-record gaging stations for each year of record. For 26 partial-record gaging stations in the New Jersey Coastal Plain, an annual base flow was calculated by correlating manual low-flow measurements made at the partial-record stations with same-day flows at index stations by using the MOVE.1 program (Colarullo and others, 2018). Annual base-flow estimates for each of the 47 basins were averaged to obtain a period-of-record average base flow for each basin.

Map with various colored lines showing updip limits, continuous and partial record
                           surface water basins, large purple box showing model grid outline.
Figure 8.

Map showing location of surface-water drainage basins with base-flow observations, New Jersey Coastal Plain.

Base-flow observations determined by using the PART hydrograph-separation technique in Barlow and others (2015) of continuous-record-station data were assigned a weight of 1.0. Observations determined from partial-record-station data correlated with measurements from continuous-record index stations were assigned a weight less than 1.0 based on the percent standard error of estimate (%SEE), which is the percent of unexplained variance relative to the value of the estimate. The %SEE is a better measure of uncertainty than the standard error of estimate (SEE) because it accounts for the proportional effect encountered in natural streamflow processes, where variance increases with flow (Amy McHugh, U.S. Geological Survey, written commun., 2020). Colarullo and others (2018) defined weights as the reciprocal of %SEE to lessen the weight of predictions made using MOVE.1 regressions that show a high degree of unexplained variability in flows and increase the weight of predictions made using regressions with a small amount of unexplained flow variability. The weights of measurements made at partial-record stations were set to a number less than 1.0 by dividing 4.5 by the squared average %SEE for each partial record station; resulting weights ranged from 0.95 to 0.22.

Model Parameters

A parameter is defined in this report as a single value assigned to a variable used in the finite-difference groundwater-flow equation at one or more model cells. The 422 parameters used in the model fall into seven categories: horizontal hydraulic conductivity (184 parameters), vertical hydraulic conductivity (184 parameters), storage coefficient (42 parameters), riverbed conductance (1 parameter), streambed conductance (1 parameter), general-head-boundary conductance (8 parameters), and recharge (2 parameter). The initial parameter values for several aquifers and confining units in New Jersey, Delaware, and eastern Maryland were changed during calibration to improve the match with observations and improve model results.

The final calibrated values for horizontal and vertical hydraulic conductivity and storage coefficients for aquifers and confining units are listed in tables 6 through 9. Table 6 lists the horizontal hydraulic conductivity (Kh) for the aquifer layers, table 7 lists the horizontal hydraulic conductivity for the confining-unit layers, and table 8 lists the vertical hydraulic conductivity (Kv) for the model layers. Table 9 lists the storage coefficients for each model layer. Parameter zones in each model layer were created using these values to allow for spatial variability of Kh and Kv. Several parameter zones were added or modified during model calibration. The conversion to a fully three-dimensional model framework and varied recharge resulted in the need for additional Kh parameter zones in the outcrop areas of all layers with updip outcrops following the strike of the hydrogeologic layers, corresponding to model layers 3 through 8 and 10 through 21. During the model calibration, some zones acquired a parameter value equal to or similar to that of an adjacent zone but remained a separate zone and was not incorporated into the adjacent zone. The ratio of horizontal to vertical hydraulic conductivity (Kh:Kv) was initially set to 10:1 for aquifers and 1:1 for confining units. Kh and Kv were adjusted independently during calibration; the sensitivity of Kv is low in most aquifer parameter zones and the sensitivity of Kh is low in most confining-unit parameter zones. The zones are designated by the parameters listed in tables 6 through 9. The location and extent of the Kh zones for the aquifers and Kv zones for the confining units are shown in figures 9 through 19. The parameter names for horizontal hydraulic conductivity for the confining units are zones with prefix “CU”. The location and extent of the horizontal hydraulic conductivity parameter zones for the confining units (table 7) are not shown in figures. The location and extent of the horizontal hydraulic conductivity parameter zones for the confining units coincide with the location and extent of the vertical hydraulic conductivity parameter zones (zones with prefix “VK”) for the confining units shown in figures 9B through 18B for those model layers.

Table 6.    

Final horizontal hydraulic conductivity parameters for the aquifers in the calibrated New Jersey Coastal Plain groundwater-flow model.

[The number of significant figures used is based on the magnitude of the parameter value. Parameter names and zones for horizontal hydraulic conductivity for aquifers are shown in figures 9A19.]

Parameter name Parameter value (feet per day) Model layer Parameter name Parameter value (feet per day) Model layer
L1HLBC 173 1 L13DWNDIP 4 13
L3DWNDIP 63 3 L13OCMLRW3 100 13
L3OCMANOK 71 3 L13OCMLRW1 500 13
L3CKKD 26 3 L13OCMLRW2 239 13
L3MANOK 500 3 L13OCMONM1 2 13
L5UPDIP 5 5 L13OCMONM2 2 13
L5OCCHOP 1,000 5 L13DWNDPSW 8.64 E−03 13
L5RIO 23 5 L13MON 6 13
L5CHOPUPDP 22 5 L13BURLMID 16 13
L5UPDIPWES 5 5 L13MIDNE 500 13
L5OUTLIER 247 5 L13MIDSW 14 13
L5CHOPTANK 5.18 E−02 5 L15DWNDIP 1.72 E−01 15
L7UPDIP 5.18 E−01 7 L15DEFTHRU 261 15
L7UPNE 47 7 L15DEEGLS 8.04 E−01 15
L7MIDUP 4.32 E−02 7 L15EGLSCA1 6 15
L7MIDDWN 6 7 L15OCEGLS1 1,000 15
L7800FT 78 7 L15OCEGLS3 432 15
L7CALVERT 24 7 L15OCEGLS2 1,000 15
L7OCCALVR 333 7 L15EGLS 86 15
L7OCKRKD1 136 7 L15MIDNE 1 15
L7OCKRKD3 1,002 7 L15MIDSW 3 15
L7OCKRKD2 1,002 7 L17DWNDIP 21 17
L9UPDIP 2 9 L17CA1 66 17
L9ATLCAPE 8 9 L17BURL 69 17
L9DEPNPN 4 9 L17DEMAGOT 2 17
L9PNPN 3 9 L17GLOU 69 17
L9DEUPDIP 3 9 L17OUTCCA1 500 17
L9DWNDIP 9 9 L17SALEM 21 17
L9CUMB 2 9 L17OCMAG1 107 17
L9SALCUMB 4 9 L17MAGMID 6 17
L9OCEAN 2 9 L17OCSAL 1,002 17
L11DWNDIP 1.73 E−02 11 L17OUTCCA2 1,000 17
L11RANCOC 500 11 L19SALEM 4 19
L11VNCNNE 6.72 E−01 11 L19DNDIP 78 19
L11OCVNCN4 579 11 L19OCUPTPS 461 19
L11OCRANCO 500 11 L19UPTPSCO 7.78 E−01 19
L11OCVNCN2 500 11 L19CA2 104 19
L11OCVNCN1 500 11 L19CA1 78 19
L11OCVNCN3 579 11 L19OCCA2 500 19
L11RANDDIP 72 11 L19OUTCCA1 100 19
L11MIDDPSW 2 11 L21DWNDIP 3.06 E−02 21
L11MIDDPNE 2.59 E−01 11 L21CA2 47 21
L11RANCMD 407 11 L21PATUX 3 21
L11VNCNSW 6.72 E−01 11 L21UPDIP 52 21
L13SALGLO 13 13 L21OCPATX2 43 21
L13BURL 4 13 L21SALEM 2.64 E−01 21
L13DEMLRW 8.02 E−01 13
Table 6.    Final horizontal hydraulic conductivity parameters for the aquifers in the calibrated New Jersey Coastal Plain groundwater-flow model.

Table 7.    

Final horizontal hydraulic conductivity parameters for confining units in the calibrated New Jersey Coastal Plain groundwater-flow model.

[The number of significant figures used is based on the magnitude of the parameter value.]

Parameter name Parameter value (feet per day) Model layer Parameter name Parameter value (feet per day) Model layer
CU2 7.43 E−04 2 CU12SEVERN 1.14 E−06 12
CU2UPCHES 1,000 2 CU12OCSEV1 7.10 E−03 12
CU2DNDIP 1.00 E−06 2 CU12OCNVS3 228 12
CU4DWNDIP 1.65 E−04 4 CU12OCNVS2 308 12
CU4OCSTMAR 448 4 CU12SEVMD 1.12 E−06 12
CU4CAPEMAY 3.29 E−02 4 CU12SAGLMD 1.99 E−06 12
CU4CM 8.64 E−07 4 CU12CABUMD 1.00 E−06 12
CU4OUTLIER 247 4 CU12BURLUD 1.00 E−06 12
CU4STMARY 5.87 E−02 4 CU14DWNDIP 8.64 E−05 14
CU4ACMID 3.58 E−05 4 CU14UPDIP 4.75 E−05 14
CU4UPDIP 5 4 CU14OCR3 20 14
CU6DWNDIP 1.00 E−06 6 CU14OCR2 20 14
CU6LCHES 2.25 E−06 6 CU14OCMAT1 2.48 E−03 14
CU6CAPEALT 1.62 E−04 6 CU14MATAWN 1.00 E−02 14
CU6CM 1.12 E−05 6 CU14OCR1 105 14
CU6UPDIP 5 6 CU14UPDPSW 3.46 E−03 14
CU6COAST 1.17 E−04 6 CU14MIDNE 8.64 E−06 14
CU6OUTLIER 247 6 CU14DEUPDP 9 14
CU6OCLCHE3 72 6 CU14OCMAT2 4.26 E−03 14
CU6OCLCHE1 72 6 CU16DWNDP 1.00 E−06 16
CU6ACMID 2.97 E−04 6 CU16MDDPSW 1.00 E−06 16
CU6OCMID 3.79 E−01 6 CU16MDDPNE 1.00 E−06 16
CU6OCLCHE2 72 6 CU16NWCP 1.00 E−06 16
CU8BASEKRK 2.71 E−04 8 CU16DE 1.00 E−06 16
CU8OCCALV2 106 8 CU16SECP 1.00 E−06 16
CU8UPDIP 8.32 E−04 8 CU16OUTC2 8 16
CU8DWNDIP 2.82 E−02 8 CU16SOC 1.00 E−06 16
CU8CALVERT 1.11 E−05 8 CU16OUTC3 20 16
CU8ALTCUMB 1.00 E−06 8 CU16MATAWN 1.77 E−05 16
CU8OCCALV3 106 8 CU16OUTC1 19 16
CU8OCCALV1 106 8 CU18CP 3.67 E−04 18
CU10ATLCMB 1.05 E−06 10 CU18DNDIP 1.03 E−06 18
CU10DWNDIP 2.09 E−06 10 CU18DEOC1 2 18
CU10DEUPD1 1 10 CU18UPDIP 1.83 E−03 18
CU10OCMAN2 1,000 10 CU18UPDPSW 1.83 E−04 18
CU10NAEMOY 2.30 E−06 10 CU18DEOC2 1,000 18
CU10UPDIP 3.65 E−03 10 CU18DE 1.00 E−06 18
CU10OCMAN1 64 10 CU18DEOC3 207 18
CU10DEUPD2 1 10 CU20DWNDIP 4.84 E−02 20
CU10UPDPSW 3.65 E−03 10 CU20DWNSW 4.84 E−02 20
CU10NOCBRL 6.15 E−05 10 CU20DEOC1 8 20
CU12DWNDIP 1.76 E−03 10 CU20SAL 1.00E−06 20
CU12MON 4.80 E−02 12 CU20DEOC2 8 20
CU12OCSEV2 7.10 E−03 12 CU20DE 1.95 E−04 20
CU12SALGLO 2.25 E−03 12 CU20GLOCAM 2.76 E−01 20
CU12OCNVS1 20 12
Table 7.    Final horizontal hydraulic conductivity parameters for confining units in the calibrated New Jersey Coastal Plain groundwater-flow model.

Table 8.    

Final vertical hydraulic conductivity parameters in the calibrated New Jersey Coastal Plain groundwater-flow model.

[The number of significant figures used is based on the magnitude of the parameter value. Parameter names and zones for vertical hydraulic conductivity for confining units are shown in figures 9B18B.]

Parameter name Parameter value (feet per day) Model layer Parameter name Parameter value (feet per day) Model layer Parameter name Parameter value (feet per day) Model layer
VK1HLBC 9 1 VK10ATLCMB 1.05 E−07 10 VK15DWNDIP 1.72 E−02 15
VK2 3.54 E−04 2 VK10DWNDIP 2.09 E−07 10 VK15DEFTHR 15 15
VK2UPCHES 100 2 VK10DEUPD1 1.19 E−01 10 VK15DEEGLS 1.67 E−03 15
VK2DNDIP 1.00 E−06 2 VK10OCMAN2 3 10 VK15EGLCA1 6.48 E−01 15
VK3DWNDIP 6 3 VK10NAEMOY 8.36 E−04 10 VK15OCEGL1 500 15
VK3OCMANOK 4 3 VK10UPDIP 1.56 E−02 10 VK15OCEGL3 43 15
VK3CKKD 4 3 VK10OCMAN1 2.90 E−01 10 VK15OCEGL2 6 15
VK3MANOK 50 3 VK10DEUPD2 1.19 E−01 10 VK15EGLS 9 15
VK4DWNDIP 4.34 E−06 4 VK10NOCBRL 4.73 E−07 10 VK15MIDNE 4 15
VK4OCSTMAR 4.69 E−01 4 VK10UPDPSW 2.07 E−05 10 VK15MIDSW 5 15
VK4CAPEMAY 7.39 E−07 4 VK11DWNDIP 5.18 E−02 11 VK16NWCP 1.13 E−07 16
VK4CM 1.39 E−06 4 VK11RANCOC 2.59 E−02 11 VK16DE 1.78 E−05 16
VK4OUTLIER 100 4 VK11VNCNNE 5.50 E−01 11 VK16SECP 1.00 E−06 16
VK4STMARY 5 4 VK11OCVNC4 58 11 VK16OUTC2 1 16
VK4UPDIP 5 4 VK11OCRANC 100 11 VK16OC 1.00 E−06 16
VK4ACMID 1.27 E−06 4 VK11OCVNC2 58 11 VK16OUTC3 2 16
VK5UPDIP 5 5 VK11OCVNC1 58 11 VK16MATAWN 3.37 E−06 16
VK5OCCHOP 100 5 VK11OCVNC3 58 11 VK16OUTC1 2 16
VK5RIO 8.64 E−01 5 VK11RANDDP 4 11 VK16DWNDP 1.12 E−06 16
VK5CHOPTNK 5.18 E−03 5 VK11VNCNSW 5.95 E−02 11 VK16MDDPSW 1.00 E−06 16
VK5CHOPUPD 2 5 VK11MDDPSW 6.05 E−02 11 VK16MDDPNE 1.00 E−06 16
VK5UPDIPWE 5.37 E−01 5 VK11MDDPNE 2.59 E−02 11 VK17DWNDIP 2 17
VK5OUTLIER 34 5 VK11RANCMD 41 11 VK17CA1 8 17
VK6DWNDIP 1.00 E−06 6 VK12DWNDIP 1.73 E−03 12 VK17BURL 7 17
VK6LCHES 7.12 E−07 6 VK12MON 8.73 E−05 12 VK17DEMAGO 6.05 E−01 17
VK6CAPEALT 1.62 E−05 6 VK12OCSEV2 7.10 E−03 12 VK17GLOU 7 17
VK6CM 3.99 E−05 6 VK12SALGLO 2.25 E−03 12 VK17OCCA1 17 17
VK6UPDIP 5 6 VK12OCNVS1 9.42 E−05 12 VK17SALEM 2 17
VK6COAST 1.29 E−05 6 VK12SEVERN 1.14 E−06 12 VK17OCMAG1 26 17
VK6OUTLIER 1,000 6 VK12OCSEV1 7.10 E−03 12 VK17MAGMID 5 17
VK6OCLCHE3 7 6 VK12OCNVS3 4 12 VK17OCSAL 100 17
VK6OCLCHE1 7 6 VK12OCNVS2 4 12 VK17OCCA2 100 17
VK6OCLCHE2 7 6 VK12SEVMD 1.12 E−06 12 VK18CP 1.57 E−04 18
VK6ACMID 1.81 E−04 6 VK12SAGLMD 7.25 E−06 12 VK18DNDIP 1.03 E−06 18
VK6OCMID 6.32 E−01 6 VK12CABUMD 1.00 E−06 12 VK18DEOC1 1.87 E−01 18
VK7UPDIP 5.37 E−01 7 VK12BURLUD 1.00 E−03 12 VK18UPDIP 3.04 E−03 18
VK7UPNE 24 7 VK13SALGLO 1 13 VK18UPDPSW 1.83 E−04 18
VK7MIDUP 4.32 E−01 7 VK13MON 6.05 E−01 13 VK18DEOC2 100 18
VK7MIDDWN 8 7 VK13DEMLRW 8.04 E−02 13 VK18DE 1.00 E−06 18
VK7800FT 2 7 VK13DWNDIP 3.64 E−03 13 VK18DEOC3 21 18
VK7CALVERT 9 7 VK13OCMLR3 3 13 VK19SALEM 4.32 E−01 19
VK7OCCALVR 100 7 VK13OCMLR1 4 13 VK19DNDIP 2 19
VK7OCKRKD1 4.83 E−01 7 VK13OCMLR2 14 13 VK19OCUPTP 23 19
VK7OCKRKD3 100 7 VK13OCMON1 3.88 E−01 13 VK19UPTPSC 7.78 E−02 19
VK7OCKRKD2 100 7 VK13OCMON2 1.61 E−01 13 VK19CA2 10 19
VK8BASEKRK 2.59 E−04 8 VK13MIDNE 38 13 VK19CA1 8 19
VK8OCCALV2 11 8 VK13MIDSW 1 13 VK19OCCA1 11 19
VK8UPDIP 4.50 E−03 8 VK13DWNDSW 8.64 E−04 13 VK19OCCA2 49 19
VK8DWNDIP 1.09 E−04 8 VK13BURL 4.32 E−01 13 VK20DWNDIP 4.84 E−02 20
VK8CALVERT 1.00 E−07 8 VK13BURLMD 2 13 VK20DWNSW 4.84 E−02 20
VK8ALTCUMB 1.00 E−07 8 VK14DWNDIP 8.64 E−05 14 VK20DEOC1 8.23 E−01 20
VK8OCCALV3 11 8 VK14UPDIP 4.75 E−05 14 VK20SAL 1.00 E−06 20
VK8OCCALV1 11 8 VK14OCR3 2 14 VK20DEOC2 8.23 E−01 20
VK9UPDIP 772 9 VK14OCR2 2 14 VK20DE 1.49 E−04 20
VK9ATLCAPE 4.82 E−01 9 VK14OCMAT1 2.48 E−03 14 VK20GLOCAM 2.76 E−01 20
VK9DEPNPN 4.22 E−01 9 VK14MATAWN 1.00 E−02 14 VK21DWNDIP 3.86 E−05 21
VK9PNPN 3.46 E−01 9 VK14OCR1 20 14 VK21CA2 5 21
VK9DEUPDIP 100 9 VK14OCMAT2 4.26 E−03 14 VK21PATUX 2.59 E−01 21
VK9DWNDIP 8.64 E−01 9 VK14UPDPSW 1.90 E−06 14 VK21UPDIP 5 21
VK9CUMB 1.73 E−01 9 VK14MIDNE 1.02 E−05 14 VK21OCPAT2 6 21
VK9SALCUMB 4.32 E−01 9 VK14DEUPDP 9 14 VK21SALEM 3.67 E−01 21
VK9OCEAN 86 9
Table 8.    Final vertical hydraulic conductivity parameters in the calibrated New Jersey Coastal Plain groundwater-flow model.

Table 9.    

Final storage-coefficient parameters in the calibrated New Jersey Coastal Plain groundwater-flow model.

Confined part of model layer Unconfined part of model layer
Parameter name Parameter value (dimensionless) Model layer Parameter name Parameter value (dimensionless) Model layer
SS_Par1 4.53 E−05 1 SS_Par2OC 1.50 E−01 2
SS_Par2 1.22 E−06 2 SS_Par3OC 4.00 E−01 3
SS_Par3 1.78 E−06 3 SS_Par4OC 8.86 E−04 4
SS_Par4 2.47 E−04 4 SS_Par5OC 1.50 E−01 5
SS_Par5 1.82 E−05 5 SS_Par6OC 1.96 E−01 6
SS_Par6 3.04 E−06 6 SS_Par7OC 4.00 E−01 7
SS_Par7 1.27 E−05 7 SS_Par8OC 4.00 E−01 8
SS_Par8 5.46 E−07 8 SS_Par10OC 2.23 E−04 10
SS_Par9 2.75 E−06 9 SS_Par11OC 5.80 E−06 11
SS_Par10 1.61 E−07 10 SS_Par12OC 5.58 E−03 12
SS_Par11 6.17 E−04 11 SS_Par13OC 7.27 E−05 13
SS_Par12 5.64 E−07 12 SS_Par14OC 4.76 E−03 14
SS_Par13 4.21 E−08 13 SS_Par15OC 1.33 E−04 15
SS_Par14 1.29 E−05 14 SS_Par16OC 4.00 E−01 16
SS_Par15 2.53 E−08 15 SS_Par17OC 1.79 E−06 17
SS_Par16 4.49 E−07 16 SS_Par18OC 1.50 E−01 18
SS_Par17 1.62 E−06 17 SS_Par19oc 2.37 E−05 19
SS_Par17up 1.00 E−07 17 SS_Par19O2 1.51 E−02 19
SS_Par18 5.15 E−07 18 SS_Par20OC 1.50 E−01 20
SS_Par19 9.50 E−05 19 SS_Par21OC 9.43 E−02 21
SS_Par20 3.26 E−07 21
SS_Par21 2.88 E−07 21
Table 9.    Final storage-coefficient parameters in the calibrated New Jersey Coastal Plain groundwater-flow model.
Red highlighting horizontal hydraulic conductivity in New Jersey and surrounding areas.
Figure 9.

Maps showing A, horizontal hydraulic conductivity zones of model layer 1 and B, vertical hydraulic conductivity of model layer 2, in the New Jersey Coastal Plain groundwater-flow model.

Light blue and various shades of red showing hydraulic conductivity values and zones
                           in New Jersey and surrounding areas.
Figure 10.

Maps showing A, horizontal hydraulic conductivity zones of model layer 3 and B, vertical hydraulic conductivity zones of model layer 4, in the New Jersey Coastal Plain groundwater-flow model.

Light blue and various shades of red showing hydraulic conductivity values and zones
                           in New Jersey and surrounding areas.
Figure 11.

Maps showing A, horizontal hydraulic conductivity zones of model layer 5 and B, vertical hydraulic conductivity zones of model layer 6, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 12.

Maps showing A, horizontal hydraulic conductivity zones of model layer 7 and B, vertical hydraulic conductivity zones of model layer 8, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 13.

Maps showing A, horizontal hydraulic conductivity zones of model layer 9 and B, vertical hydraulic conductivity zones of model layer 10, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 14.

Maps showing A, horizontal hydraulic conductivity zones of model layer 11 and B, vertical hydraulic conductivity zones of model layer 12, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 15.

Maps showing A, horizontal hydraulic conductivity zones of model layer 13 and B, vertical hydraulic conductivity zones of model layer 14, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 16.

Maps showing A, horizontal hydraulic conductivity zones of model layer 15 and B, vertical hydraulic conductivity zones of model layer 16, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 17.

Maps showing A, horizontal hydraulic conductivity zones of model layer 17 and B, vertical hydraulic conductivity zones of model layer 18, in the New Jersey Coastal Plain groundwater-flow model.

Various shades of blue and red showing horizontal hydraulic conductivity values and
                           zones in New Jersey and surrounding areas.
Figure 18.

Maps showing A, horizontal hydraulic conductivity zones of model layer 19 and B, vertical hydraulic conductivity zones of model layer 20, in the New Jersey Coastal Plain groundwater-flow model.

Map with various colors showing horizontal hydraulic conductivity in feet per day.
Figure 19.

Map showing horizontal hydraulic conductivity zones of model layer 21 in the New Jersey Coastal Plain groundwater-flow model.

The final parameters for horizontal and vertical hydraulic conductivity were compared to the values shown in tables 2 and 3 of this report. Many of the parameter values of horizontal hydraulic conductivity, such as those for the confined aquifers and the Holly Beach water-bearing zone in New Jersey, fall within the ranges of the value for horizontal hydraulic conductivity shown in the tables 2 and 3 and also table 1 of Masterson and others (2016a). Some horizontal hydraulic conductivity assigned to the outcrop areas of certain model layers for aquifers were higher (greater than 500 ft/d, table 6) because those areas are thin and the higher value increased model stability in these areas but did not affect the flow farther downdip. Higher values of horizontal hydraulic conductivity were also assigned to the outcrop areas of some confining-unit model layers (table 7) for the same reason. Lower values of horizontal hydraulic conductivity (less than 1 ft/d) were assigned in downdip areas of some aquifers—for example model layer 11 which represents the Vincentown aquifer—because the model layer is not present downdip and the small horizontal hydraulic conductivity value allowed for the lateral connection between the updip and downdip parts of the model layer and (or) a vertical connection between the model layer with the overlying and underlying model layers.

The parameters for riverbed conductance for the Delaware River cells and for streambed conductance for drains cells each have a single value. Although hydrogeologic and stream-width data could be used to justify multiple parameters for these properties, these data were considered secondary to the calibration of the regional-scale Coastal Plain model. During model calibration, simulation results indicated the model was insensitive to the riverbed-conductance parameter. Additionally, the comparison of simulated to estimated base flows was reasonable when a single streambed-conductance parameter for drain cells was used; therefore, neither parameter was subdivided to represent smaller or local areas for the river or drain cells. This process is described further in the “Parameter Sensitivity and Correlation” and “Parameter Estimation and Residual Analysis” sections below.

Two parameters were used to vary the recharge rates in the model calibration. One parameter was used for New Jersey and another parameter was used for recharge in Delaware and Maryland. The recharge parameter for New Jersey is a multiplier applied equally to the recharge value estimated for each model cell by using the Soil-Water-Balance model (described in the “Recharge” section earlier in this report). The recharge parameter was increased from an initial value of 1.0 to a final value of 1.2 to provide a reasonable match between simulated and estimated base flow in the 47 New Jersey surface-water basins (fig. 8).

Parameter Sensitivity and Correlation

A sensitivity analysis of the calibrated model was conducted to evaluate the relative effects of the various parameters on the match between simulated and observed water levels and simulated and estimated base flows. The sensitivity analysis is used to identify which parameters in the model, when modified, improve model results. Parameters which are least sensitive, when modified, have little effect on model results. The sensitivity of model parameters was calculated by using the parameter-estimation software UCODE–2014 (Poeter and others, 2014). Parameter sensitivities were calculated two ways: by using base-flow observations only and by using water-level observations only. This section discusses composite scaled sensitivities (CSS) and parameter correlation coefficients. The CSS values show how sensitive water-level and base-flow observations in the model are to changes in each of the parameters.

According to composite-scaled sensitivity analysis, the parameters most sensitive to base-flow observations (excluding recharge) are primarily streambed conductance, horizontal hydraulic conductivity, and storage in zones in the outcrop areas of several model layers (fig. 20). Recharge is the dominant parameter with respect to base-flow observations, and the recharge multiplier parameter was adjusted manually such that average simulated base flow equaled the estimated base flow of the New Jersey Coastal Plain of about 16 inches per year (in/yr), a value determined from surficial-aquifer studies completed in surface-water basins in the New Jersey Coastal Plain during 1992–2003 (Watt and Johnson, 1992; Watt and others, 1994; Lacombe and Rosman, 1995; Johnson and Watt, 1996; Johnson and Charles, 1997; Charles and others, 2001; Watt and others, 2003; Gordon, 200418). The streambed-conductance parameter for the drains (DRN_Par1) is ranked the 4th most sensitive parameter to base-flow observations (fig. 20). As indicated above, there is only one parameter for streambed conductance and any change to that parameter is applied to all streams uniformly rather than to individual streams. When the streambed conductance was decreased, groundwater levels tended to increase but the location or timing of discharge to base flow did not necessarily change. The 21 most sensitive parameters (fig. 20) other than streambed conductance (DRN_Par1) are in, or laterally or vertically adjacent to, an outcrop area.

Bar graph of composite scaled sensitivities with base-flow observations.
Figure 20.

Graph showing composite scaled sensitivities of model parameters derived from base-flow observations for the New Jersey Coastal Plain groundwater-flow model.

Results of the CSS analysis indicate that the parameters most sensitive to water-level observations are recharge, and horizontal and vertical hydraulic conductivity (fig. 21); riverbed and general-head-boundary conductances are insensitive. Of the 27 most sensitive parameters (including the two recharge parameters (RCH_SWB and RchDeFlux), 19 (70 percent) are horizontal hydraulic conductivity, and 6 (22 percent) are vertical hydraulic conductivity.

Bar graph of composite scaled sensitivities with water-level observations.
Figure 21.

Graph showing composite scaled sensitivities of model parameters derived from water-level observations for the New Jersey Coastal Plain groundwater-flow model.

Six pairs of the 130 most sensitive horizontal and vertical hydraulic conductivity parameters are correlated when recharge is excluded and only head observations are included in the parameter estimation (fig. 22). Most of the correlations are between vertically adjacent zones (for example, layers 20 and 21 in Delaware (figs. 18B and 19)) or horizontal and vertical hydraulic conductivity of the same zone (for example, layer 11 in northwestern Ocean County and northeastern Burlington County (figs. 1424.)). Four of the six pairs are parameter zones in Delaware and three of the six pairs are parameter zones in outcrop areas. Although the parameters were correlated, each correlated parameter was allowed to be estimated separately. There is not sufficient water-level data in the parameter zones to resolve the negative correlation between the horizontal hydraulic conductivity and (or) vertical hydraulic conductivity of each of the correlated parameter pairs.

Bar graph with blue bars showing negative correlation.
Figure 22.

Graph showing model parameter pairs in the New Jersey Coastal Plain groundwater-flow model with greater than 90 percent correlation.

Parameter Estimation and Residual Analysis

The parameters used in the New Jersey Coastal Plain model were tested and adjusted in an iterative calibration process. Initial parameter-estimation runs with UCODE–2014 (Poeter and others, 2014) were conducted by using all observations (water-level and base-flow). Results indicated that matching water levels in the aquifers would be less likely if base-flow observations were included. The challenge of adjusting the weights of base-flow observations to avoid their having too little or too great an effect on water levels was judged to be not effective because of the single streambed parameter for drains and the relative coarseness of the representation of the stream network in the DRN package. Therefore, after some adjustments were made to sensitive parameters to obtain an initial optimal calibration with all observations, parameter-estimation runs were done by using only base-flow observations. The recharge multiplier had the greatest effect on simulated base flows; once the recharge multiplier was held constant, the model response was generally more sensitive to parameters in the outcrop areas of aquifers and confining units. Subsequent parameter-estimation runs were done by using only water-level observations. For these runs, the parameters that were sensitive to drain (base-flow) observations, but not to water-level observations, typically remained at the value determined in the previous parameter-estimation runs in which only drain observations were used.

The number of simulated and water-level observations from wells in New Jersey used in the NJCP groundwater-flow model each total 4,243. The total includes 2,361 water-level observations from 392 wells, 920 historical water-level trends for 29 wells during 1980–2013, 210 observations of vertical gradients in 33 pairs of nested wells, and 752 observations of changes in water levels over time calculated for 134 wells. When plotted against each other (fig. 23), these values are clustered around the 1:1 correlation line where the observed water levels have a range of 411 ft and more than 96 percent of the simulated water levels fall within +/− 10 percent of the total range during 1980–2013. The residuals (simulated minus observed water level) are within + or – 10 feet for 58 percent of the water-level observations for New Jersey, and within + or – 20 feet for 84 percent. The average water-level residual is 1.5 ft; the minimum; first, second, and third quartile; and maximum residuals are −127.4, −5.5, 1.6, 10.4, and 179.4 ft, respectively.

Scatter plot with blue dots showing simulated vs observed water levels.
Figure 23.

Graph showing simulated and observed water levels in wells in New Jersey used in calibration for the New Jersey Coastal Plain groundwater-flow model, 1980–2013.

Calibration of Stream Base Flows

The simulated average annual base flows and estimated annual average base flows for the period of record for 47 surface-water basins are shown in figure 24. Hydrographs of simulated average annual base flows at continuous- and partial-record gaging stations during the period of record available for 1980–2013 generally match estimated average annual base flows and are shown in figures 25 and 26, respectively. Although the annual variations are not well correlated and simulated flows are higher than estimated base flow in some basins and lower in others, the overall difference between simulated and estimated average base flow was reasonably matched for most basins (figs. 25 and 26).

Blue showing low values and purple showing high values of base flows in New Jersey
                        and surrounding areas.
Figure 24.

Map showing simulated average and estimated average annual base flow for the period of record for selected surface-water basins in the New Jersey Coastal Plain.

Scattered black boxes and yellow circles showing estimated and simulated base flows.
Figure 25.

Hydrographs of simulated and estimated annual base flow for selected continuous-record gaging stations in the New Jersey Coastal Plain, 1980–2013.

Scattered black boxes and yellow circles showing estimated and simulated base flows.
Figure 26.

Hydrographs of simulated and estimated annual base flow for selected partial-record gaging stations in the New Jersey Coastal Plain, 1980–2013.

The residual (simulated average annual base flow minus estimated annual average base flow) for each of the 47 surface-water basins in New Jersey also is shown in table 10. Simulated average annual base flows in 37 of the 47 surface-water basins (79 percent) are within 50 percent of the estimated average annual base flows and simulated average base flows in 25 basins (53 percent) are within 25 percent. Estimated average annual base flow ranges from 0.8 to 171.7 cubic feet per second (ft3/s) and simulated average base flows range from 0.0 to 149.7 (ft3/s). The largest percent residuals in the table generally are associated with the small streamflows in basins with drainage areas 20 mi2 or less (7 of 20 surface-water basins (30 percent)). In some instances, simulated average base flow is lower than estimated base flow in the surface-water basins in updip basins (fig. 24), and simulated average base flow is higher than estimated base flow in adjacent downdip basins. For example, simulated average base flow at the Toms River continuous-record gaging station (01408500) is −34.9 ft3/s (20 percent) less than estimated average annual base flow, whereas simulated average base flow at gaging stations in the adjacent basin of Crosswicks Creek (01464500) is 26.0 ft3/s (31 percent), respectively, greater than the estimated value. Simulated base flow from this model and the previous RASA model (Voronin, 2004) are compared during model calibration to ensure that the level of calibration was maintained or improved. Voronin (2004) reports that simulated average base flows at the gaging stations at Toms River (01408500) and North Branch Rancocas Creek (0146700) are within about 10 percent of estimated values for 1926–96 and 1921–96, respectively. Simulated average base flows at Toms River (01408500) and North Branch Rancocas Creek (01467000) in this report are within −20 and 8 percent of estimated values, respectively, for 1980–2013. Similarly, Voronin (2004) reports that simulated average base flow of Raccoon Creek (01477120) is within 31 percent of the estimated value for 1966–92; simulated average base flow for this station in this report, however, is within 22 percent of that estimated for 1980–2013.

Table 10.    

Simulated average and estimated average annual base flow for the period of record at selected continuous and partial-record streamflow-gaging stations in the New Jersey Coastal Plain.

[ft3/s, cubic feet per second; NJ, New Jersey]

Station number Station name Drainage area (square miles) Model station identifier Simulated average annual base flow (ft3/s) Estimated average annual base flow (ft3/s) Residual (ft3/s) Percent residual
Surface-water basin with continuous streamflow-gaging station
01407500 Swimming River near Red Bank NJ 49.2 CRSwimmin 63.5 24.2 39.3 162
01407705 Shark River near Neptune City NJ 9.96 CRShark 16.2 8.4 7.8 93
01407760 Jumping Brook near Neptune City NJ 6.46 CRJumping 6.6 4.6 2.0 43
01408029 Manasquan River near Allenwood NJ 63.3 CRManasqu 77.4 143.1 −65.7 −46
01408120 North Branch Metedeconk River near Lakewood NJ 34.9 CRNBMeted 37.8 41.5 −3.7 −9
01408500 Toms River near Toms River NJ 123 CRTomsR 136.8 171.7 −34.9 −20
01408900 Cedar Creek near Lanoka Harbor NJ 53.1 CRCedar 76.5 81.4 −4.9 −6
01409280 Westecunk Creek at Stafford Forge NJ 15.8 CRWestecu 24.4 29.5 −5.1 −17
01409810 West Branch Wading River near Jenkins NJ 84.1 CRWBWadin 98.5 108.9 −10.4 −10
01410150 East Branch Bass River near New Gretna NJ 8.11 CREBBass 14.6 14.9 −0.3 −2
01410500 Absecon Creek at Absecon NJ 17.9 CRAbsecon 20.6 13.5 7.1 53
01410000 Oswego River at Harrisville NJ 72.5 CROswego 71.9 68.5 3.4 5
01411000 Great Egg Harbor River at Folsom NJ 57.1 CRGEHR 58.2 73.2 −15.0 −21
01411300 Tuckahoe River at Head of River NJ 30.8 CRTuckaho 36.9 33.7 3.2 9
01411500 Maurice River at Norma NJ 112 CRMaurice 128.8 143.1 −14.3 −10
01412800 Cohansey River at Seeley NJ 28 CRCohanse 34.5 29.4 5.1 17
01464500 Crosswicks Creek at Extonville NJ 81.5 CRCrosswi 110.6 84.6 26.0 31
01465850 South Branch Rancocas Creek at Vincentown NJ 64.5 CRSBRanco 80.2 68.2 12.0 18
01467000 North Branch Rancocas Creek at Pemberton NJ 118 CRNBRanc 149.7 138.3 11.4 8
01467150 Cooper River at Haddonfield NJ 17 CRCooper 25.2 17.0 8.2 48
01477120 Raccoon Creek near Swedesboro NJ 26.9 CRRacoon 33.7 27.6 6.1 22
Surface-water basin with partial record streamflow-gaging station
01408150 South Branch Metedeconk River near Lakewood NJ 27.4 PRSBMeted 41.4 45.4 −4.0 −9
01408600 Wrangle Brook near Toms River NJ 19.4 PRWrangle 25.2 34.7 −9.5 −27
01408630 Davenport Branch near Toms River NJ 12.3 PRDavenpo 17.1 18.1 −1.0 −6
01409411 Nescochague Creek at Pleasant Mills NJ 43.6 PRNescoch 54.7 58.3 −3.6 −6
01409575 Landing Creek at Phila Ave at Egg Harbor City NJ 8.3 PRLanding 12.3 5.8 6.5 112
01411071 Hospitality Branch at RR bridge near Folsom NJ 44.9 PRHospita 52.2 45.9 6.3 14
01411140 Deep Run at Weymouth NJ 20.0 PRDeepRun 23.7 25.5 −1.8 −7
01411220 South River near Belcoville NJ 20.4 PRSouthR 26.6 27.4 −0.8 −3
01411700 Muddy Run at Centerton NJ 37.7 PRMuddyRu 47.2 34.7 12.5 36
01412005 Menantico Creek at Route 49 at Millville NJ 26.4 PRMenanti 29.9 35.6 −5.7 −16
01464515 Doctors Creek at Allentown NJ 17.5 PRDoctors 22.8 15.3 7.5 49
01464530 Blacks Creek at Mansfield Square NJ 19.6 PRBlacks 25.8 14.2 11.6 82
01464540 Crafts Creek at Hedding NJ 10.5 PRCrafts 9.2 5.1 4.1 80
01464587 Assiscunk Creek at Jacksonville NJ 32.4 PRAssiscu 40.8 19.8 21.0 106
01465865 Barton Run at Tuckerton Road near Medford NJ 12.8 PRBarton 17.7 12.4 5.3 43
01465880 Southwest Branch Rancocas Creek at Medford NJ 47.1 PRSWBRanc 66.2 51.1 15.1 30
01465893 Little Creek at Chairville NJ 6.4 PRLittle 3.0 5.5 −2.5 −45
01465930 Masons Creek near Springville NJ 1.1 PRMasons 1.1 0.9 0.2 22
01467010 Parkers Creek near Mount Laurel NJ 2.5 PRParkers 2.8 0.8 2.0 250
01467020 Mill Creek at Willingboro NJ 7.8 PRMills 0.0 3.0 −3.0 −100
0146702680 Swede Run at Conrow Road at Delran NJ 4.5 PRSwede 0.0 2.5 −2.5 −100
01467057 Pompeston Creek at Cinnaminson NJ 5.8 PRPompest 1.0 1.4 −0.4 −29
01467330 South Branch Big Timber Creek at Blackwood NJ 20.8 PRSBBigTi 23.3 27.4 −4.1 −15
01467359 North Branch Big Timber Creek at Glendora NJ 18.6 PRNBBigTi 39.0 29.6 9.4 32
01475020 Mantua Creek at Sewell NJ 14.5 PRMantua 17.4 17.2 0.2 2
01475100 Edwards Run near Mantua NJ 6.4 PREdwards 6.3 8.2 −1.9 −23
Table 10.    Simulated average and estimated average annual base flow for the period of record at selected continuous and partial-record streamflow-gaging stations in the New Jersey Coastal Plain.

Calibration of Water Levels

During model calibration, simulated water levels were compared to water levels measured in the seven quinquennial New Jersey Coastal Plain synoptic studies conducted in 1983, 1988, 1993, 1998, 2003, 2008, and 2013 to ensure that the cones of depression and flow directions were simulated accurately. Discussion of simulated water levels in this report is limited to 2013 groundwater-flow conditions. Maps of the 2013 potentiometric surface for the confined aquifers in the New Jersey Coastal Plain (Cauller and Gordon, 2021) and the simulated water levels for the major aquifers in the New Jersey Coastal Plain and the confined aquifers in Delaware and residuals at wells measured in those aquifers in 2013 are shown in this section. As previously mentioned in the “Model Development” section, if the aquifer does not extend across a model layer in the MODFLOW-2005 finite difference solution method, but instead pinches out, the absent area of the aquifer is represented in the model by a thin layer that allows for the hydraulic connection between underlying and overlying model layers. These aquifers include the downdip areas of the Vincentown and lower PRM aquifers, and the Englishtown aquifer system, (model layers 11, 21, and 15, respectively), and the updip areas of the Rio Grande water-bearing zone, Atlantic 800-foot sand, and Piney Point aquifer (model layers 5, 7, and 9). Therefore, the simulated water levels for the part of a model layer representing the aquifer area only are shown on the figures for these aquifers. The residuals are the simulated water level minus the measured (observed) water level at the well and the root mean square error (RMSE) statistic of the residuals for each aquifer is computed. The RMSE is the standard deviation of the residuals. Hydrographs of water levels also were used to calibrate the model. The hydrographs were compared to simulation results from Voronin (2004) to ensure that the updated model maintained or improved the level of calibration achieved with the previous model.

Calibration of Water Levels in the Holly Beach Water-Bearing Zone

Simulated water-level contours in the Holly Beach water-bearing zone (model layer 1, Cape May County) are shown in figure 27. The Holly Beach water-bearing zone is a minor aquifer in New Jersey; the observed water levels shown in Lacombe and Carleton (2002) were used to visually compare to the simulated water-level contours.

Blue lines in southern New Jersey showing simulated water-table contours; tan showing
                           aquifer limit.
Figure 27.

Map showing simulated water-level contours for the Holly Beach water-bearing zone, Cape May County, New Jersey, 2013.

Calibration of Water Levels in the Confined Cohansey Aquifer and Unconfined Kirkwood-Cohansey Aquifer System

Simulated water-level contours in the unconfined Kirkwood-Cohansey aquifer shown in the updip area of model layer 3 closely match the observed water-table contours in general direction and slope (fig. 28). The observed water-table contours are an aggregate of data from multiple surficial-aquifer studies in the following surface-water basins of the New Jersey Coastal Plain: Great Egg Harbor River Basin (Watt and Johnson, 1992); Toms River, Metedeconk River, and Kettle Creek Basins (Watt and others, 1994); upper Maurice River Basin (Lacombe and Rosman, 1995); Mullica River Basin (Johnson and Watt, 1996); Salem River and Raccoon, Oldmans, Alloway, and Stow Creek Basins (Johnson and Charles, 1997); and Rancocas, Crosswicks, Assunpink, Blacks, and Crafts Creek Basins (Watt and others, 2003). The aggregated (observed) water-table contours are shown only in the areas where a surficial-aquifer study was completed before 2004. In Cape May County, New Jersey, as well as in Sussex County, Delaware, groundwater is confined.

Map showing a model grid outline, aquifer limits, updip limits, and multicolored contour
                           lines.
Figure 28.

Map showing simulated and composite observed water-table contours for the unconfined Kirkwood-Cohansey aquifer system in New Jersey.

The updip limit of the confined Cohansey aquifer in Cape May County is shown in figure 29. The simulated water-level contours in the confined Cohansey aquifer in New Jersey and the potentiometric surface of the 2013 New Jersey Coastal Plain water-level synoptic study (Cauller and Gordon, 2021) for this aquifer also are shown in figure 29. Simulated water levels in the confined part of the Cohansey aquifer (model layer 3, Cape May County) generally match observed water levels (table 11). The average water-level residual (simulated minus observed) for 17 observation points in Cape May County for the confined Cohansey aquifer is −0.2 ft and the root mean square error (RMSE) is 3.3 ft. The simulated water levels for the Upper Chesapeake aquifer in Delaware and the average water-level residual for four wells measured in Delaware in 2013 (table 11) also are shown in figure 29.

Table 11.    

Simulated and observed water-level observations and residuals for wells open to the confined Cohansey aquifer, Cape May County, New Jersey, and the Upper Chesapeake aquifer, Delaware, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (in feet) Simulated water level (in feet) Residual (in feet)
New Jersey wells
WLOb090048_7 090048 −10.36 –8.09 2.27
WLOb090049_7 090049 −10.42 –7.54 2.88
WLOb090054_7 090054 –13.90 –11.12 2.78
WLOb090060_7 090060 –9.06 –10.98 –1.92
WLOb090080_7 090080 –1.89 –4.02 –2.13
WLOb090089_7 090089 –2.60 –3.62 –1.02
WLOb090099_7 090099 3.55 3.06 –0.49
WLOb090150_7 090150 –8.89 –5.04 3.85
WLOb090187_6 090187 –5.95 –7.15 –1.20
WLOb090210_6 090210 –7.42 8.87 –1.45
WLOb090219_5 090219 1.41 0.28 –1.13
WLOb090292_6 090292 2.71 –0.31 –3.02
WLOb090315_3 090315 4.54 –1.91 –6.45
WLOb090354_5 090354 3.10 –2.42 –5.52
WLOb090358_4 090358 –1.14 3.73 4.87
WLOb090366_4 090366 –3.52 –3.59 –0.07
WLOb090395_5 090395 –12.91 –7.91 5.00
Delaware wells
WLObOh2502_6 Oh2502 –4.80 5.54 10.34
WLObOh2503_6 Oh2503 1.00 5.55 4.55
WLObPe5403_5 Pe5403 44.95 43.68 –1.27
WLObQj4104_6 Qj4104 3.63 0.25 –3.38
Table 11.    Simulated and observed water-level observations and residuals for wells open to the confined Cohansey aquifer, Cape May County, New Jersey, and the Upper Chesapeake aquifer, Delaware, 2013.
Map showing a model grid outline, aquifer limits, updip limits, and multicolored contour
                           lines.
Figure 29.

Map showing simulated water-level contours and 2013 potentiometric surface of, and residuals for selected wells in, the confined Cohansey aquifer, Cape May County, New Jersey, and simulated water-level contours for the Upper Chesapeake aquifer, Delaware.

The hydrographs of simulated and observed water levels for well 090049 (fig. 30) are similar, with minima, maxima, and averages for the observed water levels of –21, –8, and –13 ft, respectively, for the observed water levels and –19, –7, and –11 ft, respectively, for the simulated water levels. Simulated water levels are, on average, 2 ft greater than observed. The hydrographs of simulated and observed water levels for well 090080 in south-central Cape May County (fig. 31) are similar, with minima, maxima, and averages for the observed water levels of –11, –1, and –6 ft, respectively, for the observed water levels and –18, –0.5, and –8 ft, respectively, for the simulated water levels.

Scatter plot showing simulated versus measured water levels.
Figure 30.

Hydrograph of simulated and observed water levels in well 090049 open to the confined Cohansey aquifer, New Jersey, 1980–2013.

Scatter plot showing simulated versus measured water levels.
Figure 31.

Hydrograph of simulated and observed water levels in well 090080 open to the confined Cohansey aquifer, New Jersey, 1980–2013.

Calibration of Water Levels in the Rio Grande Water-Bearing Zone

Simulated water-level contours and the potentiometric surface of the 2013 New Jersey Coastal Plain water-level synoptic study (Cauller and Gordon, 2021) in the Rio Grande water-bearing zone (model layer 5, parts of Ocean, Atlantic, and Cape May Counties) are shown in figure 32. The 2013 simulated water levels are an average 7 ft lower than the observed water levels at six wells in Cape May County and are 35 ft lower than the observed water levels at three wells in Ocean County. The simulated 2013 potentiometric contours indicate a cone of depression along the shoreline of Cape May and Atlantic Counties caused by drawdowns in the underlying Atlantic City 800-foot sand; a cone of depression also is observed in the simulated water-level contours. The differences between the simulated and observed water levels for 10 observation points in New Jersey are listed in table 12. The average water-level residual (simulated minus observed) for the 10 observation points (figs. 32 and 33) is −16 ft and the RMSE is 24 ft. The minimum; first, second, and third quartile; and maximum residuals are −45, −31, −13, 0, and 5 ft, respectively.

Table 12.    

Simulated and observed water-level observations and residuals for wells open to the Rio Grande water-bearing zone, New Jersey, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (feet) Simulated water level (feet) Residual (feet)
WLOb090071_6 090071 –23.19 –32.57 −9.38
WLOb090149_6 090149 16.14 –18.08 –34.22
WLOb090304_4 090304 –20.83 –26.83 –6.00
WLOb090305_3 090305 –27.34 –24.39 2.95
WLOb090526_3 090526 –26.40 –24.72 1.68
WLOb090629_2 090629 –32.47 –27.70 4.77
WLOb110737_4 110737 –4.39 –20.29 –15.90
WLOb290775_6 290775 –1.97 –42.86 –40.89
WLOb290813_5 290813 0.01 –44.54 –44.55
WLOb291621_3 291621 –22.80 –42.43 –19.63
Table 12.    Simulated and observed water-level observations and residuals for wells open to the Rio Grande water-bearing zone, New Jersey, 2013.
Map showing water level contours in light blue and brown, as well as aquifer and updip
                           limits in burgundy and dark blue within the model grid outline in purple.
Figure 32.

Map showing simulated water-level contours and 2013 potentiometric surface for, and residuals for selected wells in, the Rio Grande water-bearing zone, New Jersey and simulated water-level contours for the Lower Chesapeake aquifer, Delaware.

Red dots showing water level, light and dark blue lines showing lines of correspondence.
Figure 33.

Graph showing simulated and observed water levels in wells open to the Rio Grande water-bearing zone, Cape May, Atlantic, and Ocean Counties, New Jersey, 2013.

Calibration of Water Levels in the Atlantic City 800-Foot Sand

Simulated water-level contours in the Atlantic City 800-foot sand (model layer 7, parts of Cape May, Atlantic, and Ocean Counties) generally match the 2013 potentiometric surface in coastal Atlantic and Cape May Counties (fig. 34). The simulated water-level contours are similar in shape to the 2013 potentiometric surface (Cauller and Gordon, 2021), although the gradient of the simulated contours is less steep than that of the water-level contours drawn from measurements made in updip central Atlantic County. The simulated water levels are within 10 ft of observed levels in 13 of 34 wells (38 percent) in New Jersey and within 15 ft of observed levels in 24 of 34 wells (71 percent) (figs. 34 and 35 and table 13). The 2013 simulated water levels are an average 8 ft higher than the observed water levels at 13 wells in Cape May County and are an average 13 ft lower than the observed water levels at 14 wells in Atlantic County. The simulated water levels are an average 11 ft lower than observed levels at 7 wells in Ocean County. The average water-level residual (simulated minus observed) for 34 wells in Cape May, Atlantic, and Ocean Counties is −4 ft and the RMSE is 21 ft. The minimum; first, second, and third quartile; and maximum residuals are −57, −12, 3, 12, and 17 ft, respectively. The simulated water levels for the Calvert aquifer system in Delaware and the average water-level residual for one well measured in Delaware in 2013 (table 13) also are shown in figure 34.

Table 13.    

Simulated and observed water-level observations and residuals for wells open to the Atlantic City 800-foot sand, New Jersey, and the Calvert aquifer system, Delaware, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (in feet) Simulated water level (in feet) Residual (in feet)
New Jersey wells
WLOb010037_7 010037 −74.45 −56.98 17.47
WLOb010039_7 010039 −71.02 −57.85 13.17
WLOb010180_7 010180 −45.40 −57.76 −12.36
WLOb010367_7 010367 −78.10 −63.01 15.09
WLOb010578_7 010578 −65.68 −55.73 9.95
WLOb010600_7 010600 −82.67 −67.60 15.07
WLOb010650_6 010650 12.33 −45.09 −57.42
WLOb010704_6 010704 −57.58 −70.29 −12.71
WLOb010706_6 010706 −29.39 −59.66 −30.27
WLOb010967_5 010967 −57.71 −51.79 5.92
WLOb011220_4 011220 −41.90 −95.65 −53.75
WLOb011252_3 011252 14.81 −36.07 −50.88
WLOb011253_4 011253 −43.30 −66.53 −23.23
WLOb011456_2 011456 −61.83 −73.09 −11.26
WLOb090004_7 090004 −49.54 −37.53 12.01
WLOb090079_5 090079 −42.63 −31.97 10.66
WLOb090092_7 090092 −44.16 −36.36 7.80
WLOb090106_7 090106 −61.48 −45.74 15.74
WLOb090136_7 090136 −56.20 −41.81 14.39
WLOb090144_7 090144 −63.25 −52.88 10.37
WLOb090185_5 090185 −44.03 −39.73 4.30
WLOb090302_5 090302 −31.37 −27.28 4.09
WLOb090306_5 090306 −30.05 −30.08 −0.03
WLOb090337_5 090337 −32.42 −28.62 3.80
WLOb090423_5 090423 −33.28 −32.08 1.20
WLOb090479_4 090479 −43.56 −31.35 12.21
WLOb090527_2 090527 −49.98 −37.96 12.02
WLOb290112_7 290112 −28.12 −37.85 −9.73
WLOb290457_5 290457 −39.93 −41.72 −1.79
WLOb290561_7 290561 −31.73 −37.48 −5.75
WLOb290597_5 290597 −18.08 −44.521 −26.43
WLOb290598_6 290598 −23.62 −28.78 −5.16
WLOb290814_5 290814 −28.67 −52.98 −24.31
WLOb290936_6 290936 −33.37 −39.62 −6.25
Delaware well
WLObKb3102_7 Kb3102 27.67 −19.81 −47.48
Table 13.    Simulated and observed water-level observations and residuals for wells open to the Atlantic City 800-foot sand, New Jersey, and the Calvert aquifer system, Delaware, 2013.
Map showing water level contours in light blue and brown, as well as aquifer and updip
                           limits in burgundy and dark blue within the model grid outline in purple.
Figure 34.

Map showing simulated water–level contours and 2013 potentiometric surface for, and residuals for selected wells in, the Atlantic City 800-foot sand, New Jersey, and the simulated water-level contours for, and residuals for a selected well in, the Calvert aquifer system, Delaware.

Red dots showing water level, light and dark blue lines showing lines of correspondence.
Figure 35.

Graph showing simulated and observed water levels for wells open to the Atlantic City 800-foot sand, Cape May, Atlantic, and Ocean Counties, New Jersey, 2013.

The hydrographs of wells 010180 and 010578, in northeastern and southeastern Atlantic County, respectively (fig. 36), show the same trend and amplitude as the observed water levels, but the simulated water levels are, on average, 1 ft lower and 12 ft higher, respectively, than the observed. The hydrographs of simulated and observed water levels for well 010180 (fig. 36A) are similar, with minima, maxima, and average for the observed water levels of –57, –35, and –47 ft, respectively, for the observed water levels and –60, –36, and –48 ft, respectively, for the simulated water levels. The hydrographs of simulated and observed water levels for well 010158 (fig. 36B) are similar, with minima, maxima, and average for the observed water levels of –71, –49, and –60 ft, respectively, and –60, –38, and –48 ft, respectively, for the simulated water levels. The hydrograph for well 010578 (fig. 36B) shows simulated water levels that are closer to and more closely match the trend in the observed water levels than those simulated by Voronin (2004) for 1990–98.

Scatter plot showing simulated versus measured water levels in wells A, 010180 and
                           B, 010578.
Figure 36.

Hydrographs of simulated and observed water levels in wells A, 010180 and B, 010578, open to the Atlantic City 800-foot sand, New Jersey, 1980–2013.

Calibration of Water Levels in the Piney Point Aquifer

Simulated water-level contours in the Piney Point aquifer (model layer 9, parts of Monmouth, Ocean, Atlantic, Gloucester, Cumberland, and Salem Counties, New Jersey, and Kent County in Delaware) are generally similar in shape to the 2013 potentiometric surface (Cauller and Gordon, 2021) (fig. 37), including the contours defining cones of depression in Kent County in Delaware and Cumberland and Ocean Counties in New Jersey. The simulated water levels extend farther updip in Burlington and Ocean Counties than the 2013 potentiometric surface because the updip limit of the Piney Point aquifer in the new hydrogeologic framework used in this groundwater-flow model used new well information from that area which was not available at the time of the previous framework developed by Zapecza, (1989). However, the 2013 potentiometric surface used the updip limit of the Piney Point aquifer which was based on the earlier hydrogeologic framework presented in Zapecza, (1989) for that aquifer. Available Piney Point aquifer withdrawal data indicate anomalously low rates in 2004, 2008, and 2010, likely resulting from missing data for wells with reported withdrawals in other years. Therefore, the 2007 withdrawal rates were substituted for 2004, 2008, and 2010 withdrawals for wells with presumed missing data. The average water-level residuals (simulated minus observed) are shown in table 14 and figures 37 and 38 for 26 wells in New Jersey. The average water-level residual for 26 wells in Cumberland, Atlantic, Gloucester, Camden, Burlington, and Ocean Counties in New Jersey is 9 ft and the RMSE is 23 ft. The minimum; first, second, and third quartile; and maximum residuals are –22, –2, 4, 16, and 80 ft, respectively. The simulated water levels are within 10 ft of observed water levels for 14 of the 26 wells (54 percent) in New Jersey and within 20 ft of observed levels in 19 of 26 wells (73 percent) (figures 37 and 38 and table 14). The average water-level residuals also are shown in table 14 and figure 37 for four wells in Delaware.

Table 14.    

Simulated and observed water-level observations and residuals for wells open to the Piney Point aquifer, New Jersey and Delaware, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (feet) Simulated water level (feet) Residual (feet)
New Jersey wells
WLOb010270_7 010270 19.94 27.27 7.33
WLOb010700_6 010700 16.95 −3.29 −20.24
WLOb010713_6 010713 −11.40 −7.58 3.82
WLOb010834_6 010834 −40.26 −31.66 8.60
WLOb011219_5 011219 −22.84 −22.45 0.39
WLOb050407_7 050407 48.88 42.40 −6.48
WLOb050488_7 050488 41.21 23.80 −17.41
WLOb050676_7 050676 117.76 133.69 15.93
WLOb050800_7 050800 72.73 73.46 0.73
WLOb071147_2 071147 52.25 78.15 25.90
WLOb071280_2 071280 42.17 55.07 12.90
WLOb110044_7 110044 −84.94 −32.29 52.65
WLOb110092_7 110092 −81.79 −55.56 26.21
WLOb110096_7 110096 −68.59 −52.58 16.01
WLOb110163_7 110163 −62.21 −32.37 29.84
WLOb111220_3 111220 −163.51 −83.62 79.89
WLOb151592_2 151592 −6.85 11.62 18.47
WLOb290018_7 290018 −4.61 −6.59 −1.98
WLOb290023_5 290023 −37.08 −34.82 2.26
WLOb290425_5 290425 117.25 117.25 0.00
WLOb290537_7 290537 −15.57 −22.21 −6.64
WLOb290585_6 290585 9.38 9.29 −0.09
WLOb290607_7 290607 −47.39 −52.11 −4.72
WLOb290739_7 290739 13.30 18.84 5.54
WLOb291210_4 291210 −20.48 −42.62 −22.14
WLOb291579_3 291579 23.62 28.14 4.52
Delaware wells
WLObId5501_7 Id5501 −124.39 −163.76 −39.37
WLObJd1415_6 Jd1415 −140.00 −179.04 −39.04
WLObKc3101_5 Kc3101 −66.81 −61.66 5.15
WLObNc1303_6 Nc1303 −36.71 −41.91 −5.20
Table 14.    Simulated and observed water-level observations and residuals for wells open to the Piney Point aquifer, New Jersey and Delaware, 2013.
Map showing water level contours in light blue and brown, as well as aquifer and updip
                           limits in burgundy and dark blue within the model grid outline in purple.
Figure 37.

Map showing simulated water–level contours and the 2013 potentiometric surface for, and residuals for selected wells in, the Piney Point aquifer, New Jersey and Delaware.

Red dots showing water level, light and dark blue lines showing lines of correspondence.
Figure 38.

Graph showing simulated and observed water levels for selected wells open to the Piney Point aquifer, New Jersey and Delaware, 2013.

Hydrographs for two wells screened in the Piney Point aquifer are shown in figure 39. Simulated water levels in well 050676, located in an area in southeastern Burlington County where withdrawals are minimal, have the same flat trend but a smaller amplitude than observed water levels and are, on average, 16 ft greater than observed water levels (fig. 39A). The water levels simulated for this well by Voronin (2004) also show a smaller amplitude than observed water levels and are 5 to 10 ft lower than the observed water levels during 1980–98. The simulated water levels in well 110096 in Cumberland County have the same decreasing trend and amplitude as the observed water levels (associated with increased nearby withdrawals), but the simulated water levels are on average 9 ft greater than observed water levels (fig. 39B). The water levels simulated for this well by Voronin (2004) show a pattern similar to that of the simulated water levels shown in figure 39B for 1980–98.

Simulated measurements were higher than actual measurements.
Figure 39.

Hydrographs of simulated and observed water levels in wells A, 050676 and B, 110096, open to the Piney Point aquifer, New Jersey, 1980–2013.

Calibration of Water Levels in the Vincentown Aquifer

Simulated water level contours in the Vincentown aquifer (model layer 11, present primarily in Burlington, Gloucester, Monmouth, Ocean, and Salem Counties) generally match the 2013 potentiometric surface (Cauller and Gordon, 2021) (fig. 40). The differences between the simulated and observed water levels are shown in table 15 and in figures 40 and 41 for wells measured in New Jersey in 2013. The average water-level residual (simulated minus observed) for nine wells in Burlington, Gloucester, Monmouth, Ocean, and Salem Counties is −2 ft and the RMSE is 11 ft. The minimum; first, second, and third quartile; and maximum residuals are −9, −6, −2, 8, and 29 ft, respectively. The simulated water levels are within 10 ft of observed water levels for seven of the nine wells (89 percent) in New Jersey. The simulated water levels for the Rancocas aquifer in Delaware and the average water-level residual for one well measured in Delaware in 2013 (table 15) also are shown in figure 40.

Table 15.    

Simulated and observed water-level observations and residuals for wells open to the Vincentown aquifer, New Jersey, and the Rancocas aquifer, Delaware, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (in feet) Simulated water level (in feet) Residual (in feet)
New Jersey wells
WLOb051250_4 051250 101.02 96.22 −4.80
WLOb151005_5 151005 68.42 81.03 12.61
WLOb151360_4 151360 81.15 72.54 −8.61
WLOb250636_6 250636 70.70 68.32 −2.38
WLOb250717_4 250717 123.79 123.37 -0.42
WLOb250788_4 250788 28.65 19.39 -9.26
WLOb290230_5 290230 125.46 119.67 −5.79
WLOb290658_5 290658 90.93 99.51 8.58
WLOb330292_3 330292 45.19 74.31 29.12
Delaware well
WLObFc5128_4 Fc5128 11.38 11.17 −0.21
Table 15.    Simulated and observed water-level observations and residuals for wells open to the Vincentown aquifer, New Jersey, and the Rancocas aquifer, Delaware, 2013.
Map showing water level contours in light blue and brown, as well as aquifer and updip
                           limits in burgundy and dark blue within the model grid outline in purple.
Figure 40.

Map showing simulated water-level contours and the 2013 potentiometric surface for, and residuals for selected wells in, the Vincentown aquifer, New Jersey and the simulated water-level contours for, and residual for a selected well in, the Rancocas aquifer, Delaware.

Red dots showing water level, light and dark blue lines showing lines of correspondence.
Figure 41.

Graph showing simulated and observed water levels open to the Vincentown aquifer, New Jersey, 2013.

Hydrographs of simulated and observed water levels in wells 250636 and 290139 in Monmouth and Ocean Counties, respectively, show that simulated water levels are an average of 2 ft and about 3 ft lower than observed water levels, respectively, and have the same flat trend and similar amplitude (fig. 42). The hydrograph for well 290139 in Voronin (2004) is similar to the hydrograph shown in figure 42B for 1980–98, whereas the hydrograph for well 250636 (fig. 42A) in Voronin (2004) shows simulated water levels that are 20 to 25 ft lower than observed water levels during 1980–98.

Scatter plot showing simulated versus measured water levels in wells A, 250636 and
                           B, 290139.
Figure 42.

Hydrographs of simulated and observed water levels in wells, A, 250636, and B, 290139, open to the Vincentown aquifer, New Jersey, 1980–2013.

Calibration of Water Levels in the Wenonah-Mount Laurel Aquifer

Simulated water-level contours in the Wenonah-Mount Laurel aquifer (model layer 13, Monmouth, Ocean, Burlington, Camden, Gloucester, and Salem Counties) generally match the 2013 potentiometric surface (Cauller and Gordon, 2021), but the simulated water levels are lower in updip areas of Monmouth County and central Gloucester County (fig. 43). The simulated water-level contours generally are similar in shape to the 2013 potentiometric surface over much of the extent of the aquifer and the simulated water levels are within 10 or 20 ft of the water levels that were measured at the cones of depression associated with pumping centers in coastal Ocean and Monmouth Counties and western Burlington and central Camden Counties. The differences between the simulated and observed water levels for wells measured in New Jersey in 2013 are shown in table 16 and in figures 43 and 44. The average water-level residual (simulated minus observed) for 48 wells in Salem, Gloucester, Camden, Burlington, Ocean, and Monmouth Counties is −8 ft and the RMSE is 23 ft. The minimum; first, second, and third quartiles; and maximum residuals are –57, –23, –4, 7, and 55 ft, respectively. The simulated water levels are within 10 ft of observed water levels for 19 of the 48 wells (40 percent) and within 20 ft of observed water levels for 33 of the 48 wells (69 percent) in New Jersey. The simulated water levels for the Mount Laurel aquifer in Delaware and the average water-level residual for two wells measured in Delaware in 2013 (table 16) also are shown in figure 43.

Table 16.    

Simulated and observed water-level observations and residuals for wells open to the Wenonah-Mount Laurel aquifer, New Jersey, and the Mount Laurel aquifer, Delaware, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (feet) Simulated water level (feet) Residual (feet)
New Jersey wells
WLOb050257_7 050257 47.50 43.57 –3.93
WLOb050359_7 050359 26.95 24.31 –2.64
WLOb050365_7 050365 14.56 30.15 15.59
WLOb050427_7 050427 1.61 –14.93 –16.54
WLOb050695_7 050695 7.27 15.02 7.75
WLOb050720_6 050720 –16.40 1.71 18.11
WLOb050744_7 050744 2.44 18.42 15.98
WLOb051086_4 051086 34.43 44.04 9.61
WLOb051155_5 051155 23.73 25.70 1.97
WLOb051166_5 051166 99.60 74.70 –24.90
WLOb051178_5 051178 21.62 21.02 -0.60
WLOb051387_4 051387 −1.42 –12.42 –11.00
WLOb070022_6 070022 9.51 –15.53 –25.04
WLOb070118_7 070118 67.39 56.49 –10.90
WLOb070308_7 070308 46.9 13.34 –33.56
WLOb070391_7 070391 −3.95 –37.21 –33.26
WLOb070421_7 070421 85.67 40.59 –45.08
WLOb070449_6 070449 −22.10 –35.86 –13.76
WLOb070478_7 070478 −20.58 –8.31 12.27
WLOb150687_6 150687 20.66 18.97 –1.69
WLOb150910_6 150910 82.82 27.39 –55.43
WLOb150953_6 150953 55.67 38.67 –17.00
WLOb151009_6 151009 64.28 6.90 –57.38
WLOb151060_5 151060 13.46 –7.61 –21.07
WLOb151104_5 151104 80.29 56.06 –24.23
WLOb151384_3 151384 −20.06 –24.69 –4.63
WLOb250014_7 250014 −54.47 –51.47 3.04
WLOb250088_6 250088 109.86 76.62 –33.24
WLOb250095_5 250095 140.62 95.55 –45.07
WLOb250168_6 250168 −1.23 53.34 54.57
WLOb250335_7 250335 −51.05 –28.76 22.29
WLOb250353_7 250353 14.08 15.09 1.01
WLOb250391_7 250391 −62.21 –53.42 8.79
WLOb250396_7 250396 84.72 62.54 –22.18
WLOb250405_7 250405 138.12 95.72 –42.40
WLOb250486_6 250486 −61.92 –55.63 6.29
WLOb290031_7 290031 −55.36 –48.34 7.02
WLOb290049_7 290049 −66.43 –46.60 19.83
WLOb290140_7 290140 109.52 115.31 5.79
WLOb290699_7 290699 120.17 127.76 7.59
WLOb290783_7 290783 46.95 67.11 20.16
WLOb330002_7 330002 16.70 4.21 –12.49
WLOb330020_7 330020 23.06 7.16 –15.90
WLOb330050_7 330050 4.51 7.64 3.13
WLOb330252_7 330252 −1.17 –11.25 –10.08
WLOb330842_4 330842 14.32 5.03 –9.29
WLOb330904_3 330904 48.11 14.84 –33.27
WLOb330938_4 330938 2.11 1.63 –0.48
Delaware wells
WLObFc4211_4 Fc4211 −7.67 17.53 25.20
WLObGd3304_7 Gd3304 −0.67 –12.90 –12.23
Table 16.    Simulated and observed water-level observations and residuals for wells open to the Wenonah-Mount Laurel aquifer, New Jersey, and the Mount Laurel aquifer, Delaware, 2013.
Map showing water level contours in light blue and brown, as well as aquifer and updip
                           limits in burgundy and dark blue within the model grid outline in purple.
Figure 43.

Map showing simulated water-level contours and the 2013 potentiometric surface for, and residuals for selected wells in, the Wenonah-Mount Laurel aquifer, New Jersey and simulated water-level contours for, and residuals for selected wells in, the Mount Laurel aquifer, Delaware.

Red dots showing water level, light and dark blue lines showing lines of correspondence.
Figure 44.

Graph showing simulated and observed water levels for selected wells open to the Wenonah-Mount Laurel aquifer, New Jersey, 2013.

The hydrographs of selected wells screened in the Wenonah-Mount Laurel aquifer are shown in figure 45. The hydrograph for well 070478 in Camden County indicates simulated water levels have a similar trend as observed water levels but a smaller amplitude and are lower from 1980 to 1992, but higher from 1993 to 2013, with an average difference of 8 ft (fig. 45A). These trends are also observed in the hydrograph shown in Voronin (2004) for this well during 1980–98. The hydrographs of simulated water levels for wells 250353 (fig. 45B) and 250637 (fig. 45C) in Monmouth County have the same trend and amplitude as observed water levels in those wells but with an average residual of 0.3 and 10 ft lower, respectively, than the observed water levels. Water levels in well 250353 simulated by Voronin (2004) show a smaller amplitude than the observed water levels, and the simulated water levels in well 250637 are about 10 ft lower than those observed during 1980–98. The hydrograph for well 290140 in Ocean County indicates that the simulated water levels have the same trend and amplitude as the observed water levels, and are, on average, within 3 ft of the observed water levels (fig. 45D). Simulated water levels in well 330020 in Salem County have the same trend and amplitude as observed water levels but are an average of 12 ft lower than them (fig. 45E).

Scatter plot showing simulated versus measured water levels in wells.
Figure 45.

Hydrographs of simulated and observed water levels for wells A, 070478, B, 250353, C, 250637, D, 290140, and E, 330020, open to the Wenonah-Mount Laurel aquifer, New Jersey, 1980–2013.

Calibration of Water Levels in the Englishtown Aquifer System

Simulated water-level contours in the Englishtown aquifer system (model layer 15, Monmouth, Ocean, Burlington, Camden, Gloucester, and Salem Counties) generally match the 2013 potentiometric surface (Cauller and Gordon, 2021) but the simulated water levels are lower in parts of northern Camden and Burlington counties and near a pumping center in updip Monmouth County (fig. 46). The simulated contours and 2013 potentiometric surface generally have the same minima in the cones of depression in coastal southern Monmouth County and central northern Ocean County, but the cones differ in shape, partly as a result of simulated residuals and (or) of interpretation of the limited water-level data in areas used to develop the potentiometric contours. The differences between the simulated and observed water levels for wells measured in 2013 are shown in table 17 and figures 46 and 47. The average water-level residual (simulated minus observed) for 39 wells in New Jersey is −5 ft and the RMSE is 17 ft. The minimum; first, second, and third quartile; and maximum residuals are −41, −13, −1, 6, and 28 ft, respectively. The simulated water levels are within 10 ft of observed water levels for 23 of the 39 wells (59 percent) and within 15 ft of observed water levels for 29 of the 39 wells (74 percent) in New Jersey.

Table 17.    

Simulated and observed water-level observations and residuals for wells open to the Englishtown aquifer system, New Jersey, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (feet) Simulated water level (feet) Residual (feet)
WLOb050197_7 050197 24.36 16.38 –7.98
WLOb050259_7 050259 16.46 15.67 –0.79
WLOb050375_7 050375 21.41 20.61 –0.80
WLOb050437_7 050437 58.71 29.25 –29.46
WLOb050754_7 050754 27.98 22.14 –5.84
WLOb051390_4 051390 5.79 15.68 9.89
WLOb070166_7 070166 16.63 –5.17 –21.80
WLOb070672_6 070672 19.68 –15.25 –34.93
WLOb150676_6 150676 28.70 20.26 –8.44
WLOb230104_6 230104 70.93 81.88 –10.95
WLOb250016_7 250016 –52.57 –85.96 –33.39
WLOb250030_7 250030 –65.22 –81.57 –16.35
WLOb250046_6 250046 41.39 40.49 –0.90
WLOb250080_7 250080 63.76 23.13 –40.63
WLOb250096_7 250096 59.02 32.13 –26.89
WLOb250162_7 250162 –47.36 –41.91 5.45
WLOb250250_7 250250 87.67 72.98 –14.69
WLOb250408_6 250408 102.07 90.41 –11.66
WLOb250429_7 250429 –48.47 –49.54 –1.07
WLOb250441_7 250441 –48.70 –51.57 –2.87
WLOb250638_6 250638 –3.19 –9.28 –6.09
WLOb250697_6 250697 12.64 18.88 6.24
WLOb250704_6 250704 89.15 50.44 –38.71
WLOb250710_6 250710 –65.61 –58.17 7.44
WLOb250715_5 250715 2.71 7.70 4.99
WLOb250733_4 250733 36.71 28.40 –8.31
WLOb290005_7 290005 –70.05 –74.95 –4.90
WLOb290138_7 290138 60.03 45.28 –14.75
WLOb290236_7 290236 –16.64 –7.32 9.32
WLOb290433_7 290433 –70.86 –62.34 8.52
WLOb290441_7 290441 –54.46 –53.23 1.23
WLOb290450_7 290450 –65.93 –43.49 22.44
WLOb290452_6 290452 –58.70 –47.40 11.30
WLOb290503_7 290503 –67.18 –64.53 –2.65
WLOb290530_6 290530 –75.50 –76.49 –0.99
WLOb290534_7 290534 –47.42 –1970 27.72
WLOb290938_6 290938 15.12 27.40 12.28
WLOb330168_4 330168 16.30 15.89 –0.41
WLOb330581_5 330581 12.03 18.67 6.64
Table 17.    Simulated and observed water-level observations and residuals for wells open to the Englishtown aquifer system, New Jersey, 2013.
Map showing water level contours in light blue and brown, as well as aquifer and updip
                           limits in burgundy and dark blue within the model grid outline in purple.
Figure 46.

Map showing simulated water-level contours and 2013 potentiometric surface for, and residuals for selected wells in, the Englishtown aquifer system, New Jersey and simulated water-level contours for the Matawan aquifer, Delaware.

Red dots showing water level, light and dark blue lines showing lines of correspondence.
Figure 47.

Graph showing simulated and observed water-levels for selected wells open to the Englishtown aquifer system, New Jersey, 2013.

The hydrographs of selected wells screened in the Englishtown aquifer system are shown in figure 48. The hydrograph for well 050259 in Burlington County shows that simulated water levels have the same trend and amplitude as observed water levels but are an average of 1 ft lower than observed water levels (fig. 48A). In comparison, the hydrograph in Voronin (2004) for well 050259 shows simulated water levels that are 5 to 10 ft higher than the observed water levels during 1980–98. The simulated water levels in the hydrographs for wells 250429 and 250638 in Monmouth County show the same trend of the recovering observed water levels in the 1990s and the amplitudes of the simulated water levels are similar to those of the observed water levels but are an average of –4 and –12 ft, respectively, lower than them (figs. 48B and 48C). In comparison, the hydrographs in Voronin (2004) for these two wells show simulated water levels that follow a similar pattern as those simulated by the current updated model but a smaller amplitude than observed water levels during 1980–98. Simulated water levels in well 290138 in Ocean County have the same trend and amplitude as observed water levels but are an average of 18 ft lower than them (fig. 48D). The hydrograph for this well in Voronin (2004) shows simulated water levels that are 10 to 15 ft higher than observed water levels during 1980–98.

Scatter plot showing simulated versus measured water levels in wells.
Figure 48.

Hydrographs of simulated and observed water levels in wells A, 050259, B, 250429, C, 250638, and D, 290138, open to the Englishtown aquifer system, New Jersey, 1980–2013.

Calibration of Water Levels in the Upper Potomac Raritan Magothy Aquifer

Simulated water-level contours in the upper Potomac-Raritan-Magothy aquifer (model layer 17, Middlesex, Mercer, Monmouth, Ocean, Burlington, Camden, Gloucester, and Salem Counties) closely match the 2013 potentiometric surface (Cauller and Gordon, 2021) (fig. 49). The simulated contours and 2013 potentiometric surface are similar in shape and the contours correspond well in the areas of the major cones of depression in central Camden and Ocean Counties. The differences (residuals) between the simulated and observed water levels for wells measured in New Jersey in 2013 are shown in table 18 and in figures 49 and 50. The average water-level residual for 88 wells in New Jersey is 9 ft and the RMSE is 14 ft. The minimum; first, second, and third quartile; and maximum residuals for New Jersey wells are −17, 2, 24, 15, and 35 ft, respectively. The simulated water levels are within 10 ft of observed water levels for 46 of the 88 wells (52 percent) and within 15 ft of observed water levels for 64 of the 88 wells (73 percent) in New Jersey. The simulated water levels for the Magothy aquifer in Delaware and the average water-level residual for four wells measured in Delaware in 2013 (table 18) also are shown in figure 49.

Table 18.    

Simulated and observed water-level observations and residuals for wells open to the upper Potomac-Raritan-Magothy aquifer, New Jersey, and the Magothy aquifer, Delaware, 2013.

[USGS, U.S. Geological Survey; residual, simulated water level minus observed water level]

Model well number USGS well number Observed water level (feet) Simulated water level (feet) Residual (feet)
New Jersey wells
WLOb050116_7 50116 −5.52 16.75 22.27
WLOb050167_7 50167 −56.3 −37.20 19.1
WLOb050207_6 50207 −24.25 11.1 35.35
WLOb050212_7 50212 −16.55 16.13 32.68
WLOb050218_7 50218 1.22 5.03 3.81
WLOb050229_7 50229