Skip to content

Working with Solute Transport

1. Overview

The complete MIKE SHE advection-dispersion (AD) module is comprised of four independent components, each describing the transport processes in one of the parts of the hydrological cycle. Used in combination they describe solute transport in the entire hydrological cycle. The four components are:

  • Overland Transport
  • Channel Transport (MIKE 1D)
  • Unsaturated Zone Transport
  • Groundwater Transport

A number of processes relevant for simulating reactive solute transport are included in MIKE SHE including

  • Water and solute transport in macro pores,
  • Sorption of solutes described by either equilibrium sorption isotherms (Linear, Freundlich or Langmuir) or kinetic sorption isotherms, which include effects of hysteresis in the sorption process,
  • Attenuation of solutes described by an exponential decay, and
  • Plant uptake of solutes.

a. Current Limitations

The solute transport module in MIKE SHE currently does not support

  • any solute transfer via irrigation,
  • solute migration from UZ to OL.

In the first case, the solutes will remain in the source location and only the water will be transferred. This will lead to increasing concentrations at the source.

In the second case, there is no mechanism in MIKE SHE to transfer water from UZ to OL, so there is also no means to move solutes from the UZ cells onto the ground surface. This has implications for salinity modelling, as there is no way for runoff to remove surface salts that migrate upwards due to capillarity and concentrate on the ground surface due to evaporation.

2. Flow Storing Requirements

Solute transport calculations in MIKE SHE AD are based on the water fluxes from a MIKE SHE Water Movement (WM) simulation. To ensure that all the needed WM result data types are stored, you have to specify that results should be stored for an AD simulation. See Storing of Results.

The WM data should be stored frequently enough to describe the dynamics of the flow. The selected storing frequencies of flow results will usually be a compromise between limitations in disk space and resolution of the flow dynamics. The maximum computational time steps in a transport simulation are often restricted by advective and dispersive stability criterions. However, the transport time step cannot be greater than the flow storing time step in each component.

3. Storing of Results

The simulated concentration distribution in each component as well as the mass balances and fluxes will be stored in dfs2 and dfs3 files with different time steps. Besides these result files, the program also writes output to the error log, which describes errors encountered during execution and a print log which contains execution step information, statistics on the run and a mass balance (if requested).

Normally, the results from the saturated zone (species concentration in each grid) is by far the most disk consuming parameter. So, be careful with the storing time step. Mass balances, which includes time series of mass storage and fluxes between components (and sources, drains, boundaries, etc.) can be stored at smaller time steps.

When you select the time step you should also be aware of the time scale of the process. The time scale for transport processes in groundwater is usually much larger than the time scale for transport in a river.

Enter the desired time steps - notice that the unit is hours - in each of the edit fields. There are no limitations on this time step but if you select a time step less than the simulation time step, the storing time step will be the new simulation time step.

4. Simulation and Time Step Control

Simulation time steps are to some extent controlled by the user. Several possibilities for time step control exist to make the execution as fast as possible with no numerical dispersion and instabilities.

The first possibility for controlling the simulation time steps in the different components is simply to define the maximum time step in each component. Note that time steps should be given in increasing order i.e. \(dt_{RIVER} \leq dt_{OVERLAND} \leq dt_{UZ} \leq dt_{SZ}\). Also note that this is the MAXIMUM time step. That is, the actual simulation time step is controlled by the stability criterions with respect to advective and dispersive transport given below. Furthermore, time steps for transport cannot exceed the storing time step for the relevant data in the flow result file from a MIKE SHE flow simulation.

Enter the maximum allowable Courant number for each component. The Courant number is defined by \(V x dt/dx\) (velocity times time step divided by “grid size”). This number should normally not exceed 1.0 for one- and two- dimensional transport (UZ, Overland and Channel Flow) and 0.8 for three- dimensional transport (SZ). The maximum time step will be calculated accordingly.

Enter the maximum allowable dispersive Courant number for each component. The dispersive Courant number is defined by \(D x dt/dx2\) (Dispersion coefficient times time step divided by “grid size” squared). This number should normally not exceed 0.5. The maximum time step will be calculated accordingly.

The transport limits are used to avoid negative concentrations in cases with extreme gradients (e.g. close to sources) or in areas with highly irregular velocity fields. Enter the maximum allowable transport from a node or grid as a fraction of the storage in the node or grid. A recommended value for all components is 0.9, which ensures that this option is in use (the value 0 determines that this option is not in use).

a. Calibrating and verifying the model

The advection-dispersion of solutes depends largely on the simulated flows and fluxes calculated by the MIKE SHE flow model. After your first AD simulations, you will usually have to go back and improve the calibration of your flow model. Rarely, can the simulated concentrations and mass fluxes be calibrated to the measured concentrations by tuning only the solute transport model.

