Modeling of ARD
Contents
 1 Appendix 5A
 1.1 Introduction
 1.2 Approaches to Geochemical Modeling
 1.3 Role of Thermodynamic and Kinetic Data
 1.4 Scaleup Considerations
 1.5 Examples of Major Codes
 1.6 Hydrological Modeling
 1.7 Hydrogeological Modeling
 1.8 Gas Transport Modeling
 1.9 Considerations Regarding Predictive Modeling of Effluent Quality
Appendix 5A
Introduction
Geochemical Modeling
This section describes the conceptual, thermodynamic, and kinetic fundamentals of geochemical modeling and its application to prediction of mine water quality in support of mine site characterization and remediation. The emphasis in this section is on the basic processes that models attempt to represent with discussions of the usefulness and the limitations of modeling.
In principle, geochemical modeling can be applied to all mine and process facilities, including mine portal effluent, subsurface waters (wells or underground workings), waste dumps, process tailings piles, surface waters, pit lakes, and open pits. The type of modeling used depends on both the objectives and the type of source or pathway. A wide variety of codes are available for these various environments but the critical factors are the quality of their databases, the inherent assumptions, and, most importantly, the knowledge and experience of the modeler.
Three basic approaches have been used with geochemical data: forward geochemical modeling, inverse geochemical modeling, and geostatistical analyses.
Forward modeling is also known as simulating (i.e. potential reactions between rock and water are simulated from initial conditions of a known rock type and composition). Reactions are allowed to proceed in equilibrium or kinetic or combined modes. Changes in temperature and pressure can be invoked, changes in water flow rate can be invoked, and minerals can be allowed to precipitate as they reach equilibrium solubility or dissolve as they become undersaturated. Potential reactions can be simulated to see what the consequences are. This type of modeling is the least constrained. A great many assumptions are either invoked as input data or invoked as dictated by the program that may not apply to the specific system being simulated. This approach assumes the modeler has a significant amount of information on the ability of minerals to maintain equilibrium solubility or their rates of reaction.
Inverse modeling assumes a water flow path is known and that water samples have been analyzed along that flow path. Such data can then be converted into amounts of minerals dissolved or precipitated along that flow path. Several assumptions are still made regarding the choice of minerals and their relative proportions contributing to the water chemistry, but the calculations are constrained with actual data. Inverse modeling can also be done without any recourse to kinetic or thermodynamic data, in which case it represents a relatively simple mass balance calculation. When speciation and thermodynamic and kinetic properties are included for additional constraints, the possible reactions become quite limited and the modeling is much more meaningful.
Modeling of any type does not lead to a unique solution but the possibilities are more limited with greater amounts of carefully collected field data. Martin et al. (2005) summarized the benefits and limitations of geochemical modeling as follows:
Benefits
 Provide insight into potential future conditions
 Determine which variables are most important in determining future conditions
 Assess the effects of alternative approaches to ARD management
 Assess potential effects of uncertain parameters
 Establish objectives and test conditions for field and laboratory studies
 Integrate available information
Limitations
 Insufficient input data
 Modeling can be challenging and results misinterpreted
 Uncertain and variability of the results
 Difference between modeled and actual field conditions
