Abstract

Sea level rise (SLR) can negatively affect the hydrology of coastal watersheds. However, the relevant information is incomplete and insufficient in existing literature. The objective of this study is to present a modeling approach to predict long-term effects of SLR on changes of flood peak, flood stage, and groundwater table with an assumption that the historical climate would reoccur in the future. The study was conducted for a typical coastal watershed in southeast USA. The results indicate that sea level had been rising at a rate of 4.21 mm yr−1 from 1948 to 1982 but at a faster rate of 5.16 mm yr−1 from 1983 to 2013. At such SLR rates and by 2113, the groundwater table beneath the eastern part of the watershed would be raised by 0.10 to 0.29 m, while the annual mean peak discharge and flood stage at the watershed outlet would be increased by 13.84 m3 s−1 (from 3.63 to 17.47 m3 s−1) and 0.92 m (from zero to 0.92 m), respectively. The other parts of the watershed would be relatively less affected by SLR. For coastal watersheds, SLR will probably raise the groundwater table, and increase the magnitude and occurrence of peak discharge and flood stage.

INTRODUCTION

Coastal watersheds that are hydraulically connected with oceans are expected to experience major shifts in hydrologic conditions (Kivett 2014; Goodwin 2015) as a result of changing climate-induced sea level rise (SLR) (Boon et al. 2010; USDOC 2014) as well as more intense and more frequent storms (Yin et al. 2009; Rice & Hirsch 2012; USDOC 2014). Such shifts have been predicted to make coastal regions (e.g., the Chesapeake Bay stretching from northeast to central east USA) more vulnerable to flooding (Najjar et al. 2010; USEPA 2014; Walsh et al. 2014) due to back water effect, storm surge, and/or overland runoff. In response, concerted efforts (Goodwin 2015) have been made by governments, organizations, institutions, and private sectors to develop adaptation measures to sustain the resiliency of coastal regions by mitigating impacts of climate change. However, because quantitative information is scarce and/or incomplete in existing literature regarding how SLR would affect coastal hydrology, those efforts are facing a great challenge in showing quantitative benefits for increasing the resiliency.

Previous SLR-related studies can be classified into three categories. The first category (e.g., Paolisso et al. 2012) focused on policy, socioeconomics, and strategy, with the primary goal of developing action plans to adapt to climate change (Goodwin 2015). In contrast, the second category (e.g., Boon et al. 2010; Sallenger et al. 2012; Eggleston & Pope 2013; Varnell 2014) examined physical mechanisms of SLR and aimed to predict rates of SLR. The third category (e.g., Nuttle & Portnoy 1992; Hong et al. 2010; Bjerklie et al. 2012; Rice & Hirsch 2012; Urquhart et al. 2014) analyzed potential changes in the circulation of water and constituents (i.e., salt, sediment, and nutrients) resulting from SLR, primarily investigating saltwater intrusion into regional aquifers and salinity stratification of groundwater (Masterson et al. 2014).

For the 1.07 km2 Mill Creek watershed (41°56′ N, 70°03′ W), a near-shore coastal watershed located in Massachusetts, Nuttle & Portnoy (1992) analyzed effects of rising sea level on runoff and groundwater discharge to the coastal ecosystem. The results indicate that a 21 cm rise in sea level would raise the groundwater table by 10 cm, which in turn will result in 20% decrease in groundwater discharge but 70% increase in surface runoff. Those two authors attributed the runoff increase to the increased soil moisture resulting from the raised groundwater table. Recently, Masterson et al. (2014) developed a three-dimensional groundwater model using SEAWAT (Langevin et al. 2007) and then used the model to analyze impacts of SLR (up to 60 cm) on the groundwater system beneath Assateague Island, a 60 km-long barrier island on the US mid-Atlantic coast. The results indicate that a SLR of 20 cm would lead to a substantial decrease in the depth to groundwater table because of measurable increases in the extent and depth of saltwater intrusion. On the other hand, Garcia & Loáiciga (2015) presented a method for calculating the contributions of SLR and urban growth to flood risk in the 4,452 km2 Tijuana River basin that straddles the California-Baja California section of the US–Mexico border. In the method, a set of empirical regression equations were used to calculate the flood peaks of return periods (e.g., 100 year) of interest, which in turn were used as inputs into a steady-state HEC-RAS (USACE 2010) model to calculate the responding flood stages of this same return period. Hereinafter, the terms return period and frequency are used interchangeably because the latter is simply the reciprocal of the former. The downstream boundary of the HEC-RAS model was specified as either the current (1.70 m) or a predicted future (1.95, 2.20, or 2.70 m) sea level, depending on the scenario of interest. The results indicate that urbanization will likely play a principal role in increasing the expected annual flood damage in the study area for the range of sea-level rise considered.

However, very few studies (e.g., Passeri et al. 2015; Young et al. 2016) have explored how SLR would affect long-term changes of peak discharge and flood stage along streams within a coastal watershed as well as fluctuations of groundwater table beneath the watershed. As with the others in existing literature, the studies cited above did not consider effects of tidal fluctuations and interdependences among hydroclimate, surface runoff, groundwater flow, and rising sea level. The objective of this study is to present a modeling approach to predict long-term effects of SLR on coastal hydrology in terms of changes of flood peak, flood stage, and groundwater table with an assumption that the historical climate would reoccur in the future. The study was conducted in the Lynnhaven River watershed using Storm Water Management Model (SWMM). This watershed, located on the east coast of the United States, is hydraulically connected with the Atlantic Ocean and incurs frequent coastal flooding (CVB & VDEQ undated; USACE 2013).

SWMM COMPONENTS RELEVANT TO THIS STUDY

SWMM (Rossman 2010) is a dynamic rainfall-runoff simulation model used either for single event or long-term (continuous) simulation of runoff quantity and quality from primarily urban areas. The runoff component of SWMM operates on a collection of subcatchment areas on which rain falls and runoff is generated. Average monthly wind speeds are used when computing snowmelt rates under rainfall conditions. Melt rates increase with increasing wind speed. The routing portion of SWMM transports this runoff through a conveyance system of pipes, channels, storage/treatment devices, pumps, and regulators. SWMM tracks the quantity and quality of runoff generated within each subcatchment, and the flow rate, flow depth, and quality of water in each pipe and channel during a simulation period comprising multiple time steps. While SWMM offers several methods for modeling each of the components, this study modeled the three major components of runoff, groundwater flow, and channel flow using the methods briefly described in the following context. The detailed descriptions can be found in Rossman (2010).

Method for runoff component

The runoff component treats each subcatchment area as a basic computational unit. For the computational time step of tk to tk+1 (k = 1, 2, …, N-1) (in units of s), where t1 is the start time and tN is the end time of a simulation, the water balance of a subcatchment is calculated as:  
formula
(1)
where Y (m) is the average depth of overland flow; L (m) is the length of overland flow path; B (m) is the width of overland flow path; i (m s−1) is the rainfall intensity; f (m s−1) is the infiltration rate for pervious surface; e (m s−1) is the evaporation rate; and Q (m3 s−1) is the overland flow (i.e., runoff) rate.
By assuming a wide-shallow channel, Q is computed using Manning's formula (Roberson et al. 1998) as:  
formula
(2)
where Yd (m) is the surface depression storage and S (−) is the average surface slope in fraction.
f is computed using Green-Ampt (G-A) model (Viessman & Lewis 2003; Hilpert & Glantz 2013) as:  
formula
(3)
where Ksat (m s−1) is the saturated hydraulic conductivity; ψf (m) is the soil capillary suction head at wetting front; θsat (−) is the saturated soil moisture; θa (−) is the initial soil moisture; and F (m) is the cumulative infiltration.
F is computed as:  
formula
(4)
where min{ } is the minimum function.

e is estimated using the Hargreaves method (Hargreaves & Samani 1985; Wang et al. 2006). At the end of each computational time step, the soil moisture of the upper zone (Figure 1) is redistributed using a rectangular profile model presented by Morel-Seytoux (1984, 1985).

Figure 1

Illustration of the groundwater flow component (after Rossman 2010).

Figure 1

Illustration of the groundwater flow component (after Rossman 2010).

Method for groundwater flow component

SWMM models lateral groundwater flow (QL) between receiving node (i.e., junction) and lower zone (i.e., shallow unconfined aquifer) as well as rate of percolation (QD) to deep semi-confined aquifer. Herein, it is assumed that the horizontal flow in the shallow aquifer beneath a subcatchment is fully received by the junction just downstream of the subcatchment. Also, it is assumed that the process of saltwater intrusion is much slower than that of recharge from streams into the shallow aquifer (Barlow 2003), and thus SLR simply affects the groundwater table by raising stream water surface elevation.

QL (in units of m3 s−1) is computed as:  
formula
(5)
where HGW (m) is the height of saturated zone above the bottom of aquifer; HSW (m) is the height of surface water at receiving node above aquifer bottom; HCB (m) is the height of channel bottom above aquifer bottom; and A1, A2, A3, B1, and B2 are five coefficients that can be determined by model calibration.
QD (m3 s−1) is computed as:  
formula
(6)
where LGLR (m3 s−1) is the lower groundwater loss rate parameter assigned to the subcatchment's aquifer and HGS (m) is the distance from the ground surface to the aquifer bottom.

Method for flow routing component

This study uses the dynamic wave method to route flows through the stream network. This method numerically solves the one-dimensional Saint-Venant dynamic equations (Roberson et al. 1998) expressed as:  
formula
(7)
 
formula
(8)
where y is the water depth in the channel; v is the mean velocity in the channel; Sf is the energy line gradient; S0 is the channel bed slope; and g is the gravitational acceleration.
Sf is computed using the Darcy-Weisbach equation (Finnemore & Franzini 2002) expressed as:  
formula
(9)
where f is the Darcy friction factor and Rh is the hydraulic radius.

MATERIALS AND METHODS

Lynnhaven River watershed

