1.  INTRODUCTION

There is a continuing controversy regarding the nature and applicability of the kinematic wave model. Researchers and practitioners alike have reported successes and failures of the model, with papers continuing to appear in the literature describing what the model can and cannot do (Hromadka and DeVries 1988). Current areas of concern focus on the following issues: (1) Whether the kinematic wave can describe physical diffusion, and if so, under what circumstances; (2) whether the kinematic wave can eventually replace other well-established methods of surface runoff generation such as the unit hydrograph; and (3) whether the kinematic shock is as common in practice as calculations would seem to indicate.

While answers to these questions can be found in the literature, they are dispersed among various sources and not readily accessible. This difficulty appears to be fueling the current controversy (Dawdy 1990; Goldman 1990; Hromadka and DeVries 1990; Merkel 1990; Unkrich and Woolhiser 1990; Woolhiser and Goodrich 1990). Therefore, the aim of the present paper is to review the concept of kinematic wave, delineate its range of applicability, and critically examine its overall modeling philosophy. It is hoped that this review will help focus the attention of researchers and practitioners so that the controversy may be brought to a timely end.

2.  BACKGROUND

The concept of kinematic wave is well established among the existing methods to solve unsteady, one-dimensional, gradually varied open-channel flow problems. In contrast to the dynamic wave, which features a significant inertial component, a kinematic wave is one in which the inertial component is too small to be of any practical importance. In unsteady open-channel flow modeling, a first type of kinematic wave arises when the governing equations are simplified by neglecting the local inertia, convective inertia, pressure-gradient, and momentum-source terms (Lighthill and Whitham 1955). A second, less restrictive type can be formulated by neglecting the local inertia, convective inertia, and momentum-source terms, but keeping the pressure-gradient term (Hayami 1951). To avoid confusion between these two types of kinematic waves, it is common practice to refer to the first type as kinematic wave proper, and to the second as diffusion wave (Ponce and Simons 1977).

From the physical standpoint, the kinematic wave assumption amounts to substituting a uniform flow formula (such as Manning's or Chezy's) for the equation of motion. In essence, it says that as far as momentum is concerned, the flow can be considered steady. The unsteadiness of the phenomenon, however, is preserved through the rate-of-rise term in the continuity equation (Liggett 1975). The implication of the kinematic wave assumption is that unsteady open-channel flow can be visualized as a succession of steady uniform flows, with the water surface slope remaining constant at all times. This, of course, can be reconciled with reality only if the flow unsteadiness is very mild, i.e., if the changes in momentum are indeed negligible compared to the forces driving the steady component of the motion (gravity and friction).

From the mathematical standpoint, the kinematic wave assumption results in a considerable simplification of the equation of motion, reducing it to a statement of uniform flow (such as, for instance, the Manning equation). Combining this latter equation with the equation of continuity gives also to a first-order partial differential equation, referred to as the kinematic wave equation:

∂Q             ∂Q     _____  +  c  _____  =  c qL   ∂t               ∂x

(1)

in which Q = discharge; c = kinematic wave celerity; q_L = lateral inflow; x = spatial variable; and t = temporal variable. This equation is applicable to streamflow modeling as well as to channel and gutter flow. For overland flow applications, the kinematic wave equation is expressed in terms of unit-width discharge as follows:

∂q              ∂q   _____  +  c  _____  =  c i   ∂t               ∂x

(2)

in which q = unit-width discharge; i = effective rainfall intensity; and the other terms are as defined previously.

The kinematic wave celerity is defined as the slope of the rating curve, either discharge-flow area (Q = αA^β ) in the case of streamflow, or unit- width discharge-flow depth (q = ad ^m ) for overland flow. Accordingly, c = dQ/dA = β(Q/A) = βu in the case of streamflow; and likewise, c = m (q /d ) = mu for overland flow, in which u = mean velocity. In natural channels, the kinematic wave celerity is alternatively expressed as c = (1/T ) (dQ/dy ), in which T = channel top width, and y = stage.