Approaches to Geochemical Modeling
Speciation
One of the most fundamental types of geochemical modeling is speciation modeling. “Speciation” refers to the distribution of chemical components or elements among the different possible forms or species. Aqueous speciation is the distribution of chemical components among dissolved free ions, ion pairs and triplets, and other complexes. This concept is important because research has shown that a number of processes, including mineral precipitation and dissolution, biological uptake and toxicology, and sorption are all affected by speciation.
Some species, such as redox species, have to be determined analytically. This is because most geochemical modeling codes erroneously assume that redox equilibrium is maintained while, in reality, disequilibrium among redox species is the rule, not the exception. In particular, dissolved iron is usually present in high concentration in ARD and can exist as the reduced ferrous iron or as the oxidized ferric iron. For an accurate evaluation of iron speciation, chemical analysis rather than speciation modeling is required. In NMD and SD, dissolved iron is largely absent due to formation of sparinglysoluble ferrihydrite or similar iron oxyhydroxide minerals. Solid speciation is the distribution of chemical components among various solid phases. For example, iron can precipitate from ARD as goethite, jarosite, schwertmannite, or ferrihydrite. Identifying these phases would constitute solid speciation.
Aqueous speciation results are used in a variety of modeling objectives that include modeling of saturationindex calculations for masstransfer, modeling of mineral precipitation and dissolution, modeling of adsorption and desorption, and reactivetransport modeling.
Mass Transfer (precipitation, dissolution, gas transfer)
Modeling of mineral precipitation and dissolution and gastransfer reactions can take place conceptually in one of three possible systems: equilibrium state, steadystate, or transient state. The equilibrium state assumes the system under investigation is isolated from any external exchanges of energy or mass. Although an unrealistic concept, equilibrium state is actually quite practical because many reactions approximate equilibrium even though there are gradients in water pressure or temperature. For example, in many groundwaters, calcite and gypsum quickly reach their equilibrium solubility. Even with gradients in CO2 pressure or mixing with other sources of sulphate, these minerals adjust to maintain saturation and the assumption of equilibrium may be valid. In addition, even when geochemical reactions of interest do not reach equilibrium rapidly, such reactions may achieve equilibrium over the time scale of the modeling simulation (i.e. the life of a mine and beyond). Therefore, the majority of geochemical modeling can be conducted under the assumption of equilibrium conditions.
Reactive Transport (Coupled Models)
Reactivetransport models that can be applied to simulation of ARD, NMD, and SD are generally the subject of active research, although several have been applied with considerable success. The idea is to couple flow models with chemical reaction models to determine the effects of flow on reactions and vice versa, including the effects of dispersion. Such modeling is relatively straightforward for streams and rivers because the flow path is not only visible but measurable. Considerable effort has been made to develop quantitative reactiontransport models for streams affected by acid mine drainage (Kimball et al., 1994; Runkel et al., 1996). Progress in surfacewater reactivetransport modeling has now advanced to the point where it can guide remediation decisions for complex mine sites (Runkel and Kimball, 2002; Kimball et al., 2003).
Reactivetransport modeling for groundwater has also progressed substantially over the last two decades and many of the recent codes have been applied to mine sites. Three general types of coupled models can be distinguished: those that model the groundwater only, those that model the unsaturated zone only, and those that model both. The most recent overview by Mayer et al. (2003) provides the theoretical foundations for groundwater reactivetransport modeling, methods of coupling flow with reaction, the various codes that have been used in mined environments, and case studies. An excellent example of combining laboratory testing of waste rock material with field measurements and modeling of small to mediumscale test plots of actual mine wastes to predict the consequent water quality over the short term and the long term in a very sensitive environment is in progress at the Diavik mine site near Yellowknife, Northwest Territories, Canada (Blowes et al., 2007). This investigation may be one of the first to combine labscale tests, field tests, and modeling, supported by the detailed characterization of the rock and mineral composition and their weatherability.
Role of Thermodynamic and Kinetic Data
Thermodynamic and, for some models, kinetic data are part of the basic input to codes that compute reactions and simulations for waterrock interactions. For some reactions, these data are known accurately and precisely; for others they are nonexistent or poorly known. Thermodynamic measurements and evaluations are part of ongoing research. Sometimes the conclusions of a modeling study can be greatly affected by these databases and their uncertainties and sometimes not. Rarely are modeling results evaluated from the point of view of the basic data, which reflects a general lack of QA/QC common to many modeling efforts.
Scaleup Considerations
Drainage quality prediction is made challenging by a number of factors that range in scale from small to large. Smallscale factors that influence drainage quality are related to reactions at the waterrock interface in the aqueous, gas and solid phase. Information on reactive surface area and reaction rates generally is limited. On a large scale, geology, climate, mining method, mineral processing method, and waste management practices vary within and amongst operations. Variability of these largescale factors implies that it may not always be feasible to apply information from one site to another. However, advances are being made in this respect, for instance, through the use of geoenvironmental models that present unifying principles which link mine water quality to the nature of the ore deposit, climate, and type of mine waste.
Water quality prediction typically necessitates the extrapolation of laboratoryscale results to operational scale. This extrapolation must address factors such as differences in particle size, climate conditions, water and gas transport, and duration (i.e., how these variables affect drainage composition over decades, centuries or longer). Although the construction of instrumented, largescale mine waste test cells has increased significantly in recent years and is expected to yield valuable data, little information is currently available describing the effects of these variables on wellcharacterized mine wastes over extended periods of time. Use of models therefore is required to bridge the gap between laboratory results and operational conditions (USEPA, 2003).
Examples of Major Codes
Some of the more popular codes used primarily for groundwater geochemistry but also for miningaffected sites are shown in Table A51 below. More detail on geochemical modeling, modeling codes and associated uses and limitations is presented in Alpers and Nordstrom (1999), Mayer et al. (2003), and Maest and Kuipers (2005). Section XXX on hydrogeological models in this Appendix also provides additional information.
A51: Summary of Geochemical Modeling Codes
Codes 
Type 
Reference 
PHREEQC, PHAST 
USGS codes: mass transfer and reactivetransport 
Parkhurst and Appelo (1999), Parkhurst et al. (2004) 
SOLMINEQ.GW 
USGS code: mass transfer and high temperature 
Perkins et al. (1990) 
WATEQ4F 
USGS code: speciation and lowtemperature only 
Ball and Nordstrom (1991) 
MINTEQA2 
EPA supported code: speciation and mass transfer 
USEPA (1999) 
MIN3P 
Waterloo code: saturated and unsaturated flow 
Mayer et al. (2002) 
TOUGHAMD 
Quebec code: gas and energy transfer without reaction 
Lefebvre et al. (2002) 
RETRASO 
Barcelona code: unsaturated and saturated flow and reaction 
Saaltink et al. (2002) 
SULFIDOX 
ANSTO code: gas and energy transfer 
Ritchie (2003) 
CRUNCH 
Bern/LBL code: unsaturated and saturated flow and reaction 
Steefel (2000) 
Geochemist’s Workbench 
University of Illinois code: mass transfer, saturated flow 
Bethke (1994, 1996) 
EQ 3/6 
Lawrence Livermore National Laboratory code: mass transfer and reactive transport 
Wolery and Daveler (1992) 
Hydrological Modeling
Introduction
In a general sense, a hydrological model is an analog of a natural or humanmodified hydrological system. This generic definition encompasses models of surfacewater and groundwater systems. Scientists and engineers more commonly use the term hydrological model to refer to models of surfacewater systems, and consider hydrogeological models for groundwater systems as a separate subject. This section follows the latter convention, describing hydrological models in the context of surfacewater systems. Hydrogeological models and their applications are presented in Section XXX.
Hydrological models range from simple algebraic calculations to complex reactivetransport computer codes. Physical analogs, such as stream tables, can also be useful simulations of complex surfacewater systems. Hydrological models can be used to predict the fate and transport of mine drainage through a surfacewater system, providing important input to humanhealth or ecological risk assessments. Hydrological models can also be used to estimate the waterquality and waterquantity evolution of pit lakes over time. Hydrological models can be coupled with hydrogeological and geochemical models to incorporate the interaction between surface water and groundwater into the simulation and account for geochemical reactions.
Selection of an appropriate, quantitative hydrological model depends on the type of output that is required and, critically, on the conceptual model of the system being evaluated. A robust conceptual model will identify the important physical and geochemical characteristics of the fieldscale system being evaluated. Based on that identification, an appropriate hydrological model can be selected that quantitatively represents those important processes. For complex systems or to assess a range of different types of processes, multiple hydrological models can be applied to predict the fate, transport, and potential impacts of mine discharges.
Data Needs
In common with all models, the output from a hydrological model is only as reliable as the data that are used to generate the model. Typical data requirements for many hydrological models include:
 Precipitation, either local or distributed across a region
 Evaporation from surfacewater bodies such as lakes and rivers
 Potential or actual evapotranspiration from vegetated areas and bare land
 Surface slope and land cover
 Channel slope, width, depth, and roughness for calculations of stream flow or conveyance capacity
 Concentrations of chemical constituents. These may be determined from onsite monitoring programs, laboratory or fieldscale testing programs, or estimated using geochemical models