The Lynnhaven River drains the northern part of City of Virginia Beach, Virginia, and discharges into the Chesapeake Bay at the Lynnhaven Inlet (36°54′27.90″ N, 76°05′30.90″ W) (Figure 2). The 170 km2 watershed is a treasured and pivotal part of the community in Virginia Beach (USACE 2013). It is home to thousands of boaters and residents, and has become a daily part of life for them. The watershed consists of the subareas drained by three major tributaries, namely the Broad-Linkhorn Bay, the Eastern Branch, and the Western Branch, and has nearly 242 km of shoreline. The Lynnhaven River is oligohaline and subject to the action of tides (CVB & VDEQ undated). Depending on wind and tidal patterns, the upstream portions of the river system either flow north to the Chesapeake Bay or south to the North Carolina sounds (Schoenbaum 1988). The majority of the waters outside the bays are shallow with maintained channel depths of 1.8 to 3.0 m.

Figure 2

Lynnhaven River watershed and its subcatchments, with tide station, climate stations, and phreatic well, where observed data were used in this study.

Figure 2

Lynnhaven River watershed and its subcatchments, with tide station, climate stations, and phreatic well, where observed data were used in this study.

The topographic elevations of the watershed vary from 0.0 to 26.0 m above mean sea level, with a mean of 2.4 m and a standard deviation of 6.3 m. The soils are dominated by loam (Todd & Mays 2005): 17% clay (particle size of <0.004 mm), 37% silt (particle size of 0.004 to 0.062 mm), and 46% fine sand (particle size of 0.125 to 0.25 mm). The soils, with a saturated hydraulic conductivity of 3.4 to 77.5 mm h−1 (mean 19.9 mm h−1), are moderately permeable (Todd & Mays 2005). The land is covered by 9.8% water, 19.2% impervious surface (paved roads, buildings, and parking lots), and 71.0% pervious area (forests and grasses). The watershed has a temperate and humid subtropical climate, with long warm summers and relatively short mild winters. The average summer temperature is 25 °C, with a maximum daily average of 29 °C, while the average winter temperature is 6 °C, with a minimum daily average of 0.5 °C. The watershed receives a mean annual precipitation of 1,143 mm. During fall and spring, nor'easters, which are macro-scale storms along the upper East Coast of the United States and Atlantic Canada, may impact the watershed, causing localized flooding (MESO 2002).

The aquifers of the watershed include Columbia, Yorktown, and Yorktown-Eastover (USGS 2002). The deep Yorktown-Eastover confined aquifer (upper confining bed at 28.0 m or deeper below ground surface) is overlain by the Yorktown confined aquifer (upper confining bed at 8.5 to 28.0 m below ground surface), which in turn is overlain by the unconfined Columbia aquifer (bottom confining bed at 8.5 to 24.5 m below ground surface). Because of the high permeability of the soils, the Columbia aquifer can be periodically recharged by infiltrated precipitation and high tidal water, whereas its water can be efficiently evaporated and/or taken up by vegetation during sunny days as well as discharged into the bays during low tidal periods. As a result, the groundwater table in the aquifer fluctuates with time, and exhibits both inter- and intra-annual variations. In contrast, the groundwater level in the Yorktown and Yorktown-Eastover confined aquifers is not influenced much by such localized hydrologic and marine processes. Thus, this study just used data on groundwater table in the Columbia aquifer.

Table 1

Trend statistics of the mean annual and monthly tide levelsa

Time scale for computing mean tide level Statistics Z* p-value Sen's slope SSjb (mm yr−1
Annual 8.962 0.000 4.65 
January 6.075 0.000 4.24 
February 5.712 0.000 4.21 
March 6.245 0.000 4.95 
April 6.516 0.000 5.14 
May 7.716 0.000 4.71 
June 8.101 0.000 4.59 
July 8.639 0.000 4.85 
August 8.124 0.000 4.75 
September 7.694 0.000 4.73 
October 6.448 0.000 4.31 
November 6.573 0.000 4.35 
December 7.445 0.000 5.16 
Time scale for computing mean tide level Statistics Z* p-value Sen's slope SSjb (mm yr−1
Annual 8.962 0.000 4.65 
January 6.075 0.000 4.24 
February 5.712 0.000 4.21 
March 6.245 0.000 4.95 
April 6.516 0.000 5.14 
May 7.716 0.000 4.71 
June 8.101 0.000 4.59 
July 8.639 0.000 4.85 
August 8.124 0.000 4.75 
September 7.694 0.000 4.73 
October 6.448 0.000 4.31 
November 6.573 0.000 4.35 
December 7.445 0.000 5.16 

aDetails of how to compute the statistics can be found in Wang et al. (2014).

bEqual to the average change rate of mean tide level for the responding time scale (see Equation (11)).

Data and preprocessing

