Time-dependent solutions for particle-size segregation in shallow granular avalanches

J.M.N.T Gray, M Shearer, A.R Thornton


Rapid shallow granular free-surface flows develop in a wide range of industrial and geophysical flows, ranging from rotating kilns and blenders to rock-falls, snow slab-avalanches and debris-flows. Within these flows, grains of different sizes often separate out into inversely graded layers, with the large particles on top of the fines, by a process called kinetic sieving. In this paper, a recent theory is used to construct exact time-dependent two-dimensional solutions for the development of the particle-size distribution in inclined chute flows. The first problem assumes the flow is initially homogeneously mixed and is fed at the inflow with homogeneous material of the same concentration. Concentration shocks develop during the flow and the particles eventually separate out into inversely graded layers sufficiently far downstream. Sections with a monotonically decreasing shock height, between these layers, steepen and break in finite time. The second problem assumes that the material is normally graded, with the small particles on top of the coarse ones. In this case, shock waves, concentration expansions, non-centred expanding shock regions and breaking shocks develop. As the parameters are varied, nonlinearity leads to fundamental topological changes in the solution, and, in simple-shear, a logarithmic singularity prevents a steady-state solution from being attained.


1. Introduction

Particle-size segregation is notorious in bulk solids handling. Sometimes it can be used to our advantage, such as in mineral separation technology (Wills 1979), but often it presents a major technical hurdle in manufacturing processes, where grains need to be mixed together, or transported to another location. The scale of granular materials processing is vast, with applications in the bulk chemical, pharmaceutical, mining, food and agricultural industries. In many of these processes, such as in rotational mixers (Shinbrot & Muzzio 2000) and huge continuous feed rotational kilns (Davidson et al. 2000), segregation occurs in shallow granular avalanches that develop in the free-surface layer of the flow. A knowledge of the segregation within these avalanches is, therefore, vital to understanding the segregation and mixing in more complex granular flows with solid–liquid phase transitions (Williams 1968; Gray & Hutter 1997).

There are a number of mechanisms that drive particle-size segregation, but the dominant one in shallow granular free-surface flows is that of kinetic sieving (Bridgwater 1976), with the next most important effects being diffusive remixing, particle-density differences and grain-inertia effects at size ratios greater than five (Thomas 2000). The basic kinetic sieving mechanism is simple. As the grains avalanche downslope there are fluctuations in the void space between the particles. When a void opens up under a layer of grains, the small particles are more likely to fall into the gap, because they are more likely to fit into the available space. The fines, therefore, percolate to the bottom of the flow and mass conservation dictates that there is a corresponding reverse flow of large particles towards the free-surface. Kinetic sieving is so efficient in dry granular flows that in small scale experiments a layer of 100% coarse grains develops on top of a layer of 100% fines with a sharp concentration jump between them (Savage & Lun 1988; Vallance & Savage 2000). In geology, this type of particle-size distribution is termed inverse-grading (Middleton & Hampton 1976). The local size distribution can also have a subtle feedback onto the bulk flow. For instance, when the large particles are rougher than the small ones and there is strong shear through the depth of the avalanche, the large particles tend to congregate at the front and resist the motion. This leads to an instability in which the large particles are pushed to the sides to form stationary lateral levees (Pouliquen et al. 1997) that channel the flow and lead to significantly longer debris-flow run-out distances (Vallance 2000; Iverson & Vallance 2001).

Despite the critical importance of segregation in both industrial and geophysical granular flows, there has been very little theoretical progress. Savage & Lun (1988) were the first to derive a steady-state model using statistical mechanics and information entropy ideas, but it was overlooked for many years because of its apparent complexity. It also had the disadvantage that the percolation fluxes were independent of gravity, even though gravity is the fundamental driving mechanism for the kinetic sieving process. Gray & Thornton (2005) used mixture theory to formulate mass and momentum balances for the large and small particles, which provided a natural way of introducing gravity into the theory. Assuming a linear velocity drag between the particles and that the fines carry less of the overburden pressure as they fall down through the gaps, Gray and Thornton derived an expression for the segregation flux in terms of the volume fraction of small particles. The resulting time-dependent segregation equation is considerably simpler and more general than that of Savage & Lun (1988), but is still able to quantitatively reproduce the concentration jumps observed in the laboratory experiments of Savage & Lun (1988) and Vallance & Savage (2000). This is no coincidence, since the underlying mathematical structure of the two theories is the same in steady-state, even though the derivation and the interpretation of the coefficients are quite different. Gray & Thornton's (2005) model represents the simplest possible theory for modelling three-dimensional time-dependent segregation, and its basic structure is inherent in more complex models that incorporate a dense-fluid in the pore space (Thornton et al. in press) and diffusive-remixing (Dolgunin & Ukolov 1995; Gray & Chugunov submitted).

In this paper, two exact solutions of the Gray and Thornton segregation equation are constructed; these yield detailed insight into the formation and evolution of some key segregation features, such as shock waves, expansion fans and inversely graded layers, that cannot be deduced easily from numerical simulations alone. One virtue of the two exact solutions is that they are a rigorous check that numerical solutions of the partial differential equation capture accurately details of the segregation as it evolves.

