1. Introduction
A classic problem in fluid mechanics is the terminal velocity of a spherical particle in an otherwise quiescent viscous liquid. A balance of weight, buoyancy and Stokes drag, which depends on the velocity of the particle, yields the terminal velocity. While the result is affected by nearby particles, adjacent walls, etc., the essence of the problem resides in the force balance on the particle.
In this paper, we consider the analogous problem in granular shear flow and use it to determine the segregation velocity in mixtures of two particle species having different sizes. To do this, we use a force balance of the particle weight, the segregation force and the granular drag force to determine the terminal velocity, or, equivalently, the segregation velocity, of an intruder particle in a granular flow and then extend this to the segregation velocity in mixtures of two particle species. However, the construction of the equivalent ‘granular terminal-velocity’ problem requires a crucial addition compared with the fluid problem. An intruder particle in a granular system will move only if energy, e.g. in the form of vibration or shear, is added to the system. In the case of shear, which is considered here, the local shear profile in the intruder vicinity affects the forces acting on it. Furthermore, the buoyancy force on a particle in a granular medium differs somewhat from that for a particle in a fluid. The forces related to shear and buoyancy can be expressed in terms of a ‘segregation force’ that includes both kinematics- and gravity-dependent terms (Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2021). Last, the relationship between the drag force and the particle velocity in a granular flow has recently been clarified to be Stokesian in character (Tripathi & Khakhar Reference Tripathi and Khakhar2013; Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022). Thus, a force balance of the particle weight, the segregation force and the granular drag force can be used to determine the segregation velocity for an intruder particle over a wide range of granular flow conditions.
More generally and more significantly, the granular terminal velocity of an intruder particle can be connected to the problem of particle-size segregation in a flowing mixture of small and large particles with arbitrary finite concentrations. In the mixture, small particles fall through interstices between large particles to lower parts of a flowing layer, thereby forcing large particles upward and resulting in the spatial segregation of initially mixed small and large particles (Gray Reference Gray2018; Umbanhowar, Lueptow & Ottino Reference Umbanhowar, Lueptow and Ottino2019). If segregation of two particle species is viewed as one of the two species migrating relative to the other within the bulk flow (Bridgwater, Foo & Stephens Reference Bridgwater, Foo and Stephens1985; Savage & Lun Reference Savage and Lun1988), the segregation velocity of each species can be thought of as a granular terminal-velocity problem. In this case, the segregation-force and drag-force models need to be extended from the single-intruder limit to mixtures of arbitrary species concentrations. Such extensions bridging particle-level forces and continuum models of the segregation velocity are the focus of recent work (Rousseau et al. Reference Rousseau, Chassagne, Chauchat, Maurin and Frey2021; Tripathi et al. Reference Tripathi, Kumar, Nema and Khakhar2021; Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022, Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024; Sahu et al. Reference Sahu, Kumawat, Agrawal and Tripathi2023), leading to the emergence of a new physics-based approach to segregation-velocity modelling. The present paper aims to complete this new approach and demonstrate its general applicability using a wide range of granular flow configurations.
Modelling particle-size segregation in the continuum framework usually involves solving the spatial and temporal evolution of particle concentration via an advection–diffusion–segregation equation, first suggested by Bridgwater et al. (Reference Bridgwater, Foo and Stephens1985), which differs from a standard advection–diffusion formulation in that a closure relation is needed for the segregation velocity (Gray Reference Gray2018; Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019; Thornton Reference Thornton2021). Therefore, how to model the segregation velocity is the central question of most recent segregation theories. Several different approaches have been used to understand, model and predict the segregation velocity. In one approach, the segregation velocity is determined directly from discrete element method (DEM) simulations (Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019), which empirically connect the local segregation velocity with flow kinematics (shear rate), species concentration and relative particle properties (size ratio and density ratio). Despite the semi-empirical nature of this approach, it has proved effective in capturing segregation fluxes (Fan et al. Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2014; Schlick et al. Reference Schlick, Fan, Isner, Umbanhowar, Ottino and Lueptow2015a ; Jones et al. Reference Jones, Isner, Xiao, Ottino, Umbanhowar and Lueptow2018), density- or shape-induced segregation (Xiao et al. Reference Xiao, Umbanhowar, Ottino and Lueptow2016; Zhao et al. Reference Zhao, Xiao, Umbanhowar and Lueptow2018; Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2021) and fluid effects in segregation (Cui et al. Reference Cui, Zhou and Jing2022), as well as segregation involving both gravity- and shear-gradient-related driving mechanisms (Fan & Hill Reference Fan and Hill2011a ; Liu et al. Reference Liu, Singh and Henann2023; Singh, Liu & Henann Reference Singh, Liu and Henann2024). The drawback of such empirical segregation-velocity models is the lack of a universal model suitable for a wide range of particle properties and flow configurations.
The segregation velocity can also be extracted from particle-species-specific momentum equations (force balance at the continuum level) in a flowing mixture (Gray Reference Gray2018; Thornton Reference Thornton2021). In this framework, the intrinsic pressure gradient of a particle species drives segregation and is counteracted by interspecies drag and diffusive remixing (Gray & Thornton Reference Gray and Thornton2005; Gray & Chugunov Reference Gray and Chugunov2006; Bancroft & Johnson Reference Bancroft and Johnson2021), which is a departure from mixture theory (Atkin & Craine Reference Atkin and Craine1976). Key to such approaches for segregation-velocity models is the closure relation for pressure partitioning and interspecies drag. Despite progress in kinetic-theory-based segregation models for collisional granular flows (Jenkins & Mancini Reference Jenkins and Mancini1987; Jenkins & Yoon Reference Jenkins and Yoon2002; Larcher & Jenkins Reference Larcher and Jenkins2015; Neveu et al. Reference Neveu, Larcher, Delannay, Jenkins and Valance2022), developing closures from first principles remains challenging for size-bidisperse dense granular flows, where multiple, enduring frictional contacts are common. A variety of approaches and assumptions have been used to account for pressure partitioning and drag in dense granular flows (Gray & Thornton Reference Gray and Thornton2005; Fan & Hill Reference Fan and Hill2011b ; Marks, Rognon & Einav Reference Marks, Rognon and Einav2012; Gajjar & Gray Reference Gajjar and Gray2014; Hill & Tan Reference Hill and Tan2014). Experiments have also been used to directly measure the segregation velocity of single intruder particles and mixtures undergoing shear (van der Vaart et al. Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015; Trewhela, Ancey & Gray Reference Trewhela, Ancey and Gray2021). Alternatively, a variety of approaches using virtual springs tethered to a single intruder particle or groups of particles have been used to provide insight into the forces on particles (Guillard, Forterre & Pouliquen Reference Guillard, Forterre and Pouliquen2016; Bancroft & Johnson Reference Bancroft and Johnson2021; Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2021; Liu & Müller Reference Liu and Müller2021; Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022, Reference Duan, Peckham, Umbanhowar, Ottino and Lueptow2023, Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). In all cases, the challenge is to develop a general closure model for the segregation velocity that is applicable over a broad range of granular flow geometries.
In this paper, we develop such a general closure model for the segregation velocity, exploiting recently established segregation-force and drag-force models at the particle level (Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2020, Reference Jing, Ottino, Lueptow and Umbanhowar2021, Reference Jing, Ottino, Umbanhowar and Lueptow2022) along with their extensions to mixtures of arbitrary concentrations (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022, Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). The approach follows the fluid terminal-velocity analogue by using a force balance of the particle weight, the segregation force and the granular drag force to determine the segregation velocity in mixtures of two particle species. Since the segregation-force and drag-force models are characterised based on particle-resolved simulations (i.e. DEM simulations) with measurable parameters, they serve as first-principle closures for the species momentum-balance equations at the continuum level and produce segregation-velocity predictions matching simulation results for size-bidisperse granular flows across a diverse set of flow geometries.
2. Equations of motion
Several models are combined to extract the segregation velocity. In this section, we outline these models and provide details on how to combine them to estimate the segregation velocity under a wide range of flow and particle conditions.
2.1. Particle force balance in a mixture
Consider an incompressible flow of a size-disperse mixture of two species, such as the simple shear flow shown schematically in figure 1(a). The two particle species have volume concentrations
$c_i$
, where
$i=l,s$
for large and small particles, respectively, and
$c_s+c_l=1$
. We neglect vertical acceleration terms, which is reasonable for the relatively slow segregation observed in many common granular flows including heap, chute and rotating-tumbler flows (Gray Reference Gray2018; Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019), though not necessarily all flows (such as some high-speed geophysical flows).

Figure 1. (a) A DEM simulation example of large (4 mm, blue) and small (2 mm, red) spheres in a uniform shear flow with streamwise velocity
$u(z)$
, top wall velocity
$U=u(H)$
where
$H$
is the height of the top wall above the stationary bottom wall, overburden pressure
$P_0$
and downward gravity (negative
$z$
-direction), partitioned into 2.5
$d_l$
high layers (shading) for characterising depth-varying segregation velocity. Here, large particles rise while small particles sink. The segregation direction varies in the different flow configurations analysed later. (b) Force balances on a large particle and a small particle corresponding to (2.1) and species-specific vertical segregation velocities,
$w_i$
.
A force balance at the particle level in the vertical (gravitational,
$g$
) direction for an individual non-accelerating particle with mass
$m_i$
in the flowing mixture, shown in figure 1(b), includes the segregation force,
$F^S_i$
, the weight,
$m_i g$
, and the drag force,
$F^D_i$
such that (Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022)

