Search

Simulating the interaction of groundwater with surface water networks using traditional boundary packages available with MODFLOW‐USG can be challenging for complex systems. Often several package types are required as they are typically purpose built. Moreover, these packages generally do not interact with one another which complicates accounting of groundwater discharge at different points within the system. Here, we demonstrate that the connected linear network (CLN) package of MODFLOW‐USG, and advances therein in USG‐Transport, can be used to simulate groundwater interaction with a complex surface water network comprised of creeks, ponds, wetlands, and springs, in a manner that is comparable with these other packages, but with additional benefits, including explicit routing of water between the features. Article impact statement: We demonstrate the Connected Linear Network package can be used to simulate the exchange of groundwater with complex surface water networks. [ABSTRACT FROM AUTHOR]

Numerical modeling of karst water systems is a challenge due to vastly different flow regimes in the solution channels and the porous medium, and because of the complex interactions between them. Recent developments with unstructured grid models have enabled robust simulation of such systems due to full coupling of the system of equations and the interactions, within a flexible gridding framework that allows these domains to be discretized independently of each other. These capabilities have been incorporated into the MODFLOW-USG code, which has been developed and released in the public domain by the United States Geological Survey. Further developments to MODFLOW-USG have continued and have also been released in the public domain. The significant features that enhance simulations in karst environments include various optional turbulent flow formulations, solute transport capability and a density-dependent flow simulation option that can be used to model transport of nonaqueous phase liquids (NAPLs) or freshwater–seawater interactions, for example. A dual-porosity-flow and dual-porosity-transport formulation is also included to further expand model scale and flexibility. A physically based representation of flow and transport within the karst conduits (channels) and porous medium eliminates the need for various surrogate modeling approaches in karst, thus enabling solutions to some of the common problems that groundwater professionals face when dealing with groundwater resources management in coastal karst areas. A demonstrative example shows the impact of karst on the saltwater system in a coastal aquifer. [ABSTRACT FROM AUTHOR]

Density‐dependent flow and transport solutions for coastal saltwater intrusion investigations, analyses of fluid injection into deep brines, and studies of convective fingering and instabilities of denser fluids moving through less dense fluids typically formulate the groundwater flow equation in terms of pressure or equivalent freshwater head. A formulation of the flow equation in terms of hydraulic head is presented here as an alternative. The hydraulic‐head formulation can facilitate adaptation of existing constant‐density groundwater flow codes to include density‐driven flow by avoiding the need to convert between freshwater head and hydraulic head within the code and by incorporating density‐dependent terms as a compartmentalized "correction" to constant‐density calculations already performed by the code. The hydraulic‐head formulation also accommodates complexities such as unconfined groundwater flow and Newton‐Raphson solution schemes more readily than the freshwater‐head formulation. Simulation results are presented for four example problems solved using an implementation of the hydraulic‐head formulation in MODFLOW. [ABSTRACT FROM AUTHOR]

The performance of time series models is assessed using synthetic head series simulated with a numerical model that solves Richards' equation for variably saturated flow. Heads were simulated in a homogeneous unconfined aquifer between two parallel canals; measured daily precipitation and potential evaporation are specified at the land surface and root water uptake is simulated. The head response to a precipitation event is nonlinear and depends on the saturation degree and rainfall before and after the precipitation event while evaporation reduction occurs during summers. Synthetic series were generated for 27 years and three different soil types; the unsaturated zone thickness varies between 0 and >5 m. The synthetic head series were simulated with a linear and nonlinear time series model. Performance of a linear time series model with four parameters, using a scaled Gamma response, gave R2 values ranging from 0.67 to 0.96. The nonlinear time series model with five parameters simulates recharge using a root zone reservoir after which the head response to recharge is simulated with a scaled Gamma response function. The nonlinear time series model was able to simulate all synthetic head series very well with R2 values above 0.9 for almost all models. The head response of the nonlinear model to a single precipitation event compares well to the response of the variably saturated groundwater model. The provided scripts may be used to simulate synthetic head series for other climates or for systems with additional complexity to assess the performance of other data‐driven models. [ABSTRACT FROM AUTHOR]