2. The segregation equations and the bulk flow

Gray & Thornton's (2005) theory introduces a volume fraction ϕ of small particles per unit mixture volume that lies in the rangeEmbedded Image(2.1)with ϕ=1 corresponding to 100% small particles and ϕ=0 corresponding to 100% coarse particles. The large particles are of the same density and occupy a volume fraction 1−ϕ. This simple approach implicitly assumes that the interstitial pore space is incorporated into the bulk density of the large and small particles. A more complex theory, which accounts for the buoyancy effects of an interstitial fluid of a different density to the particles (Thornton et al. in press) leads to a similar formulation. In this section, we outline the derivation of the segregation equation, referring to Gray & Thornton (2005) for the details.

Within the avalanche, the bulk pressure p is assumed to be hydrostatic through the depth h, which is measured normal to the chute, i.e.Embedded Image(2.2)where ρ is the bulk density, g is the constant of gravitational acceleration and ζ is the inclination angle of the chute. The bulk flow is assumed to be incompressible,Embedded Image(2.3)where x, y, z are coordinates in the downslope, cross-slope and normal directions to the chute, with unit coordinate vectors i, j, k, respectively, and velocity Embedded Image. In shallow granular free-surface flow models, the bulk depth-averaged downslope and cross-slope velocities Embedded Image are computed as a function of space and time, using prescribed velocity profiles. These profiles introduce shape factors into the momentum transport terms, which are equal to unity for plug flow and 4/3 for simple shear. The bulk three-dimensional velocity field u(x, t) can be reconstructed from Embedded Image with a knowledge of the shape factors together with the incompressibility relation and the no normal velocity condition at the base to determine the bulk normal velocity,Embedded Image(2.4)Assumptions (2.2) and (2.3) are common to almost all models for granular free-surface flows from the hydraulic-type avalanches theories (e.g. Grigorian et al. 1967; Eglit 1983; Gray et al. 2003), to the dry Mohr–Coulomb (e.g. Savage & Hutter 1989; Gray et al. 1999) and water saturated debris-flow models (Iverson 1997). All of these theories can be used to compute the bulk velocity field, but, for simplicity it is prescribed in this paper.

The key part of Gray & Thornton's (2005) model are the equations for the percolation of the large and small particles. These are derived from the normal component of the constituent momentum balances by assuming slow flow, a linear velocity dependent drag and a nonlinear pressure sharing law. This is based on the observation that the large particles carry more of the overburden pressure while the small particles are falling through the interstitial gaps. Full details of the derivation are given in Gray & Thornton (2005) and the resulting non-dimensional normal percolation velocities wl and ws of the large and small particles areEmbedded Image(2.5)Embedded Image(2.6)The first equation implies that the large particles move upwards relative to the bulk, until there are no more small particles left. The second equation shows that the fines percolate down until there are no more big particles. As such it represents the simplest possible model for segregation. The non-dimensional segregation number Sr is defined asEmbedded Image(2.7)where c is an inter-particle drag coefficient, the non-dimensional factor B determines the magnitude of the pressure perturbations that drive the flow and the avalanche has thickness H, length L and downslope velocity magnitude U. In general, the non-dimensional segregation number Sr may depend on a number of additional factors such as the local shear, the particle-size difference and the particle roughness. In this paper, Sr is assumed to be a constant that can be calibrated from experiment.

The large particle percolation equation (2.5) simply implies that the large particles move upwards, relative to the bulk, until they separate out into a pure coarse phase. Equation (2.6) similarly expresses the property that small particles filter down until they separate out. Substituting the small particle percolation velocity (2.6) into the small particle mass balance equationEmbedded Image(2.8)together with the assumption that the downslope and cross-slope velocity components of the small particles are the same as those of the bulk flowEmbedded Image(2.9)yields a non-dimensional segregation equation for the volume fraction of small particlesEmbedded Image(2.10)The first four terms together with the bulk incompressibility (2.3) simply express the fact that the volume fraction is advected passively with the flow. The final term on the left-hand side is responsible for the separation of the particles into inversely graded layers.

Gray & Thornton (2005) and Thornton et al. (in press) showed that sharp concentration jumps or shocks appear in steady-state flows. A generalized time-dependent jump condition across such a shock can easily be derived from the divergence structure of (2.10). Let z=z(x, y, t) be a shock surface, with space–time normal N=(zt, zx, zy,−1). Observe that equation (2.10) has the divergence form Div G=0, where Div is the space–time divergence andEmbedded Image(2.11)Consequently, the normal component of G is continuous across the shock: N. Embedded Image, where the jump bracket Embedded Image denotes the jump in a piecewise continuous function f(x, y, z, t) across the shock. Thus,Embedded Image(2.12)Dividing by Embedded Image, and assuming that the velocity (u, v, w) is continuous, we obtain the equationEmbedded Image(2.13)In the exact solutions derived in this paper, each shock z=zγ(x, y, t) is given an integer subscript γ=1, 2, 3, … to uniquely identify it.

