A KINEMATIC MODEL FOR THE VARIATION IN POROSITY IN DRY GRANULAR FLOWS

(1)

A KINEMATIC MODEL FOR THE VARIATION IN POROSITY IN DRY GRANULAR FLOWS

Chung Fang^1*

ABSTRACT

A gradient-flow constitutive theory for dry granular matters, in which the variation in porosity is modeled by using a kinematic evolution equation, is proposed on the basis of thermodynamic consistency. The model is applied to study the behavior of dry granular flows down an inclined moving plate, of which the results are compared with experimental outcomes. Results show that the present model is able to describe the velocity and porosity fields of dry granular flows with slow to moderate speed. For very rapid flows like avalanches, the variation in porosity cannot be described by using a gradient-flow theory, and should be simulated by using high-order models of the variation in porosity.

Key Words: Granular matter, Gradient-flow theory, Gravity-driven flow

INTRODUCTION

Dry granular flows are conventionally treated as elasto-visco-plastic granular continua with microstructural effects (Duran, 2000; Hutter and Rajagopal, 1994; Savage, 1993; Wang and Hutter, 2001). Among the various microstructural effects is the variation in porosity ε is the most important one (Fang et al, 2006a,b; Kirchner, 2002; Svendsen and Hutter, 1995; Svendsen et al, 1999). It indicates the effects of the redistribution of the pore space on the behavior of granular matters at macroscale. Such an effect is conventionally represented in literature by another quantity, namely the volume fraction ν defined as the solid volume divided by the total volume of a representative volume element (RVE) viz.,

1 . (1)

ε = − ν

Instead of the porosity, we should use the concept of volume fraction to discuss the effects of the pore space in the forthcoming analysis.

The microstructural effects induced by the variation of the volume fraction is commonly accounted for by introducing a scalar balance law (Wang and Hutter, 2001). However, different authors do not unanimously agree upon the form of this scalar balance equation, and the typical proposals are summarized in the following two categories:

1 Department of Civil Engineering, National Cheng Kung University, No.1, University Road, Tainan City 701, Taiwan (*Corresponding Author; Tel: +886-6-275-7575 # 32142, Fax: +886-6-235-8542, Email:

cfang@mail.ncku.edu.tw)

(2)

( )

div , (Svendsen and Hutter,1995), div div , (Wilmanski, 1996),

div div , (Bluhm et al., 1995),

div , (Goodman and Cowin, 1971), div , (Fang et al., 2006a,b).

I f

f f

II f

f ν + ν =

ν + ν = +

γ + γ = +

γναν = + γν γν ν = + γν

&

&&

&

l v

v h

h h

(2)

In the category (I) a “kinematic” equation is proposed in which the first time derivative of ν is balanced by its flux and production for a RVE. Profound models in this category are given in (2)1-3, in which γ is the true mass density of the grains, ν again the volume fraction, v the velocity and h and f the flux and production associated with the variation of ν (or γ), respectively, which are considered constitutive quantities. While in (2)1,2 the variation of ν is directly accounted for by proposing a balance equation for its first time derivative, it is taken into consideration in (2)3

implicitly by proposing a balance equation for the first time derivative of γ, which is not really different from (2)2 in view of the decomposition of the bulk density ρ = γν. Since in this category only the first time derivative of ν is employed, no power of working associated with h and f emerges, and the energy equation is thus not altered and possesses its traditional form.

In the category (II) a “dynamic” equation is postulated in which the second time derivative of ν is employed and balanced by its flux and production. Typical models are e.g. (2)4 and (2)5, where α is the equilibrated inertia, h and f the equilibrated stress vector (the associated flux) and the equilibrated intrinsic body force (the associated production) of the variation of ν, respectively, and the internal length. In fact, the equation (2)5 is the revised version of (2)4: it removed the dimensional inconsistencies in (2)4, in particular in the quantity α, by introducing an internal length . Four different considerations regarding the roles played by _land the corresponding thermodynamic analyses were performed and discussed (Fang et al., 2006a); the numerical simulations of Benchmark flow problems were also carried out (Fang et al., 2006c). Since in this category the variation in ν is modeled by using its second time derivative, additional powers of working associated with its flux and production emerge and should be taken into account in the energy equation.