This study used six national spatial data layers, namely 10-m National Elevation Dataset (NED), National Hydrography Dataset (NHD), Soil Survey Geographic database (SSURGO), National Land Cover Database 2011 (NLCD2011), Percent Developed Imperviousness (PDI), and National Wetland Dataset (NWD). The NED and NHD data were downloaded from the U.S. Geological Survey (USGS) website (https://lta.cr.usgs.gov/get_data); the NCLD2011 and PDI data from the Multi-Resolution Land Characteristics Consortium (MRLC) website (http://www.mrlc.gov/nlcd11_data.php); the NWD data from the National Fish and Wildlife Service (NFWS) website (http://www.fws.gov/wetlands/Data/Data-Download.html); and the SSURGO data from the U.S. Department of Agriculture (USDA) Natural Resources Service (NRCS) website (http://www.nrcs.usda.gov/wps/portal/nrcs/main/soils/survey/geo). These six data layers were projected in ArcMap™ 10 to have a common coordinate system of Universal Transverse Mercator (UTM) Zone 18 (Karney 2011).

The 10-m NED was developed by merging the highest-resolution, best-quality elevation data available across the USA into a seamless raster format (USGS 2014a), while the NHD is a comprehensive set of digital spatial data that contains information about surface water features such as lakes, ponds, streams, rivers, springs, and wells (USGS 2014b). The NLCD2011 is the most recent national land cover product (MRLC 2014). It provides, for the first time, the capability to assess wall-to-wall, spatially explicit, national land cover changes and trends across the United States from 2001 to 2011. As with the two previous NLCD land cover products, NLCD2011 keeps the same 16-class land cover classification scheme that has been applied consistently across the United States at a spatial resolution of 30 m. NLCD 2011 is based primarily on a decision-tree classification of circa 2011 Landsat satellite data. The PDI data were extracted from the NLCD2011 using the algorithm presented by Xian et al. (2011).

The NWD data are graphic representations of the type, size and location of the wetlands and deepwater habitats in the United States (USFWS 2014). These maps have been prepared from the analysis of high altitude imagery in conjunction with collateral data sources and field work. The maps represent reconnaissance level information on the location, type, size of wetlands habitats such that they are accurate at the nominal scale of the 1:24,000 base map for the contiguous United States and the 1:63,360 base map for Alaska. A margin of error is inherent in the use of imagery; thus, detailed, on-the-ground inspection of any particular site may result in revision of the wetland boundaries or classification established through image analysis. The data were prepared by using the Cowardin et al. (1979) definition of wetland. This definition is the national standard for wetland mapping, monitoring, and data reporting, as detailed in FGDC (1996).

The SSURGO was designed by the USDA-NRCS (2005) to provide farm-level spatial resolution for farm and ranch, landowner/user, township, county, or parish natural resources planning and management. Data for SSURGO are collected and archived in the USGS 7.5-minute topographic quadrangle units and distributed as a complete coverage for a soil survey area that usually consists of 10 or more quadrangle units. The adjoining 7.5-minute units are matched within the survey areas. Wang & Melesse (2006) found that SSURGO does not assign a classification for soils perennially covered by water and suggested it would be reasonable to assume that these areas have the same soil classifications as adjacent areas.

In addition, data on hourly rainfall and daily minimum and maximum temperatures, observed at Norfolk International Airport (COOP ID#: 446139; 36°54′12″ N, 76°11′32″ W) by the National Oceanic and Atmospheric Administration (NOAA) National Climate Data Center (NCDC), were downloaded from the NCDC website (http://www.ncdc.noaa.gov/cdo-web/results) for a record period of 5 August 1948 to 30 June 2013. The records for daily temperatures are complete, but the record for hourly rainfall has missing values from January 1952 to July 1953. This study filled the missing values using responding observed daily data by presuming a NRCS Type III distribution (Viessman & Lewis 2003). For a given day, the missing values were estimated as the multiplications of the daily value and the NRCS 24-h distribution coefficients. Also, data on daily average wind speed, observed at Norfolk NAS (ID#: USW00013750; 36°49′01″ N, 76°02′00″ W) by NOAA-NCDC, were downloaded from the same NCDC website for a record period of 1 September 1986 to 31 August 2014 and used to compute annual average monthly wind speeds. For a given month, the observed daily values were arithmetically averaged across the record period to get the annual average wind speed of this month.

Further, data on hourly sea level, observed at Sewells Point (ID#: 8638610; 36°56′48″ N, 76°19′48″ W) by NOAA, were downloaded from the NOAA's Tides & Currents website (http://tidesandcurrents.noaa.gov/waterlevels.html?id=8638610) for a record period of 5 August 1948 to 30 September 2014. Because this observational point is very close to the Lynnhaven Inlet, the time series data were used as the lower boundary condition of the SWMM model in this study. Moreover, data on daily average phreatic level in the Columbia unconfined aquifer, observed in monitoring well 61D6Sow124 (ID#: 365327076080501; 36°53′27″ N, 76°08′05″ W), were downloaded from the USGS website (http://nwis.waterdata.usgs.gov/nwis/gwlevels?site_no=365327076080501&agency_cd=USGS&format=html) for a record period of April 1978 to November 2013. The data were used to calibrate the groundwater component of the SWMM model.

The Lynnhaven River watershed had no long-term continuous streamflow gauges. Instead, flow measurements were done by City of Virginia Beach at two locations for three storms that occurred in 2011. The first location (36 °52′09″ N, 76 °01′59″ W) controls a drainage area of 77.5 ha and measured flows resulting from the storms of December 7 and 16, while the second location (36 °52′17″ N, 76 °02′23″ W) controls a drainage area of 29.6 ha and measured flows resulting from the storm of December 21. The measured flows were independently transferred to the outlet of subcatchment 8 (Figure 2) that includes the drainage areas controlled by these two observational locations using Equation (10), and the transferred flows were considered to be the responding runoffs from this entire subcatchment. This transfer approach has been commonly used to estimate flows at a location of interest based on observations at a measurement location (Wang et al. 2010). The measurement location usually controls a smaller drainage area than the location of interest. The transferred flows from the storm event of December 7 were used to calibrate the overland runoff and infiltration components of the SWMM model, whereas the transferred flows from the other two storm events were used to validate these two components.  
formula
(10)
where Qsub,t is the transferred flow or runoff at the outlet of the inclusive subcatchment at time t; QL,t is the measured/observed flow at location L at time t; DAsub is the drainage area of the subcatchment; DAL is the drainage area upstream of location L.

Trend analysis of sea level

For each of the calendar years from 1949 to 2013, the observed hourly sea levels at Sewells Point were arithmetically averaged to get the mean annual sea level of this year. Similarly, for a given month of a year, the observed hourly sea levels were arithmetically averaged to get the mean monthly sea level of this month of this year. As a result, 13 time series, one for mean annual and 12 for mean monthly sea levels, were derived. Each time series consists of 65 values, each of which responds to one of the 65 calendar years.

For the annual time series, the distribution-free cumulative sum (CUSUM) technique (McGilchrist & Woodyer 1975; Wang et al. 2014) was used to detect whether, and at which year, a significant abrupt change (i.e., upward-to-downward or downward-to-upward) occurred at a significance level of α = 0.05. In addition, the time series were plotted to visually examine the temporal variations of mean sea levels at annual and monthly time scales, respectively. Further, for each of the time series, the modified Mann-Kendall trend test technique (Mann 1945; Hirsch et al. 1982; Hamed & Rao 1998) was applied to determine whether a statistically significant (α = 0.05) temporal trend (either downward or upward) existed, while the sequential Mann-Kendall method (Taubenheim 1989; Sneyers 1990) was used to determine the year when a trend started, and the Sen's slope (Q) (Sen 1968) was computed and used to be the change rate of sea level. If the modified Mann-Kendall test statistic of Z* is positive and greater than the standard normal 1-α = 95th percentile, a statistically significant upward trend was indicated. In contrast, if Z* is negative and smaller than the standard normal 95th percentile, a statistically significant downward trend was detected. Otherwise, it was judged that no significant trend existed. The details of the test technique and method can be found in Wang et al. (2014).

Model set up

The model set up (Figure 2) consists of: (1) definition of subcatchments, conduits, and aquifers; (2) parameterization of the methods for the components described above; and (3) calibration and validation using historical data on stream flow and groundwater table.

The NED was used in the Hydrology Extension of ArcMap™ 10 to delineate the watershed and its subcatchments as well as the drainage pathways. The delineation process was iterated by manually adjusting the flow accumulation threshold (Wang et al. 2008, 2010) above which a pathway would be formed until the delineated drainage pathways were visually compatible with the stream network presented by the NHD. In addition, for each of the delineated subcatchments, its characteristics (Gironάs et al. 2009; Rossman 2010), namely drainage area, width of overland flow path (i.e., width of the overland transverse section that is wetted by flow), and average surface slope, were estimated in terms of the elevation values presented by the NED. Further, for each of the delineated drainage pathways (i.e., conduits), its length was automatically computed by the Hydrology Extension, while a representative channel cross section was extracted from the NED using the River Bathymetry Toolkit (RBT), an ArcMap™ Add-In developed by U.S. Forest Service (USFS) (McKean et al. 2009). The RBT uses topographic information to define pairs of channel bottom elevation and horizontal distance over the cross section. Herein, the representative cross section was positioned at the halfway point of the drainage pathway from its mouth. Also, the inlet of a conduit was reasonably assumed to be at the outlet of the subcatchment drained by the conduit. Subsequently, the subcatchment map was overlaid with the PDI to compute the percentage of imperviousness for each of the delineated subcatchments. For a subcatchment that includes wetlands, the length and cross section of its drainage pathway, which were previously determined from the NED, were adjusted up to reflect the characteristics of the wetlands by one hydrologic equivalent wetland (HEW) (Wang et al. 2008). The HEWs were defined based on the real wetlands presented by the NWD.

Further, the subcatchment map was overlaid with the SSURGO to estimate the aquifer properties, namely porosity (n), wilting point (θr), field capacity (θfc), saturated hydraulic conductivity (Ksat), average slope of soil tension (ψ) versus soil moisture content (θ) curve, and average slope of log(Ksat) versus θ curve. With this regard, for a given subcatchment, the values of percentage of sand within the soil map units that are included in this subcatchment were area-weight averaged to get the percentage of sand within this subcatchment, while the values of percentage of clay within the soil map units that are included in this subcatchment were area-weight averaged to get the percentage of clay within this subcatchment. Similarly, the soil organic matter content and salinity of this subcatchment were computed as the area-weighted averages of the corresponding values of the soil map units presented by the SSURGO. In turn, the computed percentage of sand, percentage of clay, organic matter content, and salinity were input into the Saxton-Rawls model (Saxton et al. 1986; Saxton & Rawls 2006) to estimate the aforementioned seven properties of the aquifer beneath this subcatchment. Moreover, the bottom elevation of the aquifer was determined from the geological profiles presented by Smith & Harlow (2002) and Heywood & Pope (2009).

The historical (August 5, 1948 to December 29, 2013) data on hourly rainfall, daily minimum and maximum temperatures, mean monthly wind speed, and hourly sea level, were used to drive the SWMM model. Within the simulation period, the first five years (August 5, 1948 to December 31, 1953) were taken as the ‘warm-up’ period to stabilize the initial values of the model states (e.g., the initial soil moistures of the unsaturated upper zones of the aquifers), whereas the remaining years (January 1, 1954 to December 29, 2013) were used for model evaluation (i.e., calibration and validation). Given that the groundwater data had a record period of April 1978 to November 2013, the data from 1978 to 1989 were used for model calibration while the data from 1990 to 2013 were used for model validation. As stated above, one of the three storm events was randomly selected for model calibration, whereas the other two storm events were used for model validation. The calibration was done by manually adjusting selected parameters (discussed below in the section ‘Calibrated model and parameter sensitivity’ in the following context) until predicted values of runoff and groundwater table were comparable with the responding observations. The model performance was visually judged in terms of plots showing predicted versus observed groundwater table or transferred runoff. Moreover, a sensitivity analysis was conducted by manually varying some of the calibration parameters within their responding ranges. When a parameter was varied, the other parameters were fixed to their calibrated values. That is, the parameters were varied one at a time and thus the sensitivity analysis did not take into account any interactions among them. Nevertheless, because of the limited available data, this study used SWMM as a ‘screening’ model to predict the ‘relative changes’ rather than ‘absolute results’ from climate change and SLR. Herein, the basic rationale is that physically based models such as SWMM can be assumed to represent the modeled watershed. Although the model parameters cannot be sufficiently adjusted because of insufficient available data, the model may give us the ‘trustable’ changes for two sets of inputs of climate and sea level. Many government agencies in the USA and throughout the world use SWMM as a ‘screening’ model without calibration/validation because of lack of data (Rossman 2010).

Simulation approach and formulation of scenarios

The lower boundary of the calibrated model was the observed hourly sea levels at Sewells Point. When the model was used to assess effects of SLR on hydrology of the study watershed, the hourly sea levels at the same station for any of the six scenarios were used to replace the observed hourly sea levels as lower boundary. The six scenarios were formulated to reflect possibly expected mean sea levels in 2023, 2033, 2043, 2063, 2093, and 2113.

For a given scenario, the sea levels were predicted as:  
formula
(11)
where s is the reference year for a scenario of interest; i (=1948, 1949, 1950, …, 2013) is the observation year; SLj,k is the sea level at time k of month j; and SSj is the Sen's slope (i.e., rising rate) of sea level in month j.

RESULTS

Trends and increasing rates of SLR

At the annual time scale, the mean sea (i.e., tide) level at Sewells Point experienced an abrupt increase around 1983, indicated by the minimum value (=− 27) of CUSUM statistics Vk outside the 95% confidence limits (Figure 3(a)). The abrupt increase can also be noticed by examining the plot showing the mean tide level versus year (Figure 3(b)): the mean tide level in 1983 was 82 mm higher than that in 1982, but it was 23 mm lower than that in 1984. In terms of the Sen's slopes, before 1983, the mean tide level had been rising by an average rate of 3.9 mm yr−1, while after 1983 the mean tide level was rising by a faster rate of 5.8 mm yr−1. Averaged across the entire record period of 1950 to 2013, the mean tide level was rising by 4.65 mm yr−1, as revealed by the significant Sen's slope at α = 0.05 (Table 1). The abrupt increase might be the inception of the accelerated rising rate after 1983.

Figure 3

Plots showing (a) CUSUM (the distribution-free cumulative sum) statistics Vk and (b) observed mean annual tide level.

Figure 3

Plots showing (a) CUSUM (the distribution-free cumulative sum) statistics Vk and (b) observed mean annual tide level.

At the monthly time scales, the mean tide levels were also rising at an accelerating rate (Figure 4). The rising rates for five months (January, February, June, October, and November) were slower than the rising rate of mean annual tide level (≤4.21 mm yr−1), whereas the rising rates for the other seven months were faster than the rising rate of mean annual tide level (≥4.71 mm yr−1) (Table 1). For all months, the rising rates were significant (p-value = 0.000+) at the significance level of α = 0.05, indicating that rising trends were persistent regardless of seasons and time scales.

Figure 4

Plot showing observed mean monthly tide levels.

Figure 4

Plot showing observed mean monthly tide levels.

Simulated scenarios

This study selected six reference years, 2023, 2033, 2043, 2063, 2093, and 2113, to formulate six scenarios (Table 2). For a given month of a reference year, the rise of mean monthly sea level was computed as the multiplication of the rising rate (i.e., Sen's slope) of this month (Table 1) and the number of years from 2013, while the time series of hourly sea level for this reference year or scenario was generated in terms of Equation (11). The scenarios were predicted to have a mean annual sea level of 0.070 to 0.476 m, which is equivalent to a rise of 0.052 to 0.459 m from the mean annual sea level of 0.017 m in 2013. However, as expected, the time series for the scenarios would have a similar temporal variation pattern to that of the historical observations, around the responding monthly mean sea levels (Figure 5).

Table 2

Summary of the scenarios simulated in this study

Scenario Reference year Mean annual tide level (m) Rise since 2013a (m) 
2023 0.070 0.052 
II 2033 0.102 0.085 
III 2043 0.148 0.131 
IV 2063 0.242 0.225 
2093 0.383 0.365 
VI 2113 0.476 0.459 
Scenario Reference year Mean annual tide level (m) Rise since 2013a (m) 
2023 0.070 0.052 
II 2033 0.102 0.085 
III 2043 0.148 0.131 
IV 2063 0.242 0.225 
2093 0.383 0.365 
VI 2113 0.476 0.459 

aThe mean tide level in 2013 is 0.017 m.

Figure 5

Plot showing the mean annual sea levels in history and for the six scenarios analyzed.

Figure 5

Plot showing the mean annual sea levels in history and for the six scenarios analyzed.

Calibrated model and parameter sensitivity

As a result of the delineation, the Lynnhaven River watershed was subdivided into 25 subcatchments (Figure 6) with a drainage area of 0.017 (at the watershed outlet) to 1,441 ha, a width of 19 to 2,074 m, and an average topographic gradient of 0.0 to 3.1% (Table 3). The three lowest subcatchments (i.e., subcatchment 1, 2, and 3) are fully covered by water, while the other subcatchments have a land cover of 0.4 to 50.3% water, 3.7 to 43.6% imperviousness, 0.0 to 67.8% wetland, and 0.0 to 83.9% pervious area. For modeling purposes, for a given subcatchment, the summation of the percentage of water and that of imperviousness was treated as the percentage of impervious area of this subcatchment.

Table 3

Characteristics of the delineated subcatchments of the Lynnhaven River watershed

Subcatchment Area (ha) Width (m) Slope (%) Imperviousness (%) Water (%) Wetland (ha) Downstream subcatchment 
2.395 91.75 0.00 0.00 100.00 0.00 
0.017 18.60 0.00 0.00 100.00 0.00 
12.648 181.69 0.88 0.00 100.00 0.00 
144.896 356.20 2.18 36.49 12.00 7.11 
449.566 665.61 1.52 18.11 11.40 16.65 
327.312 766.11 1.99 10.53 40.23 8.82 
1,441.012 886.40 2.41 5.67 3.68 976.95 
484.017 880.91 1.74 2.36 10.88 32.22 
416.160 883.12 1.60 4.82 39.90 88.83 
10 600.686 723.16 3.12 19.34 5.40 63.63 
11 865.242 1,456.56 1.83 8.16 23.10 19.08 
12 463.857 1,167.92 1.83 14.83 6.26 37.53 11 
13 165.774 697.35 1.51 3.68 50.27 2.25 
14 1,410.555 1,598.15 2.18 13.31 12.32 57.60 
15 646.394 1,054.17 1.60 24.37 12.53 32.67 13 
16 658.377 1,330.71 1.38 22.10 1.93 58.77 14 
17 626.121 1,059.65 2.15 15.26 5.89 61.20 11 
18 403.996 931.01 0.97 33.38 3.42 7.56 19 
19 57.630 298.43 2.22 13.67 21.32 3.15 14 
20 356.706 822.61 1.36 27.28 3.68 21.06 17 
21 1,439.154 2,073.71 1.49 30.38 1.93 93.78 17 
22 296.915 724.35 0.87 39.54 2.64 13.50 19 
23 668.051 1,582.56 0.49 43.57 0.44 18.99 22 
24 1,143.051 1,206.52 1.54 15.09 11.99 134.19 13 
25 1,092.009 1,680.93 0.26 30.58 2.90 20.16 22 
Minimum 0.017 18.60 0.00 0.00 0.44 0.00 
Maximum 1,441.012 2,073.71 3.12 43.57 100.00 976.95 22 
Subcatchment Area (ha) Width (m) Slope (%) Imperviousness (%) Water (%) Wetland (ha) Downstream subcatchment 
2.395 91.75 0.00 0.00 100.00 0.00 
0.017 18.60 0.00 0.00 100.00 0.00 
12.648 181.69 0.88 0.00 100.00 0.00 
144.896 356.20 2.18 36.49 12.00 7.11 
449.566 665.61 1.52 18.11 11.40 16.65 
327.312 766.11 1.99 10.53 40.23 8.82 
1,441.012 886.40 2.41 5.67 3.68 976.95 
484.017 880.91 1.74 2.36 10.88 32.22 
416.160 883.12 1.60 4.82 39.90 88.83 
10 600.686 723.16 3.12 19.34 5.40 63.63 
11 865.242 1,456.56 1.83 8.16 23.10 19.08 
12 463.857 1,167.92 1.83 14.83 6.26 37.53 11 
13 165.774 697.35 1.51 3.68 50.27 2.25 
14 1,410.555 1,598.15 2.18 13.31 12.32 57.60 
15 646.394 1,054.17 1.60 24.37 12.53 32.67 13 
16 658.377 1,330.71 1.38 22.10 1.93 58.77 14 
17 626.121 1,059.65 2.15 15.26 5.89 61.20 11 
18 403.996 931.01 0.97 33.38 3.42 7.56 19 
19 57.630 298.43 2.22 13.67 21.32 3.15 14 
20 356.706 822.61 1.36 27.28 3.68 21.06 17 
21 1,439.154 2,073.71 1.49 30.38 1.93 93.78 17 
22 296.915 724.35 0.87 39.54 2.64 13.50 19 
23 668.051 1,582.56 0.49 43.57 0.44 18.99 22 
24 1,143.051 1,206.52 1.54 15.09 11.99 134.19 13 
25 1,092.009 1,680.93 0.26 30.58 2.90 20.16 22 
Minimum 0.017 18.60 0.00 0.00 0.44 0.00 
Maximum 1,441.012 2,073.71 3.12 43.57 100.00 976.95 22 
Figure 6

Schematic of the SWMM model, with delineated subcatchments and streams of the Lynnhaven River watershed.

Figure 6

Schematic of the SWMM model, with delineated subcatchments and streams of the Lynnhaven River watershed.

In the SWMM model, the values in Table 3 were used to define the characteristics of the subcatchments, whereas Table 4 lists the calibration parameters and their adopted values. Because no multiple-station data were available, the values of the parameters related to overland runoff were assumed to be same as those of subcatchment 8, without being differentiated by subcatchment. The parameters for subcatchment 8 were manually adjusted to closely match the simulated and transferred flows at the outlet of this subcatchment. However, the parameters of the Green-Ampt (i.e., G-A) infiltration model and groundwater flow were varied subcatchment by subcatchment to better represent the spatial heterogeneities of soils and aquifers. These parameters for subcatchment 14 were manually adjusted to closely match the simulated values of groundwater table beneath this subcatchment and the responding observations at the adjacent monitoring well 61D6Sow124 (Figure 2). Herein, it was assumed that the groundwater table gradient between the well location and subcatchment 14 was negligible because of the short (<2.5 km) geographic distance.

Table 4

List of the calibrated model parameters

Parameter Definition Range/Initial value Adopted value 
Overland runoff (Equations (1) and (2)) 
 N-Imperv Manning's n for impervious area 0.012–0.02 0.015 
 N-Perv Manning's n for pervious area 0.015–0.035 0.02 
 Dstore-Imperv (mm) Depth of depression storage on impervious area 1.2–3.5 
 Dstore-Perv (mm) Depth of depression storage on pervious area 2.5–6.5 
 %Zero-Imperv (%) Percentage of impervious area with no depression storage 0–10 
 Subarea Routing Choice of internal routing between pervious and impervious sub-areas Outlet, Impervious, Pervious Pervious 
 Percent Routed (%) Percentage of runoff routed between sub-areas 0–20 
Infiltration (Equations (3) and (4))a 
 Suction Head (mm) Soil capillary suction head Estimatedb 83–362 
 Conductivity (mm hr−1Soil saturated hydraulic conductivity Estimatedc 3–78 
Channel routing (Equations (7)–(9)) 
 Length (m) Conduit length Delineatedd 222–11,494 
 Roughness Manning's n for conduit 0.015–0.035 0.03 
Groundwater (Equations (5) and (6)) 
 A1 Coefficient (m3-B1 s−1 ha−1Groundwater influence multiplier 0 to Ksat/(3.6 × 105)e 6.8 × 10−6 to 4.8 × 10−4 
 B1 Exponent Groundwater influence exponent 0–2 
 A2 Coefficient (m3-B2 s−1 ha−1Tailwater influence multiplier 0 to Ksat/(3.6 × 105)e 6.8 × 10−6 to 4.8 × 10−4 
 B2 Exponent Tailwater influence exponent 0–2 
Parameter Definition Range/Initial value Adopted value 
Overland runoff (Equations (1) and (2)) 
 N-Imperv Manning's n for impervious area 0.012–0.02 0.015 
 N-Perv Manning's n for pervious area 0.015–0.035 0.02 
 Dstore-Imperv (mm) Depth of depression storage on impervious area 1.2–3.5 
 Dstore-Perv (mm) Depth of depression storage on pervious area 2.5–6.5 
 %Zero-Imperv (%) Percentage of impervious area with no depression storage 0–10 
 Subarea Routing Choice of internal routing between pervious and impervious sub-areas Outlet, Impervious, Pervious Pervious 
 Percent Routed (%) Percentage of runoff routed between sub-areas 0–20 
Infiltration (Equations (3) and (4))a 
 Suction Head (mm) Soil capillary suction head Estimatedb 83–362 
 Conductivity (mm hr−1Soil saturated hydraulic conductivity Estimatedc 3–78 
Channel routing (Equations (7)–(9)) 
 Length (m) Conduit length Delineatedd 222–11,494 
 Roughness Manning's n for conduit 0.015–0.035 0.03 
Groundwater (Equations (5) and (6)) 
 A1 Coefficient (m3-B1 s−1 ha−1Groundwater influence multiplier 0 to Ksat/(3.6 × 105)e 6.8 × 10−6 to 4.8 × 10−4 
 B1 Exponent Groundwater influence exponent 0–2 
 A2 Coefficient (m3-B2 s−1 ha−1Tailwater influence multiplier 0 to Ksat/(3.6 × 105)e 6.8 × 10−6 to 4.8 × 10−4 
 B2 Exponent Tailwater influence exponent 0–2 

aG-A Model refers to Green-Ampt Model (Viessman & Lewis 2003).

bFor a given subcatchment, in terms of its dominant soil classification (USDA-SCS 1987), the soil bubbling pressure (ψb) and pore-size index (λ) were determined based on Assouline (2005) and Charbeneau (2006). The suction head was estimated as .

cFor a given subcatchment, it was estimated using the Saxton-Rawls model (Saxton et al. 1986; Saxton & Rawls 2006).

dFor a given subcatchment, it was automatically delineated by the Hydrology Extension of ArcMap™ 10.

eKsat is the soil saturated hydraulic conductivity, while 3.6 × 105 is a factor to convert mm hr−1 to m2 s−1 ha−1.

Visually, the model well captured the rising, primary peak, and recession of the transferred flow hydrograph of the calibration event (Figure 7(a)). The secondary peak of the transferred flow hydrograph could be a measurement error and did not actually occur because the drizzle after the primary peak had a very small intensity (<0.02 mm hr−1) and short duration (<1 hr) in terms of the rainfall hyetograph. For the two validation storm events, the transferred flow hydrographs were fairly reproduced by the model (Figure 7(b) and 7(c)). In addition, the variation pattern and magnitude of the simulated groundwater table beneath subcatchment 14 were compatible with those of the observed groundwater table in the monitoring well adjacent to the subcatchment (Figure 8, coefficient of determination R2 = 0.62 with a near-one slope of 1.09), indicating that the model successfully traced the dynamics of groundwater flow and its interactions with percolation from the overland and surface water bodies (e.g., streams and wetlands). The model underestimated the groundwater table values of greater than 4.6 m possibly because such high groundwater tables are close to the average ground surface elevation (5.13 m above mean sea level) of S14 and might be unrealistic for this subcatchment to be well handled by the model. Thus, it was judged that the SWMM model was good enough to be used for screening effects of SLR on the watershed hydrology, which had a focus of relative changes of peak discharge and flood stage rather than true values of these hydrologic variables.

Figure 7

Plots showing simulated versus transferred runoff (subcatchment 8 in Figure 6) for three selected storm events used for (a) model calibration, and (b) and (c) model validation. The hourly rainfall (input) was averaged by SWMM to the computational time interval of 15 min.

Figure 7

Plots showing simulated versus transferred runoff (subcatchment 8 in Figure 6) for three selected storm events used for (a) model calibration, and (b) and (c) model validation. The hourly rainfall (input) was averaged by SWMM to the computational time interval of 15 min.

Figure 8

Plots showing the simulated versus observed groundwater table (above mean sea level) of subcatchment 14 for the model: (a) evaluation period (1978 to 2013); (b) calibration period (1978 to 1989); and (c) validation period (1990 to 2013).

Figure 8

Plots showing the simulated versus observed groundwater table (above mean sea level) of subcatchment 14 for the model: (a) evaluation period (1978 to 2013); (b) calibration period (1978 to 1989); and (c) validation period (1990 to 2013).

Moreover, the sensitivity analysis indicated that any errors in the model parameters may have more influences on simulated peak discharge than runoff volume, and that the predicted groundwater table can be offset by −2 to 4 cm (Tables 5 and 6). As expected, the two pervious-related parameters (Table 4), namely N-Perv and Dstore-Perv, were least sensitive to the simulation, whereas the two impervious-related parameters, namely N-Imperv and Dstore-Imperv, were relatively sensitive to the simulation, in particular the simulated peak. The variations of these four parameters affected the magnitude, but had no influence on the temporal fluctuation, of simulated groundwater levels. In contrast, the parameter ‘Percent Routed’ was sensitive to simulated runoff, peak, and groundwater and its temporal fluctuation. This indicated importance of the topology (i.e., spatial connectivity) between pervious and impervious areas within the study watershed.

Table 5

Observed and simulated runoff, peak, and groundwater table

  Model parameter Storm event Ia
 
Storm event IIb
 
Storm event IIIc
 
Groundwater tabled
 
Runoff (m3 s−1Peak (m3 s−1Runoff (m3 s−1Peak (m3 s−1Runoff (m3 s−1Peak (m3 s−1Mean (m) Std. (m) 
Observed  0.06 0.23 0.05 0.20 0.09 0.59 4.11 0.43 
Simulated Adoptede 0.07 0.23 0.06 0.20 0.10 0.59 4.51 0.40 
  Model parameter Storm event Ia
 
Storm event IIb
 
Storm event IIIc
 
Groundwater tabled
 
Runoff (m3 s−1Peak (m3 s−1Runoff (m3 s−1Peak (m3 s−1Runoff (m3 s−1Peak (m3 s−1Mean (m) Std. (m) 
Observed  0.06 0.23 0.05 0.20 0.09 0.59 4.11 0.43 
Simulated Adoptede 0.07 0.23 0.06 0.20 0.10 0.59 4.51 0.40 

aFrom 12/07/2011 19:00 to 12/08/2011 0:30 at the outlet of S8 (Figure 6).

bFrom 12/16/2011 19:30 to 12/17/2011 8:00 at the outlet of S8 (Figure 6).

cFrom 12/21/2011 6:00 to 12/22/2011 0:15 at the outlet of S8 (Figure 6).

dFrom 4/20/1978 to 11/19/2013 (not continuous) underneath S14 (Figure 6).

eSee the column ‘Adopted value’ in Table 4.

Table 6

Percentage changes of the simulated runoff, peak, and groundwater table in response to percentage changes of the adopted parameter valuesa

Model parameter
 
Percentage change from the responding simulated value in Table 5 
 
Symbol Change (%) Storm event Ib
 
Storm event IIc
 
Storm event IIId
 
Groundwater tablee
 
Runoff Peak Runoff Peak Runoff Peak Mean Std. 
N-Imperv −20 14.3 13.0  5.0  3.4 −0.2  
20  −8.7  −5.0  −3.4   
33  −13.0  −5.0  −6.8 0.2  
N-Perv −25 14.3        
50 14.3        
75 14.3        
Dstore-Imperv −50 71.4 34.8 33.3 25.0   0.2  
−17 28.6 13.0 16.7 5.0     
17 −14.3 −13.0 −16.7 −5.0   −0.2  
Dstore-Perv −30 14.3      −0.2  
20 14.3        
30 14.3        
Percent Routed −60 14.3 4.3  5.0  3.4 −0.4 2.5 
100  −4.3  −5.0 −10.0 −5.1 0.4 −2.5 
200  −8.7 −16.7 −10.0 −10.0 −10.2 0.9 −2.5 
Model parameter
 
Percentage change from the responding simulated value in Table 5 
 
Symbol Change (%) Storm event Ib
 
Storm event IIc
 
Storm event IIId
 
Groundwater tablee
 
Runoff Peak Runoff Peak Runoff Peak Mean Std. 
N-Imperv −20 14.3 13.0  5.0  3.4 −0.2  
20  −8.7  −5.0  −3.4   
33  −13.0  −5.0  −6.8 0.2  
N-Perv −25 14.3        
50 14.3        
75 14.3        
Dstore-Imperv −50 71.4 34.8 33.3 25.0   0.2  
−17 28.6 13.0 16.7 5.0     
17 −14.3 −13.0 −16.7 −5.0   −0.2  
Dstore-Perv −30 14.3      −0.2  
20 14.3        
30 14.3        
Percent Routed −60 14.3 4.3  5.0  3.4 −0.4 2.5 
100  −4.3  −5.0 −10.0 −5.1 0.4 −2.5 
200  −8.7 −16.7 −10.0 −10.0 −10.2 0.9 −2.5 

aSee Table 4 for definitions of N-Imperv, N-Perv, Dstore-Imperv, Dstore-Perv, and Percent Routed. The blank cells signify zero changes.

bFrom 12/07/2011 19:00 to 12/08/2011 0:30 at the outlet of S8 (Figure 6).

cFrom 12/16/2011 19:30 to 12/17/2011 8:00 at the outlet of S8 (Figure 6).

dFrom 12/21/2011 6:00 to 12/22/2011 0:15 at the outlet of S8 (Figure 6).

eFrom 4/20/1978 to 11/19/2013 (not continuous) underneath S14 (Figure 6).

Simulated effects of SLR

As expected, rising sea level will likely cause increases of peak discharge and flood stage along the streams within the study watershed, with a larger increase at a downstream junction than an upstream one (Figures 8 and 9). For each scenario, which is signified by the number of years after 2013, the model was run for 60 years, generating 60 peak discharges and 60 peak flood stages at each outlet. Spanning the 60 assessment years, the mean of annual peak discharges at the watershed outlet was predicted to increase by 2.58 m3 s−1 after 10 years (hereinafter the reference year was 2013) to 13.84 m3 s−1 after 100 years, whereas the mean of annual peak discharges at the other junctions would increase by up to 13.93 m3 s−1 after 100 years (Figure 9(a)). In correspondence, the mean of annual peak flood stages at the junctions was predicted to increase by 0.015 to 0.92 m (Figure 10(a)), depending on the year and junction of interest. Similarly, the minimums and maximums of annual peak discharges and flood stages will be gradually increasing in the future years (Figures 9(b) and 9(c), and Figures 10(b) and (10(c))). In addition, larger annual peak discharges and flood stages were predicted to occur more often, as indicated by the decreased coefficients of skew (Cs) (Figures 9(d) and 10(d)), and those larger values would become less scattered around the responding means, as indicated by the decreased coefficients of variation (Cv) (Figures 9(e) and 10(e)). This implies that in the future the watershed will probably incur more frequent flooding with a larger magnitude even if the historical climate would reoccur in the future with no trending changes.

Figure 9

Plots showing the predicted (a) mean, (b) minimum, (c) maximum, (d) coefficient of skew (Cs), and (e) coefficient of variation (Cv), of peak discharges at the model junctions shown in Figure 6.

Figure 9

Plots showing the predicted (a) mean, (b) minimum, (c) maximum, (d) coefficient of skew (Cs), and (e) coefficient of variation (Cv), of peak discharges at the model junctions shown in Figure 6.

Figure 10

Plots showing the predicted (a) mean, (b) minimum, (c) maximum, (d) coefficient of variation (Cv), and (e) coefficient of skew (Cs), of peak stages at the model junctions shown in Figure 6.

Figure 10

Plots showing the predicted (a) mean, (b) minimum, (c) maximum, (d) coefficient of variation (Cv), and (e) coefficient of skew (Cs), of peak stages at the model junctions shown in Figure 6.

Historically, the groundwater table beneath the western part of the watershed was higher than that beneath the eastern part (Table 7), but the groundwater flow direction exhibited an inconsistent spatial pattern (Figure 2). Herein, the western part refers to the area drained by Western Branch, while the eastern part refers to the area drained by Eastern Branch and Broad-Linkhorn Bay. Excluding subcatchments 7, 9, 10, and 21, beneath which the mean annual groundwater table was higher than 4.65 m, the spatially averaged mean annual groundwater table beneath the western part ranged from 2.13 to 4.54 m, with a mean of 3.65 m, while the spatially averaged mean annual groundwater table beneath the eastern part ranged from 1.55 to 4.44 m, with a mean of 3.05 m. Such a difference of groundwater table might result in a west-to-east hydraulic gradient of 0.1‰. However, across the western part, the groundwater table beneath either an upper (e.g., 23 or 25) or a lower (e.g., 5 or 14) subcatchment was higher than that beneath a middle subcatchment (e.g., 16 or 19), whereas across the eastern part, the groundwater table had an overall upstream-to-downstream (i.e., southwest-to-northwest) hydraulic gradient. In addition to these overall spatial patterns of groundwater table, subcatchments 17 and 19 were predicted to be groundwater sinks, with net inflows from the shallow aquifers beneath their adjacent subcatchments. In the future, SLR was predicted to continuously raise groundwater table beneath the eastern part resulting from sea water intrusion into the shallow aquifer, but SLR was predicted to have a negligible influence on groundwater table beneath the western part. By 2113 (scenario VI), the increase of groundwater table was predicted to range from 0.10 to 0.29 m, depending on the subcatchment of interest (Tables 5 and 6). Based on the model predictions, on average, every 20 cm rise of sea level would cause a 10 cm increase of groundwater table. As a result of SLR, both the west-to-east hydraulic gradient across the watershed and the southeast-to-northwest hydraulic gradient across the eastern part will likely be reduced, making the regional shallow groundwater more stagnant with smaller discharges into the Chesapeake Bay.

Table 7

Summary statistics of the simulated historical groundwater table and its changes resulting from SLRa

Subcatchment Historical groundwater table (m)
 
Rise of mean groundwater table from the historical mean (m)
 
Mean Min Max Std Skew 10 years later 20 years later 30 years later 50 years later 80 years later 100 years later 
S4 1.55 0.00 2.56 0.33 0.07 0.04 0.06 0.08 0.12 0.18 0.23 
S5 4.40 1.81 4.97 0.48 −1.82 0.00 0.00 0.00 0.00 0.00 0.00 
S6 1.68 0.67 4.69 0.52 1.87 0.05 0.07 0.09 0.14 0.23 0.29 
S7 6.96 4.09 7.37 0.50 −2.66 0.00 0.00 0.00 0.00 0.00 0.00 
S8 2.40 0.93 6.10 0.74 1.74 0.05 0.07 0.09 0.14 0.23 0.29 
S9 4.80 2.21 5.70 0.60 −1.42 0.02 0.03 0.04 0.06 0.09 0.12 
S10 5.06 3.75 8.87 0.81 1.65 0.05 0.07 0.09 0.13 0.21 0.27 
S11 3.48 0.87 4.13 0.52 −1.58 0.02 0.02 0.03 0.05 0.08 0.10 
S12 3.49 1.25 4.30 0.48 −0.90 0.00 0.00 0.00 0.00 0.00 0.00 
S13 3.60 0.73 3.91 0.49 −3.71 0.01 0.01 0.01 0.02 0.03 0.03 
S14 4.43 1.68 5.05 0.54 −1.70 0.01 0.02 0.02 0.04 0.06 0.08 
S15 2.72 1.34 5.89 0.68 1.51 0.05 0.07 0.09 0.14 0.23 0.29 
S16 3.07 1.27 4.87 0.67 0.21 0.04 0.06 0.07 0.12 0.19 0.24 
S17 2.78 0.44 3.57 0.52 −0.92 0.03 0.04 0.05 0.07 0.11 0.15 
S18 4.54 2.56 6.10 0.65 −0.15 0.00 0.00 0.00 0.00 0.02 0.04 
S19 2.13 0.02 2.88 0.42 −0.48 0.03 0.04 0.05 0.07 0.12 0.15 
S20 4.32 1.38 4.53 0.41 −4.67 0.00 0.01 0.01 0.01 0.02 0.02 
S21 4.66 2.74 6.46 0.73 0.02 0.04 0.05 0.07 0.11 0.18 0.23 
S22 3.69 0.78 3.88 0.41 −4.94 0.00 0.00 0.00 0.01 0.01 0.01 
S23 2.57 1.68 5.19 0.46 1.97 0.00 0.00 0.00 0.00 0.00 0.00 
S24 4.44 2.21 5.87 0.71 −0.31 0.04 0.05 0.07 0.10 0.16 0.20 
S25 4.35 1.65 4.90 0.50 −2.02 0.00 0.00 0.00 0.00 0.00 0.00 
Subcatchment Historical groundwater table (m)
 
Rise of mean groundwater table from the historical mean (m)
 
Mean Min Max Std Skew 10 years later 20 years later 30 years later 50 years later 80 years later 100 years later 
S4 1.55 0.00 2.56 0.33 0.07 0.04 0.06 0.08 0.12 0.18 0.23 
S5 4.40 1.81 4.97 0.48 −1.82 0.00 0.00 0.00 0.00 0.00 0.00 
S6 1.68 0.67 4.69 0.52 1.87 0.05 0.07 0.09 0.14 0.23 0.29 
S7 6.96 4.09 7.37 0.50 −2.66 0.00 0.00 0.00 0.00 0.00 0.00 
S8 2.40 0.93 6.10 0.74 1.74 0.05 0.07 0.09 0.14 0.23 0.29 
S9 4.80 2.21 5.70 0.60 −1.42 0.02 0.03 0.04 0.06 0.09 0.12 
S10 5.06 3.75 8.87 0.81 1.65 0.05 0.07 0.09 0.13 0.21 0.27 
S11 3.48 0.87 4.13 0.52 −1.58 0.02 0.02 0.03 0.05 0.08 0.10 
S12 3.49 1.25 4.30 0.48 −0.90 0.00 0.00 0.00 0.00 0.00 0.00 
S13 3.60 0.73 3.91 0.49 −3.71 0.01 0.01 0.01 0.02 0.03 0.03 
S14 4.43 1.68 5.05 0.54 −1.70 0.01 0.02 0.02 0.04 0.06 0.08 
S15 2.72 1.34 5.89 0.68 1.51 0.05 0.07 0.09 0.14 0.23 0.29 
S16 3.07 1.27 4.87 0.67 0.21 0.04 0.06 0.07 0.12 0.19 0.24 
S17 2.78 0.44 3.57 0.52 −0.92 0.03 0.04 0.05 0.07 0.11 0.15 
S18 4.54 2.56 6.10 0.65 −0.15 0.00 0.00 0.00 0.00 0.02 0.04 
S19 2.13 0.02 2.88 0.42 −0.48 0.03 0.04 0.05 0.07 0.12 0.15 
S20 4.32 1.38 4.53 0.41 −4.67 0.00 0.01 0.01 0.01 0.02 0.02 
S21 4.66 2.74 6.46 0.73 0.02 0.04 0.05 0.07 0.11 0.18 0.23 
S22 3.69 0.78 3.88 0.41 −4.94 0.00 0.00 0.00 0.01 0.01 0.01 
S23 2.57 1.68 5.19 0.46 1.97 0.00 0.00 0.00 0.00 0.00 0.00 
S24 4.44 2.21 5.87 0.71 −0.31 0.04 0.05 0.07 0.10 0.16 0.20 
S25 4.35 1.65 4.90 0.50 −2.02 0.00 0.00 0.00 0.00 0.00 0.00 

aThe groundwater table beneath the subcatchments in bold would not be noticeably affected by SLR.

DISCUSSION

The observations at Sewells Point indicated that sea level has been steadily rising at 4.21 to 5.16 mm yr−1 (depending on the month of interest) over the past 65 years. This is compatible with the value (4.52 ± 0.66 mm yr−1) reported by Boon et al. (2010) but higher than the global average SLR rate of 1.3 to 2.3 mm yr−1 over 1961 to 2003 reported by IPCC (2007). In addition, the observations revealed that a significant abrupt increase occurred around 1983. The rising rate after 1983 became much faster than that before this year (5.8 versus 3.9 mm yr−1 at annual time scale). At global scale, an abrupt increase of rising rate occurred around 1993, with a post-1993 rate of 2.4 to 3.8 mm yr−1 (IPCC 2007).

As with previous studies (e.g., Nuttle & Portnoy 1992; Masterson et al. 2014), this study found that the groundwater table will increase in response to SLR. On average, every 20 cm SLR would raise groundwater table by about 10 cm. In addition, as a result of SLR, both larger peak discharges and higher flood stages would occur more often, which is consistent with the findings of Nuttle & Portnoy (1992) and Ashofteh et al. (2016). The larger peak discharges can mainly be attributed to the fact that the increased groundwater table will likely elevate average levels of soil moisture, whereas the higher flood stages can be attributed to the back-water effects of rising sea level as well as the larger discharges. Further, for coastal watersheds, including the Lynnhaven River watershed, the effects of SLR on hydrology can be exacerbated by increasing rainfall intensity resulting from climate change (Muschinski & Katz 2013).

CONCLUSIONS

This study examined the long-term temporal trends of sea level, and set up and used a SWMM model to predict effects of the trends on hydrology of the Lynnhaven River watershed, a mid-Atlantic coastal watershed that is hydraulically connected with the Chesapeake Bay. The results indicated that the sea level has been consistently rising at an annual average rate of 4.65 mm yr−1 over the past 65 years; with an abrupt increase around 1983. At monthly time scale, the rising rates varied from 4.21 to 5.16 mm yr−1. As a result of SLR, larger peak discharges and flood stages were predicted to occur more often and have a longer duration. This implies that in the future the study watershed will likely incur more frequent flooding with a larger magnitude. By 2113, the mean peak discharge and flood stage at the lower part of the watershed would be increased by 13.84 m3 s−1 and 0.92 m, respectively. These predicted increases are expected to become larger if the probable increase of rainfall intensity (Kunkel et al. 2013) were taken into account. Moreover, SLR will likely make more saltwater intrude into the shallow aquifer beneath the eastern part than that beneath the western part of the study watershed. Nevertheless, the groundwater table will likely be elevated in the future. It was predicted that every 20 cm rise of sea level would raise the groundwater table by about 10 cm.

ACKNOWLEDGEMENTS

This study was financially supported by Old Dominion University (ODU) Resilience Collaborative and the Office of Research at ODU, ODU Multidisciplinary Seed Funding Program (MSF), ODU Faculty Proposal Preparation Program (FP3), and international collaboration agreement (# 5CEC3) between ODU and Inner Mongolia Agricultural University (IMAU). Moreover, thanks are extended to Ms Amanda Kennedy, a master's student and a graduate research assistant at ODU, for her great efforts in collecting and providing the flow data.

REFERENCES

REFERENCES
Ashofteh
,
P. S.
,
Bozorg-Haddad
,
O.
,
Loáiciga
,
H. A.
&
Marino
,
M. A.
2016
Evaluation of the impacts of climate variability and human activity on streamflow at the basin scale
.
Journal of Irrigation and Drainage Engineering
142
(
8
),
04016028
.
Barlow
,
P. M.
2003
Ground Water in Freshwater-Saltwater Environments of the Atlantic Coast
.
Circular 1262
.
U.S. Geological Survey
,
Reston, VA
.
Bjerklie
,
D. M.
,
Mullaney
,
J. R.
,
Stone
,
J. R.
,
Skinner
,
B. J.
&
Ramlow
,
M. A.
2012
Preliminary Investigation of the Effects of Sea-Level Rise on Groundwater Levels in New Haven, Connecticut
.
Open-File Report 2012-1025
.
U.S. Geological Survey
,
Reston, VA
.
Boon
,
J. D.
,
Brubaker
,
J. M.
&
Forrest
,
D. R.
2010
Chesapeake Bay Land Subsidence and Sea Level Change: An Evaluation of Past and Present Trends and Future Outlook
.
A Report to the U.S. Army Corps of Engineers Norfolk District. Available at: www.vims.edu/GreyLit/VIMS/sramsoe425.pdf
(
accessed 15 May 2015
).
Charbeneau
,
R. J.
2006
Groundwater Hydraulics and Pollutant Transport
.
Waveland Press, Inc.
,
Long Grove, IL
.
Cowardin
,
L. M.
,
Golet
,
F. C.
&
LaRoe
,
E. T.
1979
Classification of Wetlands and Deepwater Habitats of the United States
.
FWS/OBS-79/31
.
U.S. Department of the Interior
,
Washington, DC
.
CVB (City of Virginia Beach, Virginia) & VDEQ (Virginia Department of Environmental Quality)
undated
Lynnhaven River Watershed, Virginia Beach, Virginia: Application for Federal No Discharge Zone Designation
. (
accessed 18 January 2015
).
Eggleston
,
J.
&
Pope
,
J.
2013
Land Subsidence and Relative Sea-Level Rise in the Southern Chesapeake Bay Region
.
U.S. Geological Survey Circular 1392
.
U.S. Geological Survey
,
Reston, VA
.
FGDC (Federal Geographic Data Committee)
1996
Classification of Wetlands and Deepwater Habitats of the United States
.
FGDC–STD-004-2013 (2nd edn)
.
U.S. Geological Survey
,
Reston, VA
.
Finnemore
,
E. J.
&
Franzini
,
J. B.
2002
Fluid Mechanics with Engineering Applications
,
10th edn
.
The McGraw-Hill Companies
,
New York, NY
.
Garcia
,
E. S.
&
Loáiciga
,
H. A.
2015
Sea-level rise and flooding in coastal riverine flood plains
.
Hydrological Sciences Journal
59
(
1
),
204
220
.
http://dx.doi.org/10.1080/02626667.2013.798660
.
Gironás
,
J.
,
Roesner
,
L. A.
,
Davis
,
J.
&
Rossman
,
L. A.
2009
Storm Water Management Model Applications Manual
.
U.S. Environmental Protection Agency
,
Cincinnati, OH
.
Goodwin
,
L.
2015
Climate Resiliency Outcomes Management Strategy (Draft)
. (
accessed 15 May 2015
).
Hamed
,
K. H.
&
Rao
,
A. R.
1998
A modified Mann-Kendall trend test for autocorrelated data
.
Journal of Hydrology
204
,
182
196
.
Hargreaves
,
G. H.
&
Samani
,
Z. A.
1985
Reference crop evapotranspiration from temperature
.
Applied Engineering in Agriculture
1
(
2
),
96
99
.
Heywood
,
C. E.
&
Pope
,
J. P.
2009
Simulation of Groundwater Flow in the Coastal Plain Aquifer System of Virginia
.
Scientific Investigations Report 2009-5039
.
U.S. Geological Survey
,
Reston, VA
.
Hilpert
,
M.
&
Glantz
,
R.
2013
Exploring the parameter space of the Green–Ampt model
.
Advances in Water Resources
53
,
225
230
.
Hirsch
,
R. M.
,
Slack
,
J. R.
&
Smith
,
R. A.
1982
Techniques for trend assessment for monthly water quality data
.
Water Resources Research
18
,
107
121
.
Hong
,
B.
,
Panday
,
N.
,
Shen
,
J.
,
Wang
,
H. V.
,
Gong
,
W.
&
Soehl
,
A.
2010
Modeling water exchange between Baltimore Harbor and Chesapeake Bay using artificial tracers: seasonal variations
.
Marine Environmental Research
70
(
1
),
102
119
.
IPCC (Intergovernmental Panel on Climate Change)
2007
Summary for policymakers
. In:
Climate Change 2007: The Physical Science Basis
(
Qin
,
S.D.
,
Manning
,
M.
,
Chen
,
Z.
,
Marquis
,
M.
,
Averyt
,
K. B.
,
Tignor
,
M.
&
Miller
,
H. L.
, eds).
Cambridge University Press
,
Cambridge
,
UK
.
Karney
,
C. F. F.
2011
Transverse Mercator with an accuracy of a few nanometers
.
Journal of Geodesy
85
(
8
),
475
485
.
Kivett
,
J.
2014
Sea level rise adaptation efforts at the regional scale: Water management
. In:
The Sixth Annual Southeast Florida Climate Leadership Summit
.
South Florida Water Management District
,
West Palm Beach, FL
.
Kunkel
,
K. E.
,
Stevens
,
L. E.
,
Stevens
,
S. E.
,
Sun
,
L.
,
Janssen
,
E.
,
Wuebbles
,
D.
,
Rennells
,
J.
,
DeGaetano
,
A.
&
Dobson
,
J. G.
2013
Regional Climate Trends and Scenarios for the U.S. National Climate Assessment, Part 1: Climate of the Northeast U.S
.
NOAA Technical Report NESDIS 142-1
.
U.S. Department of Commerce
,
Washington, DC
.
Langevin
,
C. D.
,
Thorne
,
D. J.
,
Dausman
,
A.
,
Sukop
,
M.
&
Guo
,
W.
2007
SEAWAT (ver. 4): A Computer Program for Simulation of Multi-Species Solute and Heat Transport
.
U.S. Geological Survey
,
Denver, CO
.
Mann
,
H. B.
1945
Nonparametric tests against trend
.
Econometrica
13
,
245
259
.
Masterson
,
J. P.
,
Fienen
,
M. N.
,
Thieler
,
E. R.
,
Gesch
,
D. B.
,
Gutierrez
,
B. T.
&
Plant
,
N. G.
2014
Effects of sea-level rise on barrier island groundwater system dynamics: ecohydrological implications
.
Ecohydrology
7
,
1064
1071
.
McGilchrist
,
C. A.
&
Woodyer
,
K. D.
1975
Note on a distribution-free CUSUM technique
.
Technometrics
17
(
3
),
321
325
.
McKean
,
J.
,
Nagel
,
D.
,
Tonina
,
D.
,
Bailey
,
P.
,
Wright
,
C. W.
,
Bohn
,
C.
&
Nayegandhi
,
A.
2009
Remote sensing of channels and riparian zones with a narrow-beam aquatic-terrestrial lidar
.
Remote Sensing
1
,
1065
1096
.
doi:10.3390/rs1041065
.
MESO (Multi-Community Environmental Storm Observatory)
2002
NOR'EASTERS: Comprehension, Preparation and Survival
. (
accessed 19 January 2015
).
Morel-Seytoux
,
H. J.
1984
Some recent developments in physically based rainfall-runoff modeling
. In:
Frontiers in Hydrology
.
Water Resources Publications
,
Littleton, CO
.
Morel-Seytoux
,
H. J.
1985
Conjunctive use of surface and ground waters
. In:
Artificial Recharge of Groundwater
(
Asano
,
T.
, ed.).
Butterworth Publishers
,
Boston, MA
.
MRLC (Multi-Resolution Land Characteristics Consortium)
2014
National Land Cover Database 2011 (NLCD 2011)
. (
accessed 20 January 2015
).
Muschinski
,
T.
&
Katz
,
J. I.
2013
Trends in hourly rainfall statistics in the United States under a warming climate
.
Nature Climate Change
3
,
577
580
.
doi:10.10.38/NCLIMATE1828
.
Najjar
,
R. G.
,
Pyke
,
C. R.
,
Adams
,
M. B.
,
Breitburg
,
D.
,
Hershner
,
C.
,
Kemp
,
M.
,
Howarth
,
R.
,
Mulholland
,
M.
,
Paolisso
,
M.
,
Secor
,
D.
,
Sellner
,
K.
,
Wardrop
,
D.
&
Wood
,
R.
2010
Potential climate-change impacts on the Chesapeake Bay
.
Estuarine, Coastal and Shelf Science
86
,
1
20
.
Nuttle
,
W. K.
&
Portnoy
,
J. W.
1992
Effect of rising sea level on runoff and groundwater discharge to coastal ecosystems
.
Estuarine, Coastal and Shelf Science
34
,
203
212
.
Paolisso
,
M.
,
Douglas
,
E.
,
Enrici
,
A.
,
Kirshen
,
P.
,
Watson
,
C.
&
Ruth
,
M.
2012
Climate change, justice, and adaptation among African American communities in the Chesapeake Bay region
.
Weather, Climate, and Society
4
(
1
),
34
47
.
Passeri
,
D. L.
,
Hagen
,
S. C.
,
Medeiros
,
S. C.
,
Bilskie
,
M. V.
,
Alizad
,
K.
&
Wang
,
D.
2015
The dynamic effects of sea level rise on low-gradient coastal landscapes: a review
.
Earth's Future
3
,
59
181
.
doi:10.1002/2015EF000298
.
Rice
,
K. C.
&
Hirsch
,
R. M.
2012
Spatial and Temporal Trends in Runoff at Long-Term Streamgages Within and Near the Chesapeake Bay Watershed
.
Scientific Investigation Report 2012-515
.
U.S. Geological Survey
,
Reston, VA
.
Roberson
,
J. A.
,
Cassidy
,
J. J.
&
Chaudhry
,
M. H.
1998
Hydraulic Engineering
,
2nd edn
.
John Wiley & Sons, Inc.
,
New York, NY
.
Rossman
,
L. A.
2010
Storm Water Management Model User's Manual (Version 5.0)
.
U.S. Environmental Protection Agency
,
Cincinnati, OH
.
Sallenger
,
A. H.
,
Doran
,
K. S.
&
Howd
,
P. A.
2012
Hotspot of accelerated sea-level rise on the Atlantic coast of North America
.
Nature Climate Change
2
,
884
888
.
Saxton
,
K. E.
&
Rawls
,
W. J.
2006
Soil water characteristic estimates by texture and organic matter for hydrologic solutions
.
Soil Science Society of America Journal
70
,
1569
1578
.
Saxton
,
K. E.
,
Rawls
,
W. J.
,
Romberger
,
J. S.
&
Papendick
,
R. I.
1986
Estimating generalized soil-water characteristics from texture
.
Transactions of the ASAE
50
,
1031
1035
.
Schoenbaum
,
T. J.
1988
Islands, Capes, and Sounds: The North Carolina Coast
.
John F. Blair
,
Winston-Salem, NC
.
Sen
,
P. K.
1968
Robustness of some nonparametric procedures in linear models
.
The Annals of Mathematical Statistics
39
(
6
),
1913
1922
.
Smith
,
B. S.
&
Harlow
,
G. E.
Jr
2002
Conceptual Hydrogeologic Framework of the Shallow Aquifer System at Virginia Beach, Virginia
.
Water-Resources Investigations Report 01-4262
.
U.S. Geological Survey
,
Denver, CO
.
Sneyers
,
R.
1990
Technical Note No 143 on the Statistical Analysis of Series of Observations
.
World Meteorological Organization
,
Geneva
.
Taubenheim
,
J.
1989
An easy procedure for detecting a discontinuity in a digital time series
.
Meteorologische Zeitschrift
39
(
6
),
344
347
.
Todd
,
D. K.
&
Mays
,
L. W.
2005
Groundwater Hydrology
,
3rd edn
.
John Wiley & Sons, Inc.
,
Hoboken, NJ
.
Urquhart
,
E. A.
,
Zaitchik
,
B. F.
,
Waugh
,
D. W.
,
Guikema
,
S. D.
&
Del Castillo
,
C. E.
2014
Uncertainty in model predictions of Vibrio vulnificus response to climate variability and change: a Chesapeake Bay case study
.
PLoS ONE
9
(
5
),
e98256
.
doi:10.1371/journal.pone.0098256
.
USACE (U.S. Army Corps of Engineers)
2010
Hydraulic Reference Manual of HEC-RAS River Analysis System (ver. 4.1)
.
Institute for Water Resources and Hydrologic Engineering Center
,
Davis, CA
.
USACE (U.S. Army Corps of Engineers)
2013
Draft Feasibility Report and Integrated Environmental Assessment: Lynnhaven River Basin Ecosystem Restoration
.
Main Report
.
U.S. Army Corps of Engineers Norfolk District
,
Norfolk, VA
.
USDA-NRCS (U.S. Department of Agriculture Natural Resources Conservation Service)
2005
Description of SSURGO Database
. (
accessed 20 January 2015
).
USDA-SCS (U.S. Department of Agriculture Soil Conservation Service)
1987
Soil Mechanics Level I Module 3: USDA Textural Soil Classification
.
Study Guide
.
U.S. Department of Agriculture
,
Washington, DC
.
USDOC (U.S. Department of Commerce)
2014
Sea Level Rise and Nuisance Flood Frequency Changes Around the United States
.
Technical Report NOS CO-OPS 073.58
.
Center for Operational Oceanographic Products and Services at National Oceanic and Atmospheric Administration
,
Silver Spring, MD
. (
accessed 15 May 2015
).
USEPA (U.S. Environmental Protection Agency)
2014
Climate Impacts on Coastal Areas
. (
accessed 15 May 2015
).
USFWS (U.S. Fish and Wildlife Service)
2014
National Wetlands Inventory
. (
accessed 20 January 2015
).
USGS
2014a
About NED
. (
accessed 20 January 2015
).
USGS
2014b
Hydrography
.
Available at: http://nhd.usgs.gov
(
accessed 20 January 2015
).
USGS (U.S. Geological Survey)
2002
Conceptual Hydrogeologic Framework of the Shallow Aquifer System at Virginia Beach, Virginia
.
Water-Resources Investigations Report 01-4262
.
U.S. Geological Survey and U.S. Department of the Interior
,
Richmond, VA
.
Varnell
,
L.
2014
Shoreline energy and sea level dynamics in Lower Chesapeake Bay: history and patterns
.
Estuaries and Coasts
37
(
2
),
508
523
.
doi:10.1007/s12237-013-9672-6
.
Viessman
,
W.
Jr
&
Lewis
,
G. L.
2003
Introduction to Hydrology
,
5th edn
.
Pearson Education, Inc.
,
Upper Saddle River, NJ
.
Walsh
,
J.
,
Wuebbles
,
D.
,
Hayhoe
,
K.
,
Kossin
,
J.
,
Kunkel
,
K.
,
Stephens
,
G.
,
Thorne
,
P.
,
Vose
,
R.
,
Wehner
,
M.
,
Willis
,
J.
,
Anderson
,
D.
,
Doney
,
S.
,
Feely
,
R.
,
Hennon
,
P.
,
Kharin
,
V.
,
Knutson
,
T.
,
Landerer
,
F.
,
Lenton
,
T.
,
Kennedy
,
J.
,
Somerville
,
R.
2014
Chapter 2: Our changing climate
. In:
Climate Change Impacts in the United States: The Third National Climate Assessment
(
Melillo
,
J. M.
,
Richmond
,
T. C.
&
Yohe
,
G. W.
, eds).
U.S. Global Change Research Program
,
Washington, DC
, pp.
19
67
.
Wang
,
X.
&
Melesse
,
A. M.
2006
Effects of STATSGO and SSURGO as inputs on SWAT model's snowmelt simulation
.
Journal of the American Water Resources Association
42
(
5
),
1217
1236
.
Wang
,
X.
,
Shang
,
S.
,
Qu
,
Z.
,
Liu
,
T.
,
Melesse
,
A. M.
&
Yang
,
W.
2010
Simulated wetland conservation-restoration effects on water quantity and quality at watershed scale
.
Journal of Environmental Management
91
,
1511
1525
.
Wang
,
X.
,
Yang
,
X.
,
Liu
,
T.
,
Li
,
F.
,
Gao
,
R.
,
Duan
,
L.
&
Luo
,
Y.
2014
Trend and extreme occurrence of precipitation in a mid-latitude Eurasian steppe watershed at various time scales
.
Hydrological Processes
28
,
5547
5560
.
Xian
,
G.
,
Homer
,
C.
,
Dewitz
,
J.
,
Fry
,
J.
,
Hossain
,
N.
&
Wickham
,
J.
2011
The change of impervious surface area between 2001 and 2006 in the conterminous United States
.
Photogrammetric Engineering and Remote Sensing
77
(
8
),
758
762
.
Yin
,
J.
,
Schlesinger
,
M. E.
&
Stouffer
,
R. J.
2009
Model projections of rapid sea-level rise on the northeast coast of the United States
.
Nature Geoscience
2
,
262
266
.
Young
,
R. S.
,
Tormey
,
B. R.
&
Peek
,
K. M.
2016
Coastal Hazard and Sea-Level Rise Vulnerability Assessment
.
Program for the Study of Developed Shorelines Western Carolina University
. (
accessed 30 January 2017
).