3.1. Steady fully developed flows and phase diagram

We consider steady flows over an erodible bed, i.e. flows whose properties are independent of time regardless of the fluctuations. To achieve such a flow, in addition to frictional sidewalls, a sufficient number of grains has to be used to obtain a flow rate greater than the critical one reported in Taberlet et al. (Reference Taberlet, Richard, Valance, Losert, Pasini, Jenkins and Delannay2003). At $t=0$ the angle $\theta$ is instantaneously decreased from $70^{\circ }$ to the desired value. After a transient corresponding to a global slow down of the system, the number of grains flowing out of the cell per unit time, $\partial _t N_{out}$, reaches a plateau, indicated in the figure by a dashed horizontal line. We have verified that this plateau corresponds to the occurrence of a SFD flow throughout the simulation cell. The value of the flow rate obtained at the plateau corresponds to the SFD flow rate. We have checked that this flow rate does not depend on the number of grains as long as the latter is large enough to obtain an erodible bed (Richard et al. Reference Richard, Valance, Métayer, Crassous, Louge and Delannay2010). In this case, if the number of grains is increased, the erodible bed becomes thicker but the flow zone and thus the flow rate do not change. In our simulations the creeping region corresponds to at least 10 grain layers.

Within the control parameters investigated so far ($20^{\circ } \leq \theta \leq 70^{\circ }$ and $5\leq W/d\leq 30$), we can obtain different regimes. If the angle is large enough, a SFD flow is reached after a transient, while at low angles the system remains completely jammed, i.e. without flow. This is illustrated in figure 3, where the phase diagram is reported. In addition, we have indicated the lower limit of the existence of steady flows: $\theta _1(W)$. Note also that, in vicinity of the latter, the temporal evolution of the flow rate may be intermittent, leading to a significant uncertainty in the determination of $\theta _1$. Indeed, an intermittent behaviour may be observed for a long period of time before the system comes to rest. The study of this intermittency is outside of the scope of this paper and would require an intensive statistical analysis, as done in Lemieux & Durian (Reference Lemieux and Durian2000). It is also important to emphasize here that the presence of flat frictional walls allows for the existence of SFD flows at much higher angles than those observed for unconfined flows (Silbert et al. Reference Silbert, Ertaş, Grest, Halsey, Levine and Plimpton2001). Besides, even at very steep angles the flow reaches a steady state (Taberlet et al. Reference Taberlet, Richard, Valance, Losert, Pasini, Jenkins and Delannay2003; Brodu et al. Reference Brodu, Richard and Delannay2013, Reference Brodu, Delannay, Valance and Richard2015). These results point out the crucial importance of sidewalls in the behaviour of granular flows. Additionally, as already been reported in the literature (Liu, Jaeger & Nagel Reference Liu, Jaeger and Nagel1991; Courrech du Pont et al. Reference Courrech du Pont, Gondret, Perrin and Rabaud2003; Metayer et al. Reference Metayer, Richard, Faisant and Delannay2010), one can note that the lower limit angle $\theta _1$ for steady states depends on the lateral width $W$: the wider the gap, the smaller $\theta _1$. This variation is the most important for $W<20d$ due to sidewall effects.

Figure 3. For a given gap between the sidewalls, a SFD flow over an erodible bed can only be observed above a given inclination angle. Below the latter, the system can experience an intermittent surface flow or be at rest. This is illustrated in the phase diagram reported in the parameter space $(W/d,\theta )$. The dashed line represents the angle $\theta _1(W)$ above which a SFD flow can be achieved.

3.2. Flow rate $Q$ vs $\theta$ and $W$

The flow rate per unit width, $Q$, is an important macroscopic quantity characterizing the ability of the system to flow. It is defined as $Q={\rm \Delta} m/({\rm \Delta} T W \rho _p)$, where ${\rm \Delta} m$ is the mass flowing out of the cell during the time ${\rm \Delta} T$ and $\rho _p$ the grain density. We have plotted in figure 4(a) the dimensionless flow rate per unit width, ${{Q}^{*}} = Q / d\sqrt {gd}$, vs the tangent of the angle of inclination $\theta$ for different gap widths $W$.

Figure 4. (a) For any gap between sidewalls $W$, the flow rate per unit width, ${{Q}^{*}}$, is found to increase with the angle of the free surface $\theta$. For a given $\theta$, increasing $W$ leads to an increase of ${{Q}^{*}}$. (b) The flow rate per unit width can be rescaled by $(W/d)^{5/2}$. The evolution of the rescaled flow rate vs $\tan \theta$ can be fitted by a quadratic law.

For a given angle of inclination, the flow rate per unit width ${{Q}^{*}}$ is an increasing function of the gap width $W$ as a consequence of the declining influence of the wall friction (Taberlet et al. Reference Taberlet, Richard, Valance, Losert, Pasini, Jenkins and Delannay2003; Jop et al. Reference Jop, Forterre and Pouliquen2005). For a given gap width, a one-to-one relation, which increases monotonically, exists between the flow rate and the angle. It is a consequence of sidewall friction and the increase of $\theta$ with ${{Q}^{*}}$ is more pronounced for low values of $W$ and is expected to be very weak for very large channels (Jop et al. Reference Jop, Forterre and Pouliquen2005). Our results are compatible with experimental ones (Khakhar et al. Reference Khakhar, Orpe, Andresén and Ottino2001; Courrech du Pont et al. Reference Courrech du Pont, Gondret, Perrin and Rabaud2003; Taberlet et al. Reference Taberlet, Richard, Valance, Losert, Pasini, Jenkins and Delannay2003; Jop et al. Reference Jop, Forterre and Pouliquen2005; Metayer et al. Reference Metayer, Richard, Faisant and Delannay2010) that report an increase of the angle of the flow with the flow rate. One can note again that the lower limit angle $\theta _1$ below which the flow stops is a decreasing function of the gap width (see § 3.1). The numerical evolution of ${{Q}^{*}}$ with $W$ is compatible with a $W^{5/2}$-scaling (figure 4b) consistent with that obtained by Jop et al. (Reference Jop, Forterre and Pouliquen2005) for much lower angles and larger gap widths. Also, its evolution with $\tan \theta$ can be fitted by a parabolic law $Q\propto (\tan \theta -\tan \theta _1)^{2}$.

4.1. Volume fraction profiles and definition of the flowing height

