Advances in River Sediment Research – Fukuoka et al. (eds)© 2013 Taylor & Francis Group, London, ISBN 978-1-138-00062-9

Quasi 3D numerical simulation for flow and bed variation
with various sand waves
T. Uchida & S. Fukuoka

Research and Development Initiative, Chuo University, Japan

ABSTRACT: This paper presents a quasi-3D numerical computation method for sand waves varying
with hydraulic conditions. To calculate depth scale flow structures interacting with sand waves, we developed the general Bottom Velocity Computation method (general BVC method) based on a depth integrated
model, in which bottom velocity and pressure acting on bed surface are calculated semi-­directly without
taking the assumption of the shallow water flow, such as the hydrostatic pressure distribution. The flow
separation behind the dune crest, which has been found to play an important role of the dune formation,is evaluated by the production terms of the depth integrated vorticity equations. ­Non-­equilibrium sediment transport is calculated by using momentum equation of the sediment motion with the lag ­distance.
We applied the model to experimental results on sand waves in a narrow channel and discussed the performance of the model.
1

INTRODUCTION

It is important for understanding dynamics offlow and sediment transport in rivers to evaluate
bed resistance with sand waves. The bed resistanceis separable into two components. One is the surface resistance due to sediment particles and theother is the form resistance due to sand waves.
Generally, compared to the objective flow scale,the sediment particle size is so small that the flow
around the particles is simplified as the rough walllow. On the other hand, the scale of the sand wave
and surrounding flow would be comparable to theobjective scale. Those are in an interacting system.
Although many researchers have recognized theneed to evaluate dynamics of flow with bed variation in the integrated system of sediment transport, sand waves and bed resistance (e.g., Kennedy,
1963; Ashida & Michiue, 1972; Task Committee of
JSCE, 1973; Fukuoka et al., 1982), most of previous bed variation analyses have assumed the sandwaves as the bed resistance for simplification (e.g.,
Fukuoka, 2005; Wu, 2008). It is considered for thereason that 2D model without the ability to calculate vertical distributions of velocity and pressure
intensity cannot evaluate the interactions among
sand waves, sediment transport and bed resistance. There are a lot of studies on turbulence andvortex motions in the vicinity of sand dunes (e.g.,Best, 2005). Recently, some researchers succeededin simulating dune formations in a narrow channel by vertical two-dimensional turbulence model(Giri & Shimizu, 2006; Niemann et al., 2011).

For analogous phenomena of the bed variation
with sand waves, local scouring just downstream
from a groundsill in a laboratory channel could
be reproduced by the vertical 2D model (Uchida
et al., 2004). However, it would be unrealistic toapply these models to 3D flows and bed morphology in rivers because of extremely high computational cost.
It has been commonly known that the type ofbed forms can be classified by 1D parameters suchas Froude number. The fact would implicate thatthe sand wave structures developed to the depth
scale are controlled by the depth scale dynamics offlows. For the calculations of flood flows and bed
variations in rivers, the turbulence 3D model withthe resolution capability of turbulent motions near
the bed is not always indispensable, because theeffects of sand waves on large scale phenomena of
sediment transport and flow resistance are important to be evaluated rather than the small scale
phenomena of the turbulence motions. ­However,few researchers have focused on the method forthe above purpose. Onda & Hosoda (2004) proposed the depth integrated model based on one
dimensional Boussinesq equation to computeprocesses of sand waves developing and decaying due to changing in the channel slope. They
added the effects of acceleration-deceleration near
the bottom to the irrotational flow conditions forevaluating bed shear stress on sand waves. In their
method, it would be difficult to evaluate vorticity
supply due to the flow separation behind the dune
crest, which has been found to play an important221role of the dune formation. For the condition ofa steep channel, the flow separations at hydraulicjumps should also be taken into account.
Uchida & Fukuoka (2012a, b) have developed
a new depth integrated model for reasonable bed
variation analysis, the general Bottom Velocity
Computation (BVC) method, to evaluate velocity
and pressure acting on bed surface semi-directly
without calculating those vertical distributions andassuming the shallow water flow. They have applied
the method to 3D local flow structures with horseshoe vortex around a cylinder (Uchida & Fukuoka,
2012a) and rapidly varied flows over structures(Uchida & Fukuoka, 2012b). The objective of thisstudy is to develop the general BVC method witha sediment transport model for sand waves varying
with hydraulic conditions.
2 calculation method
2.1

