Joseph D. Hughes

Scripting MODFLOW model development using Python and FloPy

Graphical user interfaces (GUIs) are commonly used to construct and postprocess numerical groundwater flow and transport models. Scripting model development with the programming language Python is presented here as an alternative approach. One advantage of Python is that there are many packages available to facilitate the model development process, including packages for plotting, array manipulation, optimization, and data analysis. For MODFLOW-based models, the FloPy package was developed by the authors to construct model input files, run the model, and read and plot simulation results. Use of Python with the available scientific packages and FloPy facilitates data exploration, alternative model evaluations, and model analyses that can be difficult to perform with GUIs. Furthermore, Python scripts are a complete, transparent, and repeatable record of the modeling process. The approach is introduced with a simple FloPy example to create and postprocess a MODFLOW model. A more complicated capture-fraction analysis with a real-world model is presented to demonstrate the types of analyses that can be performed using Python and FloPy.

Quantifying the predictive consequences of model error with linear subspace analysis

All computer models are simplified and imperfect simulators of complex natural systems. The discrepancy arising from simplification induces bias in model predictions, which may be amplified by the process of model calibration. This paper presents a new method to identify and quantify the predictive consequences of calibrating a simplified computer model. The method is based on linear theory, and it scales efficiently to the large numbers of parameters and observations characteristic of groundwater and petroleum reservoir models. The method is applied to a range of predictions made with a synthetic integrated surface-water/groundwater model with thousands of parameters. Several different observation processing strategies and parameterization/regularization approaches are examined in detail, including use of the Karhunen-Loeve parameter transformation. Predictive bias arising from model error is shown to be prediction specific and often invisible to the modeler. The amount of calibration-induced bias is influenced by several factors, including how expert knowledge is applied in the design of parameterization schemes, the number of parameters adjusted during calibration, how observations and model-generated counterparts are processed, and the level of fit with observations achieved through calibration. Failure to properly implement any of these factors in a prediction-specific manner may increase the potential for predictive bias in ways that are not visible to the calibration and uncertainty analysis process.

Feedback of land subsidence on the movement and conjunctive use of water resources

The dependency of surface- or groundwater flows and aquifer hydraulic properties on dewatering-induced layer deformation is not available in the USGSs groundwater model MODFLOW. A new integrated hydrologic model, MODFLOW-OWHM, formulates this dependency by coupling mesh deformation with aquifer transmissivity and storage and by linking land subsidence/uplift with deformation-dependent flows that also depend on aquifer head and other flow terms. In a test example, flows most affected were stream seepage and evapotranspiration from groundwater (ETgw). Deformation feedback also had an indirect effect on conjunctive surface- and groundwater use components: Changed stream seepage and streamflows influenced surface-water deliveries and returnflows. Changed ETgw affected irrigation demand, which jointly with altered surface-water supplies resulted in changed supplemental groundwater requirements and pumping and changed return runoff. This modeling feature will improve the impact assessment of dewatering-induced land subsidence/uplift (following irrigation pumping or coal-seam gas extraction) on surface receptors, inter-basin transfers, and surface-infrastructure integrity. We develop a method to simulate deformation-dependent flows for MODFLOW.We demonstrate the significance of linking subsidence to conjunctive water use.The linkages affect flows across the landscape, surface water, and groundwater.Linked flows are relevant to resource issues that include conjunctive water use.

MODFLOW-Based Coupled Surface Water Routing and Groundwater-Flow Simulation.

In this paper, we present a flexible approach for simulating one- and two-dimensional routing of surface water using a numerical surface water routing (SWR) code implicitly coupled to the groundwater-flow process in MODFLOW. Surface water routing in SWR can be simulated using a diffusive-wave approximation of the Saint-Venant equations and/or a simplified level-pool approach. SWR can account for surface water flow controlled by backwater conditions caused by small water-surface gradients or surface water control structures. A number of typical surface water control structures, such as culverts, weirs, and gates, can be represented, and it is possible to implement operational rules to manage surface water stages and streamflow. The nonlinear system of surface water flow equations formulated in SWR is solved by using Newton methods and direct or iterative solvers. SWR was tested by simulating the (1) Lal axisymmetric overland flow, (2) V-catchment, and (3) modified Pinder-Sauer problems. Simulated results for these problems compare well with other published results and indicate that SWR provides accurate results for surface water-only and coupled surface water/groundwater problems. Results for an application of SWR and MODFLOW to the Snapper Creek area of Miami-Dade County, Florida, USA are also presented and demonstrate the value of coupled surface water and groundwater simulation in managed, low-relief coastal settings.