Unidirectional dense granular flows down bumpy inclines are characterized by an almost constant volume fraction throughout the flow height (Silbert et al. Reference Silbert, Ertaş, Grest, Halsey, Levine and Plimpton2001; Bi et al. Reference Bi, Delannay, Richard, Taberlet and Valance2005, Reference Bi, Delannay, Richard and Valance2006; Delannay et al. Reference Delannay, Louge, Richard, Taberlet and Valance2007; Brodu et al. Reference Brodu, Richard and Delannay2013). Figure 5(a) shows vertical profiles of the volume fraction, $\nu$, averaged across sidewalls for various inclination angles ($30^{\circ } \le \theta \le 60^{\circ }$) and lateral widths ($5d \le W \le 30d$). Contrary to dense flows sheared over the whole depth (hereafter referred to as fully mobilized flows) an erodible bed forms at the bottom of the system and, consequently, the volume fraction is, in most cases, far from being constant and varies from $0$ to approximately $0.6$. Only for low flow angles (e.g. $\theta =30^{\circ }$ for $W/d=30$) do the volume fraction profiles tend toward a step function that is a constant volume fraction. Remarkably, the vertical profiles of the volume fraction, $\nu$, for different inclinations and widths collapse onto a single curve (see figure 5b) which can be fitted by the following equation (Richard et al. Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008):

(4.1)\begin{equation} \nu(z)=\frac{\nu_0}{2} \left[ 1 +\tanh \left( \frac{2z}{l}\right) \right], \end{equation}

with the origin such that $\nu = \nu _0/2$ at $z=0$. The characteristic length of variation of the volume fraction $l$ depends on both the gap width and flow angle and, in most cases, it is greater than the grain size. It should also be pointed out that the validity of (4.1) is questionable at low values of the characteristic length $l$ due to the possible intermittency of the flow. The length $l$ scales with $W$ and increases with inclination and (1.2) can be rewritten as

(4.2)\begin{equation} \frac{l}{W}=\frac{1}{\mu_{w,l}}(\tan\theta - \mu_{b,l}). \end{equation}

We have verified this equation for $W/d = 5,10,20$ and $30$ and obtained $\mu _{b,l} = 0.32$ and $\mu _{w,l}= 0.43$ (see figure 5c). Note that the characteristic length $l$ used here is twice the length $l_\nu$ used by Richard et al. (Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008).

Figure 5. The vertical profiles of the volume fraction for various inclination angles $\theta$ and widths $W$ (a) collapse onto a master curve $\nu (z) = ({\nu _0}/{2})[ 1 +\tanh ( {2z}/{l}) ]$ – dashed line in panel (b). The characteristic length $l$ scales with $W$ and increases with the angle of inclination $\theta$ (c). The dashed line represents the best linear fit for all the data. At low angle a slight deviation from the linear relation is found (inset of (c)). The angle for which this deviation is visible is close to $\theta _c=40^{\circ }$ and $l/W\approx 1$. For this reason, the linear fit reported is done for angles above $\theta =40^{\circ }$.

The variation of $l/W$ with $\tan \theta$ deviates from an affine law for $\theta < \theta _c \approx 40^{\circ }$; note that $\theta _c$ is the angle for which $l/W\approx 1$. As will be shown below, $\theta _c$ corresponds to the boundary between a dense flow regime (regime I) and a regime for which a strong variation of the volume fraction is observed within the flow (regime II). It is not easy to say whether the value of $\theta _c$ depends on $W$. There is likely a dependency, but it seems rather weak, with a small increase of $\theta _c$ (only a few degrees) when $W$ decreases. A similar but more pronounced effect was observed for $\theta _1$ (see figure 3).

4.2. Contact number

Another interesting property of a granular flow is the average number of contacts per grain, $N_c$. We have illustrated its variation throughout the flow and its link with the volume fraction in figures 6(a) and 6(b), respectively. As expected, it increases with increasing depth from zero in the gaseous zone to $N_c\sim 4$ within the creeping zone, which corresponds to the number of contacts obtained for an iso-static packing, i.e. a packing for which all the contacts are perfectly frictional. We have explained in § 2 that the initial state of the system might influence the reached steady state. The protocol used in the present work indeed leads to relatively loose (and fragile Berzi, Jenkins & Richard Reference Berzi, Jenkins and Richard2019) systems but others might possibly lead to denser systems and thus to a number of contacts larger than $4$ in the creeping zone. It should be pointed out that these quantities ($N_c$ and $\nu$) vary with the same length scale $l$ but they are shifted vertically: in the example reported in figure 6 the number of contacts reaches the half of its maximal value at $z/d\approx 22$ (depth for which $\nu$ tends toward $0.6$) whereas $\nu =\nu _0/2$ at $z/d=0$ (depth for which $N_c$ is below $0.1$). Thus, in the flowing zone ($z\in [-l/2,l/2]$), the flow is mainly collisional. When we plot the number of contacts as a function of the volume fraction, two different exponential regimes can be clearly identified: one for $\nu \lesssim 0.53$ for which the average number of contacts is very small (below $0.1$) and the other one for $\nu \gtrsim 0.53$ where the average number of contacts increases rapidly with $\nu$. As $N_c$ and $\nu$ vary over the same length scale $l$, the relationship between these quantities is independent of $\theta$ and $W$: the regime change corresponds to the transition between the buffer zone and the flowing zone which occurs at $z=l/2$, i.e. at $\nu (z=l/2)=\nu _0/2[ 1 + \tanh (1)]\approx 0.528$. Note also that this curve could depend on the grain stiffness $k^{n}$ used to model the normal force between grains (see § 2). It has indeed been shown that this quantity plays an important role on the number of contacts within granular flows, which increases with decreasing stiffness (Ji & Shen Reference Ji and Shen2008; Vescovi & Luding Reference Vescovi and Luding2016; Vescovi, Berzi & di Prisco Reference Vescovi, Berzi and di Prisco2018). The volume fraction can also be influenced by the grain stiffness (Vescovi & Luding Reference Vescovi and Luding2016) especially in dense zones for which the overlap between contacting grains might have a non-negligible volume. A full study of the effect of grain stiffness on our results is outside the scope of the paper. Yet, to position our study with respect to that of Vescovi & Luding (Reference Vescovi and Luding2016) we have estimated the dimensionless parameter $K^{*}=k^{n}/\rho _p d^{3} \dot \gamma$, $\rho _p$ being the grain density. Since $\dot \gamma$ varies within our system, so does $K^{*}$. Yet the order of magnitude of its lowest value is $10^{6}$. For such a value, Vescovi & Luding (Reference Vescovi and Luding2016) have found an effect of the grains’ stiffness on the volume fraction only if the latter exceed $0.625$, which is never the case in our simulations.