Frame work of the BVC method without
assumptions of the shallow water flow

It is important for bed variation analyses with sandwaves to appropriately evaluate the longitudinal
distribution of sediment transport rate affectedby bed forms (Kennedy, 1963). The general BVC
method presented in this paper (Fig. 1) makes possible to calculate bottom velocity distribution controlling sediment motions without the ­assumptionof the shallow water flow such as hydrostatic pressure distribution. The bottom velocity varies with
bed shear stress acting on bed surface and pressure distribution by the bed forms. The former
effects are evaluated in the present model mainlyby the changes in depth averaged velocity due to
bed shear stress and vertical velocity distribution
due to vorticity supply from the bed surface. To
consider the latter effects, the spatial variation
term in vertical velocity in Equation (1) and nonhydrostatic pressure terms in Equation (5) and (11)are important. And the bottom pressure intensity
is necessary to evaluate the resistance of the sandwaves. Those terms have been neglected in many
previous depth integrated models with the shallow
water assumption. The general BVC method without assuming the shallow water flows (Uchida &Fukuoka, 2012a) is necessary to calculate bed variation with sand waves.
The important equation of the BVC method isthe bottom velocity equations (1), which is derived
by depth-integrating horizontal vorticity.
ubi = usi − ε ij 3Ω j h −

∂z
∂z
∂Wh
+ ws s − wb b
∂xi
∂xi
∂xi

(1)

where, i, j = 1(x), 2(y), eij3: Levi-Civita symbol,
ubi: bottom velocity, usi: water surface velocity,Ωj: Depth Averaged (DA) horizontal vorticity,
h: water depth, W: DA vertical velocity, zs: waterlevel, zb: bed level, ws, wb: vertical velocity on water
surface and ­bottom, respectively. The bottompressure intensity is also important and given by
Eq. (2) integrated the vertical momentum equation with respect to the vertical direction for water
depth, neglecting unsteady and horizontal shear
stress terms.
dpb ∂hWUj τ bj ∂zb
+
=
∂x j
ρ
ρ ∂x j(2)where, dp: pressure deviation from hydrostaticpressure distribution (p = ρg (zs − z) + dp), dpb: dp
on bottom, Uj: DA horizontal velocity, τbj: bed
shear stress. By using Eq. (1) and Eq. (2), theBVC method allows semi-direct computations of­bottom ­velocity and bottom pressure intensity
without computing those vertical distributions.
2.2

Figure 1. Framework of the BVC method.

Governing equations of the BVC method

To evaluate bottom velocity and bottom pressure
intensity by Eqs. (1) and (2), the BVC method
solves equations for the following unknown­quantities: water depth h, DA horizontal velocity
Ui, DA turbulence kinetic energy k, DA horizontal­vorticity Ωi, horizontal velocity on water surface usi,222

DA vertical velocity W. These equations are derived
by using cubic vertical velocity distribution:

(

)

(

ui = U i + ∆ui 12η3 − 12η2 + 1 − δ ui 4η3 − 3η2) (3)where, ∆ui = usi − Ui, δui = usi − ubi, η = (zs − z)/h.
And the linear distribution of dp is assumed(dp = η/dpb).
Water depth h and DA horizontal velocity Ui
are solved by Depth Integrated (DI) continuity equation (4) and DI horizontal momentum
equations (5):∂h ∂U j h
=0
+
∂x j
∂t

(4)

∂hτ ij∂z τ∂U i h ∂U iU j h= − gh s − bi +
+
∂xi
∂x j
∂t
ρ ρ∂x j ∂hdp0 dpb ∂zb −+ρ ∂x 
 ρ∂x(5)where, τij = vSij + ui ′ u j ′ , v = vm+vt, vm: kinematic­viscosity coefficient, vt: kinematic eddy viscosity
