Grid cell distortion and MODFLOW's integrated finite-difference numerical solution

The ground water flow model MODFLOW inherently implements a nongeneralized integrated finite-difference (IFD) numerical scheme. The IFD numerical scheme allows for construction of finite-difference model grids with curvilinear (piecewise linear) rows. The resulting grid comprises model cells in the shape of trapezoids and is distorted in comparison to a traditional MODFLOW finite-difference grid. A version of MODFLOW-88 (herein referred to as MODFLOW IFD) with the code adapted to make the one-dimensional DELR and DELC arrays two dimensional, so that equivalent conductance between distorted grid cells can be calculated, is described. MODFLOW IFD is used to inspect the sensitivity of the numerical head and velocity solutions to the level of distortion in trapezoidal grid cells within a converging radial flow domain. A test problem designed for the analysis implements a grid oriented such that flow is parallel to columns with converging widths. The sensitivity analysis demonstrates MODFLOW IFD's capacity to numerically derive a head solution and resulting intercell volumetric flow when the internal calculation of equivalent conductance accounts for the distortion of the grid cells. The sensitivity of the velocity solution to grid cell distortion indicates criteria for distorted grid design. In the radial flow test problem described, the numerical head solution is not sensitive to grid cell distortion. The accuracy of the velocity solution is sensitive to cell distortion with error <1% if the angle between the nonparallel sides of trapezoidal cells is <12.5 degrees. The error of the velocity solution is related to the degree to which the spatial discretization of a curve is approximated with piecewise linear segments. Curvilinear finite-difference grid construction adds versatility to spatial discretization of the flow domain. MODFLOW-88's inherent IFD numerical scheme and the test problem results imply that more recent versions of MODFLOW 2000, with minor modifications, have the potential to make use of a curvilinear grid.     

An efficient tool for accelerating the numerical solution of the stochastic subsurface flow problem using neural networks

 An efficient numerical solution for the two-dimensional groundwater flow problem using artificial neural networks (ANNs) is presented. Under stationary velocity conditions with unidirectional mean flow, the conductivity realizations and the head gradients, obtained by a traditional finite difference solution to the flow equation, are given as input-output pairs to train a neural network. The ANN is trained successfully and a certain level of recognition of the relationship between input conductivity patterns and output head gradients is achieved. The trained network produced velocity realizations that are physically plausible without solving the flow equation for each of the conductivity realizations. This is achieved in a small fraction of the time necessary for solving the flow equations. The prediction accuracy of the ANN reaches 97.5% for the longitudinal head gradient and 94.7% for the transverse gradient. Head-gradient and velocity statistics in terms of the first two moments are obtained with a very high accuracy. The cross covariances between head gradients and the fluctuating log-conductivity (log-K) and between velocity and log-K obtained with the ANN approach match very closely those obtained by a traditional numerical solution. The same is true for the velocity components auto-covariances. The results are also extended to transport simulations with very good accuracy. Spatial moments (up to the fourth) of mean-concentration plumes obtained using ANNs are in very good agreement with the traditional Monte Carlo simulations. Furthermore, the concentration second moment (concentration variance) is very close between the two approaches. Considering the fact that higher moments of concentration need more computational effort in numerical simulations, the advantage of the presented approach in saving long computational times is evident. Another advantage of the ANNs approach is the ability to generalize a trained network to conductivity distributions different from those used in training. However, the accuracy of the approach in cases with higher conductivity variances is being investigated.     

Bromide transport under sprinkler and flood irrigation for no-till soil condition

Understanding the influence of irrigation methods on solute transport is essential to properly manage chemical use in agricultural soils. In this study, we compare the transport of a conservative solute (bromide) under sprinkler and flood irrigations on a sandy clay loam (mixed Ustollic Haplargid) under no-till condition. After spraying 148.8 kg/ha of KBr on the surface, ≈25 cm of irrigation water was applied in six increments over two months as flood irrigation on one plot and as sprinkler irrigation on another plot. The net applied water (NAW=irrigation+precipitation−evaporation) was similar for both plots, which allowed the comparison of the Br profiles for the two types of irrigation. Water content and Br concentration were sampled at 5, 19, 34, and 68 days after chemical application.The recovered mass of Br and the location of center of mass were comparable for the two types of irrigation. The spread around the center of mass, however, was higher for the flood-irrigated plot. On the flood-irrigated plot, more mass leached below the depth of 90 cm, with the differences being statistically significant. The velocity of the Br center of mass was consistently 10%–20% larger than the piston displacement velocity. Dispersion and velocity coefficients varied substantially between sampling time. A recent quasi-steady solution of the convection–dispersion equation [M.H. Nachabe, L.R. Ahuja, Quasi-analytical solution for predicting the redistribution of surface-applied chemicals. Trans. ASAE 39(5) (1996) 1659–1664], which accounts for variable flow and dispersion, simulates the Br profiles fairly well.     