Abstract: A groundwater flow model is typically used to provide the flow field for conducting groundwater solute transport simulations. The advection term of the mass conserved formulation for groundwater transport assumes that the flow field is perfectly balanced and that all water flowing into a numerical grid cell is exactly balanced by outflows after accounting for sources/sinks or internal storage. However, in many complicated regional or site‐scale models, there may be localized flow balance errors that may be difficult to eliminate through tighter flow convergence tolerances due to simulation time constraints or numerical limits on convergence tolerances. Thus, if water is erroneously gained or lost within a grid cell during the flow computation, the solutes within it will also be numerically affected in the associated transport simulation. Transport solutions neglect this error in groundwater flow as the transport equations that are solved assume no error in flow. This flow imbalance error can however have consequences on the transport solution ranging from unnoticeable errors in the resulting concentrations to spurious oscillations that can grow in time and hinder further solution. An approach has been suggested here, to explicitly handle these flow imbalances during mass conserved advective transport computations and report them in the corresponding transport mass balance output, as corrections that are needed to handle errors originating in the flow solution. Example problems are provided to explain the concepts and demonstrate the impacts. [ABSTRACT FROM AUTHOR]

Guha, Hillol and Sorab Panday, 2012. Impact of Sea Level Rise on Groundwater Salinity in a Coastal Community of South Florida. Journal of the American Water Resources Association (JAWRA) 48(3): 510-529. DOI: 10.1111/j.1752-1688.2011.00630.x Abstract: Freshwater resources of coastal communities in the United States and world over are threatened by the rate of sea level rise. According to recent estimates by various governmental agencies and climate researchers, the global sea level rise is likely to be between 0.6 and 2.1 m by the year 2100. South Florida is a coastal community and much of its coastline is subject to sea level rise and potential impacts to wetlands and the water resources of the area. To understand what the impact of sea level rise would cause to the groundwater level and salinity intrusion, an integrated groundwater and surface water model was developed for North Miami-Dade and Broward Counties of South Florida. The model was calibrated against daily groundwater heads, base flows in canals, and chloride concentrations for a period of one year and six months. Three separate sensitivity analyses were conducted by increasing the sea level by 0.6, 0.9, and 1.22 m. Results of the simulations shows increase of groundwater heads in some areas from 4 to 15%; whereas the average relative chloride concentrations increased significantly by 100-600% in some wells. The increase in groundwater elevations and chloride concentrations varies from location of the wells and its proximity to the coast. The model results indicate that even a 0.6 m increase in sea level (which is the conservative estimate) is likely to impair the vital freshwater resources in many parts of South Florida. [ABSTRACT FROM AUTHOR]

Abstract: As an alternative to grid refinement, the concept of a ghost node, which was developed for nested grid applications, has been extended towards improving sub-grid scale accuracy of flow to conduits, wells, rivers or other boundary features that interact with a finite-difference groundwater flow model. The formulation is presented for correcting the regular finite-difference groundwater flow equations for confined and unconfined cases, with or without Newton Raphson linearization of the nonlinearities, to include the Ghost Node Correction (GNC) for location displacement. The correction may be applied on the right-hand side vector for a symmetric finite-difference Picard implementation, or on the left-hand side matrix for an implicit but asymmetric implementation. The finite-difference matrix connectivity structure may be maintained for an implicit implementation by only selecting contributing nodes that are a part of the finite-difference connectivity. Proof of concept example problems are provided to demonstrate the improved accuracy that may be achieved through sub-grid scale corrections using the GNC schemes. [Copyright &y& Elsevier]

The article discusses the uncertainty on the future of groundwater modeling. It explores works on the future of hydrology, hydrogeology, and trends in groundwater modeling software development. It asserts that the future of groundwater modeling will depend largely on the types of problems relative to water supply and water quality. It also notes the use of models in connection with optimization methods to better manage resources.

