Section 9 Advection Dispersion

9.1 Overview

The TUFLOW AD Module is an add-on Module and extension of the TUFLOW hydrodynamic engines. By including passive transport fields for tracers, the advection, dispersion and fate of water quality constituents may be simulated. Examples of constituents include salinity, sediment and/or catchment inflow pollutant concentrations. Both dissolved and particulate constituents can be simulated. TUFLOW AD takes depth and velocity fields computed by the TUFLOW Classic or HPC solver and uses this information, together with initial and boundary conditions, to simulate the transport of the specified constituents. In particular, the effects of both longitudinal dispersion and transverse turbulent diffusion are included. TUFLOW AD is specifically oriented towards analyses of systems including coastal waters, estuaries, rivers, floodplains and urban areas. Up to twenty individual constituents can be simulated within TUFLOW AD. For models with linked 1D-2D domains, transport through 1D elements linked by SX connections is supported, but transport through 1D elements linked by HX connections are not.

9.1.1 2D Depth Averaged Equation of Motion

The equation of motion of a passive tracer \(\phi\) (as a concentration), in 2D depth averaged conservative formulation is:

\[\begin{equation} \frac{\partial(h \phi)}{\partial t} + \frac{\partial(h u \phi)}{\partial x} + \frac{\partial(h v \phi)}{\partial y} - \frac{\partial}{\partial x} \left( h D \frac{\partial \phi}{\partial x} \right) - \frac{\partial}{\partial y} \left( h D \frac{\partial \phi}{\partial y} \right) = h S_\phi \tag{9.1} \end{equation}\]

Where:

  • \(\phi\) = Tracer concentration as mass (or mols) per unit volume
  • \(h\) = Water depth
  • \(u\) and \(v\) = Depth averaged velocity components in the x and x directions
  • \(x\) and \(y\) = 2D spatial coordinates
  • \(t\) = Time
  • \(D\) = Isotropic dispersion plus turbulent diffusion coefficients values in the x and y directions
  • \(S_\phi\) = Tracer source as mass (or mols) per unit volume per unit time

The first term in the equation can be expanded using the chain rule and then substituting the continuity equation to get:

\[\begin{equation} \frac{\partial(h\phi)}{\partial t} = h \frac{\partial \phi}{\partial t} + \phi \frac{\partial h}{\partial t} = h \frac{\partial \phi}{\partial t} - \phi \left( \frac{\partial (h u)}{\partial x} + \frac{\partial (h v)}{\partial y} \right) \tag{9.2} \end{equation}\]

Combining Equations (9.1) and (9.2) we get:

\[\begin{equation} \frac{\partial \phi}{\partial t} + u \frac{\partial \phi}{\partial x} + v \frac{\partial \phi}{\partial y} - \frac{1}{h} \left[ \frac{\partial}{\partial x} \left( h D \frac{\partial \phi}{\partial x} \right) + \frac{\partial}{\partial y} \left( h D \frac{\partial \phi}{\partial y} \right) \right] = S_\phi \tag{9.3} \end{equation}\]

9.1.2 Solution Method (TUFLOW Classic)

The 2D advection dispersion algorithm in TUFLOW Classic is based on the ULTIMATE QUICKEST method of Leonard (1991), Leonard & Niknafs (1991) and Leonard et al. (1993). The equation of motion is based on Equation (9.3), but reformulated to make use of the offset grid and the TVD interpolation scheme (i.e. it is better to first compute \(\phi\) at the cell faces and then compute \(\partial(u\phi)/\partial x\) rather than \(\partial\phi/\partial x\)):

\[\begin{equation} \frac{\partial \phi}{\partial t} + \frac{\partial (u \phi)}{\partial x} + \frac{\partial (v \phi)}{\partial y} - \frac{1}{h} \left[ \frac{\partial}{\partial x} \left(D_{x}\frac{\partial\phi}{\partial x}\right) + \frac{\partial}{\partial y} \left(D_{y}\frac{\partial\phi}{\partial y}\right) \right] = S_\phi + \phi \left( \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} \right) \tag{9.4} \end{equation}\]

TUFLOW AD applies the dispersion formulation described by Falconer et al. (2005). This formulation computes dispersion in the X and Y directions (\(D_x\) and \(D_y\) respectively, to suit the Cartesian computational grid) from user specification of longitudinal and transverse dispersion coefficients \(K_L\) and \(K_T\), respectively. Specifically, \(D_x\) and \(D_y\) are computed dynamically at each grid cell and timestep as follows:

X Direction dispersion plus turbulent diffusion:

\[\begin{equation} D_x = max \left( \frac{(K_L u^2 + K_T v^2) h \sqrt g}{|U| C}, D_{w} \right) \tag{9.5} \end{equation}\]