Much like in the analysis of terminal velocity in a fluid, the segregation velocity in granular flows appears in the drag force,
$F^D_i$
, as we will show shortly. Thus, with an appropriate model for the segregation force,
$F^S_i$
, and a model for the dependence of the drag force,
$F^D_i$
, on the segregation velocity, we can use (2.1) to calculate the segregation velocity. The key is having appropriate models for
$F^S_i$
and
$F^D_i$
on an individual particle in the flowing mixture, which are described shortly.
Using force balance at the particle level, i.e. (2.1), differs somewhat from the continuum description of segregation using the mixture theory framework. Within this framework, the momentum balance for each species along the segregation direction (negative
$z$
-direction) in a simple shear flow scenario can be expressed as (Gray & Thornton Reference Gray and Thornton2005; Tunuguntla, Weinhart & Thornton Reference Tunuguntla, Weinhart and Thornton2017)

Here,
$ -\partial p_i / \partial z= n_i F^S_i$
is the partial pressure gradient, where
$p_i$
is the partial pressure of species
$i$
, and
$z$
is in the direction of gravity, and
$\rho _i$
is the density of species
$i$
. The interspecies momentum exchange is
$\beta _i=n_i F^D_i$
, where
$n_i = c_i \phi /V_i$
represents the particle number density,
$\phi$
is the bulk solid volume fraction, and
$V_i$
denotes the individual particle volume of species
$i$
. Combined with the bulk pressure gradient
$\partial p/\partial z = -\phi \rho g$
, where it is assumed that both species have the same density
$\rho$
, the ratio of the pressure contribution of species
$i$
to the bulk pressure
$p$
, or normal stress fraction, is
$f_i=p_i/p= n_i F^S_i/\phi \rho g=c_i F^S_i/m_i g$
in the simplified case where
$F^S_i$
remains constant with depth. Prior studies adopting the momentum-based approach have often assumed a quadratic-dependence or a particle-size or volume-weighted dependence of
$f_i$
on
$c_i$
(Marks et al. Reference Marks, Rognon and Einav2012; Tunuguntla, Bokhove & Thornton Reference Tunuguntla, Bokhove and Thornton2014; Gray & Ancey Reference Gray and Ancey2015; van der Vaart et al. Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015; Rousseau et al. Reference Rousseau, Chassagne, Chauchat, Maurin and Frey2021; Trewhela et al. Reference Trewhela, Ancey and Gray2021), along with a linear drag model (Gray & Thornton Reference Gray and Thornton2005). While these approaches have been proposed for species concentration profiles in specific scenarios, the approach we use here matches direct DEM measurements of
$f_i$
(Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022), and using (2.1) allows us to consider the segregation velocity of not only mixtures but also intruders for a wide range of granular flow conditions and geometries.
2.2. Segregation force,
$F^S_i$
Determining the segregation force
$F^S_i$
at finite concentration starts with the segregation force on a single intruder,
$F^S_{i,0}$
(subscript
$0$
indicates the single-intruder limit of species
$i$
,
$c_i \to 0$
). Here,
$F^S_{i,0}$
can be modelled with two additive terms, one related to gravity and the other to flow kinematics. Inspired by the observations of Fan & Hill (Reference Fan and Hill2011a
,Reference Fan and Hill
b
), the flow kinematics term was initially scaled with the shear stress gradient (Guillard et al. Reference Guillard, Forterre and Pouliquen2016) but was later linked to the shear rate gradient (Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2021; Singh et al. Reference Singh, Liu and Henann2024). The segregation force can be expressed as (Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2021)

where superscripts
$g$
and
$k$
indicate gravity- and kinematics-related mechanisms, respectively,
$V_i$
is the intruder particle volume,
$\dot \gamma$
is the local shear rate, and
$\rho$
is the density of both the intruder and the bed particles. The gravity term is buoyancy-like, and the kinematic term depends on the shear rate gradient in the flow. The empirical dimensionless functions
$f^g(R_d)$
and
$f^k(R_d)$
depend on the intruder-to-bed-particle-diameter ratio
$R_d=d_i/d_j$
(Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2021),


where
$R^g_1=0.92$
,
$R^g_2=2.94$
,
$c^g_1=1.43$
,
$c^g_2=3.55$
,
$f^k_\infty =0.19$
,
$R^k_1=0.59$
,
$R^k_2=5.48$
and
$c^k_2=3.63$
are fitting parameters appropriate for a range of flow conditions. In applying these functions over a range of concentrations, we need to consider both a large intruder particle surrounded by a bed of small particles and a small intruder particle surrounded by a bed of large particles corresponding to intruder-to-bed particle-size ratios of
$d_l/d_s$
and
$d_s/d_l$
, respectively.
To predict the segregation force
$F^S_i$
in mixtures at arbitrary non-zero concentrations, the intruder-segregation-force model in (2.3) is extended using semi-empirical relations for mixtures (rather than an intruder particle) based on DEM simulations (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022). For large particles in a bidisperse mixture of concentration
$c_l$
(Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024),

and for small particles,

where
$\theta$
is the angle between gravity and the segregation direction, which is generally normal to the flow.
Equations (2.5a
) and (2.5b
) can be rewritten in terms of the net forces (
$T_i=F^S_i-m_ig\cos \theta$
) that balance the interspecies drag (2.1) such that


where
$T_{i,0}=F_{i,0}-m_ig\cos \theta$
, and
$F_{i,0}$
and
$F_{i}$
denote the segregation forces for a single intruder and for mixtures, respectively. The force balance in (2.1) can be rewritten for a mixture as

What remains to be specified is an expression for the drag force in a mixture,
$F^D_i$
, including its dependence on the segregation velocity, which is described next.
2.3. Drag force,
$F^D_i$
As with the approach for the segregation force, we start with the drag force on a single intruder particle within a monodisperse flow of bed particles,
$F^D_{i,0}$
. A Stokes drag formulation can be employed to express the drag force (Tripathi & Khakhar Reference Tripathi and Khakhar2013; Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022; He et al. Reference He, Zhang, Ottino, Umbanhowar and Lueptow2025) in terms of the intruder segregation velocity relative to the local bulk flow velocity in the segregation direction,
$w_{i,0}$
, as

where
$C^D_{i,0}$
is the drag coefficient for a single intruder, and
$\eta$
is the effective bulk granular viscosity calculated from the
$\mu (I)$
rheology as described shortly. For Stokes drag on a spherical particle at low Reynolds number in a viscous fluid,
$C^D_{i,0}=3$
. However, the value of
$C^D_{i,0}$
for an intruder in a granular flow is not as simply specified. For a single spherical intruder particle,
$C^D_{i,0} \approx 2.1$
for
$1\leqslant R_d \leqslant 5$
, but the precise value depends on the flow conditions (Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022), specifically, the size-bidisperse mixture inertial number (Rognon et al. Reference Rognon, Roux, Naaïm and Chevoir2007; Tripathi & Khakhar Reference Tripathi and Khakhar2011),

For large intruders with
$R_d\geqslant 1$
,

where
$k_1=2$
,
$k_2=7$
,
$k_3=2.6$
,
$s_1=0.57$
and
$s_2=0.1$
are fitting parameters determined across a wide range of flow conditions (
$0.6\leqslant R_d \leqslant 5$
,
$1\leqslant R_{\rho } \leqslant 20$
and
$I \lesssim 1$
), and
$R_\rho$
is the intruder-to-bed-particle density ratio (Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022). The granular viscosity
$\eta$
is estimated from the
$\mu (I)$
rheology (GDR-MiDi 2004) as

where the effective friction coefficient is

where
$\tau$
is the shear stress and
$\mu _s$
,
$\mu _2$
and
$I_c$
are granular material specific parameters. The flows simulated in this study generally follow the
$\mu (I)$
rheology for dense flows with the usual caveats for collisional flow at larger
$I$
than considered here and for quasi-static flow as
$I\rightarrow 0$
, see Appendix A.
The drag coefficient,
$C^D_{i,0}$
, in (2.10) is for a single intruder particle in an otherwise homogeneous bed of the other species. Since we are interested in the segregation velocity of particles in a mixture, it is necessary to determine the dependence of the drag force in a mixture,
$F^D_i,$
on species concentration. Previous approaches to determine the mixture
$C^D_i$
have focused predominantly on density-bidisperse mixtures where
$R_d=1$
(Tripathi & Khakhar Reference Tripathi and Khakhar2013; Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2020; Bancroft & Johnson Reference Bancroft and Johnson2021) due to the simplicity of estimating the segregation force in terms of buoyancy. Reported values of
$C^D_i$
at
$R_d=1$
based on this approach range from 1.7 to 3.7 depending on the volume fraction, and
$C^D_i$
is independent of
$c_i$
for density-bidisperse mixtures. However, the concentration dependence of
$C^D_i$
in size-bidisperse mixtures (
$R_d\ne 1$
and density ratio
$R_\rho =1$
) has not been considered explicitly, although Bancroft & Johnson (Reference Bancroft and Johnson2021) mention it in passing. We use simulations of controlled shear flow later in this paper (§ 4) to demonstrate that the drag coefficient is nearly independent of the mixture concentration for the conditions we consider here. For now, in order to proceed with the analysis, we simply assume that
$C^D_i \approx C^D_{i,0}$
. Hence,

where the mixture drag coefficient,
$C^D_i$
, has been substituted for the intruder drag coefficient,
$C^D_{i,0}$
, and the species-specific segregation velocity in the mixture,
$w_i$
, for the intruder segregation velocity,
$w_{i,0}$
, in (2.8).
2.4. Segregation velocity
The species-specific segregation velocity,
$w_i$
, relative to the local bulk flow velocity in the segregation direction is now easily calculated by substituting the mixture drag force (2.12) and the segregation-force model (2.6) into the force balance (2.7) and solving for the segregation velocities of the large-particle species,
$w_l$
, and small-particle species,
$w_s$
:

and