A physically-based, spatially-distributed model is presented for simulation of surface/subsurface flow and the interactions between these domains. The model is designed for practical application to a wide variety of hydrologic evaluations, at various scales of simulation. The system is represented by the three-dimensional saturated–unsaturated flow equation for the subsurface, coupled with the diffusion wave equation for areal overland flow, both of which are coupled with the diffusion wave equation for flow through a network of streams and channels, including hydraulic structures. Ground surface unevenness at the grid scale is incorporated via the concept of detention storage, and thick vegetation or urban features are included via an obstruction storage exclusion term. Evapotranspiration from the surface and subsurface are modeled using land cover and climatic factors to define the complete water budget using a physically-based formulation. The system of equations is discretized using a fully implicit procedure, with the Newton–Raphson method to handle non-linearities efficiently. Robustness, stability and accuracy of solution are obtained for a wide variety of cases including dry systems and large surface/subsurface interaction fluxes. Adaptive time-stepping schemes and under-relaxation formulas further alleviate the computational burden. Verification and application examples demonstrate the need for a rigorous, fully-coupled solution to the set of equations, for complete hydrologic-cycle analysis. [Copyright &y& Elsevier]

A sub-time stepping method is described for computational fluid dynamics problems that utilize implicit-type time marching procedures to resolve transients. In this method, small time-step sizes are used in portions of a domain where interest and activity is high, with larger time-step sizes being applied in other locations. The sub-time step is an integral portion of a larger time step––i.e., multiple sub-time steps over a sub-timed part of a domain add up to the time interval of the full time step used over the remainder of the domain. The technique is particularly suitable for extensive simulations where large portions of a domain are temporally over-discretized. The principles underlying implementation of the implicit sub-timing procedure, the computational effort in relation to conventional implicit time-stepping methods, and an analysis of the effects of sub-timing on the matrix structure are presented. Feasibility and applicability of the implicit sub-timing method is demonstrated through three proof-of-the-concept example problems. These examples include implementation of implicit sub-timing to one-dimensional as well as two-dimensional flow and solute transport problems. [Copyright &y& Elsevier]

Construction and maintenance of engineering structures in arctic regions demands an understanding of and the ability to cope with environmental problems produced by permafrost. Where thawing of the permafrost cannot be prevented, the range, amount, and extent of degradation are factors that must be taken into account in the design of engineering structures. Predictions of these factors are facilitated by the use of a mathematical model designed to predict the physical system under consideration. This study addresses the problem by developing a predictive model to describe permafrost thaw subsidence based upon a complete formulation of the problem. This problem is one of multiphase (air, ice, water, and solid) transport and deformation in a porous medium. The conservation of mass equations for liquid water, ice with phase change, and deforming soil solids are developed. Darcy's law is extended to include moisture movement due to thermal gradients. The quasi-static equilibrium equations and stress-strain relations that assume a perfectly elastic solid matrix are also employed. Simultaneous variations of moisture retention and phase-composition curves with temperature and pore pressure are also incorporated in the model. The energy-conservation equation, which includes terms due to viscous dissipation and compression effects, along with the other equations expressed, is found to satisfactorily represent the physical processes occurring in thawing permafrost soils. Approximations to simulate the thawing of a one-dimensional soil column are then developed from the general model. [ABSTRACT FROM AUTHOR]

In this paper, we develop a numerical solution to equations derived in the first part of this study to simulate the thaw subsidence of an unsaturated frozen soil due to surface heating. Solutions obtained in terms of pore pressure, temperature, degree of unfrozen water content, ice content, and porosity provide the spatial and temporal variations of these quantities in a thawing soil column. A one-dimensional form of governing equations is discretized by employing a Crank-Nicolson finite-difference scheme. These highly nonlinear equations are then solved by an iterative procedure at every time step for predescribed positive surface temperatures that produce thawing. The numerical solutions are obtained for a 1-m silt column. Results show that, above the frozen region, unfrozen water saturation increases due to low hydraulic conductivity, which prevents the melting ice from draining, thus causing accumulation of water above it. It has also been observed that most of the settlement at the surface is due to thaw consolidation rather than settlement due to drainage. The former is completed in seven days of a 12-day simulation period. Furthermore, a sensitivity analysis was carried out to determine the response of the simulation to deviations in the values of certain parameters. It is found that the settlement estimations are very sensitive to the thaw-settlement parameter and phase-composition curve, thus requiring accurate determination of these parameters for a reliable simulation of the thaw-consolidation process. [ABSTRACT FROM AUTHOR]