Simple quantitative models of surfacewater flow such as the United States Natural Resource Conservation Service (formerly known as the Soil Conservation Service [SCS]) curvenumber method (SCS, 1972) may only require a few of before listed data elements. More detailed models, for instance those that incorporate reactive transport (e.g., Runkel and Kimball, 2002), may require additional information regarding the kinetics of reactions considered in the simulation.
Governmental agencies in many countries collect regional precipitation and evaporation data that may be used for hydrological models. Precipitation data are commonly collected with the greatest frequency through meteorological measuring stations. Evaporation data, such as pan evaporation measurements, generally are collected with less frequency. Some mine sites also collect these types of data on a local scale that can be used to refine the regional data sets.
Care must be taken if combining different types of data from different locations. The locations should be similar in terms of latitude, elevation, overall climatic zone, and cloud cover for the combined data set to be reliable. If this is not the case, statistical methods have been developed to estimate precipitation at a special site from a known precipitation network
WaterBalance and MixingCell Models
Waterbalance models apply the principle of conservation of mass to quantitatively track inflows and outflows from the various components of a conceptual model. Mass and concentration of ARDrelated constituents can be incorporated into this approach through mixingcell models. The hydrologic elements of a conceptual model, such as surfacewater reservoirs, open pits, and groundwater basins, can be represented as a series of simulated reservoirs. The connections between the reservoirs, such as the creeks or groundwater flow paths, can be represented by quantitative estimates of capacity or flow. Concentrations of individual constituents can be tracked along with water quantity to calculate the transfer of chemical mass and mathematically mixed in the model to evaluate changes in concentration over time in the reservoirs.
Waterbalance and mixingcell models can be implemented in standard spreadsheets. More complex waterbalance or mixingcell models, incorporating additional physical or chemical processes, can be addressed by using dynamic system simulators such as GoldSim^{TM} or STELLA^{TM}.
Rainfall Runoff Models
Appendix A of USEPA (2003) describes the basic approaches to modeling runoff processes based on precipitation inputs. Runoff can be thought of as the excess precipitation after processes such as infiltration and surface abstraction are evaluated. The most commonly applied model to estimate the volume of runoff is the SCS curvenumber method (SCS, 1972). The SCS curvenumber method involves estimating the vegetation and landcover characteristics of a watershed or mine facility, looking up the resulting curve number, and then applying that number along with precipitation information to develop the runoff volume for a storm event.
The unithydrograph method of runoff determination may be more appropriate for many mine sites. The method is also described in SCS (1972). A hydrograph relating runoff to precipitation is developed for a unit precipitation volume over an area, for example 1 inch or 1 centimeter of rainfall. The unit hydrograph is then used to estimate runoff from storms of greater or lesser intensity.
Water quality in wellmixed rivers and streams can be predicted using a code such as QUAL2K developed by the USEPA (Chapra et al., 2007). QUAL2K represents a modernized version of QUAL2E (Brown and Barnwell, 1987). QUAL2K is programmed in the Visual Basic for Applications language and executed within the Microsoft Excel spreadsheet environment. The program can simulate 1dimensional flow, changes in water quality along the flow path, and chemical interactions with bed sediments.
Distributedparameter rainfallrunoff models are more appropriate for larger watersheds with heterogeneous flow characteristics. SWAT2000 (Neitsch et al., 2002) is a distributedparameter model developed by the Agricultural Research Service of the U.S. Department of Agriculture to simulate runoff and water quality in large, complex watersheds. SWAT2000 and similar models break a complex watershed into hydrologic subunits, each with a uniform set of characteristics. Flow and water quality are calculated for each subunit, then aggregated to provide predictions at a complex watershed scale.
Pit Lake Modeling
Pit lake formation and the evolution of water quality can be simulated using a waterbalance approach or with complex numerical codes. Water balance models can be used to quantify the inflows to the pit lake as the pit fills after mining and dewatering ceases. Potential inflows include direct precipitation over the surface area of the lake, runoff entering the pit lake from the surrounding watershed, and groundwater inflow through the walls and floor of the pit. Outflows may include direct evaporation from the lake surface, groundwater outflow, and potentially surfacewater discharges if a spill elevation is reached.
A chemical composition can be assigned to each inflow and outflow to extend the waterbalance model to include ARDrelated impacts. For example, wallwashing results can be used to estimate the mass input of chemical constituents from seepage or overland flow coming in contact with reactive portions of the pit wall. Geochemical speciation models can be used to predict the resulting chemical quality of water in the pit.
Rainfallrunoff models can be used to develop the surfacewater inflow portions of the water balance. Groundwater inflow can be estimated using simple analytical equations (Marinelli and Niccoli, 2000). The solution to drawdown in a largediameter pumping well presented by Papadopoulos and Cooper (1967) is often used to approximate the groundwater inflow to a mine pit, and can also be used to estimate recharge to the pit lake. Cimen (2001) and Aryafar (2007) present additional analytical solutions that can be useful in pitlake studies.
Complex numerical models can also be used to estimate the groundwater inflow to a pit lake. SEEPW (Ref) and FEFLOW (WASY, XXXX) are finiteelement, variablysaturated flow models that have been applied to this problem. MODFLOW2005 (Ref), including the LAK package, is a modular, 3dimensional, finitedifference model that can be used to simulate the groundwater components of pitlake evolution. Complex models such as these, however, require more data for parameterization and calibration than the simpler approaches. Selection of more complex simulation approaches should only be made if the conceptual model and project needs require the additional computational burden.
An alternative to geochemical models for the prediction of pitlake quality is a code such as CEQUALW2 developed by the U.S. Army Corps of Engineers Waterway Experiment Station (Cole and Buchak, 1995). CEQUALW2 is suitable for applications to rivers, lakes, reservoirs, and estuaries.
Watershed Models
Watershed models are used to simulate the hydrologic cycle, including surface water, groundwater, and the interactions between the two, at the basin or watershed scale. Watershed models can be used to predict ARD impacts on downstream users and the evolution of ARDrelated water quality through a flow system. Furman (2008) summarizes the mathematics and computational tools used to simulate coupled surface and subsurface flow processes.
Watershed models can be dataintense and numerically complex. The most widely used watershed models are:
 MIKE SHE, developed by the Danish Hydraulic Institute (DHI) in Denmark
 HECHMS, developed by the U.S. Army Corps of Engineers Hydrologic Engineering Center
 WMS (Watershed Modeling System), a graphical interface developed by Environmental Modeling Systems, Inc. for a number of modules including HECHMS, CEQUALW2 and other codes
