Advances in Water Resources 53 (2013) 12–22
Contents lists available at SciVerse ScienceDirect
Advances in Water Resources journal homepage: www.elsevier.com/locate/advwatres
Net beach change in the swash zone: A numerical investigation Fangfang Zhu ⇑, Nicholas Dodd Infrastructure and Geomatics Division, Faculty of Engineering, University of Nottingham, Nottingham NG7 2RD, England, UK
a r t i c l e
i n f o
Article history: Received 8 June 2012 Received in revised form 1 October 2012 Accepted 2 October 2012 Available online 23 October 2012 Keywords: Swash Sediment transport formulae Bed profile Bed shear stress
a b s t r a c t A range of bed-load sediment transport formulae are used to run fully coupled, morphodynamic simulations of one [1] swash cycle on an erodible plane beach. A system comprising shallow water equations and Exner equation is solved, in which sediment transport rate q is either dependent only on depth averaged velocity (u), or on u and water depth (h). The results are in agreement with equivalent uncoupled results [2] in that all sediment transport formulae considered applied to the event of [1] yield net erosion in the whole of the swash. Consistent with [3], however, full coupling yields significantly less erosion for all the q ¼ qðuÞ formulae compared to the equivalent uncoupled results. It is shown that differences between uncoupled and coupled approaches (for most formulae) accumulate over the course of the swash event. The main reason for the reduced net erosion is the smaller maximum inundation. It is also shown that including a dependence on h in the bed load sediment transport formula for fully coupled simulations can result in net deposition in the upper swash. Bed shear stress described by a Chezy law is further included in fully coupled simulations to examine net beach change. Much reduced maximum inundation and net offshore sediment transport are predicted both for q ¼ qðuÞ and q ¼ qðh; uÞ. It is shown that although the net sediment flux at the base of the swash under one [1] swash event is still offshore, deposition in the middle or upper swash may be predicted when bed shear stress is included, particularly when the drag coefficient in the backwash is reduced compared to that in the uprush, consistent with some in-situ measurements. The implication is that bed shear stress must be included not just to obtain correct quantitative beach change, but also to obtain correct qualitative beach behaviour. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction The swash zone is a very dynamic and complex regime of the nearshore, in which the beach face changes rapidly. On steeper beaches bore-driven swash is the dominant swash motion and significant progress in understanding bore behaviour has been made, e.g. [4–6]. These motions also have significant implications for planform evolution in the swash and its modelling [7,8]. Peregrine and Williams [1], henceforth PW01, generalized the shallow water solution in [6], for the shoreline motion and its near region, to the whole swash, leading to a dam-break problem on a sloping bed. Work by Guard and Baldock [9] has indicated that the PW01 solution is only a special case of the swash event in [6], as it neglects the momentum behind the bore. However, subsequent work by Pritchard [10] suggests that modifying the hydrodynamic boundary conditions does not make great qualitative difference to net sediment transport. Pritchard and Hogg [2] looked in detail at this event, in particular using depth-averaged velocity ^ from the PW01 analytical solution to ^ ) and flow depth (h) (u ⇑ Corresponding author. E-mail address:
[email protected] (F. Zhu). 0309-1708/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.advwatres.2012.10.002
provide predictions of instantaneous equilibrium (steady state) sediment transport and also bed-load transport over the swash event. Then, by integrating over the swash duration (at each location) they showed that a whole range of sediment transport formulae, when interpreted as instantaneous sediment fluxes that adjust immediately to change in flow, all yield qualitatively similar net offshore transport at every point in the swash. For net deposition to result in some regions of the swash (which it must sometimes if beaches are not continually to be eroded) it was necessary either to introduce a settling and entrainment lag, or to assume that sediment is pre-suspended in the surf zone. Subsequently, Kelly and Dodd [3] considered a fully coupled model in which the flow equations and a bed evolution (Exner) ^ ¼ Au ^ instanta^ 3 (q equation with sediment transport formula q neous sediment flux and A bed mobility parameter) are solved simultaneously, and the bed change thus allowed to feed back onto the swash hydrodynamics. Kelly and Dodd [3] showed that those same PW01 initial conditions on an initially plane but now erodible bed also lead to net erosion for the whole of the PW01 swash, but significantly less than that predicted by the uncoupled model [2]. In related work, qualitatively different bed profiles were obtained by Briganti et al. [11] for the equivalent coupled, (initially)
13
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
^ ¼ Au ^ 3 formula, and flat mobile-bed dam-break problem for the q another equilibrium formula considered by Pritchard and Hogg ^u ^ ¼ Ah ^ 3 (A also bed mobility parameter). These results point [2], q to the possible importance of considering different sediment transport formulae when full morphodynamic coupling is implemented, to the extent that different equilibrium sediment transport and bed-load transport formulae may yield qualitatively different results from uncoupled predictions and from each other. A further important physical mechanism affecting swash zone dynamics is bed shear stress. Its effect is particularly significant in the region of shallow water flow near the shoreline ([12,13]). To what extent the bed shear stress influences sediment transport in the swash, and whether a qualitatively different bed profile could be predicted when bed shear stress is included are of great interest and importance for the swash simulation. Therefore in this work we examine a range of bed-load and equilibrium sediment transport formulae similar to those examined by Pritchard and Hogg [2] under the PW01 swash event in fully coupled simulations on a mobile, initially plane sloping bed. The sediment transport formulae examined are consistent with beaches where bed load transport is dominant, or on which suspended load transport occurs but with negligible scour and settling lag. The purpose is to examine the range of predictions for the different formulae, all of which are formulae in use in engineering models or extensions thereof, in fully coupled simulations, these being a true representation of morphodynamics on real beaches. Further, differences between fully coupled and uncoupled approaches are considered, as well as what this means for numerical modelling in the swash region (and also modelling of tsunami inundation). Importantly, a particular focus is whether the same event can yield net sediment movement of opposing signs depending on which sediment transport formula is used. Accordingly, we also examine the general dependence of net beach change in the swash via the sediment transport formula q, on powers of u and h. Lastly, bed shear stress, in the form of a Chezy drag law, is included in fully coupled simulations for the first time, to investigate the effects of bed friction on the beach face evolution, for different sediment transport formulae. The effect of varying drag coefficient between uprush and backwash is also examined. In the next section we present the equations examined. We then state the sediment transport formulae to be examined in Section 3. In Section 4 the PW01 event simulations are presented. Finally, conclusions are arrived at.
Fig. 1. Initial conditions for PW01 swash.
^t
x¼
^x ; h0
B¼
b ^ q B and q ¼ ; q0 h0
t¼
1=2 h0 g 1=2
h¼
;
^ h ; h0
u¼
^ u ðgh0 Þ1=2
; ð4Þ
where h0 is a length scale, and q0 represents a sediment flux scale. Substituting (4) into the governing Eqs. (1) and (2) gives:
ht þ uhx þ hux ¼ 0;
ð5Þ
ut þ uux þ hx þ Bx ¼ 0:
ð6Þ
Assuming q ¼ qðh; uÞ and substituting (4) into (3) gives:
Bt þ rqh hx þ rqu ux ¼ 0; where
ð7Þ
r ¼ g1=2nqh03=2 . 0
These equations can be written in vector form:
! !! U t þ Að U Þ U x ¼ 0
ð8Þ
with
2 3 h ! 6 7 U ¼ 4 u 5; B
2
u
! 6 Að U Þ ¼ 4 1
h u
0
3
7 1 5:
rqh rqu 0
The eigenvalues of A are the roots of the polynomial, 2. Mathematical model
k3 2uk2 þ ðu2 rqu hÞk þ rðuqu hqh Þ ¼ 0;
ð9Þ
2.1. Governing equations The system of equations governing 1D shallow water flow with a mobile bed is:
^^ þ u ^^ þ h ^u ^h ^^x ¼ 0; h x t ^^ þ g B b x^ ¼ 0; ^ ^ ^ u^t þ uu^x þ g h x b ^ þ nq ^^x ¼ 0; B t
ð1Þ ð2Þ ð3Þ
^ represents water depth (m), u ^ is a depth-averaged horizonwhere h b is the bed level (m), q ^ is sediment flux (m2 s1 ), tal velocity (ms1 ), B 1 ^ ^ which is, in general, a function of h and u; n ¼ 1p with p being bed porosity, and g is acceleration due to gravity (ms2 ). In Fig. 1 we illustrate the situation being considered. 2.2. Non-dimensionalization To make the results more intercomparable, we non-dimensionalize all variables. Dimensionless variables are:
the roots of which may be denoted k1 ; k2 and k3 , such that k1 6 k3 6 k2 . When r ! 0, one root tends to 0, and the other two tend to the corresponding hydrodynamic characteristic speeds. 2.3. Riemann equations and numerical solution The specified time interval method of characteristics (STI MOC) which has been successfully used in [3,14] is employed. The resulting Riemann (characteristic) equations for a general formula qðh; uÞ are:
RðkÞ ¼ kk
du kk þ rqh dh dB þ þ ¼ 0; dt kk u dt dt
k ¼ 1; 2; 3;
ð10Þ
which hold along the characteristics dx ¼ kk , with k defining the dt characteristic family. This method is used here to achieve very high accuracy in the vicinity of discontinuities (e.g. bores and hydraulic jumps), at which locations usual engineering codes lose some accuracy; see [15]. See [3] for further details of the numerical method.
14
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
3. Different sediment transport formulae Six of the more commonly used bed load formulae are consid^ ¼ Au ^ 3 , the Grass formula [16,17,2,14]; the Bagnold formula, ered: q see [17,18]; the Meyer-Peter–Müller formula, see [17]; the Van Rijn ^u ^ ¼ Ah ^ 3 , which formula, see [17]; the Bailard formula, e.g. [2] and q was considered by Pritchard and Hogg [2], and which we refer to as the PH formula. Note that two of these formulae (Van Rijn (IV) and Bailard (V)) have been slightly modified (thresholds have been excluded) for the purpose of examining the effects of different powers of u on sediment transport; the names of formulae, however, are not changed. Additionally, note that formulae II and III reduce to I if ue ! 0. Note that some of these formulae (Van Rijn and Bailard) have been slightly modified (thresholds have been excluded) for the purpose of examining the effects of different powers of u on sediment transport; the names of formulae, however, are not changed. These give rise to six models for the swash zone simulations. The expressions of dimensional and dimensionless q and the corresponding variables of different sediment transport formulae are summarized in Table 1, wherein it should be noted that the six different formulae have been numbered I to VI, with the different dimensional coefficients also numbered correspondingly. ^ e is the value of velocity at the threshold of motion, and ue Here u the corresponding dimensionless value. It should be noted that for the Bagnold and Meyer-Peter–Müller formulae, if juj < ue , then q ¼ 0, thus qh ¼ 0 and qu ¼ 0. When juj < ue , it is found that one of the roots of the polynomial (9) is 0, and the other two are identical to the hydrodynamic characteristics. Note also that there is a discontinuity in qu for the Bagnold formula at juj ¼ ue . 4. Swash simulation The identical swash event (i.e. identical PW01 initial conditions; subsequently flow will differ as the bed changes) is now examined with the above six different sediment transport formulae. Although each formula yields the same dimensionless form, (5)–(7), the nondimensional bed mobility parameter r is calculated differently in general (see Table 1) so that for different formulae identical values of r correspond in general to different sediment transport rates and mobilities. Therefore, the approach of the present work is to examine the same swash event and to normalise results so they are intercomparable. Note that such an approach reveals the differences due to the different dependences expressed in q ¼ qðu; hÞ rather than differences due to different values of Ai (see Table 1). Note also that the uncoupled bed changes are derived according to the analytical solution for the dam-break problem over a sloping Table 1 Expressions for six sediment transport formulae and corresponding variables.
^ q
(I) Grass
(II) Bagnold
(III) Meyer-Peter–Müller
^3 A1 u
^ 2e Þ ^ ðu ^2 u A2 u
^2 u ^ 2e Þ3=2 A3 ðu
3
2
u2e Þ 3=2
q
u
uðu
q0
A1 ðgh0 Þ3=2 0 3u2 nA1 g
A2 ðgh0 Þ 0 3u2 u2e nA2 g
A3 ðgh0 Þ3=2 0
(IV) Van Rijn
(V) Bailard
(VI) ^u ^ ¼ Ah ^ 3 (PH) q
^ q
^ ju ^ j2:4 A4 u
^ ju ^ j3 A5 u
^u ^3 A6 h
q
ujuj2:4
ujuj3
qh qu
r
q0 qh qu
r
A4 ðgh0 Þ 0
1:7
2:4
3:4juj
3jujðu2 u2e Þ1=2 nA3 g
hu 2
3
A6 h0 ðgh0 Þ3=2 u3
A5 ðgh0 Þ 0 4juj
0:2
nA4 g 1:2 h0
ðu2 u2e Þ3=2
3
Rt fixed bed [1]. From the bed evolution equation, ½Bttdi ¼ tid rqx dx (ti inundation time and td denudation time), and ti ; t d and the uncoupled DB ¼ ½Bttdi can be calculated with the known h and u values; see also [3]. 4.1. Initial and boundary conditions 4.1.1. Initial conditions The initial conditions are shown in Fig. 1. The beach is of a uniform slope a ¼ 0:1, so the initial bed elevation Bðx; t ¼ 0Þ ¼ ax. The dam is located at x ¼ 0, with water on the left side only. The water depth on the left side of the dam is hðx < 0; t ¼ 0Þ ¼ 1, and uðx < 0; t ¼ 0Þ ¼ 0. The dam is assumed to collapse at t ¼ 0. In this wet-dry dam-break problem, the discontinuity in time of velocity, i.e. between the initial conditions and the solution just after collapse (i.e. for 0 < t ¼ t 1), prohibits the use of the initial conditions for t ¼ 0 in the MOC STI solver [3]. The usual approach to circumvent this problem has been to initialise with the corresponding (initially flat bed) Riemann solution at t ¼ t [3]. However, each sediment transport formula in general corresponds to a Riemann solution, the derivation of which can be complicated. Therefore, it is desirable instead to use the (analytical) PW01 solution (here at t ¼ t ¼ 0:01; Kelly and Dodd [3] found that t ¼ 0:1 was sufficient in the Riemann problem) ) Bðx; t Þ ¼ Bðx; t ¼ 0Þ ¼ ax. Comparisons of solutions using this approach and that obtained using the Riemann wave initial data for the I and IV reveal very close agreement: see Appendix A. 4.1.2. Boundary conditions In the present swash simulation, the seaward boundary is chosen at a point far enough away to ensure that the flow is uninfluenced by the reflected wave from the land so that the water depth h and bed level B at that point remain unchanged through the computation time. However, as the flow is dominated by gravity, u at that point decreases as time increases. Here, the seaward boundary is chosen at x ¼ 250, where hð250; tÞ ¼ 1 and Bð250; tÞ ¼ 250a. As uð250; 0Þ ¼ 0, then uð250; tÞ ¼ at. The landward boundary when t > 0 is the wet-dry boundary. 4.2. Wet-dry boundary treatment pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi From (9), k1 ¼ u=2 u2 =4 þ rqu ; k2 ¼ u=2 þ u2 =4 þ rqu , and k3 ¼ u at the shoreline. Therefore, when qu > 0, k1 < u (hence k1 is also smaller than 0) and k2 > u, so that the k2 characteristic extends from the flow interior to the shoreline, and the k3 characteristic is identical to the shoreline. In the uprush for formulae I–V, we can combine the Riemann equation along the k2 characteristic with the shock relation between bed levels on both sides of the sediment bore at the tip, see [3] so that we have two equations for two unknowns at the tip (us and Bs ). The shoreline position, xs , also unknown initially, is earlier found from the shoreline relation, which can be iterated to a desired accuracy. In the backwash, as the bed level on the right side of the sediment bore at the tip is unknown, the shock relation between bed levels is not used directly to determine dependent variables of the flow region; instead for formulae I–V extrapolation of the u and B values in the flow region is used to obtain us and Bs in the backwash, with the shock relation used to obtain B on the dry beach as the shoreline retreats. Extrapolation is used to solve the wet-dry boundary for formula VI both in uprush and backwash. 4.3. Determination of the dimensionless bed mobility parameter
2
1=2
nA5 g 3=2 h0
3hu nA6 gh0
Because of the differing dimensions of Ai (and the different dependencies of q on u and h) choosing a single r value for all
15
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
Table 2 Normalisation based on identical net sediment flux at x ¼ 0 for different sediment formulae. Values of r, net offshore sediment flux at x ¼ 0; bed variation at x ¼ 0 and maximum inundation are also illustrated. (I) Grass
(II) Bagnold
(III) Meyer-Peter–Müller
(IV) Van Rijn
(V) Bailard
(VI) PH
Net offshore sediment flux at x ¼ 0
0.0100 0.4457
0.01176 0.4461
0.0126 0.4457
0.00888 0.4458
0.00693 0.4460
0.6200 0.4460
DB at x ¼ 0 Max. inundation
0.0601 17.33
0.0604 17.13
0.0615 17.06
0.0639 17.12
0.0707 16.96
0.1984 26.34
r
formulae yields very different amounts of sediment movement (dimensional and dimensionless), so in principle each formula ought to be calibrated to measured data. Here we choose r for I, roughly based on the previously calibrated value from [3], and then calculate all other values by normalising simulations for the other formulae such that we obtain the same net sediment flux under the PW01 event at x ¼ 0 (see Table 2). Note that other normalisations are possible, e.g. identical bed change at x ¼ 0, but yield a similar overall picture. Note also that we take r ¼ 0:01 for the formula I, which is lower than that of [3] (r ¼ 0:0654); this is to ensure the ‘‘equivalent’’ r for VI does not lead to unrealistically large local bed changes for the large resulting bed mobilities. The value of ue in formulae II and III is set to 0.45, which is higher than the range of values in [17], to examine the effect of the sediment motion threshold. Note that we do not relate r to specific grain sizes. The purpose here is to examine lesser or greater beach mobilities; these can, with a large degree of uncertainty, be related back to grain sizes. 4.4. Beachface evolution For all six simulations, at the moment of dam collapse there is a sudden acceleration in water velocity at the tip. After the collapse, the flow is dominated by gravity, and the velocity gradually decreases in time for all x. Maximum inundation (and run-up) is achieved where uðxs Þ ¼ 0, and in the backwash the tip subsequently recedes and offshore velocity increases. The sediment is moved onshore during run-up and offshore during backwash, with associated bed deformation. This leads to a convergence of k3 characteristics in the lower swash, resulting in a k3 shock in all six simulations. Overall, the flow is qualitatively similar for all simulations; we show the evolution of h and u obtained using I in Fig. 2(a) and (b), which is qualitatively similar for formulae I–V. Bed change is also similar if equations I–V are used, hence we show the dimensionless bed change for I and VI only, in Fig. 2(c) and (d). The most striking difference between (c) and (d) is the much larger inundation for formula VI (see Table 2), which results from a larger onshore maximum velocity for that formula at the moment of dam collapse (for VI about 1.47 times that for I). It is also found that the magnitude of instantaneous sediment deposition / erosion in the upper swash for VI is considerably smaller than that for I–V, which is because the water depth h in the upper swash for VI is considerably reduced. For all formulae sediment is initially moved onshore, causing deposition in the upper swash in the uprush with sediment coming from the base of the swash, only for the backwash subsequently to transport that and more sediment offshore. This is the result of the asymmetry in the PW01 swash velocity, see [2]. Also shown (Fig. 2(e) and (f)) are dimensionless bed change (DB) contours in the equivalent uncoupled simulations for I and VI. The pattern of erosion and deposition for formula I is very similar in coupled and uncoupled simulations (cf. Fig. 2(c) and (e)); uncoupled and coupled DB also show some similarities for VI (cf. Fig. 2(d) and (f)). Only at the base of the swash (x 0) do we see
substantial differences, in part due to the dam collapse and the aforementioned morphodynamic backwash bore formation in both coupled simulations. Additionally, we show h and u for VI in Fig. 3, which is particularly useful for comparison with equivalent simulations including bed shear stress; see Section 4.8. The shoreline trajectories of all formulae (I–VI) in fully coupled simulations, as well as the uncoupled one, are shown in Fig. 4. Fig. 4 also shows the shoreline trajectories in simulations including bed shear stress with constant C D ¼ 1 103 in uprush and backwash for I and VI, and they will be analyzed in Section 4.8. The significantly increased inundation and swash period for the formula VI (also with respect to the PW01 solution) are apparent in Fig. 4. However, those for q ¼ qðuÞ (I–V) are very close to each other and considerably reduced from the PW01 solution. This is consistent with the formation of a sediment bore at the shoreline for I– V, which does not occur for VI. This bore transports a large amount of sediment, even though h ! 0 there, and this is one of the main objections to use of formulae I–V for making engineering predictions in vanishing water depth. The existence (non-existence) of the bore can also be seen in Fig. 2(c) and (d) in the discontinuity (continuity) in bed contours when moving from wet to dry regions. 4.5. Final beach profile The final bed changes when the shoreline has receded back to the initial shoreline position during one single swash event are collected in Fig. 5 for all formulae. The bed changes of the equivalent uncoupled simulations are also shown. Final profiles for all q ¼ qðuÞ formulae indicate net offshore sediment flux, consistent with [2] with, generally, net erosion increasing offshore, but also with less net erosion than for uncoupled simulations (Fig. 5) see [3] over the whole swash. The highest power in u, yields the largest (smallest) erosion in the lower (upper) swash and smallest maximum inundation (as can be observed from the coloured squares, indicating the maximum runup of each coupled simulation). This is because j u j in the lower swash is generally larger than that in the upper swash, therefore an increase in the power of u gives an increase (decrease) in the sediment transport in the lower (upper) swash, which results in more (less) erosion in the lower (upper) swash. Profiles are generally similar, however. Note that this switch between upper and lower swash can be observed in both coupled and uncoupled simulations, indicating that only at the base of the swash, where the dam collapses and bores are formed, are the basic morphodynamics changed by full coupling. This implies that some success may be obtained with uncoupled models of beach change, as long as probable inaccuracies are understood. For a quasi-coupled approach to swash modelling, see [19] (especially Fig. A.12 of that paper). For both formulae II and III (those with thresholds) there is a region in which the bed level does not change (when j u j< ue ), and the width of this region under the PW01 swash event is u2e =ð2aÞ ¼ 1. The presence of ue reduces the sediment movement (thus Bagnold (II) and Meyer-Peter–Müller (III) formulae have a correspondingly larger r than the Grass formula (I) to ensure the
16
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
3
Fig. 2. Contour plots of flow and bed change for formulae I (q ¼ u3 ) (r ¼ 0:01) and VI (q ¼ hu ) (r ¼ 0:62). (a) h for I; (b) u for I. Note that space and time axes are normalised by maximum inundation and swash period respectively. (c) DB for I; (d) DB for VI; (e) DB for I (uncoupled); (f) DB for VI (uncoupled).
same sediment flux). However, the profiles differ little from that for I. 3 For formula VI (q ¼ hu ), as expected there is relatively little bed change in the upper swash, with significantly more (net erosion) in the lower swash, where h increases. Note that for VI more erosion is predicted in the mid- to lower swash by the fully coupled simulation than by the uncoupled one: see Fig. 5.
4.6. Net erosion and deposition: influence of power of u and h The zero contour in the dry bed region near the maximum inundation in Fig. 2(d) seems to indicate that net deposition occurs in the extreme upper swash using VI. In fact (see also Appendix B) these beach change values are smaller than the truncation error in the scheme so that deposition cannot be predicted definitively.
17
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
3
Fig. 3. Contour plots of flow for formula VI (q ¼ hu ) (r ¼ 0:62). (a) h; (b) u.
45 40 35 3
30
q=u q=u(u2−u2e)
25
q=(u2−u2)3/2
20
q=u|u|
t
e
2.4
q=u|u|3 15
q=hu3 uncoupled q=u3 with shear stress q=hu3 with shear stress
10 5 0
0
5
10
15
20
25
30
35
x
Fig. 4. Dimensionless shoreline trajectories for formulae (I–VI) along with that for PW01 (i.e. uncoupled). The equivalent shoreline trajectories in simulations with bed shear stress (constant C D ¼ 1 103 in uprush and backwash) for formulae I and VI are also illustrated.
It is therefore instructive to examine the dependence of q on h. We do so in Fig. 6, in which equivalent final profiles for normalised 0:8 1:2 simulations for three other formulae, q ¼ h u3 ; q ¼ h u3 , and 1:3 3 q ¼ h u , as well as those for a selection of powers of u (i.e. those without thresholds) are shown. In general, the relative deposition (erosion) in the upper (lower) swash increases as the power of either h or u increases; of these two factors h is the dominant one. Importantly, a power of h that is > 1 results in net deposition in the upper swash, and, in fact, the development of a small berm or swash bar there. Note that this berm is a feature of the upper swash and is not connected with the formation of a backwash bore, as found for the more generic backwash bore bars by Zhu et al. [15]. It looks qualitatively similar to the pattern of implied deposition obtained by Pritchard and Hogg [2] when a settling lag was introduced. It can be seen that (Fig. 6) the formulae qðh; uÞ predict significantly more erosion than qðuÞ in the lower swash based on the same net flux normalisation. This is due to the much larger r values we imposed on these formulae in order to achieve the same net flux.
0
4.7. Uncoupled and coupled modelling
Bed change Δ B
−0.02 −0.04
q=u3 2 2 q=u(u −ue)
−0.06
q=(u2−u2)3/2
−0.08
q=u|u|
−0.1
q=u|u|
e
In Fig. 7 following [2] we illustrate differences in instantaneous sediment transport (rq) between coupled (in black) and uncoupled
2.4 3
0
3
q=hu −0.12
−0.05
3
q=u 2 2 q=u(u −u )
−0.14
e 2 3/2 e
2
−0.1
−0.16
q=u|u|2.4 −0.18
3
q=u|u| 3
q=hu
−0.2 0
5
10
15
20
25
x
Fig. 5. Dimensionless change in bed level relative to the initially plane beach after one single swash. Solid lines: fully coupled simulations. Dashed line: uncoupled simulations. Indicates position of maximum inundation for the coupled events (for the PW01 (uncoupled) event it is at x ¼ 20). Colours indicate formula – see legend.
Bed change Δ B
q=(u −u )
−0.15
0.01
−0.2 0 −0.25
q=u3 2.4
−0.01
−0.3
q=u|u| 4
6
8
10
12
14
16
q=u|u|3
18
0.8 3
q=h
−0.35
u
3
q=hu
1.2 3
q=h
−0.4
u
q=h1.3u3 −0.45
0
5
10
15
20
25
30
x
The values are, however, indicative of a very small net beach change in this region. It appears, therefore, that formula VI is at the boundary between net erosion and deposition.
Fig. 6. Dimensionless change in bed level relative to the initially plane beach after one single swash for different powers of h and u. Indicates position of maximum inundation of simulations represented by solid lines, and indicates that of simulations represented by dashed lines. Colours indicate formula – see legend.
18
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
0.1 0.1 0.08 0.06
σy
σy
0.05
0
0.04 0.02 0
−0.05
−0.02 −0.04
−0.1
0
10
20
30
40
t (q=u3, x=0.1,5,10,15)
−0.06
0
10
20
30
40
t (q=hu3,x=0.1,8,16,24)
Fig. 7. Dimensionless instantaneous qðxÞ at a number of locations in the swash. Present simulations (black) and uncoupled (PW01) ones (red). The simulations with bed shear stress (constant C D ¼ 1 103 in uprush and backwash) are represented by blue dashed lines. (Left: I. Right: VI). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this paper.)
(in red) approaches for I and VI. The instantaneous sediment flux at different positions in simulations including bed shear stress with constant C D ¼ 1 103 in uprush and backwash for I and VI are also illustrated in Fig. 7, and they will be discussed in Section 4.8. For I (left panel) differences are apparent even for small times, due to the different speeds of the uprush tip (related to different maximum inundation). Note, however, two points: (i) Instantaneous rq values for coupled and uncoupled calculations are similar where both are non-zero; the major discrepancy is the times of inundation (t i ) and denudation (td ), and the discrepancy in td is particularly large. (ii) In the backwash rq in the coupled simulations shows slightly less offshore transport than its uncoupled equivalent. What this tells us is that the velocity, from which q is derived, is similar, but that the smaller initial velocity, which can be obtained by the Riemann solution for the wet-dry dam-break problem over a flat mobile bed ([14]), results in a less exaggerated net erosion. For VI there is no such clear correspondence. In the top panels of Fig. 8, also following [2], we plot the aforeRt mentioned net sediment flux (Q ðxÞ ¼ r tid qðx; tÞ dt) in the uprush, backwash, and over the whole event, comparing it to that obtained from the uncoupled approach. For I (indicative of I–V recall) differences between the approaches do appear to accumulate over the swash event, with relatively little discrepancy in the uprush, but significant disagreement in the backwash, which results in the main differences (see also the left bottom panel in Fig. 8, in which the difference in net sediment flux between coupled and uncoupled approaches is shown). For VI significant differences exist over all phases of the swash. 4.8. Simulations with bed shear stress 4.8.1. Numerical construction We describe the bed shear stress by the Chezy law s ¼ C D u2 ([17]) with s dimensionless bed shear stress and C D drag coefficient. When bed shear stress is included in the coupled system, the momentum Eq. (6) is changed and it becomes:
ut þ uux þ hx þ Bx ¼
C D juju ; h
ð11Þ
(5) and (7) remain the same. Consequently, the Riemann Eq. (10) become:
RðkÞ ¼ kk
du kk þ rqh dh dB C D juju þ þ ¼ kk ; dt h kk u dt dt
k ¼ 1; 2; 3;
ð12Þ
where kk is obtained from (9). The swash flow is computed by solving (12) numerically [14]. Note that the term C Dhjuju in (11) ! 1 at the shoreline, i.e. as h ! 0, in the uprush. In the flat fixed bed dam-break solution with bed shear stress [12,13], this term is balanced by hx , which implies that hx ! 1 at the wave tip. For the uprush of swash flow over a sloping beach, this balance is also appropriate, see also [20]. However, in the backwash C Dhjuju cannot be balanced by hx , as C Dhjuju ! 1 when u is negative, while hx cannot ! 1. Therefore u must be 0 such that (11) also holds at the shoreline in the backwash. This implies that there is no receding shoreline, see [20], for a detailed asymptotic analysis. Due to these considerations the numerical implementations for the shoreline solution in the uprush and backwash are different. Linear extrapolation is employed to solve the tip in the uprush, and according to the analysis of the shoreline in the backwash, us ¼ 0 is imposed in the numerical solution. Note that the shock conditions for the system with bed shear stress are identical to those without shear stress [3,15].
4.8.2. Constant C D in uprush and backwash Fully coupled simulations with sediment transport formulae I 3 (q ¼ u3 ) and VI (q ¼ hu ), representing qðuÞ and qðh; uÞ respectively, are considered with bed shear stress included. The numerical simulations and comparisons confirm the analysis in the uprush hx ! 1 near the shoreline, causing the swash flow in the uprush to be much deeper than that with no bed shear stress. The water velocity is considerably reduced in the uprush, and the maximum inundations and instantaneous sediment fluxes are accordingly reduced for both I and VI although with deeper flow near the shoreline, see Figs. 4 and 7. The backwash flow becomes very thin and slow due to the resistance of bed shear stress, which results in 3 much less offshore sediment transport, especially for q ¼ hu . As the shoreline does not recede, the swash zone is wet during the backwash and there is no denudation time unless a finite minimum allowable depth is chosen; however, the sediment flux gradually approaches 0 as the flow depth and velocity approach 0, see Fig. 7. The bed changes at t ¼ 40 with different C D values for formulae I (r ¼ 0:01) and VI (r ¼ 0:62) are shown in Fig. 9. Note that the bed in the swash zone is still wet at t ¼ 40 as when C D > 0, the shoreline does not recede seawards (us ¼ 0) in the backwash. Note also that we retain r values from Table 2 here. In reality C D and r will
19
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
0.4 0.1 0.2
0
0
Q
Q
−0.1 −0.2
−0.2
−0.3 −0.4
−0.4
−0.5 0
5
10
15
20
0
5
x (q=u3)
20
25
20
25
0.1 uprush backwash net
0.05
Q−Qun
0.1
Q−Qun
15
x (q=hu3)
0.15
0.05
0 −0.05
0
−0.05
10
−0.1
0
5
10
15
20
0
5
3
10
15 3
x (q=u )
x (q=hu )
Fig. 8. Top: Dimensionless Q ðxÞ (net sediment flux) for uprush and backwash (dashed lines) and the whole swash event (solid lines) for present simulations (black) and uncoupled (PW01) ones (red). Note that the uprush (backwash) here is defined as the duration until (after) maximum inundation is achieved. Bottom: difference between Q from present simulations and that from uncoupled approach. (Left: I. Right: VI). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this paper.) 0.01
be related, but as we do not know the relation it is appropriate to test different C D values for the given (normalised) r. It is seen from Fig. 9 that the maximum inundations for both formulae are much reduced from those with no shear stress, i.e. C D ¼ 0 (see also Fig. 4); as C D increases, the maximum inundation decreases. Fig. 4 also shows the non-receding shoreline in the backwash for both formulae. Note that for VI, the maximum inundation is significantly reduced when bed shear stress is included even with a small C D : from x ¼ 26:34 (off scale on figure) for C D ¼ 0 to x ¼ 14:61 for C D ¼ 2 104 (Fig. 9). One important feature of the non-frictional simulation with VI is the wide and thin
0 −0.01 −0.02
ΔB
−0.03 −0.04 −0.05
C =0 D
C =2×10−4
−0.06
D
CD=1×10−3
−0.07
−2
CD=1×10
−2
−0.08 −0.09
CD=2×10 0
2
4
6
8
10
12
14
16
18
20
x 0.02
0.02
0.01
−0.02
0
−0.04
−0.01
−0.06
−0.02
ΔB
ΔB
0
−0.08 −0.1
C =2×10 , C
=C /2
−0.04
C =1×10−3, C
=C /2
−0.05
CD=1×10 , CDb=CD/2 −4 C =2×10 , C =C
−4
D
Db
D
Db
D D
−2
−0.12 C =0 D
−0.14
C =2×10−4
−0.16
C =1×10−3
D
−0.06
D
−0.07
C =1×10 , C D
CD=1×10 2
4
6
8
10
12
14
16
18
−0.08
20
0
2
4
6
8
10
12
D
=C
Db
14
D
=C
Db
−2
−2
0
Db
C =1×10−3, C D
D
−0.18 −0.2
−0.03
D
16
x
x 3
3
Fig. 9. Bed changes for formulae I (q ¼ u ) and VI (q ¼ hu ) with bed shear stress at t ¼ 40 (constant C D in uprush and backwash). Top: I and bottom: VI.
Fig. 10. Comparison of bed changes for formula I (q ¼ u3 ) with varying C D and r values and constant values in the uprush and backwash (C Db ¼ C D =2 and rb ¼ 0:0035) at t ¼ 40. Solid line: varying C D and r; dashed line: constant C D and r.
20
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22 0.05
0
ΔB
−0.05 −4
C =2×10 , C D
=C /2
Db
−3
C =1×10 , C
−0.1
D
D
=C /2
Db
D
−2
CD=1×10 , CDb=CD/2 −4 CD=2×10 , CDb=CD
−0.15
−3
C =1×10 , C
=C
C =1×10−2, C
=C
D D
−0.2
0
5
10
Db Db
D D
15
x 3
Fig. 11. Comparison of bed changes for formula VI (q ¼ hu ) with varying C D and r values and constant values in the uprush and backwash (C Db ¼ C D =2 and rb ¼ 0:217) at t ¼ 40. Solid line: varying C D and r; dashed line: constant C D and r.
swash flow behind the shoreline, see Fig. 3(a), while the simulation with bed shear stress yields a much deeper tip flow. The bed shear stress confines the wide and thin tip flow to a narrower region, and therefore the flow is much deeper and the maximum inundation is greatly reduced.
The net sediment flux at x ¼ 0 is still offshore for both formulae for all tested non-zero C D values, however, the magnitude of net erosion is considerably reduced. For formula I, net erosion is predicted all over the swash zone to a reduced extent for C D ¼ 0; 2 104 ; 1 103 and 1 102 , see Fig. 9. However, net erosion in the upper swash approaches 0 as C D increases and net deposition is predicted in the region around x ¼ 3 for C D ¼ 2 102 . Bed shear stress reduces the velocity in the uprush and backwash; that in the backwash is to a greater extent, such that deposition occurs for large C D values. Note that deposition occurs in the upper swash for formula VI for all tested C D values. The sediment flux in this formula is dependent on h and u, and the greatly reduced h and u in the upper swash during the backwash period result in significantly less offshore sediment transport in this region. The combined effect of reduced onshore and offshore sediment transport is net deposition in the upper swash and erosion in the lower swash for VI, and the deposition in the upper swash is much more than that occured in the simulation with formula I. 4.8.3. Varying C D and r in uprush and backwash So far, we have assumed that the drag coefficient C D has the same value in the uprush and backwash. However, in-situ measurements, e.g. [21,22] indicate that if bed shear stress is represented by a drag law the drag coefficient in the uprush is roughly
Fig. A.12. Comparison using PW01 solution as initial data (black) with that using Riemann solution over flat mobile bed (red) for formula I (q ¼ u3 ) (r ¼ 0:01) and IV (q ¼ ujuj2:4 ) (r ¼ 0:00888). (a) h, (b) u, (c) DB, all for I; (d) DB for IV. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this paper.)
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
twice that in the backwash. Thus, it is instructive to examine the effect of a varying drag coefficient in the uprush and backwash on the bed change pattern. We set the drag coefficient in the backwash as: C Db ¼ 0:5C D . Note that the bed mobility is also changed because a smaller C D implies that a smaller shear stress exerted on the bed which will, accordingly, be less mobile when subject 3=2 r as q / s3=2 , which exto the same u; here we take rb ¼ CCDDb cludes any threshold of movement [17]. When C D is decreased in the backwash, the water recedes more quickly due to the reduced resistance due to reduced bed shear stress. Thus, more sediment is moved seawards, and more erosion could be predicted compared to a constant C D value in uprush and backwash. In contrast a decreased r in the backwash results in a less mobile bed, and therefore less offshore sediment movement, which, over a swash event, has the potential to cause deposition. We examine the combined effect of both decreased C D and r values in the backwash by the simulations for formulae I (q ¼ u3 ) 3 and VI (q ¼ hu ) with varying C D and r values. The final bed changes for varying C D and r values for formula I are shown in Fig. 10. The combined effect of varying C D and r values is considerably less erosion, and net deposition is predicted in the upper swash for all simulations and even in the middle swash when C D ¼ 0:01. The comparison of results for formula VI in Fig. 11 also shows considerably more deposition in the upper and middle swash when C D and r are decreased in the backwash. Note that when C D ¼ 0:01, the net sediment flux at x ¼ 0 is onshore; note that the ratio rb =r ¼ 0:35 is larger than that suggested by Pritchard and Hogg [2] (< 0:3) for net onshore flux to occur in the uncoupled simulation with no bed shear stress. This further suggests that the inclusion of bed shear stress is important for the swash simulation. 5. Conclusions The numerical simulations (based on the STI MOC method) of the fully coupled PW01 event with six different sediment transport formulae demonstrate that the results of [2] are robust in that net erosion over the whole swash is predicted, consistent with the coupled results. Generally, less erosion is predicted using the fully coupled approach. This is primarily because the smaller initial swash tip velocity (for sediment transport formulae q ¼ qðuÞ) results in a smaller swash period and therefore a less exaggerated velocity asymmetry. Inclusion of a threshold of movement or differing powers of u makes little difference to this. Note that the values of r (bed mobility) used here for q ¼ qðuÞ (0.007–0.013) equate to significantly mobile beaches. Inclusion of water depth in the expression for sediment transport, thought to be appropriate in very small water depths ([23]), leads to correspondingly small sediment transport and therefore also small beach change in the upper swash. Because for q ¼ qðh; uÞ the water depth is now a controlling factor, and because the PW01 swash asymmetry in h yields larger depths in the uprush, (also apparent in fully coupled simulations: see Fig. 3(a)), the resulting net bed change in the upper swash depends on two competing factors: h and u asymmetry in the PW01 event (both of which are diminished in importance in the upper swash). The result is little net change in the upper swash. Further, if we examine this h dependence by taking n q / h for n > 1, the result is net deposition in this region. Note that use of such a formula, even for large depths, implies that equilibrium concentration is proportional to depth, albeit weakly. The differences between sediment formulae I–V are small, but between I–V and VI are large. In general, therefore, differences between I–V are negligible compared to those due to differing swash asymmetries (based on [10], Fig. 11 (lowest panels, for no pre-suspended sediment – there is a settling lag here, unlike in our
21
simulations, but the E ¼ 0:3 case is closest to the present work)). Differences for VI are significant, and probably more so than those 3 in [10], but it is also apparent that the q ¼ hu formula is more appropriate for swash flows with small depths. Note also, however, that some sediment is entrained in the present simulations seaward of x ¼ 0 (i.e. the surf zone for our purposes). Furthermore, although pre-suspended load is likely to be one of the dominant effects, an event without pre-suspended load has the potential to entrain a lot of sediment and, depending on phasing of the next swash event and the grain size, this could be available as pre-suspended load for the next event. The amount of that sediment will depend on q. The influence of bed shear stress, however, is very important; its inclusion can be the difference between net erosion and deposition. Note also that this finding is not just significant for the swash itself, but also for swash modelling, which may or may not be a true representation of actual swash and beach evolution. Note that the crude Chezy law bed shear stress used here is also used by almost all engineering models in use at present, and gives a reasonable representation of bed shear stress in the swash. Simulations with bed shear stress, which becomes more effective in the flow with shallow depth, can yield qualitatively different net bed change. Bed shear stress greatly reduces the velocity, especially that in the backwash, resulting generally in less sediment transport and particularly less offshore sediment transport. Although the net sediment flux at x ¼ 0 is still offshore with constant C D in uprush and backwash, the beach change profile inside the swash zone is considerably changed for both formulae I 3 (q ¼ u3 ) and VI (q ¼ hu ), and is also closely related to the drag coefficient C D . Deposition occurs for VI for all tested non-zero C D values in the upper swash, and deposition occurs for I in the middle swash only when a relatively large C D value is imposed. Note also that formula VI can be interpreted as the flux of suspended load in equilibrium state. Although the transport of suspended sediment under unsteady swash flow is usually not in an equilibrium state, the greatly altered beach change pattern for formula VI when bed shear stress is included, indicates the possibly great effect of bed shear stress on the transport of suspended load. Simulations of decreasing C D and r values in the backwash predicts much less offshore sediment transport for both formulae I 3 (q ¼ u3 ) and VI (q ¼ hu ), and net deposition is generally predicted in the upper swash for all simulations. Net onshore sediment flux is 3 predicted when C D ¼ 0:01 for formula VI (q ¼ hu ) when decreased C D and r values are implemented in the backwash, and the ratio rb =r ¼ 0:35 is larger than that suggested by [2]. This is very indicative for the great importance of different physical mechanisms in uprush and backwash and also the inclusion of bed shear stress for simulations in the swash zone. The present work is limited and idealised in the assumptions of depth-averaged flow and the exclusion of a number of effects that, depending on circumstances can be important: in/exfiltration, detailed boundary layer, downslope effect (avalanching), air entrainment, turbulence, non-hydrostatic effects and settling lag of suspended sediment. Furthermore, the bed mobilities are not directly calibrated from field data but based on the same sediment flux normalisation. Finally, use of the PW01 solution as an initial condition for a swash event yields almost identical results to the mathematically more correct Riemann solution initial condition. We therefore recommend such an approach in numerical simulations as being much simpler and valid for all q. Acknowledgements The authors would like to express their gratitude to the China Scholarship Council, and to the International Office of The University of Nottingham for providing financial support.
22
F. Zhu, N. Dodd / Advances in Water Resources 53 (2013) 12–22
Appendix A. Verification The initial data used for the swash simulation are from the PW01 solution at time t > 0. Comparison of the results for the Grass (I) and Van Rijn formula (IV) when using the PW01 solution as initial data with those obtained using the Riemann solution over a flat mobile bed [14] is shown in Fig. A.12. Both comparisons show very close agreement, therefore the PW01 solution at a finite time can be generally used for swash simulation over a mobile bed. Appendix B. Determining net erosion/deposition for small net bed changes To determine the significance of predictions of small values of net bed change over a swash event it is useful to examine a (non-physical) sediment transport rate q ¼ Ahu ); q / water discharge, so that
Z
td
qx dt ¼ 0
ð13Þ
ti
at any location in the swash that is wetted and then dried (where t i and t d are times of inundation and denudation respectively). This therefore provides a useful measure of the error associated with numerical solution of (10) over a swash cycle. Note that taking q ¼ hu implies that (9) degenerates to a quadratic, and the eigenffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi values are: k ¼ 0; u ð1 þ rÞh. References [1] Peregrine DH, Williams SM. Swash overtopping a truncated beach. J Fluid Mech 2001;440:391–9. [2] Pritchard D, Hogg AJ. On the transport of suspended sediment by a swash event on a plane beach. Coastal Eng 2005;52:1–23.
[3] Kelly DM, Dodd N. Beach-face evolution in the swash zone. J Fluid Mech 2010;661:316–40. [4] Whitham GB. On the propagation of shock waves through regions of nonuniform area or flow. J Fluid Mech 1958;4:337–60. [5] Keller H, Levine DA, Whitham GB. Motion of a bore over a sloping beach. J Fluid Mech 1960;7:302–16. [6] Shen MC, Meyer RE. Climb of a bore on a beach. Part 3. Run-up. J Fluid Mech 1963;16:113–25. [7] Coco G, Huntley DA, O’Hare TJ. Investigation into a self-organization model for beach cusp formation. J Geophys Res 2000;105(C9):21991–2002. [8] Bryan KR, Coco G. Detecting nonlinearity in run-up on a natural beach. Nonlin Processes Geophys 2007;14:385–93. [9] Guard PA, Baldock TE. The influence of seaward boundary conditions on swash zone hydrodynamics. Coastal Eng 2007;54:321–31. [10] Pritchard D. Sediment transport under a swash event: the effect of boundary conditions. Coastal Eng 2009;56:970–81. [11] Briganti R, Dodd N, Kelly DM, Pokrajac D. An efficient and flexible solver for the simulation of the morphodynamics of fast evolving flows on coarse sediment beaches. Int J Numer Methods Fluids 2012;69:859–77. [12] Dressler RF. Hydraulic resistance effect upon the dam-break functions. J Res Natl Bur Stand 1952;49:217–25. [13] Whitham GB. The effects of hydraulic resistance in the dam-break problem. Proc R Soc London, A 1955;227:399–407. [14] Kelly DM, Dodd N. Floating grid characteristics method for unsteady flow over a mobile bed. Comput Fluids 2009;38:899–909. [15] Zhu F, Dodd N, Briganti R. Impact of a uniform bore on an erodible beach. Coastal Eng 2012;60:326–33. [16] Grass A. Sediment transport by waves and currents. Technical Report FL29, University College London, London Centre for Marine Technology; 1981. [17] Soulsby RL. Dynamics of marine sands. London: Thomas Telford; 1997. [18] Yalin M. Mechanics of sediment transport. 2nd ed. Pergamon; 1977. [19] Postacchini M, Brocchini M, Mancinelli A, Landon M. A multi-purpose, intrawave, shallow water hydro-morphodynamic solver. Adv Water Resour 2012;38:13–26. [20] Antuono M, Soldini L, Brocchini M. On the role of the Chezy frictional term near the shoreline. Theor Comput Fluid Dyn 2012;26:105–16. [21] Masselink G, Hughes M. Field investigation of sediment transport in the swash zone. Cont Shelf Res 1998;18:1179–99. [22] Barnes MP, O’Donoghue T, Alsina JM, Baldock TE. Direct bed shear stress measurements in bore-driven swash. Coastal Eng 2009;56:853–67. [23] Baldock TE, Alsina JM. Discussion of ‘‘On the transport of suspended sediment by a swash event on a plane beach’’, by D. Pritchard and A.J. Hogg. Coastal Eng 2005;52:811–4.