Economic impacts of urban flooding in South Florida: Potential consequences of managing groundwater to prevent salt water intrusion

High-value urban zones in coastal South Florida are considered particularly vulnerable to salt water intrusion into the groundwater-based, public water supplies caused by sea level rise (SLR) in combination with the low topography, existing high water table, and permeable karst substrate. Managers in the region closely regulate water depths in the extensive South Florida canal network to control closely coupled groundwater levels and thereby reduce the risk of saltwater intrusion into the karst aquifer. Potential SLR adaptation strategies developed by local managers suggest canal and groundwater levels may have to be increased over time to prevent the increased salt water intrusion risk to groundwater resources. However, higher canal and groundwater levels cause the loss of unsaturated zone storage and lead to an increased risk of inland flooding when the recharge from rainfall exceeds the capacity of the unsaturated zone to absorb it and the water table reaches the surface. Consequently, higher canal and groundwater levels are also associated with increased risk of economic losses, especially during the annual wet seasons. To help water managers and urban planners in this region better understand this trade-off, this study models the relationships between flood insurance claims and groundwater levels in Miami-Dade County. Via regression analyses, we relate the incurred number of monthly flood claims in 16 Miami-Dade County watersheds to monthly groundwater levels over the period from 1996 to 2010. We utilize these estimated statistical relationships to further illustrate various monthly flood loss scenarios that could plausibly result, thereby providing an economic quantification of a too much water trade-off. Importantly, this understanding is the first of its kind in South Florida and is exceedingly useful for regional-scale hydro-economic optimization models analyzing trade-offs associated with high water levels.

Use of general purpose graphics processing units with MODFLOW

To evaluate the use of general-purpose graphics processing units (GPGPUs) to improve the performance of MODFLOW, an unstructured preconditioned conjugate gradient (UPCG) solver has been developed. The UPCG solver uses a compressed sparse row storage scheme and includes Jacobi, zero fill-in incomplete, and modified-incomplete lower-upper (LU) factorization, and generalized least-squares polynomial preconditioners. The UPCG solver also includes options for sequential and parallel solution on the central processing unit (CPU) using OpenMP. For simulations utilizing the GPGPU, all basic linear algebra operations are performed on the GPGPU; memory copies between the central processing unit CPU and GPCPU occur prior to the first iteration of the UPCG solver and after satisfying head and flow criteria or exceeding a maximum number of iterations. The efficiency of the UPCG solver for GPGPU and CPU solutions is benchmarked using simulations of a synthetic, heterogeneous unconfined aquifer with tens of thousands to millions of active grid cells. Testing indicates GPGPU speedups on the order of 2 to 8, relative to the standard MODFLOW preconditioned conjugate gradient (PCG) solver, can be achieved when (1) memory copies between the CPU and GPGPU are optimized, (2) the percentage of time performing memory copies between the CPU and GPGPU is small relative to the calculation time, (3) high-performance GPGPU cards are utilized, and (4) CPU-GPGPU combinations are used to execute sequential operations that are difficult to parallelize. Furthermore, UPCG solver testing indicates GPGPU speedups exceed parallel CPU speedups achieved using OpenMP on multicore CPUs for preconditioners that can be easily parallelized.

EOF Analysis of Water Level Variations for Microtidal and Mangrove-Covered Frog Creek System, West–Central Florida

Abstract Zhang, J.; Wang, P., and Hughes, J., 2012. EOF analysis of water level variations for microtidal and mangrove-covered Frog Creek system, west-central Florida. At present, little is known about the mechanisms that control water level variations for the Frog Creek system, west–central Florida, which is a mangrove-covered, microtidal, shallow estuary. Although connected with open waters, topographic and roughness constraints can alter expected variability of water levels in this area. Water level variations are especially important for the Frog Creek system, because the increases or decreases in subtidal water level change the present tidal environments. In this paper, 1-year measurements of water levels were examined using the empirical orthogonal function (EOF) method. Water level data from eight stations were first filtered by employing a Lancz6-squared filter because of its ability to eliminate diurnal tidal signals while inducing minimal attenuation in variations with periods of 2 to 3 days, which are usually caused by wind. The filtered water levels were then analyzed using the EOF method to determine relative contributions of various forcing factors. The first three EOF modes explained 64.13%, 33.25%, and 1.55% of total variances and were interpreted by calculating correlation coefficients and coherences. We attributed EOF mode 1 to the variations of coastal subtidal water levels introduced by physical processes with longer periods, such as spring-neap tidal variability and seasonal freshwater river discharge variability. Next, EOF mode 2 was related to the effects of the local N-S wind component, which introduces short-timescale variations of the subtidal water level with typical periods ranging from 2 to 6 days. The effect of river discharge on water level was captured by EOF mode 3. Overall, EOF analysis provides an insightful tool to examine water level disturbances.