Hydrogeological Modeling
Introduction
Hydrogeological models address water flow and contaminant transport below the land surface. As with hydrological models, approaches to hydrogeological simulations range from simple to complex. The universe of hydrogeological models includes physical and electrical analogs. With the advent of powerful personal computers and highlevel programming languages, these approaches are rarely used in current practice. Accordingly, this discussion of hydrogeological models will focus on quantitative, mathematical approaches to subsurface water flow and contaminant transport.
A large body of literature exists regarding hydrogeological modeling, as do a number of computer programs. Zheng and Bennett (2002) provide an excellent introduction to the topic of contaminanttransport modeling. Maest and Kuipers (2005) provide a review of hydrogeological models more directly focused on ARD prediction. Other references are provided in the discussion below.
Three basic types of hydrogeological models are available, in order from simple to more complex:
1. Analytical models of flow and contaminant transport
2. Analytic element models
3. Numerical models
As a general rule, hydrogeological models should be as simple as possible while still representing the physical system with an adequate degree of precision and accuracy. More complex models should only be selected when project needs dictate, simpler models are demonstrably not adequate, and suitable data are available for model parameterization and calibration.
Hydrogeological models are useful tools for predicting the potential generation and resulting impacts of ARD. Models can be used to fill data gaps, either in space or in time. They can also be used to test alternative conceptual models in an iterative process designed to understand the complex natural or humanmodified subsurface system.
Figure 519 in Chapter 5 of this GARD Guide presents a generalized approach to the development, calibration, and use of models, including hydrogeological models. The quantitative modeling process starts with a strong conceptual model, and the quantitative model can then be used to update the conceptual model as necessary. The majority of the effort for a hydrogeological model goes into the calibration phase of the process, sometimes also referred to as inverse modeling.
Data Needs for Model Parameterization
Basic Flow and Transport Models
As model complexity grows, the data requirements for model parameterization and calibration also increase. Basic data requirements for any groundwater flow and contaminanttransport model include:
 Saturated hydraulic conductivity
 Specific yield or storativity
 Effective porosity (for calculations of contaminant transport)
UnsaturatedZone Models
Simulating flow and transport in the unsaturated zone typically requires additional information regarding the flow characteristics of the unsaturated porous medium. Unsaturated hydraulic conductivity is a function of the saturated hydraulic conductivity and the degree of saturation of the porous medium. Additional data requirements for unsaturatedzone models include the parameters for the function describing the relationship between saturation, matric suction, and unsaturated hydraulic conductivity.
Sorption and Retardation Factors
Interaction between the aquifer matrix and dissolved constituents can be an important process for ARDrelated hydrogeological models. Many contaminanttransport models simulate this interaction through the use of a retardation factor.
The retardation factor is the ratio between the apparent velocity of the contaminant front and the pore velocity of moving groundwater (Fetter, 1993). In its simplest form, the retardation factor is calculated using a distribution coefficient appropriate for a linear adsorption isotherm. More complex forms of the retardation factor can be derived using different adsorption isotherms and assumptions.
Models incorporating retardation thus require additional data, including:
 Bulk density of the aquifer matrix
 Distribution coefficient or other parameters defining the adsorption isotherm
 Rate constants for nonequilibrium sorption models
