- Email: [email protected]
A suspended sediment flow model for high solids concentration using higher order turbulence closure
A suspended sediment flow model for high solids concentration using higher order turbulence closure B. A. DeVantier and R. Narayanaswamy Department of Civil Engineering and Mechanics, Southern Illinois University, Carbondale, IL62901-6603, USA
A predictive model of sediment laden flows using a second order closure is presented. Damping of turbulence due to the presence of suspended solids is included in the model through consideration of local damping of turbulence kinetic energy. A modified k-e model is used which utilizes damping terms to account for both sediment-induced density stratification and particle interference with eddy motions. Comparisons of model predictions for solids and velocity vertical distributions to measurement are presented for both negatively buoyant and neutrally buoyant solids with good agreement. Free surface damping of turbulence is also included in the model, and the effect upon the suspended solids concentration of this damping is demonstrated. Key Words: turbulence damping, suspended sediment, turbulence modelling.
INTRODUCTION Many models have been proposed to predict the vertical distributions of velocity and solids in sediment laden flows. Recent concerns of the transport of pollutants bound to suspended sediment I have sparked renewed interest in the improvement of the predictive capabilities such models. Models which have the capability to accurately predict turbulent transport effects for fluid momentum, suspended sediment, and dissolved contaminants, all at the same time, are highly desirable. Among the earliest models of turbulent transport in flows with suspended sediment is the mixing length model and variations on it z-~. These models utilized various forms of the universal log velocity distribution. While it has long been noted that high solids concentration can have a significant effect on this 'universal' distribution 5, there is disagreement as to the universality of the von Karman constant for the log velocity expression. Still, it is well accepted that high sediment levels can significantly modify velocity profiles and can even affect the stressstrain response of the fluid mixture to the point of nonNewtonian behaviour. While non-Newtonian effects can be masked when flow is fully turbulent as is often the case in sediment transport, the difficulty of accurately predicting the effects of sediment upon turbulent transport processes is central to the disagreement concerning universal profiles. Dissatisfaction with the adequacy of the log velocity description resulting from a simple mixing length approach has led to attempts to formulate a more general description of momentum transfer in channel flows. Parker and Coleman have investigated a wake function approach to flow description 6. Another approach is to use higher order turbulence transport equations to shift the modelling emphasis away from the level of the Paper acceptedJanuary 1988. DiscussionclosesAugust 1989. © ComputationalMechanicsPublications1989 46
Adv. Water Resources, 1989, Volume 12, March
momentum description to that of the turbulence equations. The k-e model of Launder and Spaulding 7 is one of the simplest and probably the most commonly used higher order model of turbulence. It has the advantage of 'more universal' empirical modelling constants, and thus is more likely to provide the necessary detail for an acceptable description of the effects of suspended sediment on flow velocity and vice versa. THEORY I. Basic model The k-e model is appealing for application to sediment transport description, however the basic model must be modified to describe the effects of the suspended solids. A rigorous discussion of the background and development of such a flow model has been given by Devantier and Larock 8, and will provide the basis of further discussion. The momentum equations and turbulence closure equations are similar to those presented by Rodi 9 for the k-e model for density stratified flow. For steady, vertically varying, fully developed open channel flow, the governing equations of the model may be written as
0=
1 dP Pb dx
d(uw) dx
--dU
(l)
d [v, dk]
O=-UW d - + zk
--
zj-a wg
o = < fi -uw
(2)
- < W+
zkgyzj (3)
0 = d {A(1 - A)vs} - d
{7~ww}
(4)
A suspended s e d i m e n t f l o w model: B. A . D e V a n t i e r and R . N a r a y a n a s w a m y
where flow is in the x direction, the primary variables U (horizontal mean velocity), k (turbulence kinetic energy), e (dissipation rate of k), and A (mean volume concentration of solids) vary only in the vertical direction z. Also 6 is the specific gravity difference between the sediment and water (1.65 for sand). The k ~ model constants c, (0.09), CE1 (1.51), Ce2 (1.92), Cr, 3 (0.8), O"k (1.0), and tr~ (1.3) are those specified by Rodi 9, and the overbarred terms are time averaged turbulence flux terms which will require further modelling to put them in terms of the primary variables of the model. The lower case variables in these flux terms are the time fluctuating counterparts of upper case mean quantities resulting from the classical Reynolds decomposition of time varying quantities. Properties related to the flow and transport are the bulk density of the flow suspension Pb, the net sedimentation velocity of the suspension v~, and the pressure drop driving the channel flow - d P / d x . Closure for equations (1) and (4) is accomplished by fluxgradient relations given as --dU u w = - v, dz
(5)
~v~ dA t0w . . . .
(6)
(~A dz
where v t is the eddy viscosity, and trA is a sediment Schmidt number near unity (taken as 1.2 (7)). The k-e model determines vt in terms of its primary variables as
The shear velocity for equations (10)-(12) is determined by iteratively solving the following logAaw velocity equation for u, UI.o=~ ln(K ?)
where e is a wall roughness height and K is taken to be 30 as for sand roughened pipes, and U is known from the previous iteration's estimate of the velocity field. The approach to the free surface boundary conditions is to treat it as a rigid lid and apply zero flux boundary conditions for all variables. I I . M o d e l variations
The most significant modification to be considered is a change in the way that the model accounts for damping of the turbulent motion due to the suspended solids. For the highest levels of local bed values for A, the basic model underestimates the reduction of turbulence 8. If turbulence damping is due not only to energy expended to suspend sediment, but is also due to reduction of eddy length because of the simple presence of solids, then the additional damping of turbulence could be accounted by adding a term B A to the density stratification damping term B. The mathematical form of BA is determined by drawing an analogy to the Monin-Oboukhov expression for reduction of mixing length a s 9 l/lo= 1 - c t Ri
k2 vt = - -
(7)
Two significant terms which appear several times in equations (2) and (3) are --dU Pr= -uw dz
(8)
B = - 32wg
(9)
and
(14)
where lo is the mixing length with no stratification, and R i is the Richardson number defined as Ri = - B/Pr
(15)
and 7 is an empirical constant between 3 and 10 for stably stratified atmospheric boundary layer flows. The simplest algebraic form relating mixing length to local sediment concentration would be a linear relation of the form I/l o =
is the mean shear induced production of k, and B is the buoyancy-induced production. B acts as a sink of k when the suspended solids are negatively buoyant as in the case of suspended sediment. With equations (5)-(6) substituted into equations (1)-(4) the equation set is closed, and the model may be applied if boundary conditions are supplied. The near-bed boundary conditions are specified at a slight distance Zo away from the wall. The boundary condition for velocity is actually applied as a derivative condition of the form
(13)
1-
CAA
(16)
Pr
zo/Pb = u2 = - vt dU/dz
to
(10)
where u, is the shear velocity (as defined in terms of the near-bed shear stress %), and ~c is the yon Karman constant taken as 0.42. The standard application ofk and e boundary conditions at this boundary are
where CA is a proportionality constant. This form is chosen mainly for simplicity, however there are also some reasonable physical arguments for doing so. Since the effect of B A should not become significant until A is large, it may be argued that the reduced mixing length is directly related to the reduction of volume occupied by fluid. Thus the damping is in all directions, rather than only in the direction opposing gravity as in the case of B. If the Richardson number of this volume concentration damping is taken as Rih =- -- B A / P r
then by combining equation (14) (with RiA replacing Ri) with equations (16) and (17), B A c a n be expressed as A BA=
kL =
(11) 3
e]~° = u , / ( ~ Z o )
(12)
(17)
----
CA
Pr
(18)
where ~ has been combined into CAas part of the empirical fitting constant. B h is then incorporated into the model by
A d v . Water Resources, 1989, Volume 12, M a r c h
47
substituting a combined term B' everywhere for B, where the combined term is simply
B'=B+BA
(19)
The other modifications to the model are the incorporation of damping into the boundary conditions. The application or the velocity boundary condition can be modified to include Richardson number effects so that equation (13) becomes (20) The value of RiA evaluated at z0 can be taken as zo/L, when L is taken to be the Monin-Oboukhov length of the turbulence damping due to the presence of solids. Equation (20) is most commonly expressed using zo/L, but RiA can be obtained directly in the model formulation. Itakura and Kishi 1° have successfully applied the MoninOboukhov length scale in the zo/L form for the description of sediment-induced density effects on open channel flow. For a complete discussion of the MoninOboukhov length scale, the reader is directed to Turner 11. The Richardson number (and therefore L) is negative for the turbulence damping, so that shear stress near the bed is less than that of a clear water flow. Thus the use of equation (20) will result in a slight decrease in k and e at the near-bed boundary, because the inverse solution of the equation for u, will result in smaller values for u,. The other modified boundary condition accounts for damping of turbulence near the free surface as presented by Celik and Rodi ~1
e,[H-
k3/2/H
obtain convergence for the nonlinear equation set. The equation set (equations (1)-(4)) was solved in the time steady form rather than computing transient long time solutions. For the steady, unidirectional open channel flow problem chosen there are no convective transport terms, and therefore the near-bed values of U and A could be assumed without presupposing the turbulent diffusion effects in the vertical. Thus the value of A was simply set to the measured value and U was similarly set by adjusting the roughness value e. While the importance of the ability to predict these values is not to be minimized, the emphasis in here is upon verifying effects relating to vertical turbulent diffusion of momentum and solids, and this can be accomplished more clearly if inaccuracies in the specification of these boundary conditions are not allowed to interfere. The first model test examines the effect of the free surface damping of k and e for low levels of A so that damping due to B or BA should be negligible. Fig. 1 compares predictions to cases V 14 and V 17 of the data of Vanoni 2, and Fig. 2 is run 37 from the data of Montes and Ippen 13. The dashed curves are the predictions without free surface damping and the solid curves are predictions with damping. Inclusion of the effect clearly improves the prediction, and a more S-shaped curve similar to that of equation (23) results. However the model predictions are more physically realistic than that of equation (23) in the sense that A goes to a small nonzero value rather than exactly zero.
1.0
•
I
!
•
/
t
! I
(21) t
where H is the free surface above the bed, k s is the k at z = H, and cls has been empirically determined to be 5.54. With the applications of this boundary condition vt in clear water approaches the familiar parabolic profile
I
/
•
•
Z
V 14 Z
vt tcu, ~ (H - z) =
(22)
although equation (21) produces a small, nonzero surface value for vt rather than the zero level predicted by equation (22). Similarly, when equation (22) is used to estimate the A distribution, the result is the S-shaped curve described by
A=~° =
z0j
(H - z o)
0.5 •
/
(23)
where/3 is related to the ratio of u, and vs. Thus the model should also be better to reproduce free surface concentration distribution with this modification. RESULTS The model was tested by comparing model predictions for velocity and volume concentration to measured laboratory flume flow values at steady conditions. The model equations were solved by a finite element computer algorithm incorporating Newton-Raphson iteration to
48
Adv. Water Resources, 1989, Volume 12, March
!
I
I
2 A/A a Fig. 1. Model predictions and data of sediment volume fraction versus depth
Vanoni for
A suspended sediment flow model: B. A. DeVantier and R. Narayanaswamy o O-
\ 0 o0-
Predicted Meosured
\
~
,
o--r-
o-
•
m
di o:
0
°e,~
•
c5 l
I
l
I
l
0.001
!
I
!
]
l
•!
0,01
A Fig. 2.
Model predictions and data of Montes and Ippen for sediment volume fraction versus depth
The S-shaped curves can be characterized as having 3 regions. The free surface and the bed act similarly in limiting eddy excursions toward the boundaries. This allows greater net settling in these regions than in the centre portion where mixing is best and the concentration drop is nearly log-linear with depth. The near-bed and centre region effects have been noted by McTigue t4 in the comparison with the Montes and Ippen data, however all data sets of Vanoni exhibit the three distinct regions. The difference is perhaps that the McTigue comparison is for a relatively fine, uniform sediment of smaller settling velocity. The free surface drop-off may well have been occurring, but experimental error may have masked the effect. The effect of the physical damping through BA was examined by comparing the model predictions to measurements of Einstein and Chien 4 for high levels of sediment volume fraction. The fitting constant CA was determined by varying it to obtain the 'best fit' to case S16. The plot of the velocity prediction is shown in Fig. 3. The three curves compare the model prediction with no density stratification, noted as a 'Clear' water, with buoyancy damping due to density stratification alone, noted as B, and with both B and BA effects, noted as B'. The progressive increase in velocity with increased damping effect indicates that density effects through B alone are not sufficient to properly describe the damping of momentum transport. Fig. 4 presents the clear water and sediment laden (B') curves on a semi-log plot. The
straight lines superimposed upon the curves are given to demonstrate the similarity to a plot presented by Vanoni and Nomicos 5 which showed the drag reduction effect of high suspended load. Fig. 5 compares the model prediction of A to measurements for case S16. The prediction deviates significantly from measurements above 0.2 flow depths. The poor prediction of A does not have a significant effect upon prediction of U however, because the volume concentration has dropped to low enough levels by 0.2 depths to greatly reduce the modified effect upon momentum transfer. The discrepancy is probably due to the model's simplified assumption of a single average settling velocity vs, while in the measurement there was a distribution of sediment sizes for S16. The boundary condition for k (equation (20)) was tested for the effect of the sediment modified shear velocity of equation (20). Fig. 6 is a comparison of the k predictions with equation (20) versus equation (13) for case S16 of Einstein and Chien. Use of equation (20) modified the values of k somewhat, but an interesting effect is that the lowered level ofk near the bed results in a higher level of k through most of the rest of the flow. The effect upon e of the change was similar to k and was such that it compensated to keep the vt nearly the same for either method, so that the change had practically no effect upon predictions of U and A. Thus, it appears that either form is adequate for predictive purposes. To further test the adequacy of the inclusion of solids
Adv. Water Resources, 1989, Volume 12, March
49
A suspended sediment flow model." B. A. DeVantier and R. Narayanaswamy
J Clea
1.0
I
1.0
m
I t
07 /
I I
z
Z
m
H 0.1
0.5 'B'
-// / S J6
SI6 I
I
4
I
i
I
6 U
|
8
I
I
i
4
I
I
|
6
8
U ft/s
ft/s
Fig. 3. Model predictions of velocity for Einstein and Chien case S16 with clear water, with density damping, and with density and volume fraction damping
Fig. 4. Semi-log slope plots of model predictions for Einstein and Chien case S16 for clear water, and with density and volume fraction damping
1.0
Predicted -Meosured •
Z
-
, ,/
,, , / ~ w . ~
~
H
SI6 0.5
......
.
I
I
l~161
162
,,
I
ld 3
A Fig. 5.
50
Model predictions compared to volume fraction data of Einstein and Chien case S16
Adv. Water Resources, 1989, Volume 12, March
,,J
A suspended sediment flow model: B. A. DeVantier and R. Narayanaswamy qt-
'N•,,
~
damping of turbulence through BA,the model predictions for flow with neutrally buoyant solids was tested. Elata and Ippen 15 measured the reduction in K in flows with neutrally buoyant particles. The model was used to predict the effective reduction of ~cby taking straight line semi-log slopes from plots similar to Fig. 3. Fig. 7 compares the model predictions to the data of Elata and Ippen. The agreement is quite good. It appears then, that the addition of the term B6 is justified, and the fitted value of CAis adequate.
CLEAR WAFER SHEAR VELOCITY SEDIMENT CORRECITr.D SHEAR VELOCITY
",
..... •
%'~%%
0 Z,~.. 0
SUMMARY The damping of turbulence through the presence of suspended solids and the presence of a free surface has been incorporated into a higher order turbulence closure model for sediment laden flow. The addition of the effect of solids damping through BA appears to be significant, however it should only be important at high solids concentrations. The model appears to have an advantage over to the wake model of Parker and Coleman 6 in predicting the effect of sediment upon velocity, because the wake model effects are based upon a density difference between solids and water. The reduction of sediment level
% 0
O o
I
I
I
I
I
I
I
10 -3 DIMENSIONLESS KINETIC ENERGY ( k / U 2j
I
l
10
-'
Fig. 6. Turbulence kinetic energy prediction comparison using two forms for calculatin9 shear velocity
1.00 -
\ Predicted Meosured
0.90 -
% 0.80 ! -
ib II
0.70
I[, il,
N
0.60
0.50 0.00
,
|
J
I
I
I
J
I
I
|
I
I
I
I
I
I
0.10
I
I
l
]
0.20
I
I
I
I
I
I
I
I
I
[
0.30
A Fig. 7. Comparison of model prediction of reduction of effective yon Karman constant to data of Elata and Ippen for neutrally buoyant solids
Adv. Water Resources, 1989, Volume 12, March
51
A suspended sediment flow model: B. A. DeVantier and R. Narayanaswamy d u e to free s u r f a c e d a m p i n g in t h e m o d e l s h o u l d h a v e little effect u p o n velocities, h o w e v e r s e d i m e n t level p r e d i c t i o n s are m a r g i n a l l y i m p r o v e d .
7
REFERENCES
9
1 2 3 4
5 6
52
Jaffe, P. R., Parker, F. L. and Wilson, D. J. Distribution of Toxic Substances in Rivers, J. Envir. En9r. Div., ASCE, August 1982, 108, 636-649 Vanoni, V. A. Transportation of Suspended Sediment by Water, Trans. ASCE, 1946, 111, 67 133 Zeller, J. Tonerosion, Eidgen. Techn. Hoschule, Interner Bericht, T2, T3, 1963 Einstein, H. A. and Chien, N. Effects of Heavy Sediment Near the Bed on Velocity and Sediment Distribution, U. of Cal., Berkeley, and US Army Corps of Engineers, Missouri Div., 1955, Report No. 8 Vanoni, V. A. and Nomicos, N. Resistance Properties of Sediment-Laden Streams. Trans. ASCE, 1960, 125, 114(~1152 Parker, G. and Coleman, N. L. Simple Model of SedimentLaden Flow, H. Hyd. Div., ASCE, May 1986, 35(~375
Adv. Water Resources, 1989, Volume 12, March
8
10 11 12 13 14 15
Launder, B. E. and Spaulding, D. B. The Numerical Computation of Turbulent Flows, Computational Meth. Appl. Mech. Engng., 1974, 3, 269-289 DeVantier, B. A. and Larock, B. E. Sediment Transport in Stratified Turbulent Flow, J. Hyd. Div., ASCE, December 1983, 109, 1622-1635 Rodi, W. Turbulence Models and their Application in Hydraulics, IAHR publication, 1980 Itakura, T. and Kishi, T. Open Channel Flow with Suspended Sediments, J. Hyd. Div., ASCE, August 1980, 106, 1325-1342 Turner, J. S. Buoyancy Effects in Fluids, Cambridge Univ. Press, 1973 Celik, I. and Rodi, W. Simulation of Free-Surface Effects in Turbulent Channel Floe, Physiochemical Hydrodynamics, 1984, 5, 217 227 Montes, J. S. and Ippen, A. T. Interaction of Two-Dimensional Flow Suspended Particles, R.M. Parsons Lab. Water Res. and Hydrodynamics, Mass. Inst. Tech., 1973, Report 164 McTigue, D. F. Mixture Theory for Suspended Sediment Transport, J. Hyd. Div., ASCE, June 1981, 107, 659 673 Elata, C. and Ippen, A. T. The Dynamics of Open Channel Flow with Suspensions of Neutrally Buoyant Particles, Hydrodynamics Lab., Mass. Inst. Tech., 1961, Report 45
A predictive model of sediment laden flows using a second order closure is presented. Damping of turbulence due to the presence of suspended solids is included in the model through consideration of local damping of turbulence kinetic energy. A modified k-e model is used which utilizes damping terms to account for both sediment-induced density stratification and particle interference with eddy motions. Comparisons of model predictions for solids and velocity vertical distributions to measurement are presented for both negatively buoyant and neutrally buoyant solids with good agreement. Free surface damping of turbulence is also included in the model, and the effect upon the suspended solids concentration of this damping is demonstrated. Key Words: turbulence damping, suspended sediment, turbulence modelling.
INTRODUCTION Many models have been proposed to predict the vertical distributions of velocity and solids in sediment laden flows. Recent concerns of the transport of pollutants bound to suspended sediment I have sparked renewed interest in the improvement of the predictive capabilities such models. Models which have the capability to accurately predict turbulent transport effects for fluid momentum, suspended sediment, and dissolved contaminants, all at the same time, are highly desirable. Among the earliest models of turbulent transport in flows with suspended sediment is the mixing length model and variations on it z-~. These models utilized various forms of the universal log velocity distribution. While it has long been noted that high solids concentration can have a significant effect on this 'universal' distribution 5, there is disagreement as to the universality of the von Karman constant for the log velocity expression. Still, it is well accepted that high sediment levels can significantly modify velocity profiles and can even affect the stressstrain response of the fluid mixture to the point of nonNewtonian behaviour. While non-Newtonian effects can be masked when flow is fully turbulent as is often the case in sediment transport, the difficulty of accurately predicting the effects of sediment upon turbulent transport processes is central to the disagreement concerning universal profiles. Dissatisfaction with the adequacy of the log velocity description resulting from a simple mixing length approach has led to attempts to formulate a more general description of momentum transfer in channel flows. Parker and Coleman have investigated a wake function approach to flow description 6. Another approach is to use higher order turbulence transport equations to shift the modelling emphasis away from the level of the Paper acceptedJanuary 1988. DiscussionclosesAugust 1989. © ComputationalMechanicsPublications1989 46
Adv. Water Resources, 1989, Volume 12, March
momentum description to that of the turbulence equations. The k-e model of Launder and Spaulding 7 is one of the simplest and probably the most commonly used higher order model of turbulence. It has the advantage of 'more universal' empirical modelling constants, and thus is more likely to provide the necessary detail for an acceptable description of the effects of suspended sediment on flow velocity and vice versa. THEORY I. Basic model The k-e model is appealing for application to sediment transport description, however the basic model must be modified to describe the effects of the suspended solids. A rigorous discussion of the background and development of such a flow model has been given by Devantier and Larock 8, and will provide the basis of further discussion. The momentum equations and turbulence closure equations are similar to those presented by Rodi 9 for the k-e model for density stratified flow. For steady, vertically varying, fully developed open channel flow, the governing equations of the model may be written as
0=
1 dP Pb dx
d(uw) dx
--dU
(l)
d [v, dk]
O=-UW d - + zk
--
zj-a wg
o = < fi -uw
(2)
- < W+
zkgyzj (3)
0 = d {A(1 - A)vs} - d
{7~ww}
(4)
A suspended s e d i m e n t f l o w model: B. A . D e V a n t i e r and R . N a r a y a n a s w a m y
where flow is in the x direction, the primary variables U (horizontal mean velocity), k (turbulence kinetic energy), e (dissipation rate of k), and A (mean volume concentration of solids) vary only in the vertical direction z. Also 6 is the specific gravity difference between the sediment and water (1.65 for sand). The k ~ model constants c, (0.09), CE1 (1.51), Ce2 (1.92), Cr, 3 (0.8), O"k (1.0), and tr~ (1.3) are those specified by Rodi 9, and the overbarred terms are time averaged turbulence flux terms which will require further modelling to put them in terms of the primary variables of the model. The lower case variables in these flux terms are the time fluctuating counterparts of upper case mean quantities resulting from the classical Reynolds decomposition of time varying quantities. Properties related to the flow and transport are the bulk density of the flow suspension Pb, the net sedimentation velocity of the suspension v~, and the pressure drop driving the channel flow - d P / d x . Closure for equations (1) and (4) is accomplished by fluxgradient relations given as --dU u w = - v, dz
(5)
~v~ dA t0w . . . .
(6)
(~A dz
where v t is the eddy viscosity, and trA is a sediment Schmidt number near unity (taken as 1.2 (7)). The k-e model determines vt in terms of its primary variables as
The shear velocity for equations (10)-(12) is determined by iteratively solving the following logAaw velocity equation for u, UI.o=~ ln(K ?)
where e is a wall roughness height and K is taken to be 30 as for sand roughened pipes, and U is known from the previous iteration's estimate of the velocity field. The approach to the free surface boundary conditions is to treat it as a rigid lid and apply zero flux boundary conditions for all variables. I I . M o d e l variations
The most significant modification to be considered is a change in the way that the model accounts for damping of the turbulent motion due to the suspended solids. For the highest levels of local bed values for A, the basic model underestimates the reduction of turbulence 8. If turbulence damping is due not only to energy expended to suspend sediment, but is also due to reduction of eddy length because of the simple presence of solids, then the additional damping of turbulence could be accounted by adding a term B A to the density stratification damping term B. The mathematical form of BA is determined by drawing an analogy to the Monin-Oboukhov expression for reduction of mixing length a s 9 l/lo= 1 - c t Ri
k2 vt = - -
(7)
Two significant terms which appear several times in equations (2) and (3) are --dU Pr= -uw dz
(8)
B = - 32wg
(9)
and
(14)
where lo is the mixing length with no stratification, and R i is the Richardson number defined as Ri = - B/Pr
(15)
and 7 is an empirical constant between 3 and 10 for stably stratified atmospheric boundary layer flows. The simplest algebraic form relating mixing length to local sediment concentration would be a linear relation of the form I/l o =
is the mean shear induced production of k, and B is the buoyancy-induced production. B acts as a sink of k when the suspended solids are negatively buoyant as in the case of suspended sediment. With equations (5)-(6) substituted into equations (1)-(4) the equation set is closed, and the model may be applied if boundary conditions are supplied. The near-bed boundary conditions are specified at a slight distance Zo away from the wall. The boundary condition for velocity is actually applied as a derivative condition of the form
(13)
1-
CAA
(16)
Pr
zo/Pb = u2 = - vt dU/dz
to
(10)
where u, is the shear velocity (as defined in terms of the near-bed shear stress %), and ~c is the yon Karman constant taken as 0.42. The standard application ofk and e boundary conditions at this boundary are
where CA is a proportionality constant. This form is chosen mainly for simplicity, however there are also some reasonable physical arguments for doing so. Since the effect of B A should not become significant until A is large, it may be argued that the reduced mixing length is directly related to the reduction of volume occupied by fluid. Thus the damping is in all directions, rather than only in the direction opposing gravity as in the case of B. If the Richardson number of this volume concentration damping is taken as Rih =- -- B A / P r
then by combining equation (14) (with RiA replacing Ri) with equations (16) and (17), B A c a n be expressed as A BA=
kL =
(11) 3
e]~° = u , / ( ~ Z o )
(12)
(17)
----
CA
Pr
(18)
where ~ has been combined into CAas part of the empirical fitting constant. B h is then incorporated into the model by
A d v . Water Resources, 1989, Volume 12, M a r c h
47
substituting a combined term B' everywhere for B, where the combined term is simply
B'=B+BA
(19)
The other modifications to the model are the incorporation of damping into the boundary conditions. The application or the velocity boundary condition can be modified to include Richardson number effects so that equation (13) becomes (20) The value of RiA evaluated at z0 can be taken as zo/L, when L is taken to be the Monin-Oboukhov length of the turbulence damping due to the presence of solids. Equation (20) is most commonly expressed using zo/L, but RiA can be obtained directly in the model formulation. Itakura and Kishi 1° have successfully applied the MoninOboukhov length scale in the zo/L form for the description of sediment-induced density effects on open channel flow. For a complete discussion of the MoninOboukhov length scale, the reader is directed to Turner 11. The Richardson number (and therefore L) is negative for the turbulence damping, so that shear stress near the bed is less than that of a clear water flow. Thus the use of equation (20) will result in a slight decrease in k and e at the near-bed boundary, because the inverse solution of the equation for u, will result in smaller values for u,. The other modified boundary condition accounts for damping of turbulence near the free surface as presented by Celik and Rodi ~1
e,[H-
k3/2/H
obtain convergence for the nonlinear equation set. The equation set (equations (1)-(4)) was solved in the time steady form rather than computing transient long time solutions. For the steady, unidirectional open channel flow problem chosen there are no convective transport terms, and therefore the near-bed values of U and A could be assumed without presupposing the turbulent diffusion effects in the vertical. Thus the value of A was simply set to the measured value and U was similarly set by adjusting the roughness value e. While the importance of the ability to predict these values is not to be minimized, the emphasis in here is upon verifying effects relating to vertical turbulent diffusion of momentum and solids, and this can be accomplished more clearly if inaccuracies in the specification of these boundary conditions are not allowed to interfere. The first model test examines the effect of the free surface damping of k and e for low levels of A so that damping due to B or BA should be negligible. Fig. 1 compares predictions to cases V 14 and V 17 of the data of Vanoni 2, and Fig. 2 is run 37 from the data of Montes and Ippen 13. The dashed curves are the predictions without free surface damping and the solid curves are predictions with damping. Inclusion of the effect clearly improves the prediction, and a more S-shaped curve similar to that of equation (23) results. However the model predictions are more physically realistic than that of equation (23) in the sense that A goes to a small nonzero value rather than exactly zero.
1.0
•
I
!
•
/
t
! I
(21) t
where H is the free surface above the bed, k s is the k at z = H, and cls has been empirically determined to be 5.54. With the applications of this boundary condition vt in clear water approaches the familiar parabolic profile
I
/
•
•
Z
V 14 Z
vt tcu, ~ (H - z) =
(22)
although equation (21) produces a small, nonzero surface value for vt rather than the zero level predicted by equation (22). Similarly, when equation (22) is used to estimate the A distribution, the result is the S-shaped curve described by
A=~° =
z0j
(H - z o)
0.5 •
/
(23)
where/3 is related to the ratio of u, and vs. Thus the model should also be better to reproduce free surface concentration distribution with this modification. RESULTS The model was tested by comparing model predictions for velocity and volume concentration to measured laboratory flume flow values at steady conditions. The model equations were solved by a finite element computer algorithm incorporating Newton-Raphson iteration to
48
Adv. Water Resources, 1989, Volume 12, March
!
I
I
2 A/A a Fig. 1. Model predictions and data of sediment volume fraction versus depth
Vanoni for
A suspended sediment flow model: B. A. DeVantier and R. Narayanaswamy o O-
\ 0 o0-
Predicted Meosured
\
~
,
o--r-
o-
•
m
di o:
0
°e,~
•
c5 l
I
l
I
l
0.001
!
I
!
]
l
•!
0,01
A Fig. 2.
Model predictions and data of Montes and Ippen for sediment volume fraction versus depth
The S-shaped curves can be characterized as having 3 regions. The free surface and the bed act similarly in limiting eddy excursions toward the boundaries. This allows greater net settling in these regions than in the centre portion where mixing is best and the concentration drop is nearly log-linear with depth. The near-bed and centre region effects have been noted by McTigue t4 in the comparison with the Montes and Ippen data, however all data sets of Vanoni exhibit the three distinct regions. The difference is perhaps that the McTigue comparison is for a relatively fine, uniform sediment of smaller settling velocity. The free surface drop-off may well have been occurring, but experimental error may have masked the effect. The effect of the physical damping through BA was examined by comparing the model predictions to measurements of Einstein and Chien 4 for high levels of sediment volume fraction. The fitting constant CA was determined by varying it to obtain the 'best fit' to case S16. The plot of the velocity prediction is shown in Fig. 3. The three curves compare the model prediction with no density stratification, noted as a 'Clear' water, with buoyancy damping due to density stratification alone, noted as B, and with both B and BA effects, noted as B'. The progressive increase in velocity with increased damping effect indicates that density effects through B alone are not sufficient to properly describe the damping of momentum transport. Fig. 4 presents the clear water and sediment laden (B') curves on a semi-log plot. The
straight lines superimposed upon the curves are given to demonstrate the similarity to a plot presented by Vanoni and Nomicos 5 which showed the drag reduction effect of high suspended load. Fig. 5 compares the model prediction of A to measurements for case S16. The prediction deviates significantly from measurements above 0.2 flow depths. The poor prediction of A does not have a significant effect upon prediction of U however, because the volume concentration has dropped to low enough levels by 0.2 depths to greatly reduce the modified effect upon momentum transfer. The discrepancy is probably due to the model's simplified assumption of a single average settling velocity vs, while in the measurement there was a distribution of sediment sizes for S16. The boundary condition for k (equation (20)) was tested for the effect of the sediment modified shear velocity of equation (20). Fig. 6 is a comparison of the k predictions with equation (20) versus equation (13) for case S16 of Einstein and Chien. Use of equation (20) modified the values of k somewhat, but an interesting effect is that the lowered level ofk near the bed results in a higher level of k through most of the rest of the flow. The effect upon e of the change was similar to k and was such that it compensated to keep the vt nearly the same for either method, so that the change had practically no effect upon predictions of U and A. Thus, it appears that either form is adequate for predictive purposes. To further test the adequacy of the inclusion of solids
Adv. Water Resources, 1989, Volume 12, March
49
A suspended sediment flow model." B. A. DeVantier and R. Narayanaswamy
J Clea
1.0
I
1.0
m
I t
07 /
I I
z
Z
m
H 0.1
0.5 'B'
-// / S J6
SI6 I
I
4
I
i
I
6 U
|
8
I
I
i
4
I
I
|
6
8
U ft/s
ft/s
Fig. 3. Model predictions of velocity for Einstein and Chien case S16 with clear water, with density damping, and with density and volume fraction damping
Fig. 4. Semi-log slope plots of model predictions for Einstein and Chien case S16 for clear water, and with density and volume fraction damping
1.0
Predicted -Meosured •
Z
-
, ,/
,, , / ~ w . ~
~
H
SI6 0.5
......
.
I
I
l~161
162
,,
I
ld 3
A Fig. 5.
50
Model predictions compared to volume fraction data of Einstein and Chien case S16
Adv. Water Resources, 1989, Volume 12, March
,,J
A suspended sediment flow model: B. A. DeVantier and R. Narayanaswamy qt-
'N•,,
~
damping of turbulence through BA,the model predictions for flow with neutrally buoyant solids was tested. Elata and Ippen 15 measured the reduction in K in flows with neutrally buoyant particles. The model was used to predict the effective reduction of ~cby taking straight line semi-log slopes from plots similar to Fig. 3. Fig. 7 compares the model predictions to the data of Elata and Ippen. The agreement is quite good. It appears then, that the addition of the term B6 is justified, and the fitted value of CAis adequate.
CLEAR WAFER SHEAR VELOCITY SEDIMENT CORRECITr.D SHEAR VELOCITY
",
..... •
%'~%%
0 Z,~.. 0
SUMMARY The damping of turbulence through the presence of suspended solids and the presence of a free surface has been incorporated into a higher order turbulence closure model for sediment laden flow. The addition of the effect of solids damping through BA appears to be significant, however it should only be important at high solids concentrations. The model appears to have an advantage over to the wake model of Parker and Coleman 6 in predicting the effect of sediment upon velocity, because the wake model effects are based upon a density difference between solids and water. The reduction of sediment level
% 0
O o
I
I
I
I
I
I
I
10 -3 DIMENSIONLESS KINETIC ENERGY ( k / U 2j
I
l
10
-'
Fig. 6. Turbulence kinetic energy prediction comparison using two forms for calculatin9 shear velocity
1.00 -
\ Predicted Meosured
0.90 -
% 0.80 ! -
ib II
0.70
I[, il,
N
0.60
0.50 0.00
,
|
J
I
I
I
J
I
I
|
I
I
I
I
I
I
0.10
I
I
l
]
0.20
I
I
I
I
I
I
I
I
I
[
0.30
A Fig. 7. Comparison of model prediction of reduction of effective yon Karman constant to data of Elata and Ippen for neutrally buoyant solids
Adv. Water Resources, 1989, Volume 12, March
51
A suspended sediment flow model: B. A. DeVantier and R. Narayanaswamy d u e to free s u r f a c e d a m p i n g in t h e m o d e l s h o u l d h a v e little effect u p o n velocities, h o w e v e r s e d i m e n t level p r e d i c t i o n s are m a r g i n a l l y i m p r o v e d .
7
REFERENCES
9
1 2 3 4
5 6
52
Jaffe, P. R., Parker, F. L. and Wilson, D. J. Distribution of Toxic Substances in Rivers, J. Envir. En9r. Div., ASCE, August 1982, 108, 636-649 Vanoni, V. A. Transportation of Suspended Sediment by Water, Trans. ASCE, 1946, 111, 67 133 Zeller, J. Tonerosion, Eidgen. Techn. Hoschule, Interner Bericht, T2, T3, 1963 Einstein, H. A. and Chien, N. Effects of Heavy Sediment Near the Bed on Velocity and Sediment Distribution, U. of Cal., Berkeley, and US Army Corps of Engineers, Missouri Div., 1955, Report No. 8 Vanoni, V. A. and Nomicos, N. Resistance Properties of Sediment-Laden Streams. Trans. ASCE, 1960, 125, 114(~1152 Parker, G. and Coleman, N. L. Simple Model of SedimentLaden Flow, H. Hyd. Div., ASCE, May 1986, 35(~375
Adv. Water Resources, 1989, Volume 12, March
8
10 11 12 13 14 15
Launder, B. E. and Spaulding, D. B. The Numerical Computation of Turbulent Flows, Computational Meth. Appl. Mech. Engng., 1974, 3, 269-289 DeVantier, B. A. and Larock, B. E. Sediment Transport in Stratified Turbulent Flow, J. Hyd. Div., ASCE, December 1983, 109, 1622-1635 Rodi, W. Turbulence Models and their Application in Hydraulics, IAHR publication, 1980 Itakura, T. and Kishi, T. Open Channel Flow with Suspended Sediments, J. Hyd. Div., ASCE, August 1980, 106, 1325-1342 Turner, J. S. Buoyancy Effects in Fluids, Cambridge Univ. Press, 1973 Celik, I. and Rodi, W. Simulation of Free-Surface Effects in Turbulent Channel Floe, Physiochemical Hydrodynamics, 1984, 5, 217 227 Montes, J. S. and Ippen, A. T. Interaction of Two-Dimensional Flow Suspended Particles, R.M. Parsons Lab. Water Res. and Hydrodynamics, Mass. Inst. Tech., 1973, Report 164 McTigue, D. F. Mixture Theory for Suspended Sediment Transport, J. Hyd. Div., ASCE, June 1981, 107, 659 673 Elata, C. and Ippen, A. T. The Dynamics of Open Channel Flow with Suspensions of Neutrally Buoyant Particles, Hydrodynamics Lab., Mass. Inst. Tech., 1961, Report 45