The upper and lower boundaries of the avalanche, the surface z=s(x, y, t) and the bottom z=b(x, y), are special cases of the jump conditions in which there is no normal flow of small particles across the boundary. It follows from (2.12) that the surface and basal boundary conditions areEmbedded Image(2.14)Embedded Image(2.15)The square-bracketed term on the left-hand side of these equations is zero, due to the kinematic boundary condition at the surface and base of the avalanche (e.g. Gray et al. 1999, 2003). The surface and basal boundary conditions, therefore, reduce toEmbedded Image(2.16)which is satisfied when either ϕ=0 or 1 at the surface and base.

3. Steady-uniform flows

In this paper, we restrict attention to steady uniform flow, meaning granular flow in which the layer has uniform thickness, the bulk downslope velocity depends only on the depth variable z, and the cross-slope and normal velocity components are assumed to be zero. In these circumstances, the avalanche depth sets the length scale in z, so that the avalanche has unit depth in non-dimensional variables and the velocity field isEmbedded Image(3.1)The stretching transformationEmbedded Image(3.2)can be used to rescale the segregation equation (2.10) into a convenient parameter independent formEmbedded Image(3.3)where the tildes are dropped for simplicity, and the associated jump condition isEmbedded Image(3.4)The transformation (3.2) implies that the grains take greater distances to segregate into inversely graded layers and take longer to do so for smaller Sr. All solutions are constructed in the stretched coordinate system and mapped back to the physical coordinate system using the transformation (3.2) if and when necessary.

4. Segregation in a homogeneous chute flow

For our first problem, we will consider the time-dependent version of Gray & Thornton's (2005) steady-state problem in which a steady-uniform avalanche that is initially homogeneously mixed with concentration ϕ0 is continuously fed at x=0 with material of the same concentration ϕ0. The problem is, therefore, specified by the initial and boundary conditionsEmbedded Image(4.1)Embedded Image(4.2)Embedded Image(4.3)where (4.3) is the no normal flux condition at the boundaries. Within the interior of the flow the small particles immediately start to percolate down towards the bottom of the avalanche and the large ones are pushed towards the free-surface. Except near boundaries, the characteristics imply that the small particle concentration stays at ϕ0 and is swept downstream. At the lower boundary, however, there are simply no more large particles to be pushed up and the no normal flux condition (4.3) implies that the small particles separate out into a pure phase, creating a sharp concentration jump, or shock, with the initial mixture. Assuming that ϕ+=ϕ0 and ϕ=1 on either side of the shock, (3.4) with γ=1 reduces toEmbedded Image(4.4)The method of characteristics (e.g. Sneddon 1957) expresses differentiable solutions in the form F(λ, μ)=0, where F is an arbitrary function of the characteristic variables λ(x, z, t) and μ(x, z, t), which are integrals of the characteristic equationsEmbedded Image(4.5)We find expressions for λ and μ by solving pairs of these equations. Integrating the first and third equations, we obtainEmbedded Image(4.6)The second and third equations can be integrated for general velocity fields by introducing a depth-integrated velocity coordinate (Gray & Thornton 2005),Embedded Image(4.7)to giveEmbedded Image(4.8)(Here and elsewhere, we abbreviate functional dependences; thus, ψ1=ψ(z1), and equation (4.8) may be regarded as relating z1 to x along a characteristic.) Note that the depth-integrated velocity coordinate lies in the range 0≤ψ≤1 provided the downslope velocity is suitably scaled in the non-dimensionalization process.

The full solution surface can be expressed as λ=f(μ), where f is an arbitrary function. In this case, the solution is particularly simple. Since (4.6) is independent of x and the initial shock height z1(x, 0)=0, it follows that the shock is horizontal and moves upwards linearly in time:Embedded Image(4.9)for the material that was initially in the chute. For the material that enters the chute, the shock emanates from the origin, ψ1=0, so (4.8) implies that the solution is independent of time:Embedded Image(4.10)The solution, therefore, consists of a time-dependent horizontal section that propagates upwards linearly in time and is swept downstream to reveal the steady-state solution adjacent to the inflow.

At the surface of the avalanche a similar jump develops. This time there are no more small particles available to percolate downwards and the large particles separate out to form a concentration shock with the bulk mixture. Assuming that ϕ+=0 and ϕ=ϕ0 the jump condition (3.4) impliesEmbedded Image(4.11)This can be solved in the same way as for the bottom jump. It consists of a horizontal part that propagates linearly downwards with time,Embedded Image(4.12)and is swept downstream to reveal the steady-state part,Embedded Image(4.13)

The concentration jump solutions (4.9), (4.10), (4.12) and (4.13) hold for any positive velocity field u(z) with an isolated zero at the origin. Specific solutions are constructed in this paper for the linear velocity field,Embedded Image(4.14)The parameter α allows the velocity to vary from plug-flow (α=1) to simple shear (α=0) and provides a good leading order approximation to more complex nonlinear velocity fields (e.g. Gray & Thornton 2005). The corresponding depth integrated coordinate,Embedded Image(4.15)can be inverted to give z=z(ψ), i.e.Embedded Image(4.16)The development of the upper and lower shocks is shown in figure 1 for a linear velocity field with α=0.1. The position of the transition points between the steady-state and horizontal time-dependent sections of the lower and upper shocks are (respectively)Embedded Image(4.17)Embedded Image(4.18)and are shown by the open circular markers in the left-hand panels of figure 1. In plug-flow, the two transitions propagate at the same speed, but when there is shear the upper transition x2 moves faster than the lower transition x1 in response to the larger velocities at the free-surface; thus, x1<x2 for t>0.