ReactiveTransport Models
As discussed in Section XXXX under geochemical modeling, detailed evaluation of the evolution of ARDrelated constituent concentrations over time and space in an aquifer may require the use of a reactivetransport model. These types of models allow the simulation of reactions between the dissolved constituents and the aquifer matrix and reactions between the dissolved species themselves. Relative to the more basic hydrogeological models, additional data are necessary to apply these models, including:
 Nonequilibrium rate constants describing the reactions between dissolved constituents
 Proportionality constants or functions describing the solubility controls on individual species under consideration
Steefel et al. (2005) and Mayer et al. (2003) provide overviews of reactivetransport models and their associated data requirements.
Data Collection
The field data most commonly obtained in support of hydrogeological modeling are the saturated hydraulic conductivity and storage coefficients (specific yield or storativity). Saturated hydraulic conductivity can be measured on core samples in the laboratory, by using singlewell slug tests, or by using multiplewell, longterm pumping tests.
Slug tests and pumping tests provide better estimates of saturated hydraulic conductivity at the field scale than laboratory tests. Pumping tests conducted with one or more pumping wells in combination with at least one additional observation well can also provide data regarding the storage coefficients. Butler (1998) provides an extensive description of the design and performance of slug tests. Kruseman and de Ridder (2000) describe the design and performance of pumping tests.
The relationship of unsaturated hydraulic conductivity to moisture content can be measured in the field or laboratory, and the resulting data can be fitted to a number of equations. Stephens (1995) provides a detailed description of data collection and analysis related to unsaturatedzone hydrology.
Other Data Sources
Unsaturated hydraulic conductivity characteristic curves can be estimated by several methods. RETC (van Genuchten et al., 1991) and ROSETTA (Schaap, 2003?) are programs that can be used to estimate unsaturated flow characteristics from more commonly available data. SoilVision (SoilVision Systems, XXXX) contains a database of measured unsaturated hydraulic conductivity characteristic curves in addition to a number of algorithms to calculate unsaturated flow characteristics.
Adsorptionisotherm distribution coefficients for a number of metals are tabulated in Stenge and Peterson (1989). Values are included for three different pH ranges and a range of sorbent (organics, oxides, clays) contents.
Analytical Models
Analytical models are relatively simple methods for simulating groundwater flow and contaminant transport. These models are formulated as closedform equations that can be solved directly without the use of numerical methods. Transient or steadystate solutions for groundwater flow and contaminant transport with simple retardation factors in one, two or three dimensions are available.
Because of their simplicity, data needs are relatively minor for analytical models. Homogeneous, isotropic flow conditions are typically assumed. Analytical models can be useful for screeninglevel evaluations. They can also be used for more definitive assessments of groundwater flow and contaminant transport if the assumptions are judged to be valid or insufficient data are available to warrant a more complex approach.
One useful analytical model for the prediction of ARDrelated transport is the Ogata and Banks (1961) solution to the advectiondispersion equation. Domenico and Schwartz (1990) extended that solution to include a retardation factor based on a linear adsorption isotherm. The Domenico and Schwartz (1990) model can be implemented in a spreadsheet format and adapted to a wide variety of problems.
STANMOD (Simunek et al, 2003) is a public domain set of analytical solutions to the advectiondispersion equation in one, two or three dimensions. A variety of previously published solutions, already in the public domain, are included in STANMOD.
Analytic Element Models
Analytic element modeling takes advantage of the principle of linear superposition to solve groundwater flow and contaminant transport problems in more complex systems than can be addressed by analytical models. Haitjema (1995) provides the basic theoretical framework for the analytic element method and describes its use.
GFLOW (Haitjema Software, 2007) is a groundwater flow model that implements the analyticelement method. PhreFlow (Jankovic and Barnes, 2001) is a public domain analyticelement model of 3dimensional groundwater flow and contaminant transport. WhAEM2000 (Kreamer et al., 2007) is a public domain and open source general purpose groundwater flow modeling system, with strengths in representing regional flow systems, and ground water/surface water interactions. It was initially designed to facilitate capture zone delineation.
Numerical Models
Numerical models use iterative processes to solve the equations of groundwater flow and contaminant transport in complex domains. Flow and transport under saturated, unsaturated, or variablysaturated conditions in heterogeneous, anisotropic systems with various boundary conditions can be simulated using these methods. Numerical models can also require a substantial amount of data regarding parameters and for input to the simulation.
Two numerical solution schemes, finitedifference and finiteelement, are widely used in hydrogeological models (Wang and Anderson, 1982). Finitedifference models employ a rectangular discretization scheme to divide the model domain into individual cells, within which flow characteristics such as hydraulic conductivity are assumed to be uniform. Finiteelement models employ either a triangular or rectangular discretization scheme to divide the model domain into individual elements of uniform characteristics.
As a general rule, finitedifference models are more computationally efficient for a given problem compared to finiteelement models. Finiteelement models can be fitted more closely to irregular boundaries and can handle internal boundaries such as mine pits, underground workings, or faults with less numerical instability than finitedifference models. The choice of numerical solution scheme and computer code should be driven by the conceptual model, project requirements, and available computer resources.
The MODFLOW family of computer codes (e.g., MODFLOW2005, Ref) contains examples of finitedifference models. MODFLOW, originally released by the USGS in 1988 and upgraded periodically since then, is probably the most widely used hydrogeological model in the world.
Finiteelement models are exemplified by FEFLOW (WASY, XXXX). FEFLOW is a commercially available code that can be applied to a broad range of variablysaturated flow and transport problems. Compared to MODFLOW it has more capabilities for modeling mine water problems because the original program code was derived from a mining background.
Many other finitedifference and finiteelement models suitable for application to ARDprediction problems are available. Maest and Kuipers (2005) tabulate the capabilities for a range of models.
Commonly Available Models
Table A52 summarizes the characteristics of several numerical model computer codes that are widely used and can be applied to problems of ARDrelated contaminant transport. Some models are freely available in the public domain, while others are proprietary products distributed by commercial companies. Table A52 is organized by computer program and by graphical user interface (GUI). More information on GUI use and characteristics is presented below.
UnsaturatedZone Models
Unsaturatedzone models are often used to assist with predicting the formation and transport of ARD within and through wasterock dumps and unsaturated process tailings impoundments. Commonly used unsaturatedzone models include HELP (Schroeder et al., 1994), HYDRUS (Simunek et al., 2007), UNSATH (Fayer, 2000), and VADOSE/W (GeoSlope International, 2002). HELP and UNSATH are available in the public domain, as is the 1dimensional version of HYDRUS.
FractureFlow Models
The majority of hydrogeological models are strictly valid for simulating flow and transport through continuous porous media only. However, some ARD problems occur in subsurface systems dominated by flow and transport through discrete fractures or fracture networks. Even if flow and transport are primarily through fractures, porousmedium models may be adequate if the fracture density is great and the fracture aperture is small. Some models allow a dualporosity formulation that can represent the flow through a fracture network as well as flow through the porous media between fractures.
If the assumption of flow through continuous porous media is not valid, models that represent the physics of fracture flow should be considered. Two such models are FRACMAN (Golder, 2007) and FRACTRAN/FRAC3DVS (University of Waterloo, 2004).
DensityDependent Flow and Transport
Most contaminanttransport models are based on the assumption that concentrations are relatively dilute and the density of groundwater is not significantly different from fresh water. Groundwater that is heavily impacted by ARD, however, can have sufficiently large concentrations of metals, sulfate, and other species that the density effects are significant. If the conceptual site model indicates that density effects are important, a model capable of accounting for variability in density should be selected.
SEAWAT2000 (Guo and Langevin, 2002) is one such model, developed by the USGS to simulate 3dimensional, variabledensity groundwater flow in porous media. It was developed by combining MODFLOW and MT3DMS into a single program that solves the coupled flow and solutetransport equations.
ReactiveTransport Models
Most contaminanttransport models incorporate relatively simple reactions describing interaction between dissolved constituents and the aquifer matrix. These reactions are implemented in the form of retardation factors using one of several adsorption isotherms. Interactions between dissolved constituents are typically not considered.
Table A52: Hydrogeological Models and Graphical User Interfaces for those models
GUI 
Groundwater Vistas 
Groundwater Modeling Systems 
Visual MODFLOW 
Argus ONE 
PMWIN 
Description 