Y Direction dispersion plus turbulent diffusion:

\[\begin{equation} D_y = max \left( \frac{(K_T u^2 + K_L v^2) h \sqrt g}{|U| C}, D_{w} \right) \tag{9.6} \end{equation}\]

Where:

  • \(K_L\) = User specified longitudinal dispersion coefficient
  • \(K_T\) = User specified transverse diffusion coefficient
  • \(g\) = Gravitational acceleration
  • \(|U|\) = Velocity magnitude
  • \(C\) = Chezy bed friction coefficient \(C = \frac{h^\frac{1}{6}}{n}\)
  • \(D_w\) = User specified lower bound dispersion coefficient

The value of \(D_w\) can be specified as constant or spatially varied as required.

The source terms also includes sink terms such as settling (for particulate species) and decay. The scheme also includes representation of mixing due to sub-grid-scale turbulence and vertical shear via the dispersion formulation provided in Falconer et al. (2005).

The TUFLOW AD computational procedure used is an explicit scheme based on Leonard (1991). This contrasts to the TUFLOW Classic engine, which employs an implicit scheme. As such, TUFLOW AD is generally subject to stricter stability constraints than the Classic hydrodynamic engine. As such, the TUFLOW AD calculation takes the form of three steps within each timestep.

The first step involves calculation of the Courant-Friedrichs-Lewy (CFL) condition at all wet cells, where the CFL in 1 dimension is:

CFL Condition:

\[\begin{equation} CFL = \frac{u \bullet \Delta t}{\Delta x} \tag{9.7} \end{equation}\]

Where:

  • \(u\) (or \(v\)) = Fluid velocity
  • \(\Delta x\) = Grid scale
  • \(\Delta t\) = Timestep

This condition is typically required to be less than 1.0 (additive for both X and Y directions) and has a broad physical interpretation requiring that the distance fluid is advected in one timestep (\(u∆t\)) is less than one grid cell (\(∆x\)).

The second step is the computation of a similar condition for the diffusive lengthscale (related to the Peclet number) that ensures that dispersion also does not cause instability at any grid cell. The CFL and dispersion dimensionless numbers are then added and the maximum sum at any given location within a timestep is used to compute the number of sub-stepping iterations required by TUFLOW AD to remain stable within one TUFLOW timestep.

The third step within each TUFLOW timestep is to execute the advection dispersion calculations for the required number of iterations, with a modified (smaller) \(∆t\).

The original ULTIMATE QUICKEST solution method has been enhanced and improved as applied in TUFLOW AD. For example, TUFLOW AD also employs adaptive computational stencil expansion where it identifies sharp constituent concentration gradients (Leonard & Niknafs, 1991). Where possible (i.e. away from boundaries and dry cells) and required, the ULTIMATE QUICKEST stencil is expanded from the standard third order scheme to a ninth order scheme, only on principle computational axes. Cross terms greater than third order are not included. If insufficient wet cells exist to switch to ninth order, then seventh and fifth order schemes are progressively tested (with commensurately decreasing stencils) until all required wet cells are located.

Application of the ULTIMATE limiter (Leonard, 1991) has been found to induce steady flow anisotropy when extended to multi-dimensional problems and the numerical cross terms associated with additional dimensions are included in calculations. Wu & Falconer (2000) developed a modification to the ULTIMATE limiter that reduces this anisotropy, and this has been applied within the TUFLOW AD computational engine.

9.1.3 Solution Method (TUFLOW HPC)

The TUFLOW HPC solver tracks areal density, \(b=h\phi\), of the passive tracer as the primary prognostic variable. With this change, equations (9.1) and (9.2) combine to become:

\[\begin{equation} \frac{\partial b}{\partial t} + \frac{\partial (ub)}{\partial x} + \frac{\partial (vb)}{\partial y} - \frac{\partial \phi_{bu}}{\partial x} - \frac{\partial \phi_{bv}}{\partial y} = S_b \tag{9.8} \end{equation}\]

Where:

  • \(\phi_{bu}\) and \(\phi_{by}\) are the combined dispersive and diffusive unit fluxes of tracer in the x and y directions respectively
  • \(S_b\) = source terms (areal density rate)

The diffusive fluxes of tracer across cell faces are comprised of two components. The first is dispersion (induced by shear in the vertical velocity profile) which acts in the direction of the flow, and the second is turbulent diffusion which acts in both longitudinal and transverse directions. Both calculations use the same formulation, but with different coefficients (Falconer et al., 2005):

\[\begin{equation} \left( \substack{ \phi_{bL} \\ \phi_{bT}} \right) = -\frac{|U|h\sqrt{g}}{C} \left[ \substack{K_L \:\:\: 0 \\ 0 \:\:\: K_T} \right] \left( \substack{ \frac{\partial \phi}{\partial L} \\ \frac{\partial \phi}{\partial T}} \right) \tag{9.9} \end{equation}\]