Equation 1 (and by extension, Eq. 2) is a differential equation of first order; therefore, it can describe convection but not diffusion, which is a second-order process. In practice, this means that the kinematic wave equation can describe the travel of a flood wave, but not its attenuation as it propagates downstream. Although Eq. 1 cannot describe diffusion, it is a quasi-linear equation because the kinematic wave celerity is a function of discharge. This gives kinematic waves the tendency to change in shape as they propagate. If the celerity increases with discharge, the leading face of the wave will steepen; conversely, if the celerity decreases with discharge, the leading face of the wave will flatten out. In overland flow and inbank streamflow, the tendency is for the wave to steepen; in shallow-overbank streamflow, the tendency is for the wave to flatten out.

3.  KINEMATIC WAVE SOLUTIONS

Solutions to Eq. 1 (or Eq. 2) can be attempted in a variety of ways. Analytical solutions are possible for linearized analogs of the governing equations [see for instance Lighthill and Whitham (1955) and Ponce and Simons (1977)]. These solutions describe the convection of a flow quantity (either Q or q) with the celerity c in the absence of diffusion. Numerical solutions are possible using the method of characteristics or the finite difference method. Early work on the kinematic wave used the method of characteristics. In overland flow applications, the nonlinearity (or rather the quasi-linearity) of the phenomena usually led to wave steepening and the eventual development of a kinematic shock. A kinematic shock is a kinematic wave that has steepened to the point where its rising limb has an almost vertical face, wherein the flow develops a singularity and loses its gradually varied property.

4.  NUMERICAL DIFFUSION AND DISPERSION

Although kinematic wave solutions using the method of characteristics are prone to shock development, ostensibly because of their lack of diffusion, solutions using the finite difference method exhibit a somewhat different behavior. Finite difference solutions, by virtue of their discrete nature, introduce appreciable amounts of numerical diffusion and numerical dispersion. These numerical effects interfere with the physical effects, modifying them (Abbott 1976). For instance, in overland flow applications, the numerical diffusion has the effect of counteracting the tendency of the wave to steepen, thereby arresting shock development and allowing the unsteady gradually varied flow computation to continue.

The presence of numerical diffusion and dispersion in a numerical solution using the finite difference method is at the crux of the controversy surrounding the kinematic wave model. A numerical scheme is characterized by its amplitude and phase error portraits (Leendertse 1967). The amplitude portrait describes the way in which the numerical wave amplitude approaches the physical wave amplitude, with the deviation being interpreted as numerical diffusion. The phase portrait describes the way in which the numerical wave phase approaches the physical wave phase, with the deviation being interpreted as numerical dispersion. Examples of amplitude and phase portraits for convection problems are given by Cunge (1969) and Ponce et al. (1979).

Finite difference solutions of the kinematic wave equation exhibit varying amounts of numerical diffusion and dispersion depending on the type of scheme used to discretize Eq. 1 (or Eq. 2). Fully centered schemes are of second order, exhibiting no numerical diffusion. However, these schemes may exhibit numerical dispersion for Courant numbers other than 1. (In kinematic wave theory, the Courant number is defined as the ratio of the physical celerity, i.e. the kinematic wave celerity c, to the "grid celerity" Δx /Δt, with Δx = the spatial increment, or space step, and Δt = the temporal increment, or time step.) Off-centered schemes are of first order, exhibiting varying amounts of numerical diffusion, depending on the size of Δx and Δt. Smaller increments result in smaller amounts of numerical diffusion, with the numerical diffusion vanishing as the increments are driven to zero. These schemes also exhibit variable amounts of numerical dispersion for Courant numbers other than 1.