MODFLOW 6, the U.S. Geological Survey Modular Hydrologic Model

MODFLOW is a popular open-source groundwater flow model distributed by the U.S. Geological Survey. For over 30 years, the MODFLOW program has been widely used by academics, private consultants, and government scientists to accurately, reliably, and efficiently simulate groundwater flow. With time, growing interest in surface and groundwater interactions, local refinement with nested and unstructured grids, karst groundwater flow, solute transport, and saltwater intrusion, has led to the development of numerous MODFLOW versions. Although these MODFLOW versions are often based on the core MODFLOW version (presently MODFLOW-2005), there are often incompatibilities that restrict their use with other MODFLOW versions. In many cases, development of these alternative MODFLOW versions has been challenging due to the underlying program structure, which was designed for the simulation of a single groundwater flow model using a regular MODFLOW grid consisting of layers, rows, and columns. n A new object-oriented program and underlying framework called MODFLOW 6 was developed to provide a platform for supporting multiple models and multiple types of models within the same simulation. This version of MODFLOW is labeled with a 6 because it is the sixth core version of MODFLOW to be released by the USGS. In the new design, any number of models can be included in a simulation. These models can be independent of one another with no interaction, they can exchange information with one another, or they can be tightly coupled at the matrix level by adding them to the same numerical solution. Transfer of information between models is isolated to exchange objects, which allow models to be developed and used independently of one another. Within this new framework, a regional-scale groundwater model may be coupled with multiple local-scale groundwater models. Or, a surface-water flow model could be coupled to multiple groundwater flow models. The framework naturally allows for future extensions to include the simulation of solute transport. n The first release of MODFLOW 6 contains one type of hydrologic model, the Groundwater Flow (GWF) Model. The GWF Model for MODFLOW 6 is based on a generalized control-volume finite-difference (CVFD) approach in which a cell can be hydraulically connected to any number of surrounding cells. Users can define the model grid using: n 1. a regular MODFLOW grid consisting of layers, rows, and columns, n 2. a layered grid defined by (x, y) vertex pairs, or n 3. a general unstructured grid based on concepts developed for MODFLOW-USG. n For complex problems involving water-table conditions, an optional Newton-Raphson formulation, based on the formulations in MODFLOW-NWT and MODFLOW-USG, can be activated. The GWF Model is divided into packages, as was done in previous MODFLOW versions. A package is the part of the model that deals with a single aspect of simulation. Packages included with the GWF Model include those related to internal calculations of groundwater flow (discretization, initial conditions, hydraulic conductance, and storage), stress packages (constant heads, wells, recharge, rivers, general head boundaries, drains, and evapotranspiration), and advanced stress packages (streamflow routing, lakes, multi-aquifer wells, and unsaturated zone flow). An additional package is also available for moving water available in one package into the individual features of the advanced stress packages. The GWF Model also has packages for obtaining and controlling output from the model. n A new groundwater flow formulation was developed specifically for the GWF Model in MODFLOW 6. This new flow formulation is called the XT3D option. The XT3D option extends the capabilities of MODFLOW by enabling simulation of fully three-dimensional anisotropy on regular or irregular grids in a way that properly takes into account the full, three-dimensional hydraulic conductivity tensor. It can also improve the accuracy of groundwater flow simulations in cases in which the model grid violates certain geometric requirements. Thus, the XT3D option is an alternative to the Ghost-Node Correction (GNC) Package, which was developed for MODFLOW-USG. n In addition to the many new or redesigned capabilities, the MODFLOW 6 input structure has also been redesigned. Within package input files, information is divided into blocks and informative keywords are used to label numeric data and activate options. This new input structure was designed to make it easier for users to adjust simulation options in an intuitive manner, reduce user input errors, and allow new capabilities to be added without causing problems with backward compatibility. n This first MODFLOW 6 release with the GWF Model synthesizes many of the capabilities in existing MODFLOW variants: n 1. MODFLOW-2005 the GWF Model contains revisions of the commonly used flow packages, stress packages, and advanced stress packages n 2. MODFLOW-NWT the GWF Model supports an optional Newton-Raphson approach for water table aquifers n 3. MODFLOW-USG GWF Models can be developed using regular MODFLOW grids or unstructured grids n 4. MODFLOW-LGR any number of GWF Models can be specified for a single simulation; these GWF Models are tightly coupled at the matrix level n 5. MODFLOW-CDSS multiple stress and advanced stress packages of the same type can be specified for a single GWF Model n For MODFLOW users with existing models, the MODFLOW 6 distribution includes a conversion program that will translate a MODFLOW-2005, MODFLOW-NWT, or MODFLOW-LGR2 model into the MODFLOW 6 format.