coefficient, Sij: DA strain velocity, ui ′ = ui − Ui.
In this paper, superscript bar indicates depth integrated value and evaluated by vertical velocity­distribution (3).
For the kinematic eddy viscosity coefficient, thisstudy employs one equation model:∂k∂k 1 ∂  vh ∂k 
+U j=
+ Pk − ε∂t∂x j h ∂xi  σ k ∂xi ∂hDω ij∂Ω i h= ERσ i + Pω i +
∂x j
∂t

(8)

where, ERσi: rotation term of vertical ­vorticity(ERσi = usiωσs − ubiωσb), ωσs, ωσb: rotation of usi,
ubi, respectively, Pωi: production term of vorticity
from the thin bottom vortex layer, Dωij: horizontal ­vorticity flux due to convection, rotation, dispersion and turbulence diffusion. The production
terms and horizontal vorticity flux are:Pω i = C pωvtb

ω bei − ω bihDω ij = ω j ui − ω i u j +ii

The equations for horizontal voticity aredescribed:(9)

vt ∂Ωiσ ω ∂x j(10)where, Cpω = κ /α, vtb: vt on the bottom converted
into depth averaged scale, ωbi: horizontal vorticity
on bottom, ωbei: equilibrium ωbi for ubi (Eq. 18),σω = 1. ωi is evaluated by the vertical differentiationof Eq. (3).
The equations for water surface velocity arederived by assuming very thin layer under thewater depth (δzs→0, δzs: thickness of the water
surface layer):
∂u
∂usi

 ∂zs ∂dp + Psi
+ usj si = −  g + 



∂x j
∂
z
∂t
z = zs 
 ∂xi


(6)

(11)

where, Psi is shear stress under the water surface
layer and evaluated by using Eq. (3):where,()

2Pk
k1 / 2δu = 2 Sij 2 + s ′ ij 2 + Ch  i  , vt = C µ
,
 h 
vt
ε

Psi =

Cµ = 0.09, ε = Cεk3/2/∆, Cε/∆ = 1.7/h, Ch = 9(αCε)4/Cµ3, α = κ/6, κ = 0.41, s′ij: strain velocity of deviation velocity ui′. The Eq. (4) is calculated with Uihat n + 1 step by the following equation for time
variation of water depth:

∂ 
2 ∂φ
P
 (C ∆t ) +φ −φ = 0∂x j 
∂x j 

(7)

where, ∆t: time interval of the numerical integration, f n, f n+1: a value of f at the n, n + 1 time
step in the numerical calculation, respectively, f p:
predicted value of f for the n + 1 step including n
step’s variable, φ = hn+1 − hn, φP = hP − hn, C2 = gh, hP:­predicted water depth by Eq. (4) with (uih)P, (uih)P:
predicted horizontal flux by Eq. (5) with h = hn.

2vt
12 (usei − usi ) − (3δ ui − 6 ∆ui )
h2

{

}

(12)

where, usei: equilibrium water surface velocity
usei = Ui + (δui − ∆ui)/2.
The time variation in the DI vertical velocity­calculated by Eq. (7) with φ = (Wh)n+1 − (Wh)n,φP = (Wh)P, C = k1h/∆t, k1 = 1/20. (Wh)P is calculatedby the continuity equation with horizontal velocity
field by using predicted bottom velocity (ubi)P:(Wh )P ={}∂  2 n +1(h )k2 ( ∆u j )n +1 − k1(δ u j )P 
∂x j 
 ∂z∂z + h n +1  m + U j m ∂x j 
 ∂t

n +1

(13)

where, k2 = 1/10, (∆ui)P = (usi)n+1 − (ubi)P, (ubi)P is­evaluated by Eq. (1) with (Wh)n, zm = (zs + zb)/2.

223

2.3

Evaluation method for flow separation
behind the crest on dune and hydraulic jump