Numerical diffusion arises due to the neglect of the second-order term of the corresponding Taylor series expansion of the related discrete analog. Numerical dispersion arises due to the neglect of the third-order term of the Taylor series expansion (Cunge 1969; Ponce et al. 1979). Therefore, in a typical application, numerical diffusion is usually about an order of magnitude greater than numerical dispersion. In practice, this means that as a rule, most of the error of numerical solutions can be attributed to numerical diffusion rather than to numerical dispersion. Exceptions are the cases in which the Courant number is significantly less than 1, in which case the amount of numerical dispersion may grow to the point where it compares in size with the amount of numerical diffusion.

Numerical diffusion manifest itself as the diffusion or attenuation of the calculated runoff hydryograph. Since the kinematic wave equation has no built-in physical diffusion, it follows that a finite difference solution is actually simulating physical diffusion through numerical diffusion. The fact that the latter is artificial, and intrinsically related to the grid size, can be readily demonstrated by solving the same problem several times, each time halving the spatial and temporal increments [see for instance Ponce (1986)]. Carried to a practical limit, this test leads to the eventual disappearance of the numerical diffusion in question, with the result approaching the analytical solution of the kinematic wave. The recognition of this fact prompted the USDA Soil Conservation Service to retire its convex method of flood routing in the early 1980s. The convex method showed considerable sensitivity to the grid size, with the numerical diffusion vanishing as the grid size was gradually reduced.

Numerical dispersion manifests itself as dispersion; that is, as the steepening or flattening of the rising limb of the calculated runoff hydrograph. In certain extreme cases, the numerical dispersion is responsible for the wiggles, or the (usually small) negative outflows observed at the beginning or end of the calculated runoff hydrograph (Hjelmfelt 1985). In practice, these negative outflows are small, and disappear (together with all traces of numerical dispersion) as the space and time steps are chosen such that the Courant number approaches 1.

Since numerical diffusion and dispersion are inherent in the choice of space and time steps (and their ratio relative to the wave celerity), the result of a finite difference solution of Eqs. 1 or 2 is per force grid dependent; i.e., the calculated runoff hydrograph varies with the choice of grid size. Therefore, it seems pointless to try to "calibrate" a kinematic wave model by varying a physical parameter such as Manning's n in order to match calculated results with observed data. This practice amounts to curve-fitting; at best it is good conceptual modeling, but it should not be interpreted as deterministic modeling.

5.  DIFFUSION WAVE MODELING

In practice, actual runoff hydrographs do exhibit a certain amount of diffusion. To properly simulate this diffusion it is necessary to extend kinematic wave theory to encompass the related diffusion wave theory. Following Hayami (1951) and Lighthill and Whitham (1955), the diffusion wave equation is derived by neglecting the local inertia, convective inertia, and momentum-source terms in the equation of motion, leading to the following equation for streamflow and channel flow:

∂Q             ∂Q             ∂ 2Q _____  +  c  _____  =  v  ______  +  c qL   ∂t               ∂x              ∂x 2

(3)

and a similar equation for overland flow:

∂q              ∂q              ∂ 2q _____  +  c  _____  =  v  _____  +  c i   ∂t               ∂x              ∂x 2

(4)

with v = hydraulic diffusivity, defined as follows:

Q             q               v  =  ______  =  ______          2TSo         2So

(5)

and S_o = bottom slope; and other terms are as previously defined. Eqs. 3 and 4 describe the movement of kinematic waves with a diffusion component, for short diffusion waves. Unlike their counterparts Eqs. 1 and 2, Eqs. 3 and 4 are second-order parabolic differential equations and are, therefore, able to describe physical diffusion, with the diffusion coefficient defined by Eq. 5.

The diffusion wave equation can be solved either analytically, leading to Hayami's diffusion-analogy solution for flood waves, or numerically, with the aid of a numerical scheme for parabolic equations, such as the Crank-Nicolson scheme (Crandall 1956). An alternative is to extend the finite difference solution of the kinematic wave to the realm of diffusion waves by matching physical and numerical diffusivities (Cunge 1969; Dooge 1973). The physical diffusivity is the hydraulic diffusivity given by Eq. 5. The numerical diffusivity is the numerical diffusion coefficient of the discretized kinematic wave model, i.e., the coefficient of the leading (second-order) error term. When the Muskingum scheme is used to model the kinematic wave, this method is referred to as Muskingum-Cunge model (Flood 1975; Ponce and Yevjevich 1979).

