2.1. Governing equations

Consider a mixture of large and small particles whose pore space is occupied by an interstitial fluid. Each of the constituents occupies a volume fraction $\varPhi ^{\nu }\in [0,1]$ per unit mixture volume, where the constituent letters $\nu =l$, $s$ and $f$, refer to large particles, small particles and fluid, respectively. It is also possible to define volume fractions of large and small particles per unit granular volume

(2.1a,b)\begin{equation} \phi^{l}=\frac{\varPhi^{l}}{\varPhi^{l}+\varPhi^{s}},\quad \phi^{s}=\frac{\varPhi^{s}}{\varPhi^{l}+\varPhi^{s}}. \end{equation}

This is useful, because in many situations of practical interest the solids volume fraction $\varPhi =\varPhi ^{l}+\varPhi ^{s}$ does not change very much (Silbert et al. Reference Silbert, Ertas, Grest, Halsey, Levine and Plimpton2001; GDR-MiDi 2004) and can therefore be scaled out of the equations (Thornton et al. Reference Thornton, Gray and Hogg2006). In this situation, the general form of the bidisperse segregation equations for the volume fractions (concentrations) of the large and small particles (see e.g. Gray Reference Gray2018) are

(2.2) \begin{gather} \frac{\partial{\phi^{l}}}{\partial{t}} + \boldsymbol{\nabla}\boldsymbol{\cdot} (\phi^{l} \boldsymbol{u}) - \boldsymbol{\nabla}\boldsymbol{\cdot} \left(\, f_{sl}\phi^{l}\phi^{s} \frac{\boldsymbol{g}}{|\boldsymbol{g}|}\right) = \boldsymbol{\nabla} \boldsymbol{\cdot} (\mathcal{D}_{sl}\boldsymbol{\nabla}\phi^{l}), \end{gather}

(2.3) \begin{gather} \frac{\partial{\phi^{s}}}{\partial{t}} + \boldsymbol{\nabla}\boldsymbol{\cdot} (\phi^{s} \boldsymbol{u}) + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\, f_{sl}\phi^{s}\phi^{l} \frac{\boldsymbol{g}}{|\boldsymbol{g}|}\right) = \boldsymbol{\nabla} \boldsymbol{\cdot} (\mathcal{D}_{sl}\boldsymbol{\nabla}\phi^{s}), \end{gather}

respectively, where $\boldsymbol {u}$ is the bulk granular velocity field, $f_{sl}$ is the segregation velocity magnitude, g is the gravitational acceleration vector and $\mathcal {D}_{sl}$ is the diffusivity of the large and small particles. Since by definition the concentration of large and small particles sum to unity

(2.4)\begin{equation} \phi^{l} + \phi^{s} = 1, \end{equation}

the sum of (2.2) and (2.3) implies that the bulk velocity field is incompressible

(2.5)\begin{equation} \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u} =0. \end{equation}

This is reasonable leading-order approximation that is made in the incompressible $\mu (I)$-rheology (GDR-MiDi 2004; Jop, Forterre & Pouliquen Reference Jop, Forterre and Pouliquen2006; Barker & Gray Reference Barker and Gray2017) for example (where $\mu$ is the friction and I is the inertial number). It should be noted, however, that in reality there is some compressibility and indeed compressibility is needed to make the rheology thermodynamically consistent (Goddard & Lee Reference Goddard and Lee2018) and well posed (Barker et al. Reference Barker, Schaeffer, Shearer and Gray2017; Schaeffer et al. Reference Schaeffer, Barker, Tsuji, Gremaud, Shearer and Gray2019). In this paper, however, the solids volume fraction is assumed to be equal to a constant uniform value of $\varPhi =0.6$.

2.2. Segregation fluxes and velocities

In (2.2) and (2.3) the large and small particle segregation fluxes are

(2.6)\begin{gather} \boldsymbol{F}^{l}={-}f_{sl}\phi^{l}\phi^{s} \frac{\boldsymbol{g}}{|\boldsymbol{g}|}, \end{gather}

(2.7)\begin{gather} \boldsymbol{F}^{s}={+}f_{sl}\phi^{s}\phi^{l} \frac{\boldsymbol{g}}{|\boldsymbol{g}|}, \end{gather}

respectively. These are aligned with the direction of gravitational acceleration $\boldsymbol {g}$ to reflect the fact that the downward percolation of small particles due to kinetic sieving is a gravity driven process. In order to maintain bulk incompressibility there is a net reverse flow of large particles towards the surface due to squeeze expulsion. In the absence of diffusion the segregation velocity of the large and small particles relative to the bulk flow are

(2.8)\begin{gather} \hat{\boldsymbol{u}}^{l}= \boldsymbol{u}^{l}-\boldsymbol{u} ={-}f_{sl}\phi^{s} \frac{\boldsymbol{g}}{|\boldsymbol{g}|}, \end{gather}

(2.9)\begin{gather} \hat{\boldsymbol{u}}^{s} = \boldsymbol{u}^{s}-\boldsymbol{u} ={+}f_{sl}\phi^{l} \frac{\boldsymbol{g}}{|\boldsymbol{g}|}, \end{gather}

respectively. The summation constraint (2.4) implies that the segregation fluxes (2.6)–(2.7) are zero when either of the particle species is in a pure phase. If $f_{sl}$ is independent of $\phi ^{s}$ and $\phi ^{l}$ then (2.4) implies that the large particle segregation flux is quadratic in $\phi ^{s}$ and is symmetric about $\phi ^{s}=1/2$. Equation (2.8) shows that the corresponding large particle segregation velocity is linear in $\phi ^{s}$ and hence implies that the maximum segregation velocity is equal to $f_{sl}$ and is attained when there is a single large grain at concentration $\phi ^{s}=1^{-}$ (note the superscript minus is used to show that the concentration is not 100 % small particles, i.e. there is a large intruder, but it does not significantly change the local concentration from unity). Similarly, the small particle segregation velocity is linear in $\phi ^{l}$ and the maximum segregation velocity equals $-f_{sl}$ and is attained when there is a single small grain at concentration $\phi ^{s}=0^{+}$. The quadratic large particle segregation flux function and the associated linear segregation velocities are plotted in figure 1(a,b).

Figure 1. (a) Segregation flux functions and (b) the large and small particle segregation velocities as a function of the small particle concentration $\phi ^{s}$ (adapted from Gajjar & Gray Reference Gajjar and Gray2014; Gray Reference Gray2018). These assume that $f_{sl}=A_\chi (1-\chi \phi ^{s})$, where $\chi \in [0,1]$ and $A_\chi$ is a normalization factor giving all the flux functions the same amplitude as the quadratic flux when $\chi =0$ (red lines). The blue lines show the cubic flux model of Bridgwater et al. (Reference Bridgwater, Foo and Stephens1985) for $\chi =1$ and the green lines show the intermediate cases when $\chi =0.2$, 0.4, 0.6, 0.8. For $\chi \in [1/2,1]$ these develop an inflection point at $\phi ^{s}_{{inf}}=(1+\chi )/(3\chi )$ (circular markers a) and a group of large particles rise fastest at $\phi ^{s}_{{crit}}=1/(2\chi )$ (circular markers b). Sketches showing (c) the rapid percolation of a single small intruder in a matrix of large grains, (d) the slow rise of a single large particle in a matrix of fine and (e) the faster rise of a group of large particles at intermediate concentrations for a sufficiently skewed flux function.

Using shear box experiments van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) showed that there was an underlying asymmetry to the particle-size segregation, i.e. a single small particle percolating down through a matrix of large particles moves much faster than a single large particle being squeezed up through a matrix of fines. This is illustrated schematically in figure 1(c,d) and is not captured in the simple quadratic flux model. Gajjar & Gray (Reference Gajjar and Gray2014) showed that this asymmetry could, however, be captured by skewing the peak segregation flux towards smaller concentrations of fine grains, as illustrated for a cubic flux function in figure 1(a). As the peak of the large particle segregation flux moves to the left, the maximum downward small particle segregation velocity (at $\phi ^{s}=0^{+}$) is progressively increased even though the maximum amplitude of the flux function stays the same. As the cubic flux curves become skewed, the large particle segregation velocity at concentration $\phi ^{s}=1^{-}$ is progressively diminished below that of the quadratic flux model. Moreover, as the skewness parameter $\chi$ rises above one half, an inflection point develops in the flux function and the maximum segregation velocity of large particles no longer occurs at $\phi ^{s}=1^{-}$, but at an intermediate concentration $\phi ^{s}$ within the range $[1/2,1]$. As a result a group of large particles will rise faster than an individual large intruder as illustrated schematically in figure 1(d,e).

The theory shows that studying (i) a large intruder being squeezed up through a matrix of fine particles and (ii) a small intruder percolating down through a matrix of large grains can yield a considerable amount of information about the shape of the segregation flux function. A series of single-intruder refractive-index matched shear box experiments have therefore been performed in this paper.

3.1. Experimental apparatus