FLOW 
SEEP2D 
x 
A 2D finiteelement groundwater model designed to be used on profile models such as crosssections of earth dams or levees 

MODAEM 
x 
Analytic element model for simple flow and transport computations 

MODFLOW 88 
x 
xXx 
MODFLOW is a 3D, cellcentered, finite difference, saturated flow model developed by the USGS 

MODFLOW 96 
x 
x 
xXx 
MODFLOW is a 3D, cellcentered, finite difference, saturated flow model developed by the USGS 

MODFLOW 2000 
x 
x 
x 
x 
xXx 
MODFLOW is a 3D, cellcentered, finite difference, saturated flow model developed by the USGS 

MODFLOW 2005 
x 
x 
MODFLOW is a 3D, cellcentered, finite difference, saturated flow model developed by the USGS 

FEMWATER 
x 
3D finiteelement model used to simulate densitydriven coupled flow and contaminant transport in saturated and unsaturated zones 

SOLUTE 
ART3D 
x 
A threedimensional analytic reactive transport model 

MODPATH 
x 
x 
x 
x 
A particle tracking code used with MODFLOW assuming particles are transported by advection 

PATH3D 
x 
General particle tracking program for calculating groundwater paths and travel times in steadystate or transient, 2 or 3D flow fields 

PMPATH 
xXx 
Particle tracking 

MOC3D 
x 
xXx 
3D methodofcharacteristics groundwater flow and transport model integrated with MODFLOW96 

MT3D 
x 
x 
x 
x 
xXx 
Simulation of singlespecies transport via advection, dispersion and simple chemical reactions 

MT3DMS 
x 
x 
x 
x 
xXx 
Simulation of multispecies transport by advection, dispersion, and limited chemical reactions of dissolved constituents 

PHT3D 
x 
x 
A reactive transport model coupling MT3DMS and PHREEQC 

RT3D 
x 
x 
x 
x 
An advanced multispecies reactive transport model developed by the Battelle Pacific Northwest National Laboratory 

SEAM3D 
x 
Reactive transport model to simulate complex biodegradation problems involving multiple substrates and multiple electron acceptors 

SEAWAT 2000 
x 
x 
x 
Simulation of 3D, transient, variabledensity ground water flow 

