Critical Depth Primer

Written by Christopher Goodell, P.E., D.WRE  |  WEST Consultants Copyright © The RAS Solution 2014.  All rights reserved.   Critical depth is an important hydraulic parameter when evaluating hydraulic modeling results.  As we know from our college hydraulics 101 class, water flowing at depths less than critical depth is supercritical and water flowing at depths greater than critical depth is subcritical.  Supercritical flow is characterized by relatively shallow depths and high velocities.  Subcritical flow is characterized by relatively deep depths and slower velocities.  The forces dominating the movement of supercritical flow are inertial, while the forces dominating the movement of subcritical flow are gravitational.  The flow regime (subcritical or supercritical) a particular cross section, or series of cross sections is in, can be determined by the Froude Number, , where F = Froude Number, V = Velocity, g = gravitational constant, and d = depth. A Froude Number greater than 1 indicates supercritical flow, a Froude Number less than 1 indicates subcritical.  A Froude Number = 1 is considered “critical” and possesses the minimum amount of specific energy (potential energy plus kinematic energy, per unit mass).  This is considered an unstable condition in nature and is very rare.  If you ever get a Froude Number = 1 in your results, most likely, there is a problem with the computations (i.e. in steady flow, RAS could not come up with a valid solution, so it defaults to critical depth).  For those of you who have run steady flow HEC-RAS models, you know that there are three flow regime options for computing a steady flow run:  subcritical, supercritical, and mixed flow (both sub and supercritical).  If you select RAS to run in subcritical, and somewhere in your system RAS is not able to compute a valid subcritical answer, then it defaults to critical depth and moves on.  If you select RAS to run in supercritical, and somewhere in your system RAS is not able to compute a valid supercritical answer, then it also defaults to critical depth and moves on.  If you select mixed flow, RAS will compute both a subcritical and supercritical profile and anywhere there is a valid solution for both regimes, RAS will select the one that has the higher specific force value.  If, in mixed flow, there are any cross sections that “default to critical depth”, that means there was a problem with RAS obtaining a solution.  Usually, this means your cross section spacing is too far apart, you are in an area of rapidly varied flow (with not enough cross sections), or just bad input data.  In unsteady flow, contrary to intuition, checking the “Mixed Flow” box in the unsteady flow analysis window does NOT tell RAS to evaluate both sub and supercritical solutions.  This is done anyway in unsteady flow-whether “Mixed Flow Regime” is checked or not.  Let me repeat:  HEC-RAS can compute solutions in both subcritical and supercritical in unsteady flow regardless of whether “Mixed Flow Regime” is checked or not.  Checking the “Mixed Flow Regime” box in the unsteady flow analysis window simply uses a stabilizing scheme for situations near critical depth and with large changes in velocity with respect to time (the local acceleration term in the St. Venant Equation).  This is described further in http://hecrasmodel.blogspot.com/2011/04/mixed-flow-regime-options-lpi-method.html.   It’s very helpful to understand the solution you’ve obtained after running RAS by viewing the water surface profile plot with the “critical depth” variable turned on.  This allows you to gage how close you are to critical depth, and more importantly, allows you to quickly evaluate if your solution has defaulted to critical depth anywhere, indicating a problem with the solution.  The critical depth variable can be turned on in any of the graphical plots by selecting Options…Variables.  Then check the box next to Critical Depth Elevation.    When you do this, RAS will plot critical depth, but ONLY at certain locations.  Notice the plot from the Single Bridge HEC-RAS example data set.  Critical depth (in red) is only plotted downstream of the bridge, and at a couple of cross sections upstream of the bridge.  That is because RAS will only compute critical depth if your answer is supercritical, close to critical depth, or RAS is not able to come up with a valid solution (defaulting to critical depth) and at the boundaries.  If you would like HEC-RAS to compute critical depth everywhere for you, go to the steady flow analysis window, and select Options…Critical Depth Output Option. 
Then check the box next to “Critical Always Calculated.

Stabilizing a Dynamic Unsteady HEC-RAS Model

Written by Chris Goodell
Copyright © RASModel.com. 2013.  All rights reserved.