Where:

  • \(L\) and \(T\) denote longitudinal and transverse directions

The concentration gradients in the longitudinal and transverse directions can be obtained from those in cartesian coordinates with a rotational transformation:

\[\begin{equation} \left( \substack{\frac{\partial \phi}{\partial L} \\ \frac{\partial \phi}{\partial T}} \right) = \left[ \substack{cos \theta \:\:\: -sin \theta \\ sin \theta \:\:\: cos \theta} \right] \left( \substack{\frac{\partial \phi}{\partial x} \\ \frac{\partial \phi}{\partial y}} \right) \tag{9.10} \end{equation}\]

Where

  • \(\theta\) = angle of flow vector such that \(cos \theta = \frac{u}{|U|}\) and \(sin \theta = \frac{v}{|U|}\)

And similarly the dispersive and diffusive fluxes can be transformed back from flow direction coordinates to the cartesian coordinates with the reverse transform:

\[\begin{equation} \left( \substack{ \phi_{bu} \\ \phi_{bv}} \right) = \left[ \substack{cos \theta \:\:\: sin \theta \\ -sin \theta \:\:\: cos \theta} \right] \left( \substack{ \phi_{bL} \\ \phi_{bT}} \right) \tag{9.10} \end{equation}\]

Putting all of this together, and adding an isotropic diffusion term \(D_w\) (which can be used to represent wind/wave induced diffusion) yields:

\[\begin{equation} \left( \substack{\phi_{bu} \\ \phi_{by}} \right) = -hRDR^{-1} \left(\substack{\frac{\partial \phi}{\partial x} \\ \frac{\partial \phi}{\partial y}}\right) \tag{9.11} \end{equation}\]

\[\begin{equation} RDR^{-1} = \frac{|U|h\sqrt{g}}{C_z} \left[ \substack{ K_l cos^2 \theta + K_t sin^2 \theta \:\:\: (K_l - K_t) cos \theta sin \theta \\ (K_l - K_t) cos \theta sin \theta \:\:\: K_l sin^2 \theta + K_t cos^2 \theta} \right] + \left[ \substack {D_w \\ 0} \: \substack {0 \\ D_w} \right] \tag{9.12} \end{equation}\]

Note that this form is similar to that of Equations (9.5) and (9.6) but includes the off diagonal terms that arise in the matrix. Also note that the isotropic diffusion \(D_w\) in HPC is a global constant.

9.1.4 Local Constituent Transformation

In addition to pure advection and dispersion/diffusion, constituents simulated within TUFLOW AD (in both Classic and HPC) are modified by transient boundary conditions, and optional settling (for sediment constituents) and decay processes (dissolved and degradable constituents).

TUFLOW AD boundary conditions can be set to vary in time for each constituent, and can be applied to all TUFLOW boundaries that set water levels and flows (either user-specified or computed), such as HT or QT. TUFLOW AD also supports SA inflow boundaries, where flows and concentrations are used to compute mass loads that are delivered to the model domain, mixed with ambient water and then resultant concentrations computed, prior to execution of the advection routines. Consituents are not supported in rainfall boundaries, which are assumed to be pure water.

Settling of constituents to simulate removal of particulate matter from the water column has been included in the engine as a simple linear process. Once settled, constituents do not re-enter the computational domain. The local source contribution from settling is given by:

\[\begin{equation} S_\phi = -\frac{w_s}{h} \phi \: ; \: or \: S_b = -\frac{w_s}{h} b \tag{9.13} \end{equation}\]

Where:

  • \(w_s\) is the settling velocity for the particles

TUFLOW AD also supports the decay of individual species (if positive decay rates are specified) and employs first order rate equations to do so. These equations draw on user defined decay rates:

\[\begin{equation} S_\phi = -\omega \phi \: ; \: or \: S_b = -\omega b \tag{9.14} \end{equation}\]

Where:

  • \(\omega\) is the decay rate for the constituent

9.1.5 Groundwater

The 2023-03-AA release introduces tracking of constituents through the groundwater layers when run using the TUFLOW HPC solver, including return flow to surface water when horizontal hydraulic conductivity in included in the soils file. When this feature is used, groundwater return to the surface typically occurs at the bottom of slopes or adjacent to creeks, representative of catchment baseflows in reality. Refer to Section 9.3.6 for commands for setting initial concentrations for sub-surface layers. No dispersion is modelled in the groundwater layers. The 2d_bc “GT” type groundwater boundary does not support defined inflow values for constituents.

Advection-dispersion in groundwater is not supported in TUFLOW Classic or the TUFLOW HPC Quadtree Module.

9.1.6 Limitations and Recommendations