Exact solution for two-dimensional flow to a well in an anisotropic domain

Although most current applications of the analytic element method are formulated for isotropic hydraulic conductivity, anisotropic domains can be modeled with analytic elements using the well-known coordinate transformation where one coordinate axis is scaled by the square root of the anisotropy ratio. If the standard analytic solution for steady radial flow to a well is used with this coordinate transformation, the resulting solution correctly models the far field but it does not meet the constant head boundary condition at the well radius. This could be a significant shortcoming if you are interested in the flow field close to the well or want to estimate the head at the pumping well. A new solution for two-dimensional steady flow to a well in an anisotropic domain is presented. This solution satisfies the governing equations exactly and meets the constant head boundary condition at the well radius exactly. It was derived using a conformal mapping.     

Application of multiscale finite element method in the uncertainty qualification of large-scale groundwater flow

In this article, we discuss the application of multiscale finite element method (MsFEM) to groundwater flow in heterogeneous porous media. We investigate the ability of MsFEM in qualifying the flow uncertainty. Monte Carlo simulation is employed to implement the stochastic analysis, and MsFEM is used to avoid a full resolution to the spatial variable conductivity field. Large-scale flow with high variability is investigated by inspecting the single realization as well as the probability distribution functions of head and velocity. The numerical results show that the performance of MsFEM depends on the ratio between the correlation length and the coarse element size. An accurate prediction to the velocity requires a much lower ratio than the head. The MsFEM has different convergence rates for the head and the velocity, while the convergence rates do not deteriorate as the variance grows.     

Wellbore flow‐rate solution for a constant‐head test in two‐zone finite confined aquifers

The solution describing the wellbore flow rate in a constant‐head test integrated with an optimization approach is commonly used to analyze observed wellbore flow‐rate data for estimating the hydrogeological parameters of low‐permeability aquifers. To our knowledge, the wellbore flow‐rate solution for the constant‐head test in a two‐zone finite‐extent confined aquifer has never been reported so far in the literature. This article is first to develop a mathematical model for describing the head distribution in the two‐zone aquifer. The Laplace domain solutions for the head distributions and wellbore flow rate in a two‐zone finite confined aquifer are derived using the Laplace transform, and their corresponding time domain solutions are then obtained using the Bromwich integral method and residue theorem. These new solutions are expressed in terms of an infinite series with Bessel functions and not straightforward to calculate numerically. A large‐time solution for the wellbore flow rate is therefore developed by employing the relationship of small Laplace variable versus large time variable and L'Hospital's rule. The result shows that the large‐time solution is identical to the steady‐state solution obtained after applying the Tauberian theorem into the Laplace domain solution. This large‐time solution can reduce to the Thiem equation in the case of no skin. Finally, the newly developed solution is used to investigate the effects of outer boundary distance and conductivity ratio on the wellbore flow rate. Copyright © 2011 John Wiley & Sons, Ltd.     

Effect of overpumping and irrigation stress on hydrochemistry and hydrodynamics of a Saharan oasis groundwater system

Stepwise hydrochemical and isotope-based methodology was adopted to identify mineralization processes, assess the impact of resources overexploitation and flood irrigation, and conceptualize groundwater hydrodynamics in the Djérid aquifer system, Tunisia. The study demonstrates that the main processes controlling groundwater geochemistry are dissolution of evaporates and phosphate-bearing rocks, cation exchange, mixing between high and low TDS end-members, and irrigation return flow. Interpretation of isotope data demonstrates that the deep aquifer was mostly recharged by late Pleistocene palaeowater, while the shallow aquifer is entirely recharged by return flow. The intermediate aquifer groundwater is actually a mixing of early to middle Holocene palaeowater, late Pleistocene deep aquifer palaeowater and return flow waters. The established conceptual model shows that deep and shallow groundwater leakages into the intermediate aquifer are enhanced by the presence of deep faults, the high hydraulic head of the deep aquifer, the overexploitation of the intermediate aquifer, and the long-term flood irrigation.     

The exponential decline in saturated hydraulic conductivity with depth: a novel method for exploring its effect on water flow paths and transit time distribution