The shear box apparatus is shown in figure 2(a). It consists of two polyvinyl chloride (PVC) rough lateral side plates that are able to pivot about two transverse steel rods, that are located at their centres and bolted to glass panes at the front and back of the shear box. The top ends of the PVC plates are attached to a steel bar that keeps the plates parallel and connects them to a horizontally moving plunger driven by a rotational motor. The bottom of the sheared domain consists of another rough PVC plate of length $L=4.5\ \textrm {cm}$, which sets the distance between the two rotating walls, and is free to move horizontally within slotted grooves in the transparent front and back walls. As the top of the PVC side plates are driven backwards and forwards the plates pivot about their centres and push the base plate from side to side as shown in figure 2(b). To try to prevent the formation of crystalline particle structures the three PVC plates have a random stepped roughness profile, whose average length and depth are of the scale of the background particles. For instance, in the experiment shown in figure 2 the average profile length and depth are 4 mm and 2 mm, respectively, and the background medium size is $d_s=6\ \textrm {mm}$. To give structural stability the four corners of the glass panes were screwed to steel spacer tubes to maintain the transverse separation width $W=7\ \textrm {cm}$ of the shear cell. The volume between the three PVC plates is filled to a height $h=12\ \textrm {cm}$ with a mixture of large and small borosilicate glass beads of diameter $d_l$ and $d_s$, respectively. The whole cell is then submerged in a refractive-index matched fluid to make the semi-opaque intruder easy to identify and track. All the materials used to make the box were selected to be compliant with the refractive index matching technique. The dimensions of the shear cell have been selected to impose a slightly dilated bulk shearing motion on all the grains in the cell. Earlier prototypes had a larger separation distance between the pivoting walls that allowed crystalline particle packings to develop in regions of the flow, which is undesired. A relatively small spacing between the two rotating walls appears to prevent such effects.

Figure 2. (a) A photograph of the experimental apparatus showing the power supply, motor, crank and plunger that drive the upper end of the pivoted parallel PVC plates from side to side. The two pivot points are located in the middle of the sidewalls and one of them is indicated on the photo. The base plate sits in a notched groove and is able to move freely from side to side in response to the driving motion. The volume between the pivoted PVC plates is filled with a granular material of thickness $h$ that is composed of large and small borosilicate glass beads with diameters $d_l$ and $d_s$, respectively. The whole system is submersed in a refractive index matched fluid, consisting of benzyl alcohol and ethanol, which makes the majority of particles transparent and allows the position of the semi-opaque red intruder to be easily identified. (b) An image sequence showing the rise of a large red intruder of diameter $d_l=25\ \textrm {mm}$ through a background medium of small particles of diameter $d_s=6\ \textrm {mm}$ over a number of oscillatory cycles. The small grains are just visible so it is possible to see the free surface. A supplementary movie available at https://doi.org/10.1017/jfm.2021.227 showing the operation of experimental apparatus is available online.

Figure 2(b) and the first supplementary movie show how, over a series of cycles, the large red intruder is progressively squeezed upwards until it reaches the top of the cell. To analyse this motion it is useful to define a coordinate system $Oxyz$ with the horizontal $x$-axis orientated in the direction of shear, the horizontal $y$ axis across the cell and the $z$ axis pointing vertically upwards. Assuming that the PVC plates are pivoted at a height $z=z_{pivot}$ their imposed horizontal displacement field is

(3.1)\begin{equation} x=(z-z_{pivot})\,\sin(\omega t)\tan\theta_{max}, \end{equation}

where $\omega$ is the frequency and $\theta _{max}$ is the maximum angle of displacement of the sidewalls measured from the vertical axis. Typically $\theta _{max}$ ranges between $22.5^{\circ }$ and $25^{\circ }$ in the experiments presented here. Assuming that the resulting velocity $\boldsymbol {u}$ within the granular material is spatially uniform in each $(x,y)$-plane, the resulting velocity components are

(3.2)\begin{equation} u=(z-z_{pivot})\omega\cos(\omega t)\tan\theta_{max},\quad v=0,\ w=0, \end{equation}

respectively. The shear rate $\dot \gamma =2\|\boldsymbol {D}\|$ where $\|\boldsymbol {D}\| = \sqrt {\tfrac {1}{2} \mathrm {tr}(\boldsymbol {D}^{2})}$ is the second invariant of the strain-rate tensor $\boldsymbol {D} = (\boldsymbol {\nabla } \boldsymbol {u} + (\boldsymbol {\nabla } \boldsymbol {u})^{T})/2$. For the imposed bulk velocity field (3.2) it follows that the shear rate

(3.3)\begin{equation} \dot\gamma=\left|\frac{\mathrm{d}u}{\mathrm{d}z}\right|=\omega|\cos(\omega t)|\tan\theta_{max}. \end{equation}

This varies in time, but is spatially uniform throughout the shear cell. The average shear rate over one complete cycle is

(3.4)\begin{equation} \dot\gamma_m=\frac{\omega}{2{\rm \pi}}\int_0^{{2{\rm \pi}}/{\omega}}\dot\gamma\,\mathrm{d}t =\frac{2\omega}{\rm \pi}\tan\theta_{max}. \end{equation}

The motor which drives the oscillatory shear cell has a variable voltage input, so the applied shear rate depends directly on the voltage and the amplitude of the crank connected to the plunger. Monodisperse experiments with four 6 mm opaque tracers in a 6 mm translucent media sheared at $\dot {\gamma }_{m}=1.1\ \textrm {s}^{-1}$ have been performed to check that the segregation is not caused by secondary convection in the underlying flow. A supplementary movie showing this is available online.

3.2. Refractive-index matching, image acquisition and particle tracking

The experiments rely on the refractive-index matching (RIM) technique to make the intruder visible. This is achieved by using transparent borosilicate glass beads for the background media and submerging them in a mixture of benzyl alcohol and ethanol. A summary of the material properties is given in table 1. A calibration procedure was required to get the right index match. An Atago RX 5000 $\alpha$ refractometer was therefore used to measure samples of the interstitial fluid. The initial mixture was based on the measurement tables of Chen, Lin & Tu (Reference Chen, Lin and Tu2012). Further adjustments to obtain the desired value of $n_{r}=1.4726$ were done by adding small quantities of either ethanol or benzyl alcohol. A value of $n_{r}$ in between 1.471 and 1.472 were acceptable to clearly observe the intruder, a condition that improved for a period of time due to the faster evaporation of ethanol. For the same reason, values over 1.473 were undesirable despite an initially good match. The effect of the RIM is shown in figure 3(b) where a large red intruder can clearly be seen in the submersed lower two thirds of the image, but is completely obscured by the small particles in the dry upper third of the image.

Figure 3. (a) Close-up photograph to show how the random profile structure of the pivoting and sliding walls helps to prevent the formation of crystalline structures in the particles surrounding the intruder. (b) Image showing the contrast between dry and wet conditions using the RIM technique with transparent borosilicate glass beads and a mixture of benzyl alcohol and ethanol. (c) Semi-opaque red intruder identified with the image analysis code.

Table 1. Refractive indices $n_{r}$, intrinsic densities $\rho_{\ast}$, viscosities $\eta$ and suppliers for the materials used in the refractive-index matched experiments.

The RIM technique allows the intruder to be clearly seen throughout the experiment (figure 2b). Image acquisition was performed with a Basler acA2000-165uc camera positioned in front of the shear box. Different frame rates were used for each experiment, but for most experiments recording was done at 10 frames per second. Particle identification was done using the imfindcircles routine implemented in Matlab. The single-intruder tracking was performed using the code of Crocker & Grier (Reference Crocker and Grier1996). As well as tracking the opaque intruder this code also helped to identify, and exclude, false intruders, which could occur when the intruder size was close to that of the background medium (figure 3c). Figure 4 shows two space–time plots constructed by taking the vertical line of pixels through the centre of the intruder at each frame and plotting them adjacent to one another. It is striking that for comparable size ratios and shear rates the large intruder takes much longer to be squeezed up to the free surface than the small intruder takes to percolate down to the base. It is also very significant that the trajectories of the intruder centres are both curved.

Figure 4. (a) Space–time plot of a large intruder of diameter $d_l=19\ \textrm {mm}$ segregating upwards in a medium of small particles of diameter $d_s=6\ \textrm {mm}$ subject to an average shear rate $\dot \gamma _m=0.26\ \textrm {s}^{-1}$ (taken from set 2 in table 2). The space–time plot is made by plotting the vertical line of pixels that passes through the centre of the large particle at each time. The intruder therefore appears as a band of red, while the index matched interstitial fluid and the background medium of small particles appear green. The trajectory of the centre of the large intruder as a function of time is shown with the white dashed line. (b) Space–time plot of small intruder particle of diameter $d_s=6\ \textrm {mm}$ percolating down through a matrix of large particles of diameter $d_l=14\ \textrm {mm}$ at a shear rate $\dot \gamma _m=0.34\ \textrm {s}^{-1}$ (taken from set 4 in table 2). In this experiment the background medium appears black while the intruder is orange. The intruder centre as a function of time is shown with a white dashed line.

3.3. Experimental data sets

A series of five sets of experiments were performed in order to understand the functional dependence of the segregation velocity magnitude $f_{sl}$. Sets 1 and 2 correspond to experiments with a single large intruder rising through a medium of small particles. In set 1 the diameter of the large and small particles are held constant and the shear rate is varied, while in set 2 the diameter of the small particles and the shear rate are held constant and the diameter of the large intruder is changed. These sets represent experimental conditions where $\phi =1^{-}$. Sets 3, 4 and 5 all relate to a small intruder segregating through a matrix of large grains, when $\phi =0^{+}$. In sets 3 and 4 the shear rate is varied for two different small intruder diameters, while in set 5 the shear rate is held constant and the small intruder diameter is changed. The experimental conditions of all the sets are summarized in table 2.