TUFLOW AD is designed to model dissolved and particulate constituent advection and dispersion in coastal waters, estuaries, rivers, and floodplains. This is achieved through solution of the 2D transport equation combined with sub-models for settling or decay.

Limitations to note include:

  1. Simulation of constituents through 1D culverts linked to 2D SX connections is currently only on a mass balance basis. That is, it is assumed that the concentration of a constituent exiting an SX connection is the same as that at the entrance SX connection at the same timestep. This approach conserves mass to the limit that the water inflows and outflows are approximately equal and that the transit time is small compared to the timescale at which constituent concentrations vary. In recognition of this simplification, only relatively ‘short’ 1D culverts linked to SX connections (using this timescale definition) should be simulated in the present release.
  2. The dispersion scheme adopted by TUFLOW AD (Falconer et al., 2005) is such as to allow use of literature values for longitudinal and transverse dispersion coefficients (\(K_{L}\) and \(K_{T}\)). Users adopting literature values for these coefficients (as provided in (Falconer et al., 2005)), should however do so with extreme caution, as they are known to vary widely, and by up to several orders of magnitude. It is always preferable to use monitoring data to calibrate advection dispersion models (TUFLOW AD included) and this should be done whenever and wherever possible. If no such data is available, then literature values can be used for \(K_{L}\) and \(K_{T}\), however results need to be appropriately caveated, and TUFLOW AD predictions (as for any AD model) should be seen as qualitative or indicative at best.
  3. Modelling predictions should also be cross checked with desktop calculations where possible. For example, this might include a hand calculation of expected salt masses in a given tidal system, with comparison made to TUFLOW AD outputs.
  4. TUFLOW AD allows for specification and computation of large dispersion coefficients, and with the automatic substepping implementation it should generally remain stable. However, specification of large (i.e. greater than approximately 100-500) dispersion coefficients may lead to results that are not physically real or defensible. As such, (in conjunction with 2 and 3 above) results should always be sanity checked and correlated with measurements. Relying on uncalibrated model predictions is not recommended.

9.2 The Modelling Process

9.2.1 Data Input Requirements

The minimum requirements for setting up a TUFLOW AD model are:

  1. A properly constructed and stable TUFLOW hydraulic model; and
  2. Boundary conditions for constituent concentrations (e.g. ocean salinities, catchment inflow pollutant concentrations etc.).

Initial conditions, dispersion coefficients, settling and decay rates will all be set to zero if not specified to be otherwise.

Preferable (and recommended) data requirements include:

  1. Water quality calibration information as time-series data at points. This is particularly important for dispersion coefficient calibration;
  2. Spatially variant initial conditions;
  3. Particulate matter settling rates (if any); and
  4. Dissolved species decay/transformation rates (if any).

9.2.2 Calibration and Sensitivity

Advection dispersion models are usually calibrated against water quality observations. For example, salinity recovery data can be used to calibrate and validate models, with longitudinal and transverse dispersion coefficients being the primary free variables. Dissolved and/or particulate constituents can then be simulated using the derived dispersion coefficients, and can include use of settling and/or decay rates as needed.

Ideally, models should be calibrated for conditions similar to those under investigation (e.g. a catchment inflow to an estuary) although this is not always possible, particularly when data is limited. In these situations, sensitivity analyses could be carried out by increasing and decreasing calibration variables, though this not a preferred approach due to the large variability in the literature with respect to acceptable dispersion coefficients.

9.2.3 Model Resolution

The 2D domain cell size needs to be sufficiently small to reproduce advection dispersion behaviour. It is worth noting that, in general, the larger the cell size is with respect to the scale of mixing processes, the greater potential there is for numerical dispersion to play a role in the model execution process. Even though TUFLOW AD has in-built measures to reduce these effects, it is advisable to keep this in mind, and perform at least some limited testing with your model to determine its sensitivity to cell size.

9.2.4 Computational Timestep

With TUFLOW Classic, the selection of the timestep is important in that the run time is directly proportional to the number of timesteps required to calculate model behaviour for the required time period. Notwithstanding this, with TUFLOW Classic the AD module will automatically substep with respect to the main hydraulic engine on the basis of maintaining both advective and dispersive stability (see Section 9.1.2) so the selection of timestep should be focused on ensuring hydraulic stability, as AD stability should follow, providing reasonable dispersion coefficients are set.

TUFLOW HPC, being an explicit method already, uses a self-adaptive timestep which is much smaller than Classic. Further, it uses the fourth order Runge-Kutta temporal integration scheme, and the advection dispersion calculations are performed on the same timestep and with the same scheme. When running AD models in HPC There is generally no need for the user to change the timestepping from the default method. Using high dispersion coefficients may cause TUFLOW HPC to use a smaller timestep than otherwise required for the hydraulic calculations due to the diffusion (Peclet) number control.