The Muskingum-Cunge method has the significant advantage over conventional kinematic wave models of being essentially grid independent (Ponce and Theurer 1982; Ponce 1986). Therefore, calibration requires only a fine tuning of the frictional and cross-sectional parameters such as Manning's n and flow rating exponent (either β for channel flow or m for overland flow). Then, the choice of spatial and temporal increments is based solely on considerations of problem scale.

6.  APPLICABILITY OF KINEMATIC WAVES

The question of the applicability of kinematic waves has interested researchers and practitioners alike. Although kinematic waves were originally used for describing river flows (Seddon 1900), it is in the field of overland flow that questions first arose with respect to its applicability and accuracy. Notable among these contributions is that of Woolhiser and Liggett (1967), who identified a parameter k defined as follows:

SoLo               k  =  _______          doFo2

(6)

to characterize the applicability of the kinematic wave to overland flow situations. In Eq. 6, S_o = bottom slope; L_o = length of the overland flow plane; d_o = normal depth; and F_o = Froude number based on normal flow. The parameter k has been referred to in the literature as kinematic flow number (Liggett 1975). This parameter can be used as a criterion to aid in determining whether the kinematic wave solution is sufficiently accurate when used to solve overland flow problems. According to Woolhiser and Liggett (1967), a value of k ≥ 20 indicates that the flow is strongly kinematic, and therefore suited to solution using the kinematic wave equation. More recently, however, Morris and Woolhiser (1980) have stated that for low-Froude-number flows, it is also necessary that k F_o² ≥ 5, a condition compatible with the Woolhiser and Liggett criterion (k ≥ 20) for the case of F_o ≥ 0.5.

Ponce et al. (1978) used an analytical solution of the linearized equation set (Lighthill and Whitham 1955) to develop criteria for the applicability of kinematic and diffusion waves to open channel flow. They used sinusoidal perturbations to the mean flow, with L = wavelength; and T = wave period. For kinematic waves, the Ponce et al. criterion states that for the solution to be within 95% accuracy after one period of propagation, the dimensionless wave period has to be greater than 171. The dimensionless wave period τ is defined as follows:

T Souo               τ  =  ________             do

(7)

in which u_o = reference flow mean velocity; d_o = reference flow depth; and S_o = bottom slope. For practical applications, the wave period T can be taken as twice the time-of-rise of the flood wave.

For diffusion waves, a parameter τ/F_o, in which F_o = reference flow Froude number, is shown to be a better descriptor of the overall accuracy, accounting for both amplitude and phase errors. A practical applicability criterion for diffusion waves is the following (Ponce et al. 1978):

τ                      g      1/2 _____  =  TSo ( _____ )     ≥  30    Fo                    do

(8)

in which g = gravitational acceleration; and the other terms are as previously defined.

Ponce et al. (1978) confirmed the conclusions of Lighthill and Whitham (1955), as well as those of many others, that most overland flow situations would satisfy the kinematic wave criterion, and that most flood wave propagation cases in stream channels (excluding those with significant downstream control) would satisfy the diffusion wave criterion. Only in situations with markedly strong dissipative tendencies (for instance, a dam-break flood wave), flow into large reservoirs (with substantial backwater effects), or flow reversals would it be necessary to resort to the dynamic wave to properly describe the propagation of shallow-water waves.

7.  OVERLAND FLOW VERSUS UNIT HYDROGRAPH