The flows with sand waves are characterized byflow separation behind the dune crest and wave
breaking at the hydraulic jump (e.g., Kennedy,
1963; Task Committee of JSCE, 1973; Fukuoka
et al., 1982). It is important for flow and bed variation analysis with sand waves to evaluate these flow
separations at bed and water surfaces. The present
analysis uses the criterion variation in bed andwater surface gradients for the flow ­separation.
The values of 1/10 and 1/5 are setting for flow
separations criterions for bed and water surface,respectively, referring to previous study (Uchida,2011; Uchida 2012b). To consider vorticity production due to bottom flow separation and breaking wave, following vorticity fluxes are added in the
production terms:(Bottom flow separation)
−ε u / 2 ⋅ ubj ( + ε ij 3δ zbj > 1 / 10 )Pω sbi ∆x =  ij 3 bj ′(14)
( + ε ij 3δ zbi < 1 / 10 )0(Breaking wave)
+ ε u / 2 ⋅ usj ( −ε ij 3δ zsj > 1 / 5)Pω ssi ∆x =  ij 3 sj ′( −ε ij 3δ zsj < 1 / 5)
0

(15)

where, j′ = j, δzbj, δzsj: variations in bed, water surface gradients between the two neighbor grids in
j direction, respectively. For the breaking wave,
following conditions are added. One is the discontinuity condition for the advection terms in
the water surface velocity equation. The surface
velocity at the breaking point is set to zero for theadvection term. The other is the pressure intensity
just downstream of the breaking wave point is set
to zero, if negative pressure intensity is calculated
by Eq. (2).
2.4

Evaluation method of bed shear
stress and bed tractive force

The wall low at roughness bed is important toevaluate bed tractive force for sediment transport
rate. The definition of the relationship between thebottom height for bed shear stress and the heightof the surface sediment particle should be madefor this purpose (e.g., Colombini, 2004). The logarithmic velocity distribution law cannot be applied
under the particle surfaces, because of their dragforces. The bottom height in this paper is definedas Figure 2.
The bottom level zb to evaluate bed shear
stress is defined as zb = δzb + zg (zg: top height ofthe surface particle, δzb: very thin vortex layer).

Figure 2. Definitions of heights of bottom and sediment particles.

The origin of the logarithmic law z0 is found to
be defined as the volume averaged height of thesurface particle (Colombini, 2004). However,
it is known that the origin level depends on thesediment particle conditions (e.g., Uchida et al.,2001). In this paper, the origin level z0 is definedby the top height of the surface particle zg, equivalent roughness ks and coefficient a (z0 = zg − aks).
And this paper sets a = 1 for simplicity. The thickness of the vortex layer on the bed is defined as h/δzb = e3 − 1 to satisfy the vertical velocity distribution for uniform flow:
ui =

(

)

δ ui
1 − 3η2 + U i
3(16)

Assuming very thin thickness of the vortex layer,
we can evaluate bed shear stress τb and equilibrium
vorticity ωbe on the bottom:

τb
1
1  δ z + aks = u*2 = (cbub )2 ,= 8.5 + ln  b (17)ρκ cbks
ω be = 2

u*  δ zb + aks + h lnκ h  δ zb + aks (18)

The bed surface disturbance growing into sandwaves would be created mainly by unsteady 3D
turbulence motions in the vicinity of the bed. The
additional term of the disturbance is needed forthe present method based on the depth integrated
model. Although many methods for adding the
disturbance are considered, a simple way should be
required for the application to flood simulations inrivers. This paper gives the fluctuation of bed tractive force due to turbulence motions whose scale issmaller than the water depth.

224τ * = (1 + α B )

u*2
σ gd

(19)

  x + y t + 
α B = AB cos 2π 
T 
  L

where, L and T: spatial and temporal scales of the
disturbance, respectively. Although the fluctuation
αB is affected by local turbulence conditions, thispaper gives L = h, T = L/ub, AB = 1.0, for simplicity. We confirmed that the effects of these parameter in the Eq. (19) on the developed form of sandwaves can be negligible.
2.5 Governing equation for non-equilibrium
sediment transport rate with sand waves
It is important for bed variation analysis withsand dunes to evaluate non-equilibrium sedimentmotion over them (Giri & Shimizu, 2006; Niemannet al., 2011). In this study, the non-equilibrium sediment motion is calculated by the sediment particle
momentum equation (Uchida & Fukuoka, 2011):∂qBi ∂uBj qBi= ( PuBPi − DuBDi ) + m*hB (γ ei − γ i )
+
∂x j
∂t
(20)

