1 Introduction
Particles are often transported in the form of dense currents under the influence of gravity. Their bulk flow rate is greatly enhanced if part or all of their weight is supported by a gas flow through them. When particles are poured onto a slope that is less than their angle of repose, they are held stationary by the action of contact friction and merely flow down the surface of the pile in a thin layer as more grains are successively added. When they are poured onto a slope that is steeper than the angle of repose, a thin, dense current forms in which the particles move in bulk down the slope (see, for example, Ishida, Hatano & Shirai Reference Ishida, Hatano and Shirai1980; GDR MiDi 2004). If a gas is passed vertically through the particles then the drag it exerts on the particles bears some of their net weight and hence the frictional forces decrease. Consequently, the effective angle of repose of the particles decreases as does the minimum slope angle at which bulk flow takes place (Nott & Jackson Reference Nott and Jackson1992). When the gas flow is sufficiently large for the entire weight of the particles to be supported (i.e. the particles are fluidised), then bulk frictional forces are insignificant and very mobile currents form, even on horizontal surfaces. The influence of a fluidising gas flow through particles on their mobility is exploited widely in industrial settings where it is necessary to transport bulk materials either to move them from one place to another using air slides (which can be several kilometres long), or to keep horizontal surfaces clear of particles in pieces of processing equipment such as circulating fluidised beds (Savage & Oger Reference Savage and Oger2013). There are also features in many particulate environmental flows in which there is significant upward gas flow and this enhances their speed and range (e.g. Druitt Reference Druitt, Gilbert and Sparks1998; Roche et al. Reference Roche, Gilbertson, Phillips and Sparks2004).
There have been extensive studies of the flow of particles down a slope and of some of the effects of fluidisation. A common approach to mathematically modelling these motions is based on a continuum description that couples expressions of mass conservation with expressions of the balance of momentum within each phase (e.g. Nott & Jackson Reference Nott and Jackson1992). Under this approach, the fluidised material is treated as two inter-penetrating phases that interact with each other. The models do not resolve the motion of individual particles, but rather the evolution of averaged, bulk properties, which depend upon the net effect of direct interactions between particles within the current and between the particles and their surroundings. The duration of contacts between the constituent particles has important consequences: if the contacts are sustained then they are likely to be frictional in nature; if they are instantaneous then they are collisional in nature (e.g. Campbell Reference Campbell2006). The stresses induced by instantaneous collisions between pairs of particles (i.e. in dilute and rapid granular flows) can be evaluated through the use of granular kinetic theories (Jenkins & Savage Reference Jenkins and Savage1983), in which a key dependent variable is the granular temperature,
$T$
, a measure of the variance of the instantaneous velocity field. Hydrodynamic equations of motion have then been derived for granular materials that are much like those for dense gases except there is substantial energy dissipation through inelastic collisions (see Haff Reference Haff1983; Jenkins & Savage Reference Jenkins and Savage1983; Lun et al.
Reference Lun, Savage, Jeffrey and Chepurniy1984, for example). It is possible, of course, for there to be collisions over a range of durations and these ideas do not translate to dense and slowly shearing flows where contacts are prolonged and thus, in part, frictional. The action of the interstitial fluid is a further factor that needs to be considered when modelling granular flows. For example, in contrast to the original studies of granular kinetic theories, Koch & Sangani (Reference Koch and Sangani1999) proposed that interaction with the fluid could generate agitation within the flows and that fluctuating viscous forces could be the generators of particle temperature.
There have been relatively fewer studies that report granular flows that are aerated or fluidised. An early approach was to treat the fluidised particles as a non-Newtonian fluid of power-law rheology, sometimes with a yield stress (Botterill & Bessant Reference Botterill and Bessant1973; Botterill & Abdul-Halim Reference Botterill and Abdul-Halim1979; Ishida et al. Reference Ishida, Hatano and Shirai1980; Savage & Oger Reference Savage and Oger2013). This approach can be made to work well in specific practical situations (Singh, Callcott & Rigby Reference Singh, Callcott and Rigby1978), but is entirely reliant on empirical methods to determine the effective rheology in each circumstance since such approaches do not capture the fundamental dynamics of the particle motion.
A more fundamental approach is to model the evolution of averaged properties of the inter-penetrating phases. Ogawa, Umemura & Oshima (Reference Ogawa, Umemura and Oshima1980) modelled steady, one-dimensional fully fluidised currents down slopes. They derived constitutive relations based on the collisions between a particle and its neighbours, which were represented by an imaginary spherical shell surrounding it. This resulted in a balance between collisional stresses and gravitational forces. Nott & Jackson (Reference Nott and Jackson1992) coupled a kinetic theory for collisional grain flows with a Coulomb-like model for frictional effects to predict the bulk mass flow rate of aerated grains down an inclined channel. The experiments and model featured gas flow rates up to the minimum required to fully fluidise the particles. Their mathematical model of friction in the flows followed Johnson & Jackson (Reference Johnson and Jackson1987) and Johnson, Nott & Jackson (Reference Johnson, Nott and Jackson1990) and assumed that the frictional component (dominant at high particle volume fractions,
$\unicode[STIX]{x1D719}$
) was simply added to the collisional component (dominant at low
$\unicode[STIX]{x1D719}$
). They pursued a similar approach to the interaction term between the gas and the particles adding together a contribution based on the Ergun equation (dominant at high
$\unicode[STIX]{x1D719}$
) and one from the Richardson and Zaki equation (dominant at low
$\unicode[STIX]{x1D719}$
). No contribution was included from slip between the two phases in the direction of the slope. Oger & Savage (Reference Oger and Savage2013) took a similar approach (although with some different closures of the models), again retaining a frictional term, and solved the resulting equations using the MFIX numerical code to study the dynamics of granular motion within air slides, computing the steady, fully developed velocity and granular temperature fields for flows within a channel of rectangular cross-section. Finally, Eames & Gilbertson (Reference Eames and Gilbertson2000) reported the unsteady flow of fluidised materials along horizontal surfaces. For their system, they showed that collisional stresses would be small compared with those associated with fluid drag and so when fully fluidised, the force balance set hydrostatic pressure gradient against fluid drag terms. We will show below how our work differs from their modelling framework and yet is able to reproduce features of their experimental results.
Key to furthering our understanding of the dynamics of fluidised flows is direct and detailed experimental evidence against which theoretical models can be validated. However, there are few measurements of fluidised granular currents, especially down slopes. Previous experimental studies have presented bulk properties such as total flow depth and mass flow rate (e.g. Nott & Jackson Reference Nott and Jackson1992; Eames & Gilbertson Reference Eames and Gilbertson2000). Some measurements of local properties such as velocity have been made though these have often been with instruments such as optical probes or turbine elements (e.g. Botterill & Bessant Reference Botterill and Bessant1973, Reference Botterill and Bessant1976; Ishida et al. Reference Ishida, Hatano and Shirai1980; Nott & Jackson Reference Nott and Jackson1992). Whilst providing important information, the disadvantages to these techniques are that they lack spatial resolution, are intrusive (especially in fluidised particles Rowe & Masson (Reference Rowe and Masson1981) and offer only point measurements i.e. traverses are necessary to build velocity profiles and they are therefore only suited to steady flows. More recently, particle image velocimetry (PIV) has been applied to fluidised systems such as static beds (Bokkers, van Sint Annaland & Kuipers Reference Bokkers, van Sint Annaland and Kuipers2004) and dam-break experiments over horizontal surfaces of initially fluidised, fine natural volcanic ash (Girolami et al. Reference Girolami, Roche, Druitt and Corpetti2010). PIV has the advantage of offering high spatial resolution and allows instantaneous velocity fields to be calculated. The experiments of Girolami et al. (Reference Girolami, Roche, Druitt and Corpetti2010) had a short-lived phase of quasi-constant flow following the initial release of material; however, they also experienced compaction of the solid phase during the flow because the grains were not continuously fluidised along the apparatus. This meant that even though the materials were highly expanded initially because of the very small particle size, they decelerated rapidly due to the loss of mobility associated with compaction in the terminal flow phase. As such, they are not representative of fully fluidised flows.
The aim of the present work is to understand better the dynamics of fluidised granular flows by providing further experimental evidence and proposing a new unsteady model of these flows that fully takes into account the interaction between the particles and the fluid and incorporates collisional stresses. Both of these processes play a crucial role in the dynamics of fluidised granular flows in which the gas flow bears most of the weight of the particulate layer and the particle interaction contribute significantly to the shear stresses developed by the flow. This implies that the dynamics are different from ‘dry’ granular flows in which the role of the interstitial fluid is negligible (e.g. Lun et al. Reference Lun, Savage, Jeffrey and Chepurniy1984; Forterre & Pouliquen Reference Forterre and Pouliquen2008; Woodhouse, Hogg & Sellar Reference Woodhouse, Hogg and Sellar2010).
In this work, experimental measurements were made of granular currents over a range of slope inclinations and conditions and the experimental arrangement is described in § 2. The measurements were made in an apparatus that confined the flow between two walls, which enabled the overall size and shape of currents to be measured over time. In addition, PIV was used to measure the velocity profiles of the particles within the currents, enabling their overall behaviour to be linked to their rheology. Section 3 develops the general continuum model and the equations of motion for the flowing state. This builds upon the ‘two-fluid’ approach in which the gas and grains are treated as two inter-penetrating phases (Jackson Reference Jackson2000). Fully developed flows are tackled in § 4 and compared with experimental observations. The continuum model in this section is analysed in the regime for which the properties of the flowing layer vary only with distance from the underlying boundary and the solutions are computed numerically and asymptotically in a regime where the flow thickness far exceeds the diameter of an individual grain. Unsteady and transient effects found in flows along inclined channels are investigated in § 5 and a new model developed in the ‘lubrication’ regime where the downslope length scale is much large than that perpendicular to the slope. Flows along horizontal surfaces differ their counterparts along inclines (§ 6) and measurements of their inherently unsteady motion are reproduced well by a new self-similar solution to the flow model in the lubrication regime. Finally our findings are summarised and discussed in § 7. We also include two appendices. In the first we analyse the consequences of an extended kinetic theory, following the constitutive laws of Jenkins (Reference Jenkins2007). In the second the effects of the sidewalls are analysed in the regime that the flow depth is much less than the channel width.
2 Experimental approach
2.1 Experimental set-up
Figure 1. Schematic of flows and experimental set-up. Material is introduced from the raised end of the apparatus at a constant flux,
$Q$
. A flow of fluidising gas enters the apparatus at speed
$w_{g}$
through a porous distributor plate and is constant along the entire length of the apparatus and duration of an experiment. The apparatus can be included to some angle,
$\unicode[STIX]{x1D703}$
. The resulting flow has a height profile,
$h(x,t)$
, length (front position),
$x_{f}(t)$
, and longitudinal velocity profile,
$v(x,t)$
. For flows down inclined channels, the height of the current increases near to the front (head) to a constant value which is obtained towards the rear of the flow (body).
The experimental arrangement conformed to that shown in figure 1. The apparatus was a long, narrow channel (1 cm
$\times$
100 cm, 50 cm in height) which could be inclined to some angle,
$\unicode[STIX]{x1D703}$
, to the horizontal. The bottom of the channel was a porous plastic distributor material (Vyon ‘D’) through which dry air was passed from a windbox below at a speed
$w_{g}$
, but for which the pressure drop over it was much larger than that through the granular flow. This ensured that the gas flow was evenly distributed i.e. the gas flow entering the apparatus was uniform and perpendicular to the distributor plate so that at the base of the granular layer, the gas velocity is
$\boldsymbol{u}|_{z=0}=(0,w_{g}/(1-\unicode[STIX]{x1D719}(0)))$
where
$\unicode[STIX]{x1D719}(0)$
is the particle volume fraction evaluated at the base of the flow. The particles were constrained between vertical parallel walls, so that the motion is effectively two-dimensional and the motion of particles within the current could be seen. The front wall was made from a glass sheet allowing the flows to be viewed and the other sides were made of aluminium plate. The rear plate was painted black to increase the contrast between it and the white particles. Particles entered the apparatus at one end (the uppermost when inclined) through a funnel giving a constant volumetric flux,
$Q$
, which could be changed between experiments by changing the aperture of the funnel (Nedderman Reference Nedderman1992). The flow rate
$Q$
was the flow rate of the current based on the bulk volume when the particles were at rest; so, the volume flow rate per unit width of particles
$q_{0}=\unicode[STIX]{x1D719}_{m}Q/B$
where
$\unicode[STIX]{x1D719}_{m}$
is the particle volume fraction of a static bed of particles, and the distance between the front and back of the flows is
$B=1~\text{cm}$
. It could be controlled by using funnels of different sizes, each of which could then be associated with a bulk flow rate,
$Q_{nom}$
; however, this is a nominal flow rate as the actual flow rate could vary from occasion to occasion. The apparatus had closed ends; so, to avoid ‘backing up’ when running experiments with a non-horizontal slope, particles were removed from the downslope end using a vacuum cleaner. This had no measurable effect on the height profiles obtained but allowed experiments to be run for longer. No removal of particles was necessary for the slower-moving horizontal flows. The value of
$q_{0}$
was accurate to within
$\pm 3\,\%$
.
The material used for all the experiments was approximately spherical, glass beads (Potters Ballotini) with particle diameters in the range
$250{-}425~\unicode[STIX]{x03BC}\text{m}$
and a mean diameter
$d\approx 375~\unicode[STIX]{x03BC}\text{m}$
. We measured the particle volume fraction of densely packed, static material (i.e. a maximum) as
$\unicode[STIX]{x1D719}_{m}=0.610\pm 0.005$
, which is close to the maximum value of
$0.64$
for random, close-packed, mono-sized spheres (Jaeger & Nagel Reference Jaeger and Nagel1992), and the (unfluidised) bulk density was
$1.520\pm 0.008~\text{g}~\text{cm}^{-3}$
. The powder corresponds to a class B powder according to Geldart (Reference Geldart1973), so no bubble-free expansion when fluidised is expected within a ‘static’ bed. The minimum fluidisation speed,
$u_{mf}$
, was found by independent experiment where the gas flow rate through a static bed of material was gradually increased and the resulting pressure drop through the bed measured (Davidson & Harrison Reference Davidson and Harrison1963). We found that the entire weight of the bed was supported when
$w_{g}(=u_{mf})=10.77~\text{cm}~\text{s}^{-1}$
.
2.2 Shape of currents and front position extraction
The flows, viewed from the side, were recorded using a digital video camera. Calibration was performed using an image of a block of known dimensions placed in the apparatus once the camera was set-up in position for a given experiment. Still images from the recorded experiments were analysed by transforming the RGB images to grey scale. These were then turned into binary images through thresholding. Though the threshold value was calculated automatically, the contrast between black back wall and white particles meant that the resulting binary images were robust and consistent. The upper and lower surfaces of the outline of a current were defined as the first and last white pixels when descending a column of the binary image. 95 % confidence intervals for height measurements are
$\pm 0.1~\text{cm}$
. The front position was taken as the point where the top surface met the bottom one.
2.3 Velocity measurements
PIV was used to make measurements of the velocity fields of the flows using a high-speed video camera capturing at 500 frames per second (f.p.s.) close up to a particular region of the flow. The PIV measurements required the flow to be seeded with marker particles for which we used the same-sized particles as for our other experiments but approximately one third of which were dyed black. The properties of the dyed particles (
$u_{mf}$
, angle of repose etc.) were identical to the non-dyed particles.
Two-dimensional velocity fields were calculated by processing image pairs (two consecutive frames) from the video taken by the high-speed camera using the open-source Matlab-based DPIVSoft2010 code. The software makes an initial estimation of the velocity field on a coarse grid and then uses this to translate and deform the interrogation window in the second image in keeping with the deformation of the flow field. Errors associated with image pattern distortion, as is the case when velocity gradients are large, are greatly decreased using this method (Meunier & Leweke Reference Meunier and Leweke2003). Several initial iterations were run to get a good approximation for the flow field. A final run was performed with an interrogation window of 32
$\times$
32 pixels (
${\approx}5d$
) and velocity vectors were calculated using a 50 % overlap between adjacent windows. A median filter was then applied with a limit of 0.5 to remove spurious vectors (e.g. Adrian & Westerweel Reference Adrian and Westerweel2010, p. 406).
Instantaneous velocity profiles may not be representative of the flow as a whole. In particular, the bubbles of gas that could form spontaneously in the flows often disrupted the instantaneous velocity profiles. However, flows down the steeper slopes in our experiments,
$10^{\circ }$
and
$15^{\circ }$
, reached a steady state very quickly, and for these flows an ensemble average of the flow velocity could be found by averaging over both many points in time and at several positions along the flow. The quality, and hence the accuracy, of time-averaged velocity fields has been shown to be greatly improved when the average instantaneous correlation function is used to calculate the velocity field (Meinhart, Wereley & Santiago Reference Meinhart, Wereley and Santiago2000). We therefore modified the PIV routines accordingly to produce a single time-averaged velocity field per experiment using an interval of twenty frames (0.04 s) between image pairs, and fifteen image pairs per experiment. This interval is larger than characteristic time for shear (
$(\text{d}v/\text{d}z)^{-1}\approx 0.01~\text{s}$
), so the velocity fields at successive intervals are uncorrelated. Velocity profiles were then formed from the streamwise vectors of the time-averaged velocity field lying on a depthwise transect at points separated at intervals of 1 cm and averaged to form the ensemble average velocity profile. The resolution of PIV measurements can be expressed as (Adrian Reference Adrian1991),
$$\begin{eqnarray}\unicode[STIX]{x1D70E}_{u}\approx \frac{c_{1}Md}{\unicode[STIX]{x0394}t},\end{eqnarray}$$
where
$c_{1}$
is the uncertainty of locating the centroid of the correlation peaks,
$M$
is the magnification factor of the lens and
$\unicode[STIX]{x0394}t=1/500$
s is the time step between images. For our set-up,
$M=1/2$
, and
$c_{1}\approx 10\,\%$
so that
$\unicode[STIX]{x1D70E}_{u}=O(1)~\text{cm}~\text{s}^{-1}$
.
For the steady flows it is more useful to define error based on the sum of variances of the all the
$m$
profiles used to calculate the ensemble-averaged standard deviation over the
$n$
images given by
$$\begin{eqnarray}\unicode[STIX]{x1D70E}_{ens\,ave}^{2}(z)={\displaystyle \frac{1}{m\times n}}\mathop{\sum }_{i=1}^{m\times n}\unicode[STIX]{x1D70E}_{i}^{2}(z).\end{eqnarray}$$
This average standard deviation was then used to calculate 95 % confidence intervals for the velocities.
3 Equations of motion
We investigate the motion of granular currents down an inclined surface when the particles are fluidised, as shown in figure 1. These flows are gravitationally driven, but do not accelerate unboundedly; instead the principle action of particle interactions is to contribute to the shear stresses that balance the down slope acceleration and potentially lead to steady motion. We formulate a mathematical model of the two-phase motion that couples mass conservation for each phase with expressions for the balances of momentum and we show how this formulation may be applied to steady fully developed flows that vary only with distance from the underlying boundary (§ 4), and to unsteady, relatively thin flows for which the acceleration perpendicular to the underlying boundary is negligible (§§ 5 and 6).
The mathematical model is built upon a continuum description of two inter-penetrating phases which interact with each other. These models do not resolve the motion of individual particles; rather, they allow the computation of the evolution of averaged properties. Such approaches have been employed often for confined, horizontal fluidised beds (e.g. Bokkers et al. Reference Bokkers, van Sint Annaland and Kuipers2004; Goldschmidt, Beetstra & Kuipers Reference Goldschmidt, Beetstra and Kuipers2004), but these studies differ from the dynamics of the flows analysed in this contribution where there is persistent shear flow down the inclined surface. The flows analysed here also differ in an essential way from non-fluidised granular motion down inclines since the support of the weight of the grains by the imposed gas flow significantly reduces resistive forces and increases mobility. Nevertheless, we find that steady flows are admissible and thus the motion must develop sufficient shear stresses to balance gravitational acceleration. Our model assumes that these stresses arise from particle interactions and are collisional and the particle fluctuations may be characterised by a granular temperature since friction as a bulk property is virtually eliminated by fluidisation and the viscous forces associated with interstitial gas flow are negligible. The granular temperature will be shown to be relatively small and thus the interactions generate only relatively weak shear stresses, but these are sufficient to balance the gravitational acceleration.
The collisional nature of the motion is justifiable in all but some small regions of the currents, for example close to the surface of the slope. The model captures only the relatively slow evolution of averaged quantities. In particular bubbles (i.e. volumes largely evacuated of particles that travel through fluidised particles) are not explicitly resolved. Bubbles are an important feature of deep, static fluidised beds as apart from strongly affecting the local instantaneous volume fraction of particles and they are the primary source of granular temperature in such a bed (Menon & Durian Reference Menon and Durian1997). There are several processes that might lead to the generation or suppression of bubbling, most notably including the dissipation of granular temperature through collisions, which is prone to clustering instabilities (Goldhirsch & Zanetti Reference Goldhirsch and Zanetti1993; Fullmer & Hrenya Reference Fullmer and Hrenya2017). Some studies have sought to predict the onset of bubbling in static beds through linear stability analysis (e.g. see the review by Jackson Reference Jackson2000). The flows of fluidised materials analysed here are somewhat different from these stability analyses, however, due to the persistent production of granular temperature by work done by the velocity field shear against the shear stresses, a process absent in static beds; hence, by means of a scaling analysis Eames & Gilbertson (Reference Eames and Gilbertson2000) showed that the contribution of these bubbles to the overall balance for granular temperature is likely to be negligible for this downslope motion. Furthermore, shallowness in the bed is thought to suppress bubbling (Botterill et al. Reference Botterill, van der Kolk, Elliott and Mcguigan1972; Tsimring, Ramaswamy & Sherman Reference Tsimring, Ramaswamy and Sherman1999), as is shear (Botterill & Abdul-Halim Reference Botterill and Abdul-Halim1979; Ishida et al. Reference Ishida, Hatano and Shirai1980). We therefore assume that bubbling is likely to have a limited influence on the fluidised currents. Extensive bubbling was not observed in the currents. The photograph shown later in figure 9 is typical with no apparent bubbles. While agitation was visible at the top of the currents, bubbles sufficiently large to fill the width of the bed were hardly ever seen.
Most of the theoretical developments in this study will be for two-dimensional flows and the effects of the front and back walls of the apparatus are neglected. The use of this planar set-up allows the structure of the system to be seen and measured (as described in § 2), but at the expense of it being bounded by walls not present in realistic, three-dimensional systems. Arguably, because fluidisation eliminates internal friction, a large part of the effect that the presence of these walls might also be eliminated. Here, in most of what follows, we analyse the motion in the regimes that the sidewalls play a negligible role; however in § 6.1, we also analyse the case when the sidewalls have a dominant effect on horizontal currents and in appendix B, we derive the extra, weak retardation on flows down slopes that arises from sidewall drag when the depth of the flow is much smaller than the width.
The general equations of motion for a continuum model, known as a ‘two-fluid model’, of a gas–particle system have been developed by Jackson (Reference Jackson2000). The conservation of mass in each phase is given by
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}t}(1-\unicode[STIX]{x1D719})+\unicode[STIX]{x1D735}\boldsymbol{\cdot }((1-\unicode[STIX]{x1D719})\boldsymbol{u})=0 & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D719}}{\unicode[STIX]{x2202}t}+\unicode[STIX]{x1D735}\boldsymbol{\cdot }(\unicode[STIX]{x1D719}\boldsymbol{v})=0, & \displaystyle\end{eqnarray}$$
where
$\unicode[STIX]{x1D719}$
denotes the volume fraction of solids and
$\boldsymbol{u}$
and
$\boldsymbol{v}$
the velocity field of the gas and solid phase, respectively.
Following Jackson (Reference Jackson2000), the balance of momentum for the gas is given by
$$\begin{eqnarray}(1-\unicode[STIX]{x1D719})\unicode[STIX]{x1D70C}_{g}\frac{\text{D}_{g}\boldsymbol{u}}{\text{D}t}=(1-\unicode[STIX]{x1D719})\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D64E}^{g}-\boldsymbol{F}_{D}+(1-\unicode[STIX]{x1D719})\unicode[STIX]{x1D70C}_{g}\boldsymbol{g}\end{eqnarray}$$
and for the particles
$$\begin{eqnarray}\unicode[STIX]{x1D719}\unicode[STIX]{x1D70C}_{s}\frac{\text{D}_{s}\boldsymbol{v}}{\text{D}t}=\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D64E}^{s}+\unicode[STIX]{x1D719}\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D64E}^{g}+\boldsymbol{F}_{D}+\unicode[STIX]{x1D719}\unicode[STIX]{x1D70C}_{s}\boldsymbol{g},\end{eqnarray}$$
where
$\unicode[STIX]{x1D70C}_{g}$
and
$\unicode[STIX]{x1D70C}_{s}$
are the densities of the gas and solid phase respectively,
$\unicode[STIX]{x1D64E}$
is the spatially averaged stress tensor of each phase with the superscript (
$g$
,
$s$
) denoting the gas or solid phase, respectively,
$\boldsymbol{F}_{D}$
is the drag force exerted by the particles on the fluid due to the difference in their velocities and
$\boldsymbol{g}$
denotes gravitational acceleration. The material derivatives,
$\text{D}_{g}/\text{D}t$
and
$\text{D}_{s}/\text{D}t$
denote the rate of change moving with the gas and the solid phase respectively. A number of researchers, including Ergun (Reference Ergun1952) and Jackson (Reference Jackson2000), have suggested that
$\boldsymbol{F}_{D}=\unicode[STIX]{x1D6FD}(\boldsymbol{u}-\boldsymbol{v})$
, where
$\unicode[STIX]{x1D6FD}$
is a drag coefficient. Virtual mass and particle shear forces are neglected.
These equations will be solved for the situation shown schematically in figure 1. The slope is inclined at angle,
$\unicode[STIX]{x1D703}$
, to the horizontal with the underlying boundary at
$z=0$
and the upper surface of the current at
$z=h$
, while the
$x$
-axis is aligned with the basal boundary. A mixture of solid particles and gas runs down the slope under the influence of gravity.
4 Fully developed flows
4.1 Model for fully developed flows
First, fully developed flows are investigated, in which the dependent variables are functions only of the distance from the boundary,
$z$
, and the velocity fields of the gas and solids are given by
$\boldsymbol{u}=(u(z),0,w(z))$
and
$\boldsymbol{v}=(v(z),0,0)$
, respectively. Conservation of mass for the solid phase is automatically satisfied by this form, but for the fluid phase we deduce that
$$\begin{eqnarray}(1-\unicode[STIX]{x1D719})w=w_{g},\end{eqnarray}$$
where
$w_{g}$
is the fluidising gas flux per unit area normal to the boundary.
The expressions for the balance of momentum follow those proposed by Johnson & Jackson (Reference Johnson and Jackson1987) and Agrawal et al. (Reference Agrawal, Loezos, Syamlal and Sundaresan2001) where for the gas phase down the slope
$$\begin{eqnarray}\unicode[STIX]{x1D70C}_{g}w_{g}\frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}z}=(1-\unicode[STIX]{x1D719})\unicode[STIX]{x1D70C}_{g}g\sin \unicode[STIX]{x1D703}+(1-\unicode[STIX]{x1D719})\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}z}\left(\unicode[STIX]{x1D707}_{g}\frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}z}\right)+\unicode[STIX]{x1D6FD}(v-u),\end{eqnarray}$$
where
$g=|\boldsymbol{g}|$
denotes gravitational acceleration and
$\unicode[STIX]{x1D707}_{g}$
is the gas viscosity. Perpendicular to the slope, we find
$$\begin{eqnarray}\unicode[STIX]{x1D70C}_{g}w_{g}\frac{\unicode[STIX]{x2202}w}{\unicode[STIX]{x2202}z}=-(1-\unicode[STIX]{x1D719})\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}z}-(1-\unicode[STIX]{x1D719})\unicode[STIX]{x1D70C}_{g}g\cos \unicode[STIX]{x1D703}-\unicode[STIX]{x1D6FD}w+(1-\unicode[STIX]{x1D719})\frac{4}{3}\unicode[STIX]{x1D707}_{g}\frac{\unicode[STIX]{x2202}^{2}w}{\unicode[STIX]{x2202}z^{2}},\end{eqnarray}$$
where
$p$
is the pressure within the fluid phase. In (4.2) and (4.3) we have assumed that the gas phase is incompressible and can be treated as Newtonian with constant viscosity,
$\unicode[STIX]{x1D707}_{g}$
.
For the solid phase, the balance of down slope momentum is given by
$$\begin{eqnarray}0=\unicode[STIX]{x1D719}\unicode[STIX]{x1D70C}_{s}g\sin \unicode[STIX]{x1D703}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{xz}}{\unicode[STIX]{x2202}z}-\unicode[STIX]{x1D6FD}(v-u)+\unicode[STIX]{x1D719}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}z}\left(\unicode[STIX]{x1D707}_{g}\frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}z}\right),\end{eqnarray}$$
while normal to the slope,
$$\begin{eqnarray}0=-\unicode[STIX]{x1D719}\unicode[STIX]{x1D70C}_{s}g\cos \unicode[STIX]{x1D703}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{zz}}{\unicode[STIX]{x2202}z}-\unicode[STIX]{x1D719}\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}z}+\unicode[STIX]{x1D6FD}w+\unicode[STIX]{x1D719}\frac{4}{3}\unicode[STIX]{x1D707}_{g}\frac{\unicode[STIX]{x2202}^{2}w}{\unicode[STIX]{x2202}z^{2}}.\end{eqnarray}$$
In (4.4) and (4.5),
$\unicode[STIX]{x1D70E}_{xz}$
and
$\unicode[STIX]{x1D70E}_{zz}$
are components of the solid-phase stress tensor,
$\unicode[STIX]{x1D64E}^{s}$
, and at this stage we have not yet invoked any constitutive model for these stresses. Further, from (4.4) the driving force for the current is gravity and within this framework, currents over horizontal surfaces are inherently unsteady as they decelerate.
The downslope balance of momentum (4.4) differs from previous contributions. Nott & Jackson (Reference Nott and Jackson1992) implicitly assumed that there was no relative component of velocity downslope between each of the phases and thus there was no drag force (i.e.
$\unicode[STIX]{x1D6FD}(u-v)=0$
). Eames & Gilbertson (Reference Eames and Gilbertson2000) did not consider momentum balance for the fluid phase and imposed
$u=0$
; thus within their model, the drag force
$\unicode[STIX]{x1D6FD}v$
is dominant and by assumption the shear stress associated with the solid phase is negligible. We do not invoke either of these assumptions at this stage, instead maintaining the various dynamical processes until their relative magnitudes have been fully assessed below.
Adding the normal momentum equations (4.3) and (4.5), we find that
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}z}(p-\unicode[STIX]{x1D70E}_{zz})=-(\unicode[STIX]{x1D70C}_{g}(1-\unicode[STIX]{x1D719})+\unicode[STIX]{x1D70C}_{s}\unicode[STIX]{x1D719})g\cos \unicode[STIX]{x1D703}+\frac{4}{3}\unicode[STIX]{x1D707}_{g}\frac{\unicode[STIX]{x2202}^{2}w}{\unicode[STIX]{x2202}z^{2}}-\unicode[STIX]{x1D70C}_{g}w_{g}\frac{\unicode[STIX]{x2202}w}{\unicode[STIX]{x2202}z}.\end{eqnarray}$$
When the current is homogeneous so that particle volume fraction
$\unicode[STIX]{x1D719}$
is constant, from (4.1) the vertical component of the gas velocity is also constant; so, (4.6) expresses the hydrostatic balance between the vertical gradient of the normal stress from both solid and fluid phases and the weight of the fluidised grains.
It is also insightful to eliminate the fluid pressure field between (4.3) and (4.5) to find that
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{zz}}{\unicode[STIX]{x2202}z}+\frac{\unicode[STIX]{x1D719}}{(1-\unicode[STIX]{x1D719})}\unicode[STIX]{x1D70C}_{g}w_{g}\frac{\unicode[STIX]{x2202}w}{\unicode[STIX]{x2202}z}=\unicode[STIX]{x1D719}(\unicode[STIX]{x1D70C}_{s}-\unicode[STIX]{x1D70C}_{g})g\cos \unicode[STIX]{x1D703}-\frac{\unicode[STIX]{x1D6FD}w}{1-\unicode[STIX]{x1D719}}.\end{eqnarray}$$
This expression reveals the fundamental dynamical role played by fluidisation. The slope-normal component of the inter-phase drag, incorporated into the model by
$\unicode[STIX]{x1D6FD}w$
, can balance the weight of the grains and thus it is possible for the normal stress tensor of the solid phase
$\unicode[STIX]{x1D70E}_{zz}$
to be much reduced from its non-fluidised magnitude. This in turn reduces the magnitude of the solids shear stress
$\unicode[STIX]{x1D70E}_{xz}$
and thus the mobility of the fluidised flows is greatly enhanced. Equation (4.7) is different from the classical model of a static fluidised bed because velocity gradients lead to normal stresses in the solid phases and these may contribute in a non-negligible way to the balance between weight and drag as shown below.
4.1.1 Inter-phase drag and constitutive equations
The drag on the solid phase due to the fluidising gas flow is given by
$\unicode[STIX]{x1D6FD}(\boldsymbol{u}-\boldsymbol{v})$
, where the drag coefficient,
$\unicode[STIX]{x1D6FD}$
, may be written
$$\begin{eqnarray}\unicode[STIX]{x1D6FD}=\frac{\unicode[STIX]{x1D707}_{g}}{d^{2}}f_{0}(\unicode[STIX]{x1D719})+\frac{\unicode[STIX]{x1D70C}_{f}}{d}|\boldsymbol{u}-\boldsymbol{v}|f_{0}^{\ast }(\unicode[STIX]{x1D719}),\end{eqnarray}$$
where
$f_{0}$
and
$f_{0}^{\ast }$
are given in table 1 (Ergun Reference Ergun1952). The first term on the right-hand side of the equation represents the drag associated with viscous processes, and the second term with inertial processes. For the regime of interest in this study, the inertial effects are negligible since the Reynolds number, based on gas velocity and particle size, is sufficiently small
$(Re\equiv \unicode[STIX]{x1D70C}_{g}w_{g}d/\unicode[STIX]{x1D707}_{g}<10)$
; however for completeness at this stage we maintain it in the model formulation. Other expressions for the drag coefficient have been used (e.g. Agrawal et al.
Reference Agrawal, Loezos, Syamlal and Sundaresan2001; Oger & Savage Reference Oger and Savage2013) and these could replace (4.8) within this modelling framework.
Table 1. The constitutive laws for granular kinetic theory applied to fluidised systems. Here
$e$
denotes the coefficient of restitution characterising collisions between particles;
$e_{w}$
is the coefficient of restitution between the walls and the particles;
$\unicode[STIX]{x1D719}_{m}$
is the volume fraction at maximum packing;
$\unicode[STIX]{x1D713}$
is the specularity coefficient (after Johnson & Jackson Reference Johnson and Jackson1987); and
$g_{0}$
is the radial basis function which accounts for particle packing (4.15). The coefficient
$c^{\ast }=32(1-e)(1-2e^{2})/(81-17e+30e^{2}(1-e))$
(Garzo & Dufty Reference Garzo and Dufty1999).
$f_{0}$
contributes to the inter-phase drag in the viscous regime and
$f_{0}^{\ast }$
in the inertial regime (see, for example, van der Hoef, Beetstra & Kuipers Reference van der Hoef, Beetstra and Kuipers2005).
$f_{1}$
,
$f_{2}$
model the volume fraction dependence in the collisional contributions to stresses (Garzo & Dufty Reference Garzo and Dufty1999);
$f_{3}$
,
$f_{4}$
and
$f_{4}^{\ast }$
model contributions to the granular energy balance (Garzo & Dufty Reference Garzo and Dufty1999).
$f_{5}$
to
$f_{7}$
determine boundary conditions at the base of the flow:
$f_{5}$
contributes to the boundary condition for momentum balances and
$f_{6}$
and
$f_{7}$
to that for fluctuation energy (see Johnson & Jackson Reference Johnson and Jackson1987; Johnson et al.
Reference Johnson, Nott and Jackson1990).
It was noted above that particle interactions are dynamically important because of the momentum transfer arising from particle collisions (Lun et al.
Reference Lun, Savage, Jeffrey and Chepurniy1984; Garzo & Dufty Reference Garzo and Dufty1999). Here we examine the collisional stresses and follow Nott & Jackson (Reference Nott and Jackson1992) and Agrawal et al. (Reference Agrawal, Loezos, Syamlal and Sundaresan2001) amongst others who incorporate these effects into models of fluidised and aerated flows, to write the shear and normal components of stress in terms of a granular temperature,
$T$
, which measures the fluctuations of velocity about the mean, the volume fraction of solids,
$\unicode[STIX]{x1D719}$
, and the coefficient of restitution
$e$
, which characterises dissipation in the instantaneous collisions. While the constitutive laws invoked here have been validated in some scenarios by simulation and experimentation, there remains some uncertainty about their generality. Hence we pose the model quite generally so that the constitutive relations could be updated as required.
For fully developed flows, we write
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70E}_{xz}=f_{1}(\unicode[STIX]{x1D719},e)\unicode[STIX]{x1D70C}_{s}\,\text{d}T^{1/2}\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70E}_{zz}=-f_{2}(\unicode[STIX]{x1D719},e)\unicode[STIX]{x1D70C}_{s}T, & \displaystyle\end{eqnarray}$$
where
$f_{1}$
and
$f_{2}$
are dimensionless functions given in table 1. In this study, we employ the constitutive formulae derived by Garzo & Dufty (Reference Garzo and Dufty1999) and recently used for modelling dense avalanches by Jenkins & Berzi (Reference Jenkins and Berzi2010).
Granular temperature may be generated and ‘conducted’ via the flow processes and dissipated in the collisions. Following Lun et al. (Reference Lun, Savage, Jeffrey and Chepurniy1984) and Garzo & Dufty (Reference Garzo and Dufty1999) amongst others, these effects are encompassed in the following expression of energy balance within the flow
$$\begin{eqnarray}0=-\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F9}}{\unicode[STIX]{x2202}z}+\unicode[STIX]{x1D70E}_{xz}\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}-f_{3}(\unicode[STIX]{x1D719},e)\unicode[STIX]{x1D70C}_{s}\frac{T^{3/2}}{d},\end{eqnarray}$$
where the flux of granular temperature is given by
$$\begin{eqnarray}\unicode[STIX]{x1D6F9}=-f_{4}\unicode[STIX]{x1D70C}_{s}\text{d}T^{1/2}\frac{\unicode[STIX]{x2202}T}{\unicode[STIX]{x2202}z}-f_{4}^{\ast }\unicode[STIX]{x1D70C}_{s}T^{3/2}d\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D719}}{\unicode[STIX]{x2202}z},\end{eqnarray}$$
and
$f_{3}$
,
$f_{4}$
and
$f_{4}^{\ast }$
are dimensionless functions given also in table 1. In posing this balance of granular temperature we have neglected generation and dissipation of the granular temperature mediated by viscous interactions with gas (Koch & Sangani Reference Koch and Sangani1999).
Following the formulation of Agrawal et al. (Reference Agrawal, Loezos, Syamlal and Sundaresan2001), dissipation by viscous processes is much smaller than dissipation through inelastic collisions when
$$\begin{eqnarray}\frac{\unicode[STIX]{x1D707}_{g}T}{d^{2}}\ll \frac{\unicode[STIX]{x1D70C}_{s}(1-e^{2})T^{3/2}}{d}.\end{eqnarray}$$
Furthermore, the generation of granular temperature by viscous processes is much smaller than that by granular interactions when
$$\begin{eqnarray}\frac{\unicode[STIX]{x1D707}_{g}^{2}|u-v|^{2}}{\unicode[STIX]{x1D70C}_{s}d^{3}T^{1/2}}\ll \unicode[STIX]{x1D70C}_{s}\text{d}T^{1/2}\left(\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}\right)^{2}.\end{eqnarray}$$
The constitutive laws,
$f_{1}-f_{4}^{\ast }$
, as well as those involved the boundary conditions (
$f_{5}-f_{7}$
, see § 4.1.3) feature the radial basis function,
$g_{0}(\unicode[STIX]{x1D719})$
. Various authors have suggested forms for
$g_{0}$
and we employ an expression that is close to the suggestion of Vescovi et al. (Reference Vescovi, Berzi, Richard and Brodu2014), who empirically fitted a function to match data from discrete element simulations. Importantly, the radial basis function diverges as the volume fraction approaches maximum packing (as established by Torquato (Reference Torquato1995)) and following Vescovi et al. (Reference Vescovi, Berzi, Richard and Brodu2014) we write
$$\begin{eqnarray}g_{0}={\hat{g}}\frac{2-\unicode[STIX]{x1D719}}{2(1-\unicode[STIX]{x1D719})^{3}}+(1-{\hat{g}})\frac{2}{\unicode[STIX]{x1D719}_{m}-\unicode[STIX]{x1D719}},\end{eqnarray}$$
where the weighting function is given by
$$\begin{eqnarray}{\hat{g}}=\left\{\begin{array}{@{}ll@{}}1 & \unicode[STIX]{x1D719}<\unicode[STIX]{x1D719}_{\ast },\\ \displaystyle 1-\left(\frac{\unicode[STIX]{x1D719}-\unicode[STIX]{x1D719}_{\ast }}{\unicode[STIX]{x1D719}_{m}-\unicode[STIX]{x1D719}_{\ast }}\right)^{n} & \unicode[STIX]{x1D719}_{\ast }<\unicode[STIX]{x1D719}<\unicode[STIX]{x1D719}_{m}.\end{array}\right.\end{eqnarray}$$
Thus, when
$\unicode[STIX]{x1D719}<\unicode[STIX]{x1D719}_{\ast }$
the radial basis function is given by the formula proposed by Carnahan & Starling (Reference Carnahan and Starling1969), but it exceeds this value when
$\unicode[STIX]{x1D719}_{\ast }<\unicode[STIX]{x1D719}$
and diverges as maximum packing is approached. Vescovi et al. (Reference Vescovi, Berzi, Richard and Brodu2014) suggest that
$n=2$
and that
$\unicode[STIX]{x1D719}_{\ast }=0.4$
. While this choice ensures that
$g_{0}$
and its derivative are continuous at
$\unicode[STIX]{x1D719}=\unicode[STIX]{x1D719}_{\ast }$
, the second derivative is discontinuous. This is problematic for the system of differential equations that we will integrate numerically; therefore, we employ the values
$n=3$
and
$\unicode[STIX]{x1D719}_{\ast }=0.4$
, which ensure that
$g_{0}$
is sufficiently smooth. Moreover, our expression (4.16) with these values is close to those proposed by Torquato (Reference Torquato1995) and Vescovi et al. (Reference Vescovi, Berzi, Richard and Brodu2014) and appears to match the simulation data adequately.
The energetic balance encompassed in (4.11) assumes that the particles are sufficiently agitated so that the particle diameter is the appropriate correlation length scale over which dissipated occurs (and the rate of dissipation is then given by
$\unicode[STIX]{x1D70C}_{s}f_{3}T^{3/2}/d$
). Recently, however, Jenkins (Reference Jenkins2007) has suggested that at relatively high concentrations, clusters of particles begin to form and thus the correlation length increases to
$L_{c}$
(
${>}d$
) and then the rate of dissipation is given by
$\unicode[STIX]{x1D70C}_{s}f_{3}T^{3/2}/L_{c}$
. This extended kinetic theory has been applied to unfluidised flows of grains down inclined planes by Jenkins & Berzi (Reference Jenkins and Berzi2010, Reference Jenkins and Berzi2012), where an empirical formula for
$L_{c}$
, informed by comparison with simulations and experimental measurements, is proposed in terms of the dependent flow variables. In our study there is potentially the need to include this phenomenon into the modelling framework to obtain good comparison between the predicted and measured results. However, as shown in appendix A, we find that extended kinetic theory makes negligible difference to the model predictions for the fluidised flows in our regime of interest and so we do not include it in the calculations that follow.
4.1.2 Coefficient of restitution
The dynamical effects of collisions between particles and between particles and the underlying boundary are characterised in the model by three parameters: the coefficients of restitution between the particles,
$e$
, and between the particles and the boundary,
$e_{w}$
, and the specularity coefficient
$\unicode[STIX]{x1D713}$
, which governs the dynamic interaction between the particles and the bottom surface (4.32). These parameters are relatively difficult to measure directly.
The coefficient of restitution,
$e$
, plays an important role in continuum models and in discrete particle (or element) models (DPMs), which endeavour to calculate the motion of large ensembles of particles and to resolve individual particle collisions. In continuum models, the difference of
$e$
from unity is proportional to the rate at which the collisions dissipate energy (see the definition of
$f_{3}$
in table 1), whereas in DPMs it controls the ratio of normal velocities before and after binary collisions and in these models, there are potentially additional means of energy dissipation. Often values for the coefficient of restitution are adopted without independent experimental confirmation and for DPM studies, typical values are relatively high (for example,
$e=0.90$
and
$0.97$
respectively in the studies of Goldschmidt et al. (Reference Goldschmidt, Beetstra and Kuipers2004) and van der Hoef et al. (Reference van der Hoef, van Sint Annaland, Deen and Kuipers2008)). These values are close to measured values of discrete collisions (see, for example, Kharaz, Gorham & Salman Reference Kharaz, Gorham and Salman2001). When used in kinetic theory models, commonly adopted values of
$e$
are rather lower and Jenkins & Zhang (Reference Jenkins and Zhang2002) suggest a means by which the appropriate value for kinetic theories can be derived from directly measured normal and tangential coefficients of restitution and the tangential coefficient of friction. For glass spheres of 3 mm diameter, the measured data of Foerster et al. (Reference Foerster, Louge, Chang and Allia1994) correspond to an effective coefficient of
$e=0.85$
if the method of Jenkins & Zhang (Reference Jenkins and Zhang2002) is employed and this is the value we employ in this study. We have no direct measurements of the appropriate coefficient of restitution for the collisions between the particles and underlying boundary; we choose
$e_{w}=0.75$
, but note that its magnitude has very little influence upon the computed flow profiles apart from within thin basal boundary layers.
4.1.3 Boundary conditions
The boundary conditions for this problem follow the formulation of Johnson & Jackson (Reference Johnson and Jackson1987) and Johnson et al. (Reference Johnson, Nott and Jackson1990). At the base, there is no slip for the fluid phase, a slip condition for the particle phase and a condition specifying the flux of granular temperature. These are respectively given by
$$\begin{eqnarray}u=0,\quad f_{5}d\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}=v\quad \text{and}\quad \unicode[STIX]{x1D6F9}=\unicode[STIX]{x1D70C}_{s}T^{1/2}(f_{6}v^{2}-f_{7}T)\quad \text{at }z=0.\end{eqnarray}$$
What this means physically is that solid-phase slip is allowed and stress is transmitted in the downslope direction by specularity i.e. the degree to which the angle of exit of a particle after collision with the base is different from the entry angle. This is mathematically represented by the specularity coefficient
$\unicode[STIX]{x1D713}$
$(0<\unicode[STIX]{x1D713}<1)$
. Furthermore, fluctuation energy (granular temperature) is generated at the bottom surface and potentially dissipated by inelastic collisions, encompassed through a coefficient of restitution,
$e_{w}$
.
At the top of the current
$z=h$
, there are the free-surface boundary conditions that fluid shear and normal stresses vanish, given by
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}z}=0,\quad \text{and}\quad p-\frac{4\unicode[STIX]{x1D707}}{3}\frac{\unicode[STIX]{x2202}w}{\unicode[STIX]{x2202}z}=0\quad \text{at }z=h.\end{eqnarray}$$
However, in addition, the flux of granular temperature vanishes and the solid-phase normal and shear stresses adopt small values, representing the surface as being the location where collisional behaviour ends and instead the particles follow ballistic trajectories (see Johnson et al. Reference Johnson, Nott and Jackson1990). Thus we enforce
$$\begin{eqnarray}\unicode[STIX]{x1D6F9}=0\quad \text{and}\quad (\unicode[STIX]{x1D70E}_{xz},\unicode[STIX]{x1D70E}_{zz})=\frac{\unicode[STIX]{x03C0}}{6}\unicode[STIX]{x1D70C}_{s}\left(\frac{\unicode[STIX]{x1D719}}{\unicode[STIX]{x1D719}_{m}}\right)^{2/3}gd(\sin \unicode[STIX]{x1D703},-\text{cos}\unicode[STIX]{x1D703})\quad \text{at }z=h.\end{eqnarray}$$
The boundary conditions (4.17)–(4.19) are of the same character as those employed by researchers in other flow regimes (see, for example, Jenkins & Berzi Reference Jenkins and Berzi2010) and as for the constitutive laws, the framework for analysing these flows is robust to variations in the closures used for these conditions.
4.1.4 Non-dimensionalisation of equations
We now identify typical dimensional scales for the dependent variables and assess the magnitude of the various terms in the governing equations. It is convenient to sum the downslope momentum equations of each phase (4.2) and (4.4) to eliminate the inter-phase drag so that
$$\begin{eqnarray}\unicode[STIX]{x1D70C}_{g}w_{g}\frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}z}=(\unicode[STIX]{x1D70C}_{s}\unicode[STIX]{x1D719}+\unicode[STIX]{x1D70C}_{g}(1-\unicode[STIX]{x1D719}))g\sin \unicode[STIX]{x1D703}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}z}\left(\unicode[STIX]{x1D70C}_{s}f_{1}\text{d}T^{1/2}\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}\right)+\unicode[STIX]{x1D707}_{g}\frac{\unicode[STIX]{x2202}^{2}u}{\unicode[STIX]{x2202}z^{2}}.\end{eqnarray}$$
In this expression, the key driving force is the downslope gravitational acceleration and it is this term that the other terms must balance. Since the density of the gas is much smaller than that of the solid phase and the effects of gas viscosity are negligible away from boundaries in this streamwise balance, we deduce that the dominant resistance is provided by the shear stress associated with the solid phase. Coarsely scaling the variables and assuming the volume fraction and the constitutive functions of it are of order unity,
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D70C}_{s}g\sin \unicode[STIX]{x1D703}\sim \unicode[STIX]{x1D70C}_{s}\text{d}T^{1/2}v/h^{2}. & & \displaystyle\end{eqnarray}$$
Furthermore, if the granular temperature is in local equilibrium between production and dissipation (an assumption that will be tested in the numerical solutions that follow), then from (4.11),
$$\begin{eqnarray}\displaystyle \text{d}T^{1/2}(v/h)^{2}\sim T^{3/2}/d; & & \displaystyle\end{eqnarray}$$
whence, the scaling for the velocity field is given by
$$\begin{eqnarray}v\sim \left(\frac{gh^{3}\sin \unicode[STIX]{x1D703}}{d^{2}}\right)^{1/2}.\end{eqnarray}$$
It is now convenient to introduce dimensionless variables, given by
$$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle \hat{z}=z/h,\quad \hat{u} =u\left(\frac{gh^{3}\sin \unicode[STIX]{x1D703}}{d^{2}}\right)^{-1/2},\quad \hat{v}=v\left(\frac{gh^{3}\sin \unicode[STIX]{x1D703}}{d^{2}}\right)^{-1/2}\\ \displaystyle {\hat{w}}=\frac{w}{w_{g}},\quad \hat{p}=\frac{p}{\unicode[STIX]{x1D70C}_{s}gh\cos \unicode[STIX]{x1D703}}\quad \text{and}\quad \hat{T}=\frac{T}{gh\sin \unicode[STIX]{x1D703}}.\end{array}\right\} & & \displaystyle\end{eqnarray}$$
This set of scalings for the dependent variables of fluidised flows differs from those for flows of unfluidised, collisional granular media (Woodhouse et al.
Reference Woodhouse, Hogg and Sellar2010). For non-fluidised flows, the granular agitation must provide sufficient normal stress to support the weight of the overlying layer. This would require the granular temperature to be of magnitude
$gh\cos \unicode[STIX]{x1D703}$
, which is considerably larger than the estimate deduced here (4.24) unless the motion is along relatively steep inclines (i.e. when
$\tan \unicode[STIX]{x1D703}\sim 1$
). For fluidised flows, however, granular temperature is generated by collisions but the imposed gas flow through the underlying particles provides most of the normal stress to balance the weight of the flowing layer. The granular temperature, therefore, is lower and consequently the shear stresses are lower, which in turn significantly increases the mobility of these flows. Hence these fluidised flows are characterised by relatively high flow speeds and relatively weak resistance.
The model is characterised by five dimensionless groups:
$$\begin{eqnarray}\left.\begin{array}{@{}c@{}}\displaystyle S=\tan \unicode[STIX]{x1D703},\quad R=\frac{\unicode[STIX]{x1D70C}_{g}}{\unicode[STIX]{x1D70C}_{s}},\quad \unicode[STIX]{x1D6FF}=\frac{d}{h},\\ \displaystyle W_{g}=\frac{\unicode[STIX]{x1D707}_{g}w_{g}}{\unicode[STIX]{x1D70C}_{s}d^{2}g\cos \unicode[STIX]{x1D703}},\quad \text{and}\quad St=\frac{\unicode[STIX]{x1D70C}_{s}d(g\sin \unicode[STIX]{x1D703}h)^{1/2}}{\unicode[STIX]{x1D707}_{g}}\unicode[STIX]{x1D6FF}^{2},\end{array}\right\}\end{eqnarray}$$
which represent respectively the inclination of the underlying boundary, the relative density of the gas to the solid phases, the size of the particles relative to the flow depth, the magnitude of the drag exerted by the fluidising gas flow relative to the weight of the granular layer and the reduced Stokes number, which compares particle inertia to fluid viscous effects.
It is also possible to define a particle Reynolds number,
$$\begin{eqnarray}Re=\frac{\unicode[STIX]{x1D70C}_{g}w_{g}d}{\unicode[STIX]{x1D707}_{g}}=\frac{RW_{g}St^{2}}{S\unicode[STIX]{x1D6FF}^{3}}.\end{eqnarray}$$
Notably, the Reynolds number defined in this way is independent of the inclination of the slope. The magnitude of the various model parameters for the experiments are set out in table 2. These scales may be used to show that the granular temperature dissipation and generation by viscous processes are negligible compared with direct particle interactions: (4.13) is satisfied when
$\unicode[STIX]{x1D6FF}^{2}/St\ll 1$
and (4.14) when
$W_{g}\ll 1$
.
Table 2. Range of values of the physical parameters in the experiments and the dimensionless groups derived from them as defined by (4.25).
We have five governing equations: mass conservation (4.1), downslope fluid momentum conservation (4.2), the combined normal momentum equation from which the fluid pressure has been eliminated (4.7), the downslope solids momentum equation (4.4) and the equation for the conservation of granular temperature (4.11). Non-dimensionalised, these become
$$\begin{eqnarray}\displaystyle & \displaystyle (1-\unicode[STIX]{x1D719}){\hat{w}}=1, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6FF}Re\frac{\unicode[STIX]{x2202}\hat{u} }{\unicode[STIX]{x2202}\hat{z}}=RSt(1-\unicode[STIX]{x1D719})+(1-\unicode[STIX]{x1D719})\unicode[STIX]{x1D6FF}^{2}\frac{\unicode[STIX]{x2202}^{2}\hat{u} }{\unicode[STIX]{x2202}\hat{z}^{2}}+(\hat{v}-\hat{u} )(f_{0}+f_{0}^{\ast }Re\,{\mathcal{U}}), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle -S\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{z}}(f_{2}\hat{T})+\frac{\unicode[STIX]{x1D6FF}W_{g}Re\,\unicode[STIX]{x1D719}}{(1-\unicode[STIX]{x1D719})}\frac{\unicode[STIX]{x2202}{\hat{w}}}{\unicode[STIX]{x2202}\hat{z}}=\unicode[STIX]{x1D719}(1-R)-W_{g}\frac{{\hat{w}}}{1-\unicode[STIX]{x1D719}}(f_{0}+f_{0}^{\ast }Re\,{\mathcal{U}}), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle 0=\unicode[STIX]{x1D719}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{z}}\left(f_{1}\hat{T}^{1/2}\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}}\right)-\frac{(f_{0}+f_{0}^{\ast }Re\,{\mathcal{U}})}{St}(\hat{v}-\hat{u} )+\frac{\unicode[STIX]{x1D719}\unicode[STIX]{x1D6FF}^{2}}{St}\frac{\unicode[STIX]{x2202}^{2}\hat{u} }{\unicode[STIX]{x2202}\hat{z}^{2}}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle 0=\unicode[STIX]{x1D6FF}^{2}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{z}}\left(f_{4}\hat{T}^{1/2}\frac{\unicode[STIX]{x2202}\hat{T}}{\unicode[STIX]{x2202}\hat{z}}\right)+\unicode[STIX]{x1D6FF}^{2}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{z}}\left(f_{4}^{\ast }\hat{T}^{3/2}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D719}}{\unicode[STIX]{x2202}\hat{z}}\right)+f_{1}\hat{T}^{1/2}\left(\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}}\right)^{2}-f_{3}\hat{T}^{3/2}, & \displaystyle\end{eqnarray}$$
where
${\mathcal{U}}$
measures the magnitude of the dimensionless relative velocity between the phases and is given by
${\mathcal{U}}^{2}=S^{2}(\hat{u} -\hat{v})^{2}/(W_{g}St)^{2}+{\hat{w}}^{2}$
. From (4.17), the dimensionless boundary conditions at the base
$(\hat{z}=0)$
are given by
$$\begin{eqnarray}\hat{u} =0,\quad f_{5}\unicode[STIX]{x1D6FF}\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}}=\hat{v}\quad \text{and}\quad -\!\unicode[STIX]{x1D6FF}f_{4}\frac{\unicode[STIX]{x2202}\hat{T}}{\unicode[STIX]{x2202}\hat{z}}=\frac{f_{6}}{\unicode[STIX]{x1D6FF}^{2}}\hat{v}^{2}-f_{7}\hat{T},\end{eqnarray}$$
while from (4.18), we enforce at the top surface
$(\hat{z}=1)$
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\hat{u} }{\unicode[STIX]{x2202}\hat{z}}=0,\quad \frac{\unicode[STIX]{x2202}\hat{T}}{\unicode[STIX]{x2202}\hat{z}}=0\quad \text{and}\quad \left(f_{1}\hat{T}^{1/2}\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}},f_{2}\hat{T}\right)=\frac{\unicode[STIX]{x03C0}\unicode[STIX]{x1D6FF}}{6}\left(\frac{\unicode[STIX]{x1D719}}{\unicode[STIX]{x1D719}_{m}}\right)^{2/3}(1,S^{-1}).\end{eqnarray}$$
Figure 2. The volume fraction,
$\unicode[STIX]{x1D719}(\hat{z})$
, velocity of the solid phase,
$\hat{v}(\hat{z})$
, and the granular temperature,
$\hat{T}(\hat{z})$
, as functions of the dimensionless depth within the current for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.5$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$S=0.1$
,
$St=10^{3}\unicode[STIX]{x1D6FF}^{2}$
,
$W_{g}=10^{-3}$
and (i)
$\unicode[STIX]{x1D6FF}=0.1$
, (ii)
$\unicode[STIX]{x1D6FF}=0.01$
and (iii)
$\unicode[STIX]{x1D6FF}=0.001$
. Also plotted are the asymptotic solutions (dotted lines), although these are often overlain by the full solution.
Figure 3. The volume fraction,
$\unicode[STIX]{x1D719}(\hat{z})$
, velocity of the solid phase,
$\hat{v}(\hat{z})$
, and the granular temperature,
$\hat{T}(\hat{z})$
, as functions of the dimensionless depth within the current for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.5$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$S=0.1$
,
$St=0.1$
,
$\unicode[STIX]{x1D6FF}=0.01$
, (i)
$W_{g}=0.5\times 10^{-3}$
, (ii)
$W_{g}=1.0\times 10^{-3}$
and (iii)
$W_{g}=2\times 10^{-3}$
. Also plotted are the asymptotic solutions (dotted lines).
Figure 4. The volume fraction,
$\unicode[STIX]{x1D719}(\hat{z})$
, velocity of the solid phase,
$\hat{v}(\hat{z})$
, and the granular temperature,
$\hat{T}(\hat{z})$
, as functions of the dimensionless depth within the current for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.3$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$St=0.1$
,
$\unicode[STIX]{x1D6FF}=0.01$
,
$W_{g}=10^{-3}$
, (i)
$S=0.1$
; (ii)
$S=0.2$
and (iii)
$S=0.3$
. Also plotted are the asymptotic solutions (dotted lines).
The system of governing differential equations (4.27)–(4.31) and boundary conditions (4.32)–(4.33) form a seventh-order differential boundary value problem. We use (4.27) to eliminate
${\hat{w}}$
in favour of
$1/(1-\unicode[STIX]{x1D719})$
and we also evaluate
$\text{d}^{2}\unicode[STIX]{x1D719}/\text{d}z^{2}$
in (4.31) by explicitly differentiating (4.29). The system is then integrated numerically. (For this task we employ the boundary value problem solver bvp4c in MatLab.) Example solutions for the volume fraction,
$\unicode[STIX]{x1D719}(\hat{z})$
, the granular temperature,
$\hat{T}(\hat{z})$
, and the velocity of the solid phase,
$\hat{v}(\hat{z})$
, are plotted in figures 2–4 for various values of the governing dimensionless parameters. We do not plot the gas velocity,
$\hat{u} (\hat{z})$
, because outside of thin basal boundary layers, it is indistinguishable from the solids velocity,
$\hat{v}(\hat{z})$
. This basal boundary layer exists because the while the gas phase satisfies a no-slip condition, the solid phase exhibits slip.
The general trends are that the volume fraction is approximately uniform while the granular temperature decreases with distance from the bottom boundary. Additionally there is a small slip velocity for the solid phase and the velocity shear decreases with distance from the boundary. There are some systematic deviations from these general trends, most notably in regions close the upper and lower boundaries. Since the interactions with the basal boundary are more dissipative than interactions with the constituent grains, the granular temperature decreases within a region in the vicinity of the boundary. This boundary effect is diminished as the flow becomes thicker (i.e. as
$\unicode[STIX]{x1D6FF}$
decreases, see figure 2), but is magnified where either the fluidising gas flow is increased (figure 3) and the slope is increased (figure 4). We also find that throughout the bulk of the domain, away from thin layers adjacent to the upper and lower boundaries, the production of granular temperature is in close balance with its dissipation, thus confirming the dimensional scales identified above (4.24).
It is of particular interest to evaluate the dimensionless flux of solids per unit width and the average concentration of particles, respectively given by
$$\begin{eqnarray}\hat{q}=\int _{0}^{1}\unicode[STIX]{x1D719}\hat{v}\,\text{d}\hat{z}\quad \text{and}\quad \overline{\unicode[STIX]{x1D719}}=\int _{0}^{1}\unicode[STIX]{x1D719}\,\text{d}\hat{z},\end{eqnarray}$$
and these are plotted as functions of the dimensionless parameters in figures 5–7.
From figure 5, we observe that the dimensionless volume flux per unit width,
$\hat{q}$
, and the average volume fraction,
$\overline{\unicode[STIX]{x1D719}}$
, do not vary strongly with the relative particle size,
$\unicode[STIX]{x1D6FF}$
. Thus, we deduce that boundary-related effects on the bulk characteristics are negligible for flows that are in excess of twenty particles thick. This result is of particular significance when unsteady shallow flows are analysed (§ 5).
The effects of the fluidisation velocity,
$W_{g}$
, are rather more subtle (see figure 6). Increasing
$W_{g}$
increases the normal support of the weight of the granular layer due to the gas flow and this lowers the concentration of the layer. However, the net volume flux,
$\hat{q}$
, does not vary monotonically with
$W_{g}$
. Indeed, for the parameters in figure 6, it attains a maximum at a dimensionless gas flow rate
$W_{g}=2.3\times 10^{-3}$
and at that value of
$W_{g}$
,
$\unicode[STIX]{x1D719}=0.43$
. This reflects the trade-off between the increased mobility but lower solids fraction of more dilute currents. Finally there is also relatively complex behaviour with increasing slope angle (figure 7). For the computations in this figure, as we increase the inclination, we also adjust
$W_{g}$
and
$St$
, but maintain
$Re$
constant (see (4.25) and (4.26)). We find that as the slope increases, the normal stress developed by the particle collisions increases with increasing granular temperature and this supplements the fluidising gas flow, leading to a progressively decreasing average volume fraction. However the dimensionless volume flux exhibits a more complicated dependency because while the velocity fields increase with increasing slope, the volume fraction diminishes and eventually becomes sufficiently dilute for
$\hat{q}$
to be maximised at some finite value of
$S$
. (For the parameters analysed in figure 7, the local maximum in the flux occurs at
$S=0.22$
.)
Figure 5. The dimensionless volume flux per unit width transported by the flowing layer,
$\hat{q}$
, and the depth-averaged volume fraction,
$\overline{\unicode[STIX]{x1D719}}$
, as functions of the relative particle size for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.50$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$S=0.10$
,
$St=10^{3}\unicode[STIX]{x1D6FF}^{2}$
and
$W_{g}=10^{-3}$
. Also plotted are the asymptotic solutions (dashed) and the simple approximate solutions (dotted).
Figure 6. The dimensionless volume flux per unit width transported by the flowing layer,
$\hat{q}$
, and the depth-averaged volume fraction,
$\overline{\unicode[STIX]{x1D719}}$
, as functions of the dimensionless strength of the fluidising gas flow,
$W_{g}$
, for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.50$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$S=0.10$
,
$St=0.10$
and
$\unicode[STIX]{x1D6FF}=10^{-2}$
. Also plotted are the asymptotic solutions (dashed) and the simple approximate solutions (dotted).
Figure 7. The dimensionless volume flux per unit width transported by the flowing layer,
$\hat{q}$
, and the depth-averaged volume fraction,
$\overline{\unicode[STIX]{x1D719}}$
, as functions of the slope of the underlying boundary,
$S=\tan \unicode[STIX]{x1D703}$
, for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.50$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$\unicode[STIX]{x1D6FF}=10^{-2}$
,
$W_{g}=9.95\times 10^{-4}(1+S^{2})^{1/2}$
,
$St=0.317S^{1/2}(1+S^{2})^{-1/4}$
. Also plotted are the asymptotic solutions (dashed) and the simple approximate solutions (dotted). Note the local maximum at
$S=0.21$
.
4.1.5 Asymptotic solution
In the bulk of the flow away from the boundaries, it is possible to deduce an asymptotic solution to the governing equation for the regime
$\unicode[STIX]{x1D6FF}\ll 1$
and
$R\ll 1$
. This regime will have a widespread validity as
$d\ll h$
in order to use a continuum approach, and for gas–solid flows
$\unicode[STIX]{x1D70C}_{g}\ll \unicode[STIX]{x1D70C}_{s}$
. From (4.28), we note that to leading order and away from boundaries, the downslope velocities of the two phases must be equal (
$\hat{u} =\hat{v}+O(1)$
). Furthermore, from (4.31) there is a local balance between granular temperature production and dissipation such that
$$\begin{eqnarray}f_{1}\left(\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}}\right)^{2}=f_{3}\hat{T},\end{eqnarray}$$
provided the volume fraction of particles is not too small (i.e.
$\unicode[STIX]{x1D719}\gg \unicode[STIX]{x1D6FF}$
, so that the ‘conductive’ effects of the granular temperature remain negligible).
The governing equations for the normal and perpendicular momentum balances, (4.29) and (4.30), are then given by
$$\begin{eqnarray}\displaystyle & \displaystyle -S\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{z}}(f_{2}\hat{T})=\unicode[STIX]{x1D719}-W_{g}\frac{f_{0}}{(1-\unicode[STIX]{x1D719})^{2}}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{z}}((f_{1}f_{3})^{1/2}\hat{T})=-\unicode[STIX]{x1D719}. & \displaystyle\end{eqnarray}$$
These reduced governing equations neglect shear stresses in the fluid phase, which become non-negligible as the basal boundary is approached and which allow the velocities of the two phases to differ. Also the ‘conduction’ of granular temperature is neglected because, to leading order, we find a balance between production and dissipation (4.35). When
$\unicode[STIX]{x1D6FF}\ll 1$
, the lower boundary layer corresponds to a region within which the velocity of the solid phase is small, while the upper boundary layer to a region within which the granular temperature is small. The leading-order boundary conditions are then given by
$\hat{T}(1)=0$
and the volume fraction at the base is given by
$\unicode[STIX]{x1D719}(0)=\unicode[STIX]{x1D719}_{0}$
, which is determined by substituting for
$\unicode[STIX]{x2202}\hat{v}/\unicode[STIX]{x2202}\hat{z}$
from (4.35) into (4.32),
$$\begin{eqnarray}0=\left(\frac{f_{6}\,f_{5}^{2}f_{3}}{f_{1}}-f_{7}\right)\hat{T},\end{eqnarray}$$
where the constitutive functions
$(f_{i})$
are evaluated at
$\unicode[STIX]{x1D719}=\unicode[STIX]{x1D719}_{0}$
. The basal volume fraction is thus a function of
$e$
,
$e_{w}$
,
$\unicode[STIX]{x1D713}$
and
$\unicode[STIX]{x1D719}_{m}$
.
Rearranging (4.36) and (4.37) and denoting
$f=(f_{1}f_{3})^{1/2}$
, we find that
$$\begin{eqnarray}\displaystyle & \displaystyle (f{\dot{f}}_{2}-{\dot{f}}f_{2})\hat{T}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D719}}{\unicode[STIX]{x2202}\hat{z}}=\left(-\unicode[STIX]{x1D719}+\frac{W_{g}f_{0}}{(1-\unicode[STIX]{x1D719})^{2}}\right)\frac{f}{S}+f_{2}\unicode[STIX]{x1D719}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle ({\dot{f}}f_{2}-f{\dot{f}}_{2})\frac{\unicode[STIX]{x2202}\hat{T}}{\unicode[STIX]{x2202}\hat{z}}=\left(-\unicode[STIX]{x1D719}+\frac{W_{g}f_{0}}{(1-\unicode[STIX]{x1D719})^{2}}\right)\frac{{\dot{f}}}{S}+{\dot{f}}_{2}\unicode[STIX]{x1D719}, & \displaystyle\end{eqnarray}$$
where
$\dot{~}$
denotes differentiation with respect to
$\unicode[STIX]{x1D719}$
. It is straightforward to integrate numerically these coupled first-order equations subject to the boundary conditions
$\hat{T}(1)=0$
and
$\unicode[STIX]{x1D719}(0)=\unicode[STIX]{x1D719}_{0}$
. The solutions are plotted in figures 2–4 and it is evident that these asymptotic solutions accurately reproduce the numerical solution of the complete system (very often in these figures, the asymptotic curves are indistinguishable from the numerical solution of the complete system).
There is also an even simpler approximate solution. The coupled system admits a homogeneous solution
$\unicode[STIX]{x1D719}(z)=\overline{\unicode[STIX]{x1D719}}$
when
$$\begin{eqnarray}\overline{\unicode[STIX]{x1D719}}-\frac{f_{0}W_{g}}{(1-\overline{\unicode[STIX]{x1D719}})^{2}}=\frac{S\overline{\unicode[STIX]{x1D719}}f_{2}}{(f_{1}f_{3})^{1/2}},\end{eqnarray}$$
where the constitutive functions are evaluated at
$\unicode[STIX]{x1D719}=\overline{\unicode[STIX]{x1D719}}$
. In this case, the temperature gradient is constant,
$\unicode[STIX]{x2202}\hat{T}/\unicode[STIX]{x2202}\hat{z}=-\unicode[STIX]{x1D706}$
, with
$\unicode[STIX]{x1D706}=\overline{\unicode[STIX]{x1D719}}/f$
. This solution is ‘attracting’ in the sense that trajectories in the phase plane
$(\unicode[STIX]{x1D719}(\hat{z}),\hat{T}(\hat{z}))$
approach it when
$$\begin{eqnarray}f{\dot{f}}_{2}-{\dot{f}}f_{2}<0,\end{eqnarray}$$
which in turn demands that
$\overline{\unicode[STIX]{x1D719}}>\overline{\unicode[STIX]{x1D719}}_{c}(e)$
and that
$S<S_{c}(e)$
(see figure 8
c). If these inequalities are not held then the reduced system evolves towards a state different from a uniform volume fraction
$(\unicode[STIX]{x1D719}=\overline{\unicode[STIX]{x1D719}})$
and may not admit solutions at all. Physically, when
$S>S_{c}$
, the dissipation of granular temperature, here encapsulated through a coefficient of restitution
$e$
, is insufficient to allow for a steady balance between the weight of the flowing layer, the fluidising drag and the normal stresses generated through particle interactions (a balance expressed by (4.36)). When
$e\lessapprox 0.85$
, we find that this limitation does not play for parameter values associated with the flows considered in this study and that a flow with a homogeneous volume fraction of particles provides a good representation of the more complete dynamics (see figures 2–4).
Figure 8. (a) The limiting volume fraction,
$\overline{\unicode[STIX]{x1D719}}_{c}$
and slope,
$S_{c}$
, for which the uniform volume fraction is the ‘attracting’ solution as functions of the coefficient of restitution; (b) the volume flux per unit width,
$\hat{q}$
, (c) the average volume fraction,
$\overline{\unicode[STIX]{x1D719}}$
, as functions of the slope,
$S$
, for
$W_{g}=10^{-3}$
and varying values of the coefficient of restitution and (d) the product of
$\overline{\unicode[STIX]{x1D719}}$
and the mobility factor,
$F$
, as a function of
$\overline{\unicode[STIX]{x1D719}}$
for various values of the coefficient of restitution.
When the dimensionless fluidisation velocity
$W_{g}$
and the slope
$S$
are set, the average volume fraction within the current,
$\overline{\unicode[STIX]{x1D719}}$
, can be calculated using (4.41). Figure 8(a) shows the effect of
$e$
on the curves of
$\overline{\unicode[STIX]{x1D719}}$
as a function of slope
$S$
when the fluidisation flow is constant (
$W_{g}=10^{-3}$
). The curves in this plot are continued up to the maximum value of the slope,
$S_{c}(e)$
, for which the reduced model leads to a homogeneous volume fraction and it can be seen that the slope at which this can be achieved is successively reduced as dissipation in the collisions is decreased. From figure 8(b), for a given slope,
$S$
, and fluidisation gas flow rate,
$W_{g}$
, flows with lower coefficients of restitution lead to higher dimensionless volume fluxes per unit width. This is simply rationalised: as a high coefficient of restitution implies reduced dissipation and high granular temperatures. Consequentially there are higher stresses and greater resistance to the downslope motion.
Since the granular temperature must vanish at the surface
$\hat{z}=1$
to leading order, we find for the simple approximate solution with uniform volume fraction that the granular temperature is given by
$$\begin{eqnarray}\hat{T}=\frac{\overline{\unicode[STIX]{x1D719}}}{(f_{1}f_{3})^{1/2}}(1-\hat{z}),\end{eqnarray}$$
and the velocity field of the solid phase is given by
$$\begin{eqnarray}\hat{v}=F{\textstyle \frac{2}{3}}(1-(1-\hat{z})^{3/2}),\end{eqnarray}$$
where
$F=(\overline{\unicode[STIX]{x1D719}}^{2}f_{3}/f_{1}^{3})^{1/4}$
. The scaled slip velocity at the wall
$Ff_{5}\unicode[STIX]{x1D6FF}$
can be added to (4.44), but when
$\unicode[STIX]{x1D6FF}\ll 1$
, it is negligible. The velocity profile (4.44) is similar to the dimensionless ‘Bagnold’ velocity profile up to the factor
$F$
, which controls the mobility of the flowing layer and is influenced by the fluidisation velocity. In many situations the approximate solution provides a very good representation of the solution to the complete system (see figures 2 to 4).
Also in this regime
$(\unicode[STIX]{x1D6FF}\ll 1)$
, the approximate solution yields
$$\begin{eqnarray}\hat{q}={\textstyle \frac{2}{5}}\overline{\unicode[STIX]{x1D719}}F,\end{eqnarray}$$
and this is plotted in figures 5–7, once again illustrating the utility of this asymptotic solution. The quantity
$\overline{\unicode[STIX]{x1D719}}F$
thus plays a crucial role in determining the dimensionless flux,
$\hat{q}$
and in figure 8(d), we plot its dependence on volume fraction for a range of values of
$e$
. We note that
$\overline{\unicode[STIX]{x1D719}}F$
is maximised for
$\unicode[STIX]{x1D719}\approx 0.41$
(with the precise value weakly dependent on
$e$
) and vanishes both when
$\unicode[STIX]{x1D719}$
vanishes and when it approaches maximum packing. This variation reflects the balance between fast-moving dilute flows and slow-moving concentrated flows, leading to a flux maximum at intermediate values
$(\unicode[STIX]{x1D719}\approx 0.41)$
. Finally, we note that a dimensional estimate of the depth of a fluidised current may be obtained
$$\begin{eqnarray}h=\left(\frac{5}{2}\frac{d}{\overline{\unicode[STIX]{x1D719}}F}\frac{q_{0}}{\sqrt{g\sin \unicode[STIX]{x1D703}}}\right)^{2/5},\end{eqnarray}$$
where
$q_{0}$
is the dimensional flux per unit width at the source and the effects of slip at the wall has been neglected.
4.2 Experimental measurements of fully developed flows
4.2.1 Depth of currents
A typical velocity profile is shown in figure 9, superimposed on a captured image from the recording of an experiment. There is a small slip at the lower boundary, an approximately linear increase in velocity with distance from the wall until a maximum velocity is attained and then a progressive drop to zero. There appears to be a top to the current where the particle volume fraction suddenly drops and, there,
$h_{vis}(\equiv h)$
:
$h_{vis}$
is greater than the height at which the maximum velocity is attained. Above
$h_{vis}$
particles are detected, but their velocity drops with increasing height and it has a large variance. This is consistent with there being a ballistic region into which individual particles may be projected. The height at which particle velocity drops to zero is the top of the entire current and is denoted by
$h_{max}$
. The depth of the current can fluctuate a little with time (see figures 18 and 19
a). As a result, the averaging process will occasionally include points that are above the average height of the current so that the averaged velocity at these points will be necessarily lower than in the bulk of the current and the variance will be higher. The measured heights for the different currents are summarised in table 3 compared with the height estimated from (4.46).
Figure 9. Image of a granular current with the measured velocity superimposed onto it (solid line). The image is an average-intensity composite of the images used in the PIV measurements. Note the non-zero (slip) velocity at the base (
$v_{0}$
). The dashed lines indicate the 95 % confidence interval using (2.2). Also shown are the height where the velocity profile drops to zero,
$h_{max}$
, the maximum visible height of the current,
$h_{vis}$
, and the height of the peak of the velocity profile,
$z(v=v_{max})$
. For this figure
$\unicode[STIX]{x1D703}=15^{\circ }$
,
$Q=42.5~\text{cm}^{3}~\text{s}^{-1}$
,
$w_{g}/u_{mf}=1.5$
.
Table 3. Various estimates of height in steady-state currents.
$h_{vis}$
is measured from photographs of the currents such as that in figure 9. It may be compared with the prediction,
$h$
, calculated from (4.46) with
$F$
based on
$\bar{\unicode[STIX]{x1D719}}_{est}$
(see table 4),
$e=0.85$
and
$\unicode[STIX]{x1D713}=0.50$
.
$h_{max}$
is found directly as the height at which particle velocity drops to zero. In all cases
$w_{g}=1.5u_{mf}$
.
The asymptotic solution yielded an approximate formula linking height and source volume flux (4.46), and this is shown in figure 10. In this figure, fixed values of
$\overline{\unicode[STIX]{x1D719}}$
were used because from (4.41),
$\overline{\unicode[STIX]{x1D719}}$
depends on
$\unicode[STIX]{x1D703}$
and it was not possible to find a solution for the full experimental range of
$\unicode[STIX]{x1D703}$
for a fixed value of
$e=0.85$
. The decreasing effect of
$\overline{\unicode[STIX]{x1D719}}$
on mobility as its value approaches
$\overline{\unicode[STIX]{x1D719}}_{c}=0.40$
reflects the effect it has on mobility
$\overline{\unicode[STIX]{x1D719}}F$
shown in figure 8(d). The theoretical formula contains no adjusted parameters and is an approximation to the more complete description, but it yields a reasonable quantitative representation of the relationship between the depth of the flowing layer, the source flux and channel inclination.
Figure 10. The height of the flowing layer as a function of the scaled volume flux,
$q_{0}/\sqrt{g\sin \unicode[STIX]{x1D703}}$
, for a range of channel inclinations. The model curves are for (4.46) with
$\unicode[STIX]{x1D719}_{m}=0.63$
and
$e=0.85$
.
4.2.2 Velocity profiles
The ensemble-averaged velocity profiles for
$10^{\circ }$
and
$15^{\circ }$
slopes measured half-way along the tank are shown in figure 11. The features described for figure 9 are reflected in each of the profiles. The structure of the currents posited by the model is similar in structure to the experimental measurements up to the point at which the velocity is a maximum (see figure 2). Above this point, the variance of the velocity measurements increases markedly as the velocity drops off. As described above this effect could be because the PIV is simply sampling ballistic particle trajectories.
The velocity profiles measured on inclinations of
$\unicode[STIX]{x1D703}=10^{\circ }$
, with source fluxes
$Q=33~\text{cm}^{3}~\text{s}^{-1}$
and
$Q=34~\text{cm}^{3}~\text{s}^{-1}$
are distinct from each other with the former forming a current that is deeper and much faster than the latter. It is not clear why there is such a large difference between the measured profiles. One possibility could be that the flowing layer exhibits multiple states for the same imposed flux. This behaviour is known in models of unfluidised granular flows (Woodhouse et al.
Reference Woodhouse, Hogg and Sellar2010); however, we failed to find such multiplicity of solutions in the governing equations examined in this study in the parameter regime corresponding to these experiments. The relatively fast and expanded flow with
$Q=33~\text{cm}^{3}~\text{s}^{-1}$
leads to small estimates of particle volume fraction with an excessive portion of the flow where
$h>h_{vis}$
(see below, § 4.2.3. Based on
$h_{vis}$
,
$\bar{\unicode[STIX]{x1D719}}=0.33$
; based on
$h_{max}$
,
$\bar{\unicode[STIX]{x1D719}}=0.22$
), which seem physically unlikely, and so this experimental run is not reported further.
Figure 12 shows some instantaneous velocity profiles half-way along the tank for shallower slope angles whose motion may not be steady. The velocity profiles for the
$3^{\circ }$
and
$5^{\circ }$
slopes are similar in character to the averaged profiles for steeper slopes.
Figure 11. Ensemble-averaged velocity profiles for steady flows for
$\unicode[STIX]{x1D703}=10^{\circ }$
and
$\unicode[STIX]{x1D703}=15^{\circ }$
. Error bars are 95 % confidence limits calculated from (2.2) with
$m=5$
and
$n=16$
.
Figure 12. Measured velocity profiles from experiment for individual fluidised granular flows that are potentially unsteady when they are at shallow angles. Error bars are calculated from (2.1) with values of
$c_{1}$
determined from the correlation functions calculated during the PIV analysis.
4.2.3 Particle volume fraction
The expectation from the model is that the particle volume fraction
$\unicode[STIX]{x1D719}$
would be approximately uniform within the fluidised currents. It is not possible to analyse the degree to which
$\unicode[STIX]{x1D719}$
is a function of position in the currents from our experimental set-up; however, it is evident that towards the top of the current above
$h_{vis}$
,
$\unicode[STIX]{x1D719}$
drops sharply so that the current loses its opacity. This is consistent with the decreasing particle velocity there. Some bubbles were seen in the currents, but not many and at a small number of sites, even once the current had traversed the bottom of the container, and those that were seen were small in size.
From the measured velocity profiles, it is possible to estimate the average volume fraction,
$\bar{\unicode[STIX]{x1D719}}$
, by integrating the particle velocity profiles and dividing the measured particle flow rate
$\unicode[STIX]{x1D719}_{0}Q$
by the result. The results
$\unicode[STIX]{x1D719}_{meas}$
are shown in table 4, using
$h_{vis}$
as the overall depth of the current, to give measured volume fraction of the currents (
$\bar{\unicode[STIX]{x1D719}}_{meas}$
). The measured values of particle volume fraction can be compared with estimated values,
$\unicode[STIX]{x1D719}_{est}$
, which have been calculated using (4.41).
There can be quite good agreement between
$\bar{\unicode[STIX]{x1D719}}_{est}$
and
$\bar{\unicode[STIX]{x1D719}}_{meas}$
despite several inherent uncertainties in their computation. Overall, the values of
$\bar{\unicode[STIX]{x1D719}}_{meas}$
were comparable to the values of
$\bar{\unicode[STIX]{x1D719}}_{est}$
and with the typical values of
$\unicode[STIX]{x1D719}=0.50{-}0.60$
for static fluidised beds (Epstein & Young Reference Epstein and Young1962). In between
$h_{vis}$
and
$h_{max}$
particles are present, but in practice the fall off of velocity above
$h_{vis}$
is sufficiently rapid that this makes very little difference to calculations of
$\unicode[STIX]{x1D719}_{meas}$
: if
$\unicode[STIX]{x1D719}_{meas}$
is calculated on the basis of the top of the current being
$h_{max}$
rather than
$h_{vis}$
, then its value decreases by less than 0.02, except when
$\unicode[STIX]{x1D703}=15^{\circ }$
and
$Q=11.2~\text{cm}^{3}~\text{s}^{-1}$
, when it reduces by 0.09.
Table 4. Measured estimates for
$Q$
and
$\bar{\unicode[STIX]{x1D719}}$
from integration of the velocity profiles. Values of
$\bar{\unicode[STIX]{x1D719}}_{est}$
are found from solving (4.41), and those for
$\bar{\unicode[STIX]{x1D719}}_{meas}$
from integration of the measured velocity profiles up to
$h_{vis}$
. The errors are calculated using 95 % confidence limits calculated from (2.2).
4.2.4 Slip at the wall
Specularity coefficients have not been measured for fluidised granular currents and some modellers think they should not be used at all for individual collisions (Goldschmidt et al.
Reference Goldschmidt, Beetstra and Kuipers2004). Their value is sufficiently badly defined that in their investigations of bubbling fluidised beds of glass particles Altantzis, Bates & Ghoniem (Reference Altantzis, Bates and Ghoniem2015) used values between
$10^{-4}$
and
$0.5$
and Li, Grace & Bi (Reference Li, Grace and Bi2010) from
$0$
to
$0.5$
. Our computations showed that apart from within relatively narrow layers close to the boundary, the magnitude of the specularity coefficient had relatively little effect upon the flow profiles. It is, however, possible to estimate the value of
$\unicode[STIX]{x1D713}$
from the directly measured slip velocities and gradients using the boundary condition (4.17) and the definition of
$f_{5}$
in table 1, so that in terms of dimensional variables
$$\begin{eqnarray}\left.\unicode[STIX]{x1D713}=\frac{2\sqrt{3}}{\unicode[STIX]{x03C0}}\frac{\unicode[STIX]{x1D719}_{m}}{\bar{\unicode[STIX]{x1D719}}}\frac{f_{1}}{g_{0}}\frac{d}{v(0)}\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}\right|_{z=0}.\end{eqnarray}$$
The results are shown in table 5, and it can be seen that a measured average value of
$\unicode[STIX]{x1D713}$
is 0.28. The values for
$\unicode[STIX]{x1D713}$
shown in table 5 should be treated as only indicative as the velocity gradients close to the wall are shallow and the slip velocities small, so small variations in the velocity profiles can result in significant changes in the value of
$\unicode[STIX]{x1D713}$
; however, despite these uncertainties, the value of
$\unicode[STIX]{x1D713}$
is reasonably consistent. The effect of
$\unicode[STIX]{x1D713}$
on the velocity profiles is to introduce a slip velocity proportional to
$\unicode[STIX]{x1D6FF}f_{5}$
. Even for the relatively large values of
$\unicode[STIX]{x1D713}$
estimated here, the magnitude of this dimensionless term is relatively small.
Table 5. Measured slip velocities at the wall,
$v(0)$
, and velocity gradients
$\unicode[STIX]{x2202}v/\unicode[STIX]{x2202}z|_{z=0}$
with the resulting estimate for slip length and for
$\unicode[STIX]{x1D713}$
using (4.47). The slip length is expressed in terms of particle diameters and is defined as
$v(0)/\text{d}(\unicode[STIX]{x2202}v/\unicode[STIX]{x2202}z)|_{z=0}$
.
$e=0.85$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
.
4.2.5 Scaling of the velocity profiles
The measured velocity profiles scaled as
$\hat{v}$
and
$\hat{z}$
are plotted in figure 13. With the exception of the lowest flow rate when
$\unicode[STIX]{x1D703}=15^{\circ }$
, the data collapse well in the region close to the wall. The model predicts dependence of the theoretical curve on the slope angle for the flow through its influence on
$\bar{\unicode[STIX]{x1D719}}$
and hence
$F$
, but this is not reflected in the experimental curves for which the scaling appears to eliminate the effect of
$S$
.
$F$
is also affected by the value of
$e$
. Increasing
$e$
causes a decrease in the predicted dimensionless velocity: the theoretical curves shift towards the left and the difference between the curves for
$\unicode[STIX]{x1D703}=10^{\circ }$
and
$15^{\circ }$
becomes less (see inset, figure 13).
Figure 13. Velocity profiles with measurements scaled as (4.24) compared with the model scaled velocity profile under steady state (4.44) for each angle where
$e=0.85$
. The curves correspond to
$\bar{\unicode[STIX]{x1D719}}=0.51$
when
$\unicode[STIX]{x1D703}=10^{\circ }$
and
$\bar{\unicode[STIX]{x1D719}}=0.41$
when
$\unicode[STIX]{x1D703}=15^{\circ }$
.
$W_{g}=8.09\times 10^{-4}$
when
$\unicode[STIX]{x1D703}=10^{\circ }$
and
$W_{g}=8.24\times 10^{-4}$
when
$\unicode[STIX]{x1D703}=15^{\circ }$
. The inset graph shows the effect of the value of
$e$
on the model solutions with the chain-dot curves corresponding to
$e=0.95$
and the dotted curves to
$e=0.75$
.
5 Unsteady, developing flows on slopes
The mathematical model may be extended to unsteady, developing flows of fluidised currents down slopes, but it now takes a somewhat different form because streamwise gradients can no longer be neglected. In this situation, we analyse the motion in the ‘lubrication’ regime, for which a representative streamwise length scale,
$L$
, far exceeds a representative length scale perpendicular to the boundary,
$H(H/L\ll 1)$
. This means that accelerations perpendicular to the boundary are negligible and that to leading order, the normal stresses adopt the ‘hydrostatic’ balance given by (4.6). We again assume that the flows are many particles thick,
$\unicode[STIX]{x1D6FF}\ll 1$
, the density of the gas is negligible relative to that of the solids,
$R\ll 1$
, and the effects of gas viscosity are negligible (see § 4.1.5). The leading-order dimensional momentum equations of each phase parallel with the incline then take the form,
$$\begin{eqnarray}\displaystyle & \displaystyle 0=-(1-\overline{\unicode[STIX]{x1D719}})\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}x}+\unicode[STIX]{x1D6FD}(v-u) & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}\frac{\text{D}v}{\text{D}t}=-\overline{\unicode[STIX]{x1D719}}\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}x}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{xz}}{\unicode[STIX]{x2202}z}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{xx}}{\unicode[STIX]{x2202}x}+\overline{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D70C}_{s}g\sin \unicode[STIX]{x1D703}-\unicode[STIX]{x1D6FD}(v-u), & \displaystyle\end{eqnarray}$$
where the average volume fraction in the flowing layer, determined by the balance between the fluidising gas flow and the particle weight, is given by (4.41). It is interesting to note from (5.1) that now there must be a leading-order difference between the downslope velocities of the two phases. Furthermore, since the flow is spatially and temporally evolving, we must include the inertia of the solid phase, which in (5.2) is given by the term
$\unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}\,\text{D}v/\text{D}t$
(here
$\text{D}/\text{D}t$
denotes the material derivative). Summing these two momentum balances to eliminate the inter-phase drag and assuming further that the stresses in the solid phase are isotropic
$(\unicode[STIX]{x1D70E}_{xx}=\unicode[STIX]{x1D70E}_{zz})$
and that the current is in hydrostatic balance (4.6), we deduce that
$$\begin{eqnarray}\unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}\frac{\text{D}v}{\text{D}t}=\unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}g\left(\sin \unicode[STIX]{x1D703}-\cos \unicode[STIX]{x1D703}\frac{\unicode[STIX]{x2202}h}{\unicode[STIX]{x2202}x}\right)+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{xz}}{\unicode[STIX]{x2202}z}.\end{eqnarray}$$
The granular temperature of the flow in this regime is assumed to be in local balance between its production and dissipation through collisions, and these processes dominate its advective and diffusive transport. We proceed further by adopting the appropriate dimensionless scales following the distinguished scaling identified in § 4.1.4 and embodied in the dimensionless variables of (4.24) and (4.24). However, here we non-dimensionalise the depth of the flowing the current
$h$
by a representative depth-scale
$H$
,
$({\hat{h}}=h/H)$
and, additionally
$$\begin{eqnarray}\hat{v}=\frac{v}{F(g\sin \unicode[STIX]{x1D703}H^{3}/d^{2})^{1/2}},\quad \hat{x}=\frac{x}{L}\quad \text{and}\quad \hat{t}=F\left(\frac{g\sin \unicode[STIX]{x1D703}H^{3}}{d^{2}}\right)^{1/2}\frac{t}{L}.\end{eqnarray}$$
Then, using the approximate form of the solution established in § 4.1.5, we find that the depth-integrated, dimensionless momentum equation is given by
$$\begin{eqnarray}{\mathcal{R}}\left(\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{t}}\int _{0}^{{\hat{h}}}\hat{v}\,\text{d}\hat{z}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{x}}\int _{0}^{{\hat{h}}}\hat{v}^{2}\,\text{d}\hat{z}\right)+\frac{1}{\unicode[STIX]{x1D6E5}}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{x}}\left(\frac{{\hat{h}}^{2}}{2}\right)={\hat{h}}-\left(\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}}\right)_{z=0}^{2},\end{eqnarray}$$
where
$\unicode[STIX]{x1D6E5}=L\tan \unicode[STIX]{x1D703}/H$
and
$$\begin{eqnarray}{\mathcal{R}}=\frac{F^{2}H^{3}}{d^{2}L}.\end{eqnarray}$$
In this setting, as for Kumaran (Reference Kumaran2014),
${\mathcal{R}}$
measures the relative magnitude of the inertial to resistive terms. Unlike a Reynolds number for viscous fluid flows, it features only the length scales in the problem and
$F$
, because both the inertial terms and the shear stresses are proportional to the square of velocity.
To proceed further we assume that the velocity field adopts similar dependence to (4.44) on distance from the boundary and this permits the evaluation of the integral and boundary quantities in terms of the average velocity and the depth of the layer:
$$\begin{eqnarray}\int _{0}^{{\hat{h}}}\hat{v}\,\text{d}\hat{z}={\hat{h}}\overline{v},\quad \int _{0}^{{\hat{h}}}\hat{v}^{2}\,\text{d}\hat{z}=\frac{5}{4}\overline{v}^{2}{\hat{h}}\quad \text{and}\quad \left.\frac{\unicode[STIX]{x2202}\hat{v}}{\unicode[STIX]{x2202}\hat{z}}\right|_{z=0}=\frac{5\overline{v}}{2{\hat{h}}}.\end{eqnarray}$$
To complete the model, we express conservation of mass,
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}{\hat{h}}}{\unicode[STIX]{x2202}\hat{t}}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\hat{x}}(\overline{v}{\hat{h}})=0.\end{eqnarray}$$
This system is subject to the boundary condition that we impose a sustained source of particles at the origin
$$\begin{eqnarray}\overline{\unicode[STIX]{x1D719}}{\hat{h}}\overline{v}=\hat{q}_{0}\quad \text{at }\hat{x}=0.\end{eqnarray}$$
Additionally, if the flow is supercritical then we must enforce the Froude number at the source. We impose the initial condition,
${\hat{h}}(\hat{x},0)=0$
and the current forms a front
$\hat{x}=\hat{x}_{f}(\hat{t})$
, such that
${\hat{h}}(\hat{x}_{f},t)=0$
.
Before constructing solutions, it is convenient to relate the height and streamwise length scales. We choose
$L=H/\!\tan \unicode[STIX]{x1D703}$
and thus
$\unicode[STIX]{x1D6E5}=1$
. The lubrication regime requires that streamwise lengths far exceed the thickness of the flow; since the current is expanding in the streamwise direction, this regime is inevitably entered after sufficient time. However the adopted scaling may imply that in terms of these variables, the initial evolution may not be well captured by the lubrication assumption. We further choose the dimensional height
$H$
, using (4.46) so that
$$\begin{eqnarray}H=\left(\frac{d}{\overline{\unicode[STIX]{x1D719}}F}\frac{q_{0}}{\sqrt{g\sin \unicode[STIX]{x1D703}}}\right)^{2/5}.\end{eqnarray}$$
The governing equations now entail the single dimensionless parameter,
${\mathcal{R}}=\tan \unicode[STIX]{x1D703}F^{2}H^{2}/d^{2}.$
We construct travelling wave solutions for the dimensionless height and velocity fields. We write
${\hat{h}}(\hat{x},\hat{t})\equiv {\hat{h}}(\hat{x}-c\hat{t})$
and
$\hat{v}(\hat{x},\hat{t})\equiv \hat{v}(\hat{x}-c\hat{t})$
, where
$c$
is the dimensionless wave speed which is to be determined. Conservation of mass then implies that
$\hat{v}=c$
and in particular, the front speed is given by
$$\begin{eqnarray}\hat{x}_{f}=c\hat{t}.\end{eqnarray}$$
Balance of momentum leads to
$$\begin{eqnarray}\frac{{\mathcal{R}}c^{2}}{4}{\hat{h}}^{\prime }+{\hat{h}}{\hat{h}}^{\prime }={\hat{h}}-\frac{25}{4}\frac{c^{2}}{{\hat{h}}^{2}},\end{eqnarray}$$
where a prime denotes differentiation with respect to
$\unicode[STIX]{x1D702}=\hat{x}-c\hat{t}$
. Distant from the front, the flowing layer carries a constant volume flux of material determined by the source conditions (
$c{\hat{h}}\rightarrow 1$
as
$\unicode[STIX]{x1D702}\rightarrow -\infty$
). Here we note that the travelling wave solutions do not satisfy the source condition precisely at
$\hat{x}=0$
, but instead it is satisfied as
$\unicode[STIX]{x1D702}\rightarrow -\infty$
. For these flows, we find that the current adjusts over a short distance behind the front to a uniform depth and velocity and thus the travelling wave solution provides an accurate representation of the solution for the flow. Thus, we deduce that the position of the front and the far-field depth are given by
$$\begin{eqnarray}c=\left({\textstyle \frac{2}{5}}\right)^{2/5}\quad \text{and}\quad {\hat{h}}\rightarrow {\hat{h}}_{\infty }=\left({\textstyle \frac{5}{2}}\right)^{2/5}.\end{eqnarray}$$
Figure 14. The position of the front of the fluidised current as a function of time for flows along channels of varying inclinations with varying source fluxes and fluidising gas flows.
Experimental measurements of the distance travelled by fluidised currents with time are shown in figure 14, in which the inclination of the channel, the source volume flux and the fluidising gas velocity were varied; the measured flow speeds ranged over a factor of five. The measurements for the fully fluidised currents (
$u_{g}/u_{mf}>1$
) after scaling are shown in figure 15. From (5.11),
$\hat{x}_{f}$
should be proportional to
$\hat{t}$
and this is true even when the slope angles are small. Individually, the currents display a constant speed; however, the measured speeds can be significantly different from that expected from the model (
$q_{0}/\unicode[STIX]{x1D719}h$
).
The degree of data collapse for different parameters –
$w_{g}$
, the nominal flow rate
$Q_{nom}$
, and
$\unicode[STIX]{x1D703}$
– is shown in figure 16. For all three parameters, the scaling eliminates much of the scatter, but it is not fully eliminated. The collapse of data onto different lines with the same values of
$\unicode[STIX]{x1D703}$
and
$Q_{nom}$
for
$w_{g}$
is excellent. It can be quite good for
$\unicode[STIX]{x1D703}$
, especially at low
$\hat{t}$
. For
$Q_{nom}$
the collapse of data is often incomplete. There will be some variation reflecting the difference in value of the true value of
$Q$
from the nominal value
$Q_{nom}$
. It can be seen in figure 16(a) that the effect of
$w_{g}$
is small, but significant; however, it is eliminated after scaling, as shown in figure 16(b).
A systematic omission from our model is the effect of sidewall drag and this could provide an additional resistance to motion, thus slowing the speed of propagation. In appendix B, we analyse the effects of the sidewalls when the height of the current,
$H$
, is much less than the breadth of the channel,
$B$
. We demonstrate that there is a weak retardation to the dimensionless speed proportional to
$(H/B)^{2}$
when
$H/B\ll 1$
. We analysed the speed of the flow from the data plotted in figure 14 and found no systematic dependence on
$H/B$
and thus there is no evidence that these relatively shallow currents were significantly slowed by sidewall effects.
Figure 16. The effect of different variables on the distance travelled by flows down slopes before and after scaling. The graphs from the top downwards show the effects of varying
$w_{g}$
,
$Q$
and
$\unicode[STIX]{x1D703}$
. The points in the graphs correspond to the legend in figure 15. For
$w_{g}$
there are four families of curves:
$\unicode[STIX]{x1D703}=3^{\circ }$
,
$Q_{nom}=35~\text{cm}^{3}~\text{s}^{-1}$
(
$+$
);
$\unicode[STIX]{x1D703}=5^{\circ }$
,
$Q_{nom}=15~\text{cm}^{3}~\text{s}^{-1}$
(○);
$\unicode[STIX]{x1D703}=5^{\circ }$
,
$Q_{nom}=20~\text{cm}^{3}~\text{s}^{-1}$
(▫);
$\unicode[STIX]{x1D703}=15^{\circ }$
,
$Q_{nom}=35~\text{cm}^{3}~\text{s}^{-1}$
(
$\times$
). For
$Q$
there are three families of curves corresponding to
$\unicode[STIX]{x1D703}=3^{\circ }$
(
$+$
),
$5^{\circ }$
(○),
$10^{\circ }$
(▫), with
$w_{g}=1.5u_{mf}$
. For
$\unicode[STIX]{x1D703}$
there are four families of curves: when
$w_{g}=1.5u_{mf}$
,
$Q_{nom}=15~\text{cm}^{3}~\text{s}^{-1}$
(
$+$
),
$35~\text{cm}^{3}~\text{s}^{-1}$
(○),
$60~\text{cm}^{3}~\text{s}^{-1}$
(▫), and when
$w_{g}=2.0u_{mf}$
,
$Q_{nom}=35~\text{cm}^{3}~\text{s}^{-1}$
(
$\times$
).
5.1 Time taken to establish steady uniform behaviour
The dimensionless profile of a fluidised current moving down a slope is determined from (5.12) and is implicitly given by
$$\begin{eqnarray}\int _{0}^{{\hat{h}}/{\hat{h}}_{\infty }}\frac{s^{2}(A+s)}{s^{3}-1}\,\text{d}s=\frac{\unicode[STIX]{x1D702}}{{\hat{h}}_{\infty }},\end{eqnarray}$$
where
$A={\mathcal{R}}/(4{\hat{h}}_{\infty }^{3})$
. We plot in figure 17 the height of the travelling wave of material as a function of distance from the front for various values of the inertial parameter,
${\mathcal{R}}$
, and note that the length scale over which the flow adjusts to the uniform depth,
${\hat{h}}_{\infty }$
, increases with increasing
${\mathcal{R}}$
. One measure of the streamwise length,
$\unicode[STIX]{x1D6E5}_{\unicode[STIX]{x1D716}}$
, over which the flow attains its uniform depth is given by evaluating when
${\hat{h}}(-\unicode[STIX]{x1D6E5}_{\unicode[STIX]{x1D716}})={\hat{h}}_{\infty }(1-\unicode[STIX]{x1D716})$
, which when
$\unicode[STIX]{x1D716}\ll 1$
is given by
$$\begin{eqnarray}-\!\frac{\unicode[STIX]{x1D6E5}_{\unicode[STIX]{x1D716}}}{h_{\infty }}=1+\frac{1}{3}(A+1)\log \unicode[STIX]{x1D716}+\frac{1}{6}(2A-1)\log 3-\frac{\sqrt{3}\unicode[STIX]{x03C0}}{18}+\cdots\end{eqnarray}$$
At a fixed location, it is then possible to evaluate the dimensionless time-scale over which the uniform depth is established,
$\hat{t}_{\unicode[STIX]{x1D716}}=\unicode[STIX]{x1D6E5}_{\unicode[STIX]{x1D716}}/c$
.
Figure 17. The scaled height of the current,
${\hat{h}}/{\hat{h}}_{\infty }$
, as a function of position
$\unicode[STIX]{x1D702}/{\hat{h}}_{\infty }=(\hat{x}-c\hat{t})/{\hat{h}}_{\infty }$
for parameter values (i)
$A={\mathcal{R}}/(4{\hat{h}}_{\infty }^{3})=0$
; (ii)
$A=1$
; (iii)
$A=5$
; and (iv)
$A=10$
.
Figure 18. Shape of the currents over time for fluidised flows at different slope angles. The numbers on the contours indicate time in seconds after the release of material. A spatial Gaussian filter with a kernel size of 0.5 cm has been applied to the contours in order to improve clarity. Note that the vertical scale is different in each plot and the horizontal scale is different in (a) from the other diagrams. (a)
$\unicode[STIX]{x1D703}=0^{\circ }$
,
$Q=49.1~\text{cm}^{3}~\text{s}^{-1}$
. Time interval 0.25 s between contours. Note the growth and decay of the surface waves near the origin. (b)
$\unicode[STIX]{x1D703}=3^{\circ }$
,
$Q=15.0~\text{cm}^{3}~\text{s}^{-1}$
(
$A=0.09$
). The time interval between contours is 0.5 s. The flow becomes uniform after about 4 s in this case. (c)
$\unicode[STIX]{x1D703}=5^{\circ }$
,
$Q=79.5~\text{cm}^{3}~\text{s}^{-1}$
(
$A=1.57$
). The time interval between contours is 1 s with a uniform state achieved after about 2 s. (d)
$\unicode[STIX]{x1D703}=10^{\circ }$
,
$Q=38.2~\text{cm}^{3}~\text{s}^{-1}$
(
$A=0.75$
). The time interval between contours is 0.5 s and the current reaches uniform state within 1 s.
Figure 19. The transition from unsteady to uniform behaviour of fluidised currents down slopes. In both panels,
$t^{\ast }$
is the time at which the current reaches the measuring position
$x^{\ast }=$
50 cm from the source. (a) Shows the measured depth of the flowing current at a fixed position as a function of time at a fixed point for currents on slopes between 0–15
$^{\circ }$
with the same nominal source flux (
$60~\text{cm}^{3}~\text{s}^{-1}$
). Experimental conditions are:
$\unicode[STIX]{x1D703}=3^{\circ }$
,
$Q=59.44~\text{cm}^{3}~\text{s}^{-1}$
;
$\unicode[STIX]{x1D703}=5^{\circ }$
,
$Q=58.75~\text{cm}^{3}~\text{s}^{-1}$
;
$\unicode[STIX]{x1D703}=10^{\circ }$
,
$Q=58.76~\text{cm}^{3}~\text{s}^{-1}$
and
$\unicode[STIX]{x1D703}=15^{\circ }$
,
$Q=56.38~\text{cm}^{3}~\text{s}^{-1}$
. (b) Shows the scaled height of the current,
${\hat{h}}/{\hat{h}}_{\infty }$
, as a function of scaled time after the front reaches
$x^{\ast }$
using the scales of (5.4) and using the definition of
$H$
in (5.10).
$\unicode[STIX]{x1D716}=0.05$
and when
$\unicode[STIX]{x1D703}=3^{\circ }$
,
$A=0.28$
,
$\hat{t}_{\unicode[STIX]{x1D716}}=1.38$
;
$\unicode[STIX]{x1D703}=5^{\circ }$
,
$A=0.45$
,
$\hat{t}_{\unicode[STIX]{x1D716}}=1.60$
;
$\unicode[STIX]{x1D703}=10^{\circ }$
,
$A=1.07$
,
$\hat{t}_{\unicode[STIX]{x1D716}}=2.41$
;
$\unicode[STIX]{x1D703}=15^{\circ }$
,
$A=2.56$
,
$\hat{t}_{\unicode[STIX]{x1D716}}=4.38$
.
Figure 18 shows examples of the development of the fluidised currents at different angles of channel inclination. In figure 18, the flows on the steeper slopes relatively rapidly attain a uniform state in which the current does not vary along the apparatus, whereas those on shallower slopes and with smaller sources fluxes take much longer to approach this state. Flows along horizontal channels never approach a uniform state; instead currents adopt the shape of a wedge and do not progress at constant speed. These are analysed in § 6.
Figure 19 shows the change of height with time for five flows with the same nominal flow rate half-way along the apparatus before and after scaling. It can be seen that the scaled times at which the currents achieve a constant height (and systematically with
$A$
), as would be expected, but they are an order of magnitude larger than
$\hat{t}_{\unicode[STIX]{x1D716}}$
. For the expected values of
$\hat{t}_{\unicode[STIX]{x1D716}}$
, the currents would have to achieve their constant height very quickly, almost instantly, and for the values of
$A$
corresponding to the experimental flows, from figure 17 the front of the currents would be expected to be ‘blunt nosed’, with quite steep gradients of height at the front of the current. In fact, the front of the currents (figure 18) had a relatively shallow gradient.
6 Horizontal flows
Figure 20. The position of the front of the fluidised current as a function of time for varying source fluxes and fluidising gas flows.
Figure 21. Change in shape of horizontal fluidised currents with no scaling for
$Q=49.05~\text{cm}^{3}~\text{s}^{-1}$
and
$w_{g}/u_{mf}=1.5$
. The different profiles are drawn at 0.25 s intervals. The dashed line corresponds to
$t=1~\text{s}$
and the dotted line to
$t=3~\text{s}$
.
Figure 20 shows the distance travelled by the front of a fluidised flow over horizontal surface as a function of time. It is evident that the currents do not travel at a constant speed. The shapes of the currents are shown in figure 21 and, ignoring the disturbance at the start of the currents at the point that they are poured into the system, they have an approximately triangular shape, though one with a low aspect ratio (i.e. their extent far exceeds their depth). Furthermore they flow through ‘bulk’ motion, not through the build-up of lamina arising from the constant avalanching down the current’s top surface seen for non-fluidised granular flows. They must also be scaled differently because the length and time scales introduced in (5.4) become singular when
$\unicode[STIX]{x1D703}=0$
. To this end, we introduce the characteristic height scale,
$\tilde{H}$
and define the following dimensionless variables
$$\begin{eqnarray}\tilde{h}=\frac{h}{\tilde{H}},\quad \tilde{v}=\frac{v}{F(g\tilde{H}^{3})^{1/2}/d},\quad \tilde{x}=\frac{x}{\tilde{H}},\quad \text{and}\quad \tilde{t}=\frac{Ft(g\tilde{H})^{1/2}}{d}.\end{eqnarray}$$
The depth-integrated expression of momentum balance is then given by
$$\begin{eqnarray}\tilde{{\mathcal{R}}}\left(\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\tilde{t}}(\tilde{h}\tilde{v})+\frac{5}{4}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\tilde{x}}(\tilde{h}\tilde{v}^{2})\right)+\frac{\tilde{h}}{2}\frac{\unicode[STIX]{x2202}\tilde{h}}{\unicode[STIX]{x2202}\tilde{x}}=-\frac{25}{4}\frac{\tilde{v}^{2}}{\tilde{h}^{2}},\end{eqnarray}$$
instead of (5.5), where the residual dimensionless parameter
$\tilde{{\mathcal{R}}}=(FH/d)^{2}$
is the relative magnitude of inertial to resistive forces.
Conservation of mass is given by
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\tilde{h}}{\unicode[STIX]{x2202}\tilde{t}}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\tilde{t}}(\tilde{v}\tilde{h})=0.\end{eqnarray}$$
The appropriate dimensional depth scale,
$\tilde{H}$
, is determined from the source flux,
$$\begin{eqnarray}\tilde{H}=\left(\frac{q_{0}dg^{2}}{F\overline{\unicode[STIX]{x1D719}}}\right)^{1/5},\end{eqnarray}$$
so that the boundary condition is given by
$$\begin{eqnarray}\tilde{v}\tilde{h}=1\quad \text{at }x=0.\end{eqnarray}$$
Flows over horizontal surfaces decelerate as the basal drag is no longer balanced by a sustained downslope acceleration. Thus, at sufficiently early times the flow speeds and depths are set by source conditions, and after the flow has propagated for sufficient time, the resistive forces become non-negligible and the motion enters a dynamical regime in which the drag force balance the streamwise gradients of the hydrostatic pressure and the inertial forces are negligible. Analogously to Hogg & Woods (Reference Hogg and Woods2001), simple scaling shows that this regime is fully attained when
$\tilde{t}\gg \tilde{{\mathcal{R}}}$
. In this scenario, we deduce from (6.2) that
$$\begin{eqnarray}\tilde{v}=\frac{2}{5}\left(-\tilde{h}^{3}\frac{\unicode[STIX]{x2202}\tilde{h}}{\unicode[STIX]{x2202}\tilde{x}}\right)^{1/2}\end{eqnarray}$$
and consequentially from (6.3)
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\tilde{h}}{\unicode[STIX]{x2202}\tilde{t}}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\tilde{x}}\left(\frac{2}{5}\tilde{h}^{5/2}\left(-{\displaystyle \frac{\unicode[STIX]{x2202}\tilde{h}}{\unicode[STIX]{x2202}\tilde{x}}}\right)^{1/2}\right)=0,\end{eqnarray}$$
subject to the source condition
$$\begin{eqnarray}\frac{2}{5}\tilde{h}^{5/2}\left(-{\displaystyle \frac{\unicode[STIX]{x2202}\tilde{h}}{\unicode[STIX]{x2202}\tilde{x}}}\right)^{1/2}=1.\end{eqnarray}$$
Equation (6.7) may be integrated numerically to reveal the evolution of the front position as a function of time and the variation of the depth of the current along its length; however, for these currents flowing over a horizontal surface, we may also construct a quasi-analytical similarity solution for the motion.
First, we determine the gearing between spatial and temporal scales that underpins the similarity solution for unsteady flow over a horizontal surface. To do this we scale and balance terms in the governing equation (6.7) and boundary condition (6.8). This yields
$$\begin{eqnarray}\frac{\tilde{h}}{\tilde{t}}\sim \frac{\tilde{h}^{3}}{\tilde{x}^{3/2}}\quad \text{and}\quad \frac{\tilde{h}^{3}}{\tilde{x}^{1/2}}\sim 1.\end{eqnarray}$$
Thus we deduce that
$\tilde{x}\sim \tilde{t}^{6/7}$
and
$\tilde{h}\sim \tilde{t}^{1/7}$
. We may then seek a similarity solution of the form
$$\begin{eqnarray}\displaystyle & \displaystyle \tilde{h}=K^{3/4}\tilde{t}^{1/7}{\mathcal{H}}(y), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \tilde{x}_{f}(\tilde{t})=K\tilde{t}^{6/7}, & \displaystyle\end{eqnarray}$$
where
$K$
is a dimensionless constant to be determined as part of the solution and
$y=\tilde{x}/\tilde{x}_{f}(\tilde{t})$
. On substitution in the governing equation (6.7), this gives
$$\begin{eqnarray}{\textstyle \frac{1}{7}}{\mathcal{H}}-{\textstyle \frac{6}{7}}y{\mathcal{H}}^{\prime }+{\textstyle \frac{2}{5}}[(-{\mathcal{H}}^{\prime })^{1/2}{\mathcal{H}}^{5/2}]^{\prime }=0,\end{eqnarray}$$
where a prime denotes differentiation with respect to
$y$
. This ordinary differential equation (6.12) is to be integrated subject to the boundary conditions
$$\begin{eqnarray}{\mathcal{H}}(1)=0\quad \text{and}\quad {\textstyle \frac{2}{5}}K^{2/7}(-{\mathcal{H}}^{\prime })^{1/2}{\mathcal{H}}^{5/2}=1\quad \text{at }y=0.\end{eqnarray}$$
The location
$y=1$
is a singular point of the differential equation (6.12); we therefore start the numerical integration at
$y=1-\unicode[STIX]{x1D716}$
$(\unicode[STIX]{x1D716}\ll 1)$
, noting that
$$\begin{eqnarray}{\mathcal{H}}(1-\unicode[STIX]{x1D716})=\left({\textstyle \frac{30}{7}}\right)^{1/2}\unicode[STIX]{x1D716}^{1/4}\left(1-{\textstyle \frac{1}{60}}\unicode[STIX]{x1D716}+{\textstyle \frac{43}{28\,880}}\unicode[STIX]{x1D716}^{2}+\cdots \,\right).\end{eqnarray}$$
It is then straightforward to integrate the differential equation (6.12) numerically and evaluate
${\mathcal{H}}(0)=2.038$
,
${\mathcal{H}}^{\prime }(0)=-0.479$
, and so
$K=0.753$
. The dimensional expression for distance covered by a horizontal, fluidised current with time is then
$$\begin{eqnarray}x_{f}=0.753\left(\frac{gF^{2}q_{0}^{4}}{d^{2}\unicode[STIX]{x1D719}^{4}}\right)^{1/7}t^{6/7}.\end{eqnarray}$$
We plot in figure 22 the similarity solution for the height profile along the current noting that, again, the model predicts a blunt-nosed current that advances along the channel.
Figure 22. The scaled height of the current as a function of downslope distance at various instances of time from the similarity solution for unsteady propagation along a horizontal channel.
The scaled distance against time is shown in figure 23, and the data are collapsed sufficiently for the power-law form of the curve to appear to be reasonable, though the value of the exponent is different from that predicted. However, the shape of the current predicted by the scaled model is very different from the experimental measurements, taking the form of nearly flat current with a snub nose (figure 24).
6.1 Flow within a narrow channel
One of the differences between the experimentally-realised flows over horizontal surfaces and those down slopes is that the former are significantly thicker than the latter, and so it is possible for the sidewalls to have a strong influence on their development. Here, we analyse the motion of a fluidised current as it flows within a narrow channel of width
$B$
, between sidewalls for which the streamwise extent of the flow far exceeds the depth of the current
$(L\gg H)$
, which in turn far exceeds the width of the flow
$(H\gg B)$
. For this regime, it is possible to simplify the governing equations for the unsteady motion down an incline on the basis that gradients across the flow are much greater than those in any other direction. In this scenario, the dynamical balance is somewhat different from that analysed in § 4 and the resulting governing equations for the unsteady evolution of the thickness of the flow are also different (§ 5).
Our derivation of the governing equation in a narrow channel is developed from the dimensional expressions presented in § 3 and then depth integrated to establish a shallow layer model; however, it will be shown that it leads to a similarity solution with a different gearing between the spatial and temporal variables. Here we only present the governing equations for flows along horizontal channels, but the inclusion of a channel gradient is a straightforward generalisation and could lead to travelling wave solutions analogous to § 5.
It is assumed that the solid particles are fully fluidised by an imposed gas flow and attain a state in which the volume fraction is uniform,
$\unicode[STIX]{x1D719}=\overline{\unicode[STIX]{x1D719}}$
. Since the flow is relatively thin, vertical accelerations are negligible and the motion is governed by hydrostatic balance given by
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}z}(p-\unicode[STIX]{x1D70E}_{zz})=-\unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}g.\end{eqnarray}$$
In this expression we have neglected the contribution due to the weight of the gas phase since
$\unicode[STIX]{x1D70C}_{g}/\unicode[STIX]{x1D70C}_{s}\ll 1$
. In the downslope direction, after neglecting terms proportional to the density and viscosity of the gas, the combined momentum equation of both phases to leading order is given by
$$\begin{eqnarray}\unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}\frac{\text{D}v}{\text{D}t}=-\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}(p-\unicode[STIX]{x1D70E}_{zz})+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{xy}}{\unicode[STIX]{x2202}y},\end{eqnarray}$$
where
$y$
is the distance across the channel. Here it has been assumed that the normal stresses of the solid phase are isotropic
$(\unicode[STIX]{x1D70E}_{xx}=\unicode[STIX]{x1D70E}_{zz})$
and that gradients across the flow dominate all others. To complete this model, we introduce the granular temperature, which provided the channel is much wider than the grain size, is in local equilibrium between its production and dissipation. Then we may write
$$\begin{eqnarray}f_{1}d^{2}\left(\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}y}\right)^{2}=f_{3}T,\end{eqnarray}$$
and the constitutive law for the shear stress is given by
$$\begin{eqnarray}\unicode[STIX]{x1D70E}_{xy}=f_{1}d\unicode[STIX]{x1D70C}_{s}T^{1/2}\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}y}.\end{eqnarray}$$
Finally, by eliminating the fluid pressure from the normal force balances of each phase (see (4.7)), we find that
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70E}_{zz}}{\unicode[STIX]{x2202}z}=\overline{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D70C}_{s}g-\frac{\unicode[STIX]{x1D6FD}w_{g}}{(1-\overline{\unicode[STIX]{x1D719}})^{2}}.\end{eqnarray}$$
The volume fraction, temperature and, therefore, the normal stress component,
$\unicode[STIX]{x1D70E}_{zz}$
, are independent of
$z$
to leading order, and thus we deduce that
$$\begin{eqnarray}0=\overline{\unicode[STIX]{x1D719}}\unicode[STIX]{x1D70C}_{s}g-\frac{\unicode[STIX]{x1D6FD}w_{g}}{(1-\overline{\unicode[STIX]{x1D719}})^{2}}.\end{eqnarray}$$
This expression determines the average volume fraction as a function of the fluidising gas flux and marks an important departure from the shallow layer model of § 5, because to leading order the solid stresses do not contribute to the support of the granular layer.
We progress by assuming that the velocity field of the solid phase exhibits cross-stream dependence, which is identical to that found in fully developed flows with the shear in the vertical plane. Thus we write
$$\begin{eqnarray}v=\overline{v}\frac{5}{3}\left(1-\left|1-\frac{2y}{B}\right|^{3/2}\right),\end{eqnarray}$$
and consequentially
$\unicode[STIX]{x2202}v/\unicode[STIX]{x2202}y=5\overline{v}/B$
at
$y=0$
. The dimensional governing equations then express conservation of mass and after sufficient time has passed so that inertia is negligible, a balance between the streamwise pressure gradient and the sidewall stresses is given by
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}h}{\unicode[STIX]{x2202}t}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}(\overline{v}h)=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle gh\frac{\unicode[STIX]{x2202}h}{\unicode[STIX]{x2202}x}=-\frac{50d^{2}\overline{v}^{2}h}{F^{2}B^{3}}. & \displaystyle\end{eqnarray}$$
For these flows over a horizontal surface we adopt the dimensionless variables given by
$$\begin{eqnarray}\overline{h}=h/\overline{H},\quad \overline{x}=x/\overline{H}\quad \text{and}\quad \overline{t}=(q_{0}t)/(\overline{\unicode[STIX]{x1D719}}\overline{H}^{2}B),\end{eqnarray}$$
where the height scale,
$\overline{H}$
, is determined by
$$\begin{eqnarray}\overline{H}=\frac{q_{0}d}{F\overline{\unicode[STIX]{x1D719}}(gB^{5})^{1/2}}.\end{eqnarray}$$
The dimensionless governing equation then becomes
$$\begin{eqnarray}\frac{\unicode[STIX]{x2202}\overline{h}}{\unicode[STIX]{x2202}\overline{t}}+\frac{1}{5\sqrt{2}}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\overline{x}}\left(\overline{h}\left(-\frac{\unicode[STIX]{x2202}\overline{h}}{\unicode[STIX]{x2202}\overline{x}}\right)^{1/2}\right)=0,\end{eqnarray}$$
subject to
$$\begin{eqnarray}\frac{1}{5\sqrt{2}}\overline{h}\left(-\frac{\unicode[STIX]{x2202}\overline{h}}{\unicode[STIX]{x2202}\overline{x}}\right)^{1/2}=1\quad \text{at }x=0.\end{eqnarray}$$
We may construct a similarity solutions to this governing equation (6.27) and boundary condition (6.28) by first deducing the gearing between the spatial and temporal scales. The governing equation demands
$\overline{h}/\overline{t}\sim \overline{h}^{3/2}/\overline{x}^{3/2}$
, while the boundary condition leads to
$\overline{h}^{3/2}\sim \overline{x}^{1/2}$
. Thus we deduce that
$\overline{x}\sim \overline{t}^{3/4}$
,
$\overline{h}\sim \overline{t}^{1/4}$
and that the similarity solutions of the following form may be sought
$$\begin{eqnarray}\overline{h}=C^{3}\overline{t}^{1/4}{\mathcal{H}}(y)\quad \text{and}\quad \overline{x}_{f}=C\overline{t}^{3/4},\end{eqnarray}$$
where
$y=\overline{x}/\overline{x}_{f}(\overline{t})$
and
$C$
is a constant to be determined. In dimensional variables,
$$\begin{eqnarray}x_{f}=C\left(\frac{q_{0}F^{2}gB^{2}}{\overline{\unicode[STIX]{x1D719}}d^{2}}\right)^{1/4}t^{3/4}.\end{eqnarray}$$
On substitution of
$\bar{h}$
into the governing equation (6.27), we deduce that
$$\begin{eqnarray}{\textstyle \frac{1}{4}}{\mathcal{H}}-{\textstyle \frac{3}{4}}y{\mathcal{H}}^{\prime }+{\textstyle \frac{1}{5\sqrt{2}}}[{\mathcal{H}}(-{\mathcal{H}}^{\prime })^{1/2}]^{\prime }=0,\end{eqnarray}$$
subject to
${\mathcal{H}}(1)=0$
and
$C^{4}{\mathcal{H}}(0)[-{\mathcal{H}}^{\prime }(0)]^{1/2}=5\sqrt{2}$
. The similarity differential equation (6.31) is singular at
$y=1$
and the numerical solution may be initiated from a series solution valid close to that value, given by
$$\begin{eqnarray}\displaystyle {\mathcal{H}}(1-s)={\textstyle \frac{225}{9}}(s-s^{2}/5+\cdots \,), & & \displaystyle\end{eqnarray}$$
when
$s\ll 1$
. It is then straightforward to integrate (6.31) to compute the height profile, shown in figure 25, and the dimensionless constant
$C=0.5434$
.
The scaled distance against time for the experiments is shown in figure 26. There is reasonably good collapse of the data (better than in figure 23). Furthermore, it appears to follow a power law, though the exponent of the power law is slightly small than that predicted. The measured, scaled shape of the currents is plotted in figure 27. The predicted shape is now much closer to that of the experiments, although the scaling does not collapse completely all the measured data.
Figure 25. The height of the current as a function of distance along a narrow horizontal channel at various instances of time.
Figure 27. The scaled depth of the fluidised current along a horizontal channel,
${\mathcal{H}}$
, as a function of the scaled position,
$y=\overline{x}/\overline{x}_{f}$
. The data are drawn from figure 21 scaled using (5.4). The model solution is from (6.29) and is represented by the chain-dotted line.
7 Discussion and conclusions
This investigation of fluidised granular currents reveals important distinctions in their dynamical properties from both dry granular flows and static fluidised beds. Most significantly, there is substantial and sustained shear in the velocity profiles. Consequentially particles are driven into each other and this provides a mechanism for the generation of stresses. In the regime we investigated, the inertia of individual grains remains relatively high and thus the particles interact with each other through dissipative collisions, and it these interactions that lead to the shear stress that balances the downslope gravitational acceleration. Granular flows in the absence of fluidisation must generate sufficient normal stresses to support the weight of the flowing layer and thus typically their shear stresses are also relatively large; however, fluidisation changes the balance of forces acting on the current. The fluidising gas flow provides most of the normal support to the flowing layer and thus both normal and shear stresses from the particulate phase are reduced relative to their non-fluidised counterparts, leading to flows that are much more mobile.
In this study we have formed a framework of modelling fluidised currents based on the solid-phase stresses that is generated from collisions between the particles. The degree of agitation in the system is measured through the granular temperature and constitutive laws are employed to determine the stress tensor in terms of the gradients of the velocity field, the granular temperature and the volume fraction of solids, as well as several material parameters. For the regime studied here in which the flows are many particles thick, a local balance emerges between the generation and dissipation of granular temperature. This leads to an accurate asymptotic model for the complete dynamics, in which the flowing material is essentially modelled by a nonlinear local rheology. Furthermore, this reduction leads to a Bagnold-like expression between the flow depth and the flux of particles carried by the current, with the volume fraction, determined by the fluidising gas flow, contributing to this relationship. Although this approximation is a simplified description of the more complete dynamics, it embodies the key processes of these flows: the fluidised currents are granular flows in which the fluidisation affects the normal support of the layer.
The experimental measurements provide encouraging support for the model. For example, without tuning through empirical factors, the predictions are quite close to the measured flow depths and flow speed in both the uniform steady state and the transient state as it becomes established. Additionally, when applied to flows along horizontal channels, the model is able to predict the unsteady motion to reveal both the progressive deceleration and the growth in flow depth. There are, however, some systematic features in the measurements that are not reproduced in the model. Perhaps the most significant of these is the decrease in particle velocity towards the top of the layer. This feature is absent from the model and presumably corresponds to particles in the ‘free-board’ of the fluidised layer (i.e. the region above the dense current within which the volume fraction of the particles is reduced). Such a dilute layer is subject to slightly different dynamical interactions: the role of particle collisions becomes much reduced and the particles may saltate, and are potentially intermittently suspended above the denser layer below. The model predictions are also dependent upon the material properties that characterise the collisions between the particles and the boundaries. These can be difficult to measure directly, but the specularity coefficient,
$\unicode[STIX]{x1D713}$
, and the boundary coefficient of restitution,
$e_{w}$
, only play a significant role for a relatively thin boundary layer in the dynamical regime considered in this study. It is arguable that the boundary conditions require further research to refine and sharpen their formulation.
One important feature that emerges from the modelling framework is the determination of the volume fraction of particles in the fluidised current. Here we have assumed that the flows are relatively dense and that the Ergun equation provides an appropriate representation of the volume fraction dependence of the drag due to the fluidising gas flow. Other expressions could easily be used in its place (see Nott & Jackson Reference Nott and Jackson1992; Agrawal et al.
Reference Agrawal, Loezos, Syamlal and Sundaresan2001; Oger & Savage Reference Oger and Savage2013). However, perhaps of greater significance is whether there are ‘bubbles’, or inhomogeneities in the volume fraction within the fluidised current. Patches of increased voidage locally provide paths through which the fluidising gas can more readily flow and thus it is possible for the layer to exhibit fluctuations or instabilities on relatively rapid time scales. Since the local volume fraction affects the mobility of the flowing layer, one might expect fluctuations in volume fraction and velocity to be correlated and consequentially to influence the bulk dynamics. Bubbles may also affect the particle volume fraction in the bed. The classical model of fluidised beds (Toomey & Johnstone Reference Toomey and Johnstone1952) proposes that all the gas in excess of that necessary to fluidise the particles forms bubbles so that the particle volume fraction in the bulk of the flow is insensitive to
$w_{g}$
: this is contrast with (4.41). The lack of dependence on
$\unicode[STIX]{x1D703}$
of the experimental scaled velocity profiles in figure 13 also suggests that
$\bar{\unicode[STIX]{x1D719}}$
may change less with conditions than might be expected. In contrast to static fluidised beds, the stability of fully developed fluidised flow down inclines has not been assessed (Jackson Reference Jackson2000) and this appears to be an interesting topic for future research. Indeed it is intriguing that linear shear flows of unfluidised granular materials appear to exhibit transient linearised growth, but asymptotic stability (Savage Reference Savage1992; Scmid & Kytomaa Reference Scmid and Kytomaa1994). It would be interesting to investigate whether these properties carry over to sheared fluidised motions.
Our modelling framework and experimental methods could be extended to a number of related flow problems. First, one could investigate fluidised currents that are generated by instantaneous or non-sustained releases, that are not fully fluidised and for which the fluidising gas flow is localised to the region close to the source. These flows would be largely unsteady and, in situations where the fluidisation is not maintained, would introduce additional mechanisms for generating resistive shear stresses as the contact friction begins to become important. Non-monodisperse granular materials would also be interesting to investigate because the onset of full fluidisation is dependent upon grain size and the proportion of each particle component (Formisani Reference Formisani1991). It is possible that mixtures of particles segregate according to size and generate an inhomogeneous flowing current in terms of composition and therefore, the average volume fraction (
$\bar{\unicode[STIX]{x1D719}}$
). Finally, we comment that liquid fluidised systems may pose additional challenges since it is likely that viscous forces at the particle scale are non-negligible and that collisions are strongly affected by lubrication pressure in the fluid between particles.
Although direct applications have not been the focus of our study, our model formulation could be readily applied to larger-scale flows, either in industrial contexts or in nature. There are a number of practical implications of our results. For example, the transport of granular materials when they are fluidised is likely to be more efficient on even shallowly inclined surfaces than on horizontal surfaces. In addition, the transport is unlikely to be greatly improved by an increase in the gas flow rate,
$w_{g}$
, once the granular materials are fully fluidised. This is because it only directly affects the average volume fraction,
$\bar{\unicode[STIX]{x1D719}}$
, and for practical materials
$\bar{\unicode[STIX]{x1D719}}$
can strongly depend on their characteristics (e.g. the bubble-free expansion seen in small, light Geldart (Reference Geldart1973) group A particles), as can the value of the coefficient of restitution
$e$
.
The assumption at the heart of the modelling framework is that inelastic particulate collisions generate stresses that provide the resistance to motion and this is likely to be the case for larger-scale flows. Of particular note is that for steady flow, no assumption is made for the relative importance of inertial and resisting forces and hence the resulting model is valid for flows of arbitrary scale. The results show that, at least for some granular flows, full understanding of their nature can only be reached if full account is taken of both the interactions between particles and those between particles and the interstitial fluid.
Acknowledgements
The authors thank three anonymous reviewers who helped to improve our manuscript, as well as O. Pouliquen for his handling of it. This work was funded from a grant under the UK NERC Environmental Mathematics and Statistics programme (NER/S/E/2004/12600). This research was supported also in part by the National Science Foundation under grant no. NSF PHY11-25915 and A.J.H. also acknowledges support from Max Planck Institute for the Physics of Complex Systems (Two-Phase Continuum Models for Geophysical Particle-Fluid Flows). M.A.G. carried out part of this work while holding a University Research Fellowship provided by the Institute of Advanced Study at the University of Bristol. This paper is LabEx Clervolc contribution no. 259.
Appendix A. Extended kinetic theory
Figure 28. The volume fraction,
$\unicode[STIX]{x1D719}(\hat{z})$
, velocity of the solid phase,
$\hat{v}(\hat{z})$
and the granular temperature,
$\hat{T}(\hat{z})$
, as functions of the dimensionless depth within the current for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.5$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$S=0.1$
,
$St=10^{3}\unicode[STIX]{x1D6FF}^{2}$
,
$\unicode[STIX]{x1D6FF}=0.01$
and
$W_{g}=10^{-3}$
for extended kinetic theory (solid lines) and ‘standard’ kinetic theory (dashed lines).
In this appendix we analyse the consequences for the predicted flow field of employing the extended kinetic theory proposed by Jenkins (Reference Jenkins2007) and recently used to compute unfluidised flow down inclined planes by Jenkins & Berzi (Reference Jenkins and Berzi2010, Reference Jenkins and Berzi2012) and Berzi (Reference Berzi2014). In essence, the extension to kinetic theory is based upon the realisation that at higher concentrations, particles begin to form structures in the flow that have a correlation length in excess of their own diameter. Thus the rate of dissipation is reduced and in terms of the expression of the evolution of granular temperature (4.11), the dissipation term is now given by
$\unicode[STIX]{x1D70C}_{s}f_{3}T^{3/2}/L_{c}$
. Jenkins (Reference Jenkins2007) suggested a phenomenological model for the length,
$L_{c}$
, in which its magnitude is proportional to the rate of compression that occurs along at least one axis in shear flows and inversely proportional to the agitation (the granular temperature) that can destroy these structures. Thus in dimensional form for simple shear flows
$\boldsymbol{v}=v(z)\hat{\boldsymbol{x}}$
, Jenkins & Berzi (Reference Jenkins and Berzi2010) propose
$$\begin{eqnarray}\frac{L_{c}}{d}=\max \left(1,\frac{{\hat{c}}(\unicode[STIX]{x1D719}g_{0})^{1/3}d}{2T^{1/2}}\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}z}\right),\end{eqnarray}$$
where
${\hat{c}}$
is a dimensionless constant of order unity (often
${\hat{c}}=1/2$
). Jenkins & Berzi (Reference Jenkins and Berzi2010) validate this formulation empirically for unfluidised granular flows. We are not aware of any studies that have tested formulae for fluidised flows, but we can nevertheless employ this formulation (A 1) to compute profiles of the volume fraction of particles, the velocity field and the granular temperature for typical parameter values used in this study (figure 28). For a dimensionless fluidising gas flow rate,
$W_{g}$
equal to
$10^{-3}$
and a slope
$S$
of
$0.1$
, we find negligible differences in the profiles apart from very close to the base of the flow. Moreover the dimensionless volume flux per unit width for the ‘standard’ kinetic theory
$\hat{q}=0.1117$
, while for the extended kinetic theory
$\hat{q}=0.1127$
.
For more weakly fluidised flows, there can be a significant difference between the predictions of the two theories, because in these situations the concentration of particles is higher and thus
$L_{c}/d$
exceeds unity in many parts of the flow. For example when
$W_{g}=4.1\times 10^{-4}$
and
$S=0.1$
, we find that extended kinetic theory predicts more energetic and faster-moving flows (see figure 29). For these parameter values, the dimensionless volume flux,
$\hat{q}=0.0198$
for the ‘standard’ kinetic theory, whereas
$\hat{q}=0.0284$
for the extended kinetic theory. The flows that we consider in this study are more strongly fluidised than this example and thus we find it unnecessary to include this phenomenon in our analysis in the main body of this paper because it introduces negligible difference to the computed flow.
Figure 29. The volume fraction,
$\unicode[STIX]{x1D719}(\hat{z})$
, velocity of the solid phase,
$\hat{v}(\hat{z})$
, and the granular temperature,
$\hat{T}(\hat{z})$
, as functions of the dimensionless depth within the current for parameter values
$R=10^{-3}$
,
$\unicode[STIX]{x1D713}=0.5$
,
$\unicode[STIX]{x1D719}_{m}=0.63$
,
$e=0.85$
,
$e_{w}=0.75$
,
$S=0.1$
,
$St=10^{3}\unicode[STIX]{x1D6FF}^{2}$
,
$\unicode[STIX]{x1D6FF}=0.01$
and
$W_{g}=4.1\times 10^{-4}$
for extended kinetic theory (solid lines) and ‘standard’ kinetic theory (dashed lines).
Appendix B. The effects of sidewall stresses
In this appendix we analyse the effects of sidewall resistance on the motion of shallow fluidised flows down inclined channels (see § 5) and derive the first-order correction to the prediction of the front speed for flows that are unaffected by sidewalls. We show that the reduction in front speed is proportional to
$(H/B)^{2}$
, where
$H$
is the scale depth of the current given by (5.10) and
$B$
is the channel breadth.
The downslope flows studied experimentally are realised within a channel, the width of which is usually greater than the flow depth, but not far in excess of the depth. Thus, it is feasible that sidewall stresses may play a non-negligible role in the overall dynamics and may further retard the motion. Indeed for flows along horizontal surface for which the flow depth is much greater than the width of the channel, we postulate that the sidewall stresses may even play a dominant role in the resisting the driving forces (see § 6).
We analyse the motion in a channel of width
$B$
in a regime for which the volume fraction is spatially uniform and given by
$\overline{\unicode[STIX]{x1D719}}$
. On depth and width averaging the streamwise balance of momentum (5.3), we find that the dimensional governing equation is given by
$$\begin{eqnarray}\displaystyle & & \displaystyle \unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}\left(\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}t}\int _{0}^{B}\int _{0}hv\,\text{d}z\,\text{d}y+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}\int _{0}^{B}\int _{0}^{h}v^{2}\,\text{d}z\,\text{d}y+g\cos \unicode[STIX]{x1D703}\int _{0}^{B}\int _{0}^{h}\frac{\unicode[STIX]{x2202}h}{\unicode[STIX]{x2202}x}\,\text{d}z\,\text{d}y\right)\nonumber\\ \displaystyle & & \displaystyle \quad =\unicode[STIX]{x1D70C}_{s}\overline{\unicode[STIX]{x1D719}}g\sin \unicode[STIX]{x1D703}Bh-\int _{0}^{B}\unicode[STIX]{x1D70E}_{xz}(0,y)\,\text{d}y-\int _{0}^{h}(\unicode[STIX]{x1D70E}_{xy}(z,0)-\unicode[STIX]{x1D70E}_{xy}(z,B))\,\text{d}z,\end{eqnarray}$$
where
$\unicode[STIX]{x1D70E}_{xz}(y,0)=f_{1}\unicode[STIX]{x1D70C}_{s}\text{d}T^{1/2}\unicode[STIX]{x2202}v/\unicode[STIX]{x2202}z$
denotes the basal shear stress and
$\unicode[STIX]{x1D70E}_{xy}(0,z)=-\unicode[STIX]{x1D70E}_{xy}(B,z)=f_{1}\unicode[STIX]{x1D70C}_{s}\text{d}T^{1/2}\unicode[STIX]{x2202}v/\unicode[STIX]{x2202}y$
denotes the sidewall stresses.
In general, to include these sidewall effects, even if the flow had adjusted to a local balance that was independent of the streamwise coordinate, we would have to resolve the variations of the dependent fields in the
$(y,z)$
plane (see, for example, Oger & Savage Reference Oger and Savage2013). In this subsection we take a different strategy and develop a model which is appropriate to the regime
$H/B\ll 1$
, where
$H$
is the scale depth of the flow and given by (5.10). In this regime, we treat the granular temperature and flow field as predominantly varying with the distance from the basal boundary and assume that these flow fields adopt the form established in § 4.1.5 (see Jenkins & Berzi Reference Jenkins and Berzi2010). This approach was used above to derive the depth-averaged model above, but is now generalised to include lateral gradients in order to model the sidewall stresses.
Adopting the dimensionless variables using the scales of (5.4), we estimate the velocity gradient at the sidewall
$\unicode[STIX]{x2202}\hat{v}/\unicode[STIX]{x2202}{\hat{y}}=\unicode[STIX]{x1D6FC}H\hat{v}/B$
, where
$\unicode[STIX]{x1D6FC}$
is a dimensionless constant of order unity. We may then compute the depth and width averages to deduce a governing equation for a travelling wave solution
${\hat{h}}(x,t)={\hat{h}}(x-ct)$
that features the additional stresses due to the sidewalls (cf. (5.12)); it is given by
$$\begin{eqnarray}\frac{{\mathcal{R}}c^{2}}{4}{\hat{h}}^{\prime }+{\hat{h}}{\hat{h}}^{\prime }={\hat{h}}-\frac{25}{4}\frac{c^{2}}{{\hat{h}}^{2}}-\frac{25\unicode[STIX]{x1D6FC}}{9}\left(\frac{H}{B}\right)^{2}c^{2}.\end{eqnarray}$$
The final term of (B 2) represents the extra stress due to the sidewalls. Then using the uniform conditions far from the front
$(\hat{x}-c\hat{t}\rightarrow -\infty )$
we deduce that
$$\begin{eqnarray}0=\frac{1}{c}-\frac{25c^{4}}{4}-\frac{25\unicode[STIX]{x1D6FC}H^{2}c^{2}}{9B^{2}}.\end{eqnarray}$$
Thus the speed,
$c$
, is reduced by the action of the sidewall stresses and in the regime
$\unicode[STIX]{x1D6FC}(H/B)^{2}\ll 1$
, we find that
$$\begin{eqnarray}c=\left(\frac{2}{5}\right)^{2/5}\left(1-\frac{4\unicode[STIX]{x1D6FC}}{45}\left(\frac{2}{5}\right)^{-4/5}\left(\frac{H}{B}\right)^{2}+\cdots \,\right)\end{eqnarray}$$
The value of the coefficient of the second term of the expansion is
${\sim}$
0.2.