The issue of whether an overland flow kinematic wave solution can replace (and perhaps eventually retire) the unit hydrograph as a practical method of runoff generation remains clouded in controversy. Due to the fundamental differences between these two methods, a resolution of this conflict does not appear to be forthcoming in the near future. The overland flow kinematic wave solution is a deterministic, distributed-parameter, hydraulic-data-intensive method (requiring geometric and frictional parameters), which is primarily applicable to small catchments [i.e., those less than 1 sq mi (2.5 km² )], for which the idealizations inherent in mathematical modeling can be justified on practical grounds. In other words, for the kinematic wave solution to be useful, the discretization must reflect what is actually occurring in the field. When used indiscriminately, without due regard for problem scale, there is a risk that the amount of lumping introduced may interfere with the deterministic character of the method and its ability to simulate overland flows in a distributed context.

In contrast to the overland flow kinematic wave solution, the unit hydrograph is a conceptual model of runoff generation, spatially lumped, and based exclusively on hydrologic data (streamflow measurements). Although originally derived for large catchments (Sherman 1932), the method has been found to have primary applicability to midsize catchments, i.e., those in excess of 1 sq mi (2.5 km² ) and less than 400 sq mi (1,000 km²). While these limits are somewhat arbitrary, they tend to reflect current hydrologic engineering practice. Furthermore, in the proper modeling context (i.e., with catchment subdivision), the applicability of the unit hydrograph can be extended to large catchments.

Since the overland flow kinematic wave solution is primarily applicable to small catchments, and the unit hydrograph is primarily applicable to midsize (and, by extension, to large) catchments, it seems that there should be little overlap between these two methods. In practice, however, existing computer models [for instance, the U.S. Army Corps of Engineers HEC-1 (HEC-1 1985)] provide users with a choice between these two methods to solve any given runoff problem, regardless of scale. This raises the question of which method is better, or more accurate, for a given problem, a question that has no clear answer. The methods are of such different nature and have such different data needs that they are not readily comparable. Perhaps the only defensible argument in this regard is that the kinematic wave solution should increase in accuracy as the catchment size decreases; and the unit hydrograph should increase in practicality as the catchment scale increases. Specific comparisons between the two methods are likely to lead to heated arguments, but the central issue of accuracy is not likely to be settled soon. For one thing, the overland flow kinematic wave solution is based on our currently imperfect knowledge of friction mechanisms, including the estimation of Manning's n and of the rating exponent m describing the mixed laminar-turbulent regime characterizing most overland flow situations. Likewise, the unit hydrograph would have to be verified with concurrent rainfall-runoff data, which are not readily available for the typical midsize catchment application.

Keeping in mind the question of scale, the kinematic wave solution does have the significant advantage that it can describe spatial and/or temporal rainfall and roughness variations, which the unit hydrograph method, by virtue of its being lumped, cannot do. Therefore, in situations where the scale question can be reasonably compromised, the overland flow solution should provide better detail in the simulation of flood flows, including the description of runoff concentration and diffusion. Therein lies the promise of kinematic waves and the expectation of significant improvements in the accuracy of runoff prediction.