where, qBi: bed load sediment transport rate ­vector,
uBi: particle velocity vector of bed load, PuBPi: particle momentum from bed material, DuBDi: lossof particle momentum due to particles deposition, m* = µk sgcosθ/(s + 1 + CM), CM: coefficientof added mass, µk: dynamic friction coefficient ofparticles, s: specific gravity of sediment in water,
θ: bed gradient, hB: apparent thickness of bed load
(hB = qB/uB), qB2 = qBiqBi, uB2 = uBiuBi, γi, γei: unit vector of sediment movement and that for the equilibrium condition, respectively. This paper gives
the first term in Eq. (20) and the apparent thickness hB and velocity uB of bed load as (Uchida &­Fukuoka, 2011):
PuBPi − DuBDi = qBe uBei /Le − qB uBi /L
1/ 3
usgd = 3.3 (qB * ) , uBe
B
2/3hB d = 0.3 (qB * )(21)

sgd = 3.3 (qBe* )

1/ 3

(22)


where, L: lag distance of bed load motion, Le,
qBe, uBei: equilibrium lag distance, sediment transport rate and apparent velocity of bed load corresponding to local bed tractive force, respectively,
uBe2 = uBeiuBei, qB* = qB/(sgd3)1/2, qBe* = qBe/(sgd3)1/2.
The equilibrium sediment transport rate is
calculated by the Ashida & Michiue Formura(1974). The lag distance is important parameter to­ etermine non-equilibrium sediment motion, cord
relating sediment transport rate with depositionand pick-up rates. Although some empirical formulae for the lag distance have been proposed (e.g.,
Phillips & Sutherland, 1989; Fukuoka, 2005; Wu,
2008), those applicability to sand waves calculation
has not been clarified. The lag distance is identified as the average particle step length for uniform
flow condition without varying in bed tractive
force. However, for the non-equilibrium flow condition with sand waves, the lag distance would bedecreased due to acceleration and deceleration ofbottom velocity. At present, its evaluation method
has not been established yet. This paper uses equation (23) for the lag distance, which is transformed
from the previous empirical equations to apply
bed variation analysis (Uchida & Fukuoka, 2011),
adjusting coefficient αL to reproduce experimental
results by the calculation.
Le d = α L (qBe* ) , L d = α L max (qBe* , qB * )
2/3

2 /3

(23)

As increasing αL in Eq. (23), the wave length isincreasing with decreasing in the wave height, andthe plane bed appears eventually. This study sets
αL = 7, considering that the lag distance is shorter
than the average particle step length.
The temporal variation in the bed elevation iscomputed with neglecting variations in apparent
thickness of bed load:(1 − λΒ )∂zB ∂qBi
+
=0
∂t
∂xi

(24)

where, λB: the porosity of bed material (λB = 0.4).
3 calculation results
and discussions
3.1 Calculation conditions
This study applies the present method to experiments on the sand waves conducted in the narrow
channel by Fukuoka et al. (1982). This experiment investigated flow and sediment transport
rate with various types of sand waves generated in
the acrylic channel with 0.04 m width and 8.0 m
long for diverse hydraulic conditions of bed slope,
sediment particle diameter and water discharge.
Refer to the literature (Fukuoka et al., 1982) for
details of the experimental conditions and results.
This paper presents calculation results on a part oftheir experimental cases, as shown in Table 1. The
boundary conditions of the upstream end for calculations are given by the experimental discharge.
The shear stress acting on bed and side wall are225

Table 1. Hydraulic conditions for various bed forms experiment (Fukuoka et al., 1982).
Sand wave typeCaseIq (m2/s)d (mm)

Fr

ExperimentCalculation

1
2
3

1/178
1/147
1/73.043.0 × 10-3
43.0 × 10-3
44.3 × 10-3