Recall that we assume
$C^D_i \approx C^D_{i,0}$
, independent of species concentration
$c_i$
, which we confirm in § 4.
2.5. Effect of diffusion on segregation velocity
A final consideration is the diffusive flux of species driven by collisional diffusion and its effect on the measured segregation velocity. The diffusion contribution is most clearly framed in terms of the advection–diffusion–segregation transport equation based on mass balance that has been successfully used to model segregation in flowing granular media (Bridgwater et al. Reference Bridgwater, Foo and Stephens1985; Dolgunin & Ukolov Reference Dolgunin and Ukolov1995; Gray Reference Gray2018; Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019). Within this continuum framework, the concentration of species
$i$
can be expressed as

Here,
$\boldsymbol{u}_i$
is the diffusionless velocity, and the local collisional diffusion coefficient
$D$
is a scalar, although in general it is a tensor. This approximation is accurate for flows with a single dominant shear direction (Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019). Note that the diffusionless velocity,
$\boldsymbol{u}_i$
, differs slightly from the overall velocity, which is a combined effect of both advection and diffusion. With the usual assumptions of two-dimensional flow and gradual development in the streamwise direction, (2.14) in the
$z$
-direction can be written as

or, rearranging, as

where
$w$
is the local vertical velocity of the bulk. Note that
$w=0$
in the reference frames associated with the example flows used in this study.
When the normal component of flux for species
$i$
is measured from DEM simulation, it is the entire quantity within the square brackets of (2.16) that is measured. In other words, the measured normal flux
$(w_i+w)c_i-D(\partial c_i/\partial z)$
is a combination of both the total segregation flux,
$(w_i+w)c_i$
, and the diffusion flux,
$-D(\partial c_i/\partial z)$
. To compare the segregation-velocity model predictions developed in this paper with DEM measurements of the segregation velocity, the segregation flux needs to be combined with the diffusion flux, such that the net species velocity is

which can be measured directly in situations where there is a concentration gradient.
Both experimental (Bridgwater Reference Bridgwater1980; Utter & Behringer Reference Utter and Behringer2004) and computational (Fan et al. Reference Fan, Umbanhowar, Ottino and Lueptow2015; Cai et al. Reference Cai, Xiao, Zheng and Zhao2019; Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2019) studies of dense granular flows suggest that the diffusion coefficient,
$D$
, is proportional to the product of the local shear rate and the square of the local mean particle diameter,

where
$\bar d = \sum c_i d_i$
and
$A$
is a constant with reported values in the range 0.01–0.1 (Savage & Dai Reference Savage and Dai1993; Hsiau & Shieh Reference Hsiau and Shieh1999; Utter & Behringer Reference Utter and Behringer2004; Fan et al. Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2014, Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2015; Cai et al. Reference Cai, Xiao, Zheng and Zhao2019; Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2019). In this study,
$A=0.046$
based on diffusion coefficient data measured from heap-flow simulations (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022).
With the exception of demonstrating that the drag coefficient for a mixture is similar to that for an intruder particle, as noted in § 2.3, we now have all of the relationships necessary to calculate the segregation velocity. Before addressing drag in mixtures, it is first necessary to describe the simulation approach that we use.
3. Simulations
An in-house DEM code running on CUDA-enabled NVIDIA GPUs is used to simulate size-bidisperse particle mixtures with species-specific volume concentration
$c_i$
, diameter
$d_i$
and density
$\rho _l=\rho _s=1$
g cm
$^{-3}$
. Large (
$d_l=4$
mm) and small (
$d_s$
varied to adjust the size ratio,
$R_d=d_l/d_s$
) particle species have a
$\pm 10$
% uniform size distribution to minimise layering (Staron & Phillips Reference Staron and Phillips2014) (increasing the diameter variation to
$\pm 20$
% does not alter the results). The mixture is sheared in the streamwise (
$x$
) direction (see figure 1). Boundary conditions are periodic in
$x$
and
$y$
with length
$L=35d_l$
and width
$W=10d_l$
, respectively. The height is
$H=50d_l$
in the
$z$
-direction, which is normal to the flow direction (reducing
$H$
to
$25d_l$
does not alter the results). In all cases, particles fall freely under gravity to fill the domain before flow begins. Gravity may be aligned with the
$z$
-direction, as shown in figure 1, at an angle
$\theta$
with respect to
$z$
for inclined chute flow, or parallel to the flow aligned with
$x$
for vertical chute flow. In some cases, gravity is set to zero; in all other cases we use
$g=g_0\equiv 9.81$
m s
$^{-2}$
.
The standard linear spring-dashpot model (Cundall & Strack Reference Cundall and Strack1979) is used to resolve particle–particle and particle–wall contacts of spherical particles using a friction coefficient of
$\mu =0.5$
, a restitution coefficient of 0.9, and a binary collision time of 0.15 ms. We have confirmed that our results are relatively insensitive to these values except for very low friction coefficients (
$\mu \leqslant 0.2$
) (Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2020; Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2020). From 26 000 to 150 000 particles are included in each simulation, depending on the size ratio.
The segregation velocity is measured in a variety of flow configurations, including controlled shear flows and natural uncontrolled flows. These various flow configurations are explained in more detail in a previous paper in which we consider the segregation force (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). The first flow conditions that we consider are controlled shear flows in which the velocity profile is constrained to be of a certain form. The specified velocity profile,
$u(z)$
, is achieved by applying a small streamwise stabilising force
$k_v [\,u(z)-u_p(z_p)]\,$
to each particle at each DEM simulation time step in order to maintain the desired velocity profile, where
$u_p$
and
$z_p$
are the instantaneous particle velocity and position, respectively, and
$k_v$
is a gain parameter (Lerner, Düring & Wyart Reference Lerner, Düring and Wyart2012; Clark et al. Reference Clark, Thompson, Shattuck, Ouellette and O’Hern2018; Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2018; Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2020, Reference Jing, Ottino, Lueptow and Umbanhowar2021, Reference Jing, Ottino, Umbanhowar and Lueptow2022). By prescribing a specific velocity profile, we control the shear rate and shear rate gradient, which play direct roles in determining both the segregation force (2.3) and the drag (2.8) via the viscosity (2.11), and hence influence the segregation velocity (2.13). The presence of gravity results in a pressure gradient in
$z$
, which also influences the segregation force (2.3). We consider four cases: a linear velocity profile with gravity, an exponential velocity profile without gravity, a parabolic velocity profile without gravity, and an exponential velocity profile with gravity. (For reference, each of these flows is shown schematically later in the paper as insets when the results are discussed, see figure 3
a.) A wide range of inertial numbers,
$I$
, is achieved via the variation of the pressure with depth for flows with
$g\neq 0$
and by imposing large shear rates, which can lead to wall velocities of
$u(H)=U=20$
m s
$^{-1}$
in some cases. For the flows with gravity (linear and exponential velocity profiles), a small overburden pressure
$P_0$
equal to the pressure at a depth of
$2.5 d_l$
(i.e.
$P_0=0.05\rho \phi g_0 H$
, where the bulk solid fraction
$\phi \approx 0.55$
) is imposed on the upper wall, which is free to move vertically, and which fluctuates in height by no more than
$\pm$
0.05 % after an initial rapid dilatation of the particles at flow onset. For the flows without gravity (exponential and parabolic velocity profiles), the top wall is fixed vertically. These different velocity profiles allow us to consider cases with no pressure gradient (when gravity
$g=0$
) so that the gravity-related first term of (2.3) is zero, with no shear rate gradient (linear profile) so that the kinematics-related second term of (2.3) is zero, and with combinations of the gravity and shear such that both the gravity and kinematics terms in (2.3) contribute to the segregation velocity.
In addition to the controlled shear flows, we also consider four cases where the velocity field is not directly controlled. (For reference, each of these flows are shown schematically later in the paper as insets when the results are discussed, see figure 6). The flow kinematics of these uncontrolled ‘natural flows’ are driven entirely by the combined effects of gravity and boundary conditions. The walls are rough in all cases, formed from a
$2.5d_l$
thick layer of bonded large and small particles that move collectively. For the wall-driven flows, an overburden pressure
$P_0$
equal to the pressure at a depth of
$H/2$
(i.e.
$P_0=0.5\rho \phi g_0 H$
) is imposed on the upper wall. When gravity is included for plane shear flow and inclined chute flow, it results in a pressure gradient in
$z$
. In both wall-driven flows, the upper wall moves at velocity
$u(H) = 10\,$
m s
$^{-1}$
in the
$x$
-direction and the lower wall at
$u(0) = -10\,$
m s
$^{-1}$
in the negative
$x$
-direction. Both cases show little to no slip at either wall. The vertical chute flow is driven by gravity, which is aligned parallel to the rough fixed bounding walls, resulting in a generally uniform velocity at the centre of the channel that goes to zero at the walls. In this case, there is no pressure gradient in
$z$
to drive segregation, so any segregation in
$z$
is driven by shear gradients alone. Finally, the inclined chute flow lacks an upper wall (free boundary) so that particles flow due to a streamwise component of gravity. Here the pressure gradient in the segregation direction is
$g_0 \cos {\theta }$
, where
$\theta$
is the inclination angle of the base (lower wall) relative to
$\boldsymbol{g}$
.
To consider segregation for each of these flow conditions, the simulation domain is discretised into horizontal layers, each of thickness
$2.5d_l$
(1 mm) in the
$z$
-direction, for averaging purposes (see figure 1). Decreasing the layer thickness to
$1.25d_l$
increases averaging uncertainties but does not alter the mean values of the flow fields. Within each layer, various local variables are measured including the streamwise velocity (
$u$
), pressure (
$p$
), shear rate (
$\dot \gamma$
) and species concentration (
$c_i$
). Subsequently, these flow measurements, which are averaged in the
$x \hbox {-}$
and
$y \hbox {-}$
directions but vary in the
$z \hbox {-}$
direction, are used to determine intermediate variables including the net gravity and segregation force acting on each species (
$T_i$
via (2.3)–(2.6)), the local viscosity of the mixture (
$\eta$
via (2.11)), the drag coefficient for each species (
$C^D_i$
via (2.10)) and the diffusion coefficient (
$D$
via (2.18)). Finally, these computed variables are used to predict the segregation velocity using (2.13) and, where necessary, (2.17).
The predicted segregation velocities based on the model of (2.13) are compared with the segregation velocities measured from the simulations. To characterise the segregation velocity for each layer in figure 1, we assess the average centre of mass height for each species relative to the mean height of all particles over a short measurement window, calculated as