MODFLOW 
x 
x 
Enhanced simulation capabilities and robust solution methods for handling complex saturated/unsaturated flow and transport 

MODFLOWT 
x 
Version of MODFLOW that includes modules for simulating 3D contaminant transport 

SWIFT 
x 
3D model to simulate groundwater flow, heat, brine and radionuclide transport 

UTCHEM 
x 
A multiphase flow and transport model ideally suited for pump and treat simulations 

CALIBRATION 
MODFLOW 2000 
x 
x 
Parameter inversion option built into MODFLOW 2000 

UCODE 
x 
x 
x 
Developed by the USGS, UCODE is a universal inverse modeling code to solve parameter estimation problems 

PEST 
x 
x 
x 
x 
A modelindependent, nonlinear parameter estimator to assist in data interpretation, model calibration, and predictive analysis 

Stochastic 
x 
x 
Parameter inversion using Monte Carlo or Latin Hypercube 

Modac 
x 
An inverse model that calculates a K for every cell in the model (or in selected layers) using starting heads as the calibration target 

Automated 
x 
Automated sensitivity analysis that can be used for initial calibration or to test parameter sensitivity 

OPTIMIZE 
Somos 
x 
Optimization modules to aid in optimally managing water resources 

Brute Force 
x 
Optimization based on plume containment 

MODOFC 
x 
MODOFC is designed to allow the user to create and solve optimization problems for hydraulic control in groundwater systems 

MGO 
x 
x 
Optimizes groundwater management and remedial strategy based on various physical, environmental and budgetary constraints 

GRAPHICS 
GIS 
Import, export 
Import, export 
Import, export 
Import, export 

AutoCAD 
Import, export 
Import 
Import, export 

Registered 
Import 
Import 
Import 
Import 

Surfer 
Import, export 
Import 
Import, export 

EQuIS 
Import 

Slicer 
Export 

Earth Vision 
Import, export 

EVS 
Import, export 

Tecplot 
Export 

Prop. 3D 
yes 
yes 

Local Mesh Refinement? 
yes 
yes 
no 
no 
yes 
If reactions between dissolved constituents, or precipitation and redissolution of individual constituents are important processes, reactivetransport models may be necessary to adequately represent the hydrogeological system. PHAST and PHT3D are two potential choices for this type of model.
PHAST (Parkhurst et al, 2004) simulates multicomponent, reactive solute transport in threedimensional saturated groundwater flow systems. The flow and transport calculations are based on a modified version of HST3D (Kipp, 1997) that is restricted to constant fluid density and constant temperature. The geochemical reactions are simulated with the geochemical model PHREEQC, which is embedded in PHAST.
The publiclyavailable code PHT3D (Prommer, 2002) couples MT3DMS and PHREEQC and therefore works within the MODFLOW scheme. PHT3D provides the highest level of coupling between constant density flow and fully reactivetransport codes. Because PHT3D couples MT3DMS with PHREEQC, it cannot be used simultaneously with SEAWAT2000, which also uses MT3DMS to couple MODFLOW.
Models available for fully coupled reactive flow and transport with density effects are severely limited. PHWAT (Mao et al., 2006) incorporates PHREEQC2 into SEAWAT and provides the necessary capabilities. However, the model is still in development and not available commercially. It can simulate multicomponent reactive transport with variable density groundwater flow.
Graphical User Interfaces
The raw input files for many hydrogeological models, including MODFLOW, are quite userunfriendly. Model inputs are typically via multiple text files using lineentry and array format. Large models can be quite difficult to manage. Fortunately, several GUIs have been developed that are userfriendly and simplify the process of developing, calibrating, and using hydrogeological models.
In general, GUIs provide Windowsbased interfaces that simplify pre and postprocessing for MODFLOW and other hydrogeological models. Several GUIs provide interfaces with AutoCAD, Geographic Information Systems (GIS), SURFER (Golden Software, 2002) or other graphical programs to directly input material properties and boundary conditions as well as visualize model outputs.
GUIs also provide interfaces with addon modules such as calibration and optimization routines, including UCODE and PEST. Some GUIs provide interfaces with these codes in addition to the inversemodeling routines contained within MODFLOW. Further, a suite of optimization codes can be used to evaluate a variety of hydrologic issues related to groundwater pumping, plume management, cost effectiveness, and receptor management for contaminated areas.
Local Mesh Refinement (LMR) provides the ability to create submodels within a regional model. While submodels cannot be used simultaneously with a regional model, they can be used to refine calibration or predictions within a smaller area after solving the regional model. Some GUIs provide this function while others do not
Table A53 provides a comparison of capabilities of five widely used GUIs:
 Groundwater Vistas (Environmental Simulations Inc., 2007)
 Groundwater Modeling System (GMS;Environmental Modeling Systems Inc., 2007)
 Visual MODFLOW (Schlumberger Water Services, 2007)
 Processing MODFLOW for Windows (PMWIN; Chiang, 2005)
 Argus Open Numerical Environments (ONE; Argus Holdings Ltd., 1997)