Managed aquifer recharge has become a standard water resources management practice to promote the development of locally sustainable water supplies and combat water scarcity. However, installation of injection wells for replenishment purposes in urban areas with complex hydrogeology faces many challenges, such as limited land availability, potential impacts on municipal production wells and known subsurface contamination plumes, and complex spatially variable hydraulic connections between aquifer units. To assess the feasibility and cost‐effectiveness of injecting advanced treated water (ATW) into a complex urban aquifer system, a Simulation‐Optimization (SO) model was developed to automate a systematic search for the most cost‐effective locations to install new wells for injecting various quantities of ATW, if feasible. The generalized workflow presented here uses an existing MODFLOW groundwater model—along with advanced optimization routines that are publicly available—to flexibly accommodate a multiobjective function, complex constraints, and specific project requirements. The model successfully placed wells for injection of 1 to 4 MGD of ATW in aquifers underlying the study area. The injection well placement was primarily constrained by avoiding excessive impact on environmental sites with underlying groundwater plumes. The largest costs were for well installation and piping to the wells from the existing ATW pipes. This workflow is readily adaptable to other sites with different complexities, decision variables, or constraints. [ABSTRACT FROM AUTHOR]

MODFLOW 6 is the latest in a line of six "core" versions of MODFLOW released by the U.S. Geological Survey. The MODFLOW 6 architecture supports incorporation of additional hydrologic processes, in addition to groundwater flow, and allows interaction between processes. The architecture supports multiple model instances and multiple types of models within a single simulation, a flexible approach to formulating and solving the equations that represent hydrologic processes, and recent advances in interoperability, which allow MODFLOW to be accessed and controlled by external programs. The present version of MODFLOW 6 consolidates popular capabilities available in MODFLOW variants, such as the unstructured grid support in MODFLOW‐USG, the Newton‐Raphson formulation in MODFLOW‐NWT, and the support for partitioned stress boundaries in MODFLOW‐CDSS. The flexible multi‐model capability allows users to configure MODFLOW 6 simulations to represent the local‐grid refinement (LGR) capabilities available in MODFLOW‐LGR, the multi‐species transport capabilities in MT3DMS, and the coupled variable‐density capabilities available in SEAWAT. This paper provides a new, holistic and integrated overview of simulation capabilities made possible by the MODFLOW 6 architecture, and describes how ongoing and future development can take advantage of the program architecture to integrate new capabilities in a way that is minimally invasive and automatically compatible with the existing MODFLOW 6 code. [ABSTRACT FROM AUTHOR]

Axisymmetric groundwater models are used for simulating radially symmetric conditions. Groundwater simulators built specifically to model axisymmetric conditions are most commonly used for simulating aquifer tests. Although some numerical models capable of simulating flow and solute transport that are developed in the cartesian coordinate system framework offer flexibility to simulate axisymmetric conditions, most of the numerical groundwater models, such as the MODFLOW family of codes, are based on structured grids in which axisymmetric flows cannot be directly simulated. Researchers in the past have provided methods to manipulate aquifer properties to mimic axisymmetric conditions. This study presents a methodology that takes advantage of the unstructured grids of MODFLOW‐USG to simulate axisymmetric models within the MODFLOW framework. To develop axisymmetric models, the intercell interface area arrays of MODFLOW‐USG were calculated to accurately represent coaxial cylindrical model cells. Three examples are presented to demonstrate the application of MODFLOW‐USG for axisymmetric modeling: a pumping well with delayed yield effects, a vadose zone flow model simulating an infiltration basin, and a density‐dependent saltwater intrusion problem for a circular island. Results were verified against analytical solutions and published numerical codes. Article impact statement: This study presents a new methodology to create axisymmetric models using MODFLOW‐USG. [ABSTRACT FROM AUTHOR]

