Document Actions

IJEGE-11_BS-Minatti-&-Pasculli

background image
Italian Journal of Engineering Geology and Environment - Book www.ijege.uniroma1.it © 2011 Casa Editrice Università La Sapienza
467
DOI: 10.4408/IJEGE.2011-03.B-052
SPH NUMERICAL APPROACH IN MODELLING 2D MUDDY DEBRIS FLOW
L. MINATTI
(*)
& A. PASCULLI
(**)
(*)
University of Florence, Department of Civil and Environmental Engineering, Italy
(**)
University G. d’Annunzio Chieti-Pescara, Department of Science, Italy
alii, 2001). The most widely used approaches to model
debris flows are usually related either to “continuum”
or to “granular” mechanics. A combination of them can
sometimes be used, depending on the type of debris
flow under study. The content of this paper belongs to
the former kind of approach where the Navier-Stokes
equations are considered as the master tool. In the
“continuum” framework, a quite general mathematical
model, that can be used to model a physical phenom-
enon like the one discussed in this paper, would consist
of two sets of equations: a first set for the liquid phase
and a second one for the solid phase. Interactions terms
can then be added to the equations in order to model
the forces exchanged by the two phases. A model of
this kind can be derived from the mixture theory (a
t
-
kin
et alii, 1976): an implementation of such theory
in the case of debris flows can be found in i
veRson
(1997), where the possibility of using different rheo-
logical models for each phase is suggested.
The grain size distribution plays an important role
in the physics of debris flows and recent developments
in the study of water-solid mixtures have shown that
the clay fraction can be of great importance too, as
it influences the granular interactions. In particular, it
is reported by l
aiGle
et alii (1997) that debris flow
mixtures with a high clay fraction, despite retaining
their two-phase nature, behave like viscoplastic fluids
and that a single phase model of such kind can repro-
duce their physics to a good degree of approximation.
We then restrict our analysis to the specific case of
ABSTRACT
2D muddy debris flow has been simulated accord-
ing to a dam break like problem along a slope. The two
sets of equations related to the fluid and solid phases,
as considered by the debris flow mixture theory, have
been simplified in only one set of equations, consider-
ing just one equivalent material. Then the Herschel-
Bulkley fluid constitutive equations have been select-
ed. The correct parameters of the H
eRsCHel
-b
ulkley
model have been chosen in order to correctly simulate
the behaviour of mudflows. The final mathematical
model, has been solved numerically with the smoothed
particle hydrodynamics (SPH) method. SPH is a par-
ticle mesh-free Lagrangian method, well suitable for
computing highly transitory free surface flows of com-
plex fluids in complex geometries. Finally a laboratory
experimental test has been selected for comparison.
Satisfactory results have been achieved. Nevertheless,
further parametric analyses will be carried out and fur-
ther considerations about both constitutive equations
and numerical improvements will be employed and
discussed in future papers.
K
ey
w
ords
: SPH, 2D numerical modelling, muddy debris flow
INTRODUCTION
Debris flows are characterized by a mixture of wa-
ter and poorly sorted granular material flowing under
the effect of gravity (P
ieRson
et alii. 1987; t
akaHasHi
,
1991; C
oussot
et alii, 1996; i
veRson
, 1997; C
Rosta
et
background image
L. MINATTI & A. PASCULLI
468
5th International Conference on Debris-Flow Hazards Mitigation: Mechanics, Prediction and Assessment Padua, Italy - 14-17 June 2011
where:
η is the local dynamic viscosity of the fluid;
τ
c
is the yield stress;
k is called the liquid consistency;
n is called the power law index;
D is the strain rate tensor of eqn. (4).
τ
c
, k, and n represent the constitutive law param-
eters: their values will be adjusted according to experi-
ments (see Tab. 1) in order to correctly reproduce the be-
haviour of a mudflow. The physical meaning of the yield
stress is immediate, representing the stress threshold
below which the fluid starts to behave like a rigid body.
NUMERICAL MODELLING
The physics and the complexity of debris flow phe-
nomena require the selection of an appropriate numerical
approach. The related bibliography is very extensive and,
for the sake of completeness, we mention the fields of
Computational Fluid Dynamics (CFD) and of Computa-
tional Granular Dynamics (CGD). CFD (C
HunG
, 2006),
among many others) and CGD (P
ösCHel
& s
CHwaGeR
,
2005) are considered, respectively, as numerical tools
related to continuum mechanics and granular mechan-
ics. Molecular Dynamics is a further numerical tool
used to study the interactions among particles, exploit-
ing essentially Newton law related to both their transla-
tions and rotations. Nevertheless, according to the scope
of the present paper, an appropriate numerical method,
among the ones suitable for solving the Navier-Stokes
equations, should be selected. In these regards, many ap-
proaches have been proposed: in particular, techniques
capable of handling free surfaces and violent mass fluxes
debris flows characterized by a high clay fraction
(mudflows) and simplify the problem by employing a
single equivalent phase model of a viscoplastic fluid.
GOVERNING EQUATIONS
In light of the considerations reported in the intro-
duction, we are therefore interested in the equations of
motion of an incompressible non-Newtonian fluid. In
order to develop a theoretical framework for a numeri-
cal model we resort to the basic principles governing
the motion of a continuum, namely the mass conserva-
tion and the momentum equations. They can be written
in Lagrangian form as:
where:
n
is the local velocity of the continuum
ρ is the local density of the continuum
f is the body force per unit of mass exerted on the con-
tinuum;
σ is the local total stress tensor;
The total stress tensor is usually split into two
parts: an isotropic and a deviatory one. The stress ten-
sor decomposition is indicated as follows in the paper:
where:
p is the isotropic pressure;
I is the unit tensor;
τ is the deviatory part of the total stress tensor;
The isotropic part of the tensor depends on the
pressure, while the deviatory part of the total stress
tensor is usually expressed as a function – namely the
constitutive law - of the strain rate tensor D:
The constitutive law used throughout the paper to
simulate the behaviour of a mudflow is the Herschel-
Bulkley law, which has also been recently used by
l
aiGle
et alii (2007):
(1)
(2)
(3)
(4)
(5)
(6)
Tab. 1 - Rheological characteristics and experimental pa-
rameters in l
AiGle
et alii (1994) and l
AiGle
(1997)
background image
SPH NUMERICAL APPROACH IN MODELLING 2D MUDDY DEBRIS FLOW
Italian Journal of Engineering Geology and Environment - Book www.ijege.uniroma1.it © 2011 Casa Editrice Università La Sapienza
469
method that is easy to code. Its meshless and Lagrang-
ian nature make it very attractive for solving fluid flow
problems where free surface boundary conditions and
large strain rates are involved. The computational do-
main is filled with particles carrying flow field infor-
mation (e.g. pressure, velocity, density) and capable of
moving in space.
Particles are the computational frame used in the
method to solve the flow describing PDEs, as a grid or
a mesh to calculate spatial derivatives is not needed.
We shall refer to 2D cases throughout the rest of
the paper, even though all the assumptions and results
can be extended to a 3D case with little effort. The key
idea on which the method is based is the well-known
use of a convolution integral with a Dirac delta func-
tion to reproduce a generic function f(x):
In the SPH method, the Dirac function is replaced
by a “bellshaped” kernel function w (it ‘mimics’ the
Dirac delta function), and the generic function f(x)
is reproduced with a convolution integral, which in
a discrete framework takes the form of a summation
over particles:
where:
x
i
and x
j
represent the i and j particle positions in the
given frame of reference;
∆A
j
represents the tributary area (or volume in a 3D
case) associated with particle j;
Summation is extended to all the particles located
within the support domain of particle i.
The kernel function is chosen to be non negative,
even and with a support domain Ω
x
(usually circu-
lar) whose radius is a multiple of a length h, named
smoothing length. The kernel function is zero outside
the support domain and the smoothing length serves
as a scaling parameter for its arguments. It also has
the property of converging to the Dirac function as the
smoothing length approaches to zero.
It is possible to obtain the expression for the SPH
approximation of a function gradient by writing the
convolution integral of the function with the kernel
and by using the Gauss-Green formula:
appear to be suitable for solving muddy debris flows
problems.
Finite Elements (FEM), Finite Volumes (FVM),
Finite Differences (FDM) are some of the most com-
mon used methods (C
HunG
, 2006). For all the men-
tioned approaches, grid generation is an important
step to be performed. And as far as grid generation is
concerned, the possibilities are numerous: structured
grids, unstructured grids and adaptive meshes.
Other important issues related to the numerical
solutions of fluid flows problems modelled with the
Navier-Stokes equations involve incompressibility and
convective terms. In order to avoid the latter one, a La-
grangian approach is commonly selected instead of an
Eulerian one.
Furthermore, particular care should be taken con-
cerning the numerical issues related to large deforma-
tions. This is the main reason why a frequent (time
consuming) update of the grid, aimed at lowering the
excessive mesh distortion due to large deformations,
is often necessary.
From the above discussion, it is clear that meth-
ods that avoid the numerical instability due to the
convective terms, lower the grid generation time and
that are capable of easily taking free surfaces into ac-
count, are very desirable. The “Meshless Finite Ele-
ment Methods” include many of the above necessary
discussed features (see i
delsoHn
et alii, 2003; 2004;
for an overview). The “Smoothed Particle Hydrody-
namics” (SPH), briefly described in the following par-
agraphs and selected in our approach is among these.
A similar tool, the “Particle Finite Element Method”
(PFEM) (o
ñate
et alii, 2004) is very promising as it
includes all the advantages of the FEM without some
of the difficulties of the SPH method (in particular, the
imposition of boundary conditions), even if it requires
some time, (which can still be reduced by using an
extended Delaunay method), to update grid meshes.
Finally, it is worth to mention other approaches,
reliable for granular debris flow, like the “Distinct
Elements” and the “Alternate Lagrangian Eulerian”
method (ALE) (C
Rosta
et alii, 2001).
OVERVIEW OF THE SPH METHOD
The Smoothed Particle Hydrodynamic (SPH)
method is a numerical technique that was initially de-
veloped during the 1970s to solve astrophysical prob-
lems (m
onaGHan
, 2005). It is a fully meshless particle
(7)
(8)
background image
L. MINATTI & A. PASCULLI
470
5th International Conference on Debris-Flow Hazards Mitigation: Mechanics, Prediction and Assessment Padua, Italy - 14-17 June 2011
Particles positions are evolved in time, according
to the velocity of each particle.
SPH DISCRETIZATION TECHNIQUE
MAIN EQUATIONS
The equation briefly presented in the overview
paragraph can be used to discretize the Herschel-
Bulkley fluid flow governing equations in order to
create an SPH algorithm.
There are many references where it is possible to
find details on how fluid flow governing equations can
be effectively discretized into SPH equations. Among
the many of them, we indicate m
onaGHan
(2005), l
iu