Argus ONE is an open environment for creating GUIs adapted to specific models. The USGS and others have developed interfaces within Argus ONE for a number of hydrogeological models. The other 4 GUIs are distributed as packages with their respective models included in the distribution.
Model Calibration
Calibration of a hydrogeological model is an application of inverse modeling. Model calibration is the process of selecting parameter values, inputs, and boundary conditions such that model output matches related observed data with an acceptable degree of accuracy and precision.
Calibration can be a major portion of the effort required to complete the modeling phase of a project. The level of calibration required for a particular model depends on the type and amount of data available in combination with project needs. Hill and Tiedemann (2007) present suggested guidelines for effective model calibration along with a description of the calibration process. Vrugt et al. (2008) review the state of the science with respect to inverse modeling of subsurface flow and transport properties.
A number of computer programs have been developed to automate the calibration process for particular hydrogeological models. More recently, modelindependent inversemodeling programs have been developed that can be applied to a broad range of forward models. Two such programs that have been widely accepted are UCODE (Poeter et al, 2005) and PEST (Doherty, 2004). Both UCODE and PEST have been incorporated into several GUIs to speed the modelcalibration process.
Gas Transport Modeling
Introduction
Gas transport, particularly the transport of oxygen into unsaturated wasterock piles, can be an important process affecting the generation of ARD. Principal modes of oxygen transport include diffusion and advection. Wels et al. (2003) provide a comprehensive overview of the role of gas transport in ARD generation and methods that can be used to model gas transport.
Data Needs
Data required to model gas transport are similar to the data needed for equivalent modeling of water flow and transport in the subsurface. The permeability of the porous media is an important consideration. Because permeability to gas is a function of the degree of saturation of the pore space, moisture content is also important.
Permeability and moisturecontent measurements can be made in the field or the laboratory. Measurements of moisture content are reasonably straightforward using established methodologies. Field measurements of air permeability using pneumatic pumping tests are described by Baehr and Hult (1991) and are conceptually similar to groundwater pumping tests used to determine aquifer characteristics. Stonestrom and Rubin (1989) describe laboratory airpermeability measurements.
Model Selection
Relatively few models have been developed specifically to address gas transport in the subsurface and the application to ARDrelated problems. Modeling the complete set of physical and chemical processes operating within a wasterock pile requires a multiphase code capable of simulating gas and water flow in the unsaturated zone, chemical interactions with the solid matrix, heat generation and transfer, and chemical mass transfer in the liquid and gas phases.
Several generalpurpose, multiphase simulation programs have been developed that could be applied to these types of problems. The TOUGH family of codes (Pruess et al., 2004) was developed at Lawrence Berkeley National Laboratories and has been applied to a wide range of complex, multiphase problems. TOUGHAMD (Lefebvre et al., 2001) is an adaptation of TOUGH to address ARD‑related issues. TOUGHREACT (Xu et al., 2004) was developed as a comprehensive nonisothermal multicomponent reactive fluid flow and geochemical transport simulator to investigate acidmine drainage and other problems.
Groundwater flow models can be adapted to simulate air flow using appropriate transformations of variables and parameter formulations. Massman (1989) shows how groundwater solutions can be modified for gasflow problems. This type of adaptation would not be appropriate to model the most complex multiphase, reactivetransport problems, but may be adequate to address many issues of importance to the prediction of ARD generation.
Considerations Regarding Predictive Modeling of Effluent Quality
As discussed throughout Chapter 5 of the GARD Guide and this Appendix, the principal objective of mine and process water quality prediction is to evaluate the potential for geologic materials and mine and process wastes to generate acid and contaminants and affect water resources. As an important corollary, the need for and nature of mitigation measures is determined through prediction.
The nature and sophistication of the prediction effort may vary depending on the desired outcome. A prediction exercise aimed at merely answering a “yes/no” question (for instance: will the water quality criterion for arsenic be exceeded?) requires less upfront understanding of the system being evaluated, in which case while use of relatively “crude” modeling tools may suffice. In contrast, when a more quantitative answer is required (for instance: what is the expected arsenic concentration?), the complexity of the modeling effort may be quite significant, requiring both a detailed conceptualization of the system being modeled as well as use of advanced modeling codes.
It should be noted that use of more sophisticated tools does not necessary equate to more accurate and precise modeling outcomes. According to Oreskes (2000) and Nordstrom (2004), the computational abilities of codes and advanced computers currently far exceed the ability of hydrogeologists and geochemists to represent the physical, chemical and biological properties of the system at hand or to verify the model results. In light of these difficulties, the meaning of “accuracy” and “precision” in the context of mine and process water quality modeling must be reassessed on a casebycase basis, and numeric analysis needs to be conducted to reflect the uncertainty inherent in predictive modeling. Accordingly, USEPA (2003) recommends the following should be submitted at a minimum to substantiate modeling used for regulatory purposes, regardless of the specific model/code being used:
 Description of the model, its basis, and why it is appropriate for the particular use
 Identification of all input parameters and assumptions, including discussion of parameter derivation (i.e., by measurement, calculation or assumption)
 Discussion of uncertainty
 Sensitivity analysis of important input parameters
Having said all that, despite the limitations identified throughout this chapter, modeling and prediction have significant value as management tools and for gaining an understanding of the geochemical, physical and biological systems at mine and process sites (Oreskes, 2000). There is little doubt that the understanding of geologic materials, mine and process wastes and the hydrogeochemical factors that govern mine and process water quality will continue to advance through the implementation of laboratory and field experiments. In particular those experiments that isolate one variable at a time to identify its effect on overall discharge water quality will prove of great value. Similarly, ongoing characterization and monitoring of mine and process facilities will allow development of improved scaling factors needed to extrapolate results from smallerscale tests to an operational level. Lastly, the tools required for geochemical, hydrological and hydrogeological modeling already largely exist. With an increased comprehension of the factors that govern the generation of ARD, NMD or SD, modeling and prediction efforts will become increasingly reliable.