One of the most frustrating aspects of unsteady HEC-RAS modeling can be the model stabilization process. You know, you’ve gone to great lengths collecting the best survey/topo data and solid hydrology. Then you’ve painstakingly spent hours…possibly days entered all of that data only to find that once you press the “Compute” button, the model crashes. The dreaded “Red Bar”! Sometimes you can get your simulation to complete without crashing, but the listed numerical errors are so high that you can’t with good conscience submit that as your final simulation. Either way, approaching an unsteady HEC-RAS model (especially a dynamic one) as a beginner with little experience and understanding of how to stabilize it can cause significant delays in your project and worse, completely blow up your budget. I’ve uploaded two HEC-RAS projects to the following Google Drive Site: https://drive.google.com/folderview?id=0B_s8OLJOgOi0YU92SGZQZm9raWM&usp=sharing



RawleyResUnstable.prj is an unsteady flow dam breach HEC-RAS project recently sent to me for help. Although the model ran to completion without crashing, it had unacceptably high errors.
RawleyResStable.prj is the fully stabilized version of the model with no numerical errors.

The following lists out the courses of action taken to stabilize the model. Feel free to download the “Unstable” and “Stable” models and try these techniques on your own. The links following some of these items will take you to more information about that particular technique.

1. Cross Section Spacing. The initial spacing was way too coarse. A visual check alone of the geometry schematic and profile plot should encourage you to investigate a finer cross section spacing.

Geometry Schematic

Profile Plot Samuels equation suggests anywhere from 15 ft to 50 ft spacing (depending on what bed slope you use). I interpolated to 50 ft for the entire reach. http://hecrasmodel.blogspot.com/2008/12/samuels-equation-for-cross-section.html

2. Dam Breach models typically have time steps on the order of a minute or less. This model was initially set with a computation interval of 10 minutes which is high even for the largest and “slowest” of dam breach models. I changed the time step to 10 seconds. The selection of 10 seconds was based on “gut” feel and lots of experience doing dam breach models. There are some methods for approximating good timesteps, notably the Courant Condition and Fread’s equation. Also, setting the cross section spacing/timestep ratio equal to a representative stream velocity will get you close. In this case, the cross section spacing is 50 ft, and I was able to pull some velocities (prior to the model crashing) at about 6 ft/s. That suggests a time step of 8.3 seconds. Close enough to 10 seconds, so we’ll stick with that.

3. I changed the downstream boundary normal depth slope from 0.01 to 0.001. Not sure what is downstream of the first cross section, but 0.01 is awfully steep and was setting up a very low depth at the downstream boundary (which was causing instabilities). http://hecrasmodel.blogspot.com/2010/01/downstream-boundary-normal-depth.html

Notice in the following profile plot of the downstream end of the reach how the water surface at the boundary cross section is below critical depth (the red dot). This creates an overestimation of the water depth at the next upstream cross section, which in turns creates some instability over the next several timesteps.

The proper way to handle this would be to find out what is downstream of your model and select a boundary condition that best represents those conditions. In absence of downstream conditions, the 0.001 slope for normal depth provides a reasonable solution.  Either way, this underscores the importance of moving your downstream boundary far away from your area of interest in your study reach.  That way the errors that do originate from your downstream boundary assumption will have diminished to negligible levels before impacting your area of interest. 

4. In the unsteady flow editor, the initial flow and the first time step flow should always be equal. In the original model, the initial flow is left blank (which is actually okay because RAS will use the first timestep flow if left blank). However, the first time step flow is very low at 0.14 cfs. Peak discharge of your inflow hydrograph is around 600 cfs. I put in a baseflow (minimum flow) of 5% of this, which is 30 cfs. Just to be safe, I put 30 cfs in the initial flow input box as well. It’s important to make sure that the 30 cfs baseflow does not have an impact on the peak of the breach outflow hydrograph. http://hecrasmodel.blogspot.com/2009/02/minimum-flow-requirements.html

5. The upstream end of your reservoir is very steep (14%). This one was easy to spot in the profile plot because of the very nice smooth profile plot followed by a sudden spike in the energy level (green dashed line) at the upstream end. Notice to that there is a corresponding supercritical solution at the upstream most cross section. These are both indications that more cross sections are needed here. http://hecrasmodel.blogspot.com/2010/01/another-reason-for-interpolated-cross.html