9.2.5 Example AD Models

A step-by-step descrition outlining how to set up an AD model is provided on the TUFLOW Wiki. An example model dataset is also available in the Advection Dispersion Example Models.

9.3 Data Input

9.3.1 Simulation Control File

9.3.1.1 TUFLOW AD Control File (.adcf File)

The TUFLOW AD Control File or .adcf file points to the two mandatory files required for AD model execution. It is the top of the tree for the AD model and is called into the .tcf by the AD Control File command. The .adcf file must reference:

  • One global database file, using the AD Global Database command. For example,
    AD GLOBAL DATABASE == ..\bc_dbase\2d_ad_globaldatabase_Demo1.csv
  • One boundary database file, using the AD BC Database command. For example,
    AD BC DATABASE == ..\bc_dbase\2d_ad_dbase_Demo1.csv

Appendix J lists and describes these commands and their parameters. No other commands are required in the adcf file.

9.3.1.2 GIS Input File Types and Naming Conventions

Spatial variation in initial condition and minimum dispersion coefficient is an optional data entry that can be used. If adding this feature to your model, it is recommended that the prefixes described in Table 9.1 be adhered to for 2D GIS layers. This greatly enhances the data management efficiency and importantly makes it much easier for another modeller or reviewer to quickly interpret the model. This approach is also consistent with that of TUFLOW.

Table 9.1: GIS Input Data Layers and Recommended Prefixes
2D Domain GIS Layers Suggested File Prefix Description Section
2D AD Initial Conditions 2d_ad_ic_ Layer containing polygons defining the spatial distribution of initial conditions for a given constituent. Optional. Section 9.3.5
2D AD Minimum Dispersion Coefficient 2d_ad_md_ Layer containing polygons defining the spatial distribution of minimum dispersion coefficients (\(D_w\)) for a given constituent. Optional. Section 9.3.7

9.3.2 1D Geometries

TUFLOW AD does not require any further 1D data in addition to the features already included in the TUFLOW hydraulic model.

9.3.3 Specification of Constituent Properties

Constituent properties are specified in the Global Database file, which is identified in the .adcf file using the AD Global Database command. This database file has a set structure, in much the same way as TUFLOW boundary database files, and can be created in software such as Microsoft Excel. The number of constituents simulated by TUFLOW AD is simply the number of line entries in this file (excepting the header data). Constituents can be removed from the simulation by prefacing rows in the global database file with the ‘!’ or ‘#’ character. The maximum number of constituents TUFLOW AD can simulate is 20.

The global database file must be .csv (comma delimited) formatted. The first row must contain the predefined keywords (in order) as listed in Table 9.2, separated by commas. Subsequent rows contain constituent data.

An example global database file is provided in the Advection Dispersion Example Models.

Table 9.2: Global Database Keyword Descriptions
Keyword Description
Name The name of a constituent. This might be ‘TN’ or ‘Salinity’ (without the inverted commas). The Name field is limited to 40 characters and must be alphanumeric characters only. Mandatory.
Heat Name Not currently used. Leave blank.
Decay Rate

The decay rate (k) of the constituent in units of day-1. This value is used in a first order decay calculations at each timestep, i.e.

\(C(t) = C_{0}e^{-kt}\)

Where:

\(\cdot\) \(C(t)\) is constituent concentration at time \(t\)
\(\cdot\) \(C_{0}\) is a reference concentration
\(\cdot\) \(k\) is the decay rate specified here; and
\(\cdot\) \(t\) is time

If no decay is required, then either enter 0 or leave the field blank.
Settling Rate

The settling rate of the constituent in units of m.day-1. This value is used in a simple mass balance calculation that removes the constituent from the water column based on this rate.

If no settling is required, then either enter 0 or leave the field blank.
Longitudinal Dispersion Coefficient

The value \(K_{L}\) for the constituent, as per Equation (9.5) for TUFLOW Classic and Equation (9.9) for TUFLOW HPC. Allowing such variation between constituents permits simultaneous simulation of multiple constituents with varying dispersion properties. This is useful at the model calibration stage when a range of dispersion coefficients can be tested within one simulation to ascertain the best match to monitoring data. This feature can also be used in sensitivity testing. Notwithstanding this, this value should not be varied from constituent to constituent once the AD model is calibrated.

If the field is blank or set to zero, then longitudinal dispersion is switched off.
Transverse Dispersion Coefficient

The value \(K_{T}\) for the constituent, as per Equation (9.6) for TUFLOW Classic and Equation (9.9) for TUFLOW HPC. Allowing such variation can be exploited in the same manner as described above for longitudinal dispersion. This value should not be varied from constituent to constituent once the AD model is calibrated.

If the field is blank or set to zero, then transverse dispersion is switched off.
Initial Condition

