1.  INTRODUCTION

The subject of flood routing has preoccupied hydraulic and hydrologic engineers since the early 1900s. The task is to use an appropriate calculation to follow the progression or travel of a flood wave as it moves downstream along a river channel, from basin headwaters to mouth. Two practical applications are well established: (1) flood channel design, and (2) flood event forecasting.

The subject has important implications for the safety and overall well-being of society. Developed societies have an innate tendency to build on the floodplain, and this trend is clearly at the root of the recurrent flood problem. Barring fundamental changes, it is expected that societies around the globe will continue to have to cope with floods, their analysis, design, control, and management. Therefore, to this date, flood wave modeling remains a worthwhile endeavor.

In this paper, we review the issues concerning flood routing, including the available methods and their accuracy and practicality to satisfy specific objectives. We examine the related concepts of flood routing, methodologies, design, forecasting, and modeling, in order to clarify the issue for the hydraulic/hydrologic engineer engaged in this field of work. The aim is to help establish the Muskingum-Cunge method as the method of choice for a variety of flood routing applications, given its demonstrably sound theoretical basis and the substantial body of knowledge that has accrued in the more than 50 years since its original development.

2.  FLOODS

Floods are perceived to be both good and bad. Floods are good when they move sediments out of their place of origin, typically at or near basin headwaters, and move them downstream, along established flow paths, to lower elevations, where they are subject to deposition. This is the way valleys were formed, and valleys do remain an essential geomorphic feature of the contemporary societal landscape. The deposited soils feature a diversity of nutrients, which eventually make possible all types of agriculture.

Floods are bad in situations when rivers overflow their banks, eventually spreading floodwaters in developed urban settings. Typically, floodwaters carry sediments, which invariably will have a tendency to settle in unwanted places. Therefore, the issue of how to best manage floods remains in the eye of the beholder; it will depend largely on the local situation. In urban settings, floods are generally considered bad; therefore, they are to be avoided, managed, and/or controlled. This reality gives rise to the field of flood engineering, which encompasses a variety of strategies, including flood control, forecasting, management, avoidance, and, ultimately, flood tolerance.

At this juncture (2023), the state-of-the-art regarding the correct approach to floods lies somewhere in between the seemingly contrasting disciplines of flood control and flood management. Flood control is primarily aimed at fighting floods through conventional measures such as reservoirs, canals, and levees. Flood management recognizes that Nature has the advantage of time in her hands, and that solutions that do not heed this fact are eventually destined to fail. Therefore, the matter remains unsetlled. We trust that the future will provide ample experience to help resolve the issue of how to manage floods in the best way.

3.  FLOOD ROUTING

Flood routing is the process of following, through calculation, the progression or travel of a flood wave as it moves downstream along a channel course, stream, or river. At the outset, two distinct wave properties are recognized: (1) the flood wave celerity, i.e., its rate of speed; and (2) the flood wave's rate of attenuation or dissipation, described by the logarithmic decrement (Wylie, 1966; Ponce and Simons, 1977).

The wave celerity is physically related to flow concentration, described by a differential equation of first order. The wave attenuation is physically related to flow diffusion, described by a differential equation of second order. Third-order processes (dispersion) are indeed feasible (Ponce, 2020), but at this juncture, they are not utilized in practical flood routing applications.

Hydraulic engineers have delved into flood routing issues since the early 1900s. Seddon (1900), working with floods in the Lower Mississippi river, derived an expression for the celerity of a flood wave, among the earliest attempts to focus on the science. The basic elements of Seddon's findings are detailed in Box A.

4.  MUSKINGUM METHOD

The basic elements of the Muskingum method are described in Box B.

Box B.  Elements of the Muskingum method (McCarthy, 1938; Ponce, 2014). Reach inflow I Reach outlow O Reach storage S Differential equation of storage:  I - O = dS/dt Storage relation:  S = K [ X I + ( 1 - X ) O ] K:  Routing parameter, a storage coefficient, or reach travel time X:  Routing parameter, a dimensionless weighting factor, varying in the range 0 ≤ X ≤ 0.5 Routing equation (four-point scheme):   O2 = C0 I2 + C1 I1 + C2 O1 Routing coefficients:                  ( Δt / K ) - 2X C0 = _______________________              2(1 - X) + ( Δt / K )                  ( Δt / K ) + 2X C1 = _______________________              2(1 - X) + ( Δt / K )              2(1 - X) - ( Δt / K ) C2 = _______________________              2(1 - X) + ( Δt / K )