Table 2. Diameters of the large particles $d_l$, small particles $d_s$, the size ratio $R=d_l/d_s$ and average shear rates $\dot \gamma _m$ for each experimental set 1–5. Sets 1–2 correspond to a single large intruder rising through a matrix of small particles, whereas sets 3–5 correspond to a single small particle percolating down through a matrix of large grains.

4.1. Dimensional analysis

The experimental data are complex and contain a number of functional dependencies that are hard to interpret without a clear initial hypothesis. It is therefore useful to consider what dimensional analysis can say about the functional behaviour of the segregation velocity magnitude $f_{sl}$ before presenting the results. For a dry bidisperse mixture of large and small particles, $f_{sl}$ is considered to be an output of the system. The inputs are the particle sizes $d_l$ and $d_s$, the intrinsic density of the grains $\rho _*$, the volume fraction of small particles $\phi ^{s}=1-\phi ^{l}$ per unit granular volume, the solids volume fraction $\varPhi$, the bulk shear stress $\tau$, the pressure $p$, gravity $g$ and the shear rate $\dot \gamma$. This is a total of nine variables with three primary dimensions (mass, length and time), so dimensional analysis implies that there are six independent non-dimensional quantities

(4.1a–f)\begin{equation} \mu=\frac{\tau}{p},\quad I=\frac{\dot\gamma \bar d}{\sqrt{p/\rho_*}},\quad \varPhi, \quad P = \frac{p}{\rho_* g \bar d},\quad R=\frac{d_l}{d_s},\quad \phi^{s}, \end{equation}

where $\mu$ is the friction, $I$ is the generalized inertial number in which $\bar d$ is the volume fraction weighted average particle diameter

(4.2)\begin{equation} \bar d=\phi^{l} d_l+\phi^{s} d_s, \end{equation}

$P$ is the non-dimensional pressure and $R=d_l/d_s$ is the grain-size ratio. In this analysis, the experimental system size and any higher gradients of the variables are implicitly assumed not to be of importance. For a monodisperse system in the absence of gravity, only the first three non-dimensional variables in (4.1a–f) are relevant, and (4.2) reduces to the monodisperse particle diameter. It was this type of dimensional analysis combined with observations from DPM simulations and experiments that led to the development of the incompressible $\mu (I)$-rheology (GDR-MiDi 2004; Jop et al. Reference Jop, Forterre and Pouliquen2006), which has advanced our understanding of the rheology of monodisperse granular flows.

Dimensional analysis also provides a powerful way of interpreting the segregation experiments in this paper. It is assumed that it is the combination of shear-induced gravity driven percolation (kinetic sieving) and squeeze expulsion (Middleton Reference Middleton1970; Bridgwater et al. Reference Bridgwater, Foo and Stephens1985; Savage & Lun Reference Savage and Lun1988; Gray & Thornton Reference Gray and Thornton2005; Gray Reference Gray2018) that generate the dominant mechanism for segregation in dense granular flows. In particular, dimensional analysis implies that the segregation velocity magnitude should scale as

(4.3)\begin{equation} f_{sl}\sim \dot\gamma \bar d\, \mathcal{G}(\mu, I,\varPhi, P, R, \phi^{s}), \end{equation}

where $\dot \gamma \bar d$ is chosen as a reference scale and $\mathcal {G}$ is an arbitrary function of the six non-dimensional variables defined in (4.1a–f). The five sets of experiments detailed in table 2 provide key observations that help to constrain the functional form of $\mathcal {G}$.

4.2. Fundamental observations and resulting hypothesis

This paper makes four key assumptions based on the observations (a)–(d) below:- (a) Experimental data sets 1, 3 and 4 imply that the segregation velocity magnitude $f_{sl}$ scales linearly with the shear rate $\dot \gamma$ for both large and small intruders. (b) Experimental set 2 implies that the rise velocity of large intruders scales linearly with the particle size ratio $R$, and since it necessarily shuts off when $R=1$, this implies that $f_{sl}$ scales linearly in $R-1$ in the limit $\phi ^{s}=1^{-}$. (c) Set 5 shows that in the limit $\phi ^{s}=0^{+}$ small intruders exhibit an approximately quadratic dependence on $R-1$ for large size ratios, but have the same linear dependence on $R-1$, as the large intruders, for small size ratios in the range $[1,1.5]$. (d) All of the data sets in 1–5 show that the large and small intruders move along curved trajectories (see e.g. figure 4), with slower percolation and rise rates deeper in the flow, where the lithostatic pressure is higher. These trajectories are well approximated by quadratic curves. Since the pressure is lithostatic in the shear cell, this suggests that $f_{sl}$ is proportional to $1/(\mathcal {C}+P)$, where the non-dimensional constant $\mathcal {C}$ has been introduced to prevent a singularity when $P=0$ at the free surface.

The four key experimental observations (a–d) above suggest that the segregation velocity magnitude can be written as

(4.4)\begin{equation} f_{sl}= \mathcal{B}\, \dot\gamma\bar d\, \frac{\mathcal{F}(R,\phi^{s})}{\mathcal{C}+P}, \end{equation}

where $\mathcal {B}$ and $\mathcal {C}$ are non-dimensional constants and $\mathcal {F}$ is a function of $R$ and $\phi ^{s}$. The particle size-ratio dependence is encapsulated in $\mathcal {F}$, which satisfies the limits

(4.5)\begin{gather} \mathcal{F}(R,1^{-})=R-1, \end{gather}

(4.6)\begin{gather}\mathcal{F}(R,0^{+})=R-1 +\mathcal{E}(R-1)^{2}, \end{gather}

where $\mathcal {E}$ is another non-dimensional constant. A simple functional form that captures the single large and small intruder cases (4.5) and (4.6) is

(4.7)\begin{equation} \mathcal{F}=(R-1)+\mathcal{E}\varLambda(\phi^{s})(R-1)^{2}, \end{equation}

where the function $\varLambda$ satisfies the constraints

(4.8a,b)\begin{equation} \varLambda(0)=1,\quad \text{and}\quad \varLambda(1)=0. \end{equation}

A very simple linear function for $\varLambda$ will be investigated further in § 7.

Substituting for the non-dimensional pressure $P$ from the scaling (4.1a–f) the segregation velocity magnitude can be written as

(4.9)\begin{equation} f_{sl}= \mathcal{B}\,\frac{\rho_*g\dot\gamma \bar d^{2}}{\mathcal{C}\rho_*g\bar d+p}\,\mathcal{F}(R,\phi^{s}). \end{equation}

This formula captures the key processes of gravity, shear and pressure that drive kinetic sieving and squeeze expulsion during the segregation of particles of different sizes and size ratios. The local small particle concentration $\phi ^{s}=1-\phi ^{l}$ enters through the average grain size (4.2) and the nonlinear size-ratio dependence (4.7), and automatically generates asymmetry in the segregation flux functions. Equation (4.9) neglects any dependence on the friction of the particles $\mu$ and the solids volume fraction $\varPhi$. Such dependencies may exist (see Golick & Daniels Reference Golick and Daniels2009; Guillard et al. Reference Guillard, Forterre and Pouliquen2016; Jing et al. Reference Jing, Kwok and Leung2017), but they add an extra level of complexity that goes beyond the scope of this paper. Since $I=\dot \gamma \bar d/\sqrt {p/\rho _*}$ is linear in the shear rate $\dot \gamma$, (4.9) can be reformulated to have a linear inertial number dependence instead of a linear shear-rate dependence. Importantly the experimental observation (a), above, implies that there is either a linear shear-rate dependence or a linear inertial number dependence, but to leading order, at least, other shear rate or inertial number dependencies are not permitted.

4.3. Effect of the interstitial fluid

Equations (4.4) and (4.9) take no account of the interstitial fluid present in the experiments in § 3. Thornton et al. (Reference Thornton, Gray and Hogg2006) used a three phase mixture theory, composed of large grains, small grains and an interstitial fluid, to show that the segregation velocity magnitude $f_{sl}$ was moderated by the relative density difference

(4.10)\begin{equation} \hat\rho=\frac{\rho_*-\rho_*^{f}}{\rho_*}, \end{equation}

where $\rho _*^{f}$ is the density of the fluid. If this multiplicative factor is included then $f_{sl}=0$ when $\rho _*^{f}=\rho _*$, and the particles do not segregate. This is consistent with the fluid density matched bidisperse experiment of Vallance & Savage (Reference Vallance and Savage2000) in which no significant segregation was observed. Vallance & Savage (Reference Vallance and Savage2000) also made the rather counterintuitive observation that the segregation was slightly weaker in water than in a water–ethanol mixture that was 3.7 times more viscous. The viscosity therefore does not suppress segregation as Vallance & Savage (Reference Vallance and Savage2000) had hypothesized, but, instead, it is the relative density difference $\hat \rho$ that explains their results (Thornton et al. Reference Thornton, Gray and Hogg2006). All the experiments in this paper are performed with the same interstitial fluid, and for simplicity the factor $\hat \rho$ is not explicitly included throughout most of the paper. However, in § 8 the dry DPM simulations of Thornton et al. (Reference Thornton, Weinhart, Luding and Bokhove2012) are considered, and to compare the theory it is necessary to correct the experimentally measured value of $\mathcal {B}$ to account for the absence of interstitial fluid.

