Simulation of size segregation in granular flow with material point method

Segregation is common in granular flows consisting of mixtures of particles differing in size or density. In gravity-driven flows, both gradients in total pressure (induced by gravity) and gradients in velocity fluctuation fields (often associated with shear rate gradients) work together to govern the evolution of segregation. Since the local shear rate and velocity fluctuations are dependent on the local concentration of the components, understanding the co-evolution of segregation and flow is critical for understanding and predicting flows where there can be a variety of particle sizes and densities, such as in nature and industry. Kinetic theory has proven to be a robust framework for predicting this simultaneous evolution but has a limit in its applicability to dense systems where collisions are highly correlated. In this paper, we introduce a model that captures the coevolution of these evolving dynamics for high density gravity driven granular mixtures. For the segregation dynamics we use a recently developed mixture theory (Fan & Hill 2011, New J. Phys; Hill & Tan 2014, J. Fluid Mech.) which captures the combined effects of gravity and fluctuation fields on segregation evolution in high density granular flows. For the mixture flow dynamics, we use a recently proposed viscous-elastic-plastic constitutive model, which can describe the multi-state behaviors of granular materials, i.e. the granular solid, granular liquid and granular gas mechanical states (Fei et al. 2016, Powder Technol.). The platform we use for implementing this model is a modified Material Point Method (MPM), and we use discrete element method simulations of gravity-driven flow in an inclined channel to demonstrate that this new MPM model can predict the final segregation distribution as well as flow velocity profile well. We then discuss ongoing work where we are using this platform to test the effectiveness of particular segregation models under different boundary conditions.


Introduction
When mixtures of particles with different size or density flowing on an incline under gravity, segregation would happen due to the gradients in total pressure and fluctuation fields.From different perspective of views, there are many effective models that can be used to describe segregation in granular flow.By considering the micro-scale collision possibility between particles, kinetic theory offers a robust framework to explain the mechanism of granular segregation [1,2] ; on the opposite, the percolationdiffusion model [3,4] mainly focuses on the macro-scale description of segregation process by balancing the percolation and diffusion effects in granular flow.The mixture theory, in-between the former two types of model, describes granular segregation with a macro-scale frame work but also offers some particle interaction mechanisms [5,6].In this paper, the mixture theory model proposed by Fan et al [6,7] is employed, in which particle segrega-tion is driven by the space gradient of gravity, density and granular temperature.
Granular materials may behave like solid, liquid or gas when they are under different load and/or dense conditions [8,9], thus a unified constitutive model that can represent the mechanical behaviors of granular materials in different states is of importance in simulating granular flow segregation.Energy analysis of granular materials could offer some hints to develop such a model.The internal energy of particles inside granular materials can be divided as elastic energy and irreversible heat, the elastic energy comes from the reversible elastic contacts between particles, while the later one represents the irreversible dissipation that accompany particle interactions, such as irreversible deformation of a single particle, frictional sliding and structure rearrangement between particles.A nonlinear viscous-elastic-plastic constitutive model that based on this analysis [10] is proposed to capture the dynamic of granular materials in this work.
The multi-state behaviors of granular materials may lead to the co-exist of large and small deformations in a same system, making it another challenge to simulate granular flow segregation with traditional Lagrangian and Eulerian methods.As an effective numerical method that mixed the advantage of Lagrangian and Euler schemes, material point method offers a good platform for that task.In MPM, granular media is divided as a set of material points that store all the neighbor material information around them.Meanwhile, momentum conservation equations are solved on Euler background gird nodes.During the calculation, material information is transferred between Lagrangian material points and Euler background gird nodes with the help of shape function.The Lagrangian material points make it simpler to trace the deformation and movement of materials when compared to traditional Eulerian method; meanwhile, in the end of each calculation circle, Euler background grids would be reset, accuracy loss caused by large mesh deformation in Lagrangian method is avoid [18].Thus, MPM has the benefit of flexibility in tracking evolution of densities, velocities, and stress fields for different mechanical states of granular materials 2 Models

