Difference between revisions of "Surface Water Routing:Channel Routing"

From Gsshawiki
Jump to: navigation, search
(→‎5.1.1 Explicit Channel Routing Formulation)
(→‎5.1.1 Explicit Channel Routing Formulation)
 
(16 intermediate revisions by 2 users not shown)
Line 24: Line 24:
 
Several modifications were made in the implementation of the channel routing scheme to accommodate groundwater/channel interactions.  These modifications permit continuous interaction between channel nodes and the saturated groundwater cells.  The channel routing scheme was modified to allow water to remain in the channel after channel routing ends, and for water to be present in the channel when channel routing begins.  Because groundwater may discharge to the stream at anytime, channel routing is initiated anytime a minimum amount of water is in the channel network.  If the channel routing scheme indicates there is no flow in the channel, channel routing is halted during periods outside precipitation events.  Fluxes between the stream and the groundwater are still computed and adjustments to the stream volumes are made without routing.  If groundwater discharges to the stream, channel routing will resume, but at the groundwater time step, which is typically larger than the channel routing time step.
 
Several modifications were made in the implementation of the channel routing scheme to accommodate groundwater/channel interactions.  These modifications permit continuous interaction between channel nodes and the saturated groundwater cells.  The channel routing scheme was modified to allow water to remain in the channel after channel routing ends, and for water to be present in the channel when channel routing begins.  Because groundwater may discharge to the stream at anytime, channel routing is initiated anytime a minimum amount of water is in the channel network.  If the channel routing scheme indicates there is no flow in the channel, channel routing is halted during periods outside precipitation events.  Fluxes between the stream and the groundwater are still computed and adjustments to the stream volumes are made without routing.  If groundwater discharges to the stream, channel routing will resume, but at the groundwater time step, which is typically larger than the channel routing time step.
  
−
Because GSSHA uses a finite volume representation of channel flow the standard stability criterion, Courant number < 1.0, does not strictly apply.  Maintaining stability is dependent on volume changes during each time step.  Experience with the scheme indicates that stability can be maintained with a time step limitation that keeps the maximum Courant number everywhere in the network less than 1/6.  The user can specify this value using the '''MAX_COURANT_NUM''' card, which requires a real number value.  The default is 0.2.  Groundwater and over bank fluxes can induce instability and additional controls in the channel routing scheme are added to further reduce instability.  If the channel routing scheme becomes unstable (negative depth occurs in one or more cells), despite the more restrictive control on the Courant number, the time step is reduced and the channel routing calculations are repeated.  The channel routing time step may be automatically reduced to a value as low as 1/1000 of a second.  Allowing the time step to become very small during periods of sharp transition allows a larger overall model time step to be used.  For each call of the channel routing function the overall model time step is used as the beginning channel routing time step.  The time step is only reduced when the stability controls are activated, and then only for that call of the channel routing routine.  The current GSSHA version employs a predictor corrector scheme with Picard iterations on the flow area.  This has made the scheme much more stable when reservoirs or overbank flooding are simulated.
+
Because GSSHA uses a finite volume representation of channel flow the standard stability criterion, Courant number < 1.0, does not strictly apply.  Maintaining stability is dependent on volume changes during each time step.  Experience with the scheme indicates that stability can be maintained with a time step limitation that keeps the maximum Courant number everywhere in the network less than 1/6.  The user can specify this value using the '''MAX_COURANT_NUM''' card, which requires a real number value.  The default is 0.04.  Groundwater and over bank fluxes can induce instability and additional controls in the channel routing scheme are added to further reduce instability.  If the channel routing scheme becomes unstable (negative depth occurs in one or more cells), despite the more restrictive control on the Courant number, the time step is reduced and the channel routing calculations are repeated.  The channel routing time step may be automatically reduced to a value as low as 1/1000 of a second.  Allowing the time step to become very small during periods of sharp transition allows a larger overall model time step to be used.  For each call of the channel routing function the overall model time step is used as the beginning channel routing time step.  The time step is only reduced when the stability controls are activated, and then only for that call of the channel routing routine.  The current GSSHA version employs a predictor corrector scheme with Picard iterations on the flow area.  This has made the scheme much more stable when reservoirs or overbank flooding are simulated.
  
 
Several new channel routing features have been added to GSSHA since the publication of the original manual.  These features are described in the following SWWRP Tech Note.
 
Several new channel routing features have been added to GSSHA since the publication of the original manual.  These features are described in the following SWWRP Tech Note.
Line 31: Line 31:
 
   
 
   
 
[[Image:chan_route.jpg|none|frame|Figure 5 – Explicit channel routing scheme]]
 
[[Image:chan_route.jpg|none|frame|Figure 5 – Explicit channel routing scheme]]
 +
 +