where
$N_i$
and
$N$
are the number of particles of species
$i$
and the total number of particles in the horizontal averaging layer, respectively. The segregation distance for species
$i$
is the offset of its centre of mass from its initial position,
$\bar z_i - \bar z_{i,0}$
. The segregation velocity for species
$i$
is then measured as the rate of this offset,
$(\bar z_i - \bar z_{i,0})/\Delta t$
, where the measurement window is
$\Delta t=1\,s$
, which we have shown previously is sufficiently long to provide statistically meaningful data and short enough to capture temporally local results (Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2020).
To mitigate the influence of noise and kinematic acceleration on the measurements, each simulation begins with the initial flow of mixed particles subject to equal and opposite vertical restoring forces applied to particles of each species in each layer in figure 1. This technique maintains the initial uniformly mixed distribution of the two particle species and suppresses segregation while the flow develops, similar to the approach previously employed to measure the mixture segregation force (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022, Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). The flow is allowed to develop for 2 s with the prescribed velocity and concentration profiles and the segregation being suppressed. During the subsequent 1 s measurement window, the restoring forces suppressing the segregation are deactivated while primary flow field parameters, such as streamwise velocity, pressure and particle vertical positions (
$z_k$
) are recorded at intervals of 0.01 s. The species segregation velocity is then measured as
$w_i=(\bar z_i-\bar z_{i,0})/ \Delta t$
for the ensemble of each particle species in each layer, and the other local variables (
$u$
,
$p$
,
$c_i$
) in each layer are averaged over the 1 s measurement window for use in calculating the predicted segregation velocity from (2.13). We have shown previously (Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2018; Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2020; Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022, Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024) and confirmed here that the velocity profile is unaffected by the segregation for the short duration of the measurement window. Furthermore, we have confirmed that the concentration profile in the bulk changes by less than 2 % on average over the 1 s measurement window.
4. Drag force in mixtures
Before considering the segregation velocity, which is the focus of this paper, it is necessary to address the effect of mixture concentration on the drag force, as noted in § 2.3. To extend the intruder drag model of (2.8) and (2.10) to mixtures, we perform a series of simulations using an approach motivated by Bancroft & Johnson (Reference Bancroft and Johnson2021), where opposing forces are applied to particles of each species. Specifically, equal and opposite applied forces are imposed in the segregation direction to particles of each species (negative
$z$
-direction for small particles and positive
$z$
-direction for large particles) in a homogeneous shear flow of a mixture of large and small particles like that shown in figure 1 with
$g=0$
. The applied force, which is distributed across all particles of a species, drives the particle species to segregate at a rate controlled entirely by the applied force and the drag, which balance each other (Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022). In this way, (2.1) has just two terms, the mixture segregation force,
$F^S_i$
, which is equivalent to the applied force, and the mixture drag force
$F^D_i$
(since
$g=0$
). Here,
$F^D_i$
is given by (2.12) which is in terms of the effective granular viscosity,
$\eta$
, and the species-specific segregation velocity,
$w_i$
. By tracking the average motion of all of the particles of each species,
$w_i$
for species
$i$
is obtained. By calculating the overall normal and shear stresses,
$P$
and
$\tau$
, from interparticle collisions (Luding Reference Luding2008),
$\eta$
is obtained via (2.11). Using force balance (2.1) with the applied force for
$F^S_i$
, the mixture drag coefficient
$C^D_i$
is estimated using (2.12) for
$F^D_i$
with the measured
$w_i$
and calculated
$\eta$
. In these simulations, the domain is the same as the controlled shear flows used for measuring segregation velocities. The data are recorded at intervals of 0.01 s over a 1 s window after the flow reaches a steady state.

Figure 2. (a) Large-particle drag coefficient,
$C^D_l$
, versus large-particle species concentration,
$c_l$
, in a uniformly sheared flow for size ratios of
$R_d=1.5$
at
$I\approx 0.08$
(blue crosses) and
$R_d=2$
at
$I\approx 0.12$
(black circles) for
$g=0$
. Error bars show the standard deviation of
$C^D_{l}$
over a 1 s window for
$R_d=2$
; error bars for
$R_d=1.5$
are similar but omitted for clarity. Horizontal solid black line corresponds to
$C^D_{i,0}$
for
$R_d=2$
; horizontal dashed blue line corresponds to
$C^D_{i,0}$
for
$R_d=1.5$
. (b) Comparison of
$C^D_{i,0}$
with
$C^D_i$
for varying size ratio. The single-intruder drag coefficient,
$C^D_{i,0}$
is calculated from (2.10) for large (
$i=l$
for
$R_d\geqslant 1$
) (solid black curve) and small (
$i=s$
for
$R_d\lt 1$
) (dashed black curve) intruder particles. The mixture drag coefficient,
$C^D_i$
, (red curve) is calculated from (2.10) for
$R_d\geqslant 1$
and (4.4) for
$R_d\lt 1$
. Both curves represent predictions for
$I=0.2$
. Predictions of the mixture model for
$I$
values ranging from 0 (lower bound) to 0.4 (upper bound), which are typical of dense granular flows, are indicated by the shaded band.
Figure 2(a) shows the drag coefficient of large particles measured in mixtures of particles with varying large-particle concentration,
$c_l$
, in uniform shear flows with
$R_d=1.5$
and 2, having a constant inertial number in each case. The values for
$C^D_l$
are nearly independent of
$c_l$
except at very low
$c_l$
, approaching the single-intruder limit. Furthermore, the values for
$C^D_l$
match the value for
$C^D_{l,0}$
(horizontal dashed line for
$R_d=1.5$
and solid line for
$R_d=2$
). Hence, the intruder drag model (2.8) for
$C^D_{l,0}$
provides a reasonable estimate of
$C^D_l$
for mixtures at arbitrary non-zero concentrations for the cases considered here, although further study of the dependence of drag on species concentration, size ratio and inertial number is warranted.
Unlike the intruder drag model (2.8) that treats the drag force on the intruder, whether it is small or large, as if it is in a sea of particles of the other size, the drag forces for a mixture must also satisfy volume conservation when determining the species-specific segregation velocities,
$w_i$
, for a mixture. This requires that the overall drag force for all of the large particles at a given concentration of large particles must be balanced by the drag force for all of the small particles at the corresponding concentration of small particles, while at the same time volume flux must be conserved (in the laboratory reference frame):

Consequently, there exists an implicit relation between the drag coefficients for large (
$C^D_l$
) and small (
$C^D_s$
) particles to assure that volume flux is conserved. Noting from (2.7) that
$T_i+F^D_i=0$
, expressions for
$w_i$
for each species from (2.12) can be substituted into (4.1) yielding

Using (2.6), this can be expressed as

which can be rearranged as

To satisfy this constraint for size segregation with density ratio
$R_\rho =1$
, we implement the drag model so that the large-particle drag coefficient
$C^D_l$
is estimated by (2.10) but the small-particle drag coefficient
$C^D_s$
is calculated based on the correlation (4.4). The justification for this approach is that (2.10) is valid for
$0.6\leqslant R_d \leqslant 5$
, so we only use it to estimate
$C^D_l$
for large particles (
$R_d\geqslant 1$
) and use (4.4) to find the corresponding value for
$C^D_s$
for small particles that assures volume flux conservation. Figure 2(b) compares
$C^D_i$
predictions by both the original intruder drag model (2.10) and the revised mixture drag model ((2.10) for
$R_d\geqslant 1$
and (4.4) for
$R_d\lt 1$
) for
$I=0.2$
. The two approaches overlap for
$R_d\gt 0.6$
, but the revised model (red curve) allows calculation of
$C^D_s$
for
$0.3 \lesssim R_d\lt 0.6$
, while simultaneously assuring that volume flux is conserved. Note that
$C^D_i$
is nearly independent of
$I$
for
$R_d\lt 1.5$
. Outside this range,
$C^D_i$
depends on
$I$
, with the band of values in figure 2(b) corresponding to
$0\leqslant I\leqslant 0.4$
. However, the effect of
$C^D_i$
on the drag is relatively limited compared with the influence of viscosity due to the constrained range of variation in
$C^D_i$
. Further note that the lower limit in
$R_d$
below which (4.4) should not be applied is necessary because sufficiently small particles can freely pass through interstices between large particles even if the large particles are not flowing. We estimate this lower limit to be
$R_d \approx 0.3$
based on studies of the segregation velocity of small particles in sheared beds (Gao et al. Reference Gao, Ottino, Lueptow and Umbanhowar2024).
5. Segregation-velocity results
In this section, we address the central result of this paper. That is, we combine the models for the segregation force (2.6) in mixtures and the drag force for single intruders (2.8) adapted to mixtures (2.12) via a simple mixture force balance (2.7) to predict the local segregation velocity via (2.13). We then compare the predicted segregation velocity to the result measured directly from the DEM simulated flow to demonstrate the validity of the force-based modelling approach.
5.1. Controlled shear flows