Figure 1

The evolution of the volume fraction of small particles ϕ as a function of x and z. The particles enter the chute at x=0 and are swept downstream from left to right by a linear velocity field with a large shear component given by parameter α=0.1. The white region corresponds to 100% small particles, the dark region to 100% large particles and the light grey region to the mixed region with concentration ϕ0=0.2. The circular markers correspond to the position of the transition points and the ‘⊕’ is the triple-point. An mpeg movie of the solution is included in the electronic supplementary material together with the matlab-code to make movies for general α and ϕ0.

The upper and lower horizontal shocks meet where z1=z2. From (4.9) and (4.12), this occurs at t=1, at height z=ϕ0, and over the range x>xm, where xm is given by (4.18) at t=1:Embedded Image(4.19)For t>1, the resulting shock with ϕ+=0 and ϕ=1 evolves according to the jump condition (3.4), which reduces to the scalar conservation law (resembling the inviscid Burgers equation)Embedded Image(4.20)Given the initial position and time (x0, z0, t0) this has the simple solutionEmbedded Image(4.21)It follows that when the upper and lower laterally uniform (horizontal) time-dependent shocks meet, a third laterally uniform shock is formed that is stationary and advected downstream:Embedded Image(4.22)

For plug flow (α=1), the solution is now at equilibrium, but for shear flow α<1, there is subsequent dynamic behaviour, so we now discuss only the case α<1. In the region x<xm, the upper and lower shocks do not intersect. Instead, the lower horizontal shock (4.13) continues to propagate upwards intersecting the steady branch (4.8) of the top shock at a triple-point between the large, small and mixed regions:Embedded Image(4.23)This construction continues until the two steady-state branches meet at timeEmbedded Image(4.24)at the steady-state triple-point positionEmbedded Image(4.25)The triple-point is illustrated in figure 1 using a ‘⊕’ marker. It moves upstream from the initial point of meeting xm to its final steady-state position xsteady. The time it takes to reach this point is plotted in figure 2, which shows that it takes longer to develop with decreasing α and ϕ0.

Figure 2

A contour plot of the time tsteady for the two steady-state branches of the shock to meet at the steady-state triple-point as a function of velocity profile α and the inflow concentration ϕ0.

The solution upstream of the steady triple-point has reached equilibrium, but downstream, the interface between the large and small particles continues to evolve. There are three distinct sections whose dynamics are all controlled by the conservation law equation (4.20). As can be seen from the general solution (4.21), the initial position and time (x0, z0, t0) at which the shock is generated are of crucial importance. It has already been shown that the shock is horizontal with height ϕ0 in a region propagating downstream. The remaining parts of the shock are governed by the upstream motion of the triple-point as it moves to its steady position. At a given height z∈[zm,ztriple] equations (4.23) imply that its initial position isEmbedded Image(4.26)at time t0=z/ϕ0. The dynamic part of the shock is, therefore, given byEmbedded Image(4.27)Once the triple-point reaches its steady-state position, the initial conditions do not change and a straight horizontal shock of heightEmbedded Image(4.28)is formed. The development of the three sections of the third shock is shown in the right-hand panels of figure 1. When there is shear, points that are higher up in the flow move faster downstream than those below them. As a result, the central dynamic region of z3 steepens and eventually breaks at z=zm (where the interface is steepest), at timeEmbedded Image(4.29)This indicates that the breaking-time is singular when α=1 and the third shock does not break, consistent with the exact time-dependent solutions for plug-flow constructed by Gray & Thornton (2005). The breaking-time also tends to infinity as the initial concentration approaches either zero or unity. This is because the dynamic region of z3 is confined to a very narrow lateral band in which the velocities become increasingly close to one another as ϕ0 approaches these limiting states and it therefore takes an increasingly long time to break. The full parameter dependence of the breaking-time on α and ϕ0 is illustrated in figure 3. This shows that there is a small (shaded) region of parameter space for high shear and low initial concentrations in which the third shock breaks prior to the steady triple-point forming. Numerical solutions indicate that after the shock breaks a complex solution is formed that includes dynamic expansion fans and shocks, but a detailed analysis of this is left to a future paper.

Figure 3

A contour plot of the time tbreak for the third shock to break as a function of velocity profile α and the inflow concentration ϕ0. The shaded region shows the region of parameter space where the shock breaks before the steady-state triple point forms.

5. Segregation from a normally graded initial state