The conventional Muskingum method (of McCarthy) is straightforward and comparatively simple; however, considerable care is required in order to use the method effectively. Its main drawback is the strict requirement for calibration of the routing parameters, which is at best a cumbersome procedure (Ponce, 2014). A typical routing application will necessarily require a calibration of the method's parameters, essentially limiting the predictive stage (of the routing) to flood magnitudes and durations similar to those used in the calibration. This requirement substantially limits the predictive ability of the method, jeopardizing its wider applicability to other floods and/or other reaches located within the same hydrologic system. Thus, we conclude that the conventional Muskingum method is largely ineffective for its use in extensive basin hydrologic modeling, which has become the norm in more recent times.

5.  MUSKINGUM-CUNGE METHOD

Cunge (1969) substantially improved the conventional Muskingum method, giving it a decisive hydraulic flavor. He accomplished this feat by recognizing that the conventional formulation, which related reach storage to inflow and outflow through the routing parameters K and X (see Box B), could be construed as a numerical analog of the kinematic wave equation (Lighthill and Whitham, 1955). In other words, when taken in the proper context, the original McCarthy formulation (1938) was indeed a numerical solution and, thus, subject to numerical laws. We note that the conventional Muskingum method preceded the computer age and, therefore, could not click on to the new mindset that followed the availability of digital computers, ostensibly by the early 1960s.

Cunge's analysis did not limit itself to the kinematic wave, and this remains his greatest contribution. He was able to relate the leading error (the second-order error) of the numerical analog (of the Muskingum method) to the diffusion term of the kinematic-with-diffusion wave, i.e., the diffusion wave (Ponce and Simons, 1977). This made possible the derivation of an expression for the weighting factor X, substantially enhancing the methodology by giving it a deterministic flavor, since the diffusion wave is based on fundamental principles of mechanics.

Cunge's masterful use of physical and numerical principles provided the enhanced methodology, now recognized as the Muskingum-Cunge method (Flood Studies Report, 1975; Ponce and Yevjevich, 1978), featuring the important property of grid independence, which had been unattainable up to that point. In point of fact, by relating the second-order error of the numerical analog of the Muskingum-Cunge method to the intrinsic properties of the scheme (the Courant and cell Reynolds numbers), the pot of gold at the end of the rainbow could be envisioned. Indeed, the method has amply demonstrated to be essentially grid-independent, that is, the calculated flood hydrographs turn out to be the same, regardless of grid size. The basic elements of the Muskingum-Cunge method are described in Box C.

Box C.  Elements of the Muskingum-Cunge method (Cunge, 1969; Ponce, 2014). Flood discharge Q Equilibrium unit-width discharge qo Channel (stream, bottom) slope So Flood wave celerity (kinematic wave celerity):  c Reach length (space step):  Δx Temporal interval (time step):  Δt Courant number C = c (Δt /Δx) Cell Reynolds number D = qo / (So c Δx ) Routing equation (four-point scheme):   Q j+1 n+1 = C0 Q j n+1 + C1 Q j n + C2 Q j+1 n Routing coefficients              -1 + C + D C0  =  ______________               1 + C + D               1 + C - D C1  =  ______________               1 + C + D               1 - C + D C2  =  ______________               1 + C + D

We note that the fundamental tenet of the Muskingum-Cunge method, and the reason why it works so well, is the effective matching of physical and numerical diffusivities. The physical diffusivity is the (channel) hydraulic diffuvisity ν_p, originally due to Hayami (1951): ν_p = q_o /(2S_o).

The numerical diffusivity ν_n is the diffusivity derived by Cunge (1969): ν_n = c Δx (0.5 - X). The matching of physical and numerical diffusivities leads to the predictive equation for X, the weighting factor of the Muskingum method, solely in terms of hydraulic variables, as follows: X = 0.5 { 1 - [ q_o / (S_o c Δx ) ] }. Furthermore, in terms of the routing parameter cell Reynolds number D (line 8 of Box C), the formula for X reduces to: X = 0.5 (1 - D).

It is clearly seen that the operational capability of the Muskingum-Cunge method hinges on the determination of the values of the two routing parameters: (1) Courant number C, the ratio of physical celerity (the flood wave celerity, or Seddon celerity) to "numerical celerity," Δx/Δt; and (2) cell Reynolds number D, the ratio of physical diffusivity q_o /(S_o) to numerical diffusivity (c Δx). Data needs for the application of the Muskingum-Cunge method are reviewed in Box D.