[Display omitted] • Improvement of CLN and DRT packages of MODLOW-USG for BHE numerical simulation. • Testing of new procedure against previous numerical model and analytical solution. • Possibility to reproduce typical GSHP system and Thermal Response Test operation. • Easy and expeditious implementation of large borefield in a numerical model. In recent years, among renewable energies, the geothermal resource exploitation shows a constant growth; specifically, in countries engaged in CO 2 emissions reduction and dependent on energy from abroad, the low-temperature geothermal energy (geo-exchange) for air conditioning of buildings represents a cost-effective and green solution. In closed-loop systems borehole heat exchangers (BHE) are coupled with ground-source heat pumps (GSHP) constituting the key component of the heating ventilation air-conditioning (HVAC) system. Therefore, the design of the BHE and the correct interpretation of in situ Thermal Response Tests (TRT) are essential to supply the building energy demand. To support the design, Modflow-USG Connected Linear Network (CLN) and Drain Return Flow (DRT) packages are adapted and improved to reproduce the operation of one or more BHE in aquifers and to analyze the TRT. The improvements are compared with a previously developed numerical model and two different analytical solutions (infinite line source and moving line source) by imposing a constant heat rate injection into the aquifer. The results show good agreement between the new approach and previous ones (discrepancy lower than 2% for models with highly refined grid), but the new approach is much more accurate and expeditious in both implementation and execution, also allowing for an easy numerical simulation of multiple BHE. [ABSTRACT FROM AUTHOR]

[Display omitted] • Improvement of CLN and DRT packages of MODLOW-USG for BHE numerical simulation. • Testing of new procedure against previous numerical model and analytical solution. • Possibility to reproduce typical GSHP system and Thermal Response Test operation. • Easy and expeditious implementation of large borefield in a numerical model. In recent years, among renewable energies, the geothermal resource exploitation shows a constant growth; specifically, in countries engaged in CO 2 emissions reduction and dependent on energy from abroad, the low-temperature geothermal energy (geo-exchange) for air conditioning of buildings represents a cost-effective and green solution. In closed-loop systems borehole heat exchangers (BHE) are coupled with ground-source heat pumps (GSHP) constituting the key component of the heating ventilation air-conditioning (HVAC) system. Therefore, the design of the BHE and the correct interpretation of in situ Thermal Response Tests (TRT) are essential to supply the building energy demand. To support the design, Modflow-USG Connected Linear Network (CLN) and Drain Return Flow (DRT) packages are adapted and improved to reproduce the operation of one or more BHE in aquifers and to analyze the TRT. The improvements are compared with a previously developed numerical model and two different analytical solutions (infinite line source and moving line source) by imposing a constant heat rate injection into the aquifer. The results show good agreement between the new approach and previous ones (discrepancy lower than 2% for models with highly refined grid), but the new approach is much more accurate and expeditious in both implementation and execution, also allowing for an easy numerical simulation of multiple BHE. [ABSTRACT FROM AUTHOR]

The χMD matrix solver package is incorporated into USGS groundwater modeling software, such as MODFLOW-NWT, MODFLOW-USG, and MT3D. The solver is used to solve matrices assembled through numerical discretization of the groundwater flow equation, and solute transport equations. χMD has demonstrated its higher robustness, faster execution speed, and more efficient memory usage compared to the existing solvers for many types of groundwater flow problems. χMD uses preconditioned iterative Krylov-subspace methods and consists of preconditioning and acceleration modules. Because the solver package uses a variety of preconditioning features including level-based incomplete lower-upper (ILU) factorization method with a drop tolerance scheme, users must choose optimal preconditioning parameters to improve execution speed and robustness. In order to examine how the preconditioning parameters, ILU factorization level, and drop tolerance values affect the overall performance of the matrix solver, we evaluated five different groundwater model applications using MODFLOWUSG that include different numerical complexities. For those five cases, the number of discretization nodes varied from 10,000 cells to 730,300 cells. From the analysis, we found that the preconditioning parameters greatly affect execution times and memory usage of the preconditioning and acceleration procedures. In addition, a combination of the ILU level between five to seven and the drop tolerance value between 10-2 and 10-3 usually resulted in shorter overall execution time. Our study suggests that the users can elicit higher performance and robustness of the χMD matrix solver using this combination of the parameters and enhance computational efficiency of solving groundwater and solute transport problems. [ABSTRACT FROM AUTHOR]