The second problem also leads to shock wave-breaking and is the time-dependent version of Thornton et al. (in press) normally graded steady-state flows. Initially, the avalanche is assumed to be normally graded, with a layer of small particles separated from a layer of large particles below by a concentration discontinuity along z=zr. The avalanche is fed at x=0 with the same normally graded layered material and satisfies the no normal flux conditions at the surface and the base of the avalanche, i.e.Embedded Image(5.1)Embedded Image(5.2)Embedded Image(5.3)Thornton et al. (in press) have performed numerical simulations of this problem, which indicate that there is a region adjacent to the inflow that rapidly attains steady-state, a downstream region where the flow is time-dependent and laterally uniform and a transition region between the two, where there is complex spatial and temporal behaviour. The problem will, therefore, be solved in the steady-state and laterally uniform regions before the transition is considered.

(a) The steady problem

Neglecting the time-derivatives in the segregation equation (3.3) and expanding the remaining terms yields a first-order quasi-linear equation of the formEmbedded Image(5.4)which may be solved by the method of characteristics. The small particle concentration ϕ is equal to a constant ϕλ along the characteristic curve given by the subsidiary equationEmbedded Image(5.5)This may be integrated for general positive definite velocity fields with an isolated zero at the origin, by using the depth-integrated velocity coordinates defined in (4.7) to giveEmbedded Image(5.6)where (xλ, ψλ) is the starting position in mapped coordinates. Characteristic curves on which ϕλ>1/2 propagate up through the avalanche, while curves with ϕλ<1/2 propagate down towards the base of the flow. It follows that at the inflow concentration discontinuity, which lies at z=zr, or, ψ=ψr in mapped coordinates, the characteristics diverge from one another and a centred expansion fan is formed with concentrationEmbedded Image(5.7)The characteristic curve ψ=ψrx (on which ϕ=0) lies along the lower boundary of the expansion fan and separates it from a pure region of large particles. It represents the trajectory of the first small particles percolating towards the base of the avalanche, which they reach at xbottom=ψr. Similarly, the curve ψ=ψr+x (on which ϕ=1) separates the fan from a region of small particles above and represents the front of the first large particles that are being pushed up to the free-surface. They reach z=1 at xtop=1−ψr.

Once the small particles reach the base of the flow the no normal flux condition (5.3) implies that they separate out and a concentration shock is formed between the expansion fan and the pure phase. The steady-state version of the jump condition (3.4) with ϕ+ given by (5.7) and ϕ=1 reduces toEmbedded Image(5.8)in mapped coordinates. This is a linear ODE that can easily be integrated subject to the initial condition that ψ=0 at x=ψr to give the height of the bottom shockEmbedded Image(5.9)An analogous upper shockEmbedded Image(5.10)forms as the large particles separate out into a pure phase. The upper and lower shocks meet atEmbedded Image(5.11)and a third shock is formed between the inversely graded layers of large and small particles along the lineEmbedded Image(5.12)This solution has been investigated by Thornton et al. (in press) for a range of linear and nonlinear velocity profiles and values of the concentration discontinuity zr. Readers are referred to this work for a detailed discussion of the steady solution.

(b) The unsteady laterally uniform problem

Sufficiently far downstream, the lateral uniformity of the initial conditions implies that the solution is independent of the downslope coordinate x. In this case, the expanded segregation equation (3.3) reduces toEmbedded Image(5.13)which, rather interestingly, is independent of the assumed velocity profile u(z). Equation (5.13) has exactly the same form as the mapped steady-state problem (5.4) except that x is replaced by t and ψ is replaced by z. The full solution may therefore be written down immediately. It consists of an expansion fan centred at z=zr within which the concentration isEmbedded Image(5.14)

For small values of t, the expansion fan is confined to the interval |zzr|<t, the leading edges (where ϕ=0, 1) moving with speed |2ϕ−1|=1. After the leading edges hit the horizontal boundaries z=0, 1, at times tbottom=zr and ttop=1−zr, the fan is contained in the interval z1(t)<z<z2(t) bounded by the two shocksEmbedded Image(5.15)Embedded Image(5.16)Finally, when the two shocks collide, at time Embedded Image, the fan disappears and the solution is a single stationary shockEmbedded Image(5.17)

(c) The full problem and the fan-interface

Expanding the derivatives in the segregation equation (3.3) yields a first-order quasi-linear equation,Embedded Image(5.18)that governs the full problem. The method of characteristics implies that the concentration is equal to a constant ϕλ on characteristic curves (parameterized by λ) given by the subsidiary equations,Embedded Image(5.19)The first and third can be integrated directly to giveEmbedded Image(5.20)The first and second can then be integrated using (5.20) and the definition of the linear velocity field (4.14) to show thatEmbedded Image(5.21)where Embedded Image.