The strong vertical gradient in soil and subsoil saturated hydraulic conductivity is characteristic feature of the hydrology of catchments. Despite the potential importance of these strong gradients, they have proven difficult to model using robust physically based schemes. This has hampered the testing of hypotheses about the implications of such vertical gradients for subsurface flow paths, residence times and transit time distribution. Here we present a general semi‐analytical solution for the simulation of 2D steady‐state saturated‐unsaturated flow in hillslopes with saturated hydraulic conductivity that declines exponentially with depth. The grid‐free solution satisfies mass balance exactly over the entire saturated and unsaturated zones. The new method provides continuous solutions for head, flow and velocity in both saturated and unsaturated zones without any interpolation process as is common in discrete numerical schemes. This solution efficiently generates flow pathlines and transit time distributions in hillslopes with the assumption of depth‐varying saturated hydraulic conductivity. The model outputs reveal the pronounced effect that changing the strength of the exponential decline in saturated hydraulic conductivity has on the flow pathlines, residence time and transit time distribution. This new steady‐state model may be useful to others for posing hypotheses about how different depth functions for hydraulic conductivity influence catchment hydrological response. Copyright © 2016 John Wiley & Sons, Ltd.     

Adaptive-grid simulation of groundwater flow in heterogeneous aquifers

The prediction of contaminant transport in porous media requires the computation of the flow velocity. This work presents a methodology for high-accuracy computation of flow in a heterogeneous isotropic formation, employing a dual-flow formulation and adaptive gridding. The dual equations, describing the hydraulic head and the streamfunction, are numerically solved through finite element approximations. The application of classic finite-element methods requires a rather large number of nodes to represent suitably the flow in high-contrast formations. We present a mesh-adaptive approach that enhances the accuracy of the numerical flow solution for a given computational effort. We rely on an a posteriori error estimator to identify areas where refinements of the finite element mesh are needed or unrefinements are acceptable. We also demonstrate through numerical experiments that the developed methodology efficiently enhances accuracy through successive mesh adaptation.     

Numerical method of moments for solute transport in a nonstationary flow field conditioned on hydraulic conductivity and head measurements

The unconditional stochastic studies on groundwater flow and solute transport in a nonstationary conductivity field show that the standard deviations of the hydraulic head and solute flux are very large in comparison with their mean values (Zhang et al. in Water Resour Res 36:2107–2120, 2000; Wu et al. in J Hydrol 275:208–228, 2003; Hu et al. in Adv Water Resour 26:513–531, 2003). In this study, we develop a numerical method of moments conditioning on measurements of hydraulic conductivity and head to reduce the variances of the head and the solute flux. A Lagrangian perturbation method is applied to develop the framework for solute transport in a nonstationary flow field. Since analytically derived moments equations are too complicated to solve analytically, a numerical finite difference method is implemented to obtain the solutions. Instead of using an unconditional conductivity field as an input to calculate groundwater velocity, we combine a geostatistical method and a method of moment for flow to conditionally simulate the distributions of head and velocity based on the measurements of hydraulic conductivity and head at some points. The developed theory is applied in several case studies to investigate the influences of the measurements of hydraulic conductivity and/or the hydraulic head on the variances of the predictive head and the solute flux in nonstationary flow fields. The study results show that the conditional calculation will significantly reduce the head variance. Since the hydraulic head measurement points are treated as the interior boundary (Dirichlet boundary) conditions, conditioning on both the hydraulic conductivity and the head measurements is much better than conditioning only on conductivity measurements for reduction of head variance. However, for solute flux, variance reduction by the conditional study is not so significant.     

Determining drag coefficients and their application in modelling of turbulent flow with submerged vegetation

Vegetation is a key aspect of water resources and ecology in natural rivers, floodplains and irrigation channels. The hydraulic resistance of the water flow is greatly changed when submerged vegetation is present. Three kinds of drag coefficients, i.e., the drag coefficient for an isolated cylinder, the bulk drag coefficient of an array of cylinders and the vertically distributed or local drag coefficient, have been commonly used as parameters to represent the vegetation drag force. In this paper, a comprehensive experimental study of submerged stems in an open channel flow is presented. Empirical formulae for the three drag coefficients were obtained based on our experimental results and on data from previous studies. A two-layer model was developed to solve the mean momentum equation, which was used to evaluate the vertical mean velocity profile with each of the drag coefficients. By comparing the velocity distribution model predictions and the measurement results, we found that the model with the drag coefficient for an isolated cylinder and the local drag coefficient was good fit. In addition, the model with the bulk drag coefficient gave much larger velocity values than measurements, but it could be improved by adding the bed friction effect and making choice of the depth-averaged velocity within the canopy layer.     