0.76 × 10-3
0.76 × 10-3
0.76 × 10-3

0.66
0.65
0.88

4
5
6
7

1/31.01/31.01/31.01/20.039.5 × 10-3
10.8 × 10-3
26.8 × 10-3
13.3 × 10-3

0.76 × 10-3
0.19 × 10-3
0.19 × 10-3
0.19 × 10-3

1.74
1.11
1.24
0.65

Dune
Duen
Anti-dune moving
downstream
Plane bed
Anti-dune
Anti-dune
Anti-dune

Dune
Dune
Dune with longwave length
Plane bed
Anti-dune
Anti-dune
Anti-dune

I: channel slope, q: discharge per unit width, d: sediment particle diameter, Fr : averaged Froude number, Fr = q/(gh3)0.5(Local Fr varies longitudinally especially for case 5∼7).

evaluated with equivalent roughness ks = d (d: sediment particle diameter) and Manning’s roughness
coefficient n = 0.008, respectively. The computational mesh size is set to dx = dy = 0.01 m to calculate the depth scale vortex motion. For the initialcondition, flow calculation results with flat fixed
bed are given. Table 1 shows sand wave types in
the measurements and calculations. Although
the regime of anti-dune moving downstream forCase 3 was not appear in the present calculation,
the three different types of dune, plane bed andanti-dune were reproduced by the general BVC
method with non-equilibrium sediment transportmodel. The followings are discussions on calculation results for each sand wave type.
3.2

Dunes (Case 1, Case 2)

Figure 3 shows the developing stage of dune formation for Case 2. The Froude Number for the
initial condition is larger than Fr = 0.8. We can seethat water depth increases with developing the bed
form disturbance. At the initial stage (0.5 min), thebed form disturbance grows up locally, and thensmall dunes with short wave length appear in
the whole of the channel (1.0 min). These dunes
grow up with dune coalescences moving downstream and become dunes with long wave length(2.0 min). After that, the average size of sanddunes and water depth do not change (4.0 min).
However, dune sizes are different from each otherand a dune moves downstream with dune coalescences and breakups.
Followings are discussions on the results fordeveloped dunes. Figure 4 shows computed velocity distributions of flows with developed dunes
for Case 2. We can see relatively high velocity at
the middle depth and low velocity at the water
surface on the dune crest. These velocity distributions are generated under the non-hydrostatic

Figure 3. Developing stage of dune formation for Case 2.

Figure 4. Computed velocity distributions with developed dunes for Case 2.

pressure ­distributions by the effects of the dune
shape. The calculation results describe characteristics of bottom velocity of flows with sand dunes:
low ­bottom velocity in the separation zone behindthe dune crest, accelerating downstream of the226r­ eattachment point of the separation, and highbottom velocity on the crest. The computed dune
has a mild slope with upward convex upstreamof the crest and a sharp slope downstream of it,
whose characteristics have been investigated in
the previous study (e.g., Kennedy, 1963). Figure 5
shows computed water surface and bed profiles forCase 1 and Case 2. Table 2 shows quantitative comparisons between measured and computed resultsof averaged water depth, sediment discharge perunit width, wave length, wave celerity, and wave
height of dunes. It is noticed that computed water
depth are larger than measured depths and computed wave lengths are longer than measured
wave lengths. Although computed wave celerity islower than measured one, averaged wave heights
kd = 2qB/C using calculation results are larger thanthose of measurements. These results indicate thatsand dunes in the calculation are more developed
compared with those in the experiment, calculating large dune resistance. However, resulting differences of water depth and sediment discharge
between measured and calculated results are notso large. Therefore it would be safe to say thatthe present calculation method can reproduce thequalitative and quantitative characteristics of sanddunes generated in the experiments.
3.3 Plane bed and anti-dune moving downstream(Case 3, Case 4)
Figure 6 shows computed water surface and bed
profiles for Case 3 and Case 4. Table 3 shows quantitative comparisons between measured and computed results of sand waves. In the experiment,anti-dunes moving downstream were observed forthe condition of Case 3. This kind of sand waves
also occurred in the other experimental conditions(Fukuoka et al., 1982). However, flat dunes withlong wave lengths are computed for the Case 3. The
computed wave length and wave celerity are muchlarger than those of measured results. On the other
hand, differences between measured and computed
water depth and sediment discharge seem to besmall. This result would indicate that the resistanceof the anti-dunes moving downstream is smaller
than that of the dunes. More detailed discussions
on the anti-dune moving downstream are futuretasks. For the condition of Case 4, plane beds
appear both in the experiment and ­computation.