l l

Since for dry granular slow flows the collisions between the grains are relatively insignificant, a weak variation in porosity is approved. Thus, in the present study we will use the Wilmanski`s model to study the behavior of an isothermal dry granular slow flow with incompressible grains down an inclined moving plate, of which the results are compared with experimental outcomes, to estimate the validity of the gradient-flow theory. To this end, the derived constitutive models based on the two concepts in modeling the variation in porosity, namely (2)2 and (2)5, will be outlined in the next Section. In Sect. 3 we will use these two constitutive models to study the behavior of an isothermal dry granular flow with incompressible grains down an inclined moving plate, of which the results will be compared with the experimental outcomes, follow which the validity of the two concepts in modeling the variation in porosity can be estimated. The paper will be summarized in Sec. 4.

(3)

CONSTITUTIVE MODEL

The balance equations of isotropic dry granular flows are given by (Fang et al., 2006a,b)

mass 0 div , (3)

linear momentum div ,

angular momentum ,

internal friction 0= , ,

energy 0 div ,

entropy 0 div ,

T

e r

s

= γν + γν + γν

= γν − − γν

= −

−⎡⎣ ⎤⎦−

= γν − ⋅ + − γν

= γνη + − γν − π

& &

&

v

0 v t b

0 t t

Z ΩZ Φ t D q

φ

(4) (5) (6) (7) (8)

with the abbreviations

, = − ,

⎡ ⎤

⎣Ω Z⎦ ΩZ ZΩ (9)

where t is the Cauchy stress tensor, b the specific body force, t^T the transposition of t, e the specific internal energy, D the symmetric part of the velocity gradient, known as the stretching tensor, q the heat flux, r the specific energy supply, Z an Euclidean frame- indifferent, second- rank symmetric tensor (a spatial internal variable) describing the frictional and non-conservative forces inside a RVE (Svendsen and Hutter, 1995; Pitteri, 1986), Ω any orthogonal rotation of the RVE, Φ a tensor-valued constitutive relation for the production of Z,² η the specific entropy, φ the entropy flux, s the specific entropy supply and π the entropy production. In addition, the notation A．B denotes that A．B = tr (AB^T)= tr (A^TB) for two arbitrary second-rank tensors.

Equations (3), (4) and (8) are the traditional balances of mass, linear momentum and entropy, respectively, except that the bulk density ρ is decomposed into ρ = γν. Since the material is not considered micropolar or Cosserat-type and the effects of particle rotation are excluded, the balance of angular momentum reduces to its simplest form (5), namely the symmetry of the Cauchy stress tensor. Furthermore, (7) is the traditional balance of internal energy (energy equation). To account for the effects of plasticity, the internal friction and other non-conservative forces inside a granular microcontinuum is represented by Z of which the time evolution is modeled by the assigned equation (6), in which Z&− ⎡⎣Ω,Z⎤⎦is the so-called “corrotational” objective time derivative of Z. It reduces to the Jaumann derivative of Z when Ω is chosen to W, the skew- symmetric part of velocity gradient. Since in (6) the time evolution of Z is described by its first time derivative, Z does not represent any real frictional forces on the surface of the grains. It is rather an abstract ideal dealing with such forces, which per se yields no dissipation, and a “true”

variational principle for its energy contributions does not exist. As a result, the internal energy balance (7) remains unchanged, although Z is considered an independent field quantity.

Since the number of the unknowns correspond to the number of the available equations,³ the system (3)-(8) is likely mathematically well-posed, and can can obtain the values of the

2 Here we apply the concept of microcontinuum proposed by Mindlin, 1964, and regard the RVE as a deformable “granular microcontinuum”. In this respect, t is a constitutive quantity of the granular microcontinuum, whilst Z describes the internal friction and other non-conservative forces inside the granular microcontinuum which cannot be seen from the perspective outside the microcontinuum.

3It is seen that the balance of angular momentum (6) is automatically fulfilled by the constitutive models, the balance of entropy (8) is in fact an inequality and is used in the thermodynamic analysis, and the balance of energy (7) is not relevant for isothermal flows.

(4)

unknowns by integrating the equations simultaneously, provided that the constitutive equations are prescribed. From the thermodynamic analysis based on the Mueller-Liu entropy principle and the quasi-linear theory, the constitutive equations for isothermal, isochoric dry granular flows with incompressible grains are given by (Fang, 2009)

( ) 0 ² ⁸

2

1 2 3

, 0, ,

2

( ),

E

m

D

s

f p

tr II

f

ν

∞

= = β = − θλ

⎛ ν ⎞

= −νβ + λ + εν + μ γ ⎜ ⎟ ν − ν

⎝ ⎠

+ ς + ς + ς

&

h 0

t D I

I Z Z

D (10)

with the scalar parameters ς1-ς3 are scalar parameters. In deriving (10) and (11) the Cayley- Hamilton theorem has been used. For details see Fang, 2009.

NUMERICAL SIMULATIONS WITH EXPERIMENTAL COMPARISONS

Consider a steady, fully developed, isothermal dry granular flow with incompressible grains down an inclined moving plane with the coordinate system shown in Fig. 1. The flow thickness is L with a free upper surface, and the lower solid plane is moving with a constant velocity V0 in the direction oppositional to the flow direction. The motion is driven by the interaction between the gravity and the motion of the lower solid plane.

Fig. 1. Gravity-driven dry granular slow flow down an inclined moving plane and the coordinate system

It is assumed, for the considered flow, that parallel flow prevails, namely, ( ) , 0 , 0 , ( ), ( ), sin , cosθ, 0 ,

u y y p p y b b

⎡ ⎤

=⎣ ⎦ ν = ν = = −⎡⎣ θ ⎤⎦

v b (11)

hold. Substituting (11) into (3)-(4) and (10) gives rise to the dimensionless field equations in the forms⁴

(11)

( )

2

1 1

8

2

2 2

2 1

0 s

1

0 2 1 1 cos ,

1

m s m

d du

dy dy S

d S

dy

∞

⎧ ⎛ ⎞ − ν ν ⎫

⎪ ⎪

= ⎨⎪⎩ ν − ν ⎝⎜ ⎟⎠ + Ξ − ν ν⎬⎪⎭−

⎧ − ν ν ⎫

⎪ ⎪

= ⎨⎪⎩ ν ν − − Ξ − ν ν⎬⎪⎭+ ν θ in , ν θ

(12)

4In deriving equations (12) it is noted that the considered steady flow corresponds to the stationary state in hypoplastic literature (Kirchner, 2002). In this state the equation (6) is decoupled from other balance equations and can be solved separately to obtain a solution for Z_xx, Z_xy and Z_yy by using parallel flow assumption. In addition, an explicit form for f_s can also be obtained in stationary state. For details see Fang, 2009; Bauer and Herle, 2000.

(5)

with the dimensionless parameters

3

1 2 2 2

0 0 0 0 0

2 2

1 2 3 2 3

2 3 1 2

0 0

, , , , ,

, , , ,

( )

, ,

b s

m m m

m yy

m m

yy yy xy xy

m

y y L

u bL bL t

u S S

V V ³

2

0 m

Z Z Z Z

V

∞ ∞

ν ν ν ν

= ν = ν = ν = ν =

ν ν ν ν

ν −

= = = π =

μ γ α ν α γν

ζ + ζ + ζ ζ + ζ

Ξ = Ξ =

α γν μ γ

L

(13)

where νs and νb are the minimum volume fraction and the volume fraction on the solid boundary, respectively, subjected to the non-slip boundary conditions in dimensionless forms⁵

(0) _b, (0) 1 , du( ) 0 .

u L

ν = ν = − dy = (14)

(12) and (14) define the BVP in the present study, in which S2 denotes to some extent the combined effects of gravity and flow thickness, S1 represents the influences of viscosity under a fixed value of S2, and both Ξ1 and Ξ2 denote the effects of hypoplastic-related forces. For the implementations of the numerical simulations, the values of νb, ν∞ and νs need be specified.

Following previous works (Wang and Hutter, 1999a,b), the values νb = 0.51, νm = 0.555, ν∞ = 0.644 and νs = 0.25 are chosen which are appropriate for natural angular beach sand with particle diameters from 0.318 to 0.414 mm. It follows immediately that

0.919, 1.16, 0.451. (15)

b ∞ s

ν = ν = ν =

In fact, analytical solutions of the emerging ODEs (12) with the boundary conditions (14) exist.

Integrating (12)2 with the condition (14)1 gives rise to the analytical solution of ν in an implicit form

( )

( ) ( )

( )( )

2 2

2

0 3 4 cos

1 ln 1

1 1 1

b b

m b m b

m s m

m m b m b

= ν − ν − ν − ν +S y θ

⎛ ν ν − ν ⎡ − ν ν ν⎤⎞

− Ξ − ν ν ν ⎜⎜⎝ − ν ν − ν ν + ⎢⎣ − ν ν ν ⎥⎦⎟⎟⎠, (16)

where ln(a) stands for the Napierian logarithm of a. With (16), an exact solution for u can be obtained by integrating (12)1, viz.,

( )⁴ 1 1 1 ^{1/ 2} 2

1 sin 1 ,

2 1

m s m

u _∞ ⎛C S dy − ν ν ⎞

=

∫

ν − ν ⎜⎜⎝ + θ ν

∫

− Ξ − ν ν⎟⎟⎠ ^dy+^C

(17)

where C1 and C2 are two integration constants which can be determined by using (14)2,3. In the following we will display some profiles of the velocity and porosity by changing the values of S1, S2, Ξ1 andΞ1, and compare these results with experimental outcomes. For convenience, the over- bar will no longer be used to distinguish dimensionless quantities.

5Discussions of boundary conditions for granular flows can be found e.g. in Fang, 2009.

(6)

NUMERICAL RESULTS

Since the considered flow is driven by the interaction between the gravity and the resistance (reflected by the hypoplastic and viscous effects), the numerical simulations should be performed for the variations of these key factors to illuminate their physical influences. In fact, the dimensionless parameters S1 and S2 are not independent of each other, and S1 can be expressed as S1 = AS2, where A=α0νm3L²/μ0γV02, which is essentially determined in experiments when the physical properties of the grains are known. Numerical tests have shown that only the relative magnitudes of the ν- and u-profiles are influenced by the values of A, but their tendencies remain unchanged. Thus, for simplicity, A is chosen to be a constant in the numerical simulations to match the experimental outcomes from Perng et al., 2006. In addition, from previous works it is seen that the parameters Ξ1 and Ξ2, which denote the hypoplastic effects, are of equal weightings and importance, and thus the identity Ξ1 = Ξ2 is used in the numerical simulations. As a parametric study, the numerical results for the variations of S2 and Ξ1 are shown in Figs. 2-4, in which the horizontal axes denote the values of 1-νand u, while the vertical axes represent the distance y measured from the solid plane. All calculations are carried out for θ = 15.6° to match the experimental setup from Perng et al., 2006.

Fig. 2 illustrates the typical profiles of 1-νand u calculated, in which S2 = 0.01 and Ξ1 = Ξ2 =0.01.

The profiles of 1-νare in fact the profiles of the dimensionless porosity by using (1). Since the detailed information about the values of, νm, νs and ν∞ of the grains used in Perng et al., 2006 is not clear, the focus is thus on the tendencies displayed by the calculated profiles of 1-νand u, but not their absolute values. It is seen that the velocity u decreases monotonically from its boundary value toward the free surface, exhibiting a rather linear profile across the fluid layer (Fig. 2(b)).

The profile of 1-ν behaves like an “exponential” function: it increases from its minimum value on the solid plane toward the free surface, and the increasing rate becomes larger when approaching the free surface (Fig. 2(a)). This indicates that during the flow most grains are confined within the regions near the moving solid plane due to the gravitational and frictional effects, while the grains near the free surface are colliding strongly with one another, resulting in a larger value of porosity near the free surface. To estimate the performance and limitations of the present model, the experimental results quoted from Perng et al., 2006 are displayed in Figs. 2(c) and 2(d) for the profiles of 1-ν and u, respectively, in which the value of Vδ- Vδ(u0 = 0) in the horizontal axis of Fig. 2(c) corresponds to the value of 1-ν, namely the dimensionless porosity, while the value of u/u0 in the horizontal axis of Fig. 2(d) is simply the dimensionless velocity u in our notation system. It is evident that the present model delivers appropriate profile for the velocity when compared with experiments, and the profile of porosity of the porosity also matches well with the experimental outcomes. This is due to the fact that for

(7)

Fig. 2. Typical calculated profiles of 1-ν and u, in which S2 = 0.01 and Ξ1 = Ξ2 =0.01. (a)-(b): calculated results; (c)- (d): experimental results quoted from Perng et al. 2006.

Fig. 3. Typical calculated profiles of 1-ν and u for the variations of S2 (=0.001, 0.005, 0.01), in which Ξ1 = Ξ2 =0.01.

(8)

Fig. 4. Typical calculated profiles of 1-ν and u for the variations of Ξ1 and Ξ2 (=0.001, 0.005, 0.01), in which S2 =0.01.

granular flows with slow to moderate velocities, frictional effects among the grains are dominate.

This results in turn in an evolution of porosity with a “kinematic” nature, and can better be taken into account by using the concept of gradient flow theory, as shown in (2)2. These conclusions hold equally when compared with other experimental results like Pudasaini and Hutter, 2007;

GDR MiDi, 2004.

Calculated profiles of 1-ν and u for the variations of S2 are shown in Fig. 3, in which Ξ1 = Ξ2

=0.01. In view of (13)8, S2 denotes the combined effects of the gravity and fluid layer thickness.

As shown in Figs. 3 (b), as S2 increases, the free surface velocities decrease correspondingly, and the velocity profiles exhibit the same tendencies as before: an almost linear profile across the fluid layer. In addition, as S2 increases, the grains gather themselves more efficiently in the regions near the moving solid plane for granular flows with slow to moderate velocities, resulting in the profiles of 1-νwith larger amplitude toward the free surface, as shown in Fig. 3(a).

Calculations have also been performed for the variations of Ξ1 and Ξ2, and the results are illustrated in Fig. 4, in which S2 = 0.01. The parameters Ξ1 and Ξ2 indicate the hypoplastic effects, namely, the effects of plastic deformations. When plastic deformations occur, the granular body becomes easier to deform continuously (to flow). For granular flows with slow to moderate velocities, plastic deformations occur in almost the entire fluid layer, and are most significant near the free surface. Such a tendency is revealed in Fig. 4(b), in which the grains near the free surface moves at larger velocities in the direction opposite to the moving solid plane as Ξ1 and Ξ2

increase. This in turn gives rise to stronger collisions among the grains in these regions, resulting in larger values of porosity, as shown in Fig. 4(a).

CONCLUDING REMARKS

In the present study, a gradient-flow theory was proposed to account for the effects induced by the variation in porosity in dry granular flows. In order to estimate the model validity, the theory was applied the study the behavior of an isothermal dry granular slow flow down an inclined moving plate, of which the results were compared with the experimental outcomes. Results show

(9)

that the present model delivers a quasi-linear velocity profile across the fluid layer with a decreasing tendency from the boundary values on the moving solid plane toward the free surface, which is in good agreement with experiments. Moreover, the obtained profile of porosity increases from its minimum value on the moving solid plane toward the free surface “quasi- exponentially”, which is also in good agreement with experimental outcomes. This is due to the facts that for granular flows with slow to moderate velocities the collisions between the grains are relatively in significant, and thus the distribution of the porosity can better be described by using a gradient-flow theory. On the contrary, for very rapid flows like avalanches, the collisions become much stronger, and under such a circumstance the distribution of the porosity should then be described by using a dynamic model.

ACKNOWLEDGEMENTS

The author is indebted to the National Science Council, Taiwan, for the financial support through the research grant NSC 96-2221-E-006-078-MY3.

REFERENCES

Bauer E, Herle I (2000), “Stationary states in hypoplasticity”, in: Constitutive modeling of granular materials, (ed. D. Kolymbas), Springer-Verlag Berlin Heidelberg New York, 167-192 Bluhm J, de Bohr R, Wilmanski K (1995), “The thermodynamic structure of the two-component

model of porous incompressible materials with true mass densities”, Mech. Res. Comm. 22, 171-180

Duran J (2000), “Sands, Powders and Grains”, Springer-Verlag Wien New York

Fang C, Wang Y, Hutter K (2006a), “A thermo-mechanical continuum theory with internal length for cohesionless granular materials. Part I. A class of constitutive models”, Continuum.

Mech. Thermodyn. 17(8), 545-576

Fang C, Wang Y, Hutter K (2006b), “Shearing flows of a dry granular material - hypoplastic constitutive theory and numerical simulations”, Int. J. Numer. Anal. Meth. Geomech. 30, 1409—1437

Fang C, Wang Y, Hutter K (2006c), „A thermo-mechanical continuum theory with internal length for cohesionless granular materials. Part II. Non-equilibrium postulates and numerical simulations of simple shear, plane Poiseuille and gravity driven problems”, Continuum. Mech.

Thermodyn. 17(8), 577-607

Fang C (2009), “Gravity-driven dry granular slow flows down an inclined moving plane: a comparative study between two concepts of the evolution of porosity”, Rheol. Acta (in press) GDR MiDi (2004), „On dense granular flows“, Eur. Phys. J. E 14, 41–365

Goodman MA, Cowin SC (1971), “Two problems in the gravity flow of granular materials”, J.

Fluid Mech. 45, 321-339

Hutter K, Rajagopal KR (1994), “On flows of granular materials”, Continuum Mech. Thermodyn.

6, 81-139

Kirchner N (2002), “Thermodynamically consistent modeling of abrasive granular materials. I:

Non-equilibrium theory”, Proc. R. Soc. Lond. A 458, 2153-2176

Mindlin RD (1964), “Microstructure in linear elasticity”, Arch. Rat. Mech. Anal. 16, 51-78 Pitteri M (1986), “Continuum equations of balance in classical statistical mechanics”, Arch. Rat.

(10)

Mech. Anal. 94, 291-305

Perng ATH, Capart H, Chou HT (2006), “Granular configurations, motions, and correlations in slow uniform flows driven by an inclined conveyor belt”, Gran. Matt. 8, 5-17

Pudasaini S, Hutter K (2007), „Avalanche Dynamics“, Springer-Verlag Berlin Heidelberg

Savage SB (1993), “Mechanics of Granular Flows”, in: Continuum mechanics in environmental sciences and geophysics (Hutter, K. ed.), Springer-Verlag Heidelberg, 467-522

Svendsen B, Hutter K (1995), “On the thermodynamics of a mixture of isotropic materials with constraints”, Int. J. Engng. Sci. 33, Nr. 14, 2021-2054

Svendsen B, Hutter K, Laloui L (1999), „Constitutive models for granular materials including quasi-static frictional behavior: toward a thermodynamic theory of plasticity”, Continuum Mech. Thermodyn. 4, 263-275

Wang Y, Hutter K (2001), “Granular material theories revisited”, in: Geomorphological 2001.

Fluid Mechanics (Balmforth, NJ and Provenzale, A. eds.), Springer-Verlag Heidelberg, 79-107 Wang Y, Hutter K (1999a), “A constitutive theory of fluid-saturated granular materials and its

application in gravitational flows”, Rheol. Acta. 38, 214-223

Wang Y, Hutter K (1999b), “Shearing flows in a Goodman-Cowin type granular material-theory and numerical results”, Particulate Science and Technology 17, 97-124

Wilmanski K (1996), “Porous Media at Finite Strains. The new model with the balance equation of porosity”, Arch. Mech. 48, Nr. 4, 591-628