Report of Smoothed Particle Hydrodynamics Method

Introduction

Smoothed-particle hydrodynamics (SPH) is a meshfree Lagrangian numerical method developed by Gingold and Monaghan and Lucy in 1977, initially used to solve astrophysical problems. In its thirty years of existence, SPH method extend to many branches of computational physics, it has been introduced to simulate the mechanics of solid mechanics and fluidflows.

Regarding the field of geomechanics, large deformation and post-failure are essential topics, and its numerical prediction is important for engineering practice and design. In the past, the most used methods to solve geomechanics problem are finite element method (FEM) and the discrete element method (DEM). FEM has difficulties on large deformation and post-failure, while DEM cannot handle small scale stimulation. However, SPH offers a good solution to deal with these issues. The feature, application, and the advantages and deficiencies of the method in comparison with the other numerical methods are discussed in this report.

Get Help With Your Essay
If you need assistance with writing your essay, our professional essay writing service is here to help!
Essay Writing Service

Main Feature of SPH

The SPH method is a Lagrangian meshfree method which can carry field variables, including mass, density, stress tensor, etc. Different from the particle in cell method (PIC), SPH does not require a grid to calculate spatial derivatives. Instead, SPH is defined by the partial differential. It ruled by the law of conservation of energy, mass and momentum, in addition with the constitutive relations of material in the problem. Through SPH, the partial differential equations for the continuum are converted into equations of motion of these particles and then solved by updated Lagrangian numerical scheme.

Besides, SPH is so-called an initial-boundary value problem due to the need of additional information. The required information consists of the initial values for all variables at all positions at time zero, and boundary values for all variables at the boundary positions. The objective of SPH is to spread the information to all positions and times of interest.

The fundamental of SPH is an interpolation method. The integral interpolant of any function A(r) is defined by
Ar=∫Ar‘Wr–r‘,hdr‘
.

The integration is over the whole space, and W is an interpolating kernel. Noted that,
∫Wr–r‘,hdr‘=1

And
limh→0⁡Wr–r‘,h=δ(r–r‘)

For numerical work, the integral interpolant is stimulated by a summation interpolant
Asr=∑bmbAbρbWr–rb,h
,

where index b is a particle label with mass
mb
, position
rb
, density
ρb
and velocity
vb
. The value of quantity A at
rb
is
Ab
. It is important to noted that we can write a differentiable interpolant function from interpolation points by using a differentiable kernel. By ordinary differentiation, derivatives of this interpolant can be obtained. Therefore, there is not necessary to use finite differences and a grid.

Noted that, there is a first golden rule of SPH: assume the kernel is a Gaussian, if we want to find a physical interpretation of an SPH equation. The second golden rule of SPH is related to rewrite formulae with the density placed inside operators. In the initial case, we write
∇∙v=[∇∙ρv–v∙∇ρ]/ρ

Hence the divergence of the velocity at particle a is
ρa(∇∙v)a=∑b(vb–va)∙∇aWab

where
∇aWab
is the gradient of
Wr–rb,h
taken with respect to the coordinates of particle a. Now, we assume the kernel is a Gaussian, the contribution from particle b to the divergence of the velocity at particle a is
2mb(va–vb)∙(ra–rb)Wab/h2

It is noted that there is a positive contribution when the particles are moving away from each other. And the vorticity of particle a is
ρa(∇×v)a=∑bmbvab×∇aWab

It is realised that the contribution of particle b to the vorticity of particle a is proportional to the relative angular momentum per unit mass of the two particles.

Main applications of the method in Geotechnical Engineering

Desiccation induced cracks in soils

The modelling of desiccation cracking in soils involves the water loss because of moisture evaporation, the behaviour of unsaturated soils undergoing thermo-hydro-mechanical processes and large soil deformation and crack development. The direct interaction between soil and atmosphere further exacerbates the complexity of this problem. Because of these intricacy, only insufficient work has been performed to quantitatively simulate large deformations and crack development in clay undergoing thermo-hydro-mechanical processes. Models for simulating shrinkage cracks in soil with atmospheric flow processes have not published. Therefore, the understanding of drying shrinkage and related cracking in soils, as well as methods to control or avoid such cracking is still ambiguous. In this case, SPH is applied to discern the property of crack initiation and propagation in clay. Moisture evaporation, heat and mass exchange in multi-phase system and prediction of overall soil behaviour is examined by SPH.

Slope stability analysis and slope failure simulation

Limit equilibrium method (LEM) and finite element method (FEM) are applied to slope stability analysis. However, slope instability is related to a discontinues of soil, which cannot be examined by either LEM or FEM. To address this problem, SPH is developed to access the slope stability and the post failure behaviour of soil. To consider the pore water pressure, a SPH formulae is employed for soil motion. In addition, the equation can be used for saturated soil. There is a good agreement between SPH and other methods regarding safety factor and critical slip surface. It is noted that SPH can deal with large deformation and post failure of soil, hence it can apply in computational geomechanics which involve large deformation and failure of geomaterials.

Soil-structure interaction modelling