Multilayer analytic element modeling of radial collector wells

A new multilayer approach is presented for the modeling of ground water flow to radial collector wells. The approach allows for the inclusion of all aspects of the unique boundary condition along the lateral arms of a collector well, including skin effect and internal friction losses due to flow in the arms. The hydraulic conductivity may differ between horizontal layers within the aquifer, and vertical anisotropy can be taken into account. The approach is based on the multilayer analytic element method, such that regional flow and local three-dimensional detail may be simulated simultaneously and accurately within one regional model. Horizontal flow inside a layer is computed analytically, while vertical flow is approximated with a standard finite-difference scheme. Results obtained with the proposed approach compare well to results obtained with three-dimensional analytic element solutions for flow in unconfined aquifers. The presented approach may be applied to predict the yield of a collector well in a regional setting and to compute the origin and residence time, and thus the quality, of water pumped by the collector well. As an example, the addition of three lateral arms to a collector well that already has three laterals is investigated. The new arms are added at an elevation of 2 m above the existing laterals. The yield increase of the collector well is computed as a function of the lengths of the three new arms.     

Modelling combined open channel flow by artificial neural networks

Combined open channel flow is encountered in many hydraulic engineering structures and processes, such as irrigation ditches and wastewater treatment facilities. Extensive experimental studies have conducted to investigate combined flow characteristics. Nevertheless, there is no simple relationship that can fully describe the velocity profiles in a turbulent flow. The artificial neural network (ANN) has great computational capability for solving various complex problems, such as function approximation. The main objective of this study is to evaluate the applicability of the ANN for simulating velocity profiles, velocity contours and estimating the discharges accordingly. The velocity profiles measured by an acoustic doppler velocimeter in the open channel of the Chihtan purification plant, Taipei, with different discharges at fixed measuring section and different depths are presented. The total number of data sets is 640 and the data sets are split into two subsets, i.e. training and validation sets. The backpropagation algorithm is used to construct the neural network. The results demonstrate that the velocity profiles can be modelled by the ANN, and the ANN constructed can nicely fit the velocity profiles and can precisely predict the discharges for the conditions investigated. Copyright © 2005 John Wiley & Sons, Ltd.     

Feasibility of Deep Nutrients Delivery into a Prince William Sound Beach for the Bioremediation of the Exxon Valdez Oil Spill

Delivery of dissolved chemicals to potentially bioremediate oil from the Exxon Valdez oil spill was investigated at Smith Island, Prince William Sound, Alaska. Two transects for tracer application were installed: one at a clean area and another at an oiled area. Tracer delivery occurred under ambient pressure through manifolds. Lithium in LiBr/sea water solution was the inert tracer, and the solution (82.6 mg/L of lithium) was released at a flow rate of 0.23 L/min for 58.5 h. While maintaining the flow rate, the solution was switched to sea water (i.e., 0.0 mg/L of lithium) for 16 h. The results show that the tracer moved landward with rising tides and seaward with falling tides. The plume got deeper moving landward and shallower moving seaward of the manifold. Thus, in situations where oil biodegradation is limited by the availability of nutrients or oxygen, applying these chemicals by this technique would allow them to reach entrapped oil seaward of the manifold from below. The seaward plume traveling speed (around 2.0 m/d) suggests that this technique is logistically feasible from a hydraulic point of view.     

A two‐dimensional analytical solution for groundwater flow in a leaky confined aquifer system near open tidal water

Groundwater in coastal areas is commonly disturbed by tidal fluctuations. A two‐dimensional analytical solution is derived to describe the groundwater fluctuation in a leaky confined aquifer system near open tidal water under the assumption that the groundwater head in the confined aquifer fluctuates in response to sea tide whereas that of the overlying unconfined aquifer remains constant. The analytical solution presented here is an extension of the solution by Sun for two‐dimensional groundwater flow in a confined aquifer and the solution by Jiao and Tang for one‐dimensional groundwater flow in a leaky confined aquifer. The analytical solution is compared with a two‐dimensional finite difference solution. On the basis of the analytical solution, the groundwater head distribution in a leaky confined aquifer in response to tidal boundaries is examined and the influence of leakage on groundwater fluctuation is discussed. Copyright © 2001 John Wiley & Sons, Ltd.     