Figure 6. (a) Variations of the number of contacts, $N_c$, and of the volume fraction $\nu$ with the flow depth. (b) The gap between sidewalls is $W=10d$ and the angle of the flow $\theta =50^{\circ }$. The grain–grain and grain–wall friction coefficients are both equal to $0.5$. The variation of $N_c$ vs $\nu$ displays a universal behaviour with a kink corresponding to $\nu \approx 0.54$ and $N_c\approx 0.06$ (b).

4.3. Velocity profiles

The vertical profiles of the streamwise velocity are reported in figure 8(a) for various gaps $W$ between the sidewalls and various flow angles $\theta$. As expected, the two latter parameters have a strong influence: the velocities are found to increase with increasing angle and/or gap. Interestingly, at large angles, the velocity profiles display an inflection point close to $z/d=0$, around which the profile is reasonably linear. In the vicinity of the free surface, two cases are observed: either a strict maximum or a plateau. For a given $W$ we observe the plateau at low angles. When the angle is larger than a given value (which decreases strongly with increasing $W$) a strict maximum is observed. In the remainder of the paper, depending on the case, we will call $V_{max}$ the velocity of the strict maximum or of the plateau.

Following Richard et al. (Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008), we rescale the velocities and the vertical positions respectively by $\sqrt {gW}(l/d)$ and $l$ (figure 8c). The scaling holds for angles larger than $\theta _c$ (regime II). Interestingly, the maximum velocity $V_{max}$, which corresponds approximately to $V(z=-l/2)$, is found to scale always (in regime I as in regime II) with the characteristic velocity $\sqrt {gW}(l/d)$ (figure 8b) which, consequently, can be considered as a characteristic velocity of the system. Figure 8(c) also points out that $l$ is not the right characteristic length for the variation of velocity with depth in regime I. We have thus to introduce a new characteristic length. Taberlet et al. (Reference Taberlet, Richard, Valance, Losert, Pasini, Jenkins and Delannay2003) defined – somewhat arbitrarily – the flow height as twice distance between $z=0$ (i.e. the position of the inflexion point) and the depth at which the tangent to the linear part of the velocity profile intersects the vertical axis ($z$ axis). Taking into account the fact that the velocity profile is not symmetric according to $z = 0$ in regime I (see also figure 12) and the characteristic behaviour of $V_{max}$, we define here the flow height $h$ as the distance between $-l/2$ (position corresponding approximately to $V_{max}$) and the depth at which the tangent to the linear part of the velocity profile intersects the $z$ axis (see figure 7). Rescaling the depth by $h$ (inset of figure 8c), leads to profiles in regime I which depend weakly on the inclination angle $\theta$, as predicted by the $\mu (I)$ rheology (Jop et al. Reference Jop, Forterre and Pouliquen2005), and are not on the master curve of regime II which is independent of $\theta$.

Figure 7. Velocity (blue open circles) and volume fraction profiles (grey line) and its fit by (4.1) (black dashed line) for $W/d=10$ and $\theta =53^{\circ }$. The height of the flow $h$ is the distance between $z=-l/2$ (orange dashed vertical line, position corresponding approximately to $V_{max}$) and the depth at which the tangent to the linear part of the velocity profile (dashed red line) intersects the $z$ axis.

Figure 8. (a) Vertical profiles of the particle velocity in the main flow direction for various inclination angles $\theta$ and widths $W$. (b) The maximal velocity (estimated by $v_x(z=l/2)$) scales with the characteristic velocity $\sqrt {gW}(l/d)$. (c) Flows in regime II collapse onto a single curve when rescaling the velocity by $\sqrt {gW}l/d$ and the depth by $l$ or by $h$ (inset of (c)). (d) In regime II, the characteristic length $l$ determined by fitting the volume fraction profiles is close to $h$, the length determined through the velocity profile. In regime I, below $\theta _c \approx 40^{\circ }$ important relative differences are observed (c).

Note that the flow obtained for $W/d=10$ and $\theta =40^{\circ }$ clearly belongs to regime II. In contrast, that obtained for $W/d=5$ and $\theta =40^{\circ }$ seems to belong to the end of regime I, it corresponds probably to a transition between regimes I and II. This suggests a slight dependence of the boundary between the two regimes (and thus of the angle $\theta _c$) on the gap between the sidewalls.