4.4. Intruder trajectories

In the shear box the gravitational acceleration vector $\boldsymbol {g}$ points downwards parallel to the $z$ axis, i.e. $\boldsymbol {g}=-g\boldsymbol {k}$, where $\boldsymbol {k}$ is the unit vector in the $z$-direction. For the bulk velocity field defined in (3.2) it follows from (2.8) that, in the absence of diffusion, the vertical velocity of the large particles satisfies

(4.11)\begin{equation} w^{l}=\frac{\mathrm{d}z^{l}}{\mathrm{d}t}=f_{sl}\phi^{s}, \end{equation}

where $f_{sl}$ is given by (4.9). For the case of a single large intruder, $\phi ^{s}=1^{-}$, and hence the average particle size $\bar d=d_s$. The large particle trajectory therefore satisfies the ordinary differential equation (ODE)

(4.12)\begin{equation} \frac{\mathrm{d}z^{l}}{\mathrm{d}t}= \mathcal{B}\, \frac{\rho_*g \dot\gamma d_s^{2}}{\mathcal{C}\rho_*gd_s+p}\,\mathcal{F}(R,1^{-}). \end{equation}

Substituting for $\mathcal {F}(R,1^{-})$ from (4.5) and assuming that the pressure in the grains is lithostatic

(4.13)\begin{equation} p=\rho_* g\varPhi (h-z), \end{equation}

the ODE (4.12) reduces to

(4.14)\begin{equation} \frac{\mathrm{d}z^{l}}{\mathrm{d}t}= \mathcal{B}\, \frac{\dot\gamma d_s^{2}(R-1)}{\mathcal{C}d_s+\varPhi (h-z^{l})}. \end{equation}

This is separable and can be integrated, subject to the initial condition that $z^{l}=z^{l}_0$ at $t=0$, to show that time is a quadratic function $\mathcal {Z}^{l}$ of the vertical coordinate $z^{l}$, i.e.

(4.15)\begin{equation} \mathcal{K}^{l} t=Cd_s(z^{l}-z^{l}_0)+\frac{\varPhi}{2}[(h-z^{l}_0)^{2}-(h-z^{l})^{2}]=\mathcal{Z}^{l}(z^{l}), \end{equation}

where the constant

(4.16)\begin{equation} \mathcal{K}^{l}= \mathcal{B}\, \dot\gamma d_s^{2}(R-1), \end{equation}

is different for each individual experiment. Note that $d_s$ and $\varPhi$, which arise in the definition of $\mathcal {Z}^{l}$ in (4.15), are the same for all the experiments in sets 1 and 2. The quadratic equation (4.15) can be solved to give an explicit formula for the trajectory of the large intruder

(4.17)\begin{equation} z^{l} = \frac{1}{\varPhi}\left[ \mathcal{C}d_s+\varPhi h-\sqrt{\mathcal{C}^{2}d_s^{2}+2\mathcal{C}d_s\varPhi(h-z^{l}_0)+\varPhi^{2}(h-z^{l}_0)^{2}-2\varPhi \mathcal{K}^{l} t} \right]. \end{equation}

An exactly analogous argument implies that at $\phi ^{s}=0^{+}$ the trajectory of a small intruder satisfies the ODE

(4.18)\begin{equation} \frac{\mathrm{d}z^{s}}{\mathrm{d}t}={-}\mathcal{B}\, \frac{\dot\gamma d_l^{2}[(R-1)+\mathcal{E}(R-1)^{2}]}{\mathcal{C}d_l+\varPhi (h-z^{s})}. \end{equation}

This is also separable and can be integrated, subject to the initial condition that $z^{s}=z^{s}_0$ at $t=0$, to show that time is a quadratic function $\mathcal {Z}^{s}$ of the vertical coordinate $z^{s}$,

(4.19)\begin{equation} \mathcal{K}^{s} t ={-}\mathcal{C} d_l (z^{s}-z^{s}_0)-\frac{\varPhi}{2}[(h-z^{s}_0)^{2}-(h-z^{s})^{2}]= \mathcal{Z}^{s}(z^{s}), \end{equation}

where the constant

(4.20)\begin{equation} \mathcal{K}^{s}= \mathcal{B}\, \dot\gamma d_l^{2} [(R-1)+\mathcal{E}(R-1)^{2}], \end{equation}

changes between experiments. In all the experiments in sets 3–5 the values of $d_l$ and $\varPhi$ in the definition of $\mathcal {Z}^{s}$ are the same. The quadratic equation (4.19) can be solved to give an explicit formula for the small intruder trajectory

(4.21)\begin{equation} z^{s} = \frac{1}{\varPhi}\left[ \mathcal{C}d_l+\varPhi h-\sqrt{\mathcal{C}^{2}d_l^{2}+2\mathcal{C}d_l\varPhi(h-z^{s}_0)+\varPhi^{2}(h-z^{s}_0)^{2}+2\varPhi \mathcal{K}^{s} t} \right]. \end{equation}

A key test of this scaling argument and the fundamental assumptions (a)–(d) in § 4.2 is whether the intruder trajectories (4.17) and (4.21) are able to collapse all of the data with the same choices of the non-dimensional constants $\mathcal {B}$, $\mathcal {C}$ and $\mathcal {E}$.

4.5. The role of particle diffusion

It is important to note that the intruder trajectories calculated in § 4.4 implicitly assume that there is no diffusion. In fact, experimental observations suggest that there is a significant amount of self-diffusion that leads to the intruders performing random walks around the mean path. Utter & Behringer (Reference Utter and Behringer2004) used a two-dimensional Couette cell to make detailed experimental measurements of the self-diffusion in a monodisperse system of disks. They concluded that the diffusion was anisotropic, but that to leading order the diffusivity was proportional to $\dot \gamma d^{2}$, where $d$ was the particle diameter. The natural generalization of this result to polydisperse systems is to replace the particle diameter with the mean particle diameter (4.2) to give

(4.22)\begin{equation} \mathcal{D}_{sl}= \mathcal{A} \dot\gamma \bar d^{2}, \end{equation}

where $\mathcal {A}$ is a non-dimensional constant. Figure 9(b) in Utter & Behringer's (Reference Utter and Behringer2004) paper shows that the tangential diffusivity implies $\mathcal {A}=0.223$, whereas the radial diffusivity implies $\mathcal {A}=0.108$. The simulations in § 7 adopt the lower value of $\mathcal {A}$, because it is based on measurements of the diffusivity made normal to the direction of shear, which is also the predominant direction for diffusion and segregation in the shear box experiments. Note, that since $\bar d=d_s$ for a large intruder and $\bar d=d_l$ for a small intruder equation (4.22) implies that there will be much more diffusion, and hence variability, for the small intruder paths than for the large intruders.

5.1. Variable shear rate for the large intruders

The experimentally measured position of a large 19 mm intruder in a matrix of small 6 mm particles as a function of time is shown in figure 5, for four different shear rates. These experiments correspond to set 1 in table 2. The intruder struggles to segregate upwards when it is close to the bottom of the cell, but once it moves past a few layers it segregates faster and faster towards the surface. As a result the intruder trajectories all show a pronounced curvature in time that are well fitted by the quadratic curves derived in § 4.4 using the global value of $\mathcal {C}$ defined in (5.1a,b). For each individual trajectory, (4.15) is fitted to the experimental data by determining the constant $\mathcal {K}^{l}$ that gives the best approximation. Typically the fits were extremely good, with a coefficient of determination lying in the range 0.96–0.97. As the shear-rate $\dot \gamma _m$ is increased the large intruder segregates to the surface faster and therefore the time to reach the surface decreases, while the fitting constant increases. The four values of $\mathcal {K}^{l}$ are plotted as a function of $\dot \gamma _m$ in the inset plot in figure 5, together with the theoretical straight line passing through the origin defined in (4.16), i.e.

(5.2)\begin{equation} \mathcal{K}^{l}=\mathcal{B} \dot\gamma_m d_s^{2} (R-1), \end{equation}

where it is assumed that $\dot \gamma =\dot \gamma _m$ and $\mathcal {B}$ is given by the global value in (5.1a,b). The points all lie close to the straight line, although for this subset of data the global fit is not quite the best fit. The data are, however, close to the global fit and is consistent with the fundamental assumption (a) in § 4.2 that the segregation velocity magnitude is linear in the shear rate.

Figure 5. Experimentally measured position of a $19\ \textrm {mm}$ large intruder segregating through a $6\ \textrm {mm}$ matrix of small particles as a function of time, for four different imposed shear rates $\dot \gamma _{m}=0.26$, $0.42$, $0.59$ and $0.77\ \textrm {s}^{-1}$. Darker red lines correspond to higher shear rates. These data correspond to set 1 in table 2. For each case, (4.15) is fitted to the experimental data to determine the constant $\mathcal {K}^{l}$ assuming that the depth $h=12\ \textrm {cm}$ and $\mathcal {C}=0.2712$ is given by the global best fit to all the large intruder data. The fitted intruder trajectory in time (4.17) is then plotted for each $\dot \gamma _m$ with dashed lines for comparison. The coefficient of determination ranges between 0.96 and 0.97 for the proposed fits. The inset shows the fitting constants $\mathcal {K}^{l}$ as a function of $\dot \gamma _m$. These points are closely approximated by the global best fit straight (red dashed) line that passes through the origin (5.2), implying that $f_{sl}$ has a linear shear-rate dependence.