At t=0, the line x=0 separates the material that was initially in the chute from the material that subsequently enters the chute along x=0. Material initially in the chute is swept downstream, and experiences segregation in a time-dependent expansion fan (5.14). On the other hand, material entering the chute is normally graded as it enters, and then begins to segregate in a steady expansion fan (5.7) centred at x=0, z=zr. The interface between the two kinds of material points, those initially in the chute and those entering the chute, is a curve that we call the fan-interface; it is shown in the left panels of figure 4. Points in the interface are located by tracking the characteristics (5.21) with initial location x=0. Thus, f(λ)≡0. Within the time-dependent characteristic fan, ϕλ is given by (5.14), and the material on the time-dependent side of the fan-interface, therefore, lies along the straight line,Embedded Image(5.22)at time t. Equating the time-dependent and steady-state fan concentrations (5.14) and (5.7) also yields (5.22). This shows that the fan-interface (5.22) creates a continuous match between the steady-state and time-dependent fans. The solution for small times, therefore, consists of steady-state fan (5.7) upstream of the interface (5.22) and a time-dependent fan (5.14) downstream of it, which expand with increasing time. Above and below the fan are pure regions of small and large particles, as shown in the left-hand panels of figure 4. The fans continue to expand until the laterally uniform time-dependent sections intersect with the surface and base of the flow at times ttop and tbottom. Interaction with the no normal flux boundary condition (5.3) implies that shocks result from reflection at the boundary.

Figure 4

A series of contour plots showing the evolution of the volume fraction of small particles ϕ as a function of the downslope coordinate x and avalanche depth z. The bulk linear shear flow is from left to right and is given by parameter α=1/2. The normally graded particles enter the chute at x=0 and have a concentration discontinuity at zr=0.4. The white circular markers correspond to the position of the transition points between steady, fully dynamic and unsteady laterally uniform behaviour and the marker with the plus sign is the triple-point. An mpeg movie of the solution is included in the electronic supplementary material together with the matlab-code to make movies for general α and zr.

(d) The fully dynamic bottom shock

The subsequent flow involves horizontal shocks (discontinuities in ϕ) propagating in from the surface and the base. The structure of the solutions is shown schematically in figure 5. As indicated in the figure, the horizontal shocks do not match up neatly with the steady-state fan. As discussed earlier, the flow downstream consists of horizontal shocks, shown as sections CD (propagating upwards) and GH (propagating down) in figure 5. Note that the structure of the solution propagating up from the base is significantly different from that of the waves descending from the surface, due to the opposite velocity gradients encountered by the waves. In this section, we focus on the lower waves, the fully dynamic bottom shock.

Figure 5

A sketch of the full solution once the upper and lower shocks have formed.