This can be set to either a single number or a path reference to a GIS layer (2d_ad_ic_). If the former is specified, then that value is applied uniformly to all wet cells at initialisation.

If a path to a GIS layer is specified, then the layer must comprise of (only) polygons with a single attribute. The attribute must be of type Float. The field name is not used by TUFLOW AD, so can be set to a label that is meaningful to the user. The polygon data, with spatially varying attributes of initial concentration is applied to wet cells at simulation initiation. Wet cells not covered by any polygon in the specified GIS layer are set to a concentration of zero.

If the field is left blank, then a concentration of zero is applied to all cells.
Minimum Not currently used. Leave blank.
Maximum Not currently used. Leave blank.
Minimum Dispersion

This sets the value of \(D_{w}\) as per Equations (9.5) and (9.6) for TUFLOW Classic and Equation (9.12) for TUFLOW HPC. It can be set to either a single number or a path reference to a GIS layer. If the former is specified, then that value is applied uniformly to all wet cells at each timestep.

If a path to a GIS layer is specified, then the layer must comprise of (only) polygons with a single attribute. The attribute must be of type Float. The field name is not used by TUFLOW AD, so can be set to a label that is meaningful to the user. The polygon data, with spatially varying values of \(D_{w}\) is applied to wet cells at each timestep. Wet cells not covered by any polygon in the specified GIS layer will have \(D_{w}\) set to zero.

This field cannot be left blank. If a minimum dispersion of zero is required, then 0.0 must be entered in this field. Errors will result if this field is left blank.

9.3.4 Boundary Conditions

TUFLOW AD uses the same approach as TUFLOW to setting up boundary conditions, in that two types of files are required:

  • Boundary condition database; and
  • Boundary condition data files (i.e. timeseries).

Like TUFLOW, TUFLOW AD also uses a comma delimited format for both these file types.

TUFLOW AD does not require specification of any geographical information regarding the location of boundary conditions. All such required data is passed from TUFLOW to TUFLOW AD for TUFLOW boundary types HT, QT, HS, HQ, QC, VC, VT, and SA.

9.3.4.1 Boundary Condition (BC) Database

A boundary condition (BC) database is created using spreadsheet software, such as Microsoft Excel. It must be .csv (comma delimited) formatted and is identified in the .adcf file (see AD BC Database). The database contains a list of files and attribute names to search for within those files. The attribute names are then used to extract the desired boundary condition data.

The TUFLOW AD BC Database is structured in the same way as a standard TUFLOW BC Database, in that it must contain a header line with subsequent rows of information. The header line must contain the following keywords, in the below order, with meanings as per Table 9.3.

“Name, Source, Column 1, Column 2, Add Col 1, Mult Col 2, Add Col 2, Column 3, Column 4”

Table 9.3: AD BC Database Keyword Descriptions
Keyword Description
Name

The name of a BC data set. It consists of two concatenated elements as follows:

  • First: This must be the same ‘Name’ attribute used in the GIS 2d_bc layer(s) specified in TUFLOW. This is third attribute in the 2d_bc layer and follows ‘Type’ and ‘Flags’. It may contain spaces, but must not contain commas. It is not case sensitive.
  • Second: This must be the name of the appropriate constituent for which the boundary condition is being specified. This name must be one of those specified in the AD Global Database ‘Name’ field.

These two names are joined by a double underscore. As such, if there are “M” 2d_bc boundary condition objects specified in TUFLOW, and “N” constituents to be simulated by TUFLOW AD, then M*N entries are required in the AD BC Database.

An example is Ocean__Salinity, where the TUFLOW 2d_bc file has a line with Name attribute “Ocean” and the AD Global Database has a variable called “Salinity” listed.

Mandatory.
Source

The file from which to extract the BC data. It must be a .csv file. Paths are relative to the AD Global Database.

Mandatory.
Column 1

The name of the first column of data (time values) in the .csv Source File.

Mandatory.
Column 2

The name of the column of data in the .csv Source File. Values in this column are always concentrations, for all allowable BC sets. These are applied to cells corresponding to the appropriate locations of all types of BCs specified in TUFLOW 2d_bc GIS layers. These specified concentrations override any other computed concentrations at the boundary locations.

The exception to the above is for SA boundaries, where the specified concentration is multiplied by the inflow corresponding to that concentration (as passed to TUFLOW AD from TUFLOW) and the mass load over the SA polygon computed. This mass load is then added to the mass of constituent within the wet computational cells within the SA polygon, and mass conservation used to compute a new resultant ambient concentration.

Mandatory.
Add Col 1 Not used. Leave Blank.
Mult Col 2 Not used. Leave Blank.
Add Col 2 Not used. Leave Blank.
Column 3 Not used. Leave Blank.
Column 4 Not used. Leave Blank.