In GSSHA version 7.1.4 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 '''CHANNEL_MOMENTUM''' card in your project file.  This option may be helpful for simulations where channel 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.  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 '''CHANNEL_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.
  
 
==5.1.2 Boundary Conditions ==
 
==5.1.2 Boundary Conditions ==
Line 396: Line 398:
 
* Rule Curve
 
* Rule Curve
  
−
All culverts, for example, lie below a weir formed by a roadway or other embankment. In GSSHA, a hydraulic structure link can consist of any combination of the above.  Other structures consist of a number of culverts, some with different dimensions, invert elevations, etc.  This allows the user to simulate complex hydraulics from multiple drainage elements.  
+
All culverts, for example, lie below a weir formed by a roadway or other embankment. In GSSHA, a hydraulic structure link can consist of any combination of the above.  Other structures consist of a number of culverts, some with different dimensions, invert elevations, etc.  This allows the user to simulate complex hydraulics from multiple drainage elements.
  
 +
Hydraulic structures are point models that exist at a channel intersection, at the last node for the upstream link and at the first node of the downstream link, which occupy the same point.  The upstream conditions are taken from the second to last node in the upstream link.  The downstream conditions are taken from the second node in the downstream link.  It is important to note that the bed elevation of the upstream condition is the bed elevation at the junction plus half the difference between the bed elevation at the junction and the bed elevation at the next to last node of the upstream link.  The same applies to the downstream condition.  To avoid numerical issues, it is advisable to the structure invert be at or above the bed elevation at the respective end.
  
 
== [[Horizontal Weir]] ==
 
== [[Horizontal Weir]] ==
Line 470: Line 473:
 
</pre>
 
</pre>
 
|}
 
|}
 +
 +
As noted above, numerical issues can be avoiding by locating you inverts above the bed elevations at the upstream and downstream condition locations, which are at the mid-point between the junction (structure) node and the next upstream or downstream node, depending on whether it is the UPINVERT or the DOWNINVERT.
  
 
== Box Culvert ==
 
== Box Culvert ==
Line 494: Line 499:
 
</pre>
 
</pre>
 
|}
 
|}
 +
 +
As noted above, numerical issues can be avoiding by locating you inverts above the bed elevations at the upstream and downstream condition locations, which are at the mid-point between the junction (structure) node and the next upstream or downstream node, depending on whether it is the UPINVERT or the DOWNINVERT.
  
 
== Rating Curve ==
 
== Rating Curve ==
Line 571: Line 578:
 
=====5.1.4.1.4.3 - Reservoir channel links=====
 
=====5.1.4.1.4.3 - Reservoir channel links=====
  
 +
Reservoirs are defined by specifying that a reservoir is present with the '''LAKE''' card in the .cif file.  The minimum, initial, and maximum reservoir elevations are specified with the '''MINWSE''', '''INITWSE''', and '''MAXWSE''', cards respectively.  The number of cells within the maximum reservoir boundary is specified with the '''NUMPTS''' card, followed by a list of the cells in the reservoir domain, defined by row and column values.  Reservoirs are dynamic features that can increase/decrease in size, take/release stream nodes, and overland flow cells as defined in Non-ortho Channels Technical Note [[media:channels.pdf|Non-ortho Tech Note]].  The resevoir will not go into any cells not listed in the '''NUMPTS''' list, but water may flow out of the reservoir back on the overland flow plane.  All reservoirs must have a hydraulic structure defined as the reservoir outlet.  Other important restrictions are that the reservoir cannot rise above the '''MAXWSE''', move into cells not defined as part of the reservoir, join with other reservoirs, or take in stream nodes not directly in the same branch as the reservoir.  Any of these conditions will cause the program to stop, crash, or provide non-sensical results.  Large mass balance errors are indicative of problems with reservoirs violating one or more of these conditions.  Reservoirs can be kept within the desired confines using embankments.  Water can be passed through an embankment holding a reservoir using lowspots along the embankment or using overland rating curves.  Water can also pass from one lake to another using these features.  Although lakes cannot exceed the maximum water level, they can go dry.  If the water level falls below the '''MINWSE''' it will stop discharging but water can still seep out of the lake and/or evaporate.
  
 +
Output from all the reservoirs can be specified with the '''LAKE_OUTPUT''' card specified in the project file.  Lake elevations and volumes, respectively, will be written for each lake, in the order they are listed in the .cif file, for each time period specified in the HYD_FREQ card, in the file specified with the '''LAKE_OUTPUT''' card.
  
 +
Within the .cif file, a standard lake input would look like the example below: 
  
  
Line 580: Line 590:
 
<pre>
 
<pre>
 
LINK        7
 
LINK        7
−
RESERVOIR            
+
LAKE            
−
RES_MINWSE       55.000000
+
MINWSE       55.000000
−
RES_INITWSE     62.100000
+
INITWSE     62.100000
−
RES_MAXWSE       95.000000
+
MAXWSE       95.000000
−
RES_NUMPTS       210
+
NUMPTS      210
 +
18  3    18  2    19  3    17  4    17  3    18  4    16  5    17  5    17  2    19  4   
 +
16  4    16  6    15  6    18  5    17  1    15  7    16  3    20  4    17  6    15  5   
 +
14  8    14  7    19  5    16  2    16  7    15  4    14  6    18  6    16  1    15  8   
 +
13  8    20  5    13  9    15  3    17  7    14  5    14  9    13  7    19  6    15  2   
 +
16  8    12  9    14  4    18  7    13  6    12  10    15  9    15  1    20  6    12  8   
 +
14  3    17  8    13  10    13  5    11  10    19  7    12  7    14  10    14  2    16  9   
 +
11  9    13  4    11  11    18  8    12  6    15  10    10  11    14  1    20  7    12  11   
 +
11  8    13  3    17  9    12  5    13  11    10  10    19  8    11  7    10  12    16  10   
 +
13  2    9  12    14  11    12  4    10  9    18  9    11  12    11  6    20  8    15  11   
 +
13  1    9  11    10  8    17  10    12  3    12  12    8  13    11  5    9  13    19  9   
 +
9  10    13  12    16  11    10  7    12  2    10  13    8  12    18  10    11  4    9  9   
 +
14  12    10  6    7  14    20  9    11  13    12  1    8  11    17  11    15  12    8  14   
 +
11  3    9  8    12  13    10  5    7  13    19  10    8  10    9  14    16  12    9  7   
 +
6  15    11  2    13  13    18  11    7  12    10  4    10  14    6  16    8  9    14  13   
 +
20  10    9  6    7  15    6  14    11  1    17  12    7  11    11  14    5  16    10  3   
 +
8  8    15  13    8  15    19  11    9  5    6  13    12  14    7  10    5  17    5  15   
 +
16  13    8  7    9  15    10  2    18  12    4  17    13  14    6  12    9  4    7  9   
 +
7  16    20  11    10  15    5  14    8  6    4  18    10  1    17  13    14  14    4  16   
 +
6  11    9  3    11  15    7  8    3  18    19  12    8  16    8  5    5  13    15  14   
 +
6  17    6  10    4  15    12  15    9  2    7  7    18  13    3  17    16  14    9  16   
 +
3  19    5  12    8  4    7  17    2  19    13  15    20  12    6  9    4  14    5  18   
 +
</pre>
 +
|}
 +
 
 +
'''Additional Options'''
 +
 
 +
'''Seepage'''
 +
 
 +
Seepage, and/or groundwater exchange, is controlled by including the '''M_LAKE''' card in the project file, followed by the seepage rate (cm/hr).  If multiple lakes are present and different seepage values are desired for each lake then the '''M_LAKE''' is included in the .cif file for each lake.  Include the '''M_LAKE''' card, with the value, before the '''NUMPTS''' card.  See example below.
 +
 
 +
'''Dynamic verses Static Lakes'''
 +
 
 +
The default is for lakes to be dynamic features in the channel network and in the overland flow plan.  During simulations, channels and overland cells overtaken by the lake are taken out of the stream/overland domain and put in the lake.  As the lake level recedes, these nodes/cells are put back.  More details are provided in the document linked above.  This feature can be turned off by including the '''STATIC''' card in the .cif for each lake.  Include the '''STATIC''' card before the '''NUMPTS''' card (See example below).  When a lake is '''STATIC''' it means that the lake is static in relation to the grid and the stream network, not that the stage, area, and volume are static.  For a '''STATIC''' lake the footprint of the lake is assumed to be the maximum lake size, defined by the cells specified with '''NUMPTS'''.  Using '''STATIC''' lakes may speed up simulations because complexity is reduced.  However, '''STATIC''' lakes may increase infiltration and lead to underestimation of backwater effects on both the overland and in the stream network.
 +
 
 +
'''Elevation/Volume/Surface Area'''
 +
 
 +
The default is to get the stage/storage/area relationship from the grid elevations.  However, this may lead to errors if the grid size is large.  In any case, the user can specify the stage/storage/area relationship by including the '''NUMRATE''' card followed by the number of points in the rating curve.  This is followed by the values of elevation (m), volume (m3), and area (m2), one line for each '''NUMRATE'''.  These data are included AFTER the '''NUMPTS''' card and the cells that specify the lake imprint on the grid.  See example below.
 +
 
 +
If all these options were selected for previous example, the input in the .cif file might look like this:
 +
 
 +
{|
 +
|-
 +
|
 +
<pre>
 +
LINK        7
 +
LAKE         
 +
MINWSE      55.000000
 +
INITWSE      62.100000
 +
MAXWSE      95.000000
 +
M_LAKE      10.000000
 +
STATIC
 +
NUMPTS       210
 
18  3    18  2    19  3    17  4    17  3    18  4    16  5    17  5    17  2    19  4     
 
18  3    18  2    19  3    17  4    17  3    18  4    16  5    17  5    17  2    19  4     
 
16  4    16  6    15  6    18  5    17  1    15  7    16  3    20  4    17  6    15  5     
 
16  4    16  6    15  6    18  5    17  1    15  7    16  3    20  4    17  6    15  5     
Line 606: Line 668:
 
6  17    6  10    4  15    12  15    9  2    7  7    18  13    3  17    16  14    9  16     
 
6  17    6  10    4  15    12  15    9  2    7  7    18  13    3  17    16  14    9  16     
 
3  19    5  12    8  4    7  17    2  19    13  15    20  12    6  9    4  14    5  18     
 
3  19    5  12    8  4    7  17    2  19    13  15    20  12    6  9    4  14    5  18     
 +
NUMRATE 4
 +
55.0 0.0 0.0
 +
60.0 10.0 100.0
 +
80.0 1000.0 10000.0
 +
95.0 10000.0 100000.0 
 
</pre>
 
</pre>
 
|}
 
|}
Line 802: Line 869:
 
LINK        7
 
LINK        7
 
RESERVOIR           
 
RESERVOIR           
−
RES_MINWSE       55.000000
+
MINWSE       55.000000
−
RES_INITWSE     62.100000
+
INITWSE     62.100000
−
RES_MAXWSE       95.000000
+
MAXWSE       95.000000
−
RES_NUMPTS       210
+
NUMPTS       210
 
18  3    18  2    19  3    17  4    17  3    18  4    16  5    17  5    17  2    19  4     
 
18  3    18  2    19  3    17  4    17  3    18  4    16  5    17  5    17  2    19  4     
 
16  4    16  6    15  6    18  5    17  1    15  7    16  3    20  4    17  6    15  5     
 
16  4    16  6    15  6    18  5    17  1    15  7    16  3    20  4    17  6    15  5     
Line 1,096: Line 1,163:
 
Extensive regions of adverse slope are rare in natural channels.  It is recommended that regions of adverse slopes be removed from the GSSHA channel network, unless field observations or obvious geologic controls indicate that they are justified.  Smoothing of the channel thalweg profile should be done to remove stream sections of adverse slope, undulating streambed elevations, zero slope, and sharp slope transitions.   
 
Extensive regions of adverse slope are rare in natural channels.  It is recommended that regions of adverse slopes be removed from the GSSHA channel network, unless field observations or obvious geologic controls indicate that they are justified.  Smoothing of the channel thalweg profile should be done to remove stream sections of adverse slope, undulating streambed elevations, zero slope, and sharp slope transitions.   
  
−
Ogden et al (1994) developed a method of smoothing stream channel profiles based on the realization that in regions of concave topography near streams, DEMs are likely to be positively biased.  For such conditions filling depressions along stream paths is probably erroneous.  A more accurate approach is to remove higher regions to connect depressions.  WMS provides tools to smooth the channel thalweg profile both manually and with an automated version of the Ogden et al. (1994) algorithm.   
+
Ogden et al (1994) developed a method of smoothing stream channel profiles based on the realization that in regions of concave topography near streams, DEMs are likely to be positively biased.  For such conditions filling depressions along stream paths is probably erroneous.  A more accurate approach is to remove higher regions to connect depressions.  WMS provides tools to smooth the channel thalweg profile both manually and with an automated version of the Ogden et al. (1994) algorithm.   
  
−
Changes to the stream thalweg elevation alone will not result in proper simulations of channel-groundwater exchange or channel-overbank exchange.  To properly simulate these processes, the correct elevations, or at least the correct difference in elevation, of the channel thalweg and overbank area are needed.
+
Changes to the stream thalweg elevation alone will not result in proper simulations of channel-groundwater exchange or channel-overbank exchange.  To properly simulate these processes, the correct elevations, or at least the correct difference in elevation, of the channel thalweg and overbank area are needed.  Upon initialization GSSHA will check in channel thalweg profile for problems and report them, if any, in the "check_thalweg" file.
  
 
==5.1.6 Losing and Gaining Streams==
 
==5.1.6 Losing and Gaining Streams==

Latest revision as of 18:55, 26 April 2019

5.1.1 Explicit Channel Routing Formulation

The 1-D channel routing scheme is depicted in Figure 5. Inter-cell flows, Qi-1/2 and Qi+1/2 (m3/s) in the longitudinal, x, direction are computed from depths, d, at the n time level using the Manning equation for the head discharge relationship:

Equation001.gif (1)

where: n is roughness coefficient, A is the area (m2), R is the hydraulic radius, and Sf is the friction slope, calculated in the x direction as:

Equation002.gif (2)

where: Sox is the land surface slope in the x direction. If negative flow occurs (flow in the upstream direction), the head in the downstream cell is used to calculate the flow as:

Equation003.gif (3)

Since the flow direction may change at any point in the stream, especially in ephemeral streams near the beginning of rainfall events, the flow direction is determined around each node and the locally upstream cell properties are used to compute the flow. This simple local determination of the upstream cells prevents crashes in channels with adverse slopes when little or no water is present in the upstream cell. This method also allows better simulations of backwater effects.

Inter-node fluxes are used to calculate the volume, V, in each node as:

Equation004.gif (4)

where: qlat (m2/s) is the amount of lateral inflow from the overland flow cells adjacent to the node, and qrecharge (m2/s) is the exchange between the groundwater and channel. These new volumes are used to compute nodal values of A, d, and wetted perimeter at the n+1 time level. Calculations proceed from the upstream boundary to the downstream boundary.

Several modifications were made in the implementation of the channel routing scheme to accommodate groundwater/channel interactions. These modifications permit continuous interaction between channel nodes and the saturated groundwater cells. The channel routing scheme was modified to allow water to remain in the channel after channel routing ends, and for water to be present in the channel when channel routing begins. Because groundwater may discharge to the stream at anytime, channel routing is initiated anytime a minimum amount of water is in the channel network. If the channel routing scheme indicates there is no flow in the channel, channel routing is halted during periods outside precipitation events. Fluxes between the stream and the groundwater are still computed and adjustments to the stream volumes are made without routing. If groundwater discharges to the stream, channel routing will resume, but at the groundwater time step, which is typically larger than the channel routing time step.

Because GSSHA uses a finite volume representation of channel flow the standard stability criterion, Courant number < 1.0, does not strictly apply. Maintaining stability is dependent on volume changes during each time step. Experience with the scheme indicates that stability can be maintained with a time step limitation that keeps the maximum Courant number everywhere in the network less than 1/6. The user can specify this value using the MAX_COURANT_NUM card, which requires a real number value. The default is 0.04. Groundwater and over bank fluxes can induce instability and additional controls in the channel routing scheme are added to further reduce instability. If the channel routing scheme becomes unstable (negative depth occurs in one or more cells), despite the more restrictive control on the Courant number, the time step is reduced and the channel routing calculations are repeated. The channel routing time step may be automatically reduced to a value as low as 1/1000 of a second. Allowing the time step to become very small during periods of sharp transition allows a larger overall model time step to be used. For each call of the channel routing function the overall model time step is used as the beginning channel routing time step. The time step is only reduced when the stability controls are activated, and then only for that call of the channel routing routine. The current GSSHA version employs a predictor corrector scheme with Picard iterations on the flow area. This has made the scheme much more stable when reservoirs or overbank flooding are simulated.

Several new channel routing features have been added to GSSHA since the publication of the original manual. These features are described in the following SWWRP Tech Note.

Non-ortho Tech Note

Figure 5 – Explicit channel routing scheme

In GSSHA version 7.1.4 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 CHANNEL_MOMENTUM card in your project file. This option may be helpful for simulations where channel 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. 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 CHANNEL_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.

5.1.2 Boundary Conditions

The upstream boundary condition in each first order link is a no flow condition. The default watershed outlet condition is normal flow, calculated using the channel slope at the watershed outlet. The downstream boundary condition can also be a specified head. The outlet boundary type is changed to a head by using the HEAD_BOUND project card. The boundary depth, m, is specified with the BOUND_DEPTH project card. When the head boundary is specified the depth at the outlet remains at the specified depth for the entire simulation period. Flows entering the outlet cell exit at the same rate. The head boundary condition is desirable when the condition at the outlet of the basin is a known head instead of normal depth. This might occur when the basin empties into a larger water body such as a river, pond, or lake, or when a hydraulic structure is near the watershed outlet. The up-gradient method of computing inter-node discharges allows the head boundary condition influence to propagate upstream.

5.1.3 Initial Conditions

The default for the explicit scheme is to start simulations from the dry bed condition. There is also the ability to save water surface profiles and flows from one simulation and use them as the initial condition of another simulation. The last set of depths and flows from the channel routing code can be saved in one file using the EXPLIC_BACKWATER project card to specify the file name to write depth (m) and discharge (m3 s-1) for every node in the stream network. Alternately, the values of the water surface and the discharge can be saved in separate files by using the WAT_SURF_PROFILE_OUT and DIS_PROFILE_OUT project cards to provide the file names for output of the surface water profile (m), and the discharge profile, m3 s-1), respectively. To begin a simulation from the saved values, the EXPLIC_HOTSTART project card is used to prescribe the name of the file containing the initial values of depth and discharge, or the WAT_SURF_PROFILE and DIS_PROFILE project cards are used to prescribe the files of initial water surface profile and discharge, respectively. The options can be used togather; both initial conditions can be read in and ending conditions can be written out. This feature is useful when breaking a simulation into multiple time periods. Initial conditions are also useful for starting groundwater simulations, as well as many other applications.

5.1.4 Describing the Stream Network

Channel routing is specified in the Project file by including the CHAN_EXPLIC or DIFFUSIVE_WAVE cards. Along with one of these two cards, the CHAN_INPUT and STREAM_CELL cards are used to specify the stream input files as described below.

The stream network in GSSHA is described with a series of links and nodes. A node is a single computational element in the stream network. A link is a channel segment comprised of two or more computational nodes. Two files are required to describe the channel network in GSSHA:

  • Channel input file (CHAN_INPUT) – ASCII text file describing link connectivity and assigning attributes to every channel node.
  • Stream cell input file (STREAM_CELL) - ASCII text file describing the channel and grid connectivity.

The steam cell file describes the topology of the stream network and contain the information needed to provide channel connectivity with the overland flow plane and the saturated groundwater grid. The CHAN_INPUT file contains the information necessary to connect the stream network and assign attributes to each of the stream nodes. WMS can create the STREAM_CELL and CHAN_INPUT file in the required format. Because the stream network is not confined by the grid, the STREAM_CELL input file is particulary complicated and difficult to build without the WMS interface. The rules that govern the creation of channel input file are listed below.

General Rules

  1. Looped reaches cannot be simulated.
  2. Link types may not be mixed within a link.
  3. Internal boundary condition links contain two nodes, one upstream, and one downstream, of the internal boundary condition.

Required Numbering Scheme

  1. First order stream links are numbered first.
  2. The first link must be numbered 1.
  3. Link numbers must always increase in the downstream direction. Upstream links must have smaller numbers than downstream links.
  4. Link numbers must change at junctions.
  5. Link numbers may not be skipped.
  6. Node numbers increase in the downstream direction.
  7. The most upstream node in a link is node 1.

Upstream and downstream links overlap one node. The first node in the downstream link is the same as the last node in the upstream link. In the CHAN_INPUT file this node must be explicitly contained in both the upstream and downstream links, such that the upstream link in the CHAN_INPUT will contain one more node computational node than actually exist. An exception to this rule is the most downstream link containing the outlet.

As an example, we'll consider the highly conceptualized stream network shown in the figure below.

Figure 6 – Conceptualized Stream Network

This stream network has 9 stream links, including a lake, link 7, with an outlet weir, link 8.

5.1.4.1 – Channel Input File

The channel input file name is specified with the CHAN_INPUT card in the project file. The channel input file consists of four parts:
  • channel input file declaration
  • channel routing constants
  • channel network connectivity
  • link (reach) information

These four sections are assembled in the listed order to produce the channel input file.

5.1.4.1.1 Channel input file declaration

The first line of the channel input file must contain the card "GSSHA_CHAN". This card tells GSSHA that it is reading a channel input file. In the CHAN_INPUT file the first line will look like: 

GSSHA_CHAN

5.1.4.1.2 Channel routing constants

The second portion of the channel input file is used to pass physical constants and simulation parameters to the model. There are three channel routing parameters that must be prescribed. The program must be told the number of links, and the maximum number of nodes in any link in the network. In our example, the number of links is 9, and the maximum number of nodes is 9 (in link 6). Recall that all links except the outlet link must have an extra node for connectivity purposes. The total number of links is called "LINKS", and the largest number of nodes in any link in the network is called "MAXNODES". The three paramets, "ALPHA" "BETA" 'THETA" are not currently used. A value of 1.0 is suggested for these parameters.

Parameter Value Units Description
ALPHA 1.000000 none Channel routing parameter, not currently used.
BETA 1.000000 none Channel routing parameter, not currently used.
THETA 1.000000 none Channel routing parameter, not currently used.
LINKS 9 none Number of links in the channel network.
MAXNODES 9 none Number of nodes in the longest link in the network.

Table 5 - Channel file input constants

This data constitutes the second portion of channel input file and is arranged into a header which must have the form (note floating point and integers): 

ALPHA       1.000000
BETA        1.000000
THETA       1.000000
LINKS       9
MAXNODES    9

5.1.4.1.3 – Channel network connectivity

The third portion of the input file describes the network topology. This is accomplished using a line-input format, one line for each link in the network. Each line contains the identifier "CONNECT" followed by three or more integer values arranged in columns, as shown in Table 6. The total number of columns depends on the number of upstream dependencies, with one for each. The links should appear in ascending numerical order.

Column Data Description
1 CONNECT
2 Link number
3 Downstream link number (0 for outlet link)
4 Number of links upstream from this link
5-N Upstream link #1 through number of upstream links

Table 6 - Connectivity information format

For our example, the stream connectivity information looks like. 

CONNECT    1    4    0
CONNECT    2    4    0
CONNECT    3    4    0
CONNECT    4    5    3    1    2    3
CONNECT    5    6    1    4
CONNECT    6    7    1    5
CONNECT    7    8    1    6
CONNECT    8    9    1    7
CONNECT    9    0    1    8

Note that only the outlet link has 0 downstream dependencies. Links 1,2, and 3 have no upstream dependencies, while link 4 has 3.

5.1.4.1.4 - Link (Reach) information

Following the channel connectivity information is channel link and node information. Parameters are defined for each link in the stream network. Links are listed in ascending order.

5.1.4.1.4.1 - Link declaration

In the CHAN_INPUT file each link is declared with the "LINK" card, which is followed by the link number. For any CHAN_INPUT file the first line after the channel connect information will be: 

LINK           1
5.1.4.1.4.2 - Link types

The input that follows the link declaration depends on the link type. Link types may be basic fluvial links, structures, or reservoirs.

5.1.4.1.4.2.1 - Basic fluvial channel links

Channel cross-sections may be either trapezoidal, or natural channels represented by breakpoint cross-sections. Different cross-section types may be mixed in the channel network. Typically, small streams in the upland areas are modeled with trapezoidal cross-sections; larger streams in the lower reaches of the catchment are modeled with break-point cross-sections. Data for the development of break-point cross-sections can be derived from field surveys of cross-section, detailed topographic maps, or high resolution digital elevation models. Coarse resolution DEMs should not be relied upon to develop channel cross-section geometry.

The hydrodynamic channel routing in GSSHA is very sensitive to channel cross-sections, more so than the longitudinal slope. More time should be spent developing accurate cross-sectional description than to getting a highly accurate longitudinal profile. Smooth transitions in channel cross-sectional properties between all connecting fluvial links often play a vital role in the success of simulations. Abrupt changes in cross-sections can lead to numerical mass conservation errors. It may be necessary to create transition links between measured break-point and trapezoidal cross-sections when adjoining links vary greatly in cross-section. The Manning’s roughness parameter used in the model is largely a calibration parameter. Physically realistic literature values may be used as initial guesses and as bound on the calibration (e.g. Chow, 1959; Barnes, 1967).

Basic fluvial links contain nodes. For a basic fluvial link the next line after the link declaration card is the channel node increment, or channel grid size. This is declared by using the "DX" card, followed by the node spacing in the current reach. The node spacing can be different in every link. For our example, link 1 has a node spacing of 107.840396 m. The input in the CHAN_INPUT file is: 

DX             107.840396

After setting the node spacing for the link, the type of link is specified. Basic channel reaches may be defined as either "TRAPEZOID" or "BREAKPOINT" for trapezoid and natural cross sections, respectively. These basic reach types may be further described as "SUBSURFACE", for surface water/groundwater interaction, and as "ERODE", for erodible channel if performing sediment transport calculations. Both "SUBSURFACE" and "ERODE" may be combined with with either "TRAPEZOID" or "BREAKPOINT" with the use of underscores. They may used singly or combined. Such as 

TRAPEZOID_ERODE
TRAPEZOID_SUBSURFACE

or, 

TRAPEZOID_SUBSURFACE_ERODE

The basic channel type TRAPEZOID or BREAKPOINT should always be listed first, followed by the TRAPEZOID and/or SUBSURFACE descriptors.

For our example, Link 1 is a trapezoidal channel. So the input is: 

TRAPEZOID

After the link type is defined, the number of nodes is specified with the "NODE" card. For our example, Link 1 has 7 nodes, so the input is: 

NODES 7

After defining the number of nodes, information for the first node, node 1, in the link is specified. Required node information is the x and y location of the node and the elevation of the node. These are specified with "X_Y" and "ELEV" cards. For our example the information for link 1 is: 

NODE 1
X_Y  1825.200000 933.400000
ELEV 112.694141

After describing the first node, the channel cross section must be defined. How the cross section is defined depends on the type of cross section. In either case, the "X_SEC" card is included next in the CHAN_INPUT file. 

X_SEC
5.1.4.1.4.2.1.1 Trapezoidal cross-section

Trapezoidal cross-sections are symmetrical, assumed to be infinitely deep, and have a constant side slope. Trapezoidal cross-sections are defined by:

  • Manning’s N - Manning’s roughness coefficient, n (dimensionless)
  • bottom width (m)
  • channel depth (m)
  • side slope, z (change in X with a change in Y of 1).

This information is entered for each trapezoidal link with a series of cards. These cards are:

MANNINGS_N
BOTTOM_WIDTH
BANKFULL_DEPTH
SIDE_SLOPE

Each card is followed by the value. For our example Link 1 is a trapezoidal channel with Manning roughness of 0.03, a bottom width of 2 m, a bankfull depth of 2 m, and a side slope of 2.0. The X_SEC input for this reach is: 

XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
5.1.4.1.4.2.1.2 Natural cross-section

"BREAKPOINT" cross sections are used to describe natural cross sections. Breakpoint cross sections are defined by the points along the channel cross section. These are entered for each link after the "X_SEC" declarations. Internal to GSSHA, these points and the specified Manning roughness will be used to develop a table of channel information for a specified number of increments between the channel bottom and the bank full depth.

The information in the CHAN_INPUT file is the "MANNINGS_N", the number of XY pairs "NPAIRS" that describe the channel cross section, the number of increments in the channel property table used in GSSHA "NUM_INTERP", followed by the series of XY pairs, specified with the "X1" card.

For our example link 9 is a "BREAKPOINT" cross section, with a roughness of 0.03, and 4 xy points that define the cross section. We are specifying that the cross section be divided into 2000 depth increments for computation purposes. The "X_SEC" input for link 9 appears as: 

XSEC
MANNINGS_N     0.030000
NPAIRS         4
NUM_INTERP     2000
X1   -4.000000 2.000000
X1   0.000000 0.000000
X1   2.000000 0.000000
X1   6.000000 2.000000
5.1.4.1.4.2.1.3 Additional Channel Properties (X_SEC qualifiers)

Trapezoidal and breakpoint cross sections represent the two basic fluvial cross section types in GSSHA. For more complex models additional channel properties can be defined after the X_SEC card. The table below describes the additional X_SEC cards.

Card Argument Description
ERODE none specifies that the link is an erodible cross section
MAX_EROSION real maximum erosion (m) of the bed in erodible channels
SUBSURFACE none specifies that the link has groundwater interaction
M_RIVER real depth of the bed layer sediments (cm) for groundwater interaction
K_RIVER real hydraulic conductivity of the bed materials (cm/hr) for groundwater interaction
5.1.4.1.4.2.1.4 Node information

Each fluvial channel link must contain two or more nodes. The nodes are the computational elements used in GSSHA to perform the actual channel routing calculations. Nodes are defined by there XY location, as well as the thalweg elevation, the elevation of the lowest point in the channel cross section.

After describing the channel cross section, XYZ information must be entered for each computational node in the link. For each node the "NODE" card is used to define the node; the x and y position are entered with the "X_Y" card; followed by the node elevation is set with the "ELEV" card. This information is entered for each node in ascending order.

In our example, link 1 has 7 nodes. Information for node 1 was previously described, information for nodes 2 through 7 are listed as follows. 

NODE 2
X_Y  1734.200000 991.266667
ELEV 105.977735
NODE 3
X_Y  1643.200000 1049.133333
ELEV 99.186421
NODE 4
X_Y  1552.200000 1107.000000
ELEV 92.323354
NODE 5
X_Y  1461.200000 1164.866667
ELEV 85.407261
NODE 6
X_Y  1370.200000 1222.733333
ELEV 78.464655
NODE 7
X_Y  1279.200000 1280.600000
ELEV 73.846634

It is required that the bed elevation of all channel inverts at each junction be equal. In our example the elevation of the downstream node of links 1, 2 and 3 are equal to the bed elevation of the upstream end of link 4 (Figure 5).

5.1.4.1.4.2 - Structure channel links

Structures are internal boundary conditions within the GSSHA stream routing network. GSSHA allows several types of structures to be defined. The 7 types of structures recognized in GSSHA are:

  • Horizontal Weir
  • Sag Vertical Curve (Parabolic) Weir
  • Circular Culvert
  • Box Culvert
  • Rating Curve
  • Scheduled Release
  • Rule Curve

All culverts, for example, lie below a weir formed by a roadway or other embankment. In GSSHA, a hydraulic structure link can consist of any combination of the above. Other structures consist of a number of culverts, some with different dimensions, invert elevations, etc. This allows the user to simulate complex hydraulics from multiple drainage elements.

Hydraulic structures are point models that exist at a channel intersection, at the last node for the upstream link and at the first node of the downstream link, which occupy the same point. The upstream conditions are taken from the second to last node in the upstream link. The downstream conditions are taken from the second node in the downstream link. It is important to note that the bed elevation of the upstream condition is the bed elevation at the junction plus half the difference between the bed elevation at the junction and the bed elevation at the next to last node of the upstream link. The same applies to the downstream condition. To avoid numerical issues, it is advisable to the structure invert be at or above the bed elevation at the respective end.

Horizontal Weir

The "weir" card is used to specify a horizontal broad- or sharp-crested weir of specified length. The figure below shows the geometry of the weir looking downstream. The CREST_LENGTH card specifies the length L shown in the figure in meters.

Horiz weir.gif

Other required inputs are the elevation of the crest above the model elevation datum CREST_LOW_ELEV, and the forward- and reverse-flow discharge coefficients, DISCHARGE_COEFF_FORWARD and DISCHARGE_COEFF_REVERSE. 

STRUCTTYPE     WEIR
CREST_LENGTH             6.200000
CREST_LOW_ELEV           61.500000
DISCHARGE_COEFF_FORWARD  0.600000
DISCHARGE_COEFF_REVERSE  0.600000

Sag Vertical-Curve Weir

The "sag_weir" card is used to specify a non-symmetric parabolic weir of the type that is commonly created by flow over a roadway crossing a stream, when the flood plain is small. Parabolic vertical curvature of roadways is very common. The geometrical parameters shown in the Figure below are based on the slopes of the two road segments that are to be joined by the parabola, the horizontal distance between the point of vertical curvature (PVC) and point of vertical tangency (PVT), and elevation of the low point on the curve.

Sag curve weir.gif

<math>Q=C_d L H^{3/2}</math>

The required inputs for sag vertical curve weirs are: CREST_LENGTH L, which is the horizontal distance from PVC to PVT, CREST_LOW_ELEV or E_inv in the figure, the forward and reverse discharge coefficients DISCHARGE_COEFF_FORWARD and DISCHARGE_COEFF_REVERSE, the steepest approach slope entered as a negative number STEEP_SLOPE, and the shallowest approach slope entered as a positive number, SHALLOW_SLOPE. 

STRUCTTYPE     SAG_WEIR
CREST_LENGTH             25.000000
CREST_LOW_ELEV           370.000000
CREST_LOW_LOC            0
DISCHARGE_COEFF_FORWARD  0.700000
DISCHARGE_COEFF_REVERSE  0.700000
STEEP_SLOPE              -0.10
SHALLOW_SLOPE            0.08

Circular Culvert

The "round_culvert" card specifies a culvert of circular cross-section. The following six flow types are considered in GSSHA for both circular and rectangular box culverts:

Culv flow types.gif

Required inputs for circular culverts are the upstream invert elevation UPINVERT (m), the downstream invert elevation DOWNINVERT (m), the inlet discharge coefficient INLET_DISCH_COEFF, the reverse flow discharge coefficient REV_FLOW_DISCH_COEFF, culvert length LENGTH (m),the diameter of the culvert DIAMETER (m), and the Manning roughness coefficient ROUGH_COEFF. The SLOPE card is not a require input, and is calculated by the code. However, if the length of the culvert is not known, the UPINVERT and SLOPE are used to calculate the DOWNINVERT elevation. 

STRUCTTYPE     ROUND_CULVERT
UPINVERT                 409.600000
DOWNINVERT               409.500000
INLET_DISCH_COEFF        0.800000
REV_FLOW_DISCH_COEFF     0.900000
SLOPE                    0.000000
LENGTH                   10.000000
ROUGH_COEFF              0.013000
DIAMETER                 1.220000

As noted above, numerical issues can be avoiding by locating you inverts above the bed elevations at the upstream and downstream condition locations, which are at the mid-point between the junction (structure) node and the next upstream or downstream node, depending on whether it is the UPINVERT or the DOWNINVERT.

Box Culvert

The "rect_culvert" card is used to create a rectangular cross-section culvert. The same flow regimes shown in the figure above are considered.

Required inputs for rectangular culverts are the upstream invert elevation UPINVERT (m), the downstream invert elevation DOWNINVERT (m), the inlet discharge coefficient INLET_DISCH_COEFF, the reverse flow discharge coefficient REV_FLOW_DISCH_COEFF, culvert length LENGTH (m),the width of the rectangular culvert WIDTH (m), the height of the rectangular culvert HEIGHT (m), and the Manning roughness coefficient ROUGH_COEFF. The SLOPE card is not a require input, and is calculated by the code. However, if the length of the culvert is not known, the UPINVERT and SLOPE are used to calculate the DOWNINVERT elevation. If the UPINVERT and DOWNINVERT cards are specified, the slope of the culvert is calculated, and any value input through the SLOPE card is not used. 

STRUCTTYPE     RECT_CULVERT
UPINVERT                 413.700000
DOWNINVERT               413.600000
INLET_DISCH_COEFF        0.800000
REV_FLOW_DISCH_COEFF     0.900000
SLOPE                    0.010000
LENGTH                   10.000000
ROUGH_COEFF              0.013000
WIDTH                    0.910000
HEIGHT                   0.670000

As noted above, numerical issues can be avoiding by locating you inverts above the bed elevations at the upstream and downstream condition locations, which are at the mid-point between the junction (structure) node and the next upstream or downstream node, depending on whether it is the UPINVERT or the DOWNINVERT.

Rating Curve

The "rating_curve" card is used to specify a series of discharge vs. water surface elevation points, and can be used to simulate the hydraulic performance of "unusual" hydraulic structures such as irregular weirs or roadway crossings.


Scheduled Release

The "scheduled_release" card is used to enter specified discharges as a function of time, such as those from a reservoir. This functionality is particularly valuable for simulations of past events, where reservoir releases are known.


Rule Curve

The "rule_curve" card is used to enter discrete discharges as a function of water level. These rule curves are often associated with the manual operation of flood control reservoirs. At different times of the year, flood control reservoirs are operated to try to maintain a certain "target" water level in the reservoir. If the water level exceeds the target by a certain amount, the reservoir operator is authorized by the rule curve to release a certain flow rate of water from the reservoir in an attempt to restore the reservoir level to the desired target level.

Composite Structure Link

Suppose that at a certain point a channel is crossed by a parabolic-shaped roadway, and three culverts, one square and two circular are used to pass the normal flow under the roadway. The situation is shown in the following figure.

Compound structure.gif

This hydraulic structure consists of four elements: sag vertical-curve weir, box culvert, and two circular culverts. Notice that as shown in the figure, each element has a different invert elevation, which is allowed.

For the sake of discussion, assume that within the GSSHA channel network, this hydraulic structure is link number 38. The channel input file for this link would appear as: 

LINK         38
STRUCTURE
NUMSTRUCTS     4
STRUCTTYPE     RECT_CULVERT
UPINVERT                 2369.790000
DOWNINVERT               2369.540000
INLET_DISCH_COEFF        0.800000
REV_FLOW_DISCH_COEFF     0.900000
SLOPE                    0.010000
LENGTH                   18.000000
ROUGH_COEFF              0.013000
WIDTH                    2.40000
HEIGHT                   1.270000
STRUCTTYPE     ROUND_CULVERT
UPINVERT                 2370.830000
DOWNINVERT               2370.080000
INLET_DISCH_COEFF        0.800000
REV_FLOW_DISCH_COEFF     0.900000
SLOPE                    0.000000
LENGTH                   19.000000
ROUGH_COEFF              0.013000
DIAMETER                 1.220000
STRUCTTYPE     ROUND_CULVERT
UPINVERT                 2370.410000
DOWNINVERT               2369.870000
INLET_DISCH_COEFF        0.800000
REV_FLOW_DISCH_COEFF     0.900000
SLOPE                    0.000000
LENGTH                   19.000000
ROUGH_COEFF              0.013000
DIAMETER                 1.220000
STRUCTTYPE     SAG_WEIR
CREST_LENGTH             45.000000
CREST_LOW_ELEV           2372.860000
CREST_LOW_LOC            0
DISCHARGE_COEFF_FORWARD  0.700000
DISCHARGE_COEFF_REVERSE  0.700000
STEEP_SLOPE              -0.05
SHALLOW_SLOPE            0.035

Note that the order in which the hydraulic elements appear in the STRUCTURE block when there is more than one structural element is not important.

5.1.4.1.4.3 - Reservoir channel links

Reservoirs are defined by specifying that a reservoir is present with the LAKE card in the .cif file. The minimum, initial, and maximum reservoir elevations are specified with the MINWSE, INITWSE, and MAXWSE, cards respectively. The number of cells within the maximum reservoir boundary is specified with the NUMPTS card, followed by a list of the cells in the reservoir domain, defined by row and column values. Reservoirs are dynamic features that can increase/decrease in size, take/release stream nodes, and overland flow cells as defined in Non-ortho Channels Technical Note Non-ortho Tech Note. The resevoir will not go into any cells not listed in the NUMPTS list, but water may flow out of the reservoir back on the overland flow plane. All reservoirs must have a hydraulic structure defined as the reservoir outlet. Other important restrictions are that the reservoir cannot rise above the MAXWSE, move into cells not defined as part of the reservoir, join with other reservoirs, or take in stream nodes not directly in the same branch as the reservoir. Any of these conditions will cause the program to stop, crash, or provide non-sensical results. Large mass balance errors are indicative of problems with reservoirs violating one or more of these conditions. Reservoirs can be kept within the desired confines using embankments. Water can be passed through an embankment holding a reservoir using lowspots along the embankment or using overland rating curves. Water can also pass from one lake to another using these features. Although lakes cannot exceed the maximum water level, they can go dry. If the water level falls below the MINWSE it will stop discharging but water can still seep out of the lake and/or evaporate.

Output from all the reservoirs can be specified with the LAKE_OUTPUT card specified in the project file. Lake elevations and volumes, respectively, will be written for each lake, in the order they are listed in the .cif file, for each time period specified in the HYD_FREQ card, in the file specified with the LAKE_OUTPUT card.

Within the .cif file, a standard lake input would look like the example below: 


LINK         7
LAKE           
MINWSE       55.000000
INITWSE      62.100000
MAXWSE       95.000000
NUMPTS       210
18  3     18  2     19  3     17  4     17  3     18  4     16  5     17  5     17  2     19  4     
16  4     16  6     15  6     18  5     17  1     15  7     16  3     20  4     17  6     15  5     
14  8     14  7     19  5     16  2     16  7     15  4     14  6     18  6     16  1     15  8     
13  8     20  5     13  9     15  3     17  7     14  5     14  9     13  7     19  6     15  2     
16  8     12  9     14  4     18  7     13  6     12  10     15  9     15  1     20  6     12  8     
14  3     17  8     13  10     13  5     11  10     19  7     12  7     14  10     14  2     16  9     
11  9     13  4     11  11     18  8     12  6     15  10     10  11     14  1     20  7     12  11     
11  8     13  3     17  9     12  5     13  11     10  10     19  8     11  7     10  12     16  10     
13  2     9  12     14  11     12  4     10  9     18  9     11  12     11  6     20  8     15  11     
13  1     9  11     10  8     17  10     12  3     12  12     8  13     11  5     9  13     19  9     
9  10     13  12     16  11     10  7     12  2     10  13     8  12     18  10     11  4     9  9     
14  12     10  6     7  14     20  9     11  13     12  1     8  11     17  11     15  12     8  14     
11  3     9  8     12  13     10  5     7  13     19  10     8  10     9  14     16  12     9  7     
6  15     11  2     13  13     18  11     7  12     10  4     10  14     6  16     8  9     14  13     
20  10     9  6     7  15     6  14     11  1     17  12     7  11     11  14     5  16     10  3     
8  8     15  13     8  15     19  11     9  5     6  13     12  14     7  10     5  17     5  15     
16  13     8  7     9  15     10  2     18  12     4  17     13  14     6  12     9  4     7  9     
7  16     20  11     10  15     5  14     8  6     4  18     10  1     17  13     14  14     4  16     
6  11     9  3     11  15     7  8     3  18     19  12     8  16     8  5     5  13     15  14     
6  17     6  10     4  15     12  15     9  2     7  7     18  13     3  17     16  14     9  16     
3  19     5  12     8  4     7  17     2  19     13  15     20  12     6  9     4  14     5  18

Additional Options

Seepage

Seepage, and/or groundwater exchange, is controlled by including the M_LAKE card in the project file, followed by the seepage rate (cm/hr). If multiple lakes are present and different seepage values are desired for each lake then the M_LAKE is included in the .cif file for each lake. Include the M_LAKE card, with the value, before the NUMPTS card. See example below.

Dynamic verses Static Lakes

The default is for lakes to be dynamic features in the channel network and in the overland flow plan. During simulations, channels and overland cells overtaken by the lake are taken out of the stream/overland domain and put in the lake. As the lake level recedes, these nodes/cells are put back. More details are provided in the document linked above. This feature can be turned off by including the STATIC card in the .cif for each lake. Include the STATIC card before the NUMPTS card (See example below). When a lake is STATIC it means that the lake is static in relation to the grid and the stream network, not that the stage, area, and volume are static. For a STATIC lake the footprint of the lake is assumed to be the maximum lake size, defined by the cells specified with NUMPTS. Using STATIC lakes may speed up simulations because complexity is reduced. However, STATIC lakes may increase infiltration and lead to underestimation of backwater effects on both the overland and in the stream network.

Elevation/Volume/Surface Area

The default is to get the stage/storage/area relationship from the grid elevations. However, this may lead to errors if the grid size is large. In any case, the user can specify the stage/storage/area relationship by including the NUMRATE card followed by the number of points in the rating curve. This is followed by the values of elevation (m), volume (m3), and area (m2), one line for each NUMRATE. These data are included AFTER the NUMPTS card and the cells that specify the lake imprint on the grid. See example below.

If all these options were selected for previous example, the input in the .cif file might look like this: 

LINK         7
LAKE           
MINWSE       55.000000
INITWSE      62.100000
MAXWSE       95.000000
M_LAKE       10.000000
STATIC
NUMPTS       210
18  3     18  2     19  3     17  4     17  3     18  4     16  5     17  5     17  2     19  4     
16  4     16  6     15  6     18  5     17  1     15  7     16  3     20  4     17  6     15  5     
14  8     14  7     19  5     16  2     16  7     15  4     14  6     18  6     16  1     15  8     
13  8     20  5     13  9     15  3     17  7     14  5     14  9     13  7     19  6     15  2     
16  8     12  9     14  4     18  7     13  6     12  10     15  9     15  1     20  6     12  8     
14  3     17  8     13  10     13  5     11  10     19  7     12  7     14  10     14  2     16  9     
11  9     13  4     11  11     18  8     12  6     15  10     10  11     14  1     20  7     12  11     
11  8     13  3     17  9     12  5     13  11     10  10     19  8     11  7     10  12     16  10     
13  2     9  12     14  11     12  4     10  9     18  9     11  12     11  6     20  8     15  11     
13  1     9  11     10  8     17  10     12  3     12  12     8  13     11  5     9  13     19  9     
9  10     13  12     16  11     10  7     12  2     10  13     8  12     18  10     11  4     9  9     
14  12     10  6     7  14     20  9     11  13     12  1     8  11     17  11     15  12     8  14     
11  3     9  8     12  13     10  5     7  13     19  10     8  10     9  14     16  12     9  7     
6  15     11  2     13  13     18  11     7  12     10  4     10  14     6  16     8  9     14  13     
20  10     9  6     7  15     6  14     11  1     17  12     7  11     11  14     5  16     10  3     
8  8     15  13     8  15     19  11     9  5     6  13     12  14     7  10     5  17     5  15     
16  13     8  7     9  15     10  2     18  12     4  17     13  14     6  12     9  4     7  9     
7  16     20  11     10  15     5  14     8  6     4  18     10  1     17  13     14  14     4  16     
6  11     9  3     11  15     7  8     3  18     19  12     8  16     8  5     5  13     15  14     
6  17     6  10     4  15     12  15     9  2     7  7     18  13     3  17     16  14     9  16     
3  19     5  12     8  4     7  17     2  19     13  15     20  12     6  9     4  14     5  18     
NUMRATE 4
55.0 0.0 0.0
60.0 10.0 100.0
80.0 1000.0 10000.0
95.0 10000.0 100000.0

5.1.4.1.5 Assembling the Channel Input File

In this section, the three portions of the channel input file are put together to form the channel input file for our example project.

Complete Example File 
GSSHA_CHAN
ALPHA       1.000000
BETA        1.000000
THETA       1.000000
LINKS       9
MAXNODES    9
CONNECT    1    4    0
CONNECT    2    4    0
CONNECT    3    4    0
CONNECT    4    5    3    1    2    3
CONNECT    5    6    1    4
CONNECT    6    7    1    5
CONNECT    7    8    1    6
CONNECT    8    9    1    7
CONNECT    9    0    1    8

LINK           1
DX             107.840396
TRAPEZOID
NODES          7
NODE 1
X_Y  1825.200000 933.400000
ELEV 112.694141
XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
NODE 2
X_Y  1734.200000 991.266667
ELEV 105.977735
NODE 3
X_Y  1643.200000 1049.133333
ELEV 99.186421
NODE 4
X_Y  1552.200000 1107.000000
ELEV 92.323354
NODE 5
X_Y  1461.200000 1164.866667
ELEV 85.407261
NODE 6
X_Y  1370.200000 1222.733333
ELEV 78.464655
NODE 7
X_Y  1279.200000 1280.600000
ELEV 73.846634

LINK           2
DX             97.619562
TRAPEZOID
NODES          7
NODE 1
X_Y  920.200000 1743.400000
ELEV 109.559260
XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
NODE 2
X_Y  980.033333 1666.266667
ELEV 102.920801
NODE 3
X_Y  1039.866667 1589.133333
ELEV 96.282341
NODE 4
X_Y  1099.700000 1512.000000
ELEV 90.355169
NODE 5
X_Y  1159.533333 1434.866667
ELEV 84.165096
NODE 6
X_Y  1219.366667 1357.733333
ELEV 77.843573
NODE 7
X_Y  1279.200000 1280.600000
ELEV 73.846634

LINK           3
DX             119.589339
TRAPEZOID
NODES          3
NODE 1
X_Y  1445.400000 1452.600000
ELEV 80.258959
XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
NODE 2
X_Y  1362.300000 1366.600000
ELEV 75.890504
NODE 3
X_Y  1279.200000 1280.600000
ELEV 73.846634

LINK           4
DX             97.528863
TRAPEZOID
NODES          6
NODE 1
X_Y  1279.200000 1280.600000
ELEV 73.846634
XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
NODE 2
X_Y  1214.524434 1207.600341
ELEV 68.100627
NODE 3
X_Y  1149.848867 1134.600682
ELEV 64.624511
NODE 4
X_Y  1085.173301 1061.601023
ELEV 61.082366
NODE 5
X_Y  1020.497734 988.601364
ELEV 57.443587
NODE 6
X_Y  955.822168 915.601705
ELEV 53.724277

LINK           5
DX             102.425766
TRAPEZOID
NODES          3
NODE 1
X_Y  955.822168 915.601705
ELEV 53.724277
XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
NODE 2
X_Y  881.349951 845.281535
ELEV 49.777956
NODE 3
X_Y  806.877734 774.961364
ELEV 45.831635

LINK           6
DX             102.425766
TRAPEZOID
NODES          9
NODE 1
X_Y  806.877734 774.961364
ELEV 45.831635
XSEC
MANNINGS_N     0.030000
BOTTOM_WIDTH   2.000000
BANKFULL_DEPTH 2.000000
SIDE_SLOPE     2.000000
NODE 2
X_Y  732.405518 704.641194
ELEV 42.106178
NODE 3
X_Y  657.933301 634.321023
ELEV 38.368474
NODE 4
X_Y  583.461084 564.000853
ELEV 34.624646
NODE 5
X_Y  508.988867 493.680682
ELEV 30.880817
NODE 6
X_Y  434.516650 423.360512
ELEV 30.234555
NODE 7
X_Y  360.044434 353.040341
ELEV 27.028553
NODE 8
X_Y  285.572217 282.720171
ELEV 21.302121
NODE 9
X_Y  211.100000 212.400000
ELEV 14.343808

LINK         7
RESERVOIR           
MINWSE       55.000000
INITWSE      62.100000
MAXWSE       95.000000
NUMPTS       210
18  3     18  2     19  3     17  4     17  3     18  4     16  5     17  5     17  2     19  4     
16  4     16  6     15  6     18  5     17  1     15  7     16  3     20  4     17  6     15  5     
14  8     14  7     19  5     16  2     16  7     15  4     14  6     18  6     16  1     15  8     
13  8     20  5     13  9     15  3     17  7     14  5     14  9     13  7     19  6     15  2     
16  8     12  9     14  4     18  7     13  6     12  10     15  9     15  1     20  6     12  8     
14  3     17  8     13  10     13  5     11  10     19  7     12  7     14  10     14  2     16  9     
11  9     13  4     11  11     18  8     12  6     15  10     10  11     14  1     20  7     12  11     
11  8     13  3     17  9     12  5     13  11     10  10     19  8     11  7     10  12     16  10     
13  2     9  12     14  11     12  4     10  9     18  9     11  12     11  6     20  8     15  11     
13  1     9  11     10  8     17  10     12  3     12  12     8  13     11  5     9  13     19  9     
9  10     13  12     16  11     10  7     12  2     10  13     8  12     18  10     11  4     9  9     
14  12     10  6     7  14     20  9     11  13     12  1     8  11     17  11     15  12     8  14     
11  3     9  8     12  13     10  5     7  13     19  10     8  10     9  14     16  12     9  7     
6  15     11  2     13  13     18  11     7  12     10  4     10  14     6  16     8  9     14  13     
20  10     9  6     7  15     6  14     11  1     17  12     7  11     11  14     5  16     10  3     
8  8     15  13     8  15     19  11     9  5     6  13     12  14     7  10     5  17     5  15     
16  13     8  7     9  15     10  2     18  12     4  17     13  14     6  12     9  4     7  9     
7  16     20  11     10  15     5  14     8  6     4  18     10  1     17  13     14  14     4  16     
6  11     9  3     11  15     7  8     3  18     19  12     8  16     8  5     5  13     15  14     
6  17     6  10     4  15     12  15     9  2     7  7     18  13     3  17     16  14     9  16     
3  19     5  12     8  4     7  17     2  19     13  15     20  12     6  9     4  14     5  18     

LINK         8
STRUCTURE
NUMSTRUCTS     1
STRUCTTYPE     WEIR
CREST_LENGTH             0.200000
CREST_LOW_ELEV           61.500000
CREST_LOW_LOC            0
DISCHARGE_COEFF_FORWARD  0.600000
DISCHARGE_COEFF_REVERSE  0.600000

LINK           9
DX             97.920932
BREAKPOINT
NODES          4
NODE 1
X_Y  211.100000 212.400000
ELEV 14.343808
XSEC
MANNINGS_N     0.030000
NPAIRS         4
NUM_INTERP     2000
X1   -4.000000 2.000000
X1   0.000000 0.000000
X1   2.000000 0.000000
X1   6.000000 2.000000
NODE 2
X_Y  140.866667 144.166667
ELEV 7.691532
NODE 3
X_Y  70.633333 75.933333
ELEV 4.305107
NODE 4
X_Y  0.400000 7.700000
ELEV 2.116070

5.1.4.2 - Grid connectivity STREAM_CELL file

The topology relationships between the channel and the overland grid are defined in the STREAM_CELL file. This file contains two types of information, the location of the channel nodes in the grid plane in terms of I (row) and J (column) grid location, and the percentage of each stream node in the listed grid cell. Our example stream network is shown over the grid in the following figure.

Figure 7 – Stream network over grid

As is seen in the figure, every node falls into a grid cell. In addition, every node length falls in one or more grid cells. A certain percentage of each computational node falls within any grid cell that contains a channel node. This information is contained in the STREAM_CELL file. The STREAM_CELL file contains the following information.

5.1.4.2.1 - STREAM_CELL file identifier

Every STREAM_CELL file begins with the card

GRIDSTREAMFILE

5.1.4.2.2 - Number of stream cells

The second line in the STREAM_CELL file is the total number of grid cells with streams "STREAMCELLS". In our example the total number of cells with at least part of stream node is 45. So the second line in our STREAM_CELL file is:

STREAMCELLS 45

5.1.4.2.3 - Link/node/grid matching information

Every cell with some part of a stream node in it is identified with the "CELLIJ" card, which gives the row and column of the cell to be described. This is followed by the number of nodes in the cell with the "NODES" card. Followed with the "LINKNODE" card for every node in the grid cell. The "LINKNODE" card specifies the link number, node number, and the percentage of the node within the grid cell. In our example, about 27% of link 1, node 1, falls in grid cell with row 11 column 19. The remaining portion falls in grid cell 11 18, along with about 15% of link 1 node 2. This information in the STREAM_CELL file looks like: 

CELLIJ        11 19
NUMNODES      1
LINKNODE      1 1 0.276923
CELLIJ        11 18
NUMNODES      2
LINKNODE      1 1 0.723077
LINKNODE      1 2 0.150922

This information is repeated for every cell containing portions of a channel node. The information is listed is ascending link/node order. For our example, the entire STREAM_CELL file looks like:

Example 
GRIDSTREAMFILE
STREAMCELLS   45
CELLIJ        11 19
NUMNODES      1
LINKNODE      1 1 0.276923
CELLIJ        11 18
NUMNODES      2
LINKNODE      1 1 0.723077
LINKNODE      1 2 0.150922
CELLIJ        10 18
NUMNODES      1
LINKNODE      1 2 0.224903
CELLIJ        10 17
NUMNODES      2
LINKNODE      1 2 0.624176
LINKNODE      1 3 0.474725
CELLIJ        10 16
NUMNODES      1
LINKNODE      1 3 0.404307
CELLIJ        9 16
NUMNODES      2
LINKNODE      1 3 0.120968
LINKNODE      1 4 0.573626
CELLIJ        9 15
NUMNODES      2
LINKNODE      1 4 0.426374
LINKNODE      1 5 0.607143
CELLIJ        8 15
NUMNODES      1
LINKNODE      1 5 0.065385
CELLIJ        8 14
NUMNODES      2
LINKNODE      1 5 0.327473
LINKNODE      1 6 0.771429
CELLIJ        8 13
NUMNODES      5
LINKNODE      1 6 0.228571
LINKNODE      2 6 0.251513
LINKNODE      3 2 0.225581
LINKNODE      4 1 1.000000
LINKNODE      4 2 0.104115
CELLIJ        3 10
NUMNODES      1
LINKNODE      2 1 0.562662
CELLIJ        4 10
NUMNODES      2
LINKNODE      2 1 0.437338
LINKNODE      2 2 0.333705
CELLIJ        4 11
NUMNODES      1
LINKNODE      2 2 0.525414
CELLIJ        5 11
NUMNODES      3
LINKNODE      2 2 0.140882
LINKNODE      2 3 1.000000
LINKNODE      2 4 0.005014
CELLIJ        5 12
NUMNODES      1
LINKNODE      2 4 0.150561
CELLIJ        6 12
NUMNODES      2
LINKNODE      2 4 0.844425
LINKNODE      2 5 0.452031
CELLIJ        7 12
NUMNODES      1
LINKNODE      2 5 0.224292
CELLIJ        7 13
NUMNODES      3
LINKNODE      2 5 0.323677
LINKNODE      2 6 0.748487
LINKNODE      3 2 0.024719
CELLIJ        6 15
NUMNODES      1
LINKNODE      3 1 0.546330
CELLIJ        6 14
NUMNODES      1
LINKNODE      3 1 0.065298
CELLIJ        7 14
NUMNODES      2
LINKNODE      3 1 0.388372
LINKNODE      3 2 0.749699
CELLIJ        9 13
NUMNODES      1
LINKNODE      4 2 0.120459
CELLIJ        9 12
NUMNODES      2
LINKNODE      4 2 0.775426
LINKNODE      4 3 0.473984
CELLIJ        10 12
NUMNODES      1
LINKNODE      4 3 0.296769
CELLIJ        10 11
NUMNODES      2
LINKNODE      4 3 0.229247
LINKNODE      4 4 0.843854
CELLIJ        11 11
NUMNODES      2
LINKNODE      4 4 0.156146
LINKNODE      4 5 0.316932
CELLIJ        11 10
NUMNODES      2
LINKNODE      4 5 0.683068
LINKNODE      5 1 0.221867
CELLIJ        12 10
NUMNODES      1
LINKNODE      5 1 0.527704
CELLIJ        12 9
NUMNODES      2
LINKNODE      5 1 0.250430
LINKNODE      5 2 0.643934
CELLIJ        13 9
NUMNODES      2
LINKNODE      5 2 0.356066
LINKNODE      6 1 0.092353
CELLIJ        13 8
NUMNODES      2
LINKNODE      6 1 0.907647
LINKNODE      6 2 0.066001
CELLIJ        14 8
NUMNODES      1
LINKNODE      6 2 0.369135
CELLIJ        14 7
NUMNODES      2
LINKNODE      6 2 0.564864
LINKNODE      6 3 0.488068
CELLIJ        15 7
NUMNODES      1
LINKNODE      6 3 0.289850
CELLIJ        15 6
NUMNODES      2
LINKNODE      6 3 0.222082
LINKNODE      6 4 0.910135
CELLIJ        16 6
NUMNODES      2
LINKNODE      6 4 0.089865
LINKNODE      6 5 0.120701
CELLIJ        16 5
NUMNODES      2
LINKNODE      6 5 0.879299
LINKNODE      6 6 0.332202
CELLIJ        17 5
NUMNODES      1
LINKNODE      6 6 0.131281
CELLIJ        17 4
NUMNODES      2
LINKNODE      6 6 0.536516
LINKNODE      6 7 0.754269
CELLIJ        18 4
NUMNODES      1
LINKNODE      6 7 0.051997
CELLIJ        18 3
NUMNODES      3
LINKNODE      6 7 0.193734
LINKNODE      6 8 1.000000
LINKNODE      9 1 0.158045
CELLIJ        18 2
NUMNODES      1
LINKNODE      9 1 0.023685
CELLIJ        19 2
NUMNODES      2
LINKNODE      9 1 0.818271
LINKNODE      9 2 0.581870
CELLIJ        19 1
NUMNODES      1
LINKNODE      9 2 0.065419
CELLIJ        20 1
NUMNODES      2
LINKNODE      9 2 0.352711
LINKNODE      9 3 1.000000


WMS builds the STREAM_CELL file from GIS like inputs developed in the WMS interface. Building this file without the WMS interface is extremely tedious. Without the interface, it would be most practical to employ a stream network that conforms to the grid on a one to one basis. That is, the stream nodes and the grids cells are the same size, and the stream nodes do not cut across grid cells.

5.1.5 Longitudinal Channel Profile Smoothing

When available, channel cross-sectional geometry and longitudinal profile data collected from extensive field surveys can be input into GSSHA. In lieu of detailed stream bed profiles, the DEM may be used as an indicator of channel slope, assuming that over large reaches, the slope of the channel thalweg is on the order of the same slope as the land surface. This assumption is not true in the case of extremely sinuous channels. Channel slopes extracted from the grid elevations are skewed by errors in the DEM and from errors in interpolating from the DEM to the resolution of the grid. These errors may result in regions of adverse channel slope, undulating streambed elevations, or sharp transitions in channel slope. These conditions will require a reduction in computational time step for the explicit diffusive wave routing method. It may be impossible to simulate channel flow with in a channel with many uncorrected DEM derived errors.

Extensive regions of adverse slope are rare in natural channels. It is recommended that regions of adverse slopes be removed from the GSSHA channel network, unless field observations or obvious geologic controls indicate that they are justified. Smoothing of the channel thalweg profile should be done to remove stream sections of adverse slope, undulating streambed elevations, zero slope, and sharp slope transitions.

Ogden et al (1994) developed a method of smoothing stream channel profiles based on the realization that in regions of concave topography near streams, DEMs are likely to be positively biased. For such conditions filling depressions along stream paths is probably erroneous. A more accurate approach is to remove higher regions to connect depressions. WMS provides tools to smooth the channel thalweg profile both manually and with an automated version of the Ogden et al. (1994) algorithm.

Changes to the stream thalweg elevation alone will not result in proper simulations of channel-groundwater exchange or channel-overbank exchange. To properly simulate these processes, the correct elevations, or at least the correct difference in elevation, of the channel thalweg and overbank area are needed. Upon initialization GSSHA will check in channel thalweg profile for problems and report them, if any, in the "check_thalweg" file.

5.1.6 Losing and Gaining Streams

The calculation of stream losses from streams in arid regions with coarse textured substrate can be specified by including the STREAM_LOSS card in the project file. When simulating either a static or moving water table, streams may be either losing or gaining streams, and the STREAM_LOSS card need not be included in the project file. In either case, the hydraulic conductivity (Krb) (cm hr-1) and thickness (Mrb) (cm) of the substrate must be specified. Uniform values can be specified in the project file with the K_RIVER and M_RIVER cards. Including these values in the project file results in every node in every stream link, regardless of specified type, having these uniform values of stream bed parameters. Alternately, the values can be specified in the CHAN_INPUT file for each link specified as as "SUBSURFACE", as described above.

Only channel links identified in the CHAN_INPUT file as "SUBSURFACE" will gain and lose. Links can be identified as "SUBSURFACE" with the "SUBSURFACE" card. Links can also be identified as "SUBSURFACE" by appending "SUBSURFACE" to the link type declaration, i.e. "TRAPEZOID_SUBSURFACE".

Any "SUBSURFACE" link without a specified "K_RIVER" and "M_RIVER" card will be assigned the uniform values of K_RIVER and M_RIVER. The default values for M_RIVER and K_RIVER are 1.0 and 0.0, respectively, such that unless values are specified in either the project file or the CHAN_INPUT file, there will be no channel/subsurface interaction.

Calculation of the flux term is discussed in Section 8.5.1. The explicit channel routing code was modified to allow water to remain in the channel after channel routing ends, and for there to be a starting volume of water in each cell at the beginning of channel routing. When simulating groundwater, channel routing continues after a storm event ends until no water is moving in the channel or only a very small amount of water is left in the stream channels. If the stream discharges to the ground; style="ater after channel routing has ceased, the volume of the stream is updated each groundwater update call. No channel routing is performed under these conditions. However, if the groundwater is discharging to the stream, channel routing will resume, but at the groundwater time step.

5.1.7 Sediment Transport in Channels

Sediment transport in the channel network can be performed as described in Section 10. Sediment transport in the channel network requires links be identified as "ERODE" in the CHAN_INPUT file. For each link identified as "ERODE", the maximum erosion (m) should be specified with the "MAX_EROSION" card in the CHAN_INPUT file. The default value is 0.0. Typically, if sediment routing is desired, all stream links should be identified as "ERODE", with either the "ERODE" card, or appending "ERODE" to the channel type,i.e. "TRAPEZOID_ERODE".


GSSHA User's Manual

5 Surface Water Routing
5.1     Channel Routing
5.2     Overland Flow Routing
5.3     Channel Boundary Conditions
5.4     Overland Boundary Conditions
5.5     Embankments
5.6     Overland/Channel Interaction
5.7     Introducing Discharge/Constituent Hydrographs
5.8     Overland Routing with Snow
5.9     Overland Routing with BMPs


Retrieved from "https://gsshawiki.com/index.php?title=Surface_Water_Routing:Channel_Routing&oldid=6357"