Potential effects of alterations to the hydrologic system on the distribution of salinity in the Biscayne aquifer in Broward County, Florida

To address concerns about the effects of water-resource management practices and rising sea level on saltwater intrusion, the U.S. Geological Survey in cooperation with the Broward County Environmental Planning and Community Resilience Division, initiated a study to examine causes of saltwater intrusion and predict the effects of future alterations to the hydrologic system on salinity distribution in eastern Broward County, Florida. A three-dimensional, variable-density solute-transport model was calibrated to conditions from 1970 to 2012, the period for which data are most complete and reliable, and was used to simulate historical conditions from 1950 to 2012. These types of models are typically difficult to calibrate by matching to observed groundwater salinities because of spatial variability in aquifer properties that are unknown, and natural and anthropogenic processes that are complex and unknown; therefore, the primary goal was to reproduce major trends and locally generalized distributions of salinity in the Biscayne aquifer. The methods used in this study are relatively new, and results will provide transferable techniques for protecting groundwater resources and maximizing groundwater availability in coastal areas. The model was used to (1) evaluate the sensitivity of the salinity distribution in groundwater to sea-level rise and groundwater pumping, and (2) simulate the potential effects of increases in pumping, variable rates of sea-level rise, movement of a salinity control structure, and use of drainage recharge wells on the future distribution of salinity in the aquifer. Results from the simulation of historical conditions indicate that the model generally represents the observed greater westward extent of elevated salinity in the central part of the intruded area relative to the northern and southernmost parts of the intruded area. Results of sensitivity testing indicate that the extent of elevated salinity is most sensitive to pumping in areas where the source of saltwater is largely offshore, from the Atlantic Ocean, and is most sensitive to sea-level rise in areas where the source of salinity is downward leakage of brackish water from canals. Simulations of future scenarios indicate that increases in pumping near the existing interface may cause the interface to advance and decreases in pumping may cause it to retreat. Climatic effects, such as periods of prolonged drought or high precipitation, may augment or counteract long-term effects of changes in pumping on aquifer salinity at well fields. With increasing rates of sea-level rise, the freshwater-saltwater interface advances progressively inland, and flow-averaged salinities at well fields near the existing interface increase commensurately. Hypothetical southeastward (downstream) re-positioning of the existing G–54 salinity-control structure may prevent the interface from moving northwestward along and near the North New River canal, but beneficial effects are localized. Implementation of freshwater recharge wells in the city of Hallandale Beach may also have only a localized freshening effect in the aquifer and little appreciable effect on the freshwater-saltwater interface or on concentrations of salinity at well fields. Model accuracy and use are limited by uncertainty in the physical properties and boundary conditions of the system, uncertainty in historical and future conditions, and generali zations made in the mathematical relationships used to describe the physical processes of groundwater flow and transport. Because of these limitations, model results should be considered in relative rather than absolute terms. Nonetheless, model results do provide useful information on the relative scale of response of the system to changes in pumping distribution, sea-level rise, and mitigation activities. 2 Potential Effects of Alterations to the Hydrologic System on the Distribution of Salinity in the Biscayne Aquifer Introduction Saltwater intrusion of the Biscayne aquifer in Broward County, Florida (fig. 1), is a challenge for water-supply management. This issue is expected to persist, given current predictions of climate change, sea-level rise, and continued population growth. The Biscayne aquifer in Broward County is particularly susceptible to saltwater intrusion because of its high permeability, low hydraulic gradient, and proximity to an unlimited source of saltwater in the Atlantic Ocean. Additional factors that contribute to saltwater intrusion include groundwater pumping, rising sea level, anthropogenic alterations to the hydraulic system, and natural climatic variations. The Biscayne aquifer is a shallow, unconfined to semiconfined aquifer composed of highly transmissive limestone situated on top of a large, flat carbonate platform that composes the Florida peninsula. In Broward County, estimated transmissivity values in the surficial aquifer system, of which the Biscayne aquifer is the most transmissive unit, are as high as 900,000 square feet per day (ft2/d; Fish, 1988). The groundwater and surface-water systems are closely connected, and the groundwater system reacts quickly and markedly to precipitation. Average annual precipitation between 1981 and 2010 was about 58 inches (in.), and generally ranged from about 2 inches per month (in/mo) in December to almost 9 in/mo in June. Land surface elevations are typically less than 20 feet (ft) in eastern Broward County (appendix fig. 1–3), and topographic relief is minimal. The Biscayne aquifer is in contact with the Atlantic Ocean offshore. The onshore potentiometric head gradient drives fresh groundwater toward the Atlantic, whereas differences in fluid density force the denser saltwater from the Atlantic Ocean into the deeper parts of the aquifer, beneath overlying freshwater-bearing parts of the aquifer (fig. 2). The Biscayne aquifer is the primary source of water supply in Broward County, and groundwater withdrawals from the aquifer have steadily increased since the early 1900s. Groundwater extraction via pumping lowers groundwater levels in the aquifer, increasingly allowing saltwater to intrude the aquifer from the east. Groundwater pumping has been cited as a substantial cause of saltwater intrusion at several locations along the Atlantic coast (Lacombe and Carleton, 2002; Monti and others, 2009; Payne, 2010; Langevin and Zygnerski, 2013). Rising sea level in the Atlantic Ocean reduces the eastward potentiometric gradient in the Biscayne aquifer, thereby elevating groundwater levels between the recharge areas and the coast. The reduced onshore potentiometric gradient allows saltwater to enter the aquifer more readily. Local mean sea level has risen approximately 0.5 ft since 1950 (fig. 3). Numerous studies indicate that sea level in the North Atlantic has been rising for thousands of years; for example, Fairbanks (1989) indicates it has risen hundreds of feet in the past 18,000 years, and Kemp and others (2011) indicates it has risen nearly 10 ft in the past 2,500 years. The local hydrologic system has been substantially altered to allow urban development in southeastern Florida, including the eastern part of Broward County. An extensive network of canals has been constructed over decades to promote drainage and route excess water to the Atlantic Ocean during major precipitation events. Drainage tends to lower the water table and reduces the seaward movement of fresh groundwater. During the dry season (October to May), the canals are used to control saltwater intrusion by providing a source of freshwater that leaks into the Biscayne aquifer and maintains water-table elevations in the aquifer. Broward County faces the possibility of increased saltwater intrusion. Population growth in the county is expected to continue, potentially increasing groundwater pumping (fig. 4), reducing groundwater levels, and increasing saltwater inflow to the aquifer. Sea-level rise has been predicted to continue, potentially reducing the eastward potentiometric gradient in the Biscayne aquifer and in tidal canals, and driving additional saltwater into the aquifer. In response to rising sea level, alterations to surface-water management are being considered to adapt to increased coastal flooding during major storm events. Such alterations are likely to affect salinity distribution in the aquifer, given the strong connection between the surface-water and groundwater systems. To address these concerns, the U.S. Geological Survey (USGS) in cooperation with the Broward County Environmental Planning and Community Resilience Division, initiated a study to examine causes of saltwater intrusion and predict the potential effects of future alterations to the hydrologic system on salinity distribution in eastern Broward County. Purpose and Scope The purpose of this report is to (1) evaluate controls on the salinity distribution from 1950 to 2012, and (2) to simulate potential effects of possible future alterations to the hydrologic system on salinity distribution in the Biscayne aquifer in the southern and central parts of eastern Broward County. The simulated area extends from the area just east of the Atlantic coastline to the area just west of the major canals separating the Everglades from the urbanized part of the county, and from the C–14 basin in the north to the C–9 east and west basins in the south (fig. 1). Historical and future hydrologic conditions affecting saltwater intrusion in the northeastern part of the county were studied previously (Langevin and Zygnerski, 2013). The study period represents the timeframe during which the greatest alterations to the hydrologic system were made and for which records are available to quantify the alterations. The study period also extends 50 years into the future, during which time proposed hypothetical changes to management of the hydrologic system might occur. Of particular interest are production well-field areas that potentially will be most affected by