An analytic benchmark test for karst-aquifer flow

A benchmark test for flow in karstic aquifers is presented in the form of an exact solution of the harmonic variations of water flux and head within a karst conduit that is imbedded within a three-dimensional porous matrix having a free surface. The variations are driven by a prescribed variation of head applied at one end of the conduit. The benchmark consists of expressions for the spring discharge as a function of time and the conduit head and flux as functions of distance along the conduit and time. These expressions contain three dimensionless parameters, permitting development of a wide range of specific benchmark tests. The expressions are particularly simple in the case of an infinitely deep aquifer. This limiting solution should provide the most severe test for two-dimensional models of karst aquifer flow. Another limiting case of interest is that in which the conduit diameter is equal to the water depth. This limiting solution should provide the easiest test for two-dimensional models.     

Hydrograph peakedness and basin area

The peakedness of a basin and its variation with drainage area were analyzed for three areas. Peakedness of a basin is calculated as mean flow as a percentage of highest flow. A fitted power curve relating ‘Peakedness index’ (PI) to drainage area for each of the three areas indicates a break point in a basin of about 300 km². This break point divides the basins into small basins which are more peaky and large basins which are less peaky. The break point is an outcome of a difference in order of magnitude between channel flow velocity from the headwater sources and hillside flow velocity. When the basin responds to hillslope flow the runoff from the head water sources has already flown about 30 km downstream.     

Integrated water resource management in a major canal command in eastern India

The present rice‐dominated cropping system in the Hirakud canal command (eastern India) is under severe threat due to imbalance between irrigation water supply and demand. The canal water supply, which is the only source of irrigation, only meets 54% of the demand at 90% probability of exceedance (PE). In order to mitigate the irrigation water deficit from canal water, groundwater is considered as a supplemental source. Quasi‐three‐dimensional groundwater flow simulation modelling was, therefore, carried out by using Visual MODFLOW to detect the change in hydraulic head due to transient pumping stresses. The simulation model was calibrated and validated satisfactorily. Sensitivity analysis of the model parameters shows that groundwater recharge is most sensitive followed by aquifer hydraulic conductivity at almost all the sites of the command area, whereas the model is comparatively less sensitive to specific storage and specific yield. Enhanced pumping scenarios showed that groundwater extraction can be increased up to 50 times of the existing pumping without causing any adverse effect to the aquifer but the aquifer does not permit to exploit water in order to fulfill the irrigation water demand even at 10% PE. Hence, it is imperative to develop an optimal land and water resources management plan of the command area. Copyright © 2011 John Wiley & Sons, Ltd.     

Experimental and Numerical Investigation of the Flow Hydraulic in Gradual Transition Open Channels

Channel expansions are common in both natural and artificial open channels. With increasing cross-sectional dimensions in an expansion, the flow decelerates. Due to separation of flow and subsequent eddy formation, a significant head loss is occurred along the transition. This study presents the results of experimental investigations on subcritical flow along the expansive transition of rectangular to trapezoidal channels. Also, a numerical simulation was developed using the finite volume method with Reynolds Stress turbulent model. Water surface profiles and velocity distributions of flow through the transition were measured experimentally and compared with the numerical results. Also, hydraulic efficiency of the transition and coefficient of energy head loss were calculated. The results show that with increasing the upstream Froude number, hydraulic efficiency of the transition and coefficient of energy head loss are decreased and increased, respectively. The results also showed the ability of numerical simulation for simulating the flow separation zones and bed shear stress along the transition for different inlet discharges and inflow Froude numbers.     

One-dimensional uniform and time varying solute dispersion along transient groundwater flow in a semi-infinite aquifer

An analytical solution for the space-time variation of contaminant concentration in one-dimensional transient groundwater flow in a homogenous semi-infinite aquifer, subjected to time-dependent source contamination, is derived. The uniform and time varying dispersion along transient groundwater flow is investigated under two conditions. First, the flow velocity distribution in the aquifer is considered as a sinusoidally varying function, and second, the flow velocity distribution is treated as an exponentially increasing function of time. It is assumed that initially the aquifer is not solute free, so the initial background concentration is considered as an exponentially decreasing function of the space variable which is tending to zero at infinity. It is assumed that dispersion is directly proportional to the square of the velocity, noting that experimental observations indicate that dispersion is directly proportional to the velocity with a power ranging from 1 to 2. The analytical solution is illustrated using an example and may help benchmark numerical codes and solutions.