Figure 3. Depth profiles (rows) of time-averaged simulation results (symbols) and predictions (dashed black curves) for the four controlled shear flows (columns) in steady state at
$R_d=2$
. (a) Streamwise mean velocity
$u$
, (b) normalised segregation force on a large particle
$\hat F^S_l=F^S_l/m_l g_0$
, (c) bulk viscosity
$\eta$
and (d) segregation velocity,
$w_i$
, for small (red) and large (blue) particles measured from the simulation (symbols) and predicted via (2.13) (curves). Dotted vertical lines in (b) indicate segregation force equal to particle weight. In all cases,
$U=20\,$
m s
$^{-1}$
,
$c_l=c_s=0.5$
and
$H\approx 0.2\,$
m.
The simulation results and predicted values for the segregation velocity for the four controlled shear flows at
$R_d=2$
are presented in figure 3. Figure 3(a) shows the prescribed and measured streamwise velocity profiles along depth,
$z$
. The close agreement between the DEM data points and the curves representing the target velocity profile demonstrates the effectiveness of the control scheme in achieving the desired velocity profiles.
The prescribed velocity functions in figure 3(a) are used to calculate the
$z$
-profiles of
$\dot \gamma$
and
$\partial \dot \gamma /\partial z$
, while the pressure is estimated based on an idealised hydrostatic pressure,
$p=P_0+\rho \phi g (H-z)$
. These results are then used to calculate the predicted mixture segregation force, shown in figure 3(b), which we normalise with the particle weight,
$\hat F^S_i=F^S_i/m_i g_0$
. The dashed curves represent
$\hat F^S_l$
calculated using (2.3)–(2.5) based on the prescribed pressure and shear states. The symbols indicate
$\hat F^S_l$
measured directly from the simulation using the extended virtual spring approach, where forces proportional to the offset of the centres of mass of the two species are applied uniformly to particles of each species to suppress segregation (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). Error bars indicate the standard deviation of DEM data averaged over the 1 s measurement window. The agreement between the DEM data and the model predictions for
$\hat F^S_l$
across all four controlled shear flows confirms the general applicability of the segregation-force model (2.5) (see also Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). Results for
$\hat F^S_s$
are similar (not shown). Note that the total segregation force equals the depth-wise component of the particle weight for steady flows, specifically,
$c_l \hat F^S_l + c_s \hat F^S_s=\cos \theta$
, where the value on the right-hand side of this equation is shown as a vertical dotted line in figure 3(b). The large difference between measured and predicted values of
$\hat F^S_l$
in the vicinity of
$z/H=0.5$
for the parabolic velocity profile occurs because
$\dot \gamma (H=0.5)=0$
, leading to a singularity in the second term in (2.3), see Duan et al. (Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). However, this does not affect the segregation velocity because of a corresponding singularity in the viscosity, as explained next.
The effective granular viscosity,
$\eta$
, measured as the ratio of shear stress
$\tau$
to shear rate
$\dot \gamma$
, (2.11a
), averaged over the 1 s window from DEM simulation data, is plotted in figure 3(c). The dashed curve represents
$\eta$
estimated from the
$\mu (I)$
rheology (2.11) for the corresponding prescribed pressure-shear state. The predicted and measured values of viscosity match well except near the flow boundaries (see Appendix A). Similar to the segregation force, the viscosity exhibits a singularity at
$\dot \gamma \approx 0$
for the parabolic velocity profile. This singularity, however, is eliminated as
$\dot \gamma \to 0$
when the segregation force, along with the drag force, are incorporated into the force and momentum balance, because both forces scale as
$1/\dot \gamma$
.
With the knowledge of segregation force (figure 3 b), viscosity (figure 3 c) and drag coefficient ((2.10) for large particles and (4.4) for small particles), (2.13) allows the calculation of segregation velocity. Moreover, the depth-dependent segregation velocity can be measured directly for each horizontal layer by quantifying the rate of centre of mass offset between the two species (3.1) in a 1 s window after the flow reaches steady state. Figure 3(d) shows both the model predictions and DEM measurements of the segregation velocity for the four controlled shear flows.
For the linear velocity profile (first column of figure 3), the segregation velocity decreases from top to bottom, despite a constant depth-independent segregation force. This decrease is due to the increase of pressure with depth, which increases the viscosity term in the drag force, thereby slowing the segregation. A similar pattern is observed for the exponential velocity profile without gravity in the second column of figure 3. Despite the constant segregation force and pressure, the shear rate decreases with depth, which increases the viscosity and, hence, the drag. For the parabolic velocity profile (third column of figure 3), the segregation force changes sign at
$\dot {\gamma }=0$
. Consequently, the segregation velocity also changes sign at this point, consistent with large particles segregating to regions of higher shear (Fan & Hill Reference Fan and Hill2011a
). While the model accurately predicts this sign change, the flow in the
$z/H\approx 0.5$
region has negligible shear and, hence, is quasi-static. The continuum framework upon which the model is based is difficult to apply in this regime, resulting in significant discrepancies between the model predictions and the simulation measurements. Finally, in the case of the exponential velocity profile with gravity (fourth column of figure 3), the segregation force increases linearly with depth. However, the viscosity experiences a more pronounced increase, resulting in a decrease in segregation velocity with depth.
Overall, the consistent agreement between the model predictions and the simulation results across all four cases in both magnitude and sign demonstrates the effectiveness of (2.13) in predicting segregation velocity based on the local pressure-shear state. These results further indicate that the model effectively captures the underlying physical mechanisms driving segregation and the associated segregation velocity. The match between the measured segregation velocity and its predicted values in figure 3(d) is not perfect. Nevertheless, given the relatively straightforward model, the simplifying assumptions used in the model, the difficulty in isolating the various forces acting on individual particles, the complication of considering mixtures (rather than an intruder), and the inherent stochastic nature of granular flows, the match between the measured and predicted segregation velocities is surprisingly good.
5.2. Varying species concentration

Figure 4. Profiles of the segregation velocity
$w_i$
for large (blue) and small (red) particles with
$R_d=2$
for the exponential velocity profile with
$g=g_0$
and bulk large-particle concentrations of (a)
$c_l=0.2$
, (b)
$c_l=0.5$
, and (c)
$c_l=0.8$
, based on the prescribed velocity profiles (solid curves) compared with DEM measurements (symbols) averaged over 1 s after the flow reaches steady state.
The analysis in § 5.1 focuses on uniformly mixed systems with equal volume fractions of small and large particles,
$c_l=c_s=0.5$
. However, the concentration dependence of the segregation velocity as described by (2.13) should hold for any species concentration within the range
$0\leqslant c_l\leqslant 1$
, where
$c_s=1-c_l$
. To verify this, we consider controlled shear flows with the exponential velocity profile and
$g=g_0$
, because for these conditions both the segregation force and the viscosity vary with depth (see the last column of figure 3). Figure 4 shows the model predictions for
$w_i$
at
$R_d = 2$
with uniform species concentrations
$c_l = 0.2$
and 0.8, along with the results for
$c_l = 0.5$
(repeated from figure 3 with a different horizontal scale). While a wider scatter of data points is observed with
$c_l=0.2$
, the predicted segregation velocity matches the measured segregation velocity reasonably well. Furthermore, a general trend of decreasing segregation velocity with depth is evident for both particle species, similar to
$c_l=0.5$
. However, the segregation velocity for large particles increases for
$c_l=0.2$
, and the segregation velocity for small particles decreases, compared with
$c_l=0.5$
. That is, for
$c_l=0.2$
, a large particle among mostly small particles segregates faster than a small particle among few large particles. In contrast, for
$c_l=0.8$
, the model underpredicts the segregation velocity, possibly due to the difficulty in accurately determining
$\hat F^S_{s}$
for small intruders in the sea of large particles, which is highly sensitive to
$R_d$
and
$c_l$
near the monodisperse limit of
$c_l=1$
(Jing et al. Reference Jing, Ottino, Lueptow and Umbanhowar2021). Nevertheless, the decreased segregation velocity with depth and the increased (decreased) segregation velocity for small (large) particles compared with
$c_l=0.5$
is clear. Here, a small particle among many large particles segregates faster than a large particle among many small particles, consistent with previously observed asymmetry in the segregation velocity (van der Vaart et al. Reference van der Vaart, Gajjar, Epely-Chauvin, Andreini, Gray and Ancey2015; Jing, Kwok & Leung Reference Jing, Kwok and Leung2017; Jones et al. Reference Jones, Isner, Xiao, Ottino, Umbanhowar and Lueptow2018). This asymmetry in the segregation velocity originates as an asymmetry in the segregation force (2.5), which in turn can be traced back to segregation being intruder-like for small particles at large
$c_l$
and intruder-like for large particles at small
$c_l$
. The species concentration range that the segregation force is intruder-like is narrower for small particles and wider for large particles (see figures 2 and 3 in Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022), which is reflected in the hyperbolic tangent dependence of the segregation force on concentration (2.5) and, consequently, in the form of the segregation velocity (2.13). It is also consistent with the kinetic sieving model of Savage & Lun (Reference Savage and Lun1988), as demonstrated by Jones et al. (Reference Jones, Isner, Xiao, Ottino, Umbanhowar and Lueptow2018) by comparing DEM segregation-velocity results to the solution of the Savage and Lun model over a range of size ratios and species concentrations (see figures 5 and 7 in Jones et al. Reference Jones, Isner, Xiao, Ottino, Umbanhowar and Lueptow2018).