9.3.4.2 BC Database Example

The Excel spreadsheet below illustrates a simple example of a TUFLOW AD BC Database, created in a Microsoft Excel worksheet that is exported as a .csv file for use by TUFLOW AD. Four line boundaries have been specified via GIS in the TUFLOW model (with Name fields “West”, “North”, “East” and “South”) and two constituents have been specified in the AD Global Database (with Name fields “Tracer_01” and “Tracer_02”). Four boundaries and two constituents thus require eight entries in the TUFLOW AD BC Database file.

The Demo1_Concs.csv file was created by saving the below as a .csv file from Excel. All values are concentrations in mg/L, with the exception of the Time column, which has units of hours.

9.3.5 Initial Conditions

Initial conditions are specified for each constituent as either a constant value or spatially varied value (see Table 9.2). The former is simply entered as a decimal (or integer) in the Initial Condition field of the AD Global Database. This value is applied to all wet cells at model initiation. The latter is applied by entering a relative (or absolute) file path to a GIS layer in the Initial Condition field of the AD Global Database. The GIS file has the attributes described in Table 9.4.

Table 9.4: 2D Initial Conditions (2d_ad_ic) Attribute Descriptions
GIS Attribute Description Type
Conc The initial condition concentration. Float

As many polygons as needed can be included in this layer. Any wet cells not covered by these polygons will be initialised to a concentration of zero. The naming convention prefix for this layer is 2d_ad_ic_. The objects must be polygons. Rectangles, round rectangles etc. are not supported by TUFLOW. The attribute name is not read by TUFLOW AD and can be anything meaningful to the user - “Conc” is used as an example above for clarity.

9.3.6 Groundwater Initial Conditions

Initial concentrations can be set within groundwater layers using the following .tgc commands:

[Set | Read GIS | Read Grid] IGW Conc Layer N1,N2,..N Tracer M1,M2,..M == [value | /path/to/file]

Where:

  • N = The groundwater layer number (reference multiple in a single command)
  • M = The tracer number (reference multiple in a single command), as ordered in the AD Global Database

If a layer number is not referenced, it is assumed to apply to layer 1. Likewise, if a tracer number is not referenced it is assumed to apply to the first tracer. Note, the order that ‘Layer’ and ‘Tracer’ appear in the command does not matter (i.e. tracer numbers could be listed before layer numbers).

9.3.7 Minimum Dispersion Coefficient

Minimum dispersion coefficients are specified for each constituent as either a constant value or spatially varied value (see Table 9.2). The former is simply entered as a decimal (or integer) in the Minimum Dispersion field of the AD Global Database. This value is applied to all wet cells at all timesteps. The latter is applied by entering a relative (or absolute) file path to a GIS layer in the Minimum Dispersion field of the AD Global Database. The GIS file has the attributes described in Table 9.5.

Table 9.5: 2D Minimum Dispersion Coefficient (2d_ad_md) Attribute Descriptions
GIS Attribute Description Type
Conc The initial condition concentration. Float

As many polygons as needed can be included in this layer. Any wet cells not covered by these polygons will be assigned a minimum dispersion coefficient of zero. The naming convention prefix for this layer is 2d_ad_md_. The objects must be polygons. Rectangles, round rectangles etc. are not supported by TUFLOW. The attribute name is not read by TUFLOW AD and can be anything meaningful to the user - “MD” is used as an example above for clarity.

9.4 Data Output

9.4.1 TUFLOW AD Log Files

TUFLOW AD log files are written to the location specified by Log Folder command in the .tcf file, or if this command is not usedthe same location as the .tcf file. The following logs are available:

  • An AD simulation log file (.adlf);
  • A CFL condition log file - output when using Write CFL == ON, for TUFLOW Classic only; and
  • A mass balance log file - output when using Write Mass == ON.

9.4.1.1 Simulation Log File

This file is named ‘TCF name.adlf’. It contains commentary of input reading, constituent specification etc. as the simulation sets itself up.

Following model initialisation, if Write Iterations == ON is used, for TUFLOW Classic only, the following is output at each timestep. For example:

Simulation time (h) 0.0125
Finished constituent Tracer_01 at AD substep iteration number 1
Finished constituent Tracer_02 at AD substep iteration number 1
Finished constituent Tracer_02 at AD substep iteration number 2
Finished constituent Tracer_02 at AD substep iteration number 3

These lines show information regarding the simulation time for each constituent. In particular, the number of AD sub-steps needed to be executed to maintain stability under CFL and Peclet conditions is reported. In the above example, constituent Tracer_01 required no substepping, whilst constituent Tracer_02 required 3 substep iterations to maintain stability. This is caused by Tracer_02 being set up with greater dispersion coefficients than Tracer_01 in the AD Global Database. Additional rows are added as required by the number of constituents simulated.