As kinematic wave solutions continue to mature, particularly with the advent of a physically meaningful description of runoff diffusion, the way will be paved for the two methods of runoff generation to complement each other. There is an urgent need to develop synthetic unit hydrographs that reach beyond established practice [the Snyder unit hydrograph, to follow Corps of Engineers' practice; or the SCS dimensionless, to follow the Soil Conservation Service, (USDA: SCS 1985)]. Acting on this perceived need, the U.S. Bureau of Reclamation has developed a set of regional synthetic unit hydrographs for use within its jurisdiction (11 western United States) (USBR: Design 1987). In an attempt to overcome the shortcomings of conventional synthetic unit hydrographs, local agencies are engaged in developing synthetic unit hydrographs of the S-graph type (Sabol 1987, 1990). It is envisioned that under the proper modeling context, the overland flow kinematic wave model may be used as a tool to develop synthetic unit hydrographs without the burden of an extensive (and expensive) network of streamflow data collection [see Overton (1970)].

A precedent for the use of models to synthesize peak flows already exists in U.S. hydrologic practice: the TR-55 method (USDA: "Urban" 1986). This SCS method of peak flow generation was developed using the hydrologic catchment model TR-20 (USDA: "Computer" 1983) to generate synthetic peak flows that take into account the catchment's concentration properties, regional temporal rainfall distribution, and event-abstraction mechanisms of infiltration and depression storage. The TR-55 method improves on the rational method, substituting modeling for empiricism and leading to better and more consistent runoff predictions.

8.  KINEMATIC SHOCK

The kinematic shock was discussed in detail by Lighthitll and Whitham (1955); and since then, numerous studies have endeavored to analyze its causes and effects. But the subject continues to mystify researchers and practitioners alike (Cunge 1969; Kibler and Woolhiser 1970). The shock arises due to the nonlinear feature of kinematic waves, which under the right set of circumstances can result in the kinematic wave steepening to the point where it becomes for all practical purposes a wall of water. (In overland flow situations, the "wall of water" would be a small discontinuity in the water surface profile.) The shock is a direct consequence of the nonlinear steepening tendency, which is abetted when the following conditions occur (Ponce and Windingland 1985).

First, the wave is kinematic as opposed to diffusion (or dynamic). Diffusion is a mechanism acting to oppose the nonlinear steepening tendency. The more diffusive (or dynamic) a wave is, the less kinematic, and therefore, the less the steepening tendency. Second, there is a low base-to-peak flow ratio. The steepening tendency is promoted when the flow is subject to large relative changes, with baseflow being only a small fraction of peak flow. In the limit, as the baseflow approaches zero, i.e., in the case of an ephemeral stream, the steepening tendency (due to this condition) is greatest. This explains the fact that the shock is a relatively more frequent occurrence in the case of flash floods in arid and semiarid regions, which are caused by intense thunderstorms concentrating large flows into ephemeral streambeds. Third, there is a hydraulically wide and sufficiently long channel. Since the wave steepening is gradual, a sufficiently long channel is needed to give the shock a chance to develop. Strong steepening tendencies may require a shorter reach; a weak tendency may never develop the shock, given the complicating spatial effect of lateral inflows. A necessary condition for shock development is that the channel be hydraulically wide, that is, one in which the wetted perimeter is nearly constant. In practice, this hydraulic condition is reflected in a canyon-type situation, with essentially vertical walls in which as the flow increases, the wetted perimeter increases very little in comparison with the increase in flow depth and area. Mathematically, this condition is represented by a rating exponent β much greater than 1 (approaching β = 3/2 for Chezy friction; or β = 5/3 for Manning friction). Conversely, in shallow-overbank flow situations (with channels rapidly expanding their top width, such that the wetted perimeter increases in the same proportion as the flow area, and consequently, the hydraulic radius remains nearly constant), the steepening tendency is counteracted and shock development is arrested. Mathematically, this condition is represented by a value of the rating exponent β ≅ 1. Fourth, there is a high-Froude-number flow. In hydraulically wide channels, high-Froude-number flows lack sufficient diffusion to effectively counteract the steepening tendency (Ponce and Simons 1977). Therefore, at high Froude numbers, at or above critical, the steepening tendency may be promoted to the point where the shock may develop because of this condition. Theoretically, as the Froude number for turbulent flow approaches the condition of neutral stability (Fn = 1.5 for Manning friction; Fn = 2 for Chezy friction), the shock becomes unstable and the theory ceases to apply. In practice, however, such high-Froude-number flows are rare in natural steams and rivers, and the instability is seldom, if ever, observed. However, in artificial and other channels of relatively smooth boundary, free-surface instability of the type discussed here leads to the well-known roll waves (Dressler 1949; Chow 1959; Jolly and Yevjevich 1971). These roll waves are often seen in steep city streets and similar situations where prevailing laminar (or mixed laminar-turbulent) flow conditions may reduce the neutral-stability Froude number to values well below that corresponding to turbulent flow (to Fn ≅ 1 for mixed laminar-turbulent flow; and to Fn = 0.5 for laminar flow).