USGS Banner

Title: Two-Dimensional Hydrodynamic Model for the ACT and ACF River Basins

Focus Categories: MOD, FL, HYDROL

Keywords: Numerical Analysis, Rivers, Fluid Mechanics, Flood Control, Hydraulics.

Duration: 3/2000 to 2/2001

Federal funds requested: $17,996

Non-Federal matching funds: $40,278

Principal Investigator:

Dr. Fotis Sotiropoulos
School of Civil and Environmental Engineering
Georgia Institute of Technology, Atlanta, GWRI.

Congressional District: 5

Critical Regional Water Problems

The ongoing tri-state (Georgia-Alabama-Florida) negotiations over the allocation of water in the Apalachicola-Chattahoochee-Flint (ACF) and the Alabama-Coosa-Talapousa (ACT) Rivers Basins are presently at a critical stage. Even after, however, a water allocation agreement has been established, there will be a critical need to implement a river management system that ensures that the negotiated water allocation formula is observed by all riparians and that the rivers serve their water users the best way possible. Such an operational management system is not presently in existence. Though the primary emphasis of the tri-state negotiations is placed on droughts, the process also offers an opportunity to apply a modern and more effective flood management system. For the State of Georgia, this is not a luxury but an imperative need in view of the deaths and property damage caused by the killer floods of 1994 and 1996. These severe flooding events underscored the need for developing an operational river management system capable of accurate rainfall estimation, reliable flood forecasting, and effective overall water management. Such a system should be based on remote and on-site data and modern advances in rainfall-runoff, river routing, and reservoir management procedures. The proposed research addresses one critical aspect of this global modeling effort.

Expected Results and Benefits

The research proposed herein will lead to the development of an important component of a regional river management and flood forecasting system: an advanced, two-dimensional numerical model for simulating unsteady flows in natural river systems. The numerical model will be based on the depth-averaged open-channel flow equations formulated in generalized curvilinear coordinates, so that it can be readily applied to river reaches of arbitrary geometrical complexity. Special emphasis will be placed to develop a computationally efficient methodology, by designing the numerical algorithm so that it can take full advantage of state-of-the-art massively parallel computational platforms.

During the duration of this research, the numerical model will be calibrated for and applied to segments of both the ACT and ACF River Basins. Simulations will be carried out for various inflow hydrographs to understand the transient response of the two basins to discharge variations due to hydropower release events and/or flooding events.    The results of the 2D model will also be employed to evaluate and further refine the predictive capabilities of simpler one-dimensional river routing models, which are used in existing operational river management systems. Therefore, the proposed research is a critical prerequisite for developing a reliable and efficient state-of-the-science model for addressing the river management and flood forecasting needs of the State of Georgia. The proposed model will also provide the computational framework for developing a future more general numerical tool capable of tackling sediment transport and water quality issues in the ACT and ACF River Basins. For that reason the model will be formulated so that it can be readily extended in the future to solve the additional transport equations needed to tackle such issues.

Nature, Scope, and Objectives of the Research

The objective of the proposed research is to develop an advanced, two-dimensional numerical model for simulating unsteady flows in the ACF and ACT River Basins. Such a tool is presently necessary in order to: i) study and understand the complex response of the two basins to transient changes in discharge (such as those likely to occur during a hydropower release event and/or a flooding event); and ii) assess the performance and refine the predictive capabilities of simpler, one-dimensional river routing models used in existing operational river management systems. The work proposed herein will include the following tasks:

Task 1: Development of the Numerical Method. A numerical model will be developed for solving the unsteady, depth-averaged, open-channel flow equations in generalized curvilinear coordinates. The model will be formulated so that it can be readily applicable to complex, natural river systems. Special care will be exercised to ensure that the numerical methodology can take full advantage of massively parallel computational platforms so that it can serve as a practical engineering tool.

Task 2: Data Collection and Grid Generation for the ACT and ACF River Basins. The purpose of this task is two-fold. First, we will contact the Army Corps of Engineers to obtain topography and any available flow data for the ACT and ACF River Basins. Second, we will employ the topography data to generate computational grids that can adequately resolve all the important features of the two basins.

