The Vadose zone provides essential resources of water and nutrients for most surface vegetation. Besides receiving inputs of water and solutes, soil minerals undergo weathering, the products of which are contained in this biogeochemically active zone. The unsaturated zone is recognised as an important ecosystem on its own with complex nutrient cycling. Bounded by the atmosphere and the phreatic groundwater surface it has a decisive influence, for example, on the fate of pollutants from localized spills and diffuse contamination.
Numerous environmental problems require detailed process knowledge of the interaction between the water fluxes and the diverse biogeochemical reactions in the unsaturated zone in order to recognise potential threats and design sustainable management strategies. Integrated water flow and reactive transport modelling, in combination with field and laboratory experiments, provides the means for investigating important issues concerning the unsaturated zone in environmental engineering and science.
Publications
Implementation of variably saturated flow into PHREEQC for the simulation of biogeochemical reactions in the vadose zone (accepted for publication in Environmental Modelling & Software; 02.10.2009)
Abstract
A software tool for the simulation of one‐dimensional unsaturated flow and solute transport together with biogeochemical reactions in the vadose zone was developed by integrating a numerical solution of Richards’ equation into the geochemical modelling framework PHREEQC. Unlike existing software, the new simulation tool is fully based on existing capa‐bilities of PHREEQC without source code modifications or coupling to external software packages. Because of the direct integration, the full set of PHREEQC’s geochemical reactions with connected databases for reaction constants is immediately accessible. Thus, unsatu‐rated flow and solute transport can be simulated together with complex solution speciation, equilibrium and kinetic mineral reactions, redox reactions, ion exchange reactions and sur‐face adsorption including diffuse double layer calculations. Liquid phase flow is solved as a result of element advection, where the Darcy flux is calculated according to Richards’ equa‐tion. For accurate representation of advection‐dominated transport, a total‐variation‐diminishing scheme was implemented. Dispersion was simulated with a standard approach; however, PHREEQC’s multicomponent transport capabilities can be used to simulate spe‐cies‐dependent diffusion. Since liquid phase saturation is recalculated after every reaction step, biogeochemical processes that modify liquid phase saturation, such as dissolution or precipitation of hydrated minerals are considered a priori. Implementation of the numerical schemes for flow and solute transport have been described, along with examples of the ex‐tensive code verification, before illustrating more advanced application examples. These in‐clude (i) the simulation of infiltration with saturation‐dependent cation exchange capacity, (ii) changes in hydraulic properties due to mineral reactions and (iii) transport of mobile organic substances with complexation of heavy metals in varying geochemical conditions.
(C) 2009 Elsevier Ltd. All rights reserved.
Effect of mineral reactions on the hydraulic properties of unsaturated soils: Model development and application
Abstract
The selective radius shift model was used to relate changes in mineral volume due to precipitation/dissolution reactions to changes in hydraulic properties affecting flow in porous media. The model accounts for (i) precipitation/dissolution taking place only in the water-filled part of the pore space and further that (ii) the amount of mineral precipitation/dissolution within a pore depends on the local pore volume. The pore bundle concept was used to connect pore-scale changes to macroscopic soil hydraulic properties. Precipitation/dissolution induces changes in the pore radii of water-filled pores and, consequently, in the effective porosity. In a time step of the numerical model, mineral reactions lead to a discontinuous pore-size distribution because only the water-filled pores are affected. The pore-size distribution is converted back to a soil moisture characteristic function to which a new water retention curve is fitted under physically plausible constraints. The model equations were derived for the commonly used van Genuchten/Mualem hydraulic properties. Together with a mixed-form solution of Richards’ equation for aqueous phase flow, the model was implemented into the geochemical modelling framework PHREEQC, thereby making available PHREEQC’s comprehensive geochemical reactions. Example applications include kinetic halite dissolution and calcite precipitation as a consequence of cation exchange. These applications showed marked changes in the soil’s hydraulic properties due to mineral precipitation/dissolution and the dependency of these changes on water contents. The simulations also revealed the strong influence of the degree of saturation on the development of the saturated hydraulic conductivity through its quadratic dependency on the van Genuchten parameter a. Furthermore, it was shown that the unsaturated hydraulic conductivity at fixed reduced water content can even increase during precipitation due to changes in the pore-size distribution.
(C) 2009 Elsevier Ltd. All rights reserved.
Reactive transport in unsaturated soil:
Comprehensive modelling of the dynamic spatial and temporal mass balance of water and chemical components
Program capabilities:
• Transient unsaturated flow and solute transport
– Constant head and constant flux boundary conditions
– Brooks-Corey and van Genuchten hydraulic models
• Solution speciation
• Kinetic and equilibrium mineral reactions
• Redox reactions
• Solid solutions
• Surface adsorption including diffuse double layer calculations
• Influence of mineral reactions on water contents
PHREEQC software:
⇒ USGS PHREEQC project page with release notes, mail archive and links to USGS interface and Batch version.
⇒ PHREEQC interface by Vincent Post.
⇒ C.A.I. Apello’s homepage with numerous PHREEQC examples.
Code verification:
The following graphs illustrate results from extensive code verification. The simulation of water fluxes and solute transport was compared to analytical solutions [1,2] HYDRUS-1D [3]. For the verification including a complex system of geochemical reactions the HP1 software [4] was used.
Comparison of PHREEQC scheme and exact solution from Sander et al. (1988).
Comparison of PHREEQC scheme and HYDRUS-1D for the Brooks and Corey model.
Comparison of PHREEQC scheme and HYDRUS-1D for the van Genuchten model.
The use of a Kirchhoff transformation on the diffusivity term in Richards’ equation results in a higher accuracy compared to standard head-based finite difference methods and therefore allows the use of a coarser grid size than e.g. with HYDRUS-1D.
Comparison of the PHREEQC scheme (solid lines) and HP1 (dashed lines) for steady state infiltration of hyperalkaline solution in Opalinus Clay after 40 (circles) and 300 (triangles) d.
The HP1 code [4,5], which is the coupling of the variable unsaturated-saturated water flow program HYDRUS-1D [3] with the geochemical software package PHREEQC [6] was used for the validation of steady state flow and dispersive solute transport undergoing a complex set of reactions. Our scheme was compared against results from the HP1 validation case 9 (ALKALINE) [4, pp. 50-55] that simulates infiltration of a hyperalkaline solution into a 5 cm-long clay core for a period of 400 d. The core material is Opalinus Clay, which is investigated in Switzerland as a potential host formation for long term disposal of high level radioactive waste [7].
The geochemical reactions involve complex aqueous speciation, equilibrium dissolution/precipitation of calcite, sepiolite and hydrotalcite, kinetically controlled reaction of kaolinite, illite, quartz, dolomite and gypsum as well as cation exchange reactions of sodium, potassium, calcium and magnesium with a monovalent exchange species.
With respect to code verification, the results show the very good agreement of HP1 and the PHREEQC model also if complex geochemical reactions are involved. This proves the correct implementation of unsaturated reactive flow and transport. The small differences between the two models most likely result from differences in the applied solute transport schemes.
Example problems:
Transient infiltration of a hyperalkaline solution in loamy sand with complex geochemical reactions after 20 (circles) and 200 (triangles) d.
This application example shows the influence of mineral reactions on water contents and therefore flow properties in a simulation of transient infiltration.
Transient infiltration of hyperalkaline solution with surface adsorption and equilibrium phase reactions after 40 (circles) and 200 (triangles) d.
The considered geochemical reactions in this column are solution speciation of the present elements hydrogen, oxygen, sodium, calcium and carbonate together with mineral dissolution of a calcite phase and surface adsorption onto two surface types with different binding characteristics. For the definition of surface species as well as reaction equations for solution speciation, phase dissolution and surface adsorption and their parameterization with equilibrium constants the MINTEQ database (version 2; distributed with PHREEQC) was used. Surface species and their reaction constants are those of hydro-ferric oxides as derived from a compilation of experimental studies by Dzombak and Morel [8].
References:
[1] Ogata A and Banks RB. A solution of the differential equation of longitudinal dispersion in porous media. U.S. Geological Survey Professional Paper 411-A, 1961.
[2] Sander GC, Parlange JY, Kuhnel V, Hogarth WL, Lockington D, and Okane JPJ. Exact nonlinear solution for constant flux infiltration. J Hydro 1988;97(3-4):341-46.
[3] Šimunek J, van Genuchten MT, and Sejna M. The HYDRUS-1D software package for simulating the movement of water, heat, and multiple solutes in variably saturated media, version 3.0. Department of Environmental Sciences, University of California Riverside 2005.
[4] Jacques D and Šimunek J. User manual of the multicomponent variably- saturated flow and transport model HP1. SCK•CEN Waste & Disposal Department 2005.
[5] Šimunek J, Jacques D, van Genuchten MT, and Mallants D. Multicomponent geochemical transport modeling using HYDRUS-1D and HP1. J Am Water Resour As 2006;42(6):1537-47.
[6] Parkhurst DL and Appelo CAJ. User’s guide to PHREEQC (version 2): A computer program for speciation, batch-reaction, one-dimensional transport, and inverse geochemical calculations. USGS 1999.
[7] Thury M and Bossart P. The Mont Terri rock laboratory, a new international research project in a mesozoic shale formation, in Switzerland. Eng Geol 1999;52(3-4):347-59.
[8] Dzombak DA and Morel FMM. Surface complexation modeling : Hydrous ferric oxide. New York: Wiley, 1990.