Figure 5. Effect of three different spatially varying concentration profiles (columns and plotted in (a) using the large particle concentration,
$c_l$
) on the segregation velocity
$w_i$
versus depth for (b) linear (
$u=Uz/H$
) and (c) exponential (
$U\textrm {e}^{k((z/H)-1)}$
) velocity profiles with
$g=g_0$
and
$R_d=2$
. In the graphs in (b) and (c), dashed curves represent model predictions using (2.13) for
$w_i$
, solid curves represent predictions corrected by the diffusion flux, i.e.
$w^{net}_i$
from (2.17), and symbols indicate measurements from DEM simulations. Note that the volume fraction,
$\phi$
, in (a) varies only weakly with
$c_l$
.
To this point, the segregation-velocity model (2.13) has been validated against different flow configurations with uniform mixture concentrations throughout the flow domain. However, the model can be extended to scenarios where particle species concentration varies with depth. In such cases, concentration gradients induce a diffusion flux that can enhance or counteract the segregation flux, as described in § 2.5. Because differentiating between these two fluxes is impractical from a measurement standpoint, the observed segregation velocity represents the net effect of both fluxes (see (2.17)).
To consider varying species concentrations, three distinct large-particle concentration profiles, shown in figure 5(a), are investigated:
$c_l$
increasing with depth,
$c_l$
decreasing with depth, and
$c_l$
decreasing in the upper half of the flow and increasing in the lower half of the flow. Although these concentration profiles are initialised to vary linearly with depth, variations in concentration occur as the flow is established, resulting in slightly nonlinear
$c_l$
profiles at the start of the 1 s measurement window. The particle volume fraction,
$\phi$
, remains approximately constant with depth in all three cases. Figure 5(b) shows the segregation velocity for an imposed linear velocity profile with
$g=g_0$
in the
$z$
-direction. Neglecting the diffusion flux leads to an underestimation of the segregation velocity (dashed curves) for both small and large particles, as both the diffusion flux and the segregation flux act in the same direction. Specifically, both segregation and diffusion contribute to the upward movement of large particles. The improved agreement between the diffusion-corrected velocity (
$w_i^{net}$
, solid curves) and the DEM measurements shows the importance of the diffusion effect in inhomogeneous-concentration flows across all three initial conditions for
$c_l$
. Results for an exponential velocity profile with
$g=g_0$
in the
$z$
-direction in figure 5(c) indicate a generally good agreement between the predicted segregation velocity and the DEM results even without accounting for the diffusion flux, except near the top boundary where the correction for diffusion flux is clearly necessary. This is because the diffusion flux scales with shear rate, which is maximal at the top boundary and decreases rapidly with increasing depth.
For all cases in figure 5, the predicted segregation velocity accounting for diffusion matches the measured segregation velocity well. It is also worth noting that the force-based segregation-velocity model accurately captures the change in segregation direction in scenarios where
$c_l$
reaches a minimum value at a mid-depth position for both linear and exponential velocity profiles. This not only highlights the model’s ability to predict segregation velocities, even with substantial variations in the concentration gradient, but also reveals that diffusion can, under certain conditions, become the dominant mechanism, overpowering the effects of size-induced segregation.
5.3. Natural flows

Figure 6. Depth profiles (rows) of time-averaged simulation results (symbols) for the four natural shear flows (columns) in steady state at
$R_d = 1.5$
. (a) Streamwise mean velocity
$u$
, (b) normalised segregation force on a large particle
$\hat F^S_l=F^S_l/m_l g_0$
, (c) bulk viscosity
$\eta$
and (d) segregation velocity
$w_i$
for small (red) and large (blue) particles measured from the simulation (symbols) and predicted by (2.13) (dashed curves) and considering diffusion (2.17) (solid curves). Dotted vertical lines in (b) indicate segregation force equal to particle weight. In all cases,
$c_l=c_s=0.5$
and
$H\approx 0.2\,$
m.
As demonstrated earlier, the segregation-velocity model (2.13) works well for arbitrary concentration fields in a variety of flows where the velocity field is artificially prescribed. We now examine four uncontrolled wall- or gravity-driven flows, illustrated in figure 6(a), in which the velocity field develops naturally via the boundary conditions and gravity-induced body forces. Note that particles in each layer of the flow are constrained vertically using the restoring force approach so that they cannot segregate until they are released when the velocity profile reaches steady state. In wall-driven plane shear flow, which results from upper and lower walls moving in opposite directions with velocity
$U=10$
m s
$^{-1}$
, there is no gravity, so the fully developed velocity profile is linear with a nearly constant inertial number of 0.2 (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024). In the three gravity-driven cases, the orientation of gravity relative to the flow direction is parameterised by
$\theta$
. For wall-driven plane shear flow, gravity acts perpendicular to the flow direction (aligned with the
$z$
-direction such that
$\theta =0$
), resulting in a nonlinear velocity profile due to the increasing pressure with depth in the flow. In vertical chute flow, gravity aligns with the
$z$
-direction (
$\theta =\pi /2$
), such that there is no pressure gradient in the
$z$
-direction,
$\partial p/\partial z=0$
, resulting in a blunt velocity profile with zero velocity at the two confining walls. Finally, for inclined chute flow with
$\theta =28^\circ$
, which exceeds the critical angle required for flow, the resulting velocity profile is nonlinear with a maximum velocity at the free surface.
Following a similar methodology to the analysis of controlled-velocity flow fields, we plot dimensionless depth profiles of
$u$
,
$\hat F^S_i$
,
$\eta$
and
$w_i$
for the natural flows with
$c_l=c_s=0.5$
, as shown in figure 6. Here we consider
$R_d=1.5$
to demonstrate results for a different size ratio than used earlier. Unlike the controlled flows where
$\dot \gamma$
,
$\partial \dot \gamma /\partial z$
,
$p$
and
$\eta$
can be determined from prescribed functions (streamwise velocity, pressure) or assumed constant (concentration), these same variables for natural flow are measured directly from DEM simulations.
In the absence of an intrinsic velocity scale for vertical and inclined chute flows, the kinematic terms are non-dimensionalised using the acceleration due to gravity at the earth’s surface (
$g_0$
) and the flow depth (
$H$
). Unlike the controlled flows in the previous sections, here there are no predicted values for
$\hat F^S_i$
and
$\eta$
because
$u$
and, hence,
$\dot \gamma$
,
$\partial \dot \gamma /\partial z$
,
$p$
and
$\eta$
depend on the naturally developed flow conditions, rather than being prescribed. Nevertheless, the predicted value for the segregation velocity,
$w_i$
, based on measured
$z$
-dependent values of
$\dot \gamma$
,
$\partial \dot \gamma /\partial z$
,
$p$
and
$\eta$
, can be compared with the measured value of
$w_i$
.
Consider first the wall-driven plane shear flow without gravity (first column of figure 6) where the streamwise velocity decreases linearly with increasing depth, resulting in a constant shear rate. This leads to a negligible segregation force across the depth, with only minor deviations observed near the walls. Furthermore, the viscosity remains relatively constant throughout the flow domain. Consequently, the segregation velocity in the first column of figure 6(c) is approximately zero over most of the depth, except close to the walls. The minor deviation from perfect symmetry about
$z/H=0.5$
near the top and bottom walls is likely a consequence of two factors. First, the initial packing procedure, which involves particles falling freely under gravity in the negative
$z$
-direction to fill the domain, may lead to a slight initial segregation of particles. Second, random variations in roughness between the top and bottom boundaries could also contribute to the asymmetry.
With gravity (second column in figure 6), the wall-driven plane shear flow velocity profile is steep near the upper moving wall but flattens in the bottom half of the flow where the pressure is higher. As a result, the segregation force increases with depth except within the region
$z/H\lt 0.3$
. In this region, a rapid decrease in segregation force is observed, accompanied by a reversal in its sign at approximately
$z/H\approx 0.2$
. This is due to the interaction of particles with the moving bottom wall, leading to complex variations in the shear rate gradient for small
$z/H$
. A similar trend also occurs for the viscosity. Consequently, the segregation velocity in the second column of figure 6(d) demonstrates a nonlinear variation with depth, with a reversal in segregation direction at
$z/H\approx 0.2$
. Nevertheless, the predicted segregation velocity matches the measured values well except where it is affected by the walls.
The vertical chute flow (third column of figure 6) has a plug-like velocity profile, resulting in a segregation force that is antisymmetric about
$z/H=0.5$
. Within the central region of the flow, specifically
$0.2\leqslant z/H\leqslant 0.8$
, the shear rate is negligible, leading to difficulty in accurately estimating the local viscosity, which is typical of the quasi-static flow regime. Despite this problem, the model accurately predicts the measured segregation velocity to be approximately zero within this quasi-static central region, as well as the reversal of segregation direction above and below this region.
The curvature of the velocity profile for the inclined chute with
$\theta =28^\circ$
(fourth column of figure 6) is opposite that of the wall-driven flow with gravity (second column of figure 6). The combined effects of shear rate and shear rate gradient result in an almost uniform segregation force except near the static bottom boundary. However, the increase in viscosity with depth leads to a gradual decrease in the segregation velocity, except near the bottom wall, as shown in the fourth column of figure 6(c). Additionally, the model overpredicts the segregation velocity near the free surface (
$z/H \geqslant 0.9$
). This discrepancy is attributed to the inherent difficulty in accurately resolving the free surface using bin-averaging, compounded by the dilute flow regime (
$I\gt 5$
) that prevails near the free surface. Such conditions deviate from the dense flow regime for which the model was developed, resulting in the observed overprediction of
$w_i$
. Nevertheless, through most of the depth of the chute flow, the predicted segregation velocity is consistent with the measured values.