5.2. Variable size ratio for the large intruders

Figure 6 shows experimental set 2 (from table 2) in which both the shear rate $\dot \gamma _m=0.26\ \textrm {s}^{-1}$ and the small particle diameter $d_s=6\ \textrm {mm}$ are held constant, while the diameter of the large intruder $d_l$ is varied. In this limit the average particle size $\bar d=d_s$, but the size ratio $R$ changes between experiments. As the size ratio is increased the particles segregate faster and the time for a particle to rise from the bottom to the top of the cell decreases. In each case the large intruders describe curved trajectories in time as they rise through the cell. This implies that near the bottom of the cell the segregation rate is considerably lower than close to the surface, where the intruder moves notably faster. The intruder trajectories are all well approximated by the quadratic curves derived in § 4.4, although there is some evidence that the curves become more linear at high size ratios. The curves are fitted using a similar procedure to § 5.1, i.e. by determining the best value of the constant $\mathcal {K}^{l}$ for each experiment using the global value of $\mathcal {C}$ defined in (5.1a,b). The fitted values of $\mathcal {K}^{l}$ are plotted as a function of the size ratio $R-1$ in the inset of figure 6. The points lie extremely close to the theoretical straight line (5.2) that passes through the origin, indicating that the segregation rate magnitude $f_{sl}$ is linear in $R-1$ even at high size ratios, as assumed in hypothesis (b) in § 4.2.

Figure 6. Experimentally measured position of a large intruder of size $d_l=10$, $12$, $19$ and $25\ \textrm {mm}$ (darker red lines correspond to larger intruders) in a matrix of 6 mm small particles as a function of time, and at the same imposed shear rate $\dot \gamma _m=0.26\ \textrm {s}^{-1}$. This corresponds to set 2 in table 2. For each case, (4.15) is fitted to the experimental data to determine the constant $\mathcal {K}^{l}$ assuming that the depth $h=12\ \textrm {cm}$ and $\mathcal {C}=0.2712$ is given by the global best fit to all the large intruder data. The fitted intruder trajectory in time (4.17) is then plotted for each size ratio $R$ with dashed lines for comparison. The coefficient of determination ranges between 0.93 and 0.98. The inset shows the measured constants $\mathcal {K}^{l}$ as a function of $R-1$. These points are well approximated by a global best fit straight (red dashed) line that passes through the origin (5.2), implying that $f_{sl}$ has a linear dependence on $R-1$.

5.3. Collapse of all the large intruder experiments

All the large intruder trajectory data from sets 1 and 2 can be collapsed by plotting the quadratic function of the height $\mathcal {Z}^{l}$, defined in (4.15), as a function of scaled time $\mathcal {B}\dot \gamma _m d_s^{2}(R-1)t$ as shown in figure 7. This collapse is based on (4.15) and (5.2), which imply that

(5.3)\begin{equation} \mathcal{Z}^{l}=\mathcal{B} \dot\gamma_m d_s^{2} (R-1)t. \end{equation}

The mapped trajectories should therefore all lie on the 45 degree line shown in figure 7. All the trajectories follow the right trend, but move slightly off the ideal straight line as they rise up through the cell. Some variability is to be expected because the theoretical large intruder trajectory completely ignores the effect of diffusion, which will generate random walks around the average behaviour captured by the theory. Some of the variability is also likely to be an artifact of the oscillatory shear cell, which moves through a maximum angle $\theta _{max}$ before reversing direction, rather than shearing in a consistent direction. This reduces the particle rearrangements, somewhat, and creates brief interludes during which the large intruder oscillates around a given level, before rising upwards again. As can be seen in figure 7 these interludes create a horizontal offset in some of the shear-rate data, before the overall rise continues along a line parallel to the 45 degree line. The inset image in figure 7 shows the values of $\mathcal {K}^{l}$ as a function of $\dot \gamma _m d_s^{2}(R-1)$. All the experimental points for variable shear rate and variable size ratio lie close to the theoretical straight line implied by the global fit with the values of $\mathcal {B}$ and $\mathcal {C}$ from (5.1a,b). This confirms that hypotheses (a), (b) and (d), made in § 4.2, are able to collapse all the large intruder data, with the same constant values of $\mathcal {B}$ and $\mathcal {C}$.

Figure 7. The large intruder trajectories for varying shear rate $\dot \gamma _m$ and size ratio $R$ are collapsed onto an approximately straight (dashed) 45 degree line by plotting the transformed height $\mathcal {Z}^{l}$, defined in (4.15), against the scaled time $\mathcal {B}\dot \gamma _m d_s^{2}(R-1)t$. The data correspond to all the cases in sets 1 and 2 in table 2. The inset plot shows the best fit values of the constants $\mathcal {K}^{l}$ for each of the experiments as a function of $\dot \gamma _m d_s^{2}(R-1)$. The red dashed line in the inset has gradient $\mathcal {B}=0.3744$, i.e. the global best fit value of $\mathcal {B}$. The colours of the points and the lines are the same as those used in figures 5 and 6, where the legends are defined.

6.1. Variable shear rate for the small intruders

Figure 8 shows the trajectories of an 8 mm small intruder in a matrix of 14 mm large particles for three different shear rates corresponding to set 3 in table 2. The three representative profiles are plotted for each shear rate. The theoretical small particle trajectory is fitted to each curve using the same value of $\mathcal {C}$ as in (5.1a,b), and the average value of $\mathcal {K}^{s}$ is then used to determine the overall fit. There is quite a lot of variability about the individual fits, as anticipated, but the fitted trajectory captures the overall behaviour of the small intruders as they percolate downwards. As the shear rate is increased the time taken for intruder to reach the bottom of the cell decreases and $\mathcal {K}^{s}$ increases. The values of $\mathcal {K}^{s}$ are plotted in the inset graph in figure 8 and show a clear linear dependence on the average shear rate $\dot \gamma _m$, consistent with the theoretical line

(6.1)\begin{equation} \mathcal{K}^{s}=\mathcal{B} \dot\gamma_m d_l^{2} [(R-1)+\mathcal{E}(R-1)^{2}], \end{equation}

implied by (4.20). The quadratic dependence on $R-1$ in (6.1) suggests that the segregation rate could be much larger for higher size ratios. To test this, and to test that there is still a linear shear-rate dependence, figure 9 shows how a 6 mm small intruder percolates down through a matrix of 14 mm large particles at various shear rates. This corresponds to a size ratio $R=2.333$, which is not that much larger than $R=1.75$ for the previous experiments in figure 8, but the time taken to reach the bottom of the cell for $\dot \gamma _m=0.34\ \textrm {s}^{-1}$ is approximately a quarter of that when the size ratio was $R=1.75$. Importantly, however, the fitted values of $\mathcal {K}^{s}$ still exhibit a linear dependence on $\dot \gamma _m$ as shown in the inset of figure 9. This observation suggests that the dramatic enhancement in the segregation rate of the small intruders at large size ratios, can be modelled through an additional quadratic size-ratio dependence in the function $\mathcal {F}=\mathcal {F}(R,\phi ^{s})$ defined in (4.7) at low small particle volume fractions.

Figure 8. Measured position of a $8\ \textrm {mm}$ small intruder segregating through a $14\ \textrm {mm}$ matrix of large particles as a function of time for three different imposed shear rates $\dot \gamma _{m}$ (darker turquoise lines correspond to higher shear-rates). Three representative experimental profiles are plotted for each experiment. These data correspond to set 3 in table 2. For each case, (4.19) is fitted to the data set and the value is then averaged over the three realizations to determine $\mathcal {K}^{s}$, assuming the same value of $\mathcal {C}$ as in (5.1a,b). The average fitted intruder trajectory (4.21) is then plotted as a function of time for each $\dot \gamma _m$ with dashed lines for comparison. The insets show the fitting constant $\mathcal {K}^{s}$ as a function of $\dot \gamma _m$. All the points are closely approximated by a straight (dot-dashed) line that passes through the origin, implying that $f_{sl}$ has a linear shear-rate dependence at moderate size ratios.

Figure 9. Measured position of a $6\ \textrm {mm}$ small intruder segregating through a $14\ \textrm {mm}$ matrix of large particles as a function of time for three different imposed shear rates $\dot \gamma _{m}$ (darker magenta lines correspond to higher shear-rates). Three representative experimental profiles are plotted for each experiment. These data correspond to set 4 in table 2. For each case, (4.19) is fitted to the data set and the value is then averaged over the three realizations to determine $\mathcal {K}^{s}$, assuming the same value of $\mathcal {C}$ as in (5.1a,b). The average fitted intruder trajectory (4.21) is then plotted as a function of time for each $\dot \gamma _m$ with dashed lines for comparison. The insets show the fitting constant $\mathcal {K}^{s}$ as a function of $\dot \gamma _m$. All the points are closely approximated by a straight (dot-dashed) line that passes through the origin, implying that $f_{sl}$ has a linear shear-rate dependence even at larger size ratios.

