Operating Guidelines for MIDAS
MIDAS (Model Investigating Density And Size-sorting) was constructed using Fortran77 (and some Fortran90) with a main program and about 27 subroutines. It is a combined physical-empicial numerical model because most of the equations are derived from first principles, but some of the equations are based on curves derived from empirical datasets. An uncoupled, steady solution of the fluid and sediment transport equations is used; the gradually varied flow equation is solved using standard backwater calculations (see Henderson, 1966). The bedload transport of each size and density of material is calculated from a modified-Bagnold equation (Bagnold, 1973; Bridge and Dominic, 1984). Critical shear stresses for entrainment of individual sizes and densities are determined from either Komar, Egiazaroff, or James relationships (Komar, 1987; Egiazaroff, 1967; James, 1990). 
Reading input values and following subroutines
There are two input files that define the inital and boundary conditions and constants, and a third file that names the input and output data files; the main program first looks at the file files for the two input file names (here called input.1 and input.2) and the root for the output graphics files (here called freeware). The main program then reads 30 values including fluid properties, timesteps, length of reach, node spacing, starting elevation, slopes, nature of initial profile, width of reach and how it changes downstream, and various options for switching on particular equations and conditions (see Options).
After the initial width values are read in, the main program calls a subroutine called ELWIDCALC where the initial profile, widths and initial and last boundary node conditions are set.
After returning to the main program the remaining values from input.1 are read in, including data file writing instructions and the grain-size distribution. If a computer-generated distribution is required, the subroutine GSDIST is called; a user-defined grain-size distribution has input intervals, average weight percentages and smallest grain size available. 
After setting up the tributary grain-size input (TRIB), and ordering the trunk stream grain-size data by size (and tributary, if turned on) and establishing the median size (D50) (DISCR and DFIFTY), the main program will check that the input file is at the beginning of the run or is a restart file (if latter, call READRSTRT). Settling velocities and entrainment values are then calculated (SETTLE and ENTRAIN), after which the gradually varied flow equation is solved (FLDTA). In FLDTA, a hydrograph can be used in the flow routine if OPTDISC=1. Shear velocities and stresses are then calculated (SVELA). Flow data is written to an output file with the suffix a appended to the output filename (e.g., RUNNAMEa). Elevation and D50 values for each node are written to output file named plot.freeware, grain-size distributions for each density at each node are written to min1.freeware (depending on how many mineral types there are), and flow depths are written to dch.freeware. 
The next steps in the main program involve bedload and suspended load transport, and turbulence-adjusted shear stress (TURB). If grain hiding is included as an option, those subroutines including equations that cover hiding are selected (ENTRAINH, DRANGEH). The grain sizes in transport are calculated (BEDLOAD) and some of the fluid and sediment properties and calculation are printed to RUNNAMEb. The main program then calls LOGDIST to calculate the logarithmic velocity distribution and vertical diffusivity profile (EXVERT). The bed profile is the updated if there is a tributary effect (TRIB), bedload and suspended load are calculated and the bed profile is again updated (SUSP1, SUSP, and ERODEP), and then finally bed profile adjustments are made if the reach is undergoing subsidence (TECTONICS). 
Finally, the new D50 for the bed is calculated (DFIFTY), data output is tabulated and written (RESULT) to output files named RUNNAMEc and -d. If a restart file is to be saved, subroutine WRITERSTRT prints to a file named resRUNNAMEe with the timestep at saving time as a suffix. 
Options
The main input file (input.1) has several options (called OPT- ) that specify whether particular functions are turned on (on=1) or which function is to be used (=1-5). Two examples are OPTURB, which is a flag that signifies whether turbulence is treated as a normal (=1) or lognormal (<>1) shear stress distribution, and OPTERO, which determines whether erosion or deposition is allowed at the head of the reach (=1) or not. All the options are defined in the variable listing below.
Variables used in MIDAS
The variables are defined in two ways; in one large common block and/or as local variables. The local variables are those variables that only occur in a subroutine; if local variables occur in more than one subroutine they are not redefined in any following subroutines (unless their meaning has changed). Finally, all variables used in MIDAS are listed alphabetically in the Variable Index section, and the reader is refered either to the alphabetically listed common block or to a subroutine.
The timestep counter file (tst) is updated every printstep (not every timestep) to increase computational efficiency. The counter can be placed at the beginning of the main loop in AMAIN if needed.
Common block variables (in common.h)
A: a constant used in the Bridge bedload transport equation.
ACTIVE: the thickness of the active layer (m).
AVDENS (JN): the average, or bulk, density of the parent material at node JN (kg/m3).
AWPM (I,J,JN): the current weight percentage of grain size I, mineral species J at node JN.
AWPOM (I,J,JN): the original weight percentages of the parent bed material.
AWPT (I,J): the weight percentage of grain size I, mineral J of the tributary material (if OPTTRIB=1) feed in at each node downstream.
B: a constant used in the suspension criterion.
BCOUNT: a counter, which keeps track of how many timesteps have elapsed since the previous renewal of the bedload feed rate.
BEDCON (I,J,JN): the concentration of grain diameter I of mineral J at node JN at a reference height (y=a) above the bed, as expressed in "dry" mass per unit volume(kg/m3).
BETA: a weighting constant used in the convection-diffusion suspended load algorithm.
BLWTAM (I,J,JN): the “dry” mass of the bedload per unit surface area (kg/m2).
BLWTM (I,J,JN):  the “dry” bedload transport rate for grain size I, mineral species J and node IN, expressed in transported mass per unit width per unit time (kg/m/s).
BLWTO (I,J): the original (at NT=0) bedload feed for grain size I, mineral species J at the head of the reach (kg/m/s).
BLWTT (I,J): the bedload tributary feed rate for grain size I, mineral J along the reach (kg/m2).
BLWTTO (J): the original bedload tributary feed along the reach for mineral J; it is converted to BLWTT (I,J) in program with AWPTT (I,J).
BREAK: a longitudinal distance downstream where the slope changes.
CON1: the constant used in the Komar-like equation.
CON2: the constant used in the Egiazaroff-like equation.
CONC (L) the suspended concentration at level L in the flow (kg/m3).
CONIM (I,J,JN):  the total, depth integrated suspended concentration of grain size I, mineral J at rate of grain size I, mineral J at node JN (kg/m/s).
CONCM (I,J,IV): the suspended concentration of grain size I, mineral J at the previous node  at level IV (kg/m3).
COUNTER: a counter equal to the number of timesteps since the last resetting of the summary arrays.
COVAR: the coefficient of variation () of the shear stress distribution.
CVELM (I,J,JN): the suspended load at node JN (kg/m/s).
D50M (JN): the mean bed grain diameter, by volume (m).
DBDX (I,J,JN): the mass flux of grain size I, mineral J at node JN between the active layer and the moving bed layer, expressed in mass per unit length (kg/m).
DCHM (JN): the flow depth at node JN (m).
DEFEND (I,J): the smallest grain diameter larger than grain size interval I of mineral J (m).
DELQ: a coefficient used to obtain downstream discharges based on cross-stream width.
DELT: the time step, input as hours (converted to seconds in program), between calculations.
DELWID: a constant used to numerically generate cross-stream widths at downstream nodes.
DELY (JN): the vertical distance between flow levels at which suspended load calculations occur at node JN (m).
DELX: the distance between nodes, usually chosen as an integral fraction of the reach length (m).
DENDM: the upper bounds of all grain size intervals, of all the mineral densities, arranged in ascending, non-repeating order (m).
DIMIN (J): the smallest grain diameter of mineral J present in the bed (m).
D50: the midpoint of the grain size interval I of the mineral density J (m).
DINIT: the initial (t=O) depth at the control section (m).
DMAX (K,J,JN): the largest grain diameter of mineral J moving as bedload under shear stress interval K (m).
DMIN (K,J,JN): the smallest grain diameter of mineral J moving as bedload under shear stress interval K (m).
EFOLDX: the e-folding length for the exponential subsidence rate.
ELEVOM (JN): the elevation of the top of the parent material at node JN (m).
ELEVM (JN): the bed elevation relative to a user defined base at node JN (m).
ELEVTOP: the bed elevation of the first point along the reach (m).
EPSY (L,JN): the vertical fluid or sediment eddy diffusion coefficient at level L in the flow at node JN (m2/s).
FLOWEL: the elevation of the free surface of the flow at the control section (m).
FOUT: string containing the root name of the output files.
G: gravitational acceleration (m/s2).
GRAD: a new slope for application to the downstream end of the reach to make a ponding effect. 
ICONT: the node number of the control section.
IMAX (K,J,JN): grain size index of the largest grain size of mineral J moving as bedload under shear stress interval K.
IMIN (K,J,JN): grain size index of the smallest particle size of mineral J moving as bedload under shear stress interval K
KR: length of the string containing the root name of the output files.
KK1: a comminution coefficient for the abrasion algorithm (1/m).
M: the number of nodes along the reach.
MANN (JN): Manning's roughness coefficient for node JN (s/m 1/3).
MEAN (J): the desired mean grain diameter of the mineral species J (m).
MSPRTHK (I,J,JN): mass per unit area of grain size I, mineral J eroded (+) or deposited (-) at node JN (kg/m3).
N: the number of vertical flow levels dividing the stream depth.
NEXT: the next timestep that a new discharge regime will be read in.
NK: the number of shear stress intervals into which the shear stress distribution is divided.
NSIGMA: the number of distinct mineral species present in the river system.
NSIZE (J): the number of grain size intervals for mineral species J.
NT: the current time step.
NTIM: the number of time steps, of length DELT, to be performed.
NTOT: the number of unique upper bounds of the grain size intervals present in all the mineral densities.
NWAVDENS (JN): the new, updated, bulk density of the bed (kg/m3).
OPTABRAD: a flag determining whether (=1) or not (=0) abrasion will be incorporated.
OPTBINP: a flag determining the bedload sediment feed condition at the head of the reach.
OPTCONT: a flag determining how water depth at the control section is calculated.
OPTDISC: a flag determining whether discharge is based on a hydrograph file (=1)
OPTDIST (J): a flag determining whether the grain size distribution for mineral species J is to be read in (user defined; OPTDIST<>l) or if a numerically-generated phinormal distribution is desired (OPTDIST=I).
OPTEL: a flag determining whether the node bed elevations are to be user specified (OPTEL=1) or numerically generated (OPTEL=2 to 5).
OPTERO: a flag which determines whether erosion and deposition at the head of the stream is allowed (OPTERO=1) or not.
OPTFUNC: a flag that, if grain hiding has been specified, signifies a Komar-like (OPTFUNC=1) or the Egiazaroff-like (OPTFUNC=2) entrainment function.
OPTHIDE: a flag that signifies whether effects of grain hiding are to be incorporated (OPTHIDE=2) or not (OPTHIDE=1).
OPTLNODE: a flag, between 0 and 5, determining the boundary condition at the last node.
OPTSALT: a flag which determines whether the Bridge equation (OPTSALT=1) or the Einstein equation (OPTSALT=2) is used to determine grain saltation height.
OPTSINP: a flag determining the suspended load sediment feed condition at the head of the reach.
OPTSUSP: a flag determining if the suspended load is calculated using a convection-diffusion algorithm (OPTSUSP=l) or using the Rouse equation (OPTSUSP=2).
OPTTECT: a flag determing whether linear subsidence (=1) or exponential subsidence (=2) are turned on. Deflection maximum at the head of the reach.
OPTTRIB: a flag determing whether sediment seed will be added downstream.
OPTURB: a flag which signifies whether a normal (OPTURB=1) or a lognormal (OPTURB<>1) shear stress distribution is desired.
OPTWID: a flag determining whether the node widths are to be user specified (OPTEL=1) or numerically generated (OPTEL=2 to 5).
PHI: a weighting constant used in the modified Preissman difference scheme solving the sediment mass balance equation.
PRINT (J); J=1,7: a series of flags used to allow the user to choose which output files are desired. Some output files can grow to enormous proportions (especially the .b   and .d output files) and do not always contain information relevant for a certain modelling purpose. Hence:
PRINT (l): turns the a output file on (= 1) or off (<> 1).
PRINT (2): turns the b output file on  (=1) or off (<> 1).
PRINT (3): turns the c output file on  (=1) or off (<> 1).
PRINT (4): turns the d output file on  (=1) or off (<> 1).
PRINT (5): turns the e output file on  (=1) or off (<> 1).
PRINT (6): turns the graphics plotfile file on (= 1) or off (<> 1).
PRINT (7):  a constant used to limit the number of times the .b and .d  files are printed. The .b and .d files will be written to every PRINT(7)*TPRINT time.
PLOTPRINT: the longitudinal distance between the prints to the plotfiles.
PROPT (K): the proportion of the time step over which shear stress interval K acts on the bed.
PROPTJ (J): the weight proportion of the mineral species J in the bed material.
PRTHICK(I,J,JN): the thickness of the amount of material of grain size I, mineral J eroded or deposited at node JN (m).
PSI: a weighting constant used in the convection-diffusion suspended load algorithm.
QCOUNT: a counter equal to the number of times steps since the new water discharge was read in.
QDSMOD: a flag signifying whether a geomorphic stream regime relationship is used (=1) or a linearly increasing function of distance (=2) or not changing downstream (=0).
QDTMOD: a flag which signifies whether (=1) or not (=0) temporally varied water discharge is desired.
QFAC: a constant used to obtain downstream discharges based on cross-stream width.
QN (JN): the water discharge at node JN (m3/s)
QO: the inital (t=0) water discharge entering at the head of the reach (m3/s).
QSTMOD: a flag indicating whether temporally varied bedload sediment discharge is desired (=1) or not (<>1).
REGF: a constant used to calculate bed elevations.
RENEW: a temporal constant, signifying the time, in hours, after which the sediment discharge is restored to its original value.
RHO: the density of the sediment-transporting fluid (kg/m3).
RLONG: the total length of the stream or reach (m).
RUN: a string describing the purpose of the run of up to 72 characters long.
RUNNAME: a 6 character string used as the root name for the main output files.
SBLWT (I,J,JN): the instantaneous bedload transport rate of grain size I, density J at node JN (kg/m/s).
SBLWTA (I,J,JN): the time-averaged bedload transport rate of grain size I, density J at node JN (kg/m/s).
SCVELI (I,J,JN): the instantaneous suspended load transport rate of grain size I, density J at node JN (kg/m/s).
SCVELIA (I,J,JN):  the time-averaged suspended load transport rate of grain size I, density J at node JN (kg/m/s).
SDBDT (I,J,JN): instantaneous mass flux rate between the moving bed layer and the active layer (kg/m2/s).
SDBDX (I,J,JN):  the time-averaged mass flux between the moving bed layer and the active layer (kg/m2/s).
SETV (I,J): the settling velocity of minimum size I of mineral J (m/s).
SIGMA (J): the density of mineral species J (kg/m3).
SHVEL (K,JN): the midpoint of shear velocity interval K at node JN (m/s).
SHVELC (I,J,JN): the critical shear velocity necessary to entrain grain size I of mineral J at node JN (m/s).
SLOPE: a constant used to numerically generate bed elevations at each node.
SLWTAM (I,J,JN): the saltation transport rate of grain size I of mineral J at node JN per unit area, as expressed in "dry" mass transported per unit bed area(kg/m2).
SLWTM (I,J,JN): the saltation transport rate of grain size I of mineral J, at node JN, as  expressed in "dry" mass transported per unit width per unit time(kg/m/s).
SPRT'HICK (I,J,JN): the thickess of the material eroded or deposited over a timestep, of grain size I, density J at node JN (m).
STDEV (J): the standard deviation of the phinormal grain size distribution desired for mineral J (m).
STFLUX (I,J,JN): the instantaneous mass flux rate between the suspension load and the active layer (kg/m2/s).
STRDATA (1,JN): the hydraulic radius of the stream at node JN (m).
STRDATA (2,JN): the Froude number at node JN.
STRDATA (3,JN): the slope of the energy grade line at node JN.
STRDATA (4,JN): the downstream distance of node JN (m).
STRDATA (5,JN): mean stream velocity at node JN (m/s).
STRDATA (6,JN): fluid discharge at node JN (m3/s).
SVELM (JN):  the temporally averaged, skin-friction component of the bed shear velocity node JN (m/s).
SVELUM (JN): the total, temporally-averaged shear velocity at node JN (m/s).
SXFLUX (I,J,JN): the time-averaged mass flux between the suspension load and the active layer (kg/m2/s).
TC (I,J,K): the critical shear necessary to entrain minimum size I of mineral J at node JN (Pa).
TDBDX (JN): the total bedload flux at node JN, expressed in mass per unit length (kg/m).
TECRATE: tectonic subsidence rate in meters/timestep.
TFILE: the time, in hours (input) or seconds (program), between prints to the restart file.
THALFB: a temporal constant, indicating the time, in hours, over which bedload sediment  feed at the head exponentially decreases by half.
THALFQ: the time, in hours, over which the water discharge at the head exponentially decreases by half.
THETAC: the critical Shields parameter necessary to entrain a grain from a homogeneous bed under fully rough conditions (dimensionless).
THET50 (JN): the critical Shields parameter necessary to entrain the mean bed grain diameter of node JN.
T'HICKM (JN): the total thickness eroded or deposited at node JN (m).
TMSPRTHK (I,J,JN): mass per unit area from tributary material of grain size I, mineral J deposited (-) at node JN (kg/m2).
TO (K,JN): the midpoint of shear stress interval K at node JN (Pa).
TOTTHICK (I,J,JN): the total thickness eroded or deposited over the length of the run, of grain size I, density J at node JN (m).
TOXM (JN): the temporally-averaged skin-friction component of the bed shear stress at node JN (Pa).
TPRINT: the time, in timesteps between prints to files, with the exception of the restart file and the -b and -d files which print every PRINT(7)*TPRINT.
TPRTHICKM (I,J,JN): thickness of material deposited of grain size I, mineral J at node JN (m).
TXFLUX (JN): the total suspension load flux at node JN, expressed in mass per unit length (kg/m).
VEL (L,JN): the longitudinal velocity at flow level L at node JN (m/s).
VISC: the viscosity of water (kg s/m2).
VISKIN: the kinematic viscosity (viscosity/density) of the sediment-transporting fluid (m2/s).
VOLCON: the volumetric concentration of the grains in the bed (porosity) equal to the volume of the grains in the bed over the total bed volume.
WEITFAC (JN): a weighting factor to balance the incoming tributary AWPT and existing trunk stream AWPM to the new AWPM.
WIDM (JN): the cross stream width at node JN (m).
WIDTOP: the cross-stream width of the first point along the reach (m).
XFLUX (I,J,JN): the mass flux of grain size I, mineral J at node JN between the active layer and the suspension layer expressed in mass per unit length (kg/m).
XPRINT: the longitudinal (x) distance over which results are written to the main files (m).
YMAX (JN): the height of the moving bed layer at node JN (m).
ZEROLOC: the downstream distance at which there is zero subsidence effect.
Loop Variables
I: the grain size index counter.
J: the mineral species counter.
IN: the node number counter.
K: the shear stress interval counter (!except in DISCR, DFIFTY and FLDTA).
L: the flow level counter.
NT: the time step counter.
Local Variables in Subroutines
AMAIN90
AWP (I,J): the weight percentage of grain size I of mineral species J if the species J were sieved separately; for every mineral species J all the AWP (I,J) must sum to 100.
DATETIME: date-time stamped character data of discharge values that are translated in FLDTA.
FI (NRUN): the name of the main input file.
F2 (NRUN): the file that contains water discharge and depth information.
F3 (NRUN): the third filename, and the root for the graphics output files.
FILENAME: the name of the file out of which the elevation or the width data are read in (only if OPTEL or OPTWID=1).
IRSTRT (NRUN): a flag set equal to 1 when NRUN is a restarted run; else it equals 0.
IRUN: index number of the run presently being processed.
MANN: Manning's roughness coefficient for the entire reach (s/m1/3)
MIN: character string set equal to the file names of the graphics files containing the weight percentages of the different mineral species at various time steps.
NRUN: the total number of runs to be done.
NTFIN: the last time step.
NTSTAR: the time step at which we start. If a restart file has not been specified NTSTAR=0, else it equals NT+1.
BEDLOAD
CONV: a conversion factor =LOG(0.5), used in converting the ‘half life’ value of the bedload feed, THALFB.
IF: grain size index number of the largest particle size of mineral species J that can be entrained as bedload by the Kth shear stress interval.
IS: grain size index number of the smallest particle size of mineral species J that can be moved as bedload by the Kth shear stress interval.
MINISIF: grain size index number of the largest suspendable grain size of mineral J, at node JN, under the shear stress interval K.
PROPVOL: the volumetric proportion of the grain I, J currently under consideration, relative to the total active layer volume.
TMP: a variable used in the bedload transport rate calculations.
VOL: a variable proportional to the total volume of the grains in the active layer by a factor equal to 1 over the total active layer mass (m3/kg).
DFIFTY
K: the discretized grain size index number.
LOWDEN: a variable equal to the last calculated value of the discretized size array.
LOWLIM: a variable equal to last considered upper bound of a grain size interval.
SUM: a variable equal to the sum of the proportional volumes (VOLTOT).
VOL: the total proportional volume of the active layer, equal to the total volume of active  layer over the total mass of the active layer (m3/kg).
VOLTOT (K):  the  proportional volume of the grain size interval between DENDM (K-1) and DENDM (K), equal to the true volume in the interval over the total mass of the active layer (m3/kg).
DISCR
K1, K2: sorting loop variables.
NRP: the running counter of the number of repeat discretized size fractions.
TEMP: a variable used to exchange a higher with a lower value during the sorting procedure.
DRANGE
AA, BB, CC: temporary variables equal to the quadratic equation constants a, b, c respectively:  ax2 + bx + c = 0.
COUNT: a variable keeping track of the total number of iterations in the Newton-Raphson solution scheme.
DFD: a variable equal to the evaluation of the derivative of the Dietrich equation with respect to the dimensionless grain diameter at a value D*=DIMDL.
DI: a first estimate of the grain diameter that, according to the Dietrich equation, has a settling velocity SETVEL.
DIMDL: the dimensionless grain diameter, to be used in the Dietrich equation.
DIMSV: the dimensionless settling velocity, to be used in the Dietrich equation.
DIN: a variable equal to the adjusted Newton-Raphson estimate of the dimensionless grain diameter that is equivalent to the declared settling velocity.
FD: a variable equal to the Dietrich function evaluation of the calculated DIMDL.
REG: the grain Reynolds number.
SETVEL: a variable equal to the settling velocity of the largest grain in suspension under the Kth shear stress interval.
DRANGEH
D50: the mean grain diameter, by volume, in the bed (m).
DI: the bisection algorithm estimate of the maximum grain size entrainable by TO (K,JN) following the Egiazaroff equation.
DX: a temporary variable used in the bisection algorithm.
F: the Egiazaroff function evaluated at the lower bound, X1.
FMID: the Egiazaroff function evaluated at the upper bound, X2.
FUN1: the Egiazaroff function.
RTBIS: a variable used as a root guess in the bisection algorithm.
SUR: a variable equal to the denominator of the Egiazaroff equation.
X1: the grain diameter at which the Egiazaroff function reaches a minimum; used as a starting guess for the bisection algorithm.
X2: a very large grain diameter, used as the second, upper bisection algorithm guess.
XACC: the accuracy of the bisection algorithm.
ELWIDCALC
No additional local variables
ENTRAIN
DFTC: a variable equal to the derivative of the second Bridge approximation function with respect to TCI, evaluated at RELN.
FTC: a variable equal to the evaluation of the second Bridge approximation function.
RELN: the natural logarithm of the grain Reynolds number.
TCI (I,J): the critical shear stress necessary to entrain grain size I, of mineral J (Pa).
TCI: a variable signifying the estimated critical shear stress necessary to entrain a grain of size I (assuming mineral species J) (Pa).
TCIN: a variable signifying the estimated critical shear stress (Pa) necessary to entrain a grain of size I (assuming mineral species J) (Pa).
THETLN: the natural logarithm of the critical Shields parameter necessary to entrain a grain out of a homogeneous bed under fully rough conditions.
ENTRAINH
D50: the mean grain diameter, by volume, in the bed (m).
ERODEP
MSACTIVE: a variable, equal to the present mass per unit area of the active layer (kg/m2).
MSDEP: the total mass per unit area deposited (kg/m2).
MSEROACT: the mass per unit area of the parent material (kg/m2).
NEWMSACT (I,J): the new mass per unit area of grain size I of mineral species J (kg/m2).
TEMPWEIT: the weight percentage of the current grain size-density interval in the layer underlying the active layer.
TOTNEWMS: the total new mass per unit area of the active layer (kg/m2)
EXVERT
Y: variable indicating the vertical distance between the free surface of the flow and the flow level at which calculations are performed (m).
FLDTA
A, A2: cross sectional area of the node (m2)
CONST, CONST2: variables equal to parts of the gradually varied flow equation.
CONV: a conversion factor equal to LOG (0.5).
COUNT: iteration step counter of the gradually varied flow equation loop.
DCONT: the downstream distance of the control section: the point along the reach at which  the depth is known or controlled (m).
DELTAD: the depth correction of the uniform, control section flow depth estimate (m).
DMIN1: the flow depth at the previous node (m).
D2: the estimated flow depth at the next node (m).
DUR: the time, in hours, up to which the water discharge QO is operative. QO is assumed to be constant over this period unless otherwise specified.
DUMMY: a variable used to read in/skip over descriptive lines in the input file.
END: the last node at which flow depth calculations must occur in the upstream or downstream direction.
FR: the Froude number of the section.
H1, H2, H21: total head estimates (m).
K: loop counter signifying whether calculations of nodes upstream (K=2) or downstream from the control section (K=l) are occurring.
NEXT: the next timestep at which a flow discharge will be read in.
NODE: the node number for which depth estimates are being made.
P, P2: wetted perimeter of the section (m).
POS: the downstream of the section (m).
RCH, RCH2: hydraulic radius of the section (m).
RLONG: the total length of the s or reach (m).
SF1, SF2: friction slope of the flow at a node.
SIGN: the sign of the depth conection of the uniform, control section flow depth estimate.
SLOP: the bed slope at the control section (m/m).
STEP: a variable indicating if calculations are going upstream (-1) or downstream (+1).
V, V2: longitudinal velocity (m/s).
YCRIT: the critical flow depth, or the depth at which the flow changes from sub- to supercritical, at the control section (m).
DH and DQ: interpolation increments for discharge and head water level for model input values.
FLDTA- hydrograph option *this subroutine may need to be adapted to suit data type*
HLEVEL: head water level of hydrograph in meters. 
QTCOUNTER: counts how many interpolations have been performed between real hydrograph in translation to model hydrograph.
YEAR; MONTH; DAY; HOUR; MINUTE; SECOND: variables that store translated numerical data from character data date time stamps typical of hydrograph data
GSDIST
Al, A2, A3, A4, A5: constants used as series coefficients in the error function approximation.
ALIMIT: lower limit of the grain size distribution (phi).
BLIMIT: upper limit of the grain size distribution (phi).
CLLENGTH: interval length between grain size intervals (phi).
CONV: conversion factor equal to LOG (0.5)
DIST: the distance of the grain size interval to the mean (standard deviations).
DISTMIN1: variable equal to the distance from the mean of ALIMIT (standard deviations).
FLAG: set when a portion of the distribution larger than the mean grain size is being treated.
FO: probability of the lower grain size interval bound.
FI: probability of the upper grain size interval bound.
H, H1: variables equal to a part of the error function equation.
P: constant used in the error function approximation.
PISQR: constant equal to the square root of 2p.
T: variable used in the error function approximation.
LOGDIST
V: cross-sectionally averaged, longitudinal velocity (m/s).
Y: variable equal to the depth beneath the free surface of the flow (m).
PLOTFILES
NUMNODES: the number of node points along the reach at which information is printed.
POS: a local variable equal to the downstream distance of node JN.
TEMP: a temporary variable.
READRSTRT
No additional local variables.
RESULT
INSTTHIC: total instantaneous thickness of erosion or deposition at node JN (m).
TOTALTHICK: total thickness eroded or deposited over the length of the run at node JN  (m).
TOTSBLWT: total instantaneous bedload transport rate at node JN (kg/m/s).
TOTSBLWTA: total time-averaged bedload transport rate at node JN (kg/m/s).
TOTSCVELI: total instantaneous suspension load transport rate at node JN (kg/m/s).
TOTSCVELIA: total time-averaged suspension transport rate at node JN (kg/m/s).
TOTSDBDT: the total instantaneous mass flux between the moving bed layer and the active layer for node JN (kg/m2/s).
TOTSDBDX: the total time-averaged mass flux between the moving bed layer and the active layer for node JN (kg/m2/s).
TOTSPRTHICK: the total thickness of the material eroded or deposited over the timestep at node JN (m).
TOTSTFLUX: the total instantaneous mass flux between the suspension load and the active layer for node JN (kg/m2/s).
TOTSXFLUX: the total time-averaged mass flux between the suspension load and the active layer for node JN (kg/m2/s).
SETTLE
RHS: the Dietrich functional evaluation of the dimensionless grain diameter.
SUSP
AV: a variable equal to the height of the moving bed layer at the current JN (m).
B1: a variable equal to bedload + saltation load for the current I, J, JN combination (kg/m/s).
BCON: variable equal to reference concentration for the I, J, JN combination (kg/m3).
BWT: variable equal to bedload transport rate for the current I, J, JN combination (kg/m/s).
BWTA: variable equal to bedload mass per unit area for the I, J, JN combination (kg/m2).
CA: a value equal to the reference concentration of the I, J, JN combination (kg/m3).
FACTOR: a multiplying factor used in the mass flux reduction loop.
HILIM: the upper limit to the correct FACTOR value.
IN: variable equal to JN-1.
IV: the depth level counter.
LOWLIM: the lower limit to the correct FACTOR value.
MSACT (JN): the mass per unit area in the active layer (kg/m2).
NEWSIGN: a variable equal to the sign of MSPRTHK for the current I, J, JN combination, AFTER the reduction loop was entered for the first time.
PHI1: weighting factor in the modified Preissman mass balance scheme for the last node.
REDUC: a flag signifying whether (=l) or not (=0) the mass flux reduction loop has been entered.
SIGNMSP: a variable equal to the sign of MSPRTHK for the current I, J, JN combination, BEFORE the reduction loop was entered for the first time.
SWT: variable equal to saltation load for the current I, J, JN combination (kg/m/s).
SWTA: variable equal to the saltation load per unit area for the I, J, JN combination (kg/m2).
W: a variable equal to the settling velocity of the current I, J, JN combination (m/s).
Y: the depth of level L beneath the free surface (m).
Z: the Rouse number for the current I, J, JN combination.
SUSP1
NN: the number of depth levels.
SVELA
AK: variable equal to the bed roughness at node JN (m).
COUNT: a variable that keeps track of the total number of iterations in the bisection solution scheme.
DR: the hydraulic radius due to skin friction (m).
DR1, DR2, DR3: upper, lower and best guess estimates of the reduced hydraulic radius (m).
F3: functional evaluation of the estimate of hydraulic radius.
RCH: hydraulic radius of the section.
SF: friction slope of the flow at a node (m/m).
V: cross-sectionally averaged longitudinal velocity (m/s).
TABLE
ERROR: a variable used to test whether or not the time average arrays need to be reset.
TIME: a variable used to test whether or not the time average arrays need to be reset.
TIMEPRES: a variable equal to the elapsed time since run start.
TIMEMINI: a variable equal to the elapsed time at the previous timestep.
TOPRINT
DUM: a dummy variable.
RLONG: the total length of the stream or reach (m).
TECTONICS
EFOLDX: e-folding distance for exponential subsidence (OPTTECT=2).
TECRATE: tectonic subsidence rate (or uplift rate in meters/timestep)
ZEROLOC: downstream pivot point for beam of subsidence.
TRIB
TMP: temporary counter used in AWPM (I,J,JN) adjustments to bed from tributary material.
TMSDEP (JN): tributary material mass per unit area deposited (kg/m2)
TMSFLUX: mass flux into trunk stream from tributary for each size I and mineral J.
TTHICKM: thickness of material coming in from tributaries.
TURB
A1, A2, A3, A4, A5: local constants used as series coefficients in the error function approximation.
DELTO: width of the shear stress intervals (Pa).
DELTOH: half interval length (Pa).
DIST: the distance of the shear stress interval midpoint to the mean (standard deviations).
DISTMIN1: variable equal to the distance from the mean to TOLO (standard deviations).
FO: probability of the lower grain size interval bound.
FI: probability of the upper grain size interval bound.
P: constant used in the error function approximation.
STDEV: the standard deviation of the (log)normal shear stress distribution.
T: variable used in the error function approximation.
TOHI: upper limit of the shear stress distribution (Pa).
TOLO: lower limit of the shear stress distribution (Pa).
TOMEAN: the mean bed shear stress at the current node (Pa).
WRFLDTA
DCONT: the downstream distance of the control section; the point along the reach at which the depth is known or controlled (m).
VV: cross-sectionally averaged longitudinal velocity at node JN (m/s).
Variable Index 
A
A: COMMON BLOCK (CB)
Al, A2, A3, A4, A5: GSDIST, TURB
A, A2: FLDTA
AA: DRANGE
ACTIVE: CB
AK: SVELA
ALIMIT: GSDIST
AV: SUSP, SUSP1
AVDENS(JN): CB
AWP(I,J): AMAIN90, GSDIST
AWPM I,J,JN): CB
AWPOM (I,J,JN): CB
AWPT (I,J): CB
B
B: CB
B1: SUSP
BB: DRANGEH
BCON: SUSP
BCOUNT: CB
BEDCON (I,J,JN): CB
BETA: CB
BLIMIT: GSDIST
BLWTM (I,J,JN), BLWTAM (I,J,JN): CB
BLWTO (I,J): CB
BLWTT (I,J): CB
BLWTTO (J): CB
BREAK: CB
BWT, BWTA: SUSP
C
CA: SUSP, SUSP1
CC:  DRANGE
CLLENGTH: GSDIST
CONV: TURB
COVAR: CB
CON1, CON2: CB
CONC (L), CONCM (I,J,JN): CB
CONIM (I,J,JN): CB
CONST, CONST2: FLDTA
CONV: BEDLOAD, FLDTA, GSDIST
COUNT: DRANGE, DRANGEH, ENTRAIN, FLDTA, SVELA
COUNTER:  CB
CVELM (I,J,JN): CB
D
D2: FLDTA
D5O: DRANGEH, ENTRAINH
D50M (JN): CB
DATETIME: AMAIN
DBDX (I,J,JN): CB
DCHM (JN): CB
DCONT: FLDTA, WRFLDTA
DELTAD: FLDTA
DELTO, DELTOH: TURB
DELWID: CB
DELX: CB
DELY (JN): CB
DENDM: CB
DELT: CB
DFD: DRANGE, DRANGEH
DFTC: ENTRAIN
DH: FLDTA
DI: DRANGE, DRANGEH
DEFEND (I,J): CB
DIMDL: DRANGE, DRANGEH, SETTLE
DIMID (I,J): CB
DIMIN (J): CB
DIMSV: DRANGE, DRANGEH
DIN: DRANGE, DRANGEH
DINIT: CB
DIST: GSDIST, TURB
DISTMIN1: GSDIST, TURB
DMAX (K,J,JN), DMIN (K,J,JN): CB
DMIN1: FLDTA
DQ: FLDTA
DR, DR1, DR2, DR3: SVELA
DUM: TOPRINT
DUMMY: FLDTA, SUSP1
DUR: FLDTA
DX: DRANGEH
E
EFOLDX: CB
ELEVM (JN): CB
ELEVOM (JN): CB
ELEVTOP: CB
END: FLDTA
EPSY (L,JN): CB
ERROR: TABLE
F
F: DRANGEH
FACTOR: SUSP
F1(NRUN), F2(NRUN), F3(NRUN): AMAIN90
FD: DRANGE, DRANGEH
FI: GSDIST, TURB
FILENAME: AMAIN90
FLAG: GSDIST, TURB
FLOWEL: CB
FMID: DRANGEH
FO: GSDIST, TURB
FOUT: CB
FR: FLDTA, WRFLDTA
FTC: ENTRAIN
FUN1: DRANGEH
G
G: CB
GRAD: CB
GRADZERO: SUSP
H1, H2, H21: FLDTA
H, HI: GSDIST, TURB
HLEVEL: FLDTA
HILIM: SUSP
I
ICONT: CB
IF: BEDLOAD
IMAX (K,J,JN): CB
IMIN (K,J,JN): CB
IN: SUSP, LOOP VARIABLE (LV) 
INSTTHIC: RESULT
IRSTRT (NRUN): AMAIN90
IRUN: AMAIN90
IS: BEDLOAD
IV: SUSP, SUSP1
J
J: LV
JN: LV
K
K: DFIFTY, DISCR, FLDTA, LV
K1, K2: DISCR
KR: CB
KK1: CB
L
L: LV
LOWDEN: DFIFTY
LOWLIM: DFIFTY, SUSP
M
M: CB
MANN (JN): CB
MANN1: AMAIN90
MEAN (J): CB
MINISIF: BEDLOAD
MINSORT: AMAIN90
MSACT (JN): SUSP
MSACTIVE: ERODEP
MSDEP: ERODEP
MSEROACT: ERODEP
MSPRTHK (I,J,JN): CB
N
N: CB
NEWMSACT (I,J): ERODEP
NEWSIGN: SUSP
NEXT: CB, FLDTA
NK: CB
NN: SUSP1
NODE: FLDTA
NRP: DISCR
NRUN: AMAIN90
NSIGMA: CB
NSIZE (J): CB
NT: LV
NTFIN, NTSTAR: AMAIN90
NTIM: CB
NT'OT: CB
NUMNODES: PLOTFILES
NWAVDENS (JN): CB
O
OPTABRAD: CB
OPTBINP: CB
OPTDISC: CB
OPTCONT: CB
OPTDIST (J): CB
OPTEL: CB
OPTERO: CB
OPTFUNC: CB
OPTHIDE: CB
OPTLNODE: CB
OPTSALT: CB
OPTSINP: CB
OPTSUSP: CB
OPTTECT: CB
OPTTRIB: CB
OPTURB: CB
OPTWID: CB
P
P: GSDIST, TURB
P, P2: FLDTA
PHI: CB
PHI1: SUSP
PISQR: GSDIST, TURB
PLOTPRINT: CB
POS: FLDTA, PLOTFILE, WRFLDTA
PRINT (J), J=1,7: CB
PROPT (K): CB
PROPTJ (J): CB
PROPVOL: BEDLOAD
PRTHICK (I,J,JN): CB
PSI: CB
Q
QCOUNT: CB
QFAC: CB
QN (JN): CB
QO: CB
QSMOD: CB
QSPMOD: CB
QTCOUNTER: FLDTA
QTMOD: CB
RCH: SVELA, WRFLDTA
RCH, RCH2: FLDTA
REDUC: SUSP
REG: DRANGE, 
RELN: ENTRAIN, ENTRAINH
RENEW: CB
RHO: CB
RHS: SETTLE
RLONG: CB, FLDTA, TOPRINT
RTBIS: DRANGEH
RUN: CB
RUNNAME: CB
SBLWT (I,J,JN), SBLWTA (I,J,JN): CB
SCVELI (I,J,JN), SCVELIA(I,J,JN): CB
SDBDT (I,J,JN), SDBDX(I,J,JN): CB
SETV (I,J): CB
SETVEL: DRANGE, DRANGEH
SF: SVELA, FLDTA
SF1, SF2: FLDTA
SHVEL (K,JN): CB
SHVELC (I,J,JN): CB
SIGMA (J): CB
SIGN: FLDTA
SIGNMSP: SUSP
SLOP: FLDTA
SLOPE: CB
SLWTAM (I,J,JN), SLWTM (I,J,JN): CB
SPRTHICK (I,J,JN): CB
STDEV: TURB
STDEV (J): CB
STEP: FLDTA, RESULT
STFLUX (I,J,JN): CB
STRDATA (Q,JN), Q=1,6: CB
SUM: DFIFTY
SUR: DRANGEH
SVELM (JN): CB
SVELUM (JN): CB
SWF, SWTA: SUSP
SXFLUX(I,J,JN): CB
T
T: GSDIST, TURB
TC (I,J,JN): CB
TC1 (I,J): ENTRAIN
TCI: ENTRAIN, ENTRAINH
TCIN: ENTRAIN
TDBDX (JN): CB
TECRATE: CB
TEMP: DISCR, PLOTFILES
T'EMPWEIT: ERODEP
TFILE: CB
THALFB: CB
THALFQ: CB
THET50 (JN): CB
THETAC: CB
THETLN: ENTRAIN, ENTRAINH
THICKM (JN): CB
TIME: AMAIN90, RESULT, TABLE
TIMEMIN1, TIMEPRES: TABLE
TMP: BEDLOAD, TRIB
TMSDEP (JN): TRIB
TMSFLUX: TRIB
TO (K,JN): CB
TOHI, TOLO: TURB
TOMEAN: TURB
TOTALTHICK: RESULT
TOTNEWMS: ERODEP
TOTSBLWT, TOTSBLWTA: RESULT
TOTSCVELI, TOTSCVELIA: RESULT
TOTSDBDT, TOTSDBDX: RESULT
TOTSPRTHICK: RESULT
TOTSTFLUX,TOTSXFLUX: RESULT
TOTTHICK (I,J,JN): CB
TOXM (JN): CB
TPRINT: CB
TTHICKM: TRIB
TXFLUX (JN): CB
V
V: LOGDIST, SVELA
V, V2: FLDTA
VEL (L,JN): CB
VISC: CB
VISKIN: CB
VOL: BEDLOAD, DFIFTY
VOLCON: CB
VOLTOT (K): DFIFTY
VV: WRFLDTA
w
W: SUSP, SUSP1
WEITFAC (JN): CB
WIDM (JN): CB
WIDTOP: CB
x
X1, X2: DRANGEH
XACC: DRANGEH
XFLUX (I,J,JN): CB
XPRINT: CB
Y
Y: EXVERT, LOGDIST, SUSP, SUSP1
YCRIT: FLDTA
YMAX (JN): CB
Z
Z: SUSP, SUSP1
ZEROLOC: CB
References
Bagnold, R. A., 1973, The nature of saltation and “bedload” transport in water, Proceedings of the Royal Society London, v. 322A, p. 473-504.
Bridge, J. S., and Dominic, D. F., 1984, Bedload grain velocities and sediment transport rates: Water 
Resources Research, v. 20, p. 476-490.
Egiazaroff, J. V., 1967, Sediment transportation mechanics - initiation of motion, Journal of the Hydraulics Division, ASCE, v. 93, p. 281-287.
Komar, P. D., 1987, Selective entrainment by a current from a bed of mixed sizes- a reanalysis: Journal of Sedimentary Petrology, v. 57, p. 203-211.
Henderson, F. M., 1966, Open Channel Flow, New York, Macmillan, p. 522.
James, C. S., 1990, Prediction of entrainment conditions for nonuniform, noncohesive sediments, Journal of Hydraulic Research, v. 28, p. 25-41.
van Niekerk, A., Vogel, K.R., Slingerland, R.L., and Bridge, J.S., 1992, Routing of heterogeneous sediments over movable bed: model development: Journal of Hydraulic Engineering, v. 118, 2, p. 246-262.
Vogel, K. R., 1991, A field and mathematical modeling study of the 2.8 billion-year-old Ventersdorp contact reef, Witwatersrand basin, South Africa, Ph.D. thesis, Pennsylvania State University, University Park.
Vogel, K. R., van Niekerk, A., Slingerland, R. L., and Bridge, J. S., 1992, Routing of heterogeneous sediments over movable bed: model verification and testing: Journal of Hydraulic Engineering, v. 118, 2, p. 263-279.