Figure 7. Segregation velocity
$w_i$
for small (red) and large (blue) particles measured from the simulation (symbols) and predicted via (2.13) (dashed curve) and after considering diffusion via (2.17) (solid curve) for chute flow inclined at
$28^\circ$
with different size ratios. In all cases,
$c_l = c_s = 0.5$
and
$H\approx 0.2\,$
m.
To further confirm the validity of our approach for calculating the segregation velocity, model predictions are compared with chute flow simulations for size ratios other than
$R_d=1.5$
. Figure 7 shows reasonably close agreement between the model predictions and DEM measurements for
$R_d =2$
,
$2.5$
and
$3$
, although the model underpredicts the segregation velocity. The overall segregation velocity is largest for
$R_d=2$
and
$2.5$
, but is smaller for
$R_d=1.5$
and
$3$
, consistent with previous findings (Alonso, Satoh & Miyanami Reference Alonso, Satoh and Miyanami1991; Félix & Thomas Reference Félix and Thomas2004; Thornton et al. Reference Thornton, Weinhart, Luding and Bokhove2012; Jones et al. Reference Jones, Isner, Xiao, Ottino, Umbanhowar and Lueptow2018). Similar to the results for
$R_d=1.5$
in figure 6, the discrepancy for
$z/H\geqslant 0.9$
is a result of dilute flow near the free surface. Nevertheless, the agreement observed between the model predictions and the simulation results, not only for various natural flow configurations but also across different size ratios, shows the broad applicability of (2.13) to a wide range of size-bidisperse granular flows at inertial numbers typical of dense flows.
6. Comparison with a previous segregation-velocity model
The broad applicability of (2.13) can be further elucidated by demonstrating its compatibility with the well-established linear segregation-velocity model, initially developed for free-surface heap flows (Fan et al. Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2014), and successfully applied to a wide range of free-surface flows including chute flow, rotating tumbler flow, steady and intermittent heap flow, three-dimensional heap flow, multidisperse species, polydisperse species and density-disperse species (Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019). The model is expressed as (Fan et al. Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2014)

It postulates a direct proportionality between the segregation velocity
$w_i$
and the product of the concentration complement,
$1 - c_i$
, the shear rate,
$\dot \gamma$
, and the particle size,
$d_s$
(or
$d_l$
), and it is motivated by the ‘kinetic sieving’ mechanism and statistical mechanics model for dense granular flows of bidisperse mixtures of spherical particles by Savage & Lun (Reference Savage and Lun1988), as noted previously (Gray & Thornton Reference Gray and Thornton2005; Jones et al. Reference Jones, Isner, Xiao, Ottino, Umbanhowar and Lueptow2018). The segregation parameter
$S$
multiplied by the small particle diameter
$d_s$
quantifies the characteristic segregation length scale. Larger values of
$S$
indicate a stronger tendency for segregation. As determined empirically, this length scale depends on various factors, including particle properties (such as size, shape, density), flow conditions (such as flow depth) and the specific segregation mechanism at play. The effectiveness of the linear segregation-velocity model stems from its ability to encapsulate the complex interplay of segregation factors in a single length scale, while still demonstrating good agreement with simulations and experiments of free-surface flows on heaps or in rotating tumblers. For free-surface flows without strongly segregated regions, the empirically determined segregation length scale
$S$
is well approximated by (Schlick et al. Reference Schlick, Fan, Umbanhowar, Ottino and Lueptow2015a
)

However, the assumption of a linear relationship between segregation velocity and shear rate in (6.1), while simple to implement, neglects the influence of pressure. Consequently, this linear segregation model is limited to free-surface flows where the effects of lithostatic pressure on its predictions are negligible.

Figure 8. Profiles of the segregation velocity
$w_i$
for large (blue) and small (red) particles at the feed zone exit of quasi-2-D heap flows for three different size ratios: (a)
$R_d=1.5$
, (b)
$R_d=2$
and (c)
$R_d=2.5$
. Curves represent predictions from (6.1) (dashed) and (2.13) (solid). Symbols represent DEM measurements averaged over 1 s after the flow becomes steady.
To demonstrate that the force-based segregation-velocity model described in this paper is consistent with this previous simpler approach for free-surface flows, model predictions from (6.1) and (2.13) are compared with data from quasi-two-dimensional (quasi-2-D) heap-flow simulations, which are detailed elsewhere (Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2021). The simulations are performed with a two-dimensional flow rate of 20 cm
$^2$
s
$^{-1}$
and a heap length of 52 cm, resulting in a measured flowing layer thickness (
$\delta$
) of approximately 1.5 cm. The large-particle diameter is 3 mm, while the small-particle diameter is varied between 1.2 and 2 mm, to provide different size ratios. The feed concentration is
$c_l=c_s=0.5$
, which is used as the input concentration for the model. To mitigate the impact of segregation on the concentration profiles, measurements of the segregation velocity are taken immediately downstream of the feed zone, where the concentration profile remains largely unaffected by segregation, and are averaged over a 1 s window.
Previous studies of quasi-2-D heap flow kinematics (Fan et al. Reference Fan, Umbanhowar, Ottino and Lueptow2013) show that the streamwise velocity just downstream of the feed zone exit is well approximated by