I added more cross sections by interpolating the steep slope at the upstream end of the reservoir (10 ft spacing).

Notice the energy spike is not completely gone, but it is much better. Refining the HTAB parameters (http://hecrasmodel.blogspot.com/2009/02/crazy-energy-grade-line.html), turning on mixed flow and bumping up n values will collectively take care of the rest of this (see numbers 8, 9, and 10 below).

7. From RS 4577 to 3791 you have your main channel defined in a small elevated side channel. The right overbank is lower in elevation than the channel.

This is problematic in RAS. Using the graphical cross section editor, I redefined the bank stations to get rid of this problem.

8. Notice now in the HTAB parameters property table, the HTAB definition begins at the elevation of the invert of the side channel that was previously defined as the main channel. There is a big gap in the computation points in the newly defined main channel. This poses problems with HEC-RAS, particularly at low discharges/stages, since RAS will have to extrapolate to obtain a solution over that range.

I readjusted the cross section HTAB parameters to reflect the new bank station definitions (I used the “Copy Invert” button for all cross sections). While in there, I maximized the resolution of the HTAB parameters by changing number of points to 100 and minimizing the increment as much as possible while still having full coverage of each cross section. http://hecrasmodel.blogspot.com/2009/02/crazy-energy-grade-line.html http://hecrasmodel.blogspot.com/2011/03/more-on-htab-parameters.html

9. In the figure below, notice all of the “red” areas in this zoomed in section of the  upstream end of the reach.

These are areas that are near critical depth. I turned on the “Mixed Flow” option in the Unsteady Flow Analysis window, which helps stabilize the near critical depths at the upstream end of the reservoir. http://hecrasmodel.blogspot.com/2011/04/mixed-flow-regime-options-lpi-method.html If you’ve been following along correcting your own copy of the model, you’ll notice that the model now finally runs to completion without crashing. There are still a few minor errors as shown in the Computation Messages.

It’s debatable whether you need to take care of these relatively small errors. Particularly for a dam breach model, where the presumed errors of the input data probably far outweigh these small numerical errors. Nevertheless, there’s a sense of pride in putting out a very robust model, completely free of numerical errors. So…let’s continue.

10. The errors in the computation window above suggest a problem around river station 6200. That happens to be the upstream end of the reservoir and the source of the error can clearly be seen with the energy grade line and critical depth turned on.

The section at the upstream end of the reservoir is very steep at about 13.5%. In fact, that’s greater than the suggested maximum slope of 10%, as stated in the HEC-RAS manuals. You can see as the pool starts to lower, the steep reach is exposed and because of the low n values, the water surface is calculated to be supercritical. Let’s assume that the geometry is correct, and the upstream end of the reservoir really is that steep. Jarrett’s Equation suggests very high n values should be used here (~0.2), based on the bed slope and hydraulic radius (after the breach when the reservoir has drained). The original model had 0.07 and 0.035 for the overbanks and main channel, respectively. I changed all the n values in this steep reach to 0.15, because I think 0.2 might be a little overkill. http://hecrasmodel.blogspot.com/2009/12/n-values-in-steep-streams.html

There…no errors! However, there is still a Warning about extrapolating above/beyond the rating curve at a bridge (R.S. 277).

11. To get rid of that warning, I increased the HTAB headwater maximum elevation at River Station 277 (the bridge) from 278 ft to 280 ft so that RAS doesn’t have to extrapolate.

There we have it. A clean solution. No errors, no warnings.

Mixed Flow Regime Options – LPI Method

Written by Aaron A. Lee | WEST Consultants
Copyright © RASModel.com. 2011. All rights reserved.
By using the Mixed Flow Regime option for Unsteady Flow Analysis, RAS can better handle transitions from subcritical to supercritical flow. This option should be utilized only after determining that a mixed flow situation exists, which requires judgment from the modeler. One application where this could be particularly useful is dam breach modeling, or any other extreme and flashy flood event. Even though a model is stable there may still be small errors in the solution (caused by max. iterations). The Local Partial Inertia (LPI) factor may eliminate or reduce these errors, particularly if they occur when the Froude number is near 1. Figure 1 shows the Unsteady Flow Analysis window with the Mixed Flow Regime option selected. This post will focus on the LPI Filter, which is enabled when Mixed Flow Regime is selected by the modeler.





Figure 1. Unsteady Flow Analysis Window


Once the Mixed Flow Regime option is selected, additional settings can be adjusted to help stabilize the model. Navigate to Options, Mixed Flow Options. This window, shown in Figure 2, allows the user to adjust two inputs for the LPI factor.





Figure 2. Mixed Flow Options Window


For the unsteady flow computation scheme, RAS accounts for a local acceleration and convective acceleration (inertial terms) through the St. Venant equation of Conservation of Momentum. The St. Venant equations, and by extension, HEC-RAS, are designed to work best in gradually varied flow. Transitions from supercritical to subcritical flow (hydraulic jump), and to a lesser extent subcritical flow to supercritical flow, are rapidly varied flow situations. These are not gradual changes, in the hydraulic sense. Near critical depth (Froude number approaching 1) the convective acceleration terms can change very rapidly over a short distance (think of a hydraulic jump) and can lead to oscillations in the solution. These oscillations tend to grow larger until the solution goes completely unstable (HEC, 2010). The LPI factor systematically reduces these inertial terms to dampen the oscillations, helping to stabilize the model. The user can influence the magnitude of reduction by varying the two inputs in Figure 2.


The first input, m, is the exponent for Froude number reduction factor. Its default value is 10 and ranges from 1 to 128. Adjusting m will change the shape of the curve on Figure 2, thus influencing the rate of reduction of the inertial terms. You can see that by making m smaller there is an earlier and more direct reduction in the inertial terms, with respect to the Froude number. Increasing m can make the model more accurate but increases the likelihood of numerical instability.


The second input, F_T, is the Froude number threshold at which the LPI factor is set to zero. In other words, if the calculated Froude number at the current cross-section is larger than F_T the inertial terms will be eliminated from the computations at that cross-section for the current computational time step. The default value is 1. Making F_T smaller will improve the stability of the model, but will also reduce the accuracy. A larger F_T can make the model more accurate, but increases the likelihood of numerical instability as the inertial terms will be more sensitive to fluctuations in Froude number.


A good place to start is to run the simulation with the default values

Technorati Tags: HEC-RAS,RAS,Mixed Flow,LPI Method,Unsteady Flow,Supercritical,Subcritical,Hydraulic Jump,Froude,Water Surface Profile,Dam Breach,Dam Break,Flood,Numerical Model,St. Venant

and see what the profile looks like. For this flume example, the model ran without reporting any maximum water surface errors. The profile for the default LPI inputs is shown in Figure 3.





Figure 3. Profile Plot, Default Values


Next, a value of 1.6 was chosen for F_T. This simulation yielded small maximum water surface errors, but had maximum iterations at various locations. The value of m was left unchanged. Figure 4 shows the results.





Figure 4. Profile Plot, Increased Froude Number Elimination Threshold


Even though the errors were small, instabilities could be seen in the downstream end. Notice the instabilities around the transitions between the flow regimes. The value of m was reduced from the default of 10 to 7 in order to improve the stability of the model. Figure 5 shows the profile for reduced m and increased F_T. The modeler should choose the largest values of m and F_T that produce a stable model. However, check the results to make sure that the output is reasonable. Notice how the transitions between flow regimes are much better defined in Figure 5 then the default setup shown in Figure 3. That’s because the default LPI parameters (m = 10 and F_T = 1) provide dampening of the results. Though Figure 3 looks very stable (and it is), Figure 5 (m = 7 and F_T = 1.6) is both stable and (by my engineering judgment) more accurate. Also, notice how the slight increase in energy (green dashed line) is less in Figure 5 versus Figure 3. An increase in the energy elevation in the direction of flow is an indication of error in most cases. Further adjustment of the LPI parameters may help to eliminate the error in the energy grade line, while still producing a stable solution.





Figure 5. Profile Plot, Increased Froude Number Elimination Threshold and Decreased Exponent ,m