et alii (2003) & t
akeda
et alii (1994).
The SPH discretization of the mass conservation
equation (1) used in the paper is as follows:
where:
m represents the particle mass;
v
ij
= v
i
- v
j
represents the difference between the inter-
acting particles velocity;
A widely used SPH discretization of the momen-
tum equation (2) is as follows:
where:
r
ij
= r
i
- r
j
represents the difference between the inter-
acting particle position;
η
ij
is a symmetrised dynamic viscosity between inter-
acting particles, such that η
ij
= η
ji
.
Usually, a symmetrised expression for the
smoothing lengths of each pair of interacting particles
in the summations is used when calculating the kernel
gradient (e.g. hij=(hi+hj)/2). If this is the case, it can
be shown that particles exchange equal and opposite
forces, making eqn. (12) capable of conserving the
linear momentum of the particles system.
In the paper, the following expression has been
used for the symmetrised viscosity:
where:
the kernel is differentiated with respect to the x’ co-
ordinate;
n represents the normal to the support domain bound-
aries (pointing outward).
The first term of the RHS of eqn. (9) is zero if the
support domain isn’t truncated by the computational
domain boundaries, as the kernel is zero on the sup-
port domain boundaries. Another case when the term
can be zero is when the support domain is truncated
by the computational domain boundaries but there ex-
ists a boundary condition forcing the function f(x) to
vanish on the boundaries (it may be the case when f(x)
represents a velocity and a no-slip condition has to be
enforced on the computational domain boundaries). If
the first term of the RHS of eqn. (9) is zero, then the
SPH approximation of f(x) gradient takes the follow-
ing form in a discrete framework:
Eqn. (10) is often used, even when the first term
of the RHS of eqn. (9) doesn’t vanish. Nevertheless,
it is possible to offset the errors induced by neglect-
ing the term by introducing special treatments of the
boundaries, as detailed later.
There are consistency conditions that the kernel
has to satisfy in order to correctly reproduce continu-
ous field functions up to a certain order of accuracy.
They are related to the kernel moments and should
be satisfied both in the continuous and in the discrete
frameworks. More details can be found in l
iu
et alii
(2003) & in l
iu
et alii (2006).
Numerical applications of SPH share these com-
mon features: The equations describing the continuum
motion are written in Lagrangian form, by using La-
grangian instead of Eulerian time derivatives.
The spatial gradients involved in the PDEs are
discretized over the particles, by using a SPH ap-
proximation.
Time integration of variables is performed par-
ticle-wise usually by using an explicit time stepping
method.
(9)
(10)
(11)
(12)
(13)
background image
SPH NUMERICAL APPROACH IN MODELLING 2D MUDDY DEBRIS FLOW
Italian Journal of Engineering Geology and Environment - Book www.ijege.uniroma1.it © 2011 Casa Editrice Università La Sapienza
471
It can be shown that density variation is propor-
tional to the square of the Mach number (m
onaGHan
(1994)). If an artificial equation of state where sound
speed is larger than the bulk velocity of flow is used,
then it is possible to keep density variations as low
as desired. The artificial equation of state used in this
paper is as follows:
where:
c is the artificial speed of sound;
ρ
0
is the reference density of the fluid at zero pressure;
There are other possible forms for the equation of
state that can be used to simulate a weakly compressible
flow [see m
onaGHan
(1994)]. Nevertheless, we found
that the choice of the equation of state doesn’t signifi-
cantly affect the results, as long as the artificial speed of
sound is chosen to be large enough.
VISCOSITY REGULARIZATION
The Herschel-Bulkley constitutive law provides
a viscosity diverging to infinite for strain rates ap-
proaching to zero. It is therefore impossible to numer-
ically reproduce this behaviour in a straightforward
manner. In the calculations, the expression of eqn. (5)
for viscosity has been regularized and the infinite val-
ue that is obtained when the strain rate approaches to
zero has been replaced with a finite value. The expres-
sion used to regularize viscosity has been proposed by
P
aPanastasiou
(1987):
The parameter B is related to the maximum vis-
cosity that is returned by the regularization when the
strain rate is zero. In this case the maximum viscosity
value is given by:
The value assigned to the parameter B has been
chosen high enough to correctly simulate the physi-
cal behaviour of the fluid but small enough to avoid
prohibitively small time steps.
BOUNDARY TREATMENT
Since the gradient formulation given by eqn. (10)
is used even for fluid particles close to the boundaries,
The term in the right hand side of eqn. (12) involv-
ing the symmetrised viscosity represents the deviatory
stress tensor divergence. Such an expression has been
proposed by C
leaRy
(1998) & m
onaGHan
(2005) and
accounts for the presence of spatial gradients in viscos-
ity too. It is possible to find a proof in Espanol et al.
(2003) that it accounts also for a bulk viscosity coeffi-
cient, which is equal to 5/3 times the dynamic shear vis-
cosity. Despite not having any bulk viscosity in the con-
stitutive law that is being used in the paper, we still think
it could be worth using the expression in eqn. (12) for
the deviatory tensor divergence for two main reasons:
It conserves the angular momentum in the parti-
cles system as discussed in m
onaGHan
(1997), thus
leading to more accurate results in the simulations;
It provides zero-th order consistency. This means
that there are no viscous interactions between particles
if the velocity field is constant in space;
The bulk viscosity induced by the expression
should not significantly affect the simulation results,
as the velocity field in the weakly compressible fluid
approximation is nearly divergence-free;
A possible alternative to eqn. (12) is the one of
directly calculating the total stress tensor from eqn.
(10) and to use then the following SPH approximation
for its divergence:
This formulation is still able to preserve the linear
momentum of the particles system and is able to ex-
actly reproduce the Herschel- Bulkley rheology within
the accuracy of the SPH approximation, without add-
ing any kind of bulk viscosity to the flow. Neverthe-
less, it lacks some good properties of the formulation
of eqn. (12) like the zero-th order consistency and the
angular momentum conservation
ARTIFICIAL COMPRESSIBILITY
The simulation of an incompressible flow requires
the solution of a Poisson equation for the pressure,
which often leads to an increase of the computational
time. Therefore, it is more practical to approximate
the uncompressible medium with a weakly compress-
ible one, thus allowing the calculation of the pressure
from the density with a stiff equation of state which
introduces an artificial compressibility in the fluid.
(14)
(15)
(16)
(17)
background image
L. MINATTI & A. PASCULLI
472
5th International Conference on Debris-Flow Hazards Mitigation: Mechanics, Prediction and Assessment Padua, Italy - 14-17 June 2011
special treatments are required in order to offset the
errors introduced by neglecting the terms in the right
hand side of eqn. (9).
There are a number of strategies to tackle bounda-
ry-related problems. R
andles
et alii (1996) used ghost
particles to treat a symmetrical surface boundary con-
dition. Ghost particles have also been used in various
manners for particle approximations near boundaries
by t
akeda
et alii (1994), m
oRRis
et alii (1997) and f
eR
-
RaRi
et alii (2009) by using point symmetry.
In this paper, boundaries have been treated by
placing a layer of particles on them, whose spacing
is three times narrower than the fluid particles one.
Boundary particles prevent fluid particles from pass-
ing through the domain boundaries by exerting a nor-
mal force on them. They also interact with the fluid
particles via SPH summations through their viscosity,
thus enforcing the no-slip condition. The expression
used for the boundary forces has been proposed by
m
onaGHan
et alii (2009):
where:
k
b
is a constant having the dimensions of a square ve-
locity, used to correctly reproduce the bulk forces ex-
erted by the boundaries on the fluid;
Summation is extended to all the boundary particles lo-
cated within the support domain of fluid particle i.
It can be shown (m
onaGHan
et alii, 2009), that the
above summation gives negligible contribution along
the direction parallel to the boundary, making boundary
forces being exerted only along the normal direction. It
is also a symmetrical formulation, thus conserving the
linear momentum of the particles system.
TIME INTEGRATION
Time integration has been performed by means
of a symplectic Verlet scheme, as in k
aJtaR
et alii
(2008). The time stepping scheme is explicit and con-
serves the linear momentum of the particles system.
The time step Δt is controlled by a C.F.L. condition
depending on the artificial speed of sound, the viscous
interactions between particles and on the interactions
with boundary particles, according to the following
equation:
where:
The minimum time step value is sought over each
couple of interacting particle;
h
ij
=(h
i
+h
j
)/2 is a symmetrised smoothing length be-
tween the pair of interacting particles;
c
ij
=(c
i
+c
j
)/2 is a symmetrised artificial speed of sound
between the pair of interacting particles;
ρ
ij
=
i
j
)/2 is a symmetrised density between the pair
of interacting particles;
dp is the initial fluid particles spacing;
β is the ratio between the boundary and the fluid par-
ticles spacing;
The C.F.L. number, indicated as CFL, has been set
equal to 0.5.
MODEL TESTING
The ability of the SPH model to correctly repro-
duce a mudflow has been tested by simulating the
experiments performed by l
aiGle
et alii (1994) and
(1997). Their experiments consisted in creating a
mudflow dam break problem in a laboratory flume, by
quickly opening a gate. The experimental setup they
used is briefly shown in the figure below:
After the opening of the gate, the material stored
behind it was released and the three ultrasonic gages,
sketched in the picture, recorded the mudflow wave
heights in time. The authors used waterclay mixtures
prepared in laboratory with different concentrations in
order to simulate mudflows. Herschel-Bulkley rheo-
logical parameters for the used mixtures were fitted
to measures, carried out with a rheometer by the au-
thors. They performed tests with four different kinds
of mixtures, named A, B, C and D. The rheological
(18)
(19)
Fig. 1 - 3D sketch of the experimental setup in l
AiGle
et
alii (1994) and l
AiGle
(1997)
background image
SPH NUMERICAL APPROACH IN MODELLING 2D MUDDY DEBRIS FLOW
Italian Journal of Engineering Geology and Environment - Book www.ijege.uniroma1.it © 2011 Casa Editrice Università La Sapienza
473
cal dissipation but increases the computational time, as
the number of particles required to fill the entire com-
putational domain gets higher. As the physics of the
flow is dominated by viscosity, the choice of the SPH
discretization for the momentum equation plays an im-
portant role. We have tried both eqn. (12) and (14) with
different particles spacing to reproduce the experiment
of Tab.1. Accuracy is improved at higher resolutions for
both equations but even though eqn. (12) seems to pro-
vide more accurate and less noisy results, it seems also
to be more dissipative than eqn. (14). The influence of
resolution and the issues regarding what would be the
most appropriate choice for the momentum equation
will be discussed further in a future paper
Here we present the results obtained by using eqn.
(14) with a particles spacing of dp = 3 mm. The fol-
lowing figure shows the initial particle disposition,
right before the opening of the gate:
In the following figures the particle disposition is
reported in two different instants of the simulation:
We then extrapolated the wave heights below the
same locations of the ultrasonic gages used in the ex-
periments of l
aiGle
et alii. (1994) and (1997) from
the numerical results of our simulation by measuring
the thickness of the layer of particles passing below.
Finally, a comparison of our numerical wave heights
characteristics of the A mixture and the experimental
parameters are summarized in the following table:
The authors modelled the mudflow as a wide chan-
nel flow initially at rest: the resulting flow conditions
were therefore completely determined by the height
of the material stored behind the gate, the flume slope
and the material characteristics. l
aiGle
et alii (1994)
and (1997) indicated two non dimensional scaling
parameters controlling the resulting flow conditions.
The values of the scaling parameters in the tests they
performed can represent a wide range of real situa-
tions at smaller scales (C
oussot
, 1994). In particular,
the material A, along with the experimental conditions
of Tab. 1 represent a realistic natural material (with ρ
0
= 2200 kg/m
3
, τ
c
= 900 Pa, k = 290 Pa·s
1/3
) in a 120
m long slope (l
aiGle
et alii, 1997). A field determina-
tion of the Herschel-Bulkley rheological parameters
of a natural debris flow can be found in C
oussot
et
alii (1998), even though the values provided for the
parameters in the paper are in a slightly different range
than the ones indicated above. l
aiGle
et alii (1994)
and (1997) developed a 1D numerical model based
on the shallow water equations and solved it with a
Godunov conservative finite volumes scheme. We
tested our SPH model by simulating the experiment of
l
aiGle
et alii (1997) using the rheological parameters
of the Herschel-Bulkley model for the material A and
the experimental conditions indicated in Tab. 1.
In our SPH simulations, particles were stored be-
hind the gate and a layer of boundary particles was
placed to simulate the closed gate at the initial time
step. Particles were initially stored in a uniform lat-
tice according to their spacing and were given a
ρ
0
=
1410 kg/m
3
density. At first, around 3000 damped time
steps, as in m
onaGHan
(1994), were performed in or-
der to build up a hydrostatic pressure distribution. Af-
ter that, the gate opening was simulated by removing
the layer of boundary particles placed on it, thus re-
leasing the mudflow. The artificial Mach number was
set to 0.2, the boundary force constant of eqn. (18)
was set to k
b
= g·H
g
(where g = gravity acceleration,
H
g
= height of material behind the gate) and the coef-
ficient B of eqn. (17) for the viscosity regularization
was set to 10 s.
Another important parameter is the initial parti-
cles spacing dp, which plays the same role of the grid
spacing in finite differences schemes. A decrease in the
spacing improves the accuracy and reduces the numeri-
Fig. 2 - Initial particles disposition (lengths are in meters
and particles are colour coded according to their
pressure in Pa)
Fig. 3 - Particles disposition after 0.15 s from gate open-
ing (lengths are in meters and particles are colour
coded according to the module of their velocity in
m/s)
background image
L. MINATTI & A. PASCULLI
474
5th International Conference on Debris-Flow Hazards Mitigation: Mechanics, Prediction and Assessment Padua, Italy - 14-17 June 2011
break experiments made by l
aiGle
et alii (1994) and
(1997). The SPH model is fully two-dimensional and
is capable of providing more information than the
one-dimensional model. Nevertheless, a preliminary
analysis shows that it is too diffusive, while still being
able to produce a good agreement with experimental
results. While the SPH model numerical diffusivity is
reduced by increasing the resolution, other kinds of
improvements haven’t been fully tested yet, but they
seem possible by using different SPH versions of the
momentum equation. In further papers, parametric
studies with comparison with actual debris flow phe-
nomena will be discussed.
with the experimental measurements of l
aiGle
et alii
(1994) and (1997) is reported in Fig. 5:
It can be observed that while the wave heights are
well caught, the SPH model is still affected by some
numerical diffusion, as the calculated wave’s velocity
is lower than in the experimental data. This can be
improved either by further increasing the resolution or
by switching to a less diffusive SPH formulation for
the momentum equation.
CONCLUDING REMARKS
A SPH model incorporating the Herschel-Bulk-
ley non Newtonian rheology has been developed
and has been tested to simulate some mudflow dam
Fig. 4 - Particles disposition after 0.60 s from gate open-
ing (lengths are in meters and particles are colour
coded according to the module of their velocity in
m/s)
Fig. 5 - Time plot of the simulated wave heights at the lo-
cations of gage 1 (red line), gage 2 (blue line) and
gage 3 (black line). Experimental measurement s
of gages are represented by dots coloured as the
respective lines
REFERENCES
a
tkin
R.J. & C
Raine
R.e. (1976) – Continuum theories of mixtures: basic theory and historical development. Q. J. Mech. Appl.
29(2): 209-244.
C
HunG
t.J. (2006) – Computational Fluid Dynamics. Cambridge University Press, New York.
C
leaRy
P.w. (1998) – Modelling confined multi-material heat and mass flows using SPH. Appl. Math. Modelling 22: 981-993.
C
oussot
P. (1994) – Steady, laminar, flow of concentrated mud suspensions in open channel. Journal of Hydraulic Research
32(4): 535-559.
C
oussot
P. & m
eunieR
m. (1996) – Recognition, classification and mechanical description of debris flows. Earth-Sci. Rev. 40:
209-227.
C
oussot
P., l
aiGle
d., a
Rattano
m., d
eGanutti
a. & m
aRCHi
l. (1998) – Direct determination of rheological characteristics of
debris flow. Journal of Hydraulic Engineering 124(8): 865-868.
C
Rosta
G.b., C
alvetti
f., i
mPosimato
s., R
oddeman
d., f
Rattini
P. & a
GliaRdi
f. (2001) – Granular Flows and Numerical Mod-
eling of Landslides. Debrisfall Assessment in Mountain Catchments for local end-users Contract N. EVG1-CT-1999-00007.
e
sPanol
P. & R
evenGa
m. (2003) – Smoothed dissipative particle dynamics. Physical Review 67: 026705.
f
eRRaRi
a., d
umbseR
m. & t
oRo
e.f. (2009) - A new 3D parallel scheme for free surface flows. Computers & Fluids 38: 1203-
1217.
i
veRson
R.m. (1997) – The physics of debris flows. Reviews of Geophysics 35(3): 245-296.
k
aJtaR
J.b. & m
onaGHan
J.J. (2008) – SPH simulations of swimming linked bodies. Journal of Computational Physics 227:
8568-8587.
i
delsoHn
s.R., o
ñate
e., C
alvo
n. & d
el
P
in
f. (2003) – The meshless finite element method. Int. J. Numer. Meth. Engng. 58:
893-912.
i
delsoHn
s.R., o
ñate
e. & d
el
P
in
f. (2004) – The Particle Finite Element Method: a powerful tool to solve incompressible
background image
SPH NUMERICAL APPROACH IN MODELLING 2D MUDDY DEBRIS FLOW
Italian Journal of Engineering Geology and Environment - Book www.ijege.uniroma1.it © 2011 Casa Editrice Università La Sapienza
475
flows with freesurfaces and breaking waves. Int. J. Numer. Meth. Engng. 61: 964-989.
l
aiGle
d. & C
oussot
P. (1994) – Modélisation numérique des écoulements de laves torrentielles. La Houille Blanche 3: 50-56.
Laigle D. & Coussot P. (1997) – Numerical modelling of mudflows. Journal of Hydraulic Engineering 123(7): 617-623.
l
aiGle
d., l
aCHamP
P. & n
aaim
m. (2007) - SPH-based numerical investigation of mudflow and other complex fluid flow interac-
tions with structures. Comput. Geosci. 11: 297-306.
l
iu
G.R. & l
iu
m.b. (2003) - Smoothed Particle Hydrodynamics a meshfree particle method. World Scientific Publishing.
Liu M.B. & Liu G.R. (2006) - Restoring particle consistency in smoothed particle hydrodynamics. Applied Numerical Math-
ematics 56: 19-36.
m
onaGHan
J.J. (1994) - Simulating free surface flows with SPH. Journal of Computational Physics 110: 399-406.
m
onaGHan
J.J. (1997) – SPH and Riemann solvers. Journal of Computational Physics 136: 298-307.
m
onaGHan
J.J. (2005) - Smoothed particle hydrodynamics. Reports on Progress in Physics 68: 1703-1759.
m
onaGHan
J.J. & k
aJtaR
J.b. (2009) – SPH particle boundary forces for arbitrary boundaries. Computer Physics Communica-
tions 180: 1811- 1820.
m
oRRis
J.P., f
ox
P.J. & z
Hu
y. (1997) - Modeling low Reynolds number incompressible flows using SPH. Journal of Computa-
tional Physics 136: 214-226.
o
ñate
e., i
delsoHn
s.R., d
el
P
ino
f. & a
ubRy
R. (2004) – The Particle Finite Element Method. An Overview. International
Journal. of Computational Methods 1(2): 267-307.
P
aPanastasiou
t.C. (1987) – Flows of materials with yield. Journal of Rheology 31(5): 385.
P
ieRson
t.C. & C
osta
J.e. (1987) – A rheologic classification of subaerial sediment-water flows. Rev. Eng. Geol., Geol.Soc.
Am., Boulder, Co., 7: 1-12
P
ösCHel
t. & s
CHwaGeR
t. (2005) – Computational Granular Dynamics. Springer Berlin Heidelberg, New York
R
andles
P.w. & l
ibeRsky
l.d. (1996) - Smoothed particle hydrodynamics: some recent improvements and applications. Com-
puter Methods in Applied Mechanics and Engineering 138: 375-408.
t
akaHasHi
t. (1991) – Debris Flow. IAHR Monograph, published for IAHR by A.A.Balkema, Rotterdam
t
akeda
H., m
iyama
s.m. & s
ekiya
m. (1994) – Numerical simulation of viscous flow by smoothed particle hydrodynamics.
Progress of Theoretical Physics 92: 939-960.
Statistics