where
$k=2.3$
. Hence, we use this velocity profile to determine
$\dot \gamma$
and
$\partial \dot \gamma /\partial z$
, and we assume the pressure at the free surface to be the pressure of a single layer of particles,
$p_0=\rho \phi g \bar d$
. The pressure profile is
$p(z)=\rho \phi g (\bar d +\delta -z) \cos \theta$
, and the feed concentrations of the two particle species are
$c_l=c_s=0.5$
.
Predictions of the linear segregation-velocity model (6.1) and the force-based model (2.13) are generally consistent with each other for all three size ratios, as shown in figure 8. Both methods also match measurements of the local segregation velocity from the simulations, except near the free surface, particularly for
$R_d=2$
and
$2.5$
. Note that in all cases, the small and large particles are assumed to be perfectly mixed with a uniform concentration of
$c_l=c_s=0.5$
at this location in the flow, which implies equal but opposite segregation velocities for the two particle species. This is indeed the case for the measured segregation velocity for
$R_d=1.5$
and 2, where the segregation velocities for both species are nearly symmetric about zero. However, for
$R_d=2.5$
, the asymmetry in the measured segregation velocities of the two particle species suggests that some segregation has occurred in the feed zone. In this case, while both (6.1) and (2.13) accurately predict the segregation velocity of the large particles, they tend to overestimate the segregation velocity of the small particles.
In all cases in figure 8, the force-based model overpredicts the segregation velocities in the near-surface region (
$z/H\gt 0.9$
). Similar to the chute flow results in figure 7, this discrepancy is likely attributed to the dilute flow with diminishing pressure near the free surface, as well as the additional challenges associated with delineating the free surface when employing a cutoff solid fraction as the defining criterion (Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2021). Nevertheless, the close agreement between the force-based model (2.13) and the linear segregation-velocity model (6.1), along with their collective agreement with the measured segregation velocity for
$z/H\lt 0.9$
, provides compelling evidence that the two models are consistent with one another and that the linear segregation model provides a reasonable simplification for heap flow and other free-surface flow scenarios.
Finally, we note that other segregation-velocity models exist, but they are limited in several ways. For one, the effect of shear rate gradients on segregation is not included in many models and scalings (Marks et al. Reference Marks, Rognon and Einav2012; Tunuguntla et al. Reference Tunuguntla, Bokhove and Thornton2014; Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2018; Chassagne et al. Reference Chassagne, Maurin, Chauchat, Gray and Frey2020; Trewhela et al. Reference Trewhela, Ancey and Gray2021), making them inapplicable to most of the flows in this study. Other approaches are more general but lack details (such as particle size and concentration dependence) that would allow them to be applied to the specific flow situations considered here (Hill & Tan Reference Hill and Tan2014; Gray & Ancey Reference Gray and Ancey2015; Bancroft & Johnson Reference Bancroft and Johnson2021; Liu et al. Reference Liu, Singh and Henann2023; Singh et al. Reference Singh, Liu and Henann2024). Hence, even though some of these segregation-velocity models and scalings have demonstrated success for specific flow conditions and geometries, it would be challenging to meaningfully evaluate their relative performance across the broad range of flow conditions considered here.
7. Conclusions
The results in this paper complete our decade-long quest to predict particle segregation in dense granular flows. This effort began using a simple approach to predict size segregation for bidisperse mixtures of dense flowing particles using an advection–diffusion approach (2.15) that includes a term related to segregation (Fan et al. Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2014; Umbanhowar et al. Reference Umbanhowar, Lueptow and Ottino2019). This approach was first proposed four decades ago (Bridgwater et al. Reference Bridgwater, Foo and Stephens1985) and has been built upon by many other researchers (Gray & Thornton Reference Gray and Thornton2005; Marks et al. Reference Marks, Rognon and Einav2012; Tunuguntla et al. Reference Tunuguntla, Bokhove and Thornton2014; Hill & Tan Reference Hill and Tan2014). However, the key to its practical implementation was a simple expression for the segregation velocity (6.1), which was motivated by the much more complicated kinetic sieving model of Savage & Lun (Reference Savage and Lun1988). Using (6.1) for the segregation velocity allows the application of the advection–diffusion–segregation model across a wide range of surface flow geometries including heap flow, tumbler flow, chute flow and hopper flow (Fan et al. Reference Fan, Schlick, Umbanhowar, Ottino and Lueptow2014; Schlick et al. Reference Schlick, Fan, Umbanhowar, Ottino and Lueptow2015b
; Deng et al. Reference Deng, Fan, Theuerkauf, Jacob, Umbanhowar and Lueptow2020). It can be applied not only to bidisperse granular materials, but also to multidisperse and polydisperse particle mixtures (Deng et al. Reference Deng, Umbanhowar, Ottino and Lueptow2018; Gao et al. Reference Gao, Ottino, Umbanhowar and Lueptow2021), density-disperse granular materials, where
$S$
depends on the particle density ratio instead of the size ratio, combined size- and density-disperse granular materials (Duan et al. Reference Duan, Umbanhowar, Ottino and Lueptow2021) and even non-spherical particles (Jones et al. Reference Jones, Ottino, Umbanhowar and Lueptow2020, Reference Jones, Ottino, Umbanhowar and Lueptow2021). And, while it is best suited for segregation in free-surface flows, adjustments can be made such that it can be applied in situations with significant overburden pressures (Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2018).
While successful in many ways, the advection–diffusion–segregation approach using (6.1) for the segregation velocity has been limited by the empirical basis (6.2) for the segregation length scale,
$S d_s$
, which is a function of the size or density ratio, as determined from DEM simulations, although it is possible to estimate
$S$
by matching experiments to predictions using the advection–diffusion–segregation equation (Fry et al. Reference Fry, Vidyapati, Hecht, Umbanhowar, Ottino and Lueptow2020a
,Reference Fry, Vidyapati, Hecht, Umbanhowar, Ottino and Lueptow
b
). Hence, we set out to determine the segregation velocity based on particle-level forces. This requires several ingredients. Foremost is an understanding of the forces acting on an individual particle in a granular flow. While seemingly straightforward, it was not until the development of a virtual spring approach for DEM simulations (Guillard et al. Reference Guillard, Forterre and Pouliquen2016) that it became possible to extract meaningful measurements of forces on an intruder particle in a dense granular flow. This allowed a better understanding of the segregation force acting on the intruder particle, which is a combination of buoyancy-like effects and gradients in the shear rate, expressed compactly in (2.3). A related approach provided the means to measure the drag force on a particle pulled through a sheared granular medium and determine that the drag on an intruder particle is Stokes-like over a range of Reynolds numbers extending over several orders of magnitude (Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022; He et al. Reference He, Zhang, Ottino, Umbanhowar and Lueptow2025), expressed in (2.8), with the drag coefficient dependent on the intruder size ratio and density ratio as well as the inertial number. An additional ingredient in the drag force formulation is the ability to determine a granular viscosity via the
$\mu (I)$
rheology (2.11b
) (GDR-MiDi 2004). The segregation force, drag force and particle weight must balance for an individual intruder particle under the assumption of negligible acceleration, which allows the determination of its segregation velocity, much like how the terminal velocity of a sphere falling in a viscous fluid can be determined from the buoyancy force, drag force and weight. However, as in a fluid suspension where nearby particles and the local flow kinematics alter the sedimentation velocity, other nearby intruder particles and the local flow conditions alter the segregation velocity for a particle in a flowing granular mixture. Using a variant of the virtual spring approach (Duan et al. Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2022, Reference Duan, Jing, Umbanhowar, Ottino and Lueptow2024), the dependence of these forces on species concentration can be determined as given by (2.5). The final ingredient corrects the segregation velocity in a mixture to account for diffusion fluxes originating from concentration gradients, given by (2.17).
As we show in this paper, properly combining all of these ingredients results in the ability to predict the segregation velocity using forces rather than calculating it from the simple relation of (6.1), which depends on an empirical relation for the segregation length scale,
$S d_s$
, and is limited to free-surface flows. This new approach for determining the segregation velocity matches the measured segregation velocity well for the full variety of flow and mixture conditions examined here.
Although determining the segregation velocity using particle-level force models allows the application of the advection–diffusion–segregation equation across a wider range of conditions, there is more work to be done. For instance, predictions of segregation from the force-based approach could be compared with previously proposed segregation kinematics models for specific flow and material combinations (e.g. Fry et al. Reference Fry, Umbanhowar, Ottino and Lueptow2018; Rousseau et al. Reference Rousseau, Chassagne, Chauchat, Maurin and Frey2021; Trewhela et al. Reference Trewhela, Ancey and Gray2021; Jing et al. Reference Jing, Ottino, Umbanhowar and Lueptow2022; Singh et al. Reference Singh, Liu and Henann2024). Extensions to polydisperse particle-size distributions, combined size and density segregation, more complex flow geometries, and non-spherical particles are needed. Of particular interest is increasing the size ratio range that can be accurately modelled. The mixture segregation force (2.5) has only been shown to be valid for
$R_d\lt 3$
and the coefficient of drag (2.10) for
$R_d\lt 5$
. The challenge is that for large
$R_d$
the nature of particle interactions changes. In fact, for
$R_d \geqslant 6.464$
small particles fall freely through the interstices between large particles even if the particles are not undergoing shear, a process called free sifting. This effect begins to influence segregation for
$R_d\gt 3$
(Gao et al. Reference Gao, Ottino, Umbanhowar and Lueptow2023, Reference Gao, Ottino, Lueptow and Umbanhowar2024), well before unhindered free sifting occurs at
$R_d = 6.464$
. Since many practical granular flows in industry and geophysics consist of mixtures with very broad size distributions, more work is clearly warranted to better understand and predict the segregation of ‘fine’ particles with
$R_d \geqslant 4$
. Other physical effects can also play important roles in segregation. Of particular interest in the chemical and pharmaceutical industries is interparticle cohesion, which can result from moisture, surface roughness, electrostatic charging, stickiness and other causes. In some cases, cohesion reduces segregation, while in other cases fine particles agglomerate into large clusters which then segregate. A very different issue is the coupling between the granular flow field and granular segregation, which we ignore here by simply assuming that the flow is relatively unaffected by segregation. Progress is being made in this area (Barker et al. Reference Barker, Rauter, Maguire, Johnson and Gray2021; Edwards et al. Reference Edwards, Rocha, Kokelaar, Johnson and Gray2023; Liu et al. Reference Liu, Singh and Henann2023; Sahu et al. Reference Sahu, Kumawat, Agrawal and Tripathi2023; Maguire et al. Reference Maguire, Barker, Rauter, Johnson and Gray2024; Singh et al. Reference Singh, Liu and Henann2024), but more work can be done. In short, while this paper arguably closes the loop on segregation prediction for non-cohesive particles with
$R\lesssim 3$
by providing a particle-level force-based approach for predicting the segregation velocity, much work remains to be done for larger size ratios, cohesive particles, coupling flow to segregation and many other practical aspects of granular segregation.
Acknowledgements
We acknowledge the use of retrieval-augmented generative AI for language editing purposes during the drafting stage; the associated tool and prompts can be found at GitHub Repository.
Funding
This material is based upon work supported by the National Science Foundation under grant no. CBET-1929265. L.J. gratefully acknowledges financial support provided by the National Natural Science Foundation of China (12472412) and the Open Research Fund Program of State Key Laboratory of Hydroscience and Engineering (sklhse-2023-B-07).
Declaration of interests
The authors report no conflict of interest.
Appendix A. Confirmation of
$\boldsymbol\mu (\boldsymbol{I})$
parameter values
Accurate prediction of the segregation velocity from the flow field necessitates, in addition to segregation force and drag models, a value for the effective granular viscosity,
$\eta$
, in the drag model (2.8). This can be estimated from the
$\mu (I)$
rheology model for dense flows via (2.11a
) and (2.11b
). The accuracy of the
$\mu (I)$
rheology for the flows considered here is confirmed in figure 9, which shows the relationship between the local effective friction coefficient (
$\mu _{\it{eff}}$
) and the local inertial number (
$I$
) for the eight controlled and natural flows included in this study. Nearly all of the data falls on the master curve predicted by the
$\mu (I)$
model with parameters derived from a separate study using simulations with different particle properties (density-bidisperse mixture) and flow geometry (periodic chute) (Tripathi & Khakhar Reference Tripathi and Khakhar2013). Notably, outlier data points correspond to locations near the solid wall boundaries or within quasi-static regions of the flow, which are expected to deviate from the
$\mu (I)$
rheology. Corrections to the
$\mu (I)$
rheology model have been proposed for quasi-static flow (as
$I\rightarrow 0$
) and for
$I \gtrsim 0.3$
(Barker et al. Reference Barker, Schaeffer, Shearer and Gray2017; Heyman et al. Reference Heyman, Delannay, Tabuteau and Valance2017; Barker et al. Reference Barker, Rauter, Maguire, Johnson and Gray2021; Lloyd et al. Reference Lloyd, Maguire, Mistry, Reynolds, Johnson and Gray2025) that involve fitting several parameters to account for deviations from the standard
$\mu (I)$
rheology in these regions. However, these deviations from the
$\mu (I)$
rheology model only occur near the walls and in quasi-static regions where other effects also influence segregation and are not captured by our framework, so we rely on the standard
$\mu (I)$
rheology to estimate the viscosity used in the expression for the drag force, (2.12) and, subsequently, the segregation velocity. Finally, the viscosities estimated from (2.11a
) compare well to the values measured directly from the simulation, as shown in figure 3(c).

Figure 9. Effective friction coefficient
$\mu _{\it{eff}}$
versus local inertial number
$I$
for the eight controlled and natural flows included in this study. The data points represent DEM simulation measurements of the ratio of shear stress to shear rate, defined as
$\mu _{\it{eff}}$
. Circles correspond to flows shown in figures 3 and 6. Outliers represent flows near the boundaries (
$+$
), where
$z/H\lt 0.1$
or
$z/H\gt 0.9$
, and those in the quasi-static regime (
$\times$
) with
$I\lt 0.03$
. The solid curve is the prediction of (2.11b
) with
$\mu _s = 0.364$
,
$\mu _2 = 0.772$
and
$I_c = 0.434$
for data from a previous study of chute flow (Tripathi & Khakhar Reference Tripathi and Khakhar2011).