Revisiting slope influence in turbulent bedload transport: consequences for vertical flow structure and transport rate scaling

TL;DR

This study revises the understanding of slope effects in turbulent bedload transport, showing that slope influences the entire transport process and the vertical flow structure, validated by simulations and leading to a new transport rate scaling law.

Contribution

It corrects the classical slope correction of the Shields number and demonstrates that slope impacts the entire transport formula and flow structure, providing a new unified scaling law.

Findings

Classical slope correction of Shields number is based on an erroneous buoyancy formulation.

Slope influences the entire transport process and flow structure.

Transport rate scales with a rescaled Shields number, collapsing data onto a master curve.

Abstract

Gravity-driven turbulent bedload transport has been extensively studied over the past century in regard to its importance for Earth surface processes such as natural riverbed morphological evolution. In the present contribution, the influence of the longitudinal channel inclination angle on gravity-driven turbulent bedload transport is studied in an idealised framework considering steady and uniform flow conditions. From an analytical analysis based on the two-phase continuous equations, it is shown that : (i) the classical slope correction of the critical Shields number is based on an erroneous formulation of the buoyancy force, (ii) the influence of the slope is not restricted to the critical Shields number but affects the whole transport formula and (iii) pressure-driven and gravity-driven turbulent bedload transport are not equivalent from the slope influence standpoint. Analysing further the granular flow driving mechanisms, the longitudinal slope is shown to not only influence the fluid bed shear stress and the resistance of the granular bed, but also to affect the fluid flow inside the granular bed - responsible for the transition from bedload transport to debris flow. The relative influence of these coupled mechanisms allows us to understand the evolution of the vertical structure of the granular flow and to predict the transport rate scaling law as a function of a rescaled Shields number. The theoretical analysis is validated with coupled fluid-discrete element simulations of idealised gravity-driven turbulent bedload transport, performed over a wide range of Shields number values, density ratios and channel inclination angles. In particular, all the data are shown to collapse onto a master curve when considering the sediment transport rate as a function of the proposed rescaled Shields number.

Tables

Table 1. Table 1: Parameters of the simulations studied and symbol correspondence.

$\alpha$	${\rho }^p/{\rho }^f-1$	${\theta }^{\ast }$	$Re$	$Fr$	$h/d$	symbol
0.0126	0.75	0.11	35000	1.21	7.3	$▲$
0.0126	0.75	0.21	84000	1.24	13.0	$▲$
0.0126	0.75	0.34	175000	1.28	20.7	$▲$
0.0126	0.75	0.43	257000	1.31	26.3	$▲$
0.0126	0.75	0.52	356000	1.35	32.1	$▲$
0.0315	0.75	0.13	16000	1.66	3.5	$▲$
0.0315	0.75	0.22	34000	1.68	5.8	$▲$
0.0315	0.75	0.31	57000	1.75	8.0	$▲$
0.0315	0.75	0.4	85000	1.8	10.2	$▲$
0.0315	0.75	0.49	116000	1.85	12.3	$▲$
0.063	0.75	0.1	6000	2.01	1.7	$▲$
0.063	0.75	0.2	15000	2.13	2.8	$▲$
0.063	0.75	0.3	26000	2.24	4.0	$▲$
0.063	0.75	0.39	40000	2.32	5.2	$▲$
0.063	0.75	0.5	64000	2.36	7.0	$▲$
0.0191	1.5	0.11	62000	1.53	9.2	$■$
0.0191	1.5	0.21	156000	1.55	16.9	$■$
0.0191	1.5	0.3	279000	1.58	24.5	$■$
0.0191	1.5	0.49	605000	1.66	39.8	$■$
0.04775	1.5	0.12	28000	2.08	4.4	$■$
0.04775	1.5	0.21	58000	2.11	7.1	$■$
0.04775	1.5	0.28	94000	2.09	9.8	$■$
0.04775	1.5	0.41	159000	2.23	13.4	$■$
0.04775	1.5	0.52	234000	2.3	17.0	$■$
0.0955	1.5	0.11	11000	2.55	2.1	$■$
0.0955	1.5	0.2	27000	2.65	3.7	$■$
0.0955	1.5	0.3	49000	2.76	5.3	$■$
0.0955	1.5	0.39	72000	2.85	6.7	$■$
0.0955	1.5	0.49	102000	2.98	8.2	$■$
0.14325	1.5	0.11	8000	2.71	1.6	$■$
0.14325	1.5	0.21	19000	3.02	2.6	$■$
0.14325	1.5	0.31	36000	3.37	3.8	$■$
0.14325	1.5	0.39	53000	3.72	4.6	$■$
0.14325	1.5	0.48	79000	4.1	5.6	$■$
0.0253	3.0	0.11	139000	1.83	14.0	$\bullet$
0.0253	3.0	0.21	346000	1.85	25.5	$\bullet$
0.0253	3.0	0.31	615000	1.87	37.1	$\bullet$
0.0253	3.0	0.4	941000	1.91	48.6	$\bullet$
0.0253	3.0	0.5	1324000	1.95	60.2	$\bullet$
0.065	3.0	0.1	44000	2.5	5.3	$\bullet$
0.065	3.0	0.2	113000	2.52	9.8	$\bullet$
0.065	3.0	0.3	204000	2.55	14.4	$\bullet$
0.065	3.0	0.4	316000	2.62	19.0	$\bullet$
0.065	3.0	0.5	448000	2.68	23.6	$\bullet$
0.129	3.0	0.11	22000	3.08	2.9	$\bullet$
0.129	3.0	0.2	53000	3.13	5.2	$\bullet$
0.129	3.0	0.3	96000	3.22	7.5	$\bullet$
0.129	3.0	0.4	149000	3.32	9.9	$\bullet$
0.129	3.0	0.5	216000	3.45	12.3	$\bullet$
0.194	3.0	0.1	13000	3.26	2.0	$\bullet$
0.194	3.0	0.2	32000	3.37	3.5	$\bullet$
0.194	3.0	0.3	61000	3.73	5.1	$\bullet$
0.194	3.0	0.39	101000	4.15	6.6	$\bullet$
0.194	3.0	0.5	156000	4.55	8.3	$\bullet$

Equations

$$\frac{\pi }{8}{\rho }^f{d}^2{C}_D{u}_{\ast }^2-{\mu }_s\left({\rho }^{p}g\frac{\pi }{6}{d}^3-(f_b{)}_z\right)=0,$$

$${\theta }_c^0=\frac{{\rho }^f{u}_{\ast }^2}{({\rho }^p-{\rho }^f)gd}=\frac{4{\mu }_s}{3C_D}.$$

$$\frac{\pi }{8}{\rho }^f{d}^2{C}_D{u}_{\ast }^2-{\mu }_s\left({\rho }^{p}g\cos \alpha \frac{\pi }{6}{d}^3-(f_b{)}_z\right)+{\rho }^{p}g\sin \alpha \frac{\pi }{6}{d}^3+(f_b{)}_x=0.$$