As mentioned before, a peak in the velocity profile is observed for large angles and it should be pointed out that the scaling does not hold post-peak (i.e. for $z$ lower than the peak's position). This is the case, for instance, for ($\theta =60^{\circ },\ W/d=10$) and for ($\theta =45^{\circ }, W/d=20$). At the top of the flow, we observe that grains are mainly submitted to collisions with sidewalls due to the very low volume fraction. The latter being immobile, the corresponding grains are more strongly slowed down than the grains located in the bulk which collide mainly with grains that move in the main flow direction. Also the grains located at the top of the flow slow down the grains located directly below them that, in return, also slow down the grains located below them. This global slow down ceases when the probability for a grain to collide with sidewalls is not significant. Thus, the behaviour of the granular flow in the vicinity of its top is complex and strongly influenced by the gap between sidewalls.

Figure 8(d) compares the two characteristic lengths of the flow: the characteristic length of the velocity profile $h$, and the characteristic length of the volume fraction profile $l$. At large angles, in regime II, $h/W$ and $l/W$ look identical. In regime I, below $\theta _c \approx 40^{\circ }$, $h/W$ becomes clearly larger than $l/W$, but the two lengths remain correlated. Thus, in regime I, $l$ does not correspond anymore to the flowing friction height. Since $h$ characterizes the friction height for both regime I and regime II, we have reported in figure 9(a) the evolution of the tangent of the flow angle vs the ratio $h/W$. The SSH equation (1.1) is rather well verified. In regime I, the values of $l/W$ and $h/W$ are significantly different. Consequently, that of $\mu _{b,h}$ is different from that of $\mu _{b,l}$ (respectively $0.22$ and $0.32$) and that of $\mu _{w,h}$ from that of $\mu _{w,l}$ (respectively $0.51$ and $0.43$).

Figure 9. (a) The tangent of the angle of the flow depends linearly on the ratio $h/W,$ where $W$ is the gap between the sidewalls and $h$ the depth of the flow determined by fitting the velocity profile (see text). (b) The dimensionless flow rate per unit width ${{Q}^{*}}$ scales with $(W/d)^{5/2}$ at any angle. Its dependence with $h/W$ is more complex: it is compatible with a $(h/W)^{7/2}$ scaling for angles below $\theta _c\approx 40^{\circ }$ (regime I) and a $(h/W)^{2}$ scaling above (regime II). (c) Similarly, the maximal value of the velocity $V_{max}$ scales with $(W/d)^{3/2}$ at any angle and with $(h/W)^{5/2}$ ($\theta < \theta _c$, regime I) and $h/W$ ($\theta >\theta _c$, regime II). (d) The difference of scaling for the two regimes can also be confirmed by plotting $V_{max} d / (W\sqrt {gW})$ vs ${{Q}^{*}}/(W/d)^{5/2}$: the former quantity scales with the latter to the power $5/7$ and $1/2$ respectively for $\theta <\theta _c$ and $\theta >\theta _c$.

Interestingly, when $l/W\to 0$, $h/W\to (h/W)_{min} > 0$. For a given $W$ there is thus a finite minimum value of $h$: $h_{min}(W)$, and this minimum value increases proportionally to $W$: $h_{min}(W) = W (h/W)_{min}$. Note that the limit $l/W \to 0$ corresponds in our case to unconfined geometries ($W\to +\infty$), for which SSH flows have never been observed. The equation $h_{min}(W) = W (h/W)_{min}$ is consistent with the latter observation since $W\to \infty$ leads to $h_{min}\to \infty$. In other words, in unconfined geometries, SSH flows are only observed for infinitely high systems. When $l/W\to 0$, the SSH equations (1.1) and (1.2) lead to the relation: $\mu _{b,l}=\mu _{b,h}+\mu _{b,w} (h/W)_{min}$, from which we obtain $(h/W)_{min}=0.33$.

It should be pointed out that $\tan \theta _1 > \mu _{b,l} > \mu _{b,h}$. All these quantities are probably close to each other when $W \to \infty$. Note that $\theta _1$ is obtained directly from simulations in which the system transits from a dynamical state to a static one, whereas the two others are obtained from fits of numerical data obtained from the simulation of dynamical states. In the literature, $\mu _{b,h}$ is classically assumed to be gap width independent (Liu et al. Reference Liu, Jaeger and Nagel1991; Boltenhagen Reference Boltenhagen1999; Courrech du Pont et al. Reference Courrech du Pont, Gondret, Perrin and Rabaud2003; Metayer et al. Reference Metayer, Richard, Faisant and Delannay2010). Yet, thanks to our extensive simulations, which provide results for wide ranges of both gap widths and flow angles, a small but systematic increase of $\mu _{b,h}$ with $W$ is found. In contrast, the angle $\theta _1$ decreases with increasing $W$ (see the phase diagram reported in figure 3). This difference points out that these two quantities do not quantify the same physical mechanisms.

Following Jop et al. (Reference Jop, Forterre and Pouliquen2005) we have studied the scaling of both ${{Q}^{*}}$ (figure 9b) and $V_{max}$ (figure 9c) vs the gap $W$ between the sidewalls and the height $h$ of the flowing layer. In agreement with Jop et al. (Reference Jop, Forterre and Pouliquen2005), the dimensionless flow rate and the maximum velocity respectively scale with $(W/d)^{5/2}$ and $(W/d)^{3/2}$ whatever the angle of the flow. Yet, their dependence with $h/W$ confirms the existence of two regimes. At low angles, in regime I, in agreement with the results reported by Jop et al. (Reference Jop, Forterre and Pouliquen2005), ${{Q}^{*}}$ and $V_{max}/\sqrt {gd}$ scale respectively as $(h/W)^{7/2}$ and $(h/W)^{5/2}$. For larger angles, in regime II, ${{Q}^{*}}$ and $V_{max}$ scale as $(h/W)^{2}$ and $h/W$, respectively (Richard et al. Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008). We have seen before that $V_{max}/\sqrt {gd}$ scales as $l/W$ in regime I as in regime II, this shows that $h/W$ scales as $(l/W)^{2/5}$ in regime I, and thus, ${{Q}^{*}}$ scales as $(l/W)^{7/5}$. The difference of scaling for the two regimes can also be confirmed without using $h$ or $l$, for example, by plotting $V_{max} d / (W\sqrt {gW})$ vs ${{Q}^{*}}/(W/d)^{5/2}$ (see figure 9d). In agreement, with previous results, we indeed recover the existence of different scaling exponents for the two regimes, i.e. $V_{max} d / (W\sqrt {gW})$ scales with $[{{Q}^{*}}/(W/d)^{5/2}]^{5/7}$ in regime I and with $[{{Q}^{*}}/(W/d)^{5/2}]^{1/2}$ in regime II. The values of the latter exponents are consistent with the scaling laws obtained for $V_{max}$ and ${{Q}^{*}}$ vs $h/W$.

The velocity profiles are not linear both in the creeping and the dilute zones. Consequently the strain rate $\dot \gamma =\partial v_x / \partial z$ depends on $z$. Yet, the velocity profiles have a nearly linear part in the vicinity of $z\approx 0$, region which corresponds to a change of concavity for the shape of the profiles. In the remainder of the paper, we will estimate the strain rate at the depth for which the concavity of the vertical velocity profile changes (i.e. at the depth $z$ for which the strain rate is maximal corresponding to $z/d\approx 0$ for $\theta >\theta _c$). Figure 10(a) shows the evolution of the strain rate with the angle of the flow and for several cell widths $W/d$. For any value of $W/d$, the strain rate is found to increase with the flow angle. For $\theta >\theta _c\approx 40^{\circ }$ the increase weakens and the slope of $\dot \gamma$ vs $\tan \theta$ decreases with decreasing $W/d$, illustrating once more the existence of the two regimes. When plotting $\dot \gamma d h / \sqrt {gW} l$ vs $\tan \theta$ a reasonable collapse is obtained (figure 10b). In the creeping zone, there is experimental evidence (Komatsu et al. Reference Komatsu, Inagaki, Nakagawa and Nasuno2001; Crassous et al. Reference Crassous, Metayer, Richard and Laroche2008; Richard et al. Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008; Martínez et al. Reference Martínez, González-Lezcano, Batista-Leyva and Altshuler2016) that the mean velocity of the grains decays exponentially with depth and the characteristic decay length is approximately equal to the particle size. The dependence of this length with the flow rate (and thus the angle of the flow) and the gap width remains unclear despite some experimental attempts (Richard et al. Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008). Our numerical simulations indeed capture the exponential decay of the grain velocity in the creeping zone (figure 11) over several decades. We found that the characteristic length of the decay, $\lambda$, which is of the order of a grain size, depends on both the gap width and the inclination angle. It should be pointed out that our range of flow rates is much more important than those used in the literature (Komatsu et al. Reference Komatsu, Inagaki, Nakagawa and Nasuno2001; Crassous et al. Reference Crassous, Metayer, Richard and Laroche2008). Combining volume fraction measurements by $\gamma$-ray adsorption and imagery, Richard et al. (Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008) have shown that both the height of the flowing zone (i.e. in our case $h$) and the characteristic length $\lambda$ are linearly correlated. Our simulations capture this linear relation for both regimes and the slope obtained (${\approx }6$) is close to that reported in Richard et al. (Reference Richard, Valance, Métayer, Sanchez, Crassous, Louge and Delannay2008).

Figure 10. (a) The strain rate increases with the gap between the sidewalls and with the angle of the flow whatever the gap width. Yet, for angles greater that $\theta _c\approx 40^{\circ }$, the increase is reduced, confirming the existence of two regimes: $\theta < \theta _c$ and $\theta >\theta _c$. (b) The data can be approximately rescaled on a single master curve by rescaling the velocity by $\sqrt {gW}$.

Figure 11. In the creeping zone, the velocity decays exponentially with the depth. This is illustrated in (a) for $W/d=10$. (b) Rescaling the velocity by $\sqrt {gW}l/d$ and the depth by $h$ leads to a collapse of the velocity profiles in the creeping region for $\theta > \theta _c$ (b). The characteristic length of this decay, $\lambda$, is found to be correlated to the length $h$ (c).

4.4. Flux densities

Another quantity of interest is the flux density defined as the product of the volume fraction and the streamwise velocity: $\phi (z)=\nu (z)v_x(z)$. It quantifies the number of grains crossing a unit area whose normal is along the streamwise direction. This quantity is widely used in the context of aquatic and aeolian sediment transport but much less commonly in the physics of granular media.

Figure 12(a), which depicts the vertical profiles of this quantity, shows that in regime II, the flux densities are almost symmetric with respect to $z/d=0$ i.e. $\nu =\nu _0/2$, where they are maximum. Consequently, contributions to the flow rate of the low density ($\nu <\nu _0/2$) and of the high density region ($\nu >\nu _0/2$) are similar, the former being more rapid than the latter. We have seen that the characteristic lengths for the variation of the volume fraction ($l$) and the velocity ($h$) are close to each other in regime II. Consequently, the flux density being the product of the volume fraction by the velocity, $h$ (or equivalently $l$) is the relevant length to describe the variation of this quantity in regime II (figure 12b).

Figure 12. (a) In regime II, the flux densities are maximum at $z\approx 0$. (b) For the same regime, the vertical profiles of the flux densities collapse on a single curve if the velocity and the depth are respectively scaled by $\sqrt {g W} l/d$ and $h$.

We can approximate, in the vicinity of $z=0$, the volume fraction profile by $\nu (z)\approx \nu _0/2 [1+2z/l]$ and the velocity by $v_x(z)/ (\sqrt {gW}l/d)\approx \kappa (1-z/h)$ with $\kappa =\dot \gamma d/\sqrt {gW}\approx 0.4$ (where $\gamma$ is calculated at $z=0$). Note that the value of $\kappa$ is consistent with the values of $\dot \gamma d /\sqrt {gW}$ reported in figure 10(b). Using the two latter expressions and the approximation $h\approx l$ in regime II allows us to approximate the flux density $\phi (z)/(\nu _0 \sqrt {gW}l/d )$ by the parabola $\nu _0\kappa (1-2z^{2}/h^{2})$ in the vicinity of $z=0$. This parabola-shaped profile is characteristic of heap flows at large angles. Note that the contribution of the flowing layer ($z\in [-l/2,l/2]$) to the global flux is the most important. The contribution of the other zones is, however, far from negligible (${\approx }1/3$ of the global flux).

For angles below $\theta _c$, the velocity profile is not symmetric with respect to $z/d=0$, and the length $h$ characterizing the height of the flowing layer is larger than the length $l$ characterizing the variation of the packing fraction with depth. The maximum of the flux density is shifted towards higher volume fractions. The volume fraction variation is indeed more localized close to $z/d=0$ and tends toward a step function. We thus recover the behaviour observed for SFD dense flows on a bumpy bottom (Silbert et al. Reference Silbert, Ertaş, Grest, Halsey, Levine and Plimpton2001; Bi et al. Reference Bi, Delannay, Richard, Taberlet and Valance2005, Reference Bi, Delannay, Richard and Valance2006) for which the volume fraction is almost constant and, consequently, the flux density is maximum when the velocity is maximum.

4.5. Rotation

We have also studied the rotation of the grains within the flow which is a quantity of interest since, in confined systems, geometrical frustration of rotation consequence of geometrical constraints is frequent. This frustration of rotation has an important influence on the organization of the particles’ rotations and thus on the rheology of the system (Yang & Huang Reference Yang and Huang2016). The mean angular velocity $\omega _y$ (i.e. normal to the sidewalls) profiles are reported in figure 13(a). They are counted positive in the counter-clockwise direction. Kinetic theories (Jenkins & Richman Reference Jenkins and Richman1985; Lun & Savage Reference Lun and Savage1987) assume that the mean angular velocity can be approximated by the rotation rate across the flow depth: $\omega _y \approx - 0.5 \partial v_x / \partial z = 0.5 \dot \gamma$. In our system, this assumption is valid from the creeping zone up to $z/d\approx 0$, i.e. for $z/d>0$ and a volume fraction above $\nu _0/2$. In contrast, the agreement between the two quantities is poorer for $z/d<0$, i.e. for dilute flows. This discrepancy, which seems to be independent of $\theta$ and $W$, can be explained by the fact that, for $z/d>0$, the average number of contacts per grain is almost zero (see § 4.2). The probability of collision with the sidewalls becomes more significant as $z$ decreases. Consequently, their velocity distributions are not classical anymore and the relation between $\dot \gamma$ and $\partial v_x / \partial z$ does not hold for these grains. Due to the correspondence between $\omega _y$ and $\dot \gamma$, it is natural to represent the vertical profile of $\omega _y$ by scaling the latter quantity by $\sqrt {gW} l / hd$ and $z$ by $h$ (see figure 13b). As expected, the regimes I and II are recovered due to the link between $\dot \gamma$ and $\omega _y$.

Figure 13. (a) Vertical profiles of rotation velocity in the $y$ direction for various inclination angles $\theta$ and widths $W$. The dashed lines correspond to vorticity $- 0.5 \partial v_x / \partial z$. (b) The scaling with $d h/ (l \sqrt {gW}$) holds for regime II.

4.6. Temperature profile

The velocity fluctuations are characterized by the granular temperature, which is defined by $T = (T_{xx}+T_{yy}+T_{zz})/3$, where $T_{ij}=\langle v_iv_j \rangle - \langle v_i \rangle \langle v_j \rangle$, with $i, j = \{x, y, z\}$. The granular temperature is a measure of the agitation of the grains. One of the main applications of the concept of granular temperature is the derivation of kinetic and hydrodynamic theories for granular systems. It is thus an important quantity that deserves to be studied. Figure 14(a) reports the vertical profile of the granular temperature for several angles of inclination and $W/d=5$. The granular temperature is found to increase from the creeping zone (which behaves like a heat sink) to the flowing region. More details on the temperature can be found in § 5.4 where its properties along the transverse direction are reported. Note that the temperature profiles reported here are highly different from those obtained for SFD flows on a bumpy bottom, which behaves like a heat source (Silbert et al. Reference Silbert, Ertaş, Grest, Halsey, Levine and Plimpton2001). Consequently, in the latter case, the granular temperature is maximum in the vicinity of the bumpy bottom (Silbert et al. Reference Silbert, Ertaş, Grest, Halsey, Levine and Plimpton2001; Louge Reference Louge2003). It should be pointed out that the quantities $T_{xx}$, $T_{yy}$ and $T_{zz}$ behave similarly to $T$ (i.e. they increase from the creeping region to the free surface). Note that $T_{yy}\approx T_{zz}$ and that $T_{xx}$ is always larger than $T_{zz}$ and $T_{yy}$ ($T_{xx}/T_{yy}\approx 1$ in the creeping region, ${\approx }3$ in the flowing zone and up to $10$ close to the free surface). Similarly to what has been observed in the case of thestreamwise velocity profile $v_{x}(z)$, the temperature profiles can be rescaled using $g W l^{2}/d^{2}$ as a temperature scale, and $h$ as a depth scale (figure 14c).

Figure 14. (a) For a given channel width (here $W=5d$) the vertical profiles of fluctuations of the particle velocity increase with $\theta$ the angle of the flow. The vertical profiles plotted for various widths and angles can be rescaled using $l$ as a length scale and $h$ as a depth scale similarly to what is observed for the velocities (c). A zoom corresponding to low granular temperatures is reported in the inset of (c).

The kinetic theory for fast flows (Haff Reference Haff1986) assumes that fluctuations should depend on local relative motion and, consequently, predicts $T_{xx}\propto \dot \gamma ^{\chi }$ with $\chi =2$. This prediction has been experimentally verified for dilute and rapid chute flows (Azanza et al. Reference Azanza, Chevoir and Moucheront1999). Yet this scaling is no longer valid in the case of denser and slower flows where an exponent $\chi$ close to $1$ has been reported (Losert et al. Reference Losert, Bocquet, Lubensky and Gollub2000; Mueth Reference Mueth2003), a signature of non-local effects (Jenkins Reference Jenkins2007; Kamrin & Henann Reference Kamrin and Henann2015). In a geometry quite similar to ours, Zhang et al. (Reference Zhang, Yang, Lin, Chen and Yang2019) have recently reported a relation $T\propto \dot \gamma ^{\chi }$ with $\chi \approx 1.57$ in the creeping zone. Similarly, Artoni & Richard (Reference Artoni and Richard2015b) also reported data consistent with $\chi \approx 3/2$ for confined shear flows subjected to a constant pressure, which are quite similar to the flows studied in the present work. In the case of surface flows in a rotating cylinder (Orpe & Khakhar Reference Orpe and Khakhar2007), which are steady but not fully developed, the exponent $\chi$ increases with the rotational speed of the cylinder from $\chi =1$ to $\chi =2$. Figure 15 shows the variation of $T_{xx}$ with the strain rate. The reported relation between the latter quantities is consistent with a power law $T_{xx}\propto \dot \gamma ^{\chi }$ but with $\chi$ varying between $\chi \approx 1$ (in the bottom of the creeping zone) and $\chi \approx 2$ (in the buffer and the flowing zones). In agreement with the literature (Artoni & Richard Reference Artoni and Richard2015b; Zhang et al. Reference Zhang, Yang, Lin, Chen and Yang2019), if the whole curve (with the exception of the part corresponding to the gaseous zone) is fitted by a power law, an exponent $\chi$ close to $1.5$ is found. Remarkably, the variation range of $\chi$ is independent of both the inclination angle $\theta$ and the gap between the sidewalls $W$. This result is in agreement with that obtained in the case of flows in rotating drums (Orpe & Khakhar Reference Orpe and Khakhar2007) for which we also have $T_{xx} \propto \dot \gamma ^{\chi }$ with $\chi \in [1,2]$.

Figure 15. Variation of the dimensionless temperature $T_{xx} d / g W l^{2}$ with the strain rate ${\dot \gamma h d}/{\sqrt {gW}} {l}$ in the flowing layers. The dashed black line shows that the relationship can be approximated by a power law of exponent $\chi = 1$ at the bottom the creeping region and by $\chi =2$ in the buffer and the flowing zones. The whole curve (excepting the part corresponding to the gaseous zone) can be fitted in a satisfactory way by a power law of exponent $\chi =1.5$ (not shown).

5.1. Volume fraction

Since the middle of the preceding century (Brown & Hawksley Reference Brown and Hawksley1946; Verman & Banerjee Reference Verman and Banerjee1946), it has been well known that the presence of sidewalls influences the volume fraction of static granular materials in their vicinity (see Camenen, Descantes & Richard (Reference Camenen, Descantes and Richard2012) and references therein). More precisely, they lead to a reduction in volume fraction due to the wall-induced structure which extends within the packing. This result can be understood with a simple description incorporating two regions: (a) an effective boundary layer (which extends from each sidewall to a given length $\varDelta _\nu$) and (b) a bulk-like region. The latter region is assumed to have a volume fraction equal to that of an infinite system whereas the volume fraction of the former region is lower.

In our system, a lower volume fraction in the vicinity of the sidewalls is still expected. The creeping zone is not far from being static, the geometrical constraint induced by the latter remains important. In the flowing zone, as will be shown in § 5.4, the transverse profile of the temperature implies a depletion of grains in the vicinity of the sidewalls. We use a description similar to the static modelling with an effective boundary layer of size ($\varDelta _\nu$) depending on the depth within the packing.

We have reported in figure 16 the transverse volume fraction profiles for two flow depths, one in the flowing region ($z/d=0$) and the other one in the creeping zone, for a flow angle equal to $45^{\circ }$ with $W=10d$ (figure 16c) and $W=20d$ (figure 16d). The solid lines correspond to data averaged over a length of $0.1d$. They clearly reveal the layering induced by the flat walls. It is especially important in the creeping zone, but does not disappear completely in the flowing zone. To get free from the fluctuations induced by this layering, we filtered the data to produce smoothed volume fraction profiles corresponding to the dashed lines of the figure. The presence of sidewalls clearly leads to a decrease of the values of the smoothed volume fraction in their vicinity. This effect is stronger in the flowing zone that in the creeping zone. It should be pointed out that even a small modification of the volume fraction has important consequences for the global macroscopic behaviour of a granular system. In the creeping zone, the value of $\varDelta _\nu$ seems independent of $W$ and is of the order of 3 or $4d$. For $W=10d$ the volume fraction in the flowing region ($z/d=0$) is strongly influenced by the sidewalls: it increases continuously from the sidewalls to the centre of the cell. In contrast, for $W=20d$, also for $z/d=0$, a plateau of volume fraction seems to be reached between $y=-W/4$ and $y=W/4$.

Figure 16. Effect of the sidewalls on the volume fraction profiles for $\theta =45^{\circ }$. At any depth within the flow for $\theta =45^{\circ }$, and for any value of $W$, the volume fraction is smaller in the vicinity of the sidewalls with respect to that measured in the centre of the cell ((a) and (b) for respectively $W/d=10$ and $W/d=20$). The effect of the walls is visible on the transverse volume fraction profiles especially in the flow region ((c) and (d) for respectively $W/d=10$ and $W/d=20$). The solid lines correspond to data averaged over a length of $0.1d$ and the dashed lines correspond to the same data filtered by a Savitzky–Golay filter (window length $4d$ and polynomial degree $3$). The volume fraction can be considered as uniform only for large $W/d$ and far from the sidewalls.

We have reported in figure 16(a,b) the vertical profiles of the smoothed volume fraction in the middle of the cell and at the sidewalls for the same configurations. As expected, the volume fraction is always smaller in the vicinity of the sidewalls. This also confirms that this effect is stronger in the flowing zone.

It should be pointed out that, for $W=20d$ and $\theta =45^{\circ }$, around $z/d=-10$, an unexpected increase of the volume fraction (a bump) is observed in the centre of the cell in the vicinity of the free surface. This is due to the formation of dense clusters of a few grains. The exact characterization of these clusters is beyond the scope of the present paper.

5.2. Velocity profiles

Bumpy sidewalls are generally considered as non-sliding boundaries (Jop, Forterre & Pouliquen Reference Jop, Forterre and Pouliquen2006) with a vanishing velocity at the sidewalls. Here, the sidewalls are flat, the granular material is expected to have an average velocity that is not negligible at their contact. We have reported in figure 17(a,b) the vertical profiles at sidewalls and in the centre of the simulation cell for $\theta =45^{\circ }>\theta _c$, with $W=10d$ and $W=20d$, respectively. As expected, an important longitudinal velocity at the sidewalls is observed. We call it the sliding velocity and we denote it as $v_s(z)$ at depth $z$. Our simulations show that the transverse profiles of velocity are similar in regimes I and II. The nature of the regime has an effect on the evolution of the sliding velocity with the angle, but the behaviour of the reduced velocity at a given depth $z$, defined as $(v_x(y,z)-v_s(z))/(v_{x,max}(z)-v_s(z))$, where $v_{max}({z})$ is the maximum value of $v_x(y,z)$ at depth $z$, independent of the regime.

Figure 17. Effect of sidewalls on the velocity profiles for a flow angle equal to $45^{\circ }$ i.e. greater that $\theta _c$, the angle corresponding to the transition between regime I and regime II. For small gap width ((a) $W/d=10$) the difference between the vertical velocity profiles at sidewalls and in the centre of the simulation cell is small and increases with increasing $W$ ((b) $W/d=20$). The influence of the sidewalls on the transverse profile of the reduced velocity at a given depth $z$ (i.e. $(v_x-v_s)/(v_{x,max}-v_s)$ where $v_{x,max}$ and $v_s$ are respectively the maximum velocity and the sliding velocity at the aforementioned depth) is weak but increases from the flowing zone to the creeping one. This variation with depth also increases with increasing $W$ (c,d).

Transverse profiles of reduced velocity are reported in figure 17(c,d). They are found to broaden (i) from the creeping region to the flow region and (ii) with $W/d$ (see figure 17c,d). These results highlight the three-dimensional structure of the flows observed in our system.

The transverse velocity profiles can be fitted in a convenient way by $v_x(y,z)=v_{max}(z) - B[1 - \cosh (y/\varDelta _{v_x}) ]$ with $B=[v_{max}(z) - v_s(z)]/[1 - \cosh (W/2\varDelta _{v_x})]$. In the latter expression, $\varDelta _{v_x}$ is the length characterizing the influence of the sidewalls (Courrech du Pont et al. Reference Courrech du Pont, Gondret, Perrin and Rabaud2003; Zhang et al. Reference Zhang, Yang, Lin, Chen and Yang2019). We have determined $v_s$ and $\varDelta _{v_x}$ by fitting the transverse profiles of the velocities averaged between the depths $-l/2$ and $l/2$. The ratio $v_s/v_{max}$ is found to always increase with the flow angle (figure 18), but it increases much more moderately in regime II than in regime I. Interestingly, the effect of $W$ on this ratio is weak. This seems to indicate that a boundary condition on the rescaled sliding velocity, $v_s/v_{max}$, can be expressed, in a first approximation, as a function of $\theta$.

Figure 18. The longitudinal sliding velocity at the sidewalls, $v_s$, increases with the flow angle $\theta$. Interestingly, when rescaled by the maximum velocity, the latter increase is weakly dependent on the gap between the sidewalls $W/d$. As shown on the inset, the length characterizing the influence of the sidewalls (see text) increases with $\tan \theta$.

Interestingly, the characteristic length $\varDelta _{v_x}$ is much greater than $W/2$ for $W=5d$, indicating that the sidewalls influence the flow across the total width of the channel and, consequently, its determination is imprecise. For $W=20d$ and $W=30d$, the characteristic length is smaller than $W/2$ (inset of figure 18) indicating that the effect of the sidewalls is limited to a portion of the channel. The inset also shows that $\varDelta _{v_x}/W$ increases with $\theta$ and decreases with $W$. However, at a given $\theta$, $\varDelta _{v_x}$ increases with $W$ (not shown) confirming that, in agreement with Jop et al. (Reference Jop, Forterre and Pouliquen2005), the effect of the sidewalls is significant even at large gap width. The case $W=10d$ is intermediate between the two latter configurations. In any cases, the effect of the sidewalls cannot be neglected. The value of $\varDelta _{v_x}$ seems to be of the same order as $\varDelta _\nu$, the size of the boundary layer determined from the volume fraction profiles, in the flowing zone.

5.3. Rotation

The component of the rotation velocity normal to the bottom (i.e. $\omega _z$) is also of interest since it is potentially strongly influenced by the presence of the flat but frictional sidewalls (Yang & Huang Reference Yang and Huang2016). In figure 19, we have reported the transverse variations of $\omega _z$ at different flow depths. For the sake of clarity, we have reported in the inset of this figure the volume fraction and the number of contact profiles and indicate on them the latter depths. Like all the other quantities reported in this paper, this quantity has been space and time averaged (see the end of § 2). In the vicinity of the sidewalls and in the creeping zone, the rotation velocity $\omega _z$ is close to $2 v_s /d$. This indicates that, at these depths, the grains located in the vicinity of the sidewalls mainly roll without sliding at the contact with the sidewall. The influence of the sidewalls decreases in the flow zone, indicating that the propensity of the grains for sliding at contact increases. The correlations of the rotations along the transverse direction are complex. Counter-rotations (anti-correlated rotations) are observed close to the sidewalls (see figure 19). The number of involved layers increases with the depth within the flow. It should be pointed out that the observed behaviour is somewhat similar to that of a system with frustrated rotations (Radjai, Brendel & Roux Reference Radjai, Brendel and Roux1996; Khidas et al. Reference Khidas, Schliecker, Ammi, Messager and Delannay2000).

Figure 19. Profiles of vertical rotation velocity rescaled by $2 v_s/d$ vs $y$ for an angle $\theta =40^{\circ }$ and for several vertical positions. The gap between the sidewalls is $W/d=10$. For the sake of clarity, the inset reports the aforementioned vertical positions (dashed lines) on the vertical profiles of the number of contacts, $N_c$ and of the volume fraction, $\nu$.

5.4. Granular temperature

At the centre and at the sidewalls the temperature is found to increase from the creeping zone to the free surface, the differences between the two profiles being more important for larger gaps $W$. Note that for $W/d=20$ a decrease of granular temperature is observed at $z/d\approx -10$. This location corresponds to the kink observed on the volume fraction profile (see figure 16). As mentioned in § 5.1 it corresponds to the formation of clusters of grains. The transverse profile of the granular temperature is very informative. Outside the creeping zone, this quantity is greatest at the sidewalls and lowest in the centre. In contrast, in the creeping zone, the granular temperature gradually rises from its minimal value at the sidewalls to a maximum value at the centre of the cell. Note that similar results have been observed in a confined shear cell (Artoni et al. Reference Artoni, Soligo, Paul and Richard2018; Richard et al. Reference Richard, Artoni, Valance and Delannay2020). These results have important consequences. Depending on the vertical position, sidewalls could be either a granular heat source (in the flowing zone) or a sink (in the creeping zone). This points out the difficulty in stipulating a sidewall boundary condition on the granular temperature for theories aiming to capture the properties of granular flows involving both creeping and flowing zones. Note that, in the flowing zone, the observed high temperature induces an important pressure and a low volume fraction (see § 5.1).

A comparison between figures 20 and 17 suggests that the length of influence of the sidewalls on the granular temperature is weaker than that on the velocity. To illustrate this, we can compare the transverse profile of the temperature and the velocity for $W/d=20$ and $\theta =45^{\circ }$ (i.e. figures 20d and 17d, respectively). For both quantities, a plateau is observed at the centre of the cell, yet that observed for the temperature is significantly larger than that of the velocity. In this case, the length characterizing the effect of the sidewalls on the granular temperature, estimated by the distance between a sidewall and the beginning of the plateau, is approximately equal to the half of that measured on the streamwise velocity profile.

Figure 20. Transverse temperature profiles for $\theta =45^{\circ }$, $W/d=10$ (a) and $\theta =45^{\circ }$, $W/d=20$ (b) in the centre of the simulation cell and at the sidewalls. The corresponding transverse temperature profiles are respectively reported in (c,d). The absolute difference between the vertical profiles of temperature at the centre and at the sidewalls is weak in the case of small gap width (a), yet the transverse relative variations are not negligible (c). The vertical dashed lines reported in (a,b) correspond to $z/d=l$. In the case of large gap width, the transverse variations of the granular temperature are important, as shown by the vertical profiles at the sidewalls and at the centre (b) as well as on the transverse profiles (d). Interestingly, depending on the flow zone (flowing or creeping) the sidewalls behave like a source or a sink of granular temperature (c,d).