The lower shock has three sections: the steady-state part AB (as in the steady state solution), the time-dependent laterally uniform section CD, and a fully dynamic shock BC joining the two simpler shocks. The curve BC is governed by the shock condition,Embedded Image(5.23)Using the method of characteristics to solve (5.23) yields two subsidiary equations,Embedded Image(5.24)The second of these can be integrated to giveEmbedded Image(5.25)Using the linear velocity profile (4.14) and the inverse mapping (4.16), the first equation reduces toEmbedded Image(5.26)which can be integrated by substituting for ψ from (5.25), making the substitution Embedded Image to reduce the integral to standard form and hence obtainEmbedded Image(5.27)whereEmbedded Image(5.28)with Embedded Image, b=4(1−α) and c=4(1−α1. These characteristic curves are sketched in figure 5 using dashed arrowed lines; they emanate from the dot–dash curve AC, which is defined by the trajectory of the intersection between the fan-interface CF and the laterally uniform shock CD. The location of points on AC is thus given byEmbedded Image(5.29)Each characteristic is initiated at some time to from (xb(to), zb(to)), which is shown by the grey markers on the dot–dash section of AC. On this characteristic the initial conditions imply that the constants in (5.25), (5.27) are given byEmbedded Image(5.30)

The fully dynamic part of the shock BC is therefore constructed by parameterizing the characteristics that are initiated at time to∈[tbottom,t] and working out their current position at time t. Equation (5.27) and the constants (5.30) imply that at any given time t the downslope position x1(t, to) of the bottom shock that was initiated at time to, is given by solving for x1∈[xb(to),xb(t)] from the equationEmbedded Image(5.31)The corresponding height z1(t, to) of the shock is then given implicitly byEmbedded Image(5.32)While we cannot write explicit formulae for (x1, z1)(t, to), the curve BC is easily generated numerically from the sequence of equations (5.28)–(5.32). The results showing the development of the steady-state, fully dynamic and time-dependent laterally uniform sections of the lower shock appear in the right-hand panels of figure 4 for α=1/2 and zr=0.4.

In simple shear the behaviour of the characteristic emanating from point A at (xbottom, 0) is particularly interesting. As α→0, the mapped discontinuity height ψrzr2 and xb(tbottom)→zr2, which imply thatEmbedded Image(5.33)in (5.28). Since λ1(xb(tbottom))=−2zr, the constant c in (5.28) is negative andEmbedded Image(5.34)It follows that (5.31) has a logarithmic singularity as α→0, and, therefore, point B cannot move away from its initial position in simple-shear. As a consequence, it is not possible to generate a steady-state section of the lower shock, when α=0, which is a rather surprising result.

(e) The fully dynamic top shock

The structure of the top shock is slightly different to that of the bottom shock, as suggested in figure 5. There is a steady-state section EF, a fully dynamic curve FG and a laterally uniform time-dependent section GH, but this time the contact fan-interface intersects the steady-state shock at F rather than along the laterally uniform part.

The fully dynamic curve FG is governed by the shock condition,Embedded Image(5.35)since ϕ=1 above the shock, and ϕ is given by (5.14) below. The method of characteristics yields two subsidiary equations,Embedded Image(5.36)which are easily solved to show that the characteristic curves are described byEmbedded Image(5.37)whereEmbedded Image(5.38)The shocks are initiated along the trajectory of point F, which moves along the steady-state section of the upper shock EF as shown in figure 5. At any given time t the position of F is found by solving the simultaneous equations (5.10), (5.22). This is achieved in practice by solvingEmbedded Image(5.39)for xs, where the function z=z(ψ) is given by (4.16) and ψ2 is the steady-state solution (5.10). The expanding fan of characteristic curves emanating from EF can then be parameterized by xo∈[xtop,xs(t)]. For a given characteristic xo, the initial height and time of creation areEmbedded Image(5.40)which imply that the constants in (5.37) are given byEmbedded Image(5.41)The current position of the fully dynamic part of the top shock is, therefore,Embedded Image(5.42)The evolution of the three sections of the top shock are illustrated in the right-hand panels of figure 4 for α=1/2 and zr=0.4.

(f) Interaction of the dynamic top shock and the time-dependent bottom shock

At time t3 the upper and lower time-dependent laterally uniform shocks meet to form a third laterally uniform shock z3 given by (5.17) over their common region of intersection. The top-left panel of figure 6, continuing on from figure 4, shows the solution for α=1/2 and zr=0.4. The position of the triple-point between the three shocks (shown with a ‘⊕’ in figure 6) can be found by iterating to find the characteristic xotriple∈[xtop,xs(t)] whose corresponding shock position, given by (5.42), has the same height as the bottom shock (5.15) usingEmbedded Image(5.43)and then substituting the result back into (5.42). The fully dynamic part of the upper shock generated by characteristics xo∈[xotriple,xs(t)] continues to evolve as before. However, the section generated by xo∈[xtop,xotriple] has already intersected with the bottom shock at timeEmbedded Image(5.44)and at positionEmbedded Image(5.45)Embedded Image(5.46)At these points, the shock is a jump from ϕ=0 to 1, and consequently now form part of the third shock z3 which is governed by the conservation law (4.20), and has the simple solutionEmbedded Image(5.47)As can be seen from the first three left-hand panels of figure 6 the triple-point moves upstream from its initial position as it generates the non-uniform section of the third shock, which is swept downstream.

Figure 6

A series of contour plots showing the evolution of the volume fraction of small particles ϕ as a function of the downslope coordinate x and avalanche depth z for the same case as in figure 4 at later times. The bulk linear shear flow is from left to right and is given by parameter α=1/2. The normally graded particles enter the chute at x=0 and have a concentration discontinuity at zr=0.4. The white circular markers correspond to the position of the transition points between steady, fully dynamic and unsteady laterally uniform behaviour and the marker with the plus sign is the triple-point.

(g) Destruction of the fan-interface

At t=tfan, the steady-upper and dynamic-lower shocks come together for the first time, eliminating the time-dependent fan (5.14), destroying the fan-interface (5.22) and ceasing the interaction between the unsteady-uniform bottom shock and the top-dynamic shock. The fan-interface destruction time tfan, therefore, marks a fundamental switch in behaviour of the solution. It is solved by an iterative procedure. For an initial guess for Embedded Image, the position of the intersection between the top-steady shock (5.10) and the bottom unsteady laterally uniform shock (5.32) is found by calculatingEmbedded Image(5.48)and then iterating tfan until it lies on the fan-interface (5.22) and satisfiesEmbedded Image(5.49)The upper bound on the interval for tfan follows from the condition that the time-dependent laterally uniform bottom shock must intersect with the steady top shock before z=1. The shock configuration at tfan=2.14 is shown in figure 6 and the complete parameter dependence of the solution is plotted in figure 7. This shows that for fixed α, the fan-destruction time has a local maximum close to zr=1/2, while for fixed zr the time tfan increases with increasing shear. In plug-flow (α=1) the fan-destruction time tfan=t3, which is the same time that the upper and lower laterally uniform shocks initially meet.

Figure 7

A contour plot of the fan-interface destruction time tfan as a function of velocity profile α and the inflow discontinuity height zr.

(h) Interaction of the dynamic bottom shock and the steady top shock

For the vast majority of parameter values the dynamic bottom shock and the steady top shock first intersect for t>tfan to create another section of the third shock z3. The resulting triple-point is found iteratively by searching for the characteristic totriple∈[tbottom,tcontact], whose current bottom shock position x1(t,totriple) intersects with the steady top shock (5.10), i.e. whenEmbedded Image(5.50)The triple-point coordinates are then given by (5.31) and (5.32). The fully dynamic part of the bottom shock generated by the characteristics to∈[tbottom,totriple] continues to evolve as before. But, the part generated by the characteristics to∈[totriple,tfan] intersected with the steady top shock at pointEmbedded Image(5.51)at timeEmbedded Image(5.52)These points now generate a new section of z3. This is governed by the inviscid Burgers equation (4.20) and has the simple solutionEmbedded Image(5.53)where zo1=z(ψo1). The final point of this unsteady third shock is generated as point B reaches the steady-state triple-point (5.11) at (xsteady,ψsteady). This point lies on the characteristic with label tbottom and (5.31) implies that it reaches xsteady at timeEmbedded Image(5.54)This time is important because it marks the point at which the steady-state concentration fan becomes fully developed. The parameter dependence of tsteady on α and zr is illustrated in figure 8. In simple shear point B does not move from its initial position because of the logarithmic singularity (5.34). As a result the steady-state fan does not develop and there is a logarithmic singularity in tsteady along α=0. For plug-flow, the steady-state fan becomes fully developed at the same time the upper and lower laterally uniform shocks meet, i.e. tsteady=t3 when α=1. Between these two extremes there is a maximum in tsteady, for each α, that becomes progressively skewed to higher discontinuity heights zr as α decreases.

Figure 8

A contour plot of the time tsteady for the steady-state fan to fully develop as a function of velocity profile α and the inflow discontinuity height zr.

As we warned at the beginning of section, the bottom dynamic shock does not always intersect with the top shock after tfan. For a small range of parameter values in the bottom left-hand corner of figure 8, i.e. for α and zr small, the dynamic bottom shock develops a pronounced peak as shown in figure 9. This can intersect with the steady top shock prior to tfan, dividing the fan in two, and changing the topology of the problem. Numerical simulations of a similar cut-off phenomenon by Thornton et al. (in press) suggest a detached section of z3 forms, which expands and merges with the rest of z3 as the enclosed part of the fan is destroyed.

Figure 9

Two contour plots showing a zoomed in view of the small particle concentration for α=0.01 and zr=0.16. At high shear and low zr the dynamic bottom shock develops a peak (a) that intersects (b) with the top steady shock before t=tfan. Note, the vertical and horizontal scales have changed.

For the majority of parameter values, however, a laterally uniform steady-state shock ψ3, that was defined in (5.12), is generated for t>tsteady. The full structure and development of z3 is illustrated in the remaining panels of figure 6. It consists of four sections. A laterally uniform region at height 1−zr generated by the initial data, a section formed from the interaction of the dynamic top shock and the laterally uniform bottom shock (5.47), a section generated by the interaction of the steady-state top shock and the dynamic bottom shock (5.53), and, a steady-state region at height 1−ψr. The markers in figure 6 indicate the transitions between these regions. As in the first problem in §4 the fully dynamic parts of the third shock, z3, steepen and break in finite time tbreak, which may be either before or after tsteady. Although the subsequent breaking is complex, this region is swept downstream, and locally a steady-state is rapidly attained provided the logarithmic singularity in simple-shear is avoided. This is illustrated in the last panel of figure 6, where the local solution in the domain of interest has reached steady-state at t=4. In this fully developed region, large particles enter from the lower left and move up through the expansion fan, before separating out, while the small particles enter from the upper left, move down through the fan and separate into a pure phase. In this way, an inversely graded layer is formed, which is the basic building block for many pattern formation processes.

6. Discussion and conclusions

In this paper, two exact solutions to the segregation theory of Gray & Thornton (2005) have been derived for steady uniform chute flows. These graphically illustrate some of the key physical effects that develop during the segregation process, including: concentration shocks, concentration expansions, non-centred expanding shock regions, the creation of inversely graded layers and the steepening and ultimate breaking of monotonely decreasing sections of the shock between these layers. The non-dimensional segregation number Sr determines the strength of the segregation and may be calibrated from experiment. In steady-uniform flows, a spatial and temporal stretching transformation (3.2) can be used to reduce the segregation equation to a convenient parameterless form. The resulting solutions are dependent only on the assumed velocity profile with depth and the boundary and initial conditions. By constructing exact solutions we have, therefore, been able to map out the entire parameter space for the two problems. This not only provides important test cases for numerical shock-capturing methods to solve more general problems, but also yields implicit and explicit formulae for the key times and locations at which the behaviour of the solution fundamentally changes. Of particular interest is the fact that in simple shear the logarithmic singularity in (5.34) prevents the formation of a steady-state lower shock. This is qualitatively different to solutions with slip at the base, where a steady-state lower shock is formed in finite time. In addition, we have also identified a region of parameter space, where the nonlinearity leads to the pinching-off of a section of the expansion fan, changing the topology of the solution. Regardless of the precise nature of the initial and boundary conditions, or the topology changes that occur during the flow, particles separate out into inversely graded layers sufficiently far downstream and the complex shock breaking behaviour in monotonely decreasing regions is advected downstream.


The authors gratefully acknowledge the generous support of the Isaac Newton Institute Programme on ‘Granular and particle-laden flows’, which was held in Cambridge from September–December 2003. Nico Gray was also supported by an EPSRC Advanced Research Fellowship (GR/S50052/01 and GR/50069/01) and Anthony Thornton by an EPSRC Doctoral Training Account. Michael Shearer was supported in part by Grant DMS-0244491 from the National Science Foundation.



View Abstract