6.2. Variable size ratio for the small intruders

In order to determine the non-dimensional constant $\mathcal {E}$ in (6.1) a series of experiments were performed at the same shear rate, but with variable size ratio. These correspond to set 5 in table 2 and are shown in figure 10. Fits to the three sets of data shown for each case were performed in the same way as in § 6.1 and the values of $\mathcal {K}^{s}$ are shown in the inset of figure 10. For size ratios close to unity the values of $\mathcal {K}^{s}$ lie very close to the red dashed line, corresponding to (6.1) with $\mathcal {E}=0$, using the same coefficients $\mathcal {B}$ and $\mathcal {C}$ that were determined for the large intruders in (5.1a,b). However, as the size ratio increases above $R=1.5$ the values of $\mathcal {K}^{s}$ depart markedly from this line. The additional quadratic dependence on $R-1$ in (6.1) is, however, able to capture the dramatic increase in the segregation velocity at large size ratios. A least squares fit to the data suggests that

(6.2)\begin{equation} \mathcal{E}=2.0957. \end{equation}

The small intruders therefore segregate with a linear $R-1$ dependence for size ratios close to unity, but segregate with an $(R-1)^{2}$ dependence at larger size ratios. This is consistent with the notion that small intruders find it increasingly easy to find holes in which to percolate, as the size ratio increases towards the limit of spontaneous percolation at $R\simeq 6$. Above this threshold, spontaneous percolation implies that small particles can percolate through a matrix of large grains even in the absence of shear (Bridgwater & Ingram Reference Bridgwater and Ingram1971; Wilkinson & Edwards Reference Wilkinson and Edwards1982), so this observation makes physical sense. The result is also backed up by recent DPM simulations of segregation at low small particle concentrations. To see this, compare the inset in figure 10 with figure 6(b) of Chassagne et al. (Reference Chassagne, Maurin, Chauchat, Gray and Frey2020), which also shows a linear increase in the segregation rate for $R$ close to unity and a nonlinear increase at larger grain-size ratios.

Figure 10. Experimentally measured position of a small intruder of size $d_l=12$, $10$, $8$, $6$ and $4\ \textrm {mm}$ in a matrix of 14 mm large particles as a function of time, and at the same imposed shear-rate $\dot \gamma _m=0.34\ \textrm {s}^{-1}$. Three representative experimental profiles are plotted for each experiment. This corresponds to set 5 in table 2. For each case, (4.19) is fitted to the data set and the value is then averaged over the three realizations to determine $\mathcal {K}^{s}$, assuming the same value of $\mathcal {C}$ as in (5.1a,b). The average fitted intruder trajectory (4.21) is then plotted as a function of time for each size ratio $R$ with dashed lines for comparison. The inset shows the fitting constants $\mathcal {K}^{s}$ as a function of $R-1$. At size ratios close to unity the constants follow the same linear dependence on $R-1$ as the large grains (red dashed line), but for larger size ratios the constants $\mathcal {K}^{s}$ become much larger. This is well fitted by the quadratic size-ratio dependence in (6.1) with $\mathcal {E}=2.0957$ (dashed blue line).

6.3. Collapse of all the small intruder experiments

All the small intruder data from sets 3, 4 and 5 in table 2 and shown in figures 8–10 can now be collapsed onto a single curve using the same coefficients $\mathcal {B}$, $\mathcal {C}$ and $\mathcal {E}$. This is achieved by plotting the quadratic function of the height $\mathcal {Z}^{s}$ defined in (4.19) as a function of the scaled time $\mathcal {B} \dot \gamma _m d_l^{2} [(R-1)+\mathcal {E}(R-1)^{2}]t$ as shown in figure 11. This collapse is based on (4.19) and (4.20), which imply that

(6.3)\begin{equation} \mathcal{Z}^{s}=\mathcal{B} \dot\gamma_m d_l^{2} [(R-1)+\mathcal{E}(R-1)^{2}]t. \end{equation}

The mapped trajectories should therefore lie on the 45 degree line in figure 11. The collapse for the individual trajectories is far from perfect for the reasons discussed at the beginning of § 6, which is to be anticipated, but they do follow the right general trend for a wide range of shear rates and particle size ratios. The collapse is much more convincing when the individual fitted values of $\mathcal {K}^{s}$ are plotted against the theoretical line implied by (6.1), which is shown in the inset of figure 11. All the measured points lie very close to the theoretical straight line passing through the origin. This confirms that the small intruder data can be collapsed using the same values of the non-dimensional constants $\mathcal {B}$, $\mathcal {C}$ and $\mathcal {E}$. Moreover, the constants $\mathcal {B}$ and $\mathcal {C}$ are exactly the same as the values used to collapse the large intruder data in § 5. These results therefore validate our preliminary hypotheses (a)–(d) made in § 4.2.

Figure 11. The small intruder trajectories for varying shear rate $\dot \gamma _m$ and size ratio $R$ are collapsed onto an approximately straight (dashed) 45 degree line by plotting the transformed height $\mathcal {Z}^{s}$, defined in (4.19), against the scaled time $\mathcal {B}\dot \gamma _m d_l^{2} [(R-1)+\mathcal {E}(R-1)^{2}] t$. The data correspond to all the cases in sets 3, 4 and 5 in table 2. The inset plot shows the best fit values of the constants $\mathcal {K}^{s}$ for each of the experiments as a function of $\dot \gamma _m d_l^{2}\mathcal {F}(R,0^{+})$. The red dashed line in the inset has gradient $\mathcal {B}=0.3744$, i.e. the global best fit value of $\mathcal {B}$. The colours of the points and the lines are the same as those used in figures 8–10, where the legends are defined.

7.1. Summary of the equations and coefficients

In order to derive the governing equations, the particle-size distribution is assumed to be spatially uniform in the $x$ and $y$ directions. Substituting the gravity vector $\boldsymbol {g}=-g\boldsymbol {k}$ and the velocity field (3.2) into (2.3) implies that the small particle segregation equation reduces to

(7.1)\begin{equation} \frac{\partial{\phi^{s}}}{\partial{t}} - \frac{\partial}{\partial z}(\, f_{sl}\phi^{s}\phi^{l} ) = \frac{\partial}{\partial z}\left( \mathcal{D}_{sl}\,\frac{\partial\phi^s}{\partial z}\right), \end{equation}

where $f_{sl}$ and $\mathcal {D}_{sl}$ are the segregation velocity magnitude and diffusivity, respectively. Substituting the lithostatic pressure (4.13) into (4.9) and cancelling $\rho _*g$, implies that

(7.2)\begin{equation} f_{sl}= \frac{\mathcal{B}\dot\gamma \bar d^{2}\mathcal{F}}{\mathcal{C}\bar d+\varPhi (h-z)}, \end{equation}

where the size-ratio dependence is encapsulated in

(7.3)\begin{equation} \mathcal{F}=(R-1)+\mathcal{E}\varLambda(\phi^{s})(R-1)^{2}. \end{equation}

The function $\varLambda =\varLambda (\phi ^{s})$ satisfies the constraints that $\varLambda (0)=1$ and $\varLambda (1)=0$. It describes the effect of the size ratio at intermediate concentrations, and allows $\mathcal {F}$ to smoothly transition from the linear $R-1$ dependence of a large intruder (4.5) to the $(R-1)^{2}$ dependence (4.6) for a small intruder. The simplest possible form for $\varLambda$ is the linear law

(7.4)\begin{equation} \varLambda=\phi^{l}=1-\phi^{s}, \end{equation}

which will be investigated here. However, much more complicated functions are possible, so long as they satisfy the constraints (4.8a,b). As discussed in § 4.5 the diffusivity is assumed to take the form

(7.5)\begin{equation} \mathcal{D}_{sl}= \mathcal{A} \dot\gamma \bar d^{2}, \end{equation}

where the non-dimensional coefficient $\mathcal {A}$ is determined by the experiments of Utter & Behringer (Reference Utter and Behringer2004). The other non-dimensional coefficients $\mathcal {B}$, $\mathcal {C}$ and $\mathcal {E}$ have been determined by the shear box experiments in this paper and all of them are summarized in table 3.

Table 3. Non-dimensional coefficients $\mathcal {B}$, $\mathcal {C}$, $\mathcal {E}$ determined from the refractive-index matched shear-cell experiments. Note that it is anticipated that in the absence of an interstitial fluid the segregation rate will be enhanced (see § 4.3), which increases the assumed value of $\mathcal {B}$. The non-dimensional coefficient $\mathcal {A}$ is taken from figure 9 of Utter & Behringer (Reference Utter and Behringer2004).

7.2. Non-dimensionalization and numerical method

In order to solve the system of (7.1)–(7.5) it is useful to non-dimensionalize the equations using the same scalings as in van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015), i.e.

(7.6a,b)\begin{equation} t=T\tilde t,\quad z=h\tilde z, \end{equation}

where $T=13\ \textrm {s}$ is the period of one cycle and $h=0.087\ \textrm {m}$ is the depth of the cell. It follows that the segregation equation (7.1) can be written in the non-dimensional form