Task 3: Simulations of the ACT and ACF River Basins. The numerical model will be applied to model flows in both basins. Available data sets will be used to calibrate and validate the model. Simulations will be carried out for various inflow hydrographs, designed to simulate hydropower release events and/or flooding events, to understand the transient response of the two basins to discharge changes. The effects of grid and time-step resolution on the accuracy of the computed solutions will be carefully evaluated by carrying out extensive numerical sensitivity tests.

Task 4: Assessment and Refinement of 1-D River Routing Models. The results obtained with the 2D numerical model will be compared with those obtained from 1D models, such as the widely used Muskingum-Cunge model, to identify specific weakness of these simpler tools and recommend possible improvements.

Methods and Procedures

In the following sections, we present the depth-averaged open-channel flow equations, discuss the issue of turbulence closure, describe the numerical method we will employ to solve the governing equations, and outline the grid generation procedures we will employ to generate computational meshes for the ACT and ACF River Basins.

12.1 Governing Equations

The numerical model will rely on the unsteady, depth-averaged, open-channel flow equations. These equations are obtained by integrating the three-dimensional, incompressible, Reynolds-averaged Navier-Stokes (RANS) equations from the river bed to the water surface under the following assumptions: i) uniform velocity distribution in the depth direction; ii) hydrostatic pressure distribution; iii) small-channel bottom slope; iv) negligible wind shear at the surface; and v) negligible Coriolis acceleration. Although this averaging procedure simplifies the equations considerably, information concerning the vertical variations of the velocity distribution is lost. Nevertheless, the depth-averaged equations have been shown to yield very good predictions for flows in shallow channels and provide a far more practical, from the computational standpoint, engineering alternative to the full, three-dimensional RANS equations. We should note that we propose to depth-average the three-dimensional RANS equations, rather than the instantaneous Navier-Stokes, because we would like the model to be applicable to turbulent open-channel flows. The issue of turbulence closure for these equations is discussed in the subsequent section.

The depth-averaged open-channel flow equations (continuity and momentum) formulated in Cartesian coordinates read as follows (see Kuipers and Vreugdenhil (1973) for the details of the integration):

Continuity equation:

Continuity Equation                                                                                       (1.a)


X-Momentum Equation(1.b)


Y-Momentum Equation(1.c)

In the above equations h is the water depth, U and V are the mean depth-averaged velocity components in the x- and y-directions, respectively, g is the acceleration of gravity, zb is the channel bottom elevation, Tb is the bottom shear stress, and txx, txy, and tyy are mean depth-averaged effective stresses

Equation                                              (2.a)

Equation                              (2.b)

Equation                                                 (2.c)

where u and v are the time-averaged velocity components, u’ and v’ are the fluctuating components, and m is the viscosity of the water. In the above integrals, the first terms represent the viscous stresses, the second terms are the components of the Reynolds-stress tensor (arising from time-averaging the instantaneous Navier-Stokes equations), while the last terms are the so-called dispersion stresses (see subsequent section).

The bottom shear stresses can be calculated by using the following formulas:

Bottom Shear Stress Formula 1 ,   Bottom Shear Stress Formula 2                                         (3)

where Cf is the friction coefficient that needs to be determined as part of the model calibration (see related section below).

The above equations are presented, for the sake of simplicity, in Cartesian coordinates. In the proposed numerical model, however, these equations will be transformed into generalized curvilinear coordinates by invoking the partial transformation approach Equation, which maintains the Cartesian velocity components as the dependent variables. This step is essential in order to ensure that the model can simulate accurately arbitrarily shaped natural River channels. Due to space limitations the transformed equations are not included herein (see Molls and Chaudry (1995) for details).

12. 2 Turbulence Closure