Mixture theory segregation model
For a binary size, dense granular mixture flowing on an incline with a slope angel of θ, the segregation flux Ψ i T of size specie i (i = l, s for large and small particles, respectively) along segregation direction y (normal to the incline) in steady flow state can be expressed as [7] where f denotes the volume fraction of granular solid, v i and v are the flow velocity of specie i and the mixture along y direction.ρ is the density of granular mixture, g cos θ denotes the component of gravity acceleration g normal to the incline.c D and η are the linear drag and diffusion coefficients, and σ k yy = ρv v is the kinetic stress of granular material, with v the peculiar velocity of the mixture.Kinetic stress reflects the instant interaction between particles, compared to that, the enduring contact stress between particles is expressed as contact stress σ c yy .ψ c,i = σ c,i yy /σ c yy and ψ k,i = σ k,i yy /σ k yy are the stress partition coefficients of contact and kinetic stress, which show the weight of specie stresses in the mixture.φ i = f i / f denotes the concentration of specie i.
It shows clearly in Eq.( 1) that segregation is conducted by the gradients of the fluctuation of flow, the gravity and the density; the first term reflects the effect of shear in granular segregation explicitly.However, velocity and stress are just input parameters in this model.In reference [6,7], these parameters are measured directly from DEM simulation.To get the flow dynamic information of granular mixture, a constitutive model is required.

Non-linear viscous-elastic-plastic constitutive model
Based on the energy analysis mentioned ahead, mechanical behaviors inside granular materials are divided into the reversible elastic effect results from contacts between particles, and the irreversible dissipation caused by the frictional sliding and structure rearrangement inside granular materials.Respectively, the total stress of granular material σ i j is divided into two terms, the elastic stress σ e i j that reflects the elastic energy stored inside granular material, and the dissipative stress σ d i j that describes the irreversible energy dissipation of granular materials during the momentum transfer caused by the velocity gradient.
Elastic stress is deduced from the derivation of elastic energy density w e .A non-linear elastic energy density model first proposed by Jiang & Liu [11] and further improved by Sun et al. [12] is used, thus elastic stress σ e i j can be expressed as (3) with A, C the elastic modulus related parameters, which are constant in this model, C = 5q 2 φ A, q φ is the internal friction coefficient of granular material.I 1 = −ε e ii , I 2 = ε e,0 i j ε e,0 ji are invariants of the elastic strain tensor ε e i j e i e j , with ε e,0 i j = ε e i j + I 1 /3.For dense granular flows, we have I 1 > 0.
To reflect the dissipative nature of granular materials, we refer the so called μ(I) rheology [13,14] to express the dissipative stress σ d i j , with D i j = (γ i j −δ i j γkk /3)/γ, γij = dε i j /dt is the shear rate, γ = γij γ ji /2.I = γd/ σ e m /ρ m is the inertial number, ρ m and d denote the material density and diameter of particle, respectively.a and b are material property parameters that is related to the dissipative nature of granular materials.
According to Eqs. ( 3) and ( 4), the total stress of granular materials can be expressed as with H(I 1 ) the Heaviside function.When I 1 ≤ 0, H(I 1 ) = 0, which means particles lose contacts with each other, and the system would be in gaseous state.

Combined flow-segregation model
In the segregation model of Eq. ( 1), segregation flux is determined by the spatial distribution of stresses and their partition coefficients, but stress information is not supported by the model itself.In this work, however, stress profiles are achieved at the same time of segregation calculation by combing Eqs. ( 1) and (5).Considering that kinetic stress σ k yy represents the irreversible fluctuation of particles while σ c yy is the enduring contact between particles, their physical meaning are similar to that of σ d and σ e in Eq. ( 5), σ k yy and σ c yy in Eq. ( 1) are replaced by Eq. ( 4) and Eq. ( 3).What's more, along the segregation direction, elastic stress is balanced by gravity component in steady flow, and ψ d,i = ψ d,i (φ i ) can be expressed as a function of concentration φ i [7, 17], thus we have is a material parameter.In Eq. ( 6), φ i is the target valuable that is to be solved, if the initial profile of concentration φ i and the evolution of kinetic stress profile σ d,i are knew, the evolution of segregation in steady granular flow can be achieved.