(7.7)\begin{equation} \frac{\partial{\phi^{s}}}{\partial{\tilde t}} - \frac{\partial{}}{\partial{\tilde z}}( S_r \phi^{s}\phi^{l} ) = \frac{\partial{}}{\partial{\tilde z}} \left( D_r\frac{\partial{\phi^{s}}}{\partial{\tilde z}} \right), \end{equation}

where the non-dimensional segregation rate and non-dimensional diffusivity are

(7.8)\begin{gather} S_r= \frac{T f_{sl}}{h} = \frac{\mathcal{B}\dot\gamma T (\bar d/h)^{2}\mathcal{F}}{\mathcal{C}(\bar d/h)+\varPhi (1-\tilde z)}, \end{gather}

(7.9)\begin{gather} D_r= \frac{T\mathcal{D}_{sl}}{h^{2}}= \mathcal{A} \dot\gamma T (\bar d/h)^{2}, \end{gather}

respectively. The period equals $T=2{\rm \pi} /\omega$ by definition. It therefore follows from (3.4) that the average shear rate $\dot \gamma _m$ times the period

(7.10)\begin{equation} \dot\gamma_m T= 4\tan\theta_{max}, \end{equation}

where $\theta _{max}=30^{\circ }$ in the experiments of van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015). The mean value $\dot \gamma _m T$ is used to approximate $\dot \gamma T$ in (7.8) and (7.9). The non-dimensional segregation equation (7.7) is solved subject to no flux conditions at the surface and the base of the cell

(7.11)\begin{equation} S_r \phi^{s}\phi^{l} + D_r\frac{\partial{\phi^{s}}}{\partial{\tilde z}}=0, \quad\text{at}\ \tilde z=0,1, \end{equation}

and the initial condition

(7.12)\begin{equation} \phi^{s}=\begin{cases} 1, & \tilde z>0.5, \\ 0, & \tilde z\le 0.5, \end{cases} \end{equation}

so that all the small grains are above the large particles at the start of the experiment (van der Vaart et al. Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015). A Galerkin finite element solver (Skeel & Berzins Reference Skeel and Berzins1990) is used to generate the results. This method is conveniently coded in the pdepe routine in Matlab and has been extensively tested in previous papers (Gray & Ancey Reference Gray and Ancey2011; Wiederseiner et al. Reference Wiederseiner, Andreini, Epely-Chauvin, Moser, Monnereau, Gray and Ancey2011; Gray & Ancey Reference Gray and Ancey2015; van der Vaart et al. Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015). It should be noted, however, that this method fails when $\mathcal {C}=0$, because in this case the segregation rate becomes unbounded at the free surface.

7.3. Asymmetric segregation flux functions

It is interesting to see what the consequences the empirically derived scaling law has for the shape of segregation flux function. The non-dimensional segregation equation (7.7) implies that the non-dimensional segregation flux of large and small particles in the vertical direction are

(7.13)\begin{gather} \tilde{F}^{l} = S_r \phi^{s}\phi^{l}, \end{gather}

(7.14)\begin{gather} \tilde{F}^{s} ={-}S_r \phi^{s}\phi^{l}, \end{gather}

respectively. Rather than being a non-dimensional constant (as in Gray & Thornton Reference Gray and Thornton2005; Gray & Chugunov Reference Gray and Chugunov2006), the non-dimensional segregation rate $S_r$ is a function, whose value varies locally in the flow, and which is given by (7.8).

Figure 12(a) shows a series of flux curves for different size ratios $R$, assuming that $\hat z=1/2$, $d_s=4\ \textrm {mm}$ and the non-dimensional coefficients in table 3. All the flux curves are asymmetric and have an inflection point that lies to the right of the maximum. Qualitatively they look similar to the cubic flux functions for $\chi >1/2$ illustrated in figure 1. Crucially, however, these curves have been derived to quantitatively match experimental data and contain complex functional dependencies on the shear rate, the local small particle concentration, the average grain size, the grain-size ratio and the local pressure. In particular, as the grain-size ratio is increased the amplitude of the flux function increases and both the local maximum and the inflection point move to the left.

Figure 12. (a) The non-dimensional segregation flux $S_r\phi ^{s}\phi ^{l}$ as a function of the small particle concentration $\phi ^{s}$, for a range of grain-size ratios $R= 1.5$, 2, 2.5 and 3. The flux is evaluated with $h=87\ \textrm {mm}$, $T=13\ \textrm {s}$, $\theta _{max}=30^{\circ }$, $d_s=4\ \textrm {mm}$ and $\hat z=1/2$, using the non-dimensional coefficients summarized in table 3. The green markers show the position of the maximum and the blue markers show the inflection point. (b) Shows the corresponding non-dimensional large and small particle segregation velocities, $\tilde w^{l}$ and $\tilde w^{s}$, as a function of $\phi$. The yellow markers show the position of the maximum large particle segregation velocity, which occurs at intermediate concentrations.

The increasing skewness with grain-size ratio encodes the idea that small intruders find it increasingly easy to percolate downwards as the size ratio increases, whereas the large intruders find it harder to segregate upwards when there are more contacts with neighbouring small particles. The flux curves have this asymmetric behaviour even when $\mathcal {C}=0$ and $\mathcal {E}=0$, indicating that the scaling of the segregation rate on the mean particle diameter $\bar d$ in (4.9) is already sufficient to produce this effect. The flux curves look similar at other heights in the flow and grow in amplitude as one approaches the free surface. When $\mathcal {C}\ne 0$ there is an upper bound for the amplitude at $\hat z =1$, but when $\mathcal {C}=0$ the segregation flux becomes unbounded at the free surface, which is not desired.

7.4. Large and small particle segregation velocities

In the absence of diffusion the non-dimensional vertical velocities of the large and small particles are given by

(7.15)\begin{gather} \tilde{w}^{l} = S_r\phi^{s}, \end{gather}

(7.16)\begin{gather} \tilde{w}^{s} ={-}S_r\phi^{l}, \end{gather}

respectively. These are plotted in figure 12(b) for the same size ratios and parameters as the flux functions in figure 12(a). For a given size ratio, the downward small particle segregation velocity is an increasing function of the large particle concentration $\phi ^{l}=1-\phi ^{s}$, and therefore the maximum speed is attained at $\phi ^{s}=0^{+}$. Conversely, the large particles have a local maximum segregation velocity at an intermediate concentration $\phi ^{s}\in (0,1)$, with the segregation velocity of a single large intruder at $\phi ^{s}=1^{-}$ being at a much reduced rate, compared to either the maximum large particle, or the maximum small particle, segregation speed. Experimental measurements of the single intruder end states have been specifically used in this paper to determine the scaling law for segregation, and hence the segregation rate $S_r$. The particle velocities in these end states are therefore a close match to the experimental data.

It is interesting to note that when $\mathcal {C}=0$, the ratio of the maximum percolation velocity of a single small intruder to the rise velocity of a single large intruder at a given height $\hat z_0$, satisfies the simple relation

(7.17)\begin{equation} \left.\frac{|\tilde w^{s}(0^{+})|}{|\tilde w^{l}(1^{-})|}\right|_{\hat z=\hat z_0}=\left.\frac{S_r(0^{+})}{S_r(1^{-})}\right|_{\hat z=\hat z_0}=R^{2}\left(1+\mathcal{E}(R-1)\right). \end{equation}

The remarkable equation relates two seemingly disparate ends of the flux curve with a simple relation that is purely dependent on the particle size ratio $R$. It holds at all heights $\hat z_0$ in the flow, which implies that the trajectories of the particles are also related. In particular, it follows from (4.15) and (4.19) that the ratio of the time $t_l$ for a large intruder to reach the top of the cell, and the time $t_s$ for a small intruder to percolate down to the base, satisfies the same law

(7.18)\begin{equation} \frac{t_l}{t_s}=R^{2}\left(1+\mathcal{E}(R-1)\right). \end{equation}

Even when $\mathcal {C}=0.2712$, and there is some dependence on the average grain size and the position in the flow, (7.17) and (7.18) provide good approximations to the ratio of the single particle intruder speeds and the ratio of intruder times, respectively.

7.5. Comparison to van der Vaart's (2015) experimental data

Van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) used an earlier version of the shear cell and the RIM technique to measure the evolution of the particle-size distribution in a $50\ {:}\ 50$ mix of 8 mm and 4 mm particles that was initially in a normally graded configuration. To do this, the shear was stopped after each complete cycle and the cell was scanned with a laser while a series of photographs were taken. From these photos it was possible to determine the three-dimensional particle positions after each cycle and hence to build up a picture of the evolving particle-size distribution in time. A space–time plot of the van der Vaart et al.'s (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) small particle concentration data is shown in figure 13(a).

Figure 13. (a) Contour plot of van der Vaart et al.'s (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) small particle concentration data as a function of the non-dimensional depth $\tilde z=z/h$ and non-dimensional time $\tilde t=t/T$, for a $50\ {:}\ 50$ mix of large and small particles with $h=87\ \textrm {mm}$, $T=13\ \textrm {s}$, $\theta _{max}=30^{\circ }$, $d_s=4\ \textrm {mm}$ and $d_l=8\ \textrm {mm}$. (b) Contour plot of the corresponding numerical results using the non-dimensional coefficients summarized in table 3. (c) Experimental and numerical time series of the depth-averaged concentration of small particles, $\bar \phi ^{s}=\bar \phi ^{s}(t)$.

