Testing density-dependent groundwater models: two-dimensional steady state unstable convection in infinite,finite and inclined porous layers

This study proposes the use of several problems of unstable steady state convection with variable fluid density in a porous layer of infinite horizontal extent as two-dimensional (2-D) test cases for density-dependent groundwater flow and solute transport simulators. Unlike existing density-dependent model benchmarks, these problems have well-defined stability criteria that are determined analytically. These analytical stability indicators can be compared with numerical model results to test the ability of a code to accurately simulate buoyancy driven flow and diffusion. The basic analytical solution is for a horizontally infinite fluid-filled porous layer in which fluid density decreases with depth. The proposed test problems include unstable convection in an infinite horizontal box, in a finite horizontal box, and in an infinite inclined box. A dimensionless Rayleigh number incorporating properties of the fluid and the porous media determines the stability of the layer in each case. Testing the ability of numerical codes to match both the critical Rayleigh number at which convection occurs and the wavelength of convection cells is an addition to the benchmark problems currently in use. The proposed test problems are modelled in 2-D using the SUTRA [SUTRA––A model for saturated–unsaturated variable-density ground-water flow with solute or energy transport. US Geological Survey Water-Resources Investigations Report, 02-4231, 2002. 250 p] density-dependent groundwater flow and solute transport code. For the case of an infinite horizontal box, SUTRA results show a distinct change from stable to unstable behaviour around the theoretical critical Rayleigh number of 4π² and the simulated wavelength of unstable convection agrees with that predicted by the analytical solution. The effects of finite layer aspect ratio and inclination on stability indicators are also tested and numerical results are in excellent agreement with theoretical stability criteria and with numerical results previously reported in traditional fluid mechanics literature.     

Thermohaline convection in flowing groundwater

Geothermal activity creates destabilising temperature gradients which are significant in some aquifers. Usually, in such aquifers stabilising salinity gradients also exist. The combination of temperature and salinity distribution in the aquifer may induce various types of hydrodynamic instabilities which were identified in a previous article. The present article concerns the effect of anisotropic characteristics of the hydrodynamic dispersion on the growth of instabilities in the aquifer. Three different mechanisms may lead to instability of the flow field: (a) buoyancy forces may induce convection currents; if the difference between the convection velocity of salt, due to the hydraulic gradients, and that of heat is negligible, then this mechanism is generally most effective in planes parallel to the hydraulic velocity of the fluid (velocity due to the hydraulic gradient); (b) the difference between heat and salt effective diffusivities may lead to overstability; this mechanism is most effective in planes perpendicular to the hydraulic velocity; (c) the difference between the convection velocity of salt and that of heat may induce oscillations which are most effective in planes parallel to the hydraulic velocity. The growth of instabilities in an aquifer of unlimited length is different from their growth in an aquifer of limited length. In the latter thermohaline convection develops in planes perpendicular to the hydraulic velocity, whereas in the former it develops in planes forming an angle θ with the hydraulic gradient. The development of convection cells in the flow field is identified by numerical experiments. These experiments identify the convection cell length and the angle formed between the thermohaline convection plane and the hydraulic gradient.     

Assessment of Density‐Induced Tracer Movement in Groundwater Velocity Measurements with Point Velocity Probes (PVPs)

The concept of equivalent freshwater head was adapted to predict the conditions under which density‐driven flow would adversely impact measured groundwater velocities using point velocity probes (PVPs). Theoretically, vertical flow will result from any density contrast between the PVP tracer and the groundwater. However, laboratory testing of tracers with salinities ranging from 0 to 2000 mg NaCl/L showed that horizontal velocities could be determined with good accuracy with up to 60% of the total flow being vertical due to density effects in a gravel medium. The available data suggest that density effects are less likely to be pronounced in sandy sediments. The relative amount of vertical flow due to tracer density can be estimated from vertical and horizontal velocities measured with PVPs, or from the ratio of vertical to horizontal hydraulic gradients. The equivalent freshwater gradient produced from a given tracer salinity at 10 °C (a typical groundwater temperature at moderate latitudes) can be estimated from 7.80 × 10^?7 × (M_NaCl), where M_NaCl is the mass of NaCl added, in mg, to 1 L of site groundwater in the mixing of the tracer. Equations for other temperatures were also determined.     