Various approaches have been proposed to manage the nonlinearities associated with the unconfined flow equation and to simulate perched groundwater conditions using the MODFLOW family of codes. The approaches comprise a variety of numerical techniques to prevent dry cells from becoming inactive and to achieve a stable solution focused on formulations of the unconfined, partially-saturated, groundwater flow equation. Keeping dry cells active avoids a discontinuous head solution which in turn improves the effectiveness of parameter estimation software that relies on continuous derivatives. Most approaches implement an upstream weighting of intercell conductance and Newton-Raphson linearization to obtain robust convergence. In this study, several published approaches were implemented in a stepwise manner into MODFLOW for comparative analysis. First, a comparative analysis of the methods is presented using synthetic examples that create convergence issues or difficulty in handling perched conditions with the more common dry-cell simulation capabilities of MODFLOW. Next, a field-scale three-dimensional simulation is presented to examine the stability and performance of the discussed approaches in larger, practical, simulation settings. [ABSTRACT FROM AUTHOR]

Groundwater discharge through evaporation due to a shallow water table can be an important component of a regional scale water balance. Modeling this phenomenon in irrigated regions where soil moisture varies on short time scales is most accurately accomplished using variably saturated modeling codes. However, the computational demands of these models limit their application to field scale problems. The MODFLOW groundwater modeling code is applicable to regional scale problems and it has an evapotranspiration package that can be used to estimate this form of discharge, however, the use of time-invariant parameters in this module result in evaporation rates that are a function of water table depth only. This paper presents a calibration and validation of the previously developed MOD-HMS model code using lysimeter data. The model is then used to illustrate the dependence of bare soil evaporation rates on water table depth and soil moisture conditions. Finally, an approach for estimating the time varying parameters for the MODFLOW evapotranspiration package using a 1-D variably saturated MOD-HMS model is presented. [ABSTRACT FROM AUTHOR]

The sustainability of irrigated agriculture in many arid and semiarid areas of the world is at risk because of a combination of several interrelated factors, including lack of fresh water, lack of drainage, the presence of high water tables, and salinization of soil and groundwater resources. Nowhere in the United States are these issues more apparent than in the San Joaquin Valley of California. A solid understanding of salinization processes at regional spatial and decadal time scales is required to evaluate the sustainability of irrigated agriculture. A hydro-salinity model was developed to integrate subsurface hydrology with reactive salt transport for a 1,400-km² study area in the San Joaquin Valley. The model was used to reconstruct historical changes in salt storage by irrigated agriculture over the past 60 years. We show that patterns in soil and groundwater salinity were caused by spatial variations in soil hydrology, the change from local groundwater to snowmelt water as the main irrigation water supply, and by occasional droughts. Gypsum dissolution was a critical component of the regional salt balance. Although results show that the total salt input and output were about equal for the past 20 years, the model also predicts salinization of the deeper aquifers, thereby questioning the sustainability of irrigated agriculture. [ABSTRACT FROM AUTHOR]

[Display omitted] • Finite difference code used to assess the aquifer thermal perturbation of a GSHP. • CLN package of open source Modflow-USG adapted to simulate BHE operation in aquifer. • Validation of the BHE model using CLN package coupled to unstructured grid. • Possible implementation of Borefield in one numerical model by means of Modflow-USG. Simulating heat transfer in an aquifer with one or more vertical Borehole Heat Exchangers (BHEs) of a Ground Source Heat Pump (GSHP) system by means of a finite difference code is difficult because of the square or rectangular geometry grid and computational times thus limiting the types of evaluations that can be performed. The aim of this work is to explore through MODFLOW-USG code (public domain software) a different approach towards simulating a borefield that would be more efficient computationally, in order to enable simulations of larger domains with multiple BHEs. The Connected Linear Network (CLN) package, introduced in MODFLOW-USG, generally simulates 1-D linear computational cells in a 3-D grid, such as hydraulic pipes in subsoil, but for the first time has been adapted to reproduce vertical closed loop U-pipe of a BHE. Therefore, this work evaluates the MODFLOW-USG and CLN package capability to reproduce the yearly operation of one or more BHEs in an aquifer as a simpler and faster approach compared to a very fine finite-difference discretization. Once the CLN package was adapted, a sensitivity analysis on the grid size refinement was performed. There were several findings from this work. The results of the different numerical models were in good agreement with an already validated model, in terms of exchanged energies and aquifer thermal perturbation. Same analyses were carried out for different groundwater flow velocities and it was confirmed that the exchanged energy by a BHE increases with the groundwater flow velocity in accordance with literature studies. At last, a borefield of 7 BHEs was implemented in a numerical model in a more expeditious and efficient way and without any computational effort. [ABSTRACT FROM AUTHOR]