$${\theta }_c(\alpha )={\theta }_c^0\cos \alpha \left[1-\frac{\tan \alpha }{{\mu }_s}\right].$$

$${\mathbf{f}}_b=-{\rho }^f\frac{\pi }{6}{d}^3\mathbf{g}.$$

$${\mathbf{f}}_b=-\frac{\pi }{6}{d}^3\nabla .{\sigma }_u^f,$$

$${\mathbf{f}}_b=-\frac{\pi }{6}{d}^3{\eta }^f\nabla .[\nabla {\mathbf{u}}^{\mathbf{f}}+(\nabla {\mathbf{u}}^{\mathbf{f}}{)}^T],$$

$$\frac{\pi }{8}{\rho }^f{d}^2{C}_D{u}_{\ast }^2-{\mu }_s({\rho }^p-{\rho }^f)g\cos \alpha \frac{\pi }{6}{d}^3+{\rho }^{p}g\sin \alpha \frac{\pi }{6}{d}^3=0,$$

$${\theta }_c(\alpha )={\theta }_c^0\cos \alpha \left[1-\frac{\tan \alpha }{{\mu }_s}\frac{1}{1-{\rho }^f/{\rho }^p}\right],$$

$$\frac{\pi }{8}{\rho }^f{d}^2{C}_D{u}_{\ast }^2-{\mu }_s({\rho }^p-{\rho }^f)g\cos \alpha \frac{\pi }{6}{d}^3+{\rho }^{p}g\sin \alpha \frac{\pi }{6}{d}^3-\frac{\partial P^f}{\partial x}=0,$$

$$0=\frac{\partial {\tau }_{xz}^f}{\partial z}+\frac{\partial R_{xz}^f}{\partial z}+{\rho }^f(1-\varphi )g\sin \alpha -n{⟨{f_f^p}_x⟩}^p,$$

$$0=\frac{\partial {\tau }_{xz}^p}{\partial z}+{\rho }^p\varphi g\sin \alpha +n{⟨{f_f^p}_x⟩}^p,$$

$$0=-\frac{\partial P^f}{\partial z}+{\rho }^f(1-\varphi )g\cos \alpha -n{⟨{f_f^p}_z⟩}^p,$$

$$0=-\frac{\partial P^p}{\partial z}+{\rho }^p\varphi g\cos \alpha +n{⟨{f_f^p}_z⟩}^p,$$

$${\tau }_{xz}^p(z)={\rho }^{p}g\sin \alpha {\int }_z^{h_p}\varphi (\zeta )d\zeta +{\int }_z^{h_p}n{⟨{f_f^p}_x⟩}^p(\zeta )d\zeta .$$

$$P^p(z)=g\cos \alpha ({\rho }^p-{\rho }^f){\int }_z^{h_p}\varphi (\zeta )d\zeta .$$

$$0={\tau }_{xz}^f(h_p)-{\tau }_{xz}^f(z)+R_{xz}^f(h_p)-R_{xz}^f(z)+{\rho }^{f}g\sin \alpha {\int }_z^{h_p}[1-\varphi (\zeta )]d\zeta -{\int }_z^{h_p}n{⟨{f_f^p}_x⟩}^p(\zeta )d\zeta .$$

$${\int }_z^{h_p}n{⟨{f_f^p}_x⟩}^p(\zeta )d\zeta ={\tau }_b+{\rho }^{f}g\sin \alpha {\int }_z^{h_p}[1-\varphi (\zeta )]d\zeta .$$

$${\bar{\varphi }}_z(h_p-z)={\int }_z^{h_p}\varphi (\zeta )d\zeta ,$$

$$\frac{{\tau }_{xz}^p(z)}{{\tau }_{zz}^p(z)}=\frac{{\rho }^p}{{\rho }^p-{\rho }^f}\tan \alpha +\frac{{\tau }_b}{g\cos \alpha ({\rho }^p-{\rho }^f){\bar{\varphi }}_z(h_p-z)}+\frac{{\rho }^f}{{\rho }^p-{\rho }^f}\frac{1-{\bar{\varphi }}_z}{\overline{{\varphi }_z}}\tan \alpha ,$$

$$\frac{{\tau }_b}{g\cos \alpha ({\rho }^p-{\rho }^f){\bar{\varphi }}_z(h_p-z)}=\frac{{\theta }^{\ast }}{{\bar{\varphi }}_z}\frac{d}{h_p-z}.$$

$$\begin{gathered} \frac{{\tau }_{xz}^p(z)}{P^p(z)}=\frac{{\rho }^p}{{\rho }^p-{\rho }^f}\tan \alpha +\frac{{\rho }^f}{{\rho }^p-{\rho }^f}\frac{1-{\varphi }^{max}}{{\varphi }^{max}}\tan \alpha \\ =\frac{\tan \alpha }{{\rho }^p/{\rho }^f-1}\left[\frac{{\rho }^p}{{\rho }^f}+\frac{1-{\varphi }^{max}}{{\varphi }^{max}}\right]. \end{gathered}$$

$$\tan {\alpha }_0=\frac{{\mu }_s}{1+{\left[({\rho }^p/{\rho }^f-1){\varphi }^{max}\right]}^{-1}}.$$

$$Q_s={\int }_0^{\infty }\varphi (z){⟨v_x⟩}^p(z)dz={\int }_{h_c}^{h_p}\varphi (z){⟨v_x⟩}^p(z)dz=\overline{\varphi {⟨v_x⟩}^p}{\delta }_s.$$

$$\frac{{\delta }_s}{d}=\frac{{\theta }^{\ast }}{\bar{\varphi }\left[{\mu }_s-\left(1+[\bar{\varphi }({\rho }^p/{\rho }^f-1){]}^{-1}\right)\tan \alpha \right]}\simeq \frac{{\theta }^{\ast }}{{\varphi }^{max}{\mu }_s\left[1-\tan \alpha /\tan {\alpha }_0\right]}.$$

$$\overline{\varphi {⟨v_x⟩}^p}\sim \sqrt{({\rho }^p/{\rho }^f-1)g{\delta }_s}.$$

$$Q_s^{\ast }=\frac{Q_s}{\sqrt{({\rho }^p/{\rho }^f-1)g{d}^3}}\sim \frac{\sqrt{({\rho }^p/{\rho }^f-1)g{\delta }_s^3}}{\sqrt{({\rho }^p/{\rho }^f-1)g{d}^3}}\sim {\left(\frac{{\theta }^{\ast }}{{\varphi }^{max}{\mu }_s\left[1-\tan \alpha /\tan {\alpha }_0\right]}\right)}^{3/2}.$$

$$m\frac{d^2{\mathbf{x}}^p}{d{t}^2}={\mathbf{f}}_{ext}^p+\sum _{k\in \mathcal{N}}{\mathbf{f}}_c^{pk}={\mathbf{f}}_g^p+{\mathbf{f}}_f^p+\sum _{k\in \mathcal{N}}{\mathbf{f}}_c^{pk}$$