Coupled Process Models of Fluid Flow and Heat Transfer in Hydrothermal Systems in Three Dimensions

Geothermal fields and hydrothermal mineral deposits are manifestations of the interaction between heat transfer and fluid flow in the Earth’s crust. Understanding the factors that drive fluid flow is essential for managing geothermal energy production and for understanding the genesis of hydrothermal mineral systems. We provide an overview of fluid flow drivers with a focus on flow driven by heat and hydraulic head. We show how numerical simulations can be used to compare the effect of different flow drivers on hydrothermal mineralisation. We explore the concepts of laminar flow in porous media (Darcy’s law) and the non-dimensional Rayleigh number (Ra) for free thermal convection in the context of fluid flow in hydrothermal systems in three dimensions. We compare models of free thermal convection to hydraulic head driven flow in relation to hydrothermal copper mineralisation at Mount Isa, Australia. Free thermal convection occurs if the permeability of the fault system results in Ra above the critical threshold, whereas a vertical head gradient results in an upward flow field.     

Electrostatic radiation of plasma in the upper ionosphere in the inhomogeneous geomagnetic field

It has been indicated that the spectrum of electrostatic waves in the ionospheric plasma depends on the geophysical conditions and solar wind parameters. The wave field measurements in the frequency band 0.1–10 MHz in the topside ionosphere were used to analyze the electrostatic instabilities of the plasma electron content (the APEX satellite experiment). A change of the sign of one magnetic field component at the geomagnetic equator can result in the formation of the large-scale irregular plasma structure with a decay of the natural electrostatic oscillations and vortices in unstable plasma. The plasma particle polarization drift from the region of decay of electrostatic oscillations and vortices can cause large plasma density and temperature gradients across the geomagnetic field. New vortices can originate at these gradients. This mechanism of plasma vortex formation and decay can be important for mass and energy convection in the topside ionosphere.     

Flow of a double-diffusive stratified fluid in a differentially-rotating cylinder

Abstract

A study has been made of a basic state of axisymmetric flow, at large rotational Reynolds numbers, in a double-diffusive stratified fluid contained in a vertically-mounted, differentially-rotating cylindrical cavity. The aim is to describe the qualitative characteristics of the flow of a fluid, the density of which is stratified by two diffusive effects, i.e., temperature and salinity gradients. Attention is confined to situations in which the temperature and salinity gradients make opposing contributions to the overall density profile, the undisturbed stratification being gravitationally stable. Finite difference numerical solutions of the governing Navier-Stokes equations have been obtained using the Boussinesq approximation. The results are presented in a way that illustrates the explicit effects of double-diffusivity when the cavity aspect ratio, height/radius, is O(1). The principal non-dimensional parameters characterizing the flow field are identified. In the interior core, the primary dynamic balance is between the horizontal density gradient and the vertical shear of the prevailing azimuthal velocity. The effective stratification is seen to decrease as the double-diffusivity increases, even if the overall stratification parameter, St, is held constant. The solute field contains a very thin boundary layer structure at large Lewis numbers. The effective stratification increases with the Prandtl number. Results have been derived for extreme values of the cavity aspect ratio. For small cavity aspect ratios, the dominant dynamic ingredients are viscous diffusion and rotation. For large aspect ratios, the bulk of the flow field is determined by the rotating sidewall. In this case, the direct influence of the double-diffusivity is minor.     

Three-dimensional numerical simulations of crustal deformation and subcontinental mantle convection

