Smoothed-particle hydrodynamics explained
Smoothed-particle hydrodynamics (SPH) is a computational method used for simulating the mechanics of continuum media, such as solid mechanics and fluid flows. It was developed by Gingold and Monaghan[1] and Lucy[2] in 1977, initially for astrophysical problems. It has been used in many fields of research, including astrophysics, ballistics, volcanology, and oceanography. It is a meshfree Lagrangian method (where the co-ordinates move with the fluid), and the resolution of the method can easily be adjusted with respect to variables such as density.
Method
Advantages
- By construction, SPH is a meshfree method, which makes it ideally suited to simulate problems dominated by complex boundary dynamics, like free surface flows, or large boundary displacement.
- The lack of a mesh significantly simplifies the model implementation and its parallelization, even for many-core architectures.[3] [4]
- SPH can be easily extended to a wide variety of fields, and hybridized with some other models, as discussed in Modelling Physics.
- As discussed in section on weakly compressible SPH, the method has great conservation features.
- The computational cost of SPH simulations per number of particles is significantly less than the cost of grid-based simulations per number of cells when the metric of interest is related to fluid density (e.g., the probability density function of density fluctuations).[5] This is the case because in SPH the resolution is put where the matter is.
Limitations
- Setting boundary conditions in SPH such as inlets and outlets[6] and walls[7] is more difficult than with grid-based methods. In fact, it has been stated that "the treatment of boundary conditions is certainly one of the most difficult technical points of the SPH method".[8] This challenge is partly because in SPH the particles near the boundary change with time.[9] Nonetheless, wall boundary conditions for SPH are available. [10]
- The computational cost of SPH simulations per number of particles is significantly larger than the cost of grid-based simulations per number of cells when the metric of interest is not (directly) related to density (e.g., the kinetic-energy spectrum). Therefore, overlooking issues of parallel speedup, the simulation of constant-density flows (e.g., external aerodynamics) is more efficient with grid-based methods than with SPH.
Examples
Fluid dynamics
Smoothed-particle hydrodynamics is being increasingly used to model fluid motion as well. This is due to several benefits over traditional grid-based techniques. First, SPH guarantees conservation of mass without extra computation since the particles themselves represent mass. Second, SPH computes pressure from weighted contributions of neighboring particles rather than by solving linear systems of equations. Finally, unlike grid-based techniques, which must track fluid boundaries, SPH creates a free surface for two-phase interacting fluids directly since the particles represent the denser fluid (usually water) and empty space represents the lighter fluid (usually air). For these reasons, it is possible to simulate fluid motion using SPH in real time. However, both grid-based and SPH techniques still require the generation of renderable free surface geometry using a polygonization technique such as metaballs and marching cubes, point splatting, or 'carpet' visualization. For gas dynamics it is more appropriate to use the kernel function itself to produce a rendering of gas column density (e.g., as done in the SPLASH visualisation package).
One drawback over grid-based techniques is the need for large numbers of particles to produce simulations of equivalent resolution. In the typical implementation of both uniform grids and SPH particle techniques, many voxels or particles will be used to fill water volumes that are never rendered. However, accuracy can be significantly higher with sophisticated grid-based techniques, especially those coupled with particle methods (such as particle level sets), since it is easier to enforce the incompressibility condition in these systems. SPH for fluid simulation is being used increasingly in real-time animation and games where accuracy is not as critical as interactivity.
Recent work in SPH for fluid simulation has increased performance, accuracy, and areas of application:
- B. Solenthaler, 2009, develops Predictive-Corrective SPH (PCISPH) to allow for better incompressibility constraints[11]
- M. Ihmsen et al., 2010, introduce boundary handling and adaptive time-stepping for PCISPH for accurate rigid body interactions[12]
- K. Bodin et al., 2011, replace the standard equation of state pressure with a density constraint and apply a variational time integrator[13]
- R. Hoetzlein, 2012, develops efficient GPU-based SPH for large scenes in Fluids v.3[14]
- N. Akinci et al., 2012, introduce a versatile boundary handling and two-way SPH-rigid coupling technique that is completely based on hydrodynamic forces; the approach is applicable to different types of SPH solvers[15]
- M. Macklin et al., 2013 simulates incompressible flows inside the Position Based Dynamics framework, for bigger timesteps[16]
- N. Akinci et al., 2013, introduce a versatile surface tension and two-way fluid-solid adhesion technique that allows simulating a variety of interesting physical effects that are observed in reality[17]
- J. Kyle and E. Terrell, 2013, apply SPH to Full-Film Lubrication[18]
- A. Mahdavi and N. Talebbeydokhti, 2015, propose a hybrid algorithm for implementation of solid boundary condition and simulate flow over a sharp crested weir[19]
- S. Tavakkol et al., 2016, develop curvSPH, which makes the horizontal and vertical size of particles independent and generates uniform mass distribution along curved boundaries[20]
- W. Kostorz and A. Esmail-Yakas, 2020, propose a general, efficient and simple method for evaluating normalization factors near piecewise-planar boundaries[10]
- Colagrossi et al., 2019, study flow around a cylinder close to a free-surface and compare with other techniques[21]
Astrophysics
Smoothed-particle hydrodynamics's adaptive resolution, numerical conservation of physically conserved quantities, and ability to simulate phenomena covering many orders of magnitude make it ideal for computations in theoretical astrophysics.[22]
Simulations of galaxy formation, star formation, stellar collisions,[23] supernovae[24] and meteor impacts are some of the wide variety of astrophysical and cosmological uses of this method.
SPH is used to model hydrodynamic flows, including possible effects of gravity. Incorporating other astrophysical processes which may be important, such as radiative transfer and magnetic fields is an active area of research in the astronomical community, and has had some limited success.[25] [26]
Solid mechanics
Libersky and Petschek[27] [28] extended SPH to Solid Mechanics. The main advantage of SPH in this application is the possibility of dealing with larger local distortion than grid-based methods.This feature has been exploited in many applications in Solid Mechanics: metal forming, impact, crack growth, fracture, fragmentation, etc.
Another important advantage of meshfree methods in general, and of SPH in particular, is that mesh dependence problems are naturally avoided given the meshfree nature of the method. In particular, mesh alignment is related to problems involving cracks and it is avoided in SPH due to the isotropic support of the kernel functions. However, classical SPH formulations suffer from tensile instabilities[29] and lack of consistency.[30] Over the past years, different corrections have been introduced to improve the accuracy of the SPH solution, leading to the RKPM by Liu et al.[31] Randles and Libersky[32] and Johnson and Beissel[33] tried to solve the consistency problem in their study of impact phenomena.
Dyka et al.[34] [35] and Randles and Libersky[36] introduced the stress-point integration into SPH and Ted Belytschko et al.[37] showed that the stress-point technique removes the instability due to spurious singular modes, while tensile instabilities can be avoided by using a Lagrangian kernel. Many other recent studies can be found in the literature devoted to improve the convergence of the SPH method.
Recent improvements in understanding the convergence and stability of SPH have allowed for more widespread applications in Solid Mechanics. Other examples of applications and developments of the method include:
- Metal forming simulations.[38]
- SPH-based method SPAM (Smoothed Particle Applied Mechanics) for impact fracture in solids by William G. Hoover.[39]
- Modified SPH (SPH/MLSPH) for fracture and fragmentation.[40]
- Taylor-SPH (TSPH) for shock wave propagation in solids.[41]
- Generalized coordinate SPH (GSPH) allocates particles inhomogeneously in the Cartesian coordinate system and arranges them via mapping in a generalized coordinate system in which the particles are aligned at a uniform spacing.[42]
Numerical tools
Interpolations
The Smoothed-Particle Hydrodynamics (SPH) method works by dividing the fluid into a set of discrete moving elements
, referred to as particles. Their Lagrangian nature allows setting their position
by integration of their velocity
as:
| d\boldsymbol{r |
i}{dt}=\boldsymbol{v} |
i.
These particles interact through a kernel function with characteristic radius known as the "smoothing length", typically represented in equations by
. This means that the physical quantity of any particle can be obtained by summing the relevant properties of all the particles that lie within the range of the kernel, the latter being used as a weighting function
. This can be understood in two steps. First an arbitrary field
is written as a convolution with
:
A(\boldsymbol{r})=\intA\left(\boldsymbol{r\prime
}\right) W(| \boldsymbol-\boldsymbol |,h) \, \mathrmV \! \left(\boldsymbol\right).
The error in making the above approximation is order
. Secondly, the integral is approximated using a Riemann summation over the particles:
A(\boldsymbol{r})=\sumjVjAjW(|\boldsymbol{r}-\boldsymbol{r}j|,h),
where the summation over
includes all particles in the simulation.
is the
volume of particle
,
is the value of the quantity
for particle
and
denotes position. For example, the density
of particle
can be expressed as:
\rhoi=\rho(\boldsymbol{r}i)=\sumjmjWij,
where
denotes the particle mass and
the particle density, while
is a short notation for
W(|\boldsymbol{r}i-\boldsymbol{r}j|,h)
. The error done in approximating the integral by a discrete sum depends on
, on the particle size (i.e.
,
being the space dimension), and on the particle arrangement in space. The latter effect is still poorly known.
[43] Kernel functions commonly used include the Gaussian function, the quintic spline and the Wendland
kernel.
[44] The latter two kernels are compactly supported (unlike the Gaussian, where there is a small contribution at any finite distance away), with support proportional to
. This has the advantage of saving computational effort by not including the relatively minor contributions from distant particles.
Although the size of the smoothing length can be fixed in both space and time, this does not take advantage of the full power of SPH. By assigning each particle its own smoothing length and allowing it to vary with time, the resolution of a simulation can be made to automatically adapt itself depending on local conditions. For example, in a very dense region where many particles are close together, the smoothing length can be made relatively short, yielding high spatial resolution. Conversely, in low-density regions where individual particles are far apart and the resolution is low, the smoothing length can be increased, optimising the computation for the regions of interest.
Discretization of governing equations
For particles of constant mass, differentiating the interpolated density
with respect to time yields
=\sumjmj\left(\boldsymbol{v}i-\boldsymbol{v}j\right) ⋅ \nablaWij,
where
is the gradient of
with respect to
. Comparing this equation with the continuity equation in the Lagrangian description (using
material derivatives),
=-\rho\nabla ⋅ \boldsymbol{v},
it is apparent that its right-hand side is an approximation of
; hence one defines a discrete divergence operator as follows:
\operatorname{D}i\left\{\boldsymbol{v}j\right\}=-
\sumjmj\left(\boldsymbol{v}i-\boldsymbol{v}j\right) ⋅ \nablaWij.
This operator gives an SPH approximation of
at the particle
for a given set of particles with given masses
, positions
and velocities
.
The other important equation for a compressible inviscid fluid is the Euler equation for momentum balance:
} = -\frac\nabla p + \boldsymbol
Similarly to continuity, the task is to define a discrete gradient operator in order to write
}_i \left\ + \boldsymbol
One choice is
}_i\left\ = \rho_i \sum_j m_j \left(\frac + \frac\right) \nabla W_,
which has the property of being skew-adjoint with the divergence operator above, in the sense that
\sumiVi\boldsymbol{v}i ⋅ \operatorname{G
}_i \left\ = - \sum_i V_i p_i \operatorname_i\left\,
this being a discrete version of the continuum identity
\int\boldsymbol{v} ⋅ \operatorname{grad}p=-\intp\operatorname{div} ⋅ \boldsymbol{v}.
This property leads to nice conservation properties.[45]
Notice also that this choice leads to a symmetric divergence operator and antisymmetric gradient. Although there are several ways of discretizing the pressure gradient in the Euler equations, the above antisymmetric form is the most acknowledged one. It supports strict conservation of linear and angular momentum. This means that a force that is exerted on particle
by particle
equals the one that is exerted on particle
by particle
including the sign change of the effective direction, thanks to the antisymmetry property
.
Nevertheless, other operators have been proposed, which may perform better numerically or physically.For instance, one drawback of these operators is that while the divergence
is zero-order consistent (i.e. yields zero when applied to a constant vector field), it can be seen that the gradient
} is not. Several techniques have been proposed to circumvent this issue, leading to renormalized operators (see e.g.
[46]).
Variational principle
The above SPH governing equations can be derived from a least action principle, starting from the Lagrangian of a particle system:
l{L}=\sumjmj\left(
| 2 |
\tfrac{1}{2}\boldsymbol{v} | |
| j |
-ej+\boldsymbol{g} ⋅ \boldsymbol{r}j\right)
,
where
is the particle specific
internal energy. The
Euler–Lagrange equation of variational mechanics reads, for each particle:
} = \frac.
When applied to the above Lagrangian, it gives the following momentum equation:
mi
=
-\sumjmj
| \partialej |
\partial\boldsymbol{r |
i}+mi\boldsymbol{g}=
-\sumjmj
| \partial\rhoj |
\partial\boldsymbol{r |
i}+mi\boldsymbol{g}
where the chain rule has been used, since
depends on
, and the latter, on the position of the particles.Using the thermodynamic property
de=\left(p/\rho2\right)d\rho
we may write
mi
=
-\sumjmj
| \partial\rhoj |
\partial\boldsymbol{r |
i}+mi\boldsymbol{g},
Plugging the SPH density interpolation and differentiating explicitly
\tfrac{\partial\rhoj}{\partial\boldsymbol{r}i}
leads to
=-\sumjmj\left(
+
\right)\nablaWij+\boldsymbol{g},
which is the SPH momentum equation already mentioned, where we recognize the
} operator. This explains why linear momentum is conserved, and allows conservation of angular momentum and energy to be conserved as well.
[47] Time integration
From the work done in the 80's and 90's on numerical integration of point-like particles in large accelerators, appropriate time integrators have been developed with accurate conservation properties on the long term; they are called symplectic integrators. The most popular in the SPH literature is the leapfrog scheme, which reads for each particle
:
\begin{align}
&=
+
,\\
&=
+
\Deltat,\\
&=
+
,
\end{align}
where
is the time step, superscripts stand for time iterations while
is the particle acceleration, given by the right-hand side of the momentum equation.
Other symplectic integrators exist (see the reference textbook[48]). It is recommended to use a symplectic (even low-order) scheme instead of a high order non-symplectic scheme, to avoid error accumulation after many iterations.
Integration of density has not been studied extensively (see below for more details).
Symplectic schemes are conservative but explicit, thus their numerical stability requires stability conditions, analogous to the Courant-Friedrichs-Lewy condition (see below).
Boundary techniques
In case the SPH convolution shall be practiced close to a boundary, i.e. closer than, then the integral support is truncated. Indeed, when the convolution is affected by a boundary, the convolution shall be split in 2 integrals,
A(\boldsymbol{r})=\int\Omega(\boldsymbol{r)}A\left(\boldsymbol{r\prime
}\right) W(| \boldsymbol-\boldsymbol |,h) d\boldsymbol + \int_ A\left(\boldsymbol\right) W(| \boldsymbol-\boldsymbol |,h) d\boldsymbol,
where is the compact support ball centered at, with radius, and denotes the part of the compact support inside the computational domain, . Hence, imposing boundary conditions in SPH is completely based on approximating the second integral on the right hand side. The same can be of course applied to the differential operators computation,
\nablaA(\boldsymbol{r})=\int\Omega(\boldsymbol{r)}A\left(\boldsymbol{r\prime
}\right) \nabla W(\boldsymbol-\boldsymbol,h) d\boldsymbol + \int_ A\left(\boldsymbol\right) \nabla W(\boldsymbol-\boldsymbol,h) d\boldsymbol.
Several techniques has been introduced in the past to model boundaries in SPH.
Integral neglect
The most straightforward boundary model is neglecting the integral,
\intB(\boldsymbol{r)-\Omega(\boldsymbol{r})}A\left(\boldsymbol{r\prime
}\right) \nabla W(\boldsymbol-\boldsymbol,h) d\boldsymbol \simeq \boldsymbol,
such that just the bulk interactions are taken into account,
This is a popular approach when free-surface is considered in monophase simulations.[49]
The main benefit of this boundary condition is its obvious simplicity. However, several consistency issues shall be considered when this boundary technique is applied.[49] That's in fact a heavy limitation on its potential applications.
Fluid Extension
Probably the most popular methodology, or at least the most traditional one, to impose boundary conditions in SPH, is Fluid Extension technique. Such technique is based on populating the compact support across the boundary with so-called ghost particles, conveniently imposing their field values.[50]
Along this line, the integral neglect methodology can be considered as a particular case of fluid extensions, where the field,, vanish outside the computational domain.
The main benefit of this methodology is the simplicity, provided that the boundary contribution is computed as part of the bulk interactions. Also, this methodology has been deeply analyzed in the literature.[51] [50] [52]
On the other hand, deploying ghost particles in the truncated domain is not a trivial task, such that modelling complex boundary shapes becomes cumbersome. The 2 most popular approaches to populate the empty domain with ghost particles are Mirrored-Particles[53] and Fixed-Particles.[50]
Boundary Integral
The newest Boundary technique is the Boundary Integral methodology. In this methodology, the empty volume integral is replaced by a surface integral, and a renormalization:
\nablaAi=
\left(
VjAj\nablaWij+
SjAj\boldsymbol{n}jWij\right),
with the normal of the generic j-th boundary element. The surface term can be also solved considering a semi-analytic expression.
Modelling physics
Hydrodynamics
Weakly compressible approach
Another way to determine the density is based on the SPH smoothing operator itself. Therefore, the density is estimated from the particle distribution utilizing the SPH interpolation. To overcome undesired errors at the free surface through kernel truncation, the density formulation can again be integrated in time.[54]
The weakly compressible SPH in fluid dynamics is based on the discretization of the Navier–Stokes equations or Euler equations for compressible fluids. To close the system, an appropriate equation of state is utilized to link pressure
and density
. Generally, the so-called
Cole equation[55] (sometimes mistakenly referred to as the "
Tait equation") is used in SPH. It reads
p=
\right)\gamma-1\right)+p0,
where
is the reference density and
the
speed of sound. For water,
is commonly used. The background pressure
is added to avoid negative pressure values.
Real nearly incompressible fluids such as water are characterized by very high speeds of sound of the order
. Hence, pressure information travels fast compared to the actual bulk flow, which leads to very small Mach numbers
. The momentum equation leads to the following relation:
where
is the density change and
the velocity vector.In practice a value of c smaller than the real one is adopted to avoid time steps too small in the time integration scheme. Generally a numerical speed of sound is adopted such that density variation smaller than 1% are allowed. This is the so-called weak-compressibility assumption.This corresponds to a
Mach number smaller than 0.1, which implies:
where the maximum velocity
needs to be estimated, for e.g. by Torricelli's law or an educated guess. Since only small density variations occur, a linear equation of state can be adopted:
[56] p=
| 2\left(\rho-\rho |
c | |
| 0\right)
|
Usually the weakly-compressible schemes are affected by a high-frequency spurious noise on the pressure and density fields.[57] This phenomenon is caused by the nonlinear interaction of acoustic waves and by fact that the scheme is explicit in time and centered in space.[58]
Through the years, several techniques have been proposed to get rid of this problem. They can be classified in three different groups:
- the schemes that adopt density filters,
- the models that add a diffusive term in the continuity equation,
- the schemes that employ Riemann solvers to model the particle interaction.
Density filter technique
The schemes of the first group apply a filter directly on the density field to remove the spurious numerical noise. The most used filters are the MLS (moving least squares) and the Shepard filterwhich can be applied at each time step or every n time steps. The more frequent is the use of the filtering procedure, the more regular density and pressure fields are obtained. On the other hand, this leads to an increase of the computational costs. In long time simulations, the use of the filtering procedure may lead to the disruption of the hydrostatic pressure component and to an inconsistency between the global volume of fluid and the density field. Further, it does not ensure the enforcement of the dynamic free-surface boundary condition.
Diffusive term technique
A different way to smooth out the density and pressure field is to add a diffusive term inside the continuity equation (group 2) :
}=\sum _m_\left(_-_\right)\cdot \nabla W_ + \mathcal_i(\rho),}
The first schemes that adopted such an approach were described in Ferrari[59] and in Molteni[56] where the diffusive term was modeled as a Laplacian of the density field. A similar approach was also used in Fatehi and Manzari.[60] In Antuono et al.[61] a correction to the diffusive term of Molteni[56] was proposed to remove some inconsistencies close to the free-surface. In this case the adopted diffusive term is equivalent to a high-order differential operator on the density field.[62] The scheme is called δ-SPH and preserves all the conservation properties of the SPH without diffusion (e.g., linear and angular momenta, total energy,see [63]) along with a smooth and regular representation of the density and pressure fields.
In the third group there are those SPH schemes which employ numerical fluxes obtained through Riemann solvers to model the particle interactions.[64] [65] [66]
Riemann solver technique
For an SPH method based on Riemann solvers, an inter-particle Riemann problem is constructed along a unit vector
pointing form particle
to particle
. In this Riemann problem the initial left and right states are on particles
and
, respectively. The
and
states are
The solution of the Riemann problem results in three waves emanating from the discontinuity. Two waves, which can be shock or rarefaction wave, traveling with the smallest or largest wave speed. The middle wave is always a contact discontinuity and separates two intermediate states, denoted by
and
. By assuming that the intermediate state satisfies
and
, a linearized Riemann solver for smooth flows or with only moderately strong shocks can be written as
where
and
are inter-particle averages. With the solution of the Riemann problem, i.e.
and
, the discretization of the SPH method is
where
v\ast=U\asteij+(\overline{v
}_ - \overline\mathbf_) . This indicates that the inter-particle average velocity and pressure are simply replaced by the solution of the Riemann problem. By comparing both it can be seen that the intermediate velocity and pressure from the inter-particle averages amount to implicit dissipation, i.e. density regularization and numerical viscosity, respectively.
Since the above discretization is very dissipative a straightforward modification is to apply a limiter to decrease the implicit numerical dissipations introduced by limiting the intermediate pressure by[67]
where the limiter is defined as
Note that
ensures that there is no dissipation when the fluid is under the action of an expansion wave, i.e.
, and that the parameter
, is used to modulate dissipation when the fluid is under the action of a compression wave, i.e.
. Numerical experiments found the
is generally effective. Also note that the dissipation introduced by the intermediate velocity is not limited.
Viscosity modelling
In general, the description of hydrodynamic flows require a convenient treatment of diffusive processes to model the viscosity in the Navier–Stokes equations. It needs special consideration because it involves the Laplacian differential operator. Since the direct computation does not provide satisfactory results, several approaches to model the diffusion have been proposed.
Introduced by Monaghan and Gingold[68] the artificial viscosity was used to deal with high Mach number fluid flows. It reads
\Piij=\begin{cases}
\dfrac{-\alpha\bar{c}ij\phiij+\beta
} & \quad \boldsymbol_ \cdot \boldsymbol_ < 0\\ 0 & \quad \boldsymbol_ \cdot \boldsymbol_ \geq 0 \endHere,
is controlling a volume viscosity while
acts similar to the Neumann Richtmeyr artificial viscosity. The
is defined by
\phiij=
⋅ \boldsymbol{r}ij
},where
ηh is a small fraction of
h (e.g. 0.01
h) to prevent possible numerical infinities at close distances.
The artificial viscosity also has shown to improve the overall stability of general flow simulations. Therefore, it is applied to inviscid problems in the following form
\Piij=\alphahc
⋅ \boldsymbol{r}ij
}.It is possible to not only stabilize inviscid simulations but also to model the physical viscosity by this approach. To do so
is substituted in the equation above, where
is the number of spatial dimensions of the model. This approach introduces the bulk viscosity
.
For low Reynolds numbers the viscosity model by Morris[69] was proposed.
[\nu\Delta\boldsymbol{v}]ij=
2\nu
⋅ \nablawh,ij
} \, \boldsymbol_.
Additional physics
- Surface tension
- Heat transfer
- Turbulence
Multiphase extensions
Astrophysics
Often in astrophysics, one wishes to model self-gravity in addition to pure hydrodynamics. The particle-based nature of SPH makes it ideal to combine with a particle-based gravity solver, for instance tree gravity code, particle mesh, or particle-particle particle-mesh.
Solid mechanics and fluid-structure interaction (FSI)
Total Lagrangian formulation for solid mechanics
To discretize the governing equations of solid dynamics, a correction matrix
[70] [71] is first introduced to reproducing rigid-body rotation as
where
stands for the gradient of the kernel function evaluated at the initial reference configuration. Note that subscripts
and
are used to denote solid particles, and smoothing length
is identical to that in the discretization of fluid equations.
Using the initial configuration as the reference, the solid density is directly evaluated as
where
is the Jacobian determinant of deformation tensor
.
We can now discretize the momentum equation in the following form
where inter-particle averaged first Piola-Kirchhoff stress
} is defined as
Also
and
correspond to the fluid pressure and viscous forces acting on the solid particle
, respectively.
Fluid-structure coupling
In fluid-structure coupling, the surrounding solid structure is behaving as a moving boundary for fluid, and the no-slip boundary condition is imposed at the fluid-structure interface. The interaction forces
and
acting on a fluid particle
, due to the presence of the neighboring solid particle
, can be obtained as
[72] and
Here, the imaginary pressure
and velocity
are defined by
where
denotes the surface normal direction of the solid structure, and the imaginary particle density
is calculated through the equation of state.
Accordingly, the interaction forces
and
acting on a solid particle
are given by
and
The anti-symmetric property of the derivative of the kernel function will ensure the momentum conservation for each pair of interacting particles
and
.
Others
The discrete element method, used for simulating granular materials, is related to SPH.
Further reading
- Hoover, W. G. (2006); Smooth Particle Applied Mechanics: The State of the Art, World Scientific.
- Stellingwerf, R. F.; Wingate, C. A.; "Impact Modelling with SPH", Memorie della Societa Astronomia Italiana, Vol. 65, p. 1117 (1994).
- Amada, T.; Imura, M.; Yasumuro, Y.; Manabe, Y.; and Chihara, K. (2004); "Particle-based fluid simulation on GPU", in Proceedings of ACM Workshop on General-purpose Computing on Graphics Processors (August, 2004, Los Angeles, California).
- Desbrun, M.; and Cani, M.-P. (1996). "Smoothed Particles: a new paradigm for animating highly deformable bodies" in Proceedings of Eurographics Workshop on Computer Animation and Simulation (August 1996, Poitiers, France).
- Hegeman, K.; Carr, N. A.; and Miller, G. S. P.; "Particle-based fluid simulation on the GPU", in Proceedings of International Conference on Computational Science (Reading, UK, May 2006), Lecture Notes in Computer Science v. 3994/2006 (Springer-Verlag).
- Kelager, M. (2006) Lagrangian Fluid Dynamics Using Smoothed Particle Hydrodynamics, MSc Thesis, Univ. Copenhagen.
- Kolb, A.; and Cuntz, N. (2005); "Dynamic particle coupling for GPU-based fluid simulation", in Proceedings of the 18th Symposium on Simulation Techniques (2005) pp. 722–727.
- Liu, G. R.; and Liu, M. B.; Smoothed Particle Hydrodynamics: a meshfree particle method, Singapore: World Scientific (2003).
- Monaghan, Joseph J. (1992). "Smoothed Particle Hydrodynamics", Annual Review of Astronomy and Astrophysics (1992). 30 : 543–74.
- Muller, M.; Charypar, D.; and Gross, M.; "Particle-based Fluid Simulation for Interactive Applications", in Breen, D; and Lin, M. (eds.), Proceedings of Eurographics/SIGGRAPH Symposium on Computer Animation (2003).
- Vesterlund, M.; Simulation and Rendering of a Viscous Fluid Using Smoothed Particle Hydrodynamics, MSc Thesis, Umea University, Sweden.
- Violeau, D.; Fluid Mechanics and the SPH method, Oxford University Press (2012).
External links
Software
- Algodoo is a 2D simulation framework for education using SPH
- AQUAgpusph is the free (GPLv3) SPH of the researchers, by the researchers, for the researchers
- dive solutions is a commercial web-based SPH engineering software for CFD purposes
- DualSPHysics is a mostly open source SPH code based on SPHysics and using GPU computing. The open source components are available under the LGPL.
- FLUIDS v.1 is a simple, open source (Zlib), real-time 3D SPH implementation in C++ for liquids for CPU and GPU.
- Fluidix is a GPU-based particle simulation API available from OneZero Software
- GADGET http://www.mpa-garching.mpg.de/galform/gadget/ is a freely available (GPL) code for cosmological N-body/SPH simulations
- GPUSPH SPH simulator with viscosity (GPLv3)
- Pasimodo is a program package for particle-based simulation methods, e.g. SPH
- LAMMPS is a massively parallel, open-source classical molecular dynamics code that can perform SPH simulations
- Physics Abstraction Layer is an open source abstraction system that supports real time physics engines with SPH support
- PreonLab is a commercial engineering software developed by FIFTY2 Technology implementing an implicit SPH method
- Punto is a freely available visualisation tool for particle simulations
- pysph Open Source Framework for Smoothed Particle Hydrodynamics in Python (New BSD License)
- Py-SPHViewer Open Source python visualisation tool for Smoothed Particle Hydrodynamics simulations.
- RealFlow Commercial SPH solver for the cinema industry.
- RheoCube is a commercial SaaS product by Lorenz Research for the study and prediction of complex-fluid rheology and stability
- SimPARTIX is a commercial simulation package for SPH and Discrete element method (DEM) simulations from Fraunhofer IWM
- SPH-flow
- SPHERA
- SPHinXsys is an open source multi-physics, multi-resolution SPH library. It provides C++ APIs for physical accurate simulation and aims to model coupled industrial dynamic systems including fluid, solid, multi-body dynamics and beyond.
- SPHysics is an open source SPH implementation in Fortran
- SPLASH is an open source (GPL) visualisation tool for SPH simulations
- SYMPLER: A freeware SYMbolic ParticLE simulatoR from the University of Freiburg.
- Nauticle is a general-purpose computational tool for particle-based numerical methods.
- NDYNAMICS is a commercial fluid simulation software based on implicit SPH developed by CENTROID LAB currently used for internal/external flooding/nuclear/chemical engineering applications.
Notes and References
- Robert A. . Gingold . Joseph J. . Monaghan . Smoothed particle hydrodynamics: theory and application to non-spherical stars . Monthly Notices of the Royal Astronomical Society . 181 . 3 . 375–89 . 1977 . 10.1093/mnras/181.3.375 . 1977MNRAS.181..375G . free.
- L.B. Lucy . A numerical approach to the testing of the fission hypothesis. Astron. J.. 82 . 1013–1024 . 1977. 10.1086/112164. 1977AJ.....82.1013L.
- Smoothed particle hydrodynamics on GPUs . Takahiro Harada . Seiichi Koshizuka . Yoichiro Kawaguchi . Computer Graphics International . 2007 . 63–70.
- GPUs, a new tool of acceleration in CFD: efficiency and reliability on smoothed particle hydrodynamics methods . Alejandro Crespo . Jose M. Dominguez . Anxo Barreiro . Moncho Gomez-Gesteira . Benedict D. Rogers . PLOS ONE . 6 . 6 . e20685 . 2011 . 2011PLoSO...620685C . 10.1371/journal.pone.0020685 . 21695185 . 3113801 . free .
- Price, D. J.. Smoothed Particle Hydrodynamics: Things I wish my mother taught me. Advances in Computational Astrophysics: Methods. 453 . 249. 1111.1259. 2011. 2012ASPC..453..249P.
- Web site: The Smoothed Particle Hydrodynamics Method vs. Finite Volume Numerical Methods . 2018-08-30. 2018-03-21.
- Adami, S. and Hu, X. Y. and Adams, N. A... A generalized wall boundary condition for smoothed particle hydrodynamics. Journal of Computational Physics. 231. 21. 7057–7075 . 2012 . 10.1016/j.jcp.2012.05.005. 2012JCoPh.231.7057A.
- Shadloo, M. S. and Oger, G. and Touze, D. L... Smoothed particle hydrodynamics method for fluid flows, towards industrial applications: Motivations, current state, and challenges. Computers and Fluids. 136. 11–34. 2016. 10.1016/j.compfluid.2016.05.029.
- Fraser, K.and Kiss, L. I. and St-George, L.. A generalized wall boundary condition for smoothed particle hydrodynamics. 14th International LS-DYNA Conference. 2016.
- Kostorz. A semi-analytical boundary integral method for radial functions with application to Smoothed Particle Hydrodynamics. Journal of Computational Physics. 417. 2020. 109565. 10.1016/j.jcp.2020.109565. 2020JCoPh.41709565K. 219438340.
- Solenthaler. Predictive-Corrective Incompressible SPH. 2009.
- Imhsen. Boundary handling and adaptive time-stepping for PCISPH. 2010. Workshop on Virtual Reality Interaction and Physical Simulation VRIPHYS.
- Bodin. Constraint Fluids . 2011. IEEE Transactions on Visualization and Computer Graphics. 18 . 3 . 516–26 . 10.1109/TVCG.2011.29 . 22241284 . 14023161 .
- Hoetzlein. Fluids v.3, A Large scale, Open Source Fluid Simulator . 2012.
- Akinci. Versatile Rigid-Fluid Coupling for Incompressible SPH . 2012. ACM Transactions on Graphics. 31 . 4 . 1–8 . 10.1145/2185520.2185558 . 5669154 .
- Macklin. Position Based Fluids . 2013. ACM Transactions on Graphics. 32 . 4 . 1–12 . 10.1145/2461912.2461984 . 611962 .
- Akinci. Versatile Surface Tension and Adhesion for SPH Fluids SPH . 2013. ACM Transactions on Graphics. 32 . 6 . 1–8 . 10.1145/2508363.2508395 . 10.1.1.462.8293 . 12550964 .
- Journal of Tribology. Application of Smoothed Particle Hydrodynamics to Full-Film Lubrication. 2013.
- Mahdavi and Talebbeydokhti. Scientia Iranica, Transaction A, Civil Engineering. 22. 4. 1457–1469. A hybrid solid boundary treatment algorithm for smoothed particle hydrodynamics . 2015.
- International Journal for Numerical Methods in Fluids. Curvilinear smoothed particle hydrodynamics. International Journal for Numerical Methods in Fluids. 83. 2. 115–131. 2016. 2017IJNMF..83..115T. 10.1002/fld.4261. 124425780 .
- Colagrossi. Viscous flow past a cylinder close to a free surface: benchmarks with steady, periodic and metastable responses, solved by meshfree and mesh-based schemes. Computers and Fluids. 181. 345–363. 2019. 10.1016/j.compfluid.2019.01.007. 128143912.
- 0903.5075. Price. Daniel J. Astrophysical Smooth Particle Hydrodynamics. New Astron.rev. 53. 4–6. 78–104. 2009. 10.1016/j.newar.2009.08.007. 2009NewAR..53...78R. 129246.
- 10.1007/lrca-2015-1. SPH Methods in the Modelling of Compact Objects. Living Rev Comput Astrophys. 1. 1. 1. 2015. Rosswog. Stephan. 2015LRCA....1....1R. 1406.4224. 119119783.
- astro-ph/0512532. Price. Daniel J. Rockefeller. Gabriel. Warren. Michael S. SNSPH: A Parallel 3-D Smoothed Particle Radiation Hydrodynamics Code. Astrophys. J.. 643. 292–305. 2006. 10.1086/501493. 16733573.
- Web site: Star Formation with Radiative Transfer.
- Price . Daniel J. . 2012-02-01 . Smoothed particle hydrodynamics and magnetohydrodynamics . Journal of Computational Physics . Special Issue: Computational Plasma Physics . 231 . 3 . 759–794 . 10.1016/j.jcp.2010.12.011 . 0021-9991. 1012.1885 .
- Book: Libersky. L.D.. Petschek. A.G.. Smooth Particle Hydrodynamics with Strength of Materials, Advances in the Free Lagrange Method. Lecture Notes in Physics. 1990. 395. 248–257. 10.1007/3-540-54960-9_58. 978-3-540-54960-4.
- L.D. Libersky . A.G. Petschek . A.G. Carney . T.C. Hipp . J.R. Allahdadi . F.A. High. Strain Lagrangian hydrodynamics: a three-dimensional SPH code for dynamic material response. J. Comput. Phys.. 1993. 109. 1 . 67–75. 1993JCoPh.109...67L. 10.1006/jcph.1993.1199.
- J.W. Swegle . D.A. Hicks . S.W. Attaway. Smooth particle hydrodynamics stability analysis. J. Comput. Phys.. 1995. 116. 1 . 123–134. 10.1006/jcph.1995.1010. 1995JCoPh.116..123S .
- T. Belytschko . Y. Krongauz . J. Dolbow . C. Gerlach. On the completeness of meshfree particle methods. Int. J. Numer. Methods Eng.. 1998. 43. 5 . 785–819. 10.1002/(sici)1097-0207(19981115)43:5<785::aid-nme420>3.0.co;2-9. 1998IJNME..43..785B . 10.1.1.28.491 .
- W.K. Liu . S. Jun . Y.F. Zhang. Reproducing kernel particle methods. Int. J. Numer. Methods Eng.. 1995. 20. 8–9. 1081–1106. 10.1002/fld.1650200824. 1995IJNMF..20.1081L .
- P.W. Randles. L.D. Libersky . Recent improvements in SPH modelling of hypervelocity impact. Int. J. Impact Eng.. 1997. 20. 6–10 . 525–532. 10.1016/s0734-743x(97)87441-6.
- G.R. Johnson. S.R. Beissel . Normalized smoothing functions for SPH impact computations. Int. J. Numer. Methods Eng.. 1996. 39. 16 . 2725–2741. 10.1002/(sici)1097-0207(19960830)39:16<2725::aid-nme973>3.0.co;2-9. 1996IJNME..39.2725J .
- C.T. Dyka. R.P. Ingel . An approach for tension instability in Smoothed Particle Hydrodynamics. Comput. Struct.. 1995. 57. 4 . 573–580. 10.1016/0045-7949(95)00059-p.
- C.T. Dyka . P.W. Randles . R.P. Ingel. Stress points for tension instability in SPH. Int. J. Numer. Methods Eng.. 1997. 40. 13 . 2325–2341. 10.1002/(sici)1097-0207(19970715)40:13<2325::aid-nme161>3.0.co;2-8. 1997IJNME..40.2325D .
- P.W. Randles. L.D. Libersky . Normalized SPH with stress points. Int. J. Numer. Methods Eng.. 2000. 48. 10 . 1445–1462. 10.1002/1097-0207(20000810)48:10<1445::aid-nme831>3.0.co;2-9. 2000IJNME..48.1445R .
- T. Belytschko . Y. Guo . W.K. Liu . S.P. Xiao. A unified stability analysis of meshless particle methods. Int. J. Numer. Methods Eng.. 2000. 48. 9 . 1359–1400. 10.1002/1097-0207(20000730)48:9<1359::aid-nme829>3.0.co;2-u. 2000IJNME..48.1359B .
- J. Bonet. S. Kulasegaram . Correction and stabilization of smooth particle hydrodynamics methods with applications in metal forming simulations. Int. J. Numer. Methods Eng.. 2000. 47. 6 . 1189–1214. 10.1002/(sici)1097-0207(20000228)47:6<1189::aid-nme830>3.0.co;2-i. 2000IJNME..47.1189B .
- W. G. Hoover. C. G. Hoover . Computing in Science and Engineering. SPAM-based recipes for continuum simulations. 2001. 3. 2. 78–85. 10.1109/5992.909007 . 2001CSE.....3b..78H . free.
- T. Rabczuk. J. Eibl . L. Stempniewski . Simulation of high velocity concrete fragmentation using SPH/MLSPH. Int. J. Numer. Methods Eng.. 2003. 56. 10 . 1421–1444. 10.1002/nme.617. 2003IJNME..56.1421R . 119799557 .
- M.I. Herreros. M. Mabssout . A two-steps time discretization scheme using the SPH method for shock wave propagation. Comput. Methods Appl. Mech. Engrg.. 2011. 200. 21–22 . 1833–1845. 10.1016/j.cma.2011.02.006. 2011CMAME.200.1833H .
- S. Yashiro. T. Okabe . Smoothed particle hydrodynamics in a generalized coordinate system with a finite-deformation constitutive model. Int. J. Numer. Methods Eng.. 2015. 103. 11. 781–797. 10.1002/nme.4906. 2015IJNME.103..781Y. 2324/4476056. 121935150 . free.
- Truncation error in mesh-free particle methods . N.J. Quinlan . M. Basa . M. Lastiwka . International Journal for Numerical Methods in Engineering . 66 . 13 . 2064–2085 . 2006 . 10.1002/nme.1617 . 10379/1170 . 2006IJNME..66.2064Q . 53701216 . free .
- Piecewise polynomial, positive definite and compactly supported radial functions of minimal degree . H. Wendland . Advances in Computational Mathematics . 4 . 4 . 389–396 . 1995 . 10.1007/BF02123482 . 36452865 .
- Investigation of wall bounded flows using SPH and the unified semi-analytical wall boundary conditions . A. Mayrhofer . B.D. Rogers . D. Violeau . M. Ferrand . Computer Physics Communications . 184 . 11 . 2515–2527 . 2013 . 10.1016/j.cpc.2013.07.004 . 10.1.1.770.4985 . 2013CoPhC.184.2515M . 1304.3692 . 35008128 .
- Variational and momentum preservation aspects of Smoothed Particle Hydrodynamics formulations . J. Bonet . T.S. Lok . Computers Methods in Applied Mechanical Engineering . 180 . 1–2 . 97–115 . 1999 . 10.1016/S0045-7825(99)00051-1 . 1999CMAME.180...97B .
- Smoothed particle hydrodynamics . J.J. Monaghan . Reports on Progress in Physics . 68 . 8 . 1703–1759 . 2005 . 10.1088/0034-4885/68/8/R01 . 2005RPPh...68.1703M . 5987481 .
- Book: Geometric Numerical Integration . E. Hairer . C. Lubich . G. Wanner . 2006 . Springer . 978-3-540-30666-5 .
- Andrea Colagrossi . Matteo Antuono . David Le Touzè . Theoretical considerations on the free-surface role in the smoothed-particle-hydrodynamics model . Physical Review E . 79 . 5 . 056701 . 2009. 2009PhRvE..79e6701C . 10.1103/PhysRevE.79.056701 . 19518587 .
- Bejamin Bouscasse . Andrea Colagrossi . Salvatore Marrone . Matteo Antuono . Nonlinear water wave interaction with floating bodies in SPH . Journal of Fluids and Structures . 42 . 112–129 . 2013. 10.1016/j.jfluidstructs.2013.05.010 . 2013JFS....42..112B .
- Fabricio Macià . Matteo Antuono . Leo M González . Andrea Colagrossi . Theoretical analysis of the no-slip boundary condition enforcement in SPH methods . Progress of Theoretical Physics . 125 . 6 . 1091–1121 . 2011. 10.1143/PTP.125.1091 . 2011PThPh.125.1091M . free .
- Jose Luis Cercos-Pita . Matteo Antuono . Andrea Colagrossi . Antonio Souto . SPH energy conservation for fluid--solid interactions . Computer Methods in Applied Mechanics and Engineering . 317 . 771–791 . 2017. 10.1016/j.cma.2016.12.037 . 2017CMAME.317..771C .
- J. Campbell . R. Vignjevic . L. Libersky . A contact algorithm for smoothed particle hydrodynamics . Computer Methods in Applied Mechanics and Engineering . 184 . 1 . 49–65 . 2000. 10.1016/S0045-7825(99)00442-9 . 2000CMAME.184...49C .
- M. Ferrand, D.R. Laurence, B.D. Rogers, D. Violeau, C. Kassiotis . Unified semi-analytical wall boundary conditions for inviscid, laminar or turbulent flows in the meshless SPH method. International Journal for Numerical Methods in Fluids. Int. J. Numer. Meth. Fluids. 71. 4. 446–472. 2013. 2013IJNMF..71..446F. 10.1002/fld.3666. 124465492.
- Book: H. R. Cole . Underwater Explosions . Princeton University Press . Princeton, New Jersey . 1948.
- D. Molteni, A. Colagrossi. A simple procedure to improve the pressure evaluation in hydrodynamic context using the SPH. Computer Physics Communications. 180. 6. 861–872. 2009. 10.1016/j.cpc.2008.12.004. 2009CoPhC.180..861M.
- Colagrossi . Andrea . Landrini . Maurizio . Numerical simulation of interfacial flows by smoothed particle hydrodynamics . Journal of Computational Physics . 2003 . 191 . 2 . 448–475 . 10.1016/S0021-9991(03)00324-3 . 2003JCoPh.191..448C .
- Book: Randall J. LeVeque . Finite difference methods for ordinary and partial differential equations: steady-state and time-dependent problems. Siam. 2007.
- A new 3D parallel SPH scheme for free surface flows. A. Ferrari, M. Dumbser, E. Toro, A. Armanini. Computers & Fluids. 38. 6. 1203–1217. 2009. Elsevier. 10.1016/j.compfluid.2008.11.012.
- A remedy for numerical oscillations in weakly compressible smoothed particle hydrodynamics. Fatehi, R and Manzari, MT . International Journal for Numerical Methods in Fluids. 67. 9. 1100–1114. 2011. Wiley Online Library . 10.1002/fld.2406 . 2011IJNMF..67.1100F . 121381641 .
- Free-surface flows solved by means of SPH schemes with numerical diffusive terms. M. Antuono, A. Colagrossi, S. Marrone, D. Molteni. Computer Physics Communications. 181. 3. 532–549. 2010. Elsevier . 10.1016/j.cpc.2009.11.002. 2010CoPhC.181..532A.
- Numerical diffusive terms in weakly-compressible SPH schemes. M. Antuono, A. Colagrossi, S. Marrone. Computer Physics Communications. 183. 12. 2570–2580. 2012. Elsevier . 10.1016/j.cpc.2012.07.006. 2012CoPhC.183.2570A.
- Energy balance in the δ-SPH scheme. Antuono Matteo and Marrone S and Colagrossi A and Bouscasse B. Computer Methods in Applied Mechanics and Engineering. 289. 209–226. 2015. Elsevier . 10.1016/j.cma.2015.02.004 . 2015CMAME.289..209A .
- On particle weighted methods and smooth particle hydrodynamics. JP. Vila. Mathematical Models and Methods in Applied Sciences. 9. 2. 161–209. 1999. World Scientific . 10.1142/S0218202599000117.
- Free surface flows simulations in Pelton turbines using an hybrid SPH-ALE method. Marongiu Jean-Christophe and Leboeuf Francis and Caro Joëlle and Parkinson Etienne. Journal of Hydraulic Research. 48. S1. 40–49. 2010. Taylor & Francis. 10.1080/00221686.2010.9641244 . 121493014 .
- Book: Modelisation d'écoulements visqueux par methode SPH en vue d'application à l'hydrodynamique navale. De Leffe, Matthieu. 2011. PhD Thesis, Ecole centrale de Nantes.
- A weakly compressible SPH method based on a low-dissipation Riemann solver. Chi Zhang and Xiangyu Hu and Nikolaus Adams. Journal of Computational Physics. 335. 605–620. 2017. 10.1016/j.jcp.2017.01.027. 2017JCoPh.335..605Z .
- Joseph J. . Monaghan . Robert A. . Gingold . Shock Simulation by the Particle Method . Journal of Computational Physics . 1983. 52. 2 . 347–389. 1983JCoPh..52..374M . 10.1016/0021-9991(83)90036-0 .
- J. P. Morris. P. J. Fox. Y. Zhu. Modeling Low Reynolds Number Incompressible Flows Using SPH. Journal of Computational Physics. 1997. 136. 1. 214–226. 1997JCoPh.136..214M. 10.1006/jcph.1997.5776.
- Vignjevic, Rade . Reveles, Juan R. . Campbell, James . SPH in a total Lagrangian formalism . Computer Modeling in Engineering and Sciences . 44 . 181–198 . 2006 .
- Han, Luhui . Hu, Xiangyu . SPH modeling of fluid-structure interaction . Journal of Hydrodynamics . 30 . 62–69 . 2018 . 1 . 10.1007/s42241-018-0006-9 . 2018JHyDy..30...62H . 125369012 .
- A multi-resolution SPH method for fluid-structure interactions. Chi Zhang . Massoud Rezavand . Xiangyu Hu. Journal of Computational Physics. 110028 . 2020. 429. 10.1016/j.jcp.2020.110028. 0021-9991. 1911.13255. 208513116.