Material point method for granular flow segregation
Like many other MPM simulations of two phase flow , two sets of materials could be used for binary mixture granular flow, with each set represents a separate size specie.However, taking advantage of the rheology of binary granular mixture can make the MPM simulation simpler.Even though components of the mixture differ in size, they obey the same rheology rule, i.e.Eq. ( 5); what's more, the flow velocity profiles of species in binary mixture of granular materials are the same to that of the mixture, and the mixture velocity profile is determined by applying specie average density and particle size ρ = 2 i=1 φ i ρ i , d = 2 i=1 φ i d i in Eq. ( 4) [15,16], then the space profile of σ d yy can be determined by Eq. ( 4).Since concentration profile φ i only depends on σ d yy (Eq.( 6)), segregation evolution could be achieved.Thus, in this work, just one set of material points is used to represent the mixture, and material parameters for the mixture are averaged from binary species.Note that the ordinary MPM is space-averaged over each Euler background grid cell, while segregation leads to space gradient of concentration inside each grid cell, thus calculating segregation with the ordinary MPM may lead to numerical viscus on the spacial evolution of specie concentration.To avoid the effect of numerical viscus on segregation result, a modified MPM frame with a extra Euler grid is used, granular flow is simulated on the original MPM grids, just like the ordinary MPM, while evolution of specie concentration is solved on the extra Euler grid nodes, as shown in figure 1.

Results and Discussion
An equal-volume mixture of two size granite particles (2 mm and 3 mm) flowing on an open channel is simulated with the combine model Eqs.( 3), ( 4) and ( 6) by MPM, and results are compared with DEM simulation.In both MPM and DEM simulations, the channel is 22mm wide, with a slope angle of 26 • , periodic boundaries are used in the front and rear of the channel, and the channel length is 88 mm, bottom is set as non-slip boundary.Sidewalls are vertical to the bottom incline.In the DEM simulation, 1841 large particles and 5168 small particles are mixed randomly and then deposited on the channel by gravity, 310 particles are glued to the lower incline to form a nonslip boundary.Frictional sidewalls are used, friction coefficient between sidewall and particles is μ w = 0.3, be- tween particles is μ p = 0.4.Elastic modulus of particle is 29 GPa, particle density is 2650 kg/m 3 [7].In MPM simulation, 7200 material points are used to represent the mixture, the spacing between points is 0.25mm; other material parameters A = 1.79 × 10 9 , q φ = 0.38, a = 0.46, b = 1.0, φ = 0.6 are used for flow dynamic calculation [10], and η = 1.11 × 10 −5 m 2 s −2 , c D = 6.8s −1 s −1 , b k,s = 0.39 for segregation calculation [7].Simulation results are shown in figure 2. As we can see, both the concentration profiles of small particle (figure 2-a, 2-b) and the mixture velocity profile (figure 2-c, 2-d) of MPM simulation can match well with those of DEM; small particles are transferred from top to bottom during the flow.Due to the frictional sidewall, particles near the bottom are almost still during the flow process, they form the erodible bed and change the mixture velocity profile from Bagnold type to exponent type profile [14,20]; accordingly, there is almost no segregation near the bottom, concentration of small particle in this area keeps the original value of 0.5.Simulation results suggest that, the combine model Eq. ( 6) can represent the flow dynamics as well as segregation evolution of binary size granular mixture flow well, and material point method turns out to be a reliable numerical scheme to simulate granular flow.We also notice some deviations of velocity profile between DEM simulation and our model in figure 2-d.These deviation may be caused by the erodible bed near the bottom.Actually the erodible bed plays an important role on the rheology of granular flow [21], and it also have a significant effect on segregation results [22].At the stage, our ongoing work is focusing on the relation between erodible bed and granular segregation.By applying kinetic theory and mixture theory in granular flows with different sidewall conditions, and comparing the segregation prediction between these two types of model, more hints about the relation between erodible bed, granular mixture rheology and segregation results may emerge.

Figure 1 .
Figure 1.Calculation process of Material Point Method for binary mixture granular flow.(a) discrete granular mixture; (b) continuum granular mixture; (c) divide continuum granular mixture into material points (MP), background grid nodes (GN) are used as assists to solve field equations.A calculation cycle contains two parts, the flow calculation (the ordinary MPM [19]) (d)∼(h) and the segregation calculation (i)∼(k).(d) mapping MP information to GN; (e) solving momentum equation on GN; (f) mapping information back MP from GN; (g) updating MP information (h) resetting deformed grid.Segregation calculation begin when flow calculation is at step (f).(i) mapping specie concentration and mixture dynamic information from MP to extra Euler grid node; (j)transfer specie concentration between extra Euler grid nodes according to Eq. (6); (k) specie concentration are mapped back to MP.

Figure 2 .
Figure 2. Results for MPM (red line) and DEM (black circle) at 1s (left) and 20s (right) of the simulation.(a), (b) are concentration profile of small particles; (c), (d) are velocity profiles of the mixture.