3-D simulations of mantle convection allowing for continental crust are explored to study the effects of crustal thickening on lithosphere stability and of continents on large-scale mantle flow. Simulations begin with a crustal layer within the upper thermal boundary layer of a mantle convection roll in a 1 × 1 × 1 Cartesian domain. Convective stresses cause crust to thicken above a sheet-like mantle downwelling. For mild convective vigor an initial crustal thickness variation is required to induce 3-D lithospheric instability below the zone of crustal convergence. The amplitude of the required variation decreases with increasing convective vigor. Morphologically, instability is manifest in formation of drip-like thermals that exist within the large-scale roll associated with initial crustal thickening. A strong surface signature of the drips is their ability to cause deviations from local Airy compensation of topography. After the initial thickening phase, the crustal accumulation that forms serves as a model analog to a continent. Its presence leads to mantle flow patterns distinctly different from the steady-state roll that results in its absence. Large lateral thermal gradients are generated at its edge allowing this region to be the initiation site for continued small-scale thermal instabilities. Eventually these instabilities induce a restructuring of large-scale mantle flow, with the roll pattern being replaced by a square cell. Although preliminary and idealized, the simulations do show the fluid dynamical plausibility behind the idea that significant mantle variations can be generated along the strike of a largely 2-D mountain chain by the formation of the chain itself. The ability of a model continent to cause a change in fundamental convective planform also suggests that the effects of continental crust on mantle convection may be low-order despite the seemingly trivial volume of crust relative to mantle.     

Testing Hypotheses for the Mount Isa Copper Mineralisation with Numerical Simulations

We have applied reactive transport simulations to evaluate conceptual models of hydrothermal fluid flow related to the Mesoproterozoic Mount Isa copper mineralisation. Numerical experiments have been performed specifically to investigate whether fluid flow was driven by mechanical deformation, higher than hydrostatic fluid pressure gradients, or thermal buoyancy, and what the mechanism of ore deposition was. One distinct feature of the Mount Isa mineralising system is a region of massive silica-rich alteration that surrounds the copper ore bodies within the Urquhart shale, indicating upward flow of a cooling fluid. Hydromechanical modelling revealed that contraction and horizontal shear can produce a dilation pattern that favours upward fluid flow, whereas strike slip movement causes dilation of pre-existing vertical structures. Reactive transport models show that hydraulic head driven flow is more likely to produce a more realistic silica alteration pattern than free thermal convection, but neither process generates a flow pattern capable of precipitating copper at the appropriate location. Instead we propose that gravity driven flow of a dense oxidised basin brine led to chalcopyrite mineralisation by fluid-rock reaction.     

Groundwater–surface water interactions in a river estuary and the importance of geomorphology: Insights from hydraulic,thermal and geophysical observations

This study presents the groundwater flow and salinity dynamics along a river estuary, the Werribee River in Victoria, Australia, at local and regional scales. Along a single reach, salinity across a transverse section of the channel (~80 m long) with a point bar was monitored using time-lapse electrical resistivity (ER) through a tidal cycle. Groundwater fluxes were concurrently estimated by monitoring groundwater levels and temperature profiles. Regional porewater salinity distribution was mapped using 6-km long longitudinal ER surveys during summer and winter. The time-lapse ER across the channel revealed a static electrically resistive zone on the side of the channel with a pronounced cut bank. Upward groundwater flux and steep vertical temperature gradients with colder temperatures deeper within the sediment suggested a stable zone of fresh groundwater discharge along this cut bank area. Generally, less resistive zones were observed at the shallow portion of the inner meander bank and at the channel center. Subsurface temperatures close to surface water values, vertical head gradients indicating both upward and downward groundwater flux, and higher porewater salinity closer to that of estuary water suggest strong hyporheic circulation in these zones. The longitudinal surveys revealed higher ER values along deep and sinuous segments and low ER values in shallow and straighter reaches in both summer and winter; these patterns are consistent with the local channel-scale observations. This study highlights the interacting effects of channel morphology, broad groundwater–surface water interaction, and hyporheic exchange on porewater salinity dynamics underneath and adjacent to a river estuary.     