It is important to recognise that a transport model must be calibrated. This is true for all applications larger than the laboratory scale since model output cannot necessarily be compared directly to measured values. Measurements are mostly point measurements at a certain time whereas results often are mean values over larger volumes and longer times.

The purpose of the calibration is to tune the model so that it is able to reproduce measured conditions for a particular period in a satisfactory way. This period known as the calibration period - should be chosen long enough to include events of similar kind as the ones you are going to investigate.

A satisfactory calibration is reached when the model is able to reproduce the measured values taking the following conditions into account: - uncertainty in the measurements (time, space, equipment) - representativeness of measurements (point/average grid values) - differences between your conceptual model and nature - uncertainty in other model parameters and data (source description etc.)

In general, it is impossible to specify an exact level of divergence between measured data and computed results before the model is satisfactorily calibrated. In each application you have to consider all factors influencing your result.

After the calibration, you should verify your model by running one or more simulations for which measurements are available without changing your model parameters. If the model is able to reproduce the validation measurements you can consider your calibration to be successful. This ensures that simulations can be made for any period similar to the calibration and the verification period with satisfactory results.

5. Executing MIKE SHE WQ

In the top icon bar, there is a three-button set of icons for running your model.

.

WQ - The WQ button starts the Water Quality simulation. After you have successfully run a water movement (WM) simulation to completion, you can run a water quality simulation.

In addition to the three icon buttons, there is a Run menu. In this menu, you can check on and off all three of the above options. Finally, there is an Execute... menu sub-item that runs only the checked items above it. The Execute option can also be launched using the Alt - R - E hot-key sequence.

MIKE SHE WQ can also be launched from a batch file with or without the MZLaunch function. For more information on this, see Using Batch Files

6. Output

The output of the MIKE SHE AD is stored to several dfs0, dfs2 and dfs3 files which can be viewed and processed with the different tools available for these files in MIKE ZERO.

For each species, a concentration file is created for each hydrologic process - a dfs2 file for OL, and a dfs3 file each for UZ and SZ.

For each species, the total WQ mass balance is stored in two dfs0 files. The xx_species_AllItems.dfs0 includes all of the possible water quality mass balance items. The xx_species.dfs0 is a copy of the _AllItems.dfs0 file with all of the non-zero items removed.

The first 20 items in the _AllItems.dfs0 file define the global mass balance (see Table 38.2). There are four items: Storage, Input, Output and Error. Each of these is calculated for each of the five storage items: SZ, Immob(SZ), UZ, MP(UZ), and OL.

The item Immob(SZ) is for solutes stored in the SZ matrix porosity when the dual porosity option in SZ is turned on. In this case, water and solutes move in the fractures and solutes diffuse into the rock matrix. The fractures are then the primary porosity.

The item MP(UZ) is for solutes stored in the UZ macropores when the macropore option in UZ is turned on. In this case, water and solutes move in both the matrix and the macropores, but the volume of water in the macropores is generally much less than the volume of water in the matrix. So, the primary porosity is the UZ matrix.

The Error is calculated for each of the five items as:

\(Error = Output - Input + \Delta Storage\)

However, the mass in the SZ and UZ items includes the mass in both the primary and secondary porosities.

The Output, Input and Storage are all displayed as positive values in the dfs0 file and the WQ log file. A positive change in storage denotes an increase in mass.

Table 38.1 WQ mass balance items in the _AllItems.dfs0 file

Mass balance item Component
Storage SZ, Immob(SZ), UZ, MP(UZ), OL
Input SZ, Immob(SZ), UZ, MP(UZ), OL
Output SZ, Immob(SZ), UZ, MP(UZ), OL
Error SZ, Immob(SZ), UZ, MP(UZ), OL

The rest of the items in the _AllItems.dfs0 file are only non-zero if the item is relevant for the current WQ simulation. Table 38.2 lists all of the rest of the items in the _AllItems.dfs0 organized by the source of the solute.

Table 38.2 Available WQ mass balance items in the _AllItems.dfs0 file