This TUFLOW AD timestep information is not available for TUFLOW HPC, or needs review by the modeller, due to the fact the TUFLOW HPC hydraulic engine and advection dispersion calculations are performed on the same timestep and with the same scheme. See Section 9.2.4 for discussion regarding timestep differences relating to TUFLOW Classic and TUFLOW HPC.

9.4.1.2 CFL Log File

This file is named ‘TCF name_ADcfl.csv’. It can be output when using TUFLOW Classic in conjunction with the Write CFL == ON command. It contains commentary on CFL and Peclet numbers, maximum dispersion coefficients and number of substep iterations. The column data is as per Table 9.6:

Table 9.6: _ADcfl.csv File Columns
Column Description
time (h) The simulation time in hours.
Constituent Name The name of the constituent as specified in the AD global database.
Max_CFL_u The maximum CFL for u velocities anywhere in the computational domain at that timestep.
Max_CFL_v The maximum CFL for v velocities anywhere in the computational domain at that timestep.
Max_Peclet_u The maximum Peclet number for u dispersion anywhere in the computational domain at that timestep.
Max_Peclet_v The maximum Peclet number for v dispersion anywhere in the computational domain at that timestep.
Max_sum_u The maximum of the sum of CFL and Peclet numbers in the x (u) direction anywhere in the computational domain.
Max_sum_v The maximum of the sum of CFL and Peclet numbers in the y (v) direction anywhere in the computational domain.
Max_Disp_x The maximum dispersion coefficient in the x (u) direction anywhere in the computational domain.
Max_Disp_y The maximum dispersion coefficient in the y (v) direction anywhere in the computational domain.
Num_iterations The number of iterations required by TUFLOW AD to remain stable. This can vary from constituent to constituent if different dispersion coefficients are applied.

As per the TUFLOW AD timestep information in the ‘.adlf’ file. This log information is not available for TUFLOW HPC, or needs review by the modeller, due to the fact the TUFLOW HPC hydraulic engine and advection dispersion calculations are performed on the same timestep and with the same scheme. See Section 9.2.4 for discussion regarding timestep differences relating to TUFLOW Classic and TUFLOW HPC.

9.4.1.3 Mass Log File

This file is named ‘TCF name_ADmass.csv’. It can be output when using TUFLOW Classic or TUFLOW HPC in conjunction with the Write Mass == ON command. It contains commentary on mass conservation. Specifically, the column data is as per Table 9.7:

Table 9.7: _ADmass.csv File Columns
Column Description
time (h) The simulation time in hours.
Constituent Name 1 The total mass of that constituent in the computational domain. TUFLOW AD assumes that the constituent concentrations are specified in mg/L, and this number is then in tonnes of constituent. If the concentration is g/L, then this number should be multiplied by 1000 to be in units of tonnes.
Constituent Name 2 As per constituent 1. Column repeated until all constituents have been accounted for

9.4.2 Check Files

When using the AD module, an _AD_check file will be output. For check file information, see the _AD_check TUFLOW Wiki Page.

9.4.3 Result Files

Result files contain the computed spatial and temporal evolution of simulated constituents. The output is a temporal format (maximums are not available) and will be output in the format specified by the Map Output Format command. For example, if Map Output Format == XMDF is used, the XMDF output will contain a result type called “Conc Tracer_X”.

The AD output will automatically be produced unless the Map Output Data Types command is set to “AD OFF”. This off option can also be set on a Map Output format basis (e.g. FLT Map Output Data Types == AD OFF).

Plot Outputs are currently not supported for TUFLOW AD outputs.

References

Falconer, R. A., Lin, B., & Kashefipour, S. M. (2005). Modelling Water Quality Processes in Estuaries. In P. B. Bates, S. N. Lane, & R. I. Furguson (Eds.), Computational fluid dynamics: Applications in environmental hydraulics. John Wiley; Sons.
Leonard, B. P. (1991). The ULTIMATE Conservative Difference Scheme Applied to Unsteady One-Dimensional Advection. Computer Methods in Applied Mechanics and Engineering, 88(1), 17–74.
Leonard, B. P., MacVean, M. K., & Lock, A. P. (1993). Positivity-Preserving Numerical Schemes for Multidimensional Advection (No. 93N27091). NASA Centre for Aerospace Information.
Leonard, B. P., & Niknafs, H. S. (1991). Cost-Effective Accurate Coarse-Grid Method for Highly Convective Multidimensional Unsteady Flows (No. 91N21075). NASA Centre for Aerospace Information.
Wu, & Falconer, R. A. (2000). A Mass Conservative 3-d Numerical Model for Predicting Solute Fluxes in Estuarine Waters. Advances in Water Resources, 23(5), 531–543.