三维球壳地幔流动形式探讨

建立三维球壳模型，将地壳和地幔作为不可压缩牛顿粘性流体处理，考虑了大陆、大洋、洋脊及俯冲带的不同物性参数，利用有限元法求解流体运动的基本方程组，研究地幔流动、板块运动与大地构造之间的关系。计算表明，速度边界条件、粘度分布等是影响地幔流动形式的主要因素。若地表径向取应力边界条件，且径向速度自由，则没有封闭的对流环存在。计算得到的地表径向速度分布与全球的构造分区具有显著的对应关系：在岛弧区、造山区速度向上运动，在大洋区、边缘海及裂谷区速度向下运动。     

Numerically quantifying the relative importance of topography and buoyancy in driving groundwater flow

Both topography and buoyancy can drive groundwater flow;however,the interactions between them are still poorly understood.In this paper,the authors conduct numerical simulations of variable-density fluid flow and heat transport to quantify their relative importance.The finite element modeling experiments on a 2-D conceptual model reveal that the pattern of groundwater flow depends largely upon the relative magnitude of the flow rate due to topography alone and the flow rate due to buoyancy alone.When fluid ...     

Mantle flow with existence of plates and generation of the toroidal field

The observed plate velocities contain two types of motions. The poloidal component is related to the formation of ridges and subduction zones and the toroidal field expresses the shearing of surface plates. One very important consideration in modeling flow in the earth's mantle is the existence and motion of the lithospheric plates. The motion of plates represents a large-scale circulation with strong viscous coupling to the mantle underneath. The mantle flow probably is neither a purely free convection driven by buoyancy forces due to nonadiabatic temperature gradients in the mantle nor a forced convection generated by boundary forces, but a mixed convection that combines the effects of boundary and buoyancy forces. We present, in this paper, the mixed convection model resulting in a surface velocity field that contains both the observed poloidal and toroidal components.     

Hydrodynamic modelling for groundwater flow through permeable reactive barriers

Hydrodynamic modelling for analysis of groundwater flow through permeable reactive barriers (PRBs) is addressed in this paper. Permeable reactive barriers constitute an emerging technology for in situ remediation of groundwater contamination and have many advantages over the traditional ex situ treatment methods. The transport domains during groundwater flow through PRBs often may involve free‐flow or non‐porous sections. To model the fluid mobility efficiently in such situations, the free and porous flow zones (PRBs) must be studied in conjunction with each other. The present paper is devoted to the analysis of groundwater flow through combined free flow domains and PRBs. The free‐flow regime is modelled using the Navier–Stokes equations whereas the permeable barriers are simulated by either the Darcy or the Brinkman equation. In order to couple the governing equations of motions, well‐posed mathematical formulations of matching boundary conditions are prescribed at the interface between the free‐groundwater‐flow zones and the permeable barriers. Combination of the Navier–Stokes equations with the Brinkman equation is more straightforward owing to their analogous forms. However, the Navier–Stokes and Darcy equations are incompatible mathematically and cannot be linked directly. The problem is resolved in this paper by invoking validated hydrodynamical expressions for describing the flow behaviour at the interfaces between free‐flow and porous zones. Three schemes for the analyses of fluid flow in combined domains are applied to the case of groundwater flow through permeable reactive barriers and different model results are compared. Copyright © 2002 John Wiley & Sons, Ltd.     

Modelling Feedback of Chemical Reactions on Flow Fields in Hydrothermal Systems

Coupled reactive transport models of hydrothermal systems provide new insights and deeper understanding of the processes occurring due to fluid flow, heat transfer, solute transport, and chemical reactions. Basic concepts of species transport (diffusion, dispersion, and advection) and chemical precipitation and dissolution reactions are discussed, and five end-member types of reactive transport environments are introduced. One of these reactive transport environments, named ‘reactions within thermal gradients’, is used to demonstrate how free thermal convection can lead to redeposition of minerals and, due to the feedback of reaction on the flow field, a change of the convection pattern. The direct consequence of changing the flow field is a significant variation of the temperature distribution within the modelled area. With the example it is shown how reactive transport simulation can be applied for the detailed study of fossil and recent hydrothermal systems.     