Closure of the depth-averaged equations requires the determination of the components of the effective stress tensor. Since our emphasis is on high Reynolds-number flows, we will neglect the viscous stresses from eqns. (2). We will also neglect the depressive stresses as these are difficult to quantify and model. It is common practice in the literature to neglect these stresses, a simplification which does appear to have a significant impact on the accuracy of the model. The remaining depth-averaged Reynolds stresses will be modeled using the Boussinesq eddy-viscosity concept:

Equation                                           (4)

where i,j = 1,2, nt is the depth-averaged eddy-viscosity, k is the depth-averaged turbulence kinetic energy, and dij is Kronecker’s delta.

There are several alternatives to obtain the eddy-viscosity needed to close the governing equations. The simplest approach is to assume a constant eddy-viscosity along the entire basin. The value of this constant can be determined by calibrating the model equations with available data. Although this treatment is computationally very attractive, due to its simplicity and computational expedience, it constitutes an oversimplifying treatment of the very complex turbulent transport processes within natural rivers. A far more rigorous approach, albeit somewhat more expensive computationally, is to determine the eddy-viscosity distribution by solving two additional transport equations for the turbulence kinetic energy, k, and its rate of dissipation e, the so-called k-e model. These equations provide turbulence velocity and length scales that can be used to calculate the eddy-viscosity as follows:

Equation                                                                    (5)

where Cm is a model constant (Cm = 0.09). The depth-averaged k-e model equations can be found in Rastogi and Rodi (1978) and McGuirk and Rodi (1978). Other closure alternatives are available (see Rodi (1980) for a comprehensive review of turbulence models for hydraulic engineering problems) but their sophistication and associated computational overhead are beyond the scope of the present modeling effort, which will focus on basin-wide flow simulations.

In the proposed research, we will evaluate both turbulence closure options. Implementing the constant eddy-viscosity model is, of course, rather trivial but will, as already indicated, require considerable calibration with field data. We will also implement the k-e model and compare its performance relatively to the constant-eddy viscosity treatment.

A very important issue when modeling turbulent flows is the approach adopted to resolve the flow near solid boundaries (the channel banks, in the present case). Given the scales the proposed model is intended for, the only practical approach for modeling the near-wall flow is to employ the so-called wall-functions approach. This approach makes use of the logarithmic law of the wall to bridge the gap between the fully turbulent region and the laminar sublayer and allows the simulation of flows at very high Reynolds numbers with rather coarse computational meshes. An important advantage of the wall functions approach is that it can be readily extended, using appropriately modified expressions for the logarithmic velocity profile, to account for roughness effects, which are very important and need to be accounted for in natural rivers (see Rodi (1980) for details).

12.3 The Numerical Method

The depth-averaged continuity, momentum, and turbulence closure equations will be formulated in strong-conservation form and discretized using the finite-volume approach. Second-order accurate central-differencing will be employed to discretize all spatial derivatives. Due to the purely dispersive errors introduced by the central-differencing of the convective terms, artificial dissipation terms need to be explicitly introduced for stability. We will employ third-order, fourth-difference artificial dissipation terms scaled appropriately with the eigenvalues of the Jacobian matrices of the convective flux vectors (see Lin and Sotiropoulos 1997, for a three-dimensional implementation of various artificial dissipation schemes for incompressible flow equations). To facilitate the resolution of regions of steep gradients of the water depth, as well as to ensure that the proposed model can also be applied to simulate traveling hydraulic jumps, we will implement a non-linear artificial dissipation model by extending to depth-averaged flows techniques developed for capturing shock waves in compressible flows (see Jameson et al. 1981). More specifically, the artificial dissipation terms will consist of both first and third-order terms, each multiplied by appropriate functions designed to switch-on or off certain terms based on the local intensity of the spatial gradients of the water elevation. Thus, in regions of steep water-surface gradients the first-order terms are automatically turned on by the model to preserve the robustness and monotonicity of the numerical model. In smooth regions of the flow, however, only the third-order terms are activated, thus, preserving the second-order accuracy of the numerical method.