Figure 5. Computed water surface and bed profiles forCase 1 and Case 2.

Table 2. Measured and computed data for Case 1 andCase 2.
qB(m2/s)

L
(m)

C
(m/s)

kd
(m)

Case 1
Exp. 0.074
Cal. 0.079

1.8 × 10-5
1.4 × 10-5

0.17
0.26

5.7 × 10-3
2.4 × 10-3

6.3 × 10-3
11.6 × 10-3

Case 2Exp. 0.077
Cal. 0.079

2.0 × 10-5
1.8 × 10-5

0.25
0.30

7.8 × 10-3
3.0 × 10-3

5.1 × 10-3
12.0 × 10-3

h: averaged water depth, qB: averaged sediment discharge
per unit width, L: averaged wave length of sand dunes,
C: averaged wave celerity of sand dunes, kd: average wave
height of sand dunes (kd = 2qB/C).

h(m)

Figure 6. Computed water surface and bed profiles forCase 3 and Case 4.
Table 3.
Case 4.

Measured and computed data for Case 3 and
h
(m)

qB
(m2/s)

L
(m)

C
(m/s)

Case 3
Exp.
Cal.

0.064
0.063

4.7 × 10-5
8.9 × 10-5

0.23
0.48

7.6 × 10-3
1.7 × 10-3

Case 4
Exp.
Cal.

0.037
0.037

2.3 × 10-4
4.9 × 10-4

Plane bed
Plane bed

227

The computed water depth is almost the same as
the measured depth. And the rapidly increased
sediment discharge in the experimental result ofCase 4 compared with Case 1∼3 is well-explainedby the calculation.
3.4

Anti-dunes (Case 5∼7)

Figure 7 shows longitudinal water surface profiles
and bed topographies for Case 5∼7. Table 3 shows
data of sand waves in the experiment and the calculation for Case 5∼7. In the present calculation,
anti-dunes moving upstream are generated as theexperimental results. For the condition, generations
and disappearances of anti-dunes and hydraulicjumps are seen successively (Fig. 7 shows snap
shots when hydraulic jumps appear). Therefore,
temporal changes in the hydraulic conditions forCase 5∼7 are large. Table 4 presents the maximum
amplitudes with average values. The computed
water depths are smaller than those of measured
results. One of the reasons would be large temporal and spatial variations in water depth. We
can see large temporal variations in sediment discharge, especially 70% fluctuations for Case 7. The
averaged sediment discharges in the experiment arealmost reproduced by the calculation. The computation results could explain the wave lengths andwave celerity for anti-dune conditions.
4 conclusions

Figure 7. Computed water surface and bed profiles forCase 5∼7.
Table 4. Measured and computed data for Case 5∼7.
h*
(m)

qB*
(m2/s)

Case 5
Exp. 0.021
1.1 × 10-4
Cal. 0.016 ± 0.003 (1.9 ± 0.6)
× 10-4
Case 6
Exp. 0.036
2.7 × 10-4
Cal. 0.030 ± 0.002 (2.8 ± 0.6)
× 10-4
Case 7
Exp. 0.035
3.2 × 10-4
Cal. 0.016 ± 0.005 (2.9 ± 1.8)
× 10-4

L*
(m)

C
(m/s)