The Effect of Porous Medium Storage on Unstable Density‐Driven Solute Transport

Unstable density‐driven groundwater flow and solute transport (i.e., free convection) leads to spatiotemporal variations in pressure. Specific storage (S_o) indicates the capability of a confined aquifer to release or store groundwater associated with a pressure change. Although S_o is known to dampen pressure propagation, S_o has been implicitly assumed to have a negligible impact on the unstable free convective process in prior studies. This work explores the effect of S_o on both the classic onset criterion and the fingering process using numerical models. Results show that the classic onset criterion is applicable when S_o is smaller than 10^–1 m^–1. Results also demonstrate that S_o does not play a significant role in the free convective fingering process unless it is greater than 10^–3 m^–1. For most practical purposes in hydrogeology (large Rayleigh number and small S_o), the implicit assumption of small or zero S_o is appropriate.     

Bidirectional stream–groundwater flow in response to ephemeral and intermittent streamflow and groundwater seasonality

It is often assumed that the net groundwater flow direction is towards the channel in headwater streams in humid climates, with magnitudes dependent on flow state. However, studies that characterize stream–groundwater interactions in ephemeral and intermittent streams in humid landscapes remain sparse. Here, we examined seasonally driven stream–groundwater interactions in response to temporary streamflow on the basis of field observations of streamflow and groundwater on an adjacent hillslope. The direction of hydraulic head gradients between the stream and groundwater shifted seasonally. The stream gained water (head gradients were towards the stream) when storage state was high. During this period, streamflow was persistent. The stream lost water to the groundwater system (head gradients were away from the stream) when storage state was low. During this period, streamflow only occurred in response to precipitation events, and head gradients remained predominantly away from the stream during events. This suggested that mechanisms other than deep groundwater contributions produced run‐off when storage was low, such as surface and perched subsurface flowpaths above the water table. Analysis of the annual water balance for the study period showed that the residual between precipitation inputs and streamflow and evapotranspiration outputs, which were attributed to the loss of water to the deeper, regional groundwater system, was similar in magnitude to streamflow. This, coupled with results that showed bidirectionality in stream–groundwater head gradients, indicated that headwaters composed of temporary (e.g., ephemeral and intermittent) streams can be important focal areas for regional groundwater recharge, and both contribute to and receive water, solutes, and materials from the groundwater system.     

A quasi three-dimensional model for flow and transport in unsaturated and saturated zones: 1. Implementation of the quasi two-dimensional case

A quasi three-dimensional (QUASI 3-D) model is presented for simulating the subsurface water flow and solute transport in the unsaturated and in the saturated zones of soil. The model is based on the assumptions of vertical flow in the unsaturated zone and essentially horizontal groundwater flow. The 1-D Richards equation for the unsaturated zone is coupled at the phreatic surface with the 2-D flow equation for the saturated zone. The latter was obtained by averaging 3-D flow equation in the saturated zone over the aquifer thickness. Unlike the Boussinesq equation for a leaky-phreatic aquifer, the developed model does not contain a storage term with specific yield and a source term for natural replenishment. Instead it includes a water flux term at the phreatic surface through which the Richards equation is linked with the groundwater flow equation. The vertical water flux in the saturated zone is evaluated on the basis of the fluid mass balance equation while the horizontal fluxes, in that equation, are prescribed by Darcy law. A 3-D transport equation is used to simulate the solute migration. A numerical algorithm to solve the problem for the general quasi 3-D case was developed. The developed methodology was exemplified for the quasi 2-D cross-sectional case (QUASI2D). Simulations for three synthetic problems demonstrate good agreement between the results obtained by QUASI2D and two fully 2-D flow and transport codes (SUTRA and 2DSOIL). Yet, simulations with the QUASI2D code were several times faster than those by the SUTRA and the 2DSOIL codes.     