The discrete governing equations will be integrated in time using an explicit, multi-stage, Runge-Kutta algorithm designed to ensure second-order temporal accuracy. In spite of the fact that explicit integration algorithms typically impose more severe time-step limitations, as compared to fully implicit procedures, they are preferable because they can take full advantage of parallel computer architectures. Such algorithms do not involve matrix inversion operations and can, thus, be parallelized very effectively and in a rather straightforward manner. The resulting computer code should be able to run very efficiently on multi-processor computer architectures and, thus, the additional cost due to the increased number of time steps required to simulate a given time interval (because of the stability restriction of the explicit algorithm) can be greatly offset by the drastic reduction of CPU time per time step.

12.4 Grid Generation Procedure

The ACT and ACF River Basins are both very complex systems, each consisting of three major rivers and several tributaries. The use of curvilinear grids is, therefore, essential in order to ensure that the various topographical features of these systems are adequately resolved. In the proposed research, we will model directly only the three major rivers of each system. The effects of the tributaries on the flow will be accounted for by specifying inflow or outflow-type boundary conditions at the locations where the various tributaries intersect the main river.

Using available topography data (from the Army Corp of Engineers), we will re-construct the geometry of the riverbanks and the variation of the bed elevation, zb=zb(x,y), that appears in the depth-averaged equations. Curvilinear computational grids will be constructed using the GRIDGEN commercial grid generation software, available at the Computational Fluid Dynamics laboratory of CEE at Georgia Tech. Spline interpolation will be used to obtain the values of zb at the nodes of the computational grid.

12.5 Model Calibration

In general, there are two model constants that need to be calibrated with input from field observations: the eddy viscosity coefficient nt, and the bed friction coefficient Cf. The former coefficient requires calibration only when turbulence closure is achieved by assuming a constant eddy-viscosity. When the k-e turbulence model is employed, on the other hand, nt is directly calculated from the model equations and does not need to specified as a model constant. As we have already indicated, in the proposed research we intent to investigate both turbulence closures. For the constant eddy-viscosity model, we will compute nt by adjusting its level empirically until reasonable agreement with field observations is achieved. The predictions of the k-e model may also be used to provide a global estimate for the level of nt.

The bed friction coefficient results from the depth-averaging of the governing equations and can only be determined by calibrating the numerical model. Following Wenka et al. (1991), Cf is related to the friction coefficient l (Cf=l/8), which can be calculated from the Prandtl-Colebrook formula in terms of the Reynolds number (based on channel depth), the local channel depth, and the equivalent roughness height. Since the roughness distribution in natural River systems is unknown (and very difficult to estimate experimentally without laboratory scale studies), we propose to calculate the equivalent roughness height from the Manning friction coefficient formula. The Manning coefficient depends strongly on the specific river situation and needs to be determined by calibrating the model using available data. The procedure we will employ is to adjust Manning’s coefficient until the calculated slope of the water-surface profile for a given discharge matches closely that observed in the field.

Related Research

Depth-averaged models have been employed in the past to simulate flows in man-made open channels as well as in natural environments, such as river systems and estuaries. Since the emphasis of this proposal is on natural river basins, we present here a brief review of those previous studies that focused on shallow flows in natural environments. It should be pointed out that what follows is not intended to be an extensive review of related literature. We rather focus on few studies that underscore the progress made in modeling open channel flows with 2-D models during the past decade.

Among the most notable early studies is the work of Tingsanchali and Rodi (1986) who developed a depth-averaged k-e model for simulating suspended sediment transport in the Neckar River in Germany. Wenka et al. (1991) employed a similar k-e model to simulate flooding over a stretch of the river Rine in Germany. They reported good agreement between the predicted flow patterns and those observed in aerial photographs of the river during the flooding event.    Muin and Spaulding (1996) developed a two-dimensionsal model for simulating tidal circulation patterns in the Providence River. Their model neglects all viscous and turbulent stresses, except, of course, the stresses due to bed friction. Hu and Kot (1997) developed a depth-averaged model for studying tides in the Pearl River estuary in China. They employed the constant eddy-viscosity assumption for turbulence closure. Shankar and Cheong (1997) employed various formulations of a 2-D numerical model to study tidal motions in Singapore coastal waters. They also adopted the constant eddy-viscosity assumption.