Figure 13(b) shows a contour plot of the computed small particle concentration as a function of the non-dimensional depth and time, using the experimentally determined non-dimensional coefficients in table 3. The simulations are in very good quantitative agreement with the experimental data of van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015). In particular, they accurately capture the time at which the first small grains reach the base of the flow, the slightly later arrival of the first large grains at the top of the flow, as well as the overall time scale for the solution to approach the steady state by $\tilde t=120$. During this evolution the grains reverse position, with the large particles rising to the surface and the small grains percolating down to the base, to form a slightly diffuse inversely graded particle-size distribution. It should be emphasized that this excellent overall agreement is achieved with the value of $\mathcal {A}$ determined from the experiments of Utter & Behringer (Reference Utter and Behringer2004) and the same non-dimensional coefficients $\mathcal {B}$, $\mathcal {C}$ or $\mathcal {E}$ measured in this paper, together with the simplest possible functional form for $\varLambda$.

The slight mismatches between the simulated and experimental concentration contours in the expansion fan, between $\tilde t=0$ and $60$ time units, are most likely due to experimental error. To demonstrate this, it is useful to integrate the segregation equation (7.7) from the base up to the free surface, and apply the surface and basal no flux boundary conditions (7.11). Exchanging the order of integration and differentiation implies that

(7.19)\begin{equation} \frac{\mathrm{d}}{\mathrm{d}\tilde t}\left(\int_0^{1} \phi^{s}\,\text{d}\tilde z\right)=0, \end{equation}

and hence that the depth-averaged concentration of small particles is conserved, i.e.

(7.20)\begin{equation} \bar\phi^{s}=\int_0^{1} \phi^{s}(\tilde z,\tilde t) \,\text{d}\tilde z=\int_0^{1} \phi^{s}(\tilde z,0)\,\text{d}\tilde z=\bar\phi^{s}_0. \end{equation}

Figure 13(c) shows the depth-averaged concentration $\bar \phi ^{s}$ as a function of time for the numerical simulation and the experiment. While the simulated $\bar \phi ^{s}$ is indeed constant, and equal to the initial depth-averaged concentration of small particle $\bar \phi ^{s}_0=0.5$, the experimental value of $\bar \phi ^{s}$ is up to 5 % below this value in the first $\tilde t=0$ and $60$ time units. This is due to the error in identifying all the small particles with the refractive-index matched fluid technique used by van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015).

7.6. Properties of the steady-state solution

A first-order ODE for the steady-state concentration profile can be formulated by assuming that the concentration is independent of time. Integrating the non-dimensional segregation equation (7.7) with respect to $\tilde z$, applying the no flux condition (7.11) and cancelling the $\dot \gamma T (\bar d/h)^{2}$ dependence yields

(7.21)\begin{equation} \frac{\text{d}\phi^{s}}{\text{d}\tilde z} ={-} \frac{(\mathcal{B}/\mathcal{A})\phi^{s}\phi^{l} \mathcal{F}(R,\phi^{s})}{\mathcal{C}(\bar d/h)+\varPhi (1-\tilde z)}. \end{equation}

This equation has some interesting properties. In particular, the non-dimensional coefficient $\mathcal {C}$ is included in the theory primarily to remove the singularity in the segregation rate at the free surface in the numerical method. When $\mathcal {C}=0$, (7.21) just depends on the particle size ratio $R$, and is completely independent of $\bar d$ and $h$. The thickness does enter in the non-dimensionalization (7.6a,b), so although the non-dimensional solution will be identical for a fixed size ratio, the physical solution will be stretched vertically proportionately to $h$. As a result, a flow of double the thickness will have steady-state concentration gradients that are half as strong as the original flow for the same size ratio. This is a highly non-intuitive consequence of the combination of the scaling law (7.2) and the diffusivity (7.5), and is a bold prediction of the resulting theory.

7.7. Exact steady-state solution for the case $\mathcal {C}=0$

In the case when $\mathcal {C}=0$, it is possible to derive an exact steady-state solution to (7.21). Substituting (7.3) and (7.4) into (7.21) and writing it as an ODE for $\tilde z=\tilde z(\phi ^{s})$ implies

(7.22)\begin{equation} \frac{-1}{(1-\tilde z)}\frac{\text{d}\tilde z}{\text{d}\phi^{s}} = \frac{\mathcal{A}\varPhi}{\mathcal{B}(R-1)\phi^{s}(1-\phi^{s}) \left[1+\mathcal{E}(1-\phi^{s})(R-1)\right]}. \end{equation}

Using partial fractions the right-hand side can be integrated to give the exact solution

(7.23)\begin{equation} \tilde z = 1 - K (1-\phi^{s})^{ -\lambda_1} (1+\mathcal{E}(1-\phi^{s})(R-1))^{\lambda_2} (\phi^{s})^{\lambda_3}, \end{equation}

where $K$ is a constant of integration and the coefficients $\lambda _1$, $\lambda _2$ and $\lambda _3$ are

(7.24a–c)\begin{align} \lambda_1=\frac{\varPhi\mathcal{A}}{\mathcal{B}(R-1)},\quad \lambda_2=\frac{\varPhi\mathcal{A}\mathcal{E}}{\mathcal{B}(1+\mathcal{E}(R-1))},\quad \lambda_3=\frac{\varPhi\mathcal{A}}{\mathcal{B}(R-1)(1+\mathcal{E}(R-1))}, \end{align}

respectively. The average volume fraction of small particles in the cell is obtained by integrating the area under the curve $\phi ^{s}=\phi ^{s}(\tilde z)$ between $\tilde z$ is zero and unity. Since, this area is exactly the same as the area under the curve $\tilde z=\tilde z(\phi ^{s})>0$ in the range $\phi ^{s}=[0,1]$, it is easy to find the constant of integration for any given average concentration of small particles in the cell.

7.8. Steady-state comparison to van der Vaart's (2015) experiment

Figure 14 shows the predicted steady-state exact solution for the experiment of van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015), where $K = 0.4041$ for the initial $50\ {:}\ 50$ mix of particles. At $\hat t=221$ non-dimensional units the exact solution lies very close to the computed concentration profile using the non-zero value of $\mathcal {C}=0.2712$ determined in this paper. In fact both the curves lie very close to van der Vaart et al.'s (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) experimental data between 177 and 221 cycles, when the experiment is essentially in steady state (see figure 14). Following Wiederseiner et al. (Reference Wiederseiner, Andreini, Epely-Chauvin, Moser, Monnereau, Gray and Ancey2011) and Gray (Reference Gray2018), the Péclet number for segregation is defined as the ratio of the non-dimensional segregation rate to the rate of non-dimensional diffusion, i.e.

(7.25)\begin{equation} Pe=\frac{S_r}{D_r}=\frac{\mathcal{B}(R-1)(1+\mathcal{E}(1-\phi^{s})(R-1))}{\mathcal{A}(\mathcal{C}(\bar d/h)+\varPhi(1-\tilde z))}. \end{equation}

It provides a measure of the strength of the segregation relative to the diffusion. As opposed to some earlier theories (see e.g. Gray & Chugunov Reference Gray and Chugunov2006; Wiederseiner et al. Reference Wiederseiner, Andreini, Epely-Chauvin, Moser, Monnereau, Gray and Ancey2011), which had a constant Péclet number, the segregation model developed in this paper produces a Péclet number that varies strongly with depth. This is due to the changing local grain-size distribution and the decreasing lithostatic pressure distribution with increasing height. For the computed solution (see inset figure 14) the Péclet number starts at $Pe=5.74$ at the base of the flow, rises to $Pe=18.86$ at $\tilde z=1/2$ and then tends to a finite, but, large, value of $Pe=430$ at the top. Note that in the case when $\mathcal {C}=0$, the Péclet number will be similar over most of the range, but the singularity in pressure at the free-surface implies that $Pe\rightarrow \infty$ at $\tilde z\rightarrow 1$. Since the diffuse interface separating high concentrations of large particles from high concentrations of fines, in figure 14, lies between $\tilde z=0.4$ and $0.6$, the Péclet number $Pe=11.38$ at the lower end of this range and rises sharply to $Pe=37$ near the top. As a result, the solution is much more rounded and diffuse at the bottom, which is a defining feature of the current theory.

Figure 14. Plot of the computed small concentration $\phi ^{s}$ at $\tilde t=221$ non-dimensional time units (red line) as a function of the non-dimensional depth, for $h=87\ \textrm {mm}$, $T=13\ \textrm {s}$, $\theta _{max}=30^{\circ }$, $d_s=4\ \textrm {mm}$, $d_l=8\ \textrm {mm}$ and using the non-dimensional coefficients summarized in table 3. The blue circles show the corresponding experimental data of van der Vaart et al. (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) between $\hat t=177$ and 221 non-dimensional time units (cycles). The yellow curve shows the exact steady-state solution (when $\mathcal {C}=0$) with a constant of integration $K = 0.4041$ for a $50\ {:}\ 50$ mix of particles as in van der Vaart et al.'s (Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015) experiments. The inset image shows the variation of the Péclet number with depth for the computed solution when $\mathcal {C}\ne 0$.