From To Comment
Sources SZ, UZ, OL, River Note that sources can be specified as positive or negative.
Ext. Sources (OpenMI) SZ, UZ, OL This is non-zero only if a model is linked to MIKE SHE by OpenMI that adds mass to the component.
Ext. Input (OpenMI) SZ Drain This is non-zero only if a model is linked to MIKE SHE by OpenMI that adds mass to the component. However, this is a special case because you can add mass directly to the SZ drain without it actually becoming part of the SZ model. The mass is then added to the model at the location were the drain discharges (i.e. river link, SZ boundary, or local SZ depression)
Boundary SZ, UZ, OL, River The (Boundary to River) item is typically zero, but can be non-zero in a couple of rare cases: If you have a fixed head SZ boundary, or a fixed concentration boundary, next to a river link, then mass from this cell to the river will be added to (Boundary to River).
Precip UZ, OL Mass from precipitation is always added to either OL or UZ, even though precipitation can be added to SZ in the WM module. If you have a SZ-only simulation, then mass from precipitation is included in the (Precip to OL) and then (OL to SZ)
Decay SZ, Immob(SZ), UZ, MP(UZ), OL This is the mass that has decayed in each of the processes
Sorp/DeSorp SZ, Immob(SZ), UZ, MP(UZ), OL This is the net Sorption and Desorption to the soil matrix in each of the processes. If mass is sorbed to the soil matrix it is removed from solution and this value will be negative. If mass desorbs from the soil matrix it is added to solution and this value will be positive.
Colloid-Sorp/DeSorp SZ, Immob(SZ), UZ, MP(UZ), OL This is the net Sorption and Desorption to colloids in each of the processes. If mass is sorbed to the colloids it is removed from solution and this value will be negative. If mass desorbs from the colloids it is added to solution and this value will be positive. However, this is normally zero, because colloid transport is not available in the commercial version of MIKE SHE, only a research version.
EcoLab SZ, Immob(SZ), UZ, MP(UZ), OL This is the mass change resulting from passing solutes to and from MIKE ECO Lab - positive if MIKE ECO Lab causes the mass to increase and negative if mass decreases.
SZ UZ, MP(UZ), OL, Sinks, Sources, Ext.Sinks(OpenMI), Decay, Sorp/DeSorp, Colloid-Sorp/DeS- orp, EcoLab, Plant Uptake Note that there is no mass transfer to SZ drains because the SZ drains have not storage. Mass that discharges to SZ drains passes straight through the drain and is added to the end recipient (i.e. a river, boundary, or local depression.
Immob(SZ) UZ, Decay, Sorp/DeSorp, Colloid- Sorp/DeSorp, EcoLab
SZ Baseflow River This is a special item that includes only the mass from SZ to rivers.
SZ Drain River, SZ (Local Dep), Boundary This is a special item that divides up the SZ mass discharge to drains by the end recipient.
SZ Flow Boundary This is a special item to distinguish between SZ mass discharge to boundaries via drains and direct discharge to the boundary.
SZ(Fract) Immob(SZ) This is a special item that includes only the mass exchange between fractures and the matrix when the dual porosity option is selected in the SZ.
UZ SZ, Immob(SZ), OL, Sinks, Sources, Ext.Sinks(OpenMI), Boundary, Decay, Sorp/DeSorp, Colloid-Sorp/DeSorp, EcoLab, Plant Uptake Note that the (UZ to Boundary) item refers to mass discharge from UZ to SZ boundaries. This can arise, for example, when a UZ column discharges into an SZ cell that contains an internal boundary condition, such as a constant head. Note that the (UZ to Immob(SZ)) item will only be non-zero when the dual porosity option in SZ is turned on. In this case, as the water table increases, there will be a transfer of UZ mass to water in both the SZ fractures and SZ immobile matrix based on the ratios of their porosities. Related to the above, if macropores are active, then mass in the UZ macropores will be distributed to both the SZ matrix and the SZ fractures.
MP(UZ) SZ, Immob(SZ), UZ(Matr), Decay, Sorp/DeSorp, Colloid-Sorp/DeS- orp, EcoLab
UZ(Matr) MP(UZ) This is a special item that includes only the mass exchange between macropores and the matrix when the macropore option is selected in the UZ.
OL SZ, UZ, MP(UZ), River, Sinks, Sources, Ext.Sinks(OpenMI), Boundary, Decay, Sorp/DeSorp, Colloid-Sorp/DeSorp, EcoLab
River SZ, OL, Boundary (River to OL) is always zero because this exchange has not yet been implemented.

7. Coupling Water Quality in MIKE SHE and Rivers

Detailed information on the MIKE 1D water quality modules are found in the MIKE+ documentation.

The coupling between MIKE 1D and the rest of MIKE SHE’s hydrologic processes is relatively automatic. You must set up a MIKE 1D AD model independently of MIKE SHE using MIKE+ and specify the .m1dx river setup file in the Rivers and Lakes dialogue. This river setup file must only have the same network geometry as the WM river model setup file. It does not have to be the same river setup.

The MIKE 1D AD model can also include MIKE ECO Lab, which will allow you simulate eutrophication, etc. in the surface water.

There are a few caveats/limitations that you need to be aware of: - Species names must be identical in MIKE SHE and MIKE 1D. If they are not identical, then the solutes will be transferred to the river as an infinite sink, but will not be transported in MIKE 1D. - The overland WQ must be included if you want to simulate water quality coupled to MIKE 1D. - Recycling of WM results is not supported in MIKE 1D. This means that if you want to simulate the coupling between MIKE 1D and the rest of MIKE SHE, your WQ simulation must be continuous.