Coupled flow and salinity transport modelling in semi-arid environments: The Shashe River Valley,Botswana

Numerical groundwater modelling is used as the base for sound aquifer system analysis and water resources assessment. In many cases, particularly in semi-arid and arid regions, groundwater flow is intricately linked to salinity transport. A case in point is the Shashe River Valley in Botswana. A freshwater aquifer located around an ephemeral stream is depleted by the combined effect of transpiration and pumping. Quantitative system analysis reveals that the amount of water taken by transpiration is far more than the quantities pumped for water supply. Furthermore, the salinity distribution in and around Shashe River Valley as well as its temporal dynamics can be satisfactorily reproduced if the transpiration is modelled as a function of groundwater salinity. The location and dynamics of the saltwater–freshwater interface are highly sensitive to the parameterization of evaporative and transpirative salt enrichment. An existing numerical code for coupled flow/transport simulations (SEAWAT) was adapted to this situation. Model results were checked against a large set of field data including water levels, water chemistry, isotope data and ground and airborne geophysical data. The resulting groundwater model was able to reproduce the long-term development of the freshwater lens located in Shashe River Valley as well as the decline in piezometric heads observed over the last decade. Furthermore, the old age of the saline water surrounding the central freshwater lens could be explained.     

Impact of fracture network geometry on free convective flow patterns

The effect of fracture network geometry on free convection in fractured rock is studied using numerical simulations. We examine the structural properties of fracture networks that control the onset and strength of free convection and the patterns of density-dependent flow. Applicability of the equivalent porous medium approach (EPM) is also tested, and recommendations are given, for which situations the EPM approach is valid. To date, the structural properties of fracture networks that determine free convective flow are examined only in few, predominantly simplified regular fracture networks. We consider fracture networks containing continuous, discontinuous, orthogonal and/or inclined discrete fractures embedded in a low-permeability rock matrix. The results indicate that bulk permeability is not adequate to infer the occurrence and magnitude of free convection in fractured rock. Fracture networks can inhibit or promote convection depending on the fracture network geometry. Continuous fracture circuits are the crucial geometrical feature of fracture networks, because large continuous fracture circuits with a large vertical extent promote convection. The likelihood of continuous fracture circuits and thus of free convection increases with increasing fracture density and fracture length, but individual fracture locations may result in great deviances in strength of convection between statistically equivalent fracture networks such that prediction remains subject to large uncertainty.     

Understanding and modelling spatial drain–aquifer interactions in a low‐lying coastal aquifer—the Thurne catchment,Norfolk, UK

Seawater intrusion into fresh groundwater formations generally results inadvertently from human activities, such as over‐abstraction from coastal aquifers. This article describes the data analysis to quantify drain–aquifer interactions in a low‐lying pump‐drained coastal aquifer, which is subject to saline intrusion due to widespread land drainage, and the resulting development and application of a numerical groundwater model to understand the spatial groundwater system behaviour (including groundwater salinity fluxes). Without measured flow data in this pump‐drained catchment, a novel groundwater head‐dependent approach to hydrograph separation is described. Time‐variant and time‐invariant MODFLOW analyses are utilised to examine the flow processes. A new approach to calculate drain coefficients, which represent the extensive network of drainage ditches in the regional model, using field information, is described; the sum of the drainage coefficients are close to the values independently estimated from the head‐dependent hydrograph separation. Results show that (1) the groundwater flows into the drainage systems are well reproduced using the new drain coefficients, (2) particle tracking of fresh and saline water can explain observed spatial salinity distribution within drainage networks and (3) the modelled flow of seawater across the coast is approximately 25% greater than that discharged by the pumps, demonstrating the need for drainage management to be aware of the slow response of groundwater systems to past drainage system changes. The article demonstrates that numerical groundwater modelling can produce the improved understanding needed to inform management decisions in such complex environments. Copyright © 2010 John Wiley & Sons, Ltd.