General Questions on RT3D
- What is RT3D? What does it do?
- Why not just use a screening tool like Bioscreen or Biochlor?
- Why is RT3D finite difference? Why not use a finite element code?
- Can I do Monte Carlo or sensitivity analysis using RT3D?
- Can RT3D handle multiple phases (water, solids, gas, NAPL)? Will it in the future?
- I want to sell RT3D. / I want to modify RT3D and sell it. Is that okay?
Getting Help with RT3D
- I need help on RT3D. Where do I start?
- I have read all the documentation, looked around at other resources, but I'm at my wits end trying to solve my problem. What can I do? Who can I contact?
- Can you help me for free?
- You have to help me. There's a bug in your RT3D code.
- Are there training courses available?
- Why do I have to register to download RT3D and what do you do with my information?
- What platforms will RT3D run on? (Mac, DOS, Win32, UNIX varieties)
- How do I install RT3D? I don't see a Setup.exe file.
- What is the U.S. export control classification for RT3D?
- Why do I get an error (e.g., "FLOW-TRANSPORT LINK FILE SAVED BY VER 1 OF [LKMT] PACKAGE") about the link between MODFLOW and RT3D?
- What MODFLOW file does RT3D need in order to run?
- I want to use MODFLOW 200X. Do I need to do anything special for this to work with RT3D?
- Can I use a code other than MODFLOW?
- Must the flow be steady state?
- I get an error about "MATRIX IS SEVERELY NON-DAGONALLY DOMINANT. CHECK INPUT FILES. STOP EXECUTION." What is wrong with my Modflow simulation?
Setting Up a Simulation Model / Selecting Options
- What units can I use to get correct answers? (i.e., Can I define a grid in feet and concentrations in µg/L?)
- How do I do transient simulations?
- When should I use the finite difference (FD) advection solver versus method of characteristics (MOC, MMOC, or HMOC) or the third-order total variation diminishing (TVD) solvers?
- What should the value of the [longitudinal, transverse, vertical] dispersivity be?
- How do I represent a non-aqueous phase liquid (NAPL) source?
Reaction Modules / Kinetics
- Is bioremediation the only thing that RT3D is useful for?
- Is it sufficient to just use the pre-programmed reaction modules as is?
- What half-lives/rate constants/utilization factors do I use?
- Are first-order kinetics suitable/appropriate?
- Which reaction solver should I use?
- What values of RTOL and ATOL should I use?
- How do I use my own kinetics (i.e., in a different way than entering site-specific rate constants into the pre-programmed reaction modules)?
- Can RT3D model zones of differing microbial activity?
- What is rt3dbat1 and how do I use it?
- Why is my transport step size so small?
- My simulation run is taking too long. How can I make my run go faster? How do I make the transport step size bigger?
- My simulation crashed and said that the "Reaction solver failed." What can I do?
- I selected the option to use reaction solver #4 (solver based on a semi-implicit extrapolation method) [or #5 (Non-stiff Runge-Kutta solver)], but the simulation doesn't run and I get a message "Solver 4 and 5 are not available..." What's up with that?
Post-Processing of Output Data
- How can I get RT3D to output UCN files (MT3DMS standard files)?
- How can I get RT3D output files in ASCII text format (i.e., for import into my [plotting, spreadsheet, database] program)?
Advanced Operations & Customization
- Can you [debug, provide a reaction module, provide rate constants, etc.] for me for free?
- What compilers can I use with the RT3D code?
- What about using RT3D in a parallel processing environment (multiple processors on one machine or a cluster of machines)?
- How do I simulate recirculation of water between an injection/extraction well pair?
General Questions on RT3D
What is RT3D? What does it do?
RT3D (Reactive Transport in 3-Dimensions) is a software code for simulating three-dimensional, multi-species, reactive transport of chemical compounds in groundwater. RT3D provides an easy-to-use and flexible framework that you can apply to natural attenuation, accelerated bioremediation, or other reactive transport modeling scenarios. Predefined modules are available for common bioremediation kinetics, but you have the flexibility to add any reaction kinetics desired/suitable to represent multiple chemical species in aqueous and sorbed phases. Several graphical user interfaces to RT3D are available to organize the process of defining the required RT3D input information and to provide visualization capabilities for the RT3D output. At the same time, the availability of the Fortran 90 source code provides flexibility for advanced users, regardless of the target computing platform.
Why not just use a screening tool like Bioscreen or Biochlor?
Bioscreen (for hydrocarbons) and Biochlor (for chlorinated ethenes) are simple spreadsheet models to perform screening assessments of natural attenuation using no biodegradation or simple biodegradation models. These tools can be effective for 1) an up-front screening assessment or 2) "simple" sites that lack complexities in flow, geology, and reaction. RT3D on the other hand, can be used to model the complexities of flow, geology, and reaction for natural attenuation or accelerated in situ bioremediation scenarios. Although targeted at monitored natural attenuation, the documentation for the Scenarios Evaluation Tool Truex et al., 2006 discusses assessment of when to use a conceptual model, screening (analytical) model, or a numerical model.
- 2006 Scenarios Evaluation Tool for Chlorinated Solvent Monitored Natural Attenuation WSRC-STI-2006-00096 Savannah River National Laboratory, Washington Savannah River Company, Aiken, South Carolina
Why is RT3D finite difference? Why not use a finite element code?
The short answer is that RT3D (and the MT3D family of codes) is finite difference because they use the finite difference flow solution from MODFLOW. MODFLOW was selected because it is the U.S. Geological Survey code, it is very popular/commonly used, and it is a public domain code.
Either finite difference or finite element can be used to effectively solve the equations that represent contaminant transport. Finite difference (FD) methods are simple to use and easy to understand. FD methods are easy to implement for rectangular regions (2D and 3D), but are more difficult for complicated geometries. Many fast FD solvers and packages (e.g., FFT) exist for regular domains (rectangular and circular regions). Finite element (FE) methods work well for problems with complicated boundaries or for elliptic PDEs. FE methods require careful selection of the triangulation mesh, but this is made easier by the availability of numerous commercial software packages.
The following three references provide detailed discussions/explanations of FD and FE methods in the context of groundwater modelling.
- 1984 Comparison of Finite Difference and Finite Element Methods In: Fundamentals of Transport Phenomena in Porous Media (NATO ASI Series, Series E, Applied Sciences, Vol. 82) J. Bear and M.Y. Corapcioglu, editors Martinus Nijhoff Publishers, Boston, Massachusetts pp. 899-952
- 1983 Computational Methods in Subsurface Flow Academic Press, New York
- 1982 Introduction to Groundwater Modeling: Finite Difference and Finite Element Methods W H Freeman & Co., San Francisco, California
Can I do Monte Carlo or sensitivity analysis using RT3D?
There is no automated feature to perform a Monte Carlo or sensitivity analysis with RT3D. However, it would be straightforward for the user to set up a series of simulations for such an analysis.
Can RT3D handle multiple phases (water, solids, gas, NAPL)? Will it in the future?
RT3D is not a multi-phase code. RT3D explicitly tracks concentrations in the groundwater and implicitly tracks chemical species sorbed to the soil phase. Non-aqueous phase liquids (NAPLs) are represented as source cells, but fluid flow of the NAPL is not represented. Pre-programmed reaction modules 4 and 5 or a user-defined reaction module could be used to explicitly track contaminants not in the mobile (groundwater) phase. There are no plans to convert RT3D to a multi-phase code.
I want to sell RT3D. / I want to modify RT3D and sell it. Is that okay?
LICENSE AND COPYRIGHT INFORMATION
Like any other literary work, computer programs are protected by copyright. Unauthorized reproduction or distribution of this computer code, or any portion of it, is not permitted. Users can, however, modify the code and use it for solving a specific problem. Users are not permitted to resell or redistribute a the RT3D code or a modified version thereof (or any portion of RT3D) directly or via computer bulletin boards, web pages, or other publicly accessible archives. If you have modified the code for any research application and would like to share the code with other researchers, please contact the developers first and provide details of your work. In all your publications, please cite the RT3D manual (PNNL-11720) and the 1998 Groundwater Monitoring and Remediation journal paper by Clement et al.. Also, consider citing other published RT3D papers that are appropriate. If you have any questions or comments, please contact the developers.
Getting Help with RT3D
I need help on RT3D. Where do I start?
The first places to look are the RT3D documentation and the RT3D FAQ. The RT3D documentation (manuals, release notes/update notes, tutorials, white papers, and FAQ) and accompanying examples discuss the basis for RT3D and how to use the code in detail. The MT3DMS manual Zheng and Wang, 1999 is also a good resource for understanding RT3D because both codes share a common basis. The publications relating to RT3D provide another resource for understanding the principles and application of RT3D to reactive transport modelling. The U.S. Geological Survey (USGS) MODFLOW web site has a number of reports on the use and functionality of MODFLOW (particularly TWRI 6-A1, Open-File Report 96-485, and Open-File Report 01-82.
- 1999 MT3DMS: A Modular Three-Dimensional Multispecies Transport Model for Simulation of Advection, Dispersion, and Chemical Reactions of Contaminants in Groundwater Systems; Documentation and User's Guide SERDP-99-1 United States Army Corps of Engineers, Engineer Research and Development Center, Vicksburg, Mississippi
I have read all the documentation, looked around at other resources, but I'm at my wits end trying to solve my problem. What can I do? Who can I contact?
If you have exhausted the above-mentioned resources, you have two options for requesting help. If you are using a graphical user interface (GUI), then contact the GUI vendor for help. Otherwise, see the contact information for the RT3D team and we'll see what we can do.
Can you help me for free?
Generally speaking, No. We can answer simple questions (or point you to the resources mentioned above) and will address issues that are actual "bugs" with the code. However, like most folks, we need funding to pay for our time. We can set up a contract with a scope to fit your needs on a project. Tasks that we could perform include simply getting your model set up initially, developing custom reaction modules, microcosm studies to define reaction kinetics, customization of RT3D, doing a full natural attenuation evaluation, or designing an accelerated in situ bioremediation system.
You have to help me. There's a bug in your RT3D code.
Some 95% or more of the "bugs" that people claim to have found are actually problems with the setup of the model input files. RT3D has been extensively tested and verified, so we are quite confident in it. Certainly, as upgrades are made to the code there is a potential for errors to occur and we welcome feedback so we can remedy any such problems. We cannot address problems caused by the use of graphical user interfaces other than to pass along the information to the GUI vendor (which you could probably do better than us).
Are there training courses available?
Training notices are posted on the Solutions page. Custom training may be an option; contact us for further information. Additionally, the graphical user interface vendors and the NGWA often offer training courses, though they are typically geared towards working with their product and groundwater modelling in general (versus focusing on RT3D).
Why do I have to register to download RT3D and what do you do with my information?
Registration provides a method to insure that you understand the terms of downloading RT3D and provides us with important information. This information is used for two purposes: 1) to send notices of RT3D updates via e-mail (on the order of one e-mail per year) and 2) to collect demographics information. The latter information is important for the RT3D developers in obtaining the resources required to improve RT3D. Individual user information will not be released; statistics on collective user demographics (e.g., the number of users in California) may be released to the public.
What platforms will RT3D run on? (Mac, DOS, Win32, UNIX varieties)
Binary executables for RT3D are available for the Intel x86/Win32 platform (compiled with Intel or Compaq Visual Fortran). Standard Fortran source code is available for users to compile on any platform with a Fortran compiler.
How do I install RT3D? I don't see a Setup.exe file.
RT3D is a simple executable and doesn't require a setup program. Here are the recommended steps for "installation" of RT3D. First download the RT3D zip file from the RT3D Web Site to a convenient directory on your hard drive (e.g., c:\rt3d). Unzip this file, including the directory structure (see the help file for your favorite unzipping software to learn how to do this). The RT3D executables are in the "<path>\RT3D_and_MF_exe" directory (where <path> is the directory you selected in the first step). Add this directory to your search path; on most newer windows systems go to Control Panel -- System applet -- Advanced tab -- Environment Variables button, then edit the "Path" system variable to include the RT3D directory (i.e., add ";c:\rt3d\RT3D_and_MF_exe;" to the end of the Path string). After applying this change you can go to a command prompt (i.e., "DOS window") and type "rt3dv25.exe" to start RT3D. Of course, you need to execute this RT3D command from a directory that contains a valid set of RT3D input files and the MODFLOW flow-transport link file.
If you are using a groundwater modelling software package, then see that software's help file for instructions on how to set the RT3D executable. Most such packages either use a predefined location or have the ability to pick the directory & executable file to use.
What is the U.S. export control classification for RT3D?
The RT3D software has been classified as Export Control Classification Number EAR99 by the United States Department of Commerce, Bureau of Industry and Security. This software may be exported without a license (i.e., "NLR") PROVIDED THAT the recipient meets all regulations governing export. To meet the Export Administration Regulations, no recipient may 1) be a national of or located within designated (i.e., embargoed/terrorist) countries, 2) appear on the Denied Persons List, 3) appear on the Unverified List, 4) appear on the Entity List, 5) appear on the Specially Designated Nationals List, or 6) appear on the Debarred List. All lists are published by the Bureau of Industry and Security and can be accessed at http://www.bis.doc.gov/complianceandenforcement/liststocheck.htm.
This software is exported from the United States in accordance with the Export Administration Regulations. Diversion contrary to U.S. law is prohibited; that is, the software may not be further distributed to third or subsequent parties who are ineligible to receive exports based on the Export Administration Regulations and the associated Bureau of Industry and Security lists.
A notice stating that you are inelibible for downloading RT3D indicates that you do not satisfy the Export Administration Regulations. If you feel you have received this notice in error, please contact us and relate appropriate details (i.e., name, citizenship, etc.).
Why do I get an error (e.g., "FLOW-TRANSPORT LINK FILE SAVED BY VER 1 OF [LKMT] PACKAGE") about the link between MODFLOW and RT3D?
Several things can cause this error. The heart of the problem is that the *.hff (flow-transport link file) cannot be read by RT3D. This may occur for any of the following reasons.
- You are using an old version of Modflow. Update to a newer version of MODFLOW that uses version 2 or above of the LMT package (also known as LMT6 or LKMT).
- There is no *.hff file found. Run your flow model or move the *.hff file to the same directory as your transport model.
- You are using an incorrectly written/formatted *.hff file. If you are using a flow code other than MODFLOW, check the format of the *.hff file that is being written for consistency with the standard header options described in the USGS Open-File Report 01-82. Make sure that the *.hff file was written with a flow code compiled with the same compiler as RT3D. Head and flow files that are written with a Lahey-compiled MODFLOW cannot be read correctly by a Visual Fortran-compiled RT3D.
What MODFLOW file does RT3D need in order to run?
RT3D uses the "Flow/Transport Link" file that is generated by the LMT package (formerly the LKMT package) in MODFLOW. See the USGS Open-File Report 01-82 for details on the LMT package included with Modflow-2000. Depending on the interface used for MODFLOW and/or user preferences, the Flow/Transport link file can have any of a number of file extensions. Some common extensions used are: *.hff, *.flo, *.ftl, and *.mt3.
I want to use MODFLOW 200X. Do I need to do anything special for this to work with RT3D?
For any version of MODFLOW, you must do two things. First, compile MODFLOW and RT3D using the same compiler. The executable for MODFLOW 2000 that is distributed by the U.S. Geological Survey (http://water.usgs.gov/nrp/gwsoftware/modflow.html) is compiled with Lahey Fortran, while the RT3Dv2.5 executable is compiled with Compaq Visual Fortran. Second, in the LMT6 input file for MODFLOW 200X, you must specify the "Standard" header (or use a blank file); see U.S. Geological Survey Open-File Report 01-82 (at the MODFLOW 2000 web site) for more details.
Can I use a code other than MODFLOW?
MODFLOW is the most convenient groundwater flow code to use. However, any block centered finite difference flow code may be used. The code must be able to produce a head and flow file as required by RT3D. See the USGS Open-File Report 01-82 at the USGS MODFLOW web site for details on the format for the head and flow "link" file.
Must the flow be steady state?
No. RT3D can accommodate transient flow conditions. Just be sure to coordinate your time step and stress period information between MODFLOW and RT3D.
I get an error about "MATRIX IS SEVERELY NON-DAGONALLY DOMINANT. CHECK INPUT FILES. STOP EXECUTION." What is wrong with my Modflow simulation?
A likely cause of this error is that the top and bottom elevations overlap, resulting in cells that have a thickness of zero or a negative value. Check the layer thicknesses and make corrections as needed.
Setting Up a Simulation Model / Selecting Options
What units can I use to get correct answers? (i.e., Can I define a grid in feet and concentrations in µg/L?)
Consistency in units is required for two groups of parameters: concentration-related and grid/aquifer-related parameters. Concentrations, bulk density, and sorption coefficients (e.g., Kd) must all have consistent units such as mg/L, mg/L, and L/mg, respectively. This may mean that the bulk density and Kd values are in unconventional units and are values of a much different magnitude that one is used to working with. Additionally, the reaction constants (e.g., Monod half-saturation or inhibition coefficients) must have units that correspond to the concentration units. Whether mass concentration, mole concentration, or either are acceptable will depend on the specific reaction module being used; the built-in reaction modules use mass concentrations. Parameters related to model grid dimensions or aquifer parameters must use consistent units. For example, the hydraulic conductivity, grid cell size, and dispersivity must all be consistent with units such as ft/day, ft, and ft, respectively. The use of mixed units between these two groups (e.g., grid/aquifer parameter units in feet and concentration units in mg/L) will not cause errors in the RT3D concentration outputs (*.con or *.obs files). The mass balance output files (*.mas), however, will have "mass" numbers that are a mixture of the units from both groups and thus may not give the actual mass, rather the values would be proportional to the mass.
Note that the statement on page 33 of the RT3D Manual (PNNL-11720) and the statement on page 10 of the PNNL-15938 document suffer from the same failing -- the statements incorrectly use "mg/L" instead of "mass concentration." The intent was to inform the user that mass concentration (e.g., mg/L, µ/L, kg/m³, lb/ft³) must be used, and not mole concentration (e.g., mol/L, µmol/L, lb-mol/ft³). Using mole concentration will give incorrect results.
Some of the graphical user interface software will manage units for the user. That is, data can be entered in a variety of units and the software will eventually convert the data to the (user-specified) set of consistent units during the process of generating simulation input files.
How do I do transient simulations?
Transient reactive transport simulations can be done with either steady state or transient flow conditions. There are two things that need to be done to set up a transient simulation: define the stress periods (in the *.btn file) and set up the transient values for sources & sinks (in the *.ssm file). The number of stress periods is defined with the 4th item on the 3rd line of the *.btn file. At the end of the *.btn file, repeated sets of records A21-A23 are required (one set for each stress period). The table below summarizes these input data records (see the MT3DMS manual for details). In the case of steady state flow, the values for record A21 can be those that make sense for the transport simulation. However, when the flow model is transient, the values in record A21 (and A22, if used) must correspond to the values used in the flow simulation.
|A21||PERLEN, NSTP, TSMULT||F10, I10, F10||
Length of the stress period, number of flow time steps, and flow time step multiplier for geometric progression.
Length of flow time steps if geometric progression of time steps is not used in flow model.
This record is only used if TSMULT ≤ 0
|A23||DT0, MXSTRN, TTSMULT, TTSMAX||F10, I10, 2F10||
Transport step size, maximum number of transport steps allowed in one flow time step, transport step multiplier for geometric progression, and maximum transport step size.
The last two only apply when the GCG solver is being used; otherwise they are dummy values.
The sources & sinks must all be defined for each stress period in the *.ssm file.
When should I use the finite difference (FD) advection solver versus method of characteristics (MOC, MMOC, or HMOC) or the third-order total variation diminishing (TVD) solvers?
The choice of advection solver is based on several criteria: whether the problem is advection dominated (affecting the amount of numerical dispersion and artificial oscillation), concern over mass conservation, and run time (affected by grid spacing, stability constraints, etc.). The table below summarizes the differences between three of the advection solvers.
|Type of Advection Solver||Eulerian (finite difference)||Lagrangian (particle tracking)||higher order finite difference (with flux limiter)|
|Good for Advection-Dominated Problems?||Maybe (with fine enough grid or large enough physical dispersion)||Yes (better than TVD)||Yes|
|Implicit / Explicit Options||Explicit (upstream only) or Implicit (upstream or central-in-space) solution of advection. |
Non-advection terms are solved in the same manner as advection terms.
|Not applicable for advection terms. |
Non-advection terms can be solved with explicit or implicit FD.
|Explicit solution of advection.|
Non-advection terms can be solved with explicit or implicit FD.
|Comments||For advection-dominated problems, the upstream weighting option is more susceptible to numerical dispersion while the central-in-space weighting option is susceptible to artificial oscillation. Method is computationally efficient.||Good for advection-dominated problems. Relatively slow run times. Higher memory requirements (depending on number of particles tracked).||Not as computationally efficient as FD. Not as effective as MOC at advection-dominated problems. TVD is mass conservative, essentially oscillation free, and has lower memory requirements than MOC.|
One reference (you may find others if you do a literature search) that may also be of use regarding this topic is:
- 2001 A Comparison of Solute-Transport Solution Techniques and Their Effect on Sensitivity Analysis and Inverse Modeling Results Ground Water 392300-307
What should the value of the [longitudinal, transverse, vertical] dispersivity be?
The typical approach to setting dispersivity is to 1) select an initial value based on correlations to tabulated data (or rules of thumb) and 2) adjust the dispersivity as a fitting parameter in the model calibration. The use of tabulated data and correlations based on such data has several drawbacks. First, site-specific heterogeneity is a major influence on the dispersivity and the tabulated data (in most cases) was not measured at your site. Secondly, an examination of the tabulated data will reveal that for a give plume scale, the dispersivity can range over an order of magnitude or more. On the other hand, measurement of dispersivity at your site may be impractical because of the expense and/or time required to collect the data.
A good tabulation of dispersivity versus the field scale (plume scale) is given in Gelhar et al. (1992). A number of papers describe correlations for the longitudinal dispersivity, which are given below. The Xu and Eckstein equations are recommended. Dispersivity in the directions transverse and vertical to the plume axis is generally lower than that in the longitudinal direction. Less data is available on transverse and vertical dispersivity and what is available is mostly of low reliability (Gelhar et al., 1992). Rule of thumb type equations for transverse and vertical dispersivity are given below.
L = field scale; i.e., the length of the plume (m) alpha_L = longitudinal dispersivity (m) alpha_T = transverse dispersivity (m) alpha_Z = vertical dispersivity (m) Log = base 10 logarithm x^y = raise x to the power of y -------------------------------------------------------------------------- alpha_L = 1.2 ⋅ ((Log(L))^2.958) for no weighting of the data points (Eqn. 12b of Xu and Eckstein, 1995) alpha_L = 0.83 ⋅ ((Log(L))^2.414) for the 1:2:3 weighting scheme of low:intermediate:high reliability data points (Eqn. 14b of Xu and Eckstein, 1995) -------------------------------------------------------------------------- alpha_L = 0.0175 ⋅ ((L)^1.46) for L < 100 m, using high reliability data only (Eqn. 22 of Neuman, 1990) alpha_L = 0.32 ⋅ ((L)^0.83) for L ≥ 100 m, using high reliability data only (Eqn. 32 of Neuman, 1990) -------------------------------------------------------------------------- alpha_L = 0.177 ⋅ ((L)^0.728) using field data only (Eqn. 2.7 of Ayra, 1986) -------------------------------------------------------------------------- alpha_L = 0.1 ⋅ L for a stratified granular aquifer (Pickens and Grisak, 1981) -------------------------------------------------------------------------- alpha_T = 0.1 ⋅ alpha_L alpha_T = (1/3) ⋅ alpha_L -------------------------------------------------------------------------- alpha_Z = 0.1 ⋅ alpha_T
- 1986 Dispersion and Reservoir Heterogeneity Ph.D. Dissertation University of Texas, Austin, Texas
- 1992 A Critical Review of Data on Field-Scale Dispersion in Aquifers Water Resources Research 2871955-1974
- 1990 Universal Scaling of Hydraulic Conductivities and Dispersivities in Geologic Media Water Resources Research 2681749-1758
- 2002 Prediction of Dispersivity for Undisturbed Soil Columns from Water Retention Parameters Soil Sci. Soc. Am. J. 663696-701
- 1981 Scale-Dependent Dispersion in a Stratified Granular Aquifer Water Resources Research 1741191-1211
- 1995 Use of Weighted Least Squares Method in Evaluation of the Relationship Between Dispersivity and Field-Scale Ground Water 336905-908
- 1997 Statistical Analysis of the Relationships Between Dispersivity and Other Physical Properties of Porous Media Hydrogeology J. 544-20
There are a number of Internet resources on dispersivity, including online calculators at:
How do I represent a non-aqueous phase liquid (NAPL) source?
NAPL transport itself is not modeled because RT3D does not handle multiple fluid phases. However, NAPLs may be represented as spatially fixed source cells in several ways.
- For a fixed concentration of all species throughout the entire simulation, a constant concentration "boundary condition" (which need not actually be at the grid boundary) may be specified by setting the ICBUND value of a cell to -1. The starting concentration for that cell is then used as the constant concentration.
- For a source concentration of one or more species that is constant during a stress period, but which may vary between stress periods, define a specified concentration source (in the SSM package, ITYPE = -1) in the desired grid cells. The concentration of this source type is defined as part of the point source sink information. Species that are not source constituents and whose concentrations should be allowed to change in the specified concentration cell are indicated by setting their concentration to a negative number in the point source/sink definition. Starting concentrations are ignored for source species and are used for the species flagged as having changing concentration in the cell.
- For a source concentration of one or more species that changes over time according to a first-order decay, define a decaying source (in the SSM package, ITYPE = -2) in the desired grid cells. The concentration of this source type is defined as part of the point source sink information. Species that are not source constituents and whose concentrations are unaffected by the decay calculations in the decaying source cell are indicated by setting their concentration to a negative number in the point source/sink definition. Starting concentrations are ignored for source species in the decaying source cell.
- Another method of representing a source is to define a well in the flow simulation with a very low flow rate and a sufficiently high concentration to represent the mass influx. This is a reasonable approximation for single species transport if the amount of water injected into the aquifer does not alter the flow field. However, for multi-species simulations, additional care should be exercised to maintain the injection rate at a very low level. At discharge points (i.e., at the hypothetical injection wells) the concentrations of species other than the primary contaminant are usually assumed to be zero and an injection flow rate that is too high may introduce artificial dilution of these species. For example, if a TCE source is created using an injection well, then the concentration of TCE biotransformation byproducts (DCE and VC) will be assumed to be zero in the injected water. Maintaining the injection flow rate at a high value would lead to artificial local dilution of the DCE and VC plumes near the injection well. Maintaining low flow rates will help minimize this dilution effect.
Reaction Modules / Kinetics
Is bioremediation the only thing that RT3D is useful for?
No. RT3D is suitable for modelling any type of reactive transport in groundwater. The impetus for developing RT3D was for use in accelerated bioremediation and natural attenuation scenarios. However, RT3D is flexible enough to be used for subsurface flow and transport modelling of other types of processes (e.g., chemical oxidation). The pre-programmed reaction modules are geared towards bioremediation, although they could represent any process with similar properties (i.e., module 6 for sequential first-order decay can represent any sequential decay process with up to 4 species, not just PCE/TCE/DCE/VC). The flexibility of RT3D stems from the ability of the user to define their own reaction module in an easy-to-implement manner. The user can add their kinetic expressions of whatever desired complexity to represent the reactions of interest. Examples of they types of things the user could do with a custom reaction module include: inhibition terms could be used to limit degradation, bacterial population kinetics could be used to track the amount of reaction that will occur, non-equilibrium interactions with solid phase species could be represented, or loss of contaminant mass to the vadose zone could be included.
Is it sufficient to just use the pre-programmed reaction modules as is?
Generally speaking, No. The pre-programmed reaction modules are good models, based on information available in the literature. As with any model, however, the user needs to think about what they are doing and to understand the scope, assumptions, and limitations of the model if they want to obtain meaningful output. It's the old saying of "garbage in = garbage out." In this situation, the equations and basis for the reaction models are solidly explained in the RT3D user manual. Where the user must apply caution is in the values used for the reaction constants. Default values are often provided (depending on the graphical user interface to RT3D that you are using) based on literature information and assumptions about the species involved in the reactions. But these values may not be appropriate for the conditions at your site and may be inaccurate if you are modelling other chemical species. So rather than blindly charging ahead and using the pre-programmed reaction models as is, the user is encouraged to examine the basis of the reaction modules (i.e., read the RT3D manual and the literature citations provided), to do your own literature search, and to use site specific kinetic data (field or laboratory microcosm) where available.
What half-lives/rate constants/utilization factors do I use?
The best rate constants to use would be those derived from laboratory microcosms with site sediment or field tests at the site. If such data is not available, then a literature search may provide a range of values to use as input for your model. Care must be exercised when using literature values because conditions tend to vary widely from site to site. However, literature data can provide an initial starting point and perhaps bounding cases (depending on how much data is available). Initial degradation rates based on literature can then be revised during the model calibration process. Several example references are listed below for a starting point (two are taken from the EPA HWIR Risk Assessment web site).
- 1999 Biodegradation Rates for Fuel Hydrocarbons and Chlorinated Solvents in Groundwater Bioremediation J. 34337-362
- 1999 Anaerobic Biodegradation Rates of Organic Chemicals in Groundwater: A Summary of Field and Laboratory Studies U.S. Environmental Protection Agency, Office of Solid Waste, Washington, D.C PDF File Note: This is a follow-on report to Aronson & Howard 
- 1998a Aerobic Biodegradation of Organic Chemicals in Environmental Media: A Summary of Field and Laboratory Studies U.S. Environmental Protection Agency, Office of Research and Development, Athens, Georgia PDF File
- 1997 Anaerobic Biodegradation of Organic Chemicals in Groundwater: A Summary of Field and Laboratory Studies SRC TR-97-0223F Syracuse Research Corporation, North Syracuse, New York PDF File
- 2005 Impact of Cold Temperatures on Biodegradation Rates for Natural Attenuation of Petroleum Hydrocarbons Publication No. 809 Alberta Environment, Edmonton, Alberta PDF File
There are some default rates for the pre-programmed reaction modules in the documentation or inserted by the graphical user interface (depending on which GUI you are using). These defaults are typical values, but may not be appropriate for your specific site.
Are first-order kinetics suitable/appropriate?
Some previous studies have shown that first-order kinetics can be reasonable for modeling dechlorination mechanisms at field scales U.S. EPA, 1998b; Clement et al., 2000 and at low concentrations Schmidt et al., 1985. Yet other studies Scow et al., 1986 have determined that first-order kinetics were not appropriate. A number of factors (substrate concentration, microbial population density, etc.) can impact the type of model (first-order, Monod, logistic, etc.) that best represents biodegradation of a chemical Simkins and Alexander, 1984; Alexander, 1985; Schmidt et al., 1985. The first-order assumption should be assessed Bekins et al., 1998 as part of the determination of whether the reaction module is suitable for the intended application.
First-order kinetics express the biodegradation rates as being primarily a function of the concentration of the contaminant. Implicit to the assumption of first-order kinetics is that the microbial population is nongrowing, i.e., the substrate concentration is much less than the Monod half saturation constant. First-order kinetics do not have an explicit dependence on the type or concentration of electron donor, the type or quantity of bacteria present, nutrient limitations, or the geochemical conditions. All of these factors should be evaluated prior to selecting a first-order model to represent biodegradation. Additionally, it may be desireable to modify the first-order degradation kinetics if the contaminants themselves inhibit dechlorination. If a biostimulation or bioaugmentation technique is used at a site (increasing the quantity of subsurface bacteria) then the assumptions associated with first-order kinetics may not be valid and the user should apply a more suitable kinetic model.
- 1985 Biodegradation of Organic Chemicals Environ. Sci. Technol. 182106-111
- 1998 A Comparison of Zero-Order, First-Order, and Monod Biotransformation Models Ground Water 362261-268
- 2000 Natural Attenuation of Chlorinated Ethene Compounds: Model Development and Field-Scale Application at the Dover Site J. Contaminant Hydrology 422-4113-140
- 1985 Models for the Kinetics of Biodegradation of Organic Compounds Not Supporting Growth Appl. Environ. Microbiol. 502323-331
- 1986 Kinetics of Mineralization of Organic Compounds at Low Concentrations in Soil Appl. Environ. Microbiol. 5151028-1035
- 1984 Models for Mineralization Kinetics with the Variables of Substrate Concentration and Population Density Appl. Environ. Microbiol. 4761299-1306
- 1998b Technical Protocol for Evaluating Natural Attenuation of Chlorinated Solvents in Ground Water EPA/600/R-98/128 United States Environmental Protection Agency, Office of Research and Development, Washington D.C
Which reaction solver should I use?
There are two types of reaction solvers that can be used in RT3D. The first type is based on the Runge-Kutta class of ordinary differential equation (ODE) solvers (solvers 3 and 5). Runge-Kutta solvers are more efficient for reaction kinetics that are not stiff (i.e., well behaved - see any text on differential equations for discussion of stiffness). However, for stiff problems, the Runge-Kutta solvers must take very small time steps and may not be able to provide a solution at all. The other type of solver is designed for stiff problems and is either based on the work of Hindmarsh (solvers 1 and 2; Hindmarsh, 1983) or the work of Bader and Deuflhard (solver 4; Bader and Deuflhard, 1983). Solvers 1 and 2 use Adams methods (predictor-corrector) in the non-stiff case and Backward Differentiation Formula (BDF) methods (the Gear methods) in the stiff case, with automatic switching. The difference between solvers 1 and 2 is that the latter uses an analytical Jacobian matrix while the former computes the Jacobian numerically. The Bader and Deuflhard solver is a generalization of the Bulirsch-Stoer algorithm that uses a semi-implicit extrapolation method and (here) requires an analytical Jacobian matrix.
A good place to start is with solver 1 because it adjusts for stiff and non-stiff cases. The Runge-Kutta solvers are preferred when the problem is not stiff.
- 1983 A Semi-Implicit Mid-Point Rule for Stiff Systems of Ordinary Differential Equations Numer. Math. 413373-398
- 1983 ODEPACK, A Systematized Collection of ODE Solvers In: Scientific Computing (vol. 1 of IMACS Transactions on Scientific Computation) R. S. Stepleman et al. (eds.) North-Holland, Amsterdam pp. 55-64 Also published as LLNL Report UCRL-88007, August 1982
Numerous resources pertaining to ODEs exist both in the literature and on the Internet (e.g., http://mathworld.wolfram.com/StiffDifferentialEquation.html).
What values of RTOL and ATOL should I use?
If ISOLVER = 1 or 2, then the atol and rtol values specified for each species are used by the solvers. If ISOLVER = 3, then only the first pair of values, atol(1) and rtol(1), is used by the solver. If ISOLVER = 4 or 5, then only the first absolute error tolerance value, atol(1), is used by the solvers. However, to maintain a general structure for the input data, atol and rtol must be specified for each species, regardless of the choice of reaction solver.
For ISOLVER = 1 or ISOLVER = 2 (i.e., the Gear-based solvers), relative (rtol) and absolute (atol) error tolerances are needed for the solver to determine how accurate the solution for the species concentration, y(i), must be. The differential equation solver will attempt to control the species concentration local error, e(i), such that, for each species, the local error has (approximately) a smaller magnitude than the error weight, ewt(i):
|e(i)| ≤ ewt(i) where, ewt(i) = rtol(i) ⋅ |y(i)| + atol(i)
Here, ewt(i) is a vector of weights which must always be positive, and the values of rtol and atol should all be non-negative. Setting atol(i) = 0.0 results in a pure relative error test on that component (species). Setting rtol(i) = 0.0 results in a pure absolute error test on that component. A mixed test with non-zero rtol(i) and atol(i) corresponds roughly to a relative error test when the solution component, y(i), is much bigger than atol and to an absolute error test when the solution component is smaller than the threshold atol. For practical problems, the following rules of thumb may be used to set atol and rtol values. The atol(i) should be set to a tolerable absolute error level selected by the user (around 10-6 to 10-9). Defining m as the number of significant digits required in the solution for species y(i), set rtol(i) = 10- (m+1). Caution: Actual (global) errors may exceed the local tolerances, so choose atol(i) and rtol(i) conservatively.
For ISOLVER = 3, both relative (rtol) and absolute (atol) tolerances are required. One value for each tolerance is used with the local error and solution vectors to control local error at each step according to:
|e(i)| ≤ rtol ⋅ |y(i)| + atol
The relative and absolute local error tolerances must be non-negative. The relative error tolerance should not be smaller than about 10-8. Pure absolute error control is not permitted; the relative error must be a minimum value, which is machine dependant.
For ISOLVER = 4 or ISOLVER = 5, one global tolerance value, known as the error residue "eps", specified as atol(1), is used to control convergence. For further details about error residue refer to Press et al., 1992, pg. 711. For practical problems, set eps to a small value between 10-6 to 10-9.
- 1992 Numerical Recipes in Fortran 77: The Art of Scientific Computing Cambridge University Press, New York
How do I use my own kinetics (i.e., in a different way than entering site-specific rate constants into the pre-programmed reaction modules)?
Essentially, you need to add the set of differential equations defining the kinetics to the "rxns" subroutine (and the "jacrxns" subroutine if you will be using reaction solvers 2 or 4) and then compile the code. There are two methods for compiling user-defined reactions: build a dynamic link library (DLL) as was done with version 1 of RT3D or compile the whole RT3D program with the user-defined "rxns" and "jacrxns" subroutines. See the RT3D Version 2.5 Update Document and/or the RT3D Manual for details on how to set up a user-defined reaction module.
Can RT3D model zones of differing microbial activity?
Certainly. As an example, consider reaction module 7 (Aerobic/Anaerobic PCE/TCE Dechlorination) where we want to define spatially different zones where anoxic and aerobic activity are occurring. Since this reaction module was constructed to calculate the effects of two types of metabolism, there are corresponding aerobic and anoxic reaction rate constants. In the reaction package input file (*.rct), you would define 9 spatially variable reaction constants (nvrxndata = 9 and ncrxndata = 0). Each reaction constant would be entered as an array with one value for each grid cell. Then, in zones (grid cells) of aerobic activity you would set the aerobic rate constants (e.g., Kv2) to some positive values and the anoxic rate constants (e.g., Kp) to zero. In anoxic grid cells, you would do just the opposite with nonzero anoxic rates and zero for the aerobic rates.
What is rt3dbat1 and how do I use it?
RT3D comes with a batch reaction utility for testing kinetics. This utility, rt3dbat1, is a stripped down version of RT3D version 1.0 that does not do any transport calculations. The utility merely takes initial concentrations and runs the reaction kinetics for a specified length of time, producing an output file with time series data for all species. This utility can be used to troubleshoot user-defined kinetics or to test pre-programmed reaction modules. The reaction kinetics equations, user-defined or pre-programmed, must be set up in the "rxns" subroutine for rt3dbat1 and then the program must be compiled. Currently, the rt3dbat1 code is set up to use a dynamic link library call to rxns.dll, but the code could easily be modified to compile everything at once. Note that rt3dbat1 lacks several of the reaction solver options that are available in RT3D version 2.5.
Why is my transport step size so small?
If you are using the implicit solver (GCG package), then it is small because that's what you told RT3D to do. Set the transport step size to what you want (or what gives you an appropriately accurate solution).
If you are solving the transport terms explicitly, then there are stability constraints that limit the transport step size, even if you told RT3D to use a certain step size. The stability constraints are listed in the MT3DMS manual Zheng and Wang, 1999, page 54 for the advection, dispersion, and source/sink terms (the reaction stability constraint does not apply in RT3D). As the flow velocity, dispersion, or source/sink flux increase, the transport step size is reduced so that the solution is assured to be stable. As the retardation factor (or product of the porosity and the retardation factor) increases (i.e., as there is more equilibrium sorption), the transport step size for a stable solution will increase. Increasing the grid cell size(s) will also increase the maximum allowable transport step size. However, such modifications may not make sense for a particular model or they might decrease solution accuracy to an unacceptable level. You can look at the *.out (output) file to determine what term is limiting the transport step size and at what location. Search the file for "MAXIMUM STEPSIZE WHICH MEETS STABILITY CRITERION" (in 3 locations) and "MAXIMUM STEPSIZE DURING WHICH ANY PARTICLE CANNOT MOVE MORE THAN ONE CELL". These statements identify the maximum step size and the cell indices where the value was calculated. If the simulation is transient, the maximum step size calculations are repeated for each stress period. See the FAQ below for some ideas on how to address small time steps.
My simulation run is taking too long. How can I make my run go faster? How do I make the transport step size bigger?
There are several adjustments that can be made to affect the run time and/or transport step size.
- With RT3D versions 1 and above, you can set the transport step size directly. A transport step size of 0.0 will result in automatic determination of the step size. When the transport step size is set to any positive value, RT3D will use that value unless stability restrictions (for the explicit solver) require a smaller step size. The implicit solver (discussed next) does not have advection stability constraints.
- With version 2.x, you can use an implicit solver (the GCG package). This allows you to set the transport step size to what you want. However, one must be cautious about this. As you set the transport step size larger, your solution will tend to lose accuracy. You should run some tests with a variety of step sizes using the implicit solver and compare the results (and also compare to the explicit solution).
- With version 2.x, there are some additional reaction solver options that may speed up your simulation run (although they won't affect your transport step size). The Runge-Kutta solvers are faster than the Gear solvers. Which solver you can use will depend on how stiff your reaction kinetics are. See the RT3D Version 2.5 Update Document for details on the input file format and options for the reaction solvers.
Additionally, you can look at the model setup (particularly in light of the FAQ above) to see where adjustments may be necessary.
- Look at adjusting the grid size; larger spacing will increase the timestep at the expense of accuracy. Perhaps variable grid spacing would be appropriate.
- Examine the reaction kinetics. Is there a way to make them less complex and/or less stiff? Perhaps adjusting the scale of the concentration values to avoid concentrations of widely different magnitude. As with using a different reaction solver, you will not see a change in transport step size if you change the reaction kinetics, but you may see a change in simulation run time.
- Examine the flow solution. Are there high flux/velocity sources or sinks (e.g., an injection well) that could be adjusted? Do the boundary conditions make sense?
- Is the dispersivity too large?
- Are you using appropriate porosity, bulk density, and sorption parameters? Make sure units correspond to obtain a dimensionless retardation factor.
- Make sure you didn't do something unintentional, like printing out results at a much too frequent interval.
My simulation crashed and said that the "Reaction solver failed." What can I do?
Depending on the failure code, additional information (which differs depending on whether you are using the Gear or Runge-Kutta solver) is presented to guide you in taking corrective action. Failure of the Runge-Kutta solver is generally related to either the choice of the error tolerance or the fact that the reaction kinetics are too stiff. See the RT3D Version 2.5 Update Document for details on the error tolerance. Although the Gear solver has a wider variety of possible problems, they are addressed similarly to Runge-Kutta errors. The user should check tolerance values, confirm that there are not errors in the input data, try smaller time steps, and/or re-examine the reaction kinetics equations for errors.
I selected the option to use reaction solver #4 (solver based on a semi-implicit extrapolation method) [or #5 (Non-stiff Runge-Kutta solver)], but the simulation doesn't run and I get a message "Solver 4 and 5 are not available..." What's up with that?
RT3D reaction solvers 4 and 5 are optional solver routines that cannot be distributed via the Internet because of copyright issues. If you need access to these solvers then please contact the developers.
Post-Processing of Output Data
How can I get RT3D to output UCN files (MT3DMS standard files)?
To obtain MT3DMS standard "Unformatted Concentration" files (in addition to the RT3D *.con output files), the user must uncomment several blocks of code in the BTNRTv25.f file and then recompile RT3D. Uncomment the blocks at lines 311-320, 467-477, and 1065-1071 (line numbers may vary with versions of RT3D other than 2.5, but the blocks should be in the same general location). Additionally, the calls to the "OPENFL" subroutine in the first two blocks must be modified to have only 3 arguments by deleting the last two arguments. These three blocks are shown below with the corrections made to the call to OPENFL.
! -- Block at lines 311-320 -- IF(SAVUCN) THEN WRITE(IOUT,1046) IUCN FLNAME='RT3Dnnn.UCN' DO INDEX=1,NCOMP WRITE(FLNAME(5:7),'(I3.3)') INDEX CALL OPENFL(-(IUCN+INDEX),0,FLNAME) ENDDO ELSE WRITE(IOUT,1047) ENDIF ! -- Block at lines 467-477 -- IF(SAVUCN) THEN FLNAME='RT3D.CNF' CALL OPENFL(ICNF,0,FLNAME) WRITE(ICNF,*) NLAY,NROW,NCOL WRITE(ICNF,*) (DELR(J),J=1,NCOL) WRITE(ICNF,*) (DELC(I),I=1,NROW) WRITE(ICNF,*) ((HTOP(J,I),J=1,NCOL),I=1,NROW) WRITE(ICNF,*) (((DZ(J,I,K),J=1,NCOL),I=1,NROW),K=1,NLAY) WRITE(ICNF,*) CINACT CLOSE(ICNF) ENDIF ! -- Block at lines 1065-1071 -- IF(SAVUCN .AND. PRTOUT) THEN TEXT='CONCENTRATION' DO K=1,NLAY WRITE(IUCN+ICOMP) NTRANS,KSTP,KPER,TIME2,TEXT,NCOL,NROW,K WRITE(IUCN+ICOMP) ((CNEW(J,I,K,ICOMP),J=1,NCOL),I=1,NROW) ENDDO ENDIF
How can I get RT3D output files in ASCII text format (i.e., for import into my [plotting, spreadsheet, database] program)?
There are three ways to convert RT3D output (*.con) files into ASCII text files. You can load the *.con files into your favorite graphical user interface and then use the GUI's facilities for manipulating/saving data. You can write a small Fortran program to read the *.con files and save the data in your desired format. Or you can contact us to request a copy of our conversion utility.
Advanced Operations & Customization
Can you [debug, provide a reaction module, provide rate constants, etc.] for me for free?
No. We can provide some general guidance, but we too need funding to pay for our time. Users are expected to read all the documentation and use literature and related resources to obtain information for running their RT3D model. Bug reports are welcome and will be addressed, however a great percentage of the time such "bugs" are actually problems with the setup of the model input files. If you want to make use our expertise, we could set up a contract to perform an appropriate scope of work for you. See the information about working together to achieve solutions.
What compilers can I use with the RT3D code?
Any Fortran compiler should work to compile the RT3D code. Version 2.5 was tested with Intel & Compaq Visual Fortran on an Intel x96/Win32 platform, as well as with Fortran on DEC Ultrix and Cray platforms. Since we are unable to test all combinations of compilers and platforms, please feel free to let us know if you find that compilation options or source code must be modified for RT3D to compile for your compiler/platform.
What about using RT3D in a parallel processing environment (multiple processors on one machine or a cluster of machines)?
Investigations are ongoing into a version of RT3D capable of parallel execution on a multi-processor (e.g., dual/quad core) machine.
How do I simulate recirculation of water between an injection/extraction well pair?
True recirculation of groundwater from one extraction well to another is not available for the current version of RT3D. Battelle PNWD has a proprietary version of RT3D with this capability and can work with clients to provide this capability.