Surface Water Routing:Overland Flow Routing
5.2.1 Overland Flow Routing Formulation
Overland flow in GSSHA employs the same methods described for 1-D channel routing, except the calculations are made in two dimensions. Flow is routed in two orthogonal directions in each grid cell during each time step. The watershed boundary represents a no flow boundary for the overland flow routing and when a grid cell lies on the watershed boundary, flow is not routed across the boundary. There is an option in GSSHA to impose head boundaries along the watershed boundary and compute flows across the overland boundary. In GSSHA, Δx = Δy. Inter-cell fluxes in the x and y directions, p and q, respectively, are computed in cell ij from the depth, dij, at the nth time level using the Manning equation for the head discharge relationship in the x and y directions, respectively, as
Depths in each cell are calculated at the n+1 time level based on the flows for each cell (Julien and Saghafian, 1991):
In addition to this original formulation in the CASC2D model, two additional methods of solving the equations have been added, an alternating direction explicit scheme (ADE) and an ADE scheme with an additional predictor-corrector step (ADE-PC) (Downer 2002a; Downer et al., 2000). Both the ADE and ADE-PC methods employ the up-gradient difference technique, Equation 3, for flows in the upstream direction (Downer, 2002a). Fluxes other than inter-cell fluxes, direct evaporation (DET), infiltration, exfiltration, are accounted for before overland routing is computed.
In the ADE method, inter-cell flows are first calculated in the x direction according to Equation 6. Depths in each row are updated based on the flows in the x direction:
Inter-cell flows in the y direction are computed using the updated depths:
Depths in each column are updated based on the flows in the y direction:
With the ADE-PC method additional steps are added to improve accuracy and stability. As before, during each sweep, by rows or by columns, an estimate of heads is made based on the calculated flows, Equations 9 and 11. Next, using the updated depths, updated estimates of flow are computed at the n+1 time level
The original flows and the updated flows are then averaged to come up with an estimate of flows for the time step:
These flows are then used to update the original depths, Equations 9 and 11. This procedure is essentially the MacCormack method (MacCormack, 1969) except up-gradient differences are used in both the predictor and corrector steps. A similar method was successfully implemented by Wang and Hjelmfelt (1998).
In benchmark tests using the three methods: original explicit EXPLICIT, ADE, and ADE-PC, in an contrived watershed consisting of two converging planes (open book), the Goodwin Creek Experimental Watershed (GCEW) (Senarath et al., 2000), and Poplar Creek (Downer et al., 2002a), the ADE and ADE-PC methods ran with significantly larger time steps (Downer et al., 2000). Depending on the test case, time steps could be increased from 20% to 240% with commensurate decreases in simulation times.
The routing scheme is selected using the OVERTYPE card. The default value is ADE. The most efficient scheme to use depends on the particular watershed. The ADE-PC scheme can generally handle rougher terrain and typically requires less smoothing of the DEM, but the additional computational steps result in greater computation time, unless use of the ADE-PC scheme permits substantially greater time steps than with one of the other two methods. For smoothed DEMs or in watersheds with smoother terrain, the ADE and EXPLICIT methods usually can be employed, with a resulting savings in execution time.
To improve stability, the timestep in both the EXPLICIT and ADE schemes is variable. The time step is not variable in the ADE-PC method. The addition of a variable timestep has allows for increased overall model timestep when using either the EXPLICIT and ADE methods. These new methods have not been benchmarked against the ADE-PC method.
In GSSHA version 7.1.1 and higher, an inertial formulation of the flow equations, as described in Bates et al. (2010), is included as an option. To use this option include the OVERLAND_MOMENTUM card in your project file along with ADE specified with the OVERTYPE card . This option may be helpful for simulations where overland hydrodynamics are important, such as flooding due to a storm surge. For these cases stability may be increased due the inclusion of inertial terms and a better accounting of the effects of friction. The method may allow longer time steps and potentially reduce simulation times. In addition, in GSSHA v7.1.1 and higher, a strict overland time step for the diffusive wave formulation, as described in Bates et al. (2010) has been included. This option is specified by including the OVERLAND_STRICT_DT card in the project file. This option also only works with the ADE method. While use of this option may increase stability for high overland flows, the timestep is very restrictive compared to the normal criteria used for the ADE method. If instability arises with use of the ADE method, it may be better to reduce the overall timestep, use the inertial formulation for flow OVERLAND_MOMENTUM, or use the ADE-PC method, in lieu of the OVERLAND_STRICT_DT option. For versions 7.14 and beyond, the user can control the time step coefficient for the momentum formulation by specifying a real value between 0.0 and 1.0 after the OVERLAND_MOMENTUM card. The default value is 0.2. Larger values increase speed, lower values increase stability. The user is referred to Bates et al. (2010) for details: Bates, P. D., Horritt, M. S., and T. J. Fewtrell (2010) A simple inertial formulation of the shallow water equations for efficient two-dimensional flood inundation modelling. J. Hydro. 387 (2010) 33-45.
Special overland flow cell types that may defined and result in a special type of overland routing are wetlands described in the GSSHA primer here: http://www.gsshawiki.com/Wetlands:Wetlands
And overland lakes, which are lakes not connected to the stream network. Overland lakes are specified with the OVERLAND_LAKE_MAP card that specifies the overland lake mask. The overland lake mask is a map with zeros and integer values, like the MASK_MAP except the integer values occur where overland lakes are desired. This mask map can be created in WMS by making a coverage in WMS for overland lakes and then assigning integer values to polygons that define the lakes, and then creating an index map from that WMS coverage. The significance of an overland lake is that no overland routing occurs within the overland lake cells. Other overland processes, such as infiltration, still occur.
5.2.2 Overland Flow Hydraulic Roughness
The GSSHA model requires that Manning roughness coefficients be assigned to every cell in the watershed mask. There are three ways to specify the hydraulic roughness of the overland flow planes in GSSHA. The first method is to apply a constant value over the entire watershed through the MANNING_N project file card. The second method is to use the MAPPING_TABLE to assign roughness coefficients using tabled values referenced to an index map. The third method is to produce a GRASS ASCII map of roughness coefficient, and provide the name of this map to GSSHA using the ROUGHNESS project file card. Table 8 provides typical values of the Manning roughness coefficient for overland flow over various surfaces:
|Land Use or Cover|| Recommended
|Concrete or asphalt||0.011a|| 0.01-0.013a|
|Bare clay-loam (eroded)||0.02a||0.012-0.033a|
|Bare field – no residue||0.05a||0.006-0.16a|
|Grass and pasture||-||0.05 – 0.15a|
|Clover||-||0.08 – 0.25a|
|Small grain||-||0.1 – 0.4a|
|Row crops||-||0.07 – 0.2a|
|Grass (bluegrass sod)||0.45a|| 0.39-0.63a|
|Short grass prairie||0.15a||0.10-0.20a,d|
Table 8 - Values of overland flow roughness coefficient Notes: aEngman (1986), bDowner (2002b), cSenarath et al (2000), dHEC (1985)
These values should be considered guidelines. Manning roughness coefficients are typically assigned from literature values and then adjusted through calibration. Additional sources of literature values are Liong et al. (1989), Engman (1986), and Ree et al. (1977). If calibrated values differ significantly from published values, there may be appropriate justification.
Typical Manning roughness coefficient values for open channel flow are considerably smaller than overland flow values because of deeper flow depths in the channel. In cases where the overland flow may become very deep, such as simulating a tidal surge, the user may want to use depth varying overland roughness. GSSHA is currently formulated to calculate the overland roughness at any depth (d) using the formula
where nd is the Manning roughness at depth (d), n0 is the specified roughness in the ROUGHNESS mapping table, and B is an exponent specified in the ROUGHNESS mapping table. The default value for B is zero, which results in a static value roughness with depth. To modify this relationship, a value of B is specified for each roughness category in the ROUGHNESS table. The value is listed in the table after the specified roughness values, n0. So that the table will look like
ID DESCRIPTION1 DESCRIPTION2 ROUGH EXPONENT
1 Roughness ID 0.300000 0.5
Positive values of B result in the roughness decreasing with depth, common; negative values result in roughness increasing with depth, uncommon. A typical value is 0.5. The exponent may require calibration.
5.2.3 Runoff Retention
Natural land surfaces contain micro-topography, small depressions, that retain water prior to runoff. The water held in the grid cell, or retention storage, never becomes direct runoff and can only be removed from the land surface as infiltration or direct evaporation. In certain regions, the retention storage can be significant. Retention storage is input as a depth (mm) in each grid cell and may be optionally input to GSSHA as:
- a uniform value using the RETENTION card,
- as a table of values related to index maps by using the MAPPING_TABLE project card, and the RETEN_DEPTH card without specifying a file name or
- as an ASCII GRASS map through the use of the RETEN_DEPTH project file card with a specified file name.
5.2.4 Specifying Initial Depths on the Watershed
Initial depths on the overland flow plane at the beginning of the simulation may be specified by use of the OV_INIT_DEPTH project file card. This card is used to specify an ASCII GRASS map containing initial overland depth values, m. The initial depth can also be specified with the READ_OV_HOTSTART card. A map of the final overland flow values from a simulation can be written out using the WRITE_OV_HOTSTART card. This map can then be used to hotstart subsequent simulations. A single value of initial water surface elevation (m) can be specified using INIT_ELEV_HEAD card. The water surface elevation in every cell with an elevation less the specified elevation will be set to this value and the depths computed from the specified elevation and the cell elevation.
5.2.5 Simulations without Channel Routing
It might be desirable to perform simulations without channel routing in very small watersheds lacking a defined channel network, or in the beginning stages of developing a GSSHA model of a watershed. It is always prudent to build a GSSHA model one process at a time. Getting just the overland flow portion of the model to run is always the first step in building a model. For these reasons, GSSHA may be used without channel routing. Because the channel network normally provides the outlet point, the overland cell containing the outlet must be specified during simulations without channel routing. Normally this grid cell will have the lowest elevation in the watershed. The row and column containing the outlet grid cell are specified using the OUTROW and OUTCOL cards, while the slope of the outlet grid cell is specified using the OUTSLOPE card. If channel routing is enabled through the inclusion of CHAN_EXPLIC in the project file, the OUTROW, OUTCOL, and OUTSLOPE cards are not required, and ignored if present.
- 5 Surface Water Routing