SPH has been employed to model the soil-structure interaction. Herein, combine the Drucker Prager yield criterion with elastoplastic stress-strain relationship for the soil model. While the concrete structure is an isotopically elastic-perfectly plastic material associated with the Von-Misses yield criterion. Numerical results have shown a great work of SPH in modelling large deformation and failure of soil and concrete structure.

Seepage failure/erosion modelling

Slope or river embankments failures are result of soil erosion or soil boiling problems, which were caused by seepage flow. In the past, FEM combined with finite difference method (FDM) had been used to investigate the failure mechanism. However, FEM is limited to handle large deformation and failure of geomaterials. To cope with this, a two-stage SPH soil-water coupling model is established, based on the mixing of SPH frames for geomaterials and SPH frames for incompressible water streams.

The main advantages and deficiencies of the method in comparison with the other numerical methods used in the field of Geotechnical Engineering such as the Finite Element Method.

Advantages

In the past, the finite element method (FEM) and the discrete element method (DEM) are often used solve geotechnical problems. Nevertheless, dealing with the large deformation and post-failure is challenging for FEM because FEM is a grid-based method. Thus, FEM gets grid distortions, which cause errors in the solution or even fail to do the computation because of the negative values of Jacobian determinants at nodes of numerical integration.

DEM simulation unlike FEM, it does not have those difficulties. Yet, processor power limits DEM to small-scale simulation with a few hundred thousand of particles. DEM is unable to handle large scale problems. Furthermore, specification of DEM parameters is slightly cryptic, and there are no reliable guidelines yet.

Currently, the SPH method has been established and applied to simulate large deformation of continuum or dispersed material. Compared with other grid-based numerical methods, SPH has several advantages such as it can handle large deformation and post-failure very well due to its Lagrangian and adaptive nature; complex free surfaces are modeled naturally without any special treatments; complex geometries can be handled without any difficulties and it is relatively easy to incorporate complicated physics.

In comparison with other grid-based numerical methods, SPH has the benefits of Lagrangian and adaptive, which can deal with large deformation and post failure well. The complex free surface is naturally modelled without any special treatment. Besides, complex geometries can be applied with ease, and it is accessible to combine with complex physics.

In general, SPH is a versatile and uncomplicated method for numerical fluid dynamics, which an active research topic in geotechnical engineering. Compared with other numerical methods, SPH offers a significant dynamic range in spatial resolution and density. Besides SPH has a great conservation property for angular momentum, which is not automatically guaranteed in Eulerian codes., SPH accurately conserves the total energy when dealing with self-gravity. This is not obvious in most mesh-based methods. Last but not least, SPH is Galilean invariant, without any errors from advection, which is another benefit compared to Eulerian mesh-based method.

Deficiencies

SPH need a very massive number of particles to take complex interactions in a large domain. Compared with conventional Eulerian codes, Lagrangian code has a poor performance on the computational efficiency and parallel scalability, and is incapable of accurately examing sharp phenomena such as shocks and interphase boundaries. 

Secondly, in SPH, setting boundary conditions like inlets, outlets and walls is much challenging than with grid-based approaches. A study is concluded that “the treatment of boundary conditions is certainly one of the most difficult technical points of the SPH method”. It is tough to set boundary condition in SPH as the particles near the boundary change with time.

Finally, for instance, considering kinetic energy spectrum, if the metric of interest is indirectly associated with density, the computational cost of each particle’s SPH simulation is significantly greater than the cost of each cell’s grid-based simulation. Therefore, ignoring the problem of parallel acceleration, using a grid-based approach to simulate a constant density flow such as external aerodynamics is more efficient than using SPH.

References

Bui, H. H., Fukagawa, R., Sako, K., & Ohno, S. (2008). Lagrangian meshfree particles method (SPH) for large deformation and failure flows of geomaterial using elastic–plastic soil constitutive model. International Journal for Numerical and Analytical Methods in Geomechanics, 32(12), 1537-1570.

H.H. Bui et al. (2015). Soil cracking modelling using the mesh-free SPH method, 12th Australian New Zealand Conference on Geomechanics (ANZ12), 22-25 February, Wellington, New Zealand.

Ha H. Bui, K. Sako , R. Fukagawa (2009): Slope stability analysis and slope failure simulation by SPH, The 17th International Conference for Soil Mechanics and Geotechnical Engineering(ICSMGE), Vol.2. pp.1578-1581, Egypt, Oct., 05-09, 2009.

Ha H. Bui, R. Fukagawa, K. Sako & J.C. Wells (2011): Slope stability analysis and discontinuous slope failure simulation by elasto-plastic smoothed particle hydrodynamics (SPH), Géotechnique, Vol.61. pp.565 -574, 2011 

Ha H. Bui, K. Sako, R. Fukagawa, and J.C. Wells: SPH-based numerical solutions for large deformation of geomaterial considering soil-structure interaction, The 12th International Association for Computer Methods and Advances in Geomechanics (IACMAG), pp.570-578, Goa, India, Oct., 01-06, 2008.

Bui H.H and Nguyen G.D. (2017). A coupled fluid-solid SPH approach to modelling flow through deformable porous media, International Journal of Solids and Structures, Vol.125, pp.224-246 (doi: 10.1016/j.ijsolstr.2017.06.022), ISI-Q1 Journal.