An important conclusion that follows from the above brief review is that the enormous complexities of natural flows necessitate the development of numerical models that are specifically tailored to the idiosyncrasies of a specific site. To the best of our knowledge a two-dimensional numerical model for the ACT and ACF River Basins, such as the one proposed here, has neither been developed nor is currently under development elsewhere. Moreover, and with only exception the works of Tingsanchali and Rodi (1986) and Wenka et al. (1991), we are not aware of any other depth-averaged model that has employed a rigorous turbulence closure model, such as the k-e model proposed herein, to simulate flows in natural rivers (recall that most previous studies have assumed constant eddy viscosity).


Hu, S., and Kot, S. C. (1997), “Numerical Model of Tides in pearl River Estuary with Moving Boundary,” ASCE J. Hydr. Eng. 123(1), pp. 21-29.

Jameson, A., Schmidt, W., and Turkel, E. (1981), “Numerical Simulation of the Euler Equations by Finite Volume Methods using Runge-Kutta Time Stepping Schemes,” AIAA paper 81-1259.

Kuipers, J., and Vreugdenhil, C. B. (1973), “Calculation of Two-Dimensional Horizontal Flow,” Report S163, Part I, Delft Hydraulic Laboratory, Delft, The Netherlands.

Lin, F., and Sotiropoulos, (1997) F., "Assessment of Artificial Dissipation Models for 3-D Incompressible Flow Solutions," ASME J. Fluids Eng. 119, 331-340.

McGuirc, J. J., and Rodi, W. (1978), “A Depth-Averaged Mathematical Model for the Near Field of Side Discharges into Open-Channel Flow,” J. Fluid Mech. 86(4), pp. 761-781.

Molls, T., and Chaudhry, M. H. (1995), “Depth-Averaged Open-Channel Flow Model,” ASCE J. Hydr. Eng. 121(6), pp. 453-465.

Muslim, M., and Spaulding, M. (1996), “Two-Dimensional Boundary-Fitted Circulation Model in Spherical Coordinates,” ASCE J. Hydr. Eng. 122(9), pp. 512-521.

Rastogi, A. K., and Rodi, W. (1978), “Predictions of Heat and Mass Transfer in Open Channels,” ASCE J. Hydr. Div. 104(3), pp. 397-419.

Rodi, W. (1980), “Turbulence Models and Their Application in Hydraulics,” IAHR monograph, Delft, The Netherlands.

Shankar, J., and Cheong, H.-F. (1996), “Boundary Fitted Grid Models for Tidal Motions in Singapore Coastal Waters,” IAHR J. Hydr. Res. 35(1), pp. 3-19.

Tingsanchali, T., and Rodi, W. (1986), “Depth-Averaged Calculation of Suspended Sediment Transport in Rivers,” 3rd Int. Symposium on River Sendimentation, The Univ. of Mississipi, March 31-April 4.

Wenka, T., Valenta, P., and Rodi, W. (1991), “Depth-Averaged Calculation of Flood Flow in a River with Irregular Geometry,” proc. of XXIV IAHR Congress, Madrid, Spain, Sept. 9-13, pp. A-225 to A-232.

Deliverables and Information Transfer Plan

The deliverables of this research will be: i) a computer code for simulating unsteady flows in the ACT and ACF River Basins using depth-averaged equations; ii) a user’s manual describing the numerical method and the computer code; iii) a paper submitted to a refereed journal; and iv) a proposal submitted to an appropriate state agency requesting funds for further development of the numerical model. The deliverables will be made available to the Georgia State Department of Natural Resources.

U.S. Department of the Interior, U.S. Geological Survey
Maintained by: John Schefter
Last Updated: Wednesday October 26, 2005 1:09 PM
Privacy Statement || Disclaimer
|| Accessibility