This study developed the general BVC methodwith non-equilibrium sediment transport model tocalculate bed variation with sand waves for various
hydraulic conditions. The present method was able
to calculate three different types of sand waves:dune, plane bed and anti-dune and to describe
the magnitude of sand waves, water depth andsediment transport rate in the experiment. The
governing equations in the present method were
derived from Reynolds-averaged Navier-Stokes
equations without the shallow water assumption,using vertical velocity and pressure distributions.
Non-equilibrium sediment motion was calculated
by sediment particle momentum equations. We
did not employ special formulae and coefficients,
except for bed tractive force fluctuations (Eq. 19)
and the lag distance (Eq. 23). Because little influence of parameter variations in Eq. (19) on developed sand wave conditions are confirmed, it would
be emphasized that the present method could
explain characteristics of various kinds of sandwaves with constant parameter of Eq. (23). We
plan to develop the method to calculate flood flows
and bed variations in river with 3D sand waves.
REFERENCES

0.57-8.1 × 10-3
1.2 ± 0.1 -1.3 × 10-2
1.5
-1.2 × 10-2
1.5 ± 0.2 -1.6 × 10-2
0.78
-9.2 × 10-3
1.2 ± 0.4 -2.0 × 10-2*: Computed values indicate averaged value ± maximum
amplitudes for channel center (x = 3∼6 m).

Ashida, K. and Michiue, M. 1972. Study on hydraulic resistance and bed-load transport rate in alluvial
streams, Proc. of JSCE 206:59–69.
Best, J. 2005. The fluid dynamics of river dunes: A review
and some future research directions, J. Geophys. Res.
110 F04S02, doi:10.1029/2004JF000218.
Colombini, M. 2004. Revisiting the linear theory of sanddune formation, J. Fluid Mech. 502:1–16.
Fukuoka, S. 2005. Flood Hydraulics and River Channel
Design, Morikita publishing company.
Fukuoka, S., Okutsu, K. and Yamasaka, M. 1982.
Dynamic and kinematic features of sand waves inupper regime, Proc. of JSCE 323:77–89.

228Giri, S. and Shimizu, Y. 2006. Numerical computation of sand dune migration with free surface flow,
Water Resources Researches, Vol. 42, W10422,
doi:10.1029/2005WR004588.
Kennedy, J.F. 1963. The mechanics of dunes and antidunes in erodible-bed channel, J. Fluid Mech.
16:521–544.
Niemann, S.L., Fredsøe, J. and Jacobsen, N.G. 2011.
Sand dunes in steady flow at low Froude numbers:dune height evaluation and flow resistance, J. Hydraul.
Eng. 137(1):5–14.
Onda, S. and Hosoda, T. 2004. Numerical simulation
on development process of dunes and flow resistance, Annual Journal of Hydraulic Engineering JSCE
48:973–978.
Phillips, B.C. and Sutherland, A.J. 1989. Spatial lag
effects in bed load sediment transport, J. Hydraul.
Res. 27(1):115–133.
Task Committee on the Bed Configration and Hydraulic Resistance of Alluvial Streams, Committee on
Hydraulic and Hydraulic Engineering, JSCE. 1973.
The Bed configuration and roughness of alluvial
streams, Proc. of JSCE 210:65–91.

Uchida, T. and Fukuoka, S. 2011. A bottom velocity
computation method for estimating bed variation in achannel with submerged groins, Journal of JSCE Ser.
B1 (Hydraulic Engineering) 67(1):16–29.
Uchida, T. and Fukuoka, S. 2012a. Bottom velocity computation method by depth integrated model without
shallow water assumption, Journal of JSCE Ser.B1(Hydraulic Engineering) 68(4):I_1225-I_1230.
Uchida, T. and Fukuoka, S. 2012b A computation
method for flow over structures, Advances in River
Engineering 18:351–356.
Uchida, T., Fukuoka, S., Fukushima, T. and Tanaka, M.
2001. Two-dimensional shallow water analysis overlarge roughness elements and its application, Journalof Hydraulic Coast and Environmental Engineering
JSCE 691/II-57:93–104.
Uchida, T., Fukuoka, S. and Watanabe, A. 2004. Vertical two-dimentional analysis for local scour justdownstream from a groundsill, Journal of Hydraulic
Coast and Environmental Engineering JSCE 768/II-68:
45–54.
Wu, W. 2008. Computational river dynamics, Taylor &Francis, London.

229