$${\mathbf{f}}_b^p=-\frac{\pi d^3}{6}\nabla P^f,$$

$${\mathbf{f}}_D^p=\frac{1}{2}{\rho }^f\frac{\pi d^2}{4}{C}_D\left|\left|{⟨\mathbf{u}⟩}_{{\mathbf{x}}^{\mathbf{p}}}^f-{\mathbf{v}}^{\mathbf{p}}\right|\right|\left({⟨\mathbf{u}⟩}_{{\mathbf{x}}^{\mathbf{p}}}^f-{\mathbf{v}}^{\mathbf{p}}\right),$$

Full text

1 Introduction

Turbulent bedload transport is of major importance for the prediction of riverbed evolution and coastal processes, which represent important issues for public safety, management of water resources and environmental sustainability. In this framework, the key parameter to predict is the dimensionless sediment transport rate (Einstein, 1942), $Q_{s}^{*}=Q_{s}/\sqrt{(\rho^{p}/\rho^{f}-1)gd^{3}}$, as a function of the dimensionless fluid bed shear stress denoted as the Shields number (Shields, 1936), $\theta^{*}=\tau_{b}/[(\rho^{p}-\rho^{f})gd]$, where $Q_{s}$ is the volumetric sediment transport rate per unit width, $g$ is the acceleration of gravity, $d$ the particle’s diameter, $\rho^{p}$ and $\rho^{f}$ are the particle and fluid densities and $\tau_{b}$ is the fluid bed shear stress. Due to the inherent complexity of granular media behaviour and turbulent fluid flows, turbulent bedload transport understanding remains limited despite a century of modern research on the subject (Gilbert, 1914; Bagnold, 1956; Frey & Church, 2011; Duran et al., 2012; Aussillous et al., 2013). This is illustrated by the poor predictions provided by the classical formulas linking Einstein and Shields numbers - such as the Meyer-Peter & Müller (1948) formula - which lead to sediment transport rate up to two orders of magnitude different from what is observed in the field (Recking et al., 2013). Accordingly, the present paper focuses on the analysis of the slope influence in turbulent bedload transport, which might be one of the key aspect of the observed data dispersion.

Most applications of turbulent bedload transport involve the presence of a slope, for example in the case of a beach in coastal sediment transport, a river or a mountain stream. The slope inclination angle is expected to affect the sediment transport rate through a modification of the particle’s mobility. This is classically accounted for by considering a force balance on a single grain at the top of the granular bed close to the onset of motion (Fredsøe & Deigaard, 1992; Andreotti et al., 2013). In the zero-slope case (see figure 1), the streamwise force balance at the onset of motion reduces to an equality between the streamwise force induced by the fluid flow and the resistive sliding friction force on the granular layer below. Considering only the main fluid forces to apply to the grain, i.e. drag and buoyancy (Schmeeckle et al., 2007) the friction force can be expressed as a granular friction coefficient, $\mu_{s}$, multiplied by the buoyant weight representing the vertical force applied to the grain. Then, the force balance on a grain at the onset reads:

[TABLE]

where $C_{D}$ is the drag coefficient, $\mathbf{f}_{b}$ is the buoyancy force and $u_{*}$ is the velocity scale at the granular bed assimilated to the fluid bed friction velocity. From this balance, by expressing the buoyancy force, the critical Shields number can be written as:

[TABLE]

In the presence of a longitudinal slope inclination angle, $\alpha$, (see figure 1) two additional positive terms appear in the force balance due to the projection of the particle weight and buoyancy force along the streamwise axis, while the friction force is reduced due to the projection of the particle weight along the vertical axis:

[TABLE]

Taking the buoyancy force as $\mathbf{f}_{b}=-\rho^{f}\frac{\pi}{6} d^{3}\mathbf{g}$, leads to the following reduction of the critical Shields number with increasing slope:

[TABLE]

This expression of the modification of the critical Shields number has been first formulated by Fernandez Luque & Van Beek (1976) considering the onset of motion to be due to rolling instead of sliding. The consecutive moment balance on the grain leads to the same expression of the modified critical Shields number, with $\mu_{s}$ the tangent of the so-called pocket angle formed by the local arrangement between the particle and its neighbours (Wiberg & Smith, 1987). Performing turbulent bedload transport experiments in an inclined rectangular pressure-driven closed conduit, with variation of particle diameter ($d\in[0.9,3.3]mm$), density ratio $\rho^{p}/\rho^{f}\in[1.34,4.5]$ and channel inclination angle ($\alpha\in[0,22]^{\circ}$), Fernandez Luque & Van Beek (1976) found a relatively good agreement between the theoretical prediction and experimental results providing a fit of the pocket angle, which appeared unexpectedly large. Following this pioneering work, Chiew & Parker (1994) reproduced a similar approach, deriving equation (4) from a force balance associating $\mu_{s}$ with the granular medium repose angle and comparing the prediction to experimental data in pipe flows. Varying the particle mean diameter $d\in[0.5,3]mm$, the repose angle $\Phi\in[33,38]^{\circ}$ and the channel inclination angle $\alpha\in[0,31]^{\circ}$, they showed that the critical Shields number follows the prediction of equation (4) for downward slopes. This work has been further generalised by Seminara et al. (2002) and Dey (2003) to the combined effect of transverse and longitudinal slopes, validating the analysis with experiments in pipe flows in the latter study.

This type of approach is widespread in turbulent sediment transport and can be found in classical textbooks of sediment transport (Fredsøe & Deigaard, 1992) and granular media (Andreotti et al., 2013), as well as in the aeolian saltation community (Iversen & Rasmussen, 1994, 1999). It is also known in the morphodynamic community as the Ikeda-Coleman-Iwagaki model (Wiberg & Smith, 1987) and has been applied in turbulent bedload transport to experiments and field transport rate predictions (see e.g. Li & Cheng (1999); Wilcock & Crowe (2003); Karmaker & Dutta (2016)). In addition, it has been extended to account for the lift force (Wiberg & Smith, 1987; Chiew & Parker, 1994; Armanini & Gregoretti, 2005), the fluid viscous sub-layer at low slope (Wiberg & Smith, 1987) and the duration of a given applied fluid force in the context of turbulent fluid flows (Diplas et al., 2008; Valyrakis et al., 2010).

In the meantime, the longitudinal slope effect has been studied experimentally in gravity-driven turbulent bedload transport. Considering a large range of slopes, Smart & Jaeggi (1983); Smart (1984) and later Rickenmann (1991, 2001) have determined empirical relationships between the dimensionless sediment transport rate and a modified Shields number. The latter, classically evaluated from the water depth, $h$, was determined from the so-called mixture depth, $h_{m}=h+\delta_{s}$, including both the water depth, $h$ and the thickness of the granular layer in motion, $\delta_{s}$. As a consequence, the Shields number has a different meaning and this approach cannot be compared directly to the classical one, as the mobile layer thickness $\delta_{s}$ scales with the classical Shields number and the slope (Sumer et al., 1996; Hsu et al., 2004; Revil-Baudard & Chauchat, 2013). These dependencies are expected to modify the classical scaling laws, underlining the importance of adopting the same definition of the Shields number when comparing experimental data and theoretical predictions. Although in a different framework, these studies show that the slope not only modifies the critical Shields number but also the Shields number definition. As a consequence, the empirical law relating the dimensionless sediment transport rate to the Shields number is also modified. A similar modification of the sediment transport formula has been proposed by Cheng & Chen (2014), replacing the gravity contribution by a slope-modified gravity in both the Shields and Einstein number formulations and considering the classical Meyer-Peter & Müller (1948) formula. The obtained formula fits the existing experimental data better than the classical corrections (Cheng & Chen, 2014), but lacks a solid theoretical justification. Lastly, it is interesting to note that Damgaard et al. (1997) also proposed an empirical correction of the transport formula based on experimental data in a closed conduit.

The literature review underlines two different trends that seem to be associated with pressure-driven and gravity-driven configurations in turbulent bedload transport. On one hand, gravity-driven experiments exhibit a modification of the transport formula as a function of the slope. On the other hand, the variation of the critical Shields number with the slope seems to be well predicted from a force/moment balance on a single particle in pressure-driven configurations. However, there is a priori no theoretical justification why a variation of the slope should only affect the critical Shields number and the study of Damgaard et al. (1997) at moderate Shields number values suggests a behaviour similar to the gravity-driven configuration. This lack of characterisation together with the absence of clear theoretical bases in the literature suggests the need for further analysis.

In the present contribution, we attempt to give a better understanding of the longitudinal slope influence on turbulent bedload transport by adopting an idealised and theoretical point of view. Focusing on gravity-driven turbulent bedload transport under steady uniform flow conditions, we discuss the bases of the critical Shields number derivation and analyse the granular entrainment mechanisms in the framework of the two-phase continuous equations (section 2). This allows us to characterise the influence of the slope on the vertical flow structure and propose a re-scaling of the Shields number to account for the slope influence on the sediment transport rate. The proposed scaling is tested against fluid-discrete element method simulations (section 3) and the theoretical results are discussed more generally in the light of the numerical results, considering in particular the vertical flow structure and the difference between gravity-driven and pressure-driven configurations (section 4).

2 Theoretical analysis

2.1 Discussion on the classical critical Shields number derivation

The classical derivation of the critical Shields number reproduced in the introduction relies on the following expression of the buoyancy force applied to a particle:

[TABLE]

While this expression is classically used, it does not apply a priori to all the different configurations explored in fluid-particle flow, as stressed by Christensen (1995).

The buoyancy force is defined as the force a fluid element would undergo if it was occupying the position of the particle (Maxey & Riley, 1983). It can be derived explicitly in the Stokes flow case and leads to the following formulation (Maxey & Riley, 1983):

[TABLE]

where $\sigma^{f}_{u}$ is the undisturbed fluid velocity stress tensor, i.e. the fluid stress tensor based on the fluid velocity field undisturbed by the presence of the particles. This formulation is used in Jackson (2000) at the level of the continuous medium, including the Reynolds stresses contribution ($R^{f}$) inside the fluid stress tensor, i.e. replacing $\sigma_{u}^{f}$ by $\sigma_{u}^{f}+R^{f}$ inside equation (6). However, the formulation of the buoyancy force should not a priori include the contribution from the Reynolds stresses considering that they are associated with advection and related to the time averaging of the fluid velocity fluctuations. Indeed, considering a single particle in a turbulent Newtonian fluid flow, it experiences a time-averaged buoyancy force depending on the average fluid velocity field but independent of the fluctuations with time. This is the case because the expression of the buoyancy is linear with the second spatial derivative of the fluid velocity for a Newtonian fluid,

[TABLE]

with $\eta^{f}$ the dynamic fluid viscosity, so that the time averaging of the buoyancy force reduces to the contribution from the average fluid velocity field. From these considerations, the expression of the buoyancy force should follow the equation (6) in all the different fluid flow regimes (Lhuillier, 2011).

As a consequence, the above-mentioned derivation is not valid in uniform gravity-driven turbulent bedload transport as the undisturbed fluid equation, $0=-\nabla R^{f}-\nabla\sigma_{u}^{f}+\rho^{f}\mathbf{g}$, together with the regular expression of the buoyancy force (eq. 6), does not allow us to write the buoyancy force as equation (5).

Expressing the buoyancy directly from equation (6), it follows that $(f_{b})_{z}=\pi d^{3}/6\rho^{f}g\cos\alpha$ and $(f_{b})_{x}=-\pi d^{3}/6~{}\partial_{z}\sigma^{f}_{xz}$, by using the hydrostatic characteristics of the fluid along $z$. Considering that the viscous shear stress contribution is negligible in turbulent bedload transport $\partial_{z}\sigma^{f}_{xz}\simeq 0$ (Revil-Baudard & Chauchat, 2013; Maurin et al., 2016), the force balance performed on a particle lying on the granular bed at the onset of motion (eq. 3) reads:

[TABLE]

where only the particle weight projection along the streamwise axis is non-negligible. This leads to the following modification of the critical Shields number due to the slope effect:

[TABLE]

where the density ratio $\rho^{p}/\rho^{f}$ modifies the classical slope effect characterised by equation (4). This difference is of importance for underwater bedload transport, where the density ratio is typically of order one ($\rho^{p}\gtrsim\rho^{f}$).

Adopting the regular expression for the buoyancy force also affects the derivation of the slope influence on the critical Shields number in the pressure-driven configuration. In this case, the streamwise component of the buoyancy force is completed by the term due to the presence of a longitudinal pressure gradient:

[TABLE]

which does not enable us to recover the classical modification of the critical Shields number derived in the introduction (eq. 4). The latter equation remains only valid for laminar flows and in the case where the fluid mass column stays at a constant level perpendicular to the gravity, i.e. for underwater avalanches or coastal sediment transport.

2.2 Vertical granular flow structure and entrainment mechanisms

Instead of focusing on the critical Shields number modification in a discrete particle framework, we adopt a more general approach based on the continuous two-phase flow framework to analyse the slope influence from the onset of motion up to intense turbulent bedload transport. In order to better understand the vertical structure of the granular flow and the local mechanisms at play, let us express the shear-to-normal granular stress ratio as a function of the depth $\tau^{p}_{xz}(z)/P^{p}(z)$. The yield criterion $\tau^{p}_{xz}/P^{p}>\mu_{s}$ is characteristic of the onset of granular flow, so that the positive contributions to the shear-to-normal stress ratio at a given height $z$ reflects the local entrainment mechanism of the granular medium.

The two-phase volume-averaged equations read for steady uniform flows (Anderson & Jackson, 1967; Jackson, 2000; Chauchat, 2017):

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Where $\sigma_{ij}^{f}=-P^{f}\delta_{ij}+\tau_{ij}^{f}$ is the effective fluid stress tensor with $\delta_{ij}$ the Kronecker’s delta, $\sigma_{ij}^{p}=-P^{p}\delta_{ij}+\tau_{ij}^{p}$ is the granular stress tensor, $R^{f}_{ij}$ is the Reynolds stress tensor, $\phi$ and $\epsilon$ are the solid and fluid volume fraction respectively and $n\left<{f_{f}^{p}}_{k}\right>^{p}$ is the $k$ component of the fluid-solid momentum transfer, with $n$ the particle density and $\left<{f_{f}^{p}}_{k}\right>^{p}$ the volume-averaged fluid-particle interaction force. These equations represents the volume-averaged momentum balance over the streamwise and wall-normal directions of the fluid and solid phases and the brackets with superscript $p$ denotes the classical volume averaging over the particle phase defined in Anderson & Jackson (1967).

Integrating these equations between a position $z<h_{p}$ in the moving granular layer and $h_{p}$, the top of it, allows us to simplify the equations and to express the shear-to-normal stress ratio as a function of the vertical position $z$. Focusing on the granular equations (12) and (14), the granular stress tensor is zero at $h_{p}$ and the integrated streamwise component of the granular momentum balance leads to:

[TABLE]

Assuming that the fluid wall-normal pressure is hydrostatic, the fluid-particle interaction term $n\left<{f_{f}^{p}}_{z}\right>^{p}$ can be expressed from equation 13. Replacing it in equation (14), we obtain:

[TABLE]

Lastly, the fluid streamwise momentum balance (eq. 11) can be integrated between $z$ and $h_{p}$ in order to express the last term on the right-hand side of equation (15):

[TABLE]

It has been shown previously in turbulent bedload transport and sheet-flow simulations that the effective viscous shear stress tensor is negligible throughout the depth with respect to the Reynolds stresses and to the slope contribution (Revil-Baudard & Chauchat, 2013; Maurin et al., 2016). Furthermore, the Reynolds stresses are fully damped in the granular bed and we assume that their contribution is negligible in the moving granular layer $R_{xz}^{f}(z)\sim 0$. Also, the fluid bed shear stress used to define the Shields number is taken here as the maximum of the turbulent shear stress (Revil-Baudard & Chauchat, 2013; Maurin et al., 2016), assuming that the latter is located at the top of the granular layer, $\tau_{b}=R_{xz}^{f}(h_{p})$. Therefore, equation (17) reads:

[TABLE]

Combining equations (15), (17) and (18) and defining the average solid volume fraction between $z$ and $h_{p}$, $\bar{\phi}_{z}$ as:

[TABLE]

the shear-to-normal stress ratio $\tau_{xz}^{p}(z)/\tau_{zz}^{p}(z)$ can be expressed as:

[TABLE]

Three contributions to the granular phase shear-to-normal stress ratio can be identified. The first term on the right-hand side of equation (20) is constant within the depth and represents the slope effect on the granular phase. It is analogous to the classical configuration of a dry granular medium on an inclined plane, where the shear-to-normal stress ratio is constant within the depth and equal to the tangent of the inclination angle $\alpha$ (Andreotti et al., 2013). In bedload transport, the presence of a buoyancy force along the vertical axis leads to a modification of this term by the density ratio. The two last terms on the right-hand side correspond to the fluid contribution on the granular phase, which has been split into two contributions coming from the fluid flow above the bed and the fluid flow inside the granular layer respectively. The Shields number can be identified in the former and the term reduces to:

[TABLE]

This represents the contribution from the fluid bed shear stress to the entrainment of particles. Its surface nature is characterised by a decrease of this contribution inside the granular bed, i.e. as $z$ decreases. It is accounted for in bedload transport studies, by assuming that the Shields number is representative of the total imposed fluid shear stress. However, the contribution from the fluid flow inside the granular layer, corresponding to the third term on the right-hand side of equation (20), is usually neglected while it is seen to be of importance when increasing the slope angle. This contribution depends on the slope and characterises the effect of the slope on the fluid flow inside the granular layer. The latter affects the granular phase through the fluid-particle interaction and in consequence depends also on the specific density. In addition, the different dependence on the vertical position $z$ of the two fluid terms implies that a modification of the slope at given Shields number and specific density would induce a change in the vertical repartition of the shear-to-normal stress ratio and therefore on the vertical structure of the granular flow. This could explain why the velocity and solid volume fraction profiles are not self-similar in gravity-driven bedload transport (Armanini et al., 2005; Larcher et al., 2007; Capart & Fraccarollo, 2011; Frey, 2014; Maurin, 2015; Maurin et al., 2015).

In order to gain more insight into the physical meaning of the fluid flow contribution inside the granular bed, let us consider the simple case of a gravity-driven fluid flow through a quasi-static porous granular bed. The Shields number contribution is negligible in this configuration and the granular bed is at rest so that the average solid volume fraction is equal to the maximum packing fraction $\bar{\phi}_{z}=\phi^{max}$. This leads to a shear-to-normal stress ratio independent of the vertical position $z$:

[TABLE]

As a consequence, there exists a critical angle above which the entire granular layer will be entrained independently from its thickness. This avalanche angle, $\alpha_{0}$, corresponds to the configuration for which the shear-to-normal stress ratio exceeds the static friction coefficient of the granular medium $\mu_{s}$:

[TABLE]

This corresponds to the critical angle for the onset of debris flow as derived by Takahashi (1978) and is well known in the debris flow community (Takahashi, 2007). The derivation presented here makes the link between debris flows and gravity-driven turbulent bedload transport and evidences the importance of the slope influence through the fluid flow inside the granular layer when considering steep slope configurations. In addition, it shows that the expression of the critical angle of debris flow $\alpha_{0}$ is characteristic of the slope influence on the granular phase through both its direct impact and the fluid-induced contribution in gravity-driven bedload transport, two mechanisms that depend on the specific density $\rho^{p}/\rho^{f}-1$.

2.3 Sediment transport rate scaling

Taking advantage of the developed analysis, we look for a scaling law of the sediment transport rate accounting for the slope effect. Following Bagnold (1956), the sediment transport rate per unit width $Q_{s}$ can be expressed as a function of the solid volume fraction $\phi$ and the average particle velocity $\left<v_{x}\right>^{p}$ profiles:

[TABLE]

where $h_{c}$ is the lower bound of the mobile granular layer, $\delta_{s}=h_{p}-h_{c}$ is the mobile granular layer thickness and $\bar{\bullet}$ denotes the average over $z$ along the mobile granular layer thickness. Considering $z=h_{c}$ in equation (20), the mobile layer thickness can be expressed as:

[TABLE]

where equation (23) has been used and it has been assumed that $\bar{\phi}\simeq\phi^{max}$ as the average solid volume fraction $\bar{\phi}$ is only weakly varying in the problem.

The term $\overline{\phi\left<v_{x}\right>^{p}}$ has the dimension of a velocity, i.e. the square root of the product of an acceleration and a length scale. In the problem, the only acceleration scale is the buoyant reduced gravity $\tilde{g}=(\rho^{p}/\rho^{f}-1)g$ (Duran et al., 2012) and we have seen that the relevant length scale is the mobile granular layer thickness. Thus, the volume flux should scale as:

[TABLE]

Compiling equations (24), (25) and (26), the scaling law of the dimensionless sediment transport rate $Q_{s}^{*}$ can be written as:

[TABLE]

This approach not only allows us to recover the scaling law of the dimensionless sediment transport rate with the Shields number and the power $3/2$ in the high Shields number limit, but also allows us to express the influence of the slope through an effective granular friction coefficient, $\mu_{eff}=\mu_{s}\left[1-\tan\alpha/\tan\alpha_{0}\right]$. The proposed scaling law is consistent with the original definition of the Shields number as a ratio between the driving and resistive forces applied to the granular phase. Indeed, it represents the ratio between the tangential fluid bed shear stress $\tau_{b}$ and the tangential shear stress associated with the friction between granular layers $\tau=\mu_{eff}P^{p}=\phi\mu_{eff}(\rho^{p}-\rho^{f})gd$. The dependence of the effective friction coefficient on the distance to the critical debris flow angle, $1-\tan\alpha/\tan\alpha_{0}$, is consistent with the analysis of section 2.2 and allows us to take into account both the direct and indirect slope effects on the granular phase in gravity-driven bedload transport. It suggests that the key parameter is not the slope in itself but the ratio between the slope and the critical angle of debris flow. This is of particular importance as the latter is not constant and depends on the maximum solid volume fraction and on the specific density.

3 Numerical Analysis

In order to verify the developed analysis, numerical simulations of turbulent bedload transport are performed using an existing coupled fluid-discrete element model (Maurin et al., 2015).

3.1 Numerical model

The model is based on a three-dimensional (3-D) discrete element method (DEM) coupled with a 1-D volume-averaged fluid momentum balance. The fine resolution of the granular phase and the two-way coupling ensures an accurate description of the vertical depth profiles and the momentum conservation of the system on average. The model has been described in detail and validated experimentally in Maurin et al. (2015), so that only the main characteristics are recalled here. The typical configuration considered here is shown in figure 2.

The bi-periodic 3-D DEM is based on the explicit resolution of the particle phase, solving Newton’s equations of motion for each particle $p$ at position $x_{p}$:

[TABLE]

where the sum of the contact forces $\mathbf{f}_{c}^{pk}$ is made over the ensemble of nearest neighbours ${\cal N}$, $\mathbf{f}_{g}^{p}$ is the gravity force, $\mathbf{f}_{f}^{p}$ is the force applied by the fluid on particle $p$. Similarly, the rotation of the particles are solved from Newton’s equations of motion. The contact force between particles is determined from the particles overlap using the classical spring-dashpot contact law (Schwager & Pöschel, 2007) which defines a unique normal restitution coefficient $e_{n}$ and considers tangential friction characterised by a friction coefficient $\mu_{p}$. The latter two are taken as respectively $e_{n}=0.5$ and $\mu_{p}=0.4$, as determined from experimental comparisons (Maurin et al., 2015). The interaction with the fluid phase ($\mathbf{f}_{f}^{p}$) is mainly restricted to buoyancy, drag and lift forces in turbulent bedload transport (Nino & Garcia, 1994; Niño & García, 1998). While it is clear that the lift force plays a non-negligible role in turbulent bedload transport (Ji et al., 2013), its expression has only been derived in the limits of Stokes flow (Saffman, 1965) or an inviscid fluid (Schmeeckle et al., 2007) and does not apply to the fluid flow regimes associated with turbulent bedload transport (Schmeeckle et al., 2007). As a consequence, it has been decided to avoid including a controversial expression of the lift force, which appeared to be unnecessary to reproduce accurately turbulent bedload transport experiments (Maurin et al., 2015). Therefore, the fluid-particle interactions are restricted in the present work to the three-dimensional buoyancy ($\mathbf{f}_{b}^{p}$) and drag forces ($\mathbf{f}_{D}^{p}$) (Jackson, 2000; Maurin et al., 2015):

[TABLE]

[TABLE]

where the average fluid velocity and the fluid pressure are taken at the centre of particle $p$ and $\mathbf{v^{p}}$ represents the velocity of particle $p$. The drag coefficient $C_{D}$ depends on the particle Reynolds number $\Rey_{p}=\left|\left|\left<\mathbf{u}\right>^{f}_{\mathbf{x^{p}}}-\mathbf{v^{p}}\right|\right|d/\nu^{f}$ and takes into account hindrance effects (DallaValle, 1948; Richardson & Zaki, 1954): $C_{D}=\left(0.4+24.4/Re_{p}\right)(1-\phi)^{-3.1}$. Knowing the position and velocity of the particles, the fluid pressure and velocity fields, Newton’s equations of motion are solved for the ensemble of particles using the open-source DEM code YADE (Šmilauer et al., 2015).

The coupled fluid description solves the volume-averaged momentum balance of the fluid phase, corresponding to equations (11) and (13):

[TABLE]

[TABLE]

Due to the steady uniform character of the problem, the fluid velocity field reduces to its streamwise component and depends only on the vertical position: $\left<\mathbf{u}\right>^{f}=\left<u_{x}\right>^{f}(z)\ \mathbf{e_{x}}$, as sketched in figure 2. Therefore, the buoyancy force (eq. 29) is restricted to its wall-normal component and equation (32) leads to a hydrostatic fluid pressure distribution, the wall-normal average drag force being negligible. The solution of the streamwise fluid phase momentum balance (eq. 31) requires closures for the viscous and Reynolds stress tensors, as well as the determination of the solid volume fraction $\phi$ and the momentum transfer associated with the granular phase interactions. Considering a Newtonian fluid, the effective viscous shear stress is expressed as:

[TABLE]

with $\nu^{f}$ the clear fluid kinematic viscosity and $\left<u_{x}\right>^{f}$ the volume-averaged streamwise fluid phase velocity. The Reynolds shear stress is based on the eddy viscosity concept ($\nu^{t}$) using a mixing length formulation:

[TABLE]

where the mixing length is taken similarly to Li & Sawamoto (1995) as:

[TABLE]

with $\kappa=0.41$, the von Kármán constant. The formulation adopted allows us to recover the law of the wall (Prandtl, 1926) in a clear fluid, while completely damping the turbulence inside the granular bed at maximum packing fraction ($\phi^{max}$).

The solid volume fraction ($\phi$) and the fluid-particle interaction term ($n\left<{f_{f}^{p}}\right>^{p}$) are determined from spatial averaging of the discrete solid phase. Considering cubic boxes of finite wall-normal length scale $l_{z}$, these two terms are averaged over the whole length and width of the granular bi-periodic cell. In order to solve the important wall-normal gradients present in turbulent bedload transport, a small wall-normal weighting function length scale has been adopted (typically $l_{z}=d/30$) and this choice has been confirmed by the experimental validation (Maurin et al., 2015).

The 3-D DEM and the fluid model are solved as transient problems applying a fixed bottom boundary condition for both the fluid ($\left<\mathbf{u}\right>^{f}(z=0)=0$) and the particle phase (fixed random particles) and imposing the position of the water free surface ($d\left<u_{x}\right>^{f}/dz(z=h)=0$). In order to achieve a stable integration, the DEM time step is bounded by the propagation time of the fastest wave over a particle diameter (Maurin, 2015; Maurin et al., 2015). The fluid resolution time step corresponds to a typical characteristic evolution time scale of the granular medium and is taken much larger than the DEM one (Maurin et al., 2015): $\Delta t_{f}=10^{-2}s$ with respect to $\Delta t_{p}\sim O(10^{-4}-10^{-5})\ s$.

The model has been compared with experiments and has shown its capability to describe accurately both the sediment transport rate and the granular depth structure in turbulent bedload transport (Maurin et al., 2015).

3.2 Results

In order to study the bulk flow behaviour and to investigate the parameter space, bi-periodic three-dimensional numerical simulations are performed by varying the Shields number $\theta^{*}$ between the onset of motion $\theta^{*}_{c}$ and intense bedload transport ($\theta^{*}\sim 0.6$), the relative slope angle $\tan\alpha/\tan\alpha_{0}$ between $0.1$ and $0.75$ and the specific density $\rho^{p}/\rho^{f}-1$ between $0.75$ (lightweight plastic), $1.5$ (glass/natural sediment) and $3$ (metal). Detailed parameters of the runs are shown in table 1. In the simulations, the water free surface position and the channel slope angle are imposed before letting the system evolves under the effect of gravity. After reaching the steady state, the data are averaged over space every typical granular evolution time scale, $0.1s$, over 1000 measurements (Maurin et al., 2015). Note that the eddy turnover time varies between $0.01s$ and $0.1s$ in the configurations considered. Consistently with the theoretical derivation of section 2.2, the Shields number is evaluated by taking the fluid bed shear stress as the maximum of the time-averaged turbulent Reynolds shear stresses. The sediment transport rate is evaluated from equation (24), neglecting therefore the local instantaneous lateral and vertical contributions. The latter have been checked to be negligible with respect to the streamwise contribution.

Figure 3(a,b) shows the dimensionless sediment transport rate as a function of the Shields number for variation of the channel inclination angle and the specific density respectively. While being far from the threshold of motion, an increase of up to an order of magnitude is observed in the dimensionless sediment transport rate when the slope increases at constant Shields number. This important influence at high Shields number evidences the impact of the bed slope on the sediment transport rate. In addition, the variation of the specific density at constant slope and Shields number is also seen to affect the dimensionless sediment transport rate (see fig. 3b), as expected from the theoretical analysis presented in section 2.

Combining the variation of the slope and of the specific density, figure 4a shows a complex picture of the phenomenon, with an important dispersion at given Shields number values resulting from the coupling between the two mechanisms.

Evaluating the granular medium friction coefficient from dry inclined plane avalanche simulations ($\mu_{s}\simeq 0.4$), the results are plotted in figure 4b as a function of the modified Shields number proposed in the previous section. All the data are shown to collapse on a single curve of power law close to 3/2, with some dispersion at high modified Shields number values. This collapse shows that the rescaling characterises at first order the effect of the slope and of the specific density variations on the sediment transport rate in gravity-driven turbulent bedload transport. While the small dispersion observed at high modified Shields number values highlights the limits of the proposed scaling law, its intuitive nature, its simplicity and the fact that it encompasses the different physical mechanisms represent a clear improvement compared with previous corrections.

Additional simulations have been performed close to the onset of motion to characterise the critical Shields number. The latter has been evaluated by forcing the flow slightly above the onset of motion before letting the simulations evolve at the given Shields number values in order to characterise the cessation threshold (Ouriemi et al., 2007). The results are shown in figure 5 and suggest that the proposed rescaling is valid at first order close to the critical Shields number. Even though the data exhibit more dispersion due to the complex characterisation of the onset of motion (Clark et al., 2015), the rather good collapse observed suggests that the proposed formulation allows us to define a unique critical Shields number independent of the slope and the specific density, which is close to $\theta^{*}_{c}/[\mu_{s}(1-\tan\alpha/\tan\alpha_{0})]\simeq 0.1$ in the present simulations.

In order to illustrate the underlying mechanisms of the slope influence on the vertical granular flow structure, the solid depth profiles of volume fraction, $\phi$, average streamwise particle velocity, $<v_{x}^{p}>$ and transport rate density, $q_{s}=\phi<v_{x}^{p}>$ are shown in figure 6 for various Shields number values at constant $\tan\alpha/\tan\alpha_{0}=0.1$.

The solid volume fraction depth profiles are characterised by the presence of a possible fixed point at the interface between the granular bed and the fluid flow, around $\phi=0.3$. The Shields number increase leads to an increase of the solid volume fraction above the interface and a decrease below it, associated with an increase of the mobile granular layer thickness. The corresponding solid velocity profiles keep the same shape and exhibit a shift to higher values with increasing Shields number values, similarly to the transport rate density profiles which are broadening in the same time due to the increased mobile layer thickness.

The presence of the fixed point in the solid volume fraction profiles at a low relative slope inclination angle suggests that an appropriate non-dimensionalisation of the $z$ component could collapse all the solid volume fraction profiles obtained at different Shields number values. However, as we have seen in section 2.2, the vertical structure of the granular flow results from a competition between the slope influence on the lower part of the granular flow and the surface contribution from the dimensionless fluid bed shear stress, i.e. the Shields number. These two effects are mostly independent and they are dominant in the lower part and the upper part of the granular flow respectively. This is well observed when considering the variation of the slope influence ($\tan\alpha/\tan\alpha_{0}$) at a constant Shields number (see figure 7).

Increasing the slope induces an increasing vertical asymmetry in the transport rate profiles towards the bottom, resulting from a shift of the concentration shoulder to higher solid volume fraction and an increase of the solid velocity in the lower part of the mobile layer. The absence of a fixed point and this changing asymmetry in the solid volume fraction and transport rate density profiles show that the granular depth profiles cannot be self-similar in gravity-driven turbulent bedload transport.

4 Discussion

The present theoretical analysis and numerical simulations of gravity-driven turbulent bedload transport demonstrate that the effect of the slope on the sediment transport rate is not restricted to the critical Shields number, but influence the whole transport rate formula, in agreement with part of the literature (Smart & Jaeggi, 1983; Smart, 1984; Rickenmann, 1991; Cheng & Chen, 2014). It highlights the importance of the critical angle of debris flow in the transport rate prediction. The debris flow angle is representative of the effective resistance of the granular layer at rest and includes both the direct and indirect slope influence on the granular medium. Therefore, the analysis suggests that the key parameter in turbulent bedload transport is not the slope but the ratio between the tangent of the slope angle and the tangent of the critical angle of debris flow, $\tan\alpha/\tan\alpha_{0}$.

As a consequence, the specific density or density ratio plays a crucial role in the slope effect on the sediment transport rate. A variation of the density ratio has important consequences due to the associated modification of the critical debris flow angle and care should be taken when using model material in experiments. In particular, plastic particles are commonly used for sheet flows and bedload transport experiments (Armanini et al., 2005; Larcher et al., 2007; Capart & Fraccarollo, 2011; Ni & Capart, 2015; Revil-Baudard et al., 2015), with density of approximately $\rho^{p}/\rho^{f}=1.2$, leading to a critical debris flow angle divided by four with respect to the classical density ratio for underwater natural sediment, $\rho^{p}/\rho^{f}=2.65$. This modification strongly influences both the sediment transport rate and the vertical structure of the granular flow, preventing a direct application of the results to natural sediment transport. In particular, the absence of concentration shoulder in the solid volume fraction profile measured by Capart & Fraccarollo (2011) and Sumer et al. (1996) could be explained by the use of plastic particles in their experiments. Indeed, the resulting higher importance of the fluid contribution inside the granular layer tends to smooth out the concentration shoulder as seen in figure 7.

Lastly, we discuss the common points and differences between the gravity-driven and the pressure-driven configurations in steady uniform turbulent bedload transport. To illustrate the comparison, we consider the entrainment mechanisms of the granular medium in the pressure-driven configuration, similarly to section 2.2 for the gravity-driven case. Performing the same integration of the two-phase equation, the shear-to-normal stress ratio in the pressure-driven configuration reads (see appendix A):

[TABLE]

The two first terms on the right-hand side correspond to the direct influence of the slope on the granular phase and the Shields number dependence. They are exactly the same as for the gravity-driven case (see equation 20). The last term on the right-hand side represents the effect of the fluid flow induced by the pressure gradient inside the granular layer and does not vanish at zero slope.

Similarly to the gravity-driven configuration, increasing the slope at a constant Shields number decreases the resistance of the granular bed to entrainment and one might expect a similar consequence on the transport rate scaling. This indicates that the effect of the slope in pressure-driven configurations would not only be restricted to the critical Shields number as suggested in the literature (Fernandez Luque & Van Beek, 1976; Fredsøe & Deigaard, 1992; Chiew & Parker, 1994; Dey, 2003; Andreotti et al., 2013), but also affects the whole transport rate formula, as evidenced by the experiments of Damgaard et al. (1997). In addition, similarly to the gravity-driven case developed in this paper, the importance of the fluid flow inside the granular layer cannot be accounted for in the classical framework in term of Shields number - i.e. dimensionless fluid shear stress at the top of the granular bed - and should be considered when modelling pressure-driven bedload transport. This has already been achieved and verified experimentally in laminar bedload transport and leads to a transport rate scaling law in terms of dimensionless fluid flow rate (Aussillous et al., 2013).

Besides, the gravity-driven and pressure-driven configurations do not exhibit the same dependency on the slope. In the former, the fluid flow inside the granular bed is driven by the slope (eq. 20), while in the latter it is driven by the pressure gradient and is independent of the slope (eq. 36). As a consequence, the two configurations are different in nature and not equivalent in terms of slope dependency. Therefore, care should be taken to consider configuration where the subsurface fluid flow contribution is negligible when applying results established for pressure-driven flows (e.g. Chiew & Parker (1994); Dey (2003)) to field prediction and experimental gravity-driven configurations (e.g. Li & Cheng (1999); Wilcock & Crowe (2003); Karmaker & Dutta (2016) among others).

5 Conclusion

Analysing turbulent bedload transport from a theoretical and numerical point of view, we attempted to clarify the mechanisms and origin of the slope influence and made a step toward a better understanding of the phenomenon. In particular, it has been evidenced that the classical modification of the critical Shields number relies on an expression of the buoyancy force which is not valid for uniform turbulent bedload transport. Focusing on gravity-driven configurations in steady uniform conditions, we have evidenced the entrainment mechanisms of the granular phase and shown the neglected importance of the fluid flow inside the bed. The relative importance of the latter with respect to the surface contribution characterised by the Shields number, affects the vertical structure of the granular flow and the sediment transport rate. The proposed modification of the Shields number to account macroscopically for these mechanisms has been shown to make all the present numerical data collapse onto a single master curve. It evidences that the key parameter in gravity-driven turbulent bedload transport is not the slope but the ratio between the tangent of the slope angle and the tangent of the critical debris flow angle. The difference is of importance when considering ideal particles for which the density (e.g. plastic) and the shape (e.g. spheres) are different from natural sediments. In addition, the theoretical analysis has evidenced a difference in nature between gravity-driven and pressure-driven configurations with respect to the slope influence. This result underlines the necessity to differentiate the two configurations in the analysis of slope influence on turbulent bedload transport.

The present contribution provides a better understanding of the slope influence in idealised configurations and represents a step in the understanding and prediction of the slope influence in turbulent bedload transport. In order to go further, it would be interesting to perform precisely controlled experiments and fully resolved simulations at varying slope inclination angle, to validate the present approach and characterise the importance of the neglected mechanisms such as the submergence, the lift force and the coherent structures.

Acknowledgement

We would like to thank the three anonymous reviewers for their interesting and constructive questions. We are grateful to Ashley Dudill for English corrections. RM would like to thank Thomas Pähtz, Alain Recking and Laurent Lacaze for fruitful discussions. This research was supported by Irstea (formerly Cemagref), the French national research agency project SegSed ANR-16-CE01-0005, the labex OSUG@2020 and the French Institut National des Sciences de l’Univers program EC2CO-BIOHEFECT and EC2CO-LEFE MODSED.

Appendix A Pressure-driven vertical granular flow structure

For the pressure-driven case, only the streamwise fluid momentum balance is modified, by replacing the driving slope term with a pressure gradient:

[TABLE]

Contrary to the gravity-driven case, there exists a pressure gradient in the streamwise direction so that the buoyancy force has a streamwise component which can be expressed as $-\phi\frac{\partial P^{f}}{\partial x}$ and leads to the following reformulation:

[TABLE]

Performing the same manipulation of the equations as for the gravity driven case, the shear-to-normal granular stress ratio can be shown to read:

[TABLE]

leading to the final equation:

[TABLE]

Appendix B Detailed characteristics of the numerical runs

Table 1 presents the detailed characteristics and the associated symbols of the simulations performed far from the threshold of motion, excluding the simulations of figure 5 which are secondary. The Reynolds number, $Re=Uh/\nu^{f}$, the Froude number $Fr=U/\sqrt{g_{z}h}$ and the relative submergence $h/d$ are evaluated considering the mean fluid velocity $U$ inside water depth $h$, defined as the difference between the position of the water free surface and the maximum of the Reynolds stresses. This is consistent with the definition adopted for the Shields number, based on the maximum Reynolds stresses and equivalent to $\theta^{*}=h\sin\alpha/[(\rho^{p}/\rho^{f}-1)d]$.