Smoothed-particle hydrodynamics

Smoothed-particle hydrodynamics (SPH) is a computational method used for simulating fluid flows. It was developed by Gingold and Monaghan (1977) and Lucy (1977) initially for astrophysical problems. It has been used in many fields of research, including astrophysics, ballistics, volcanology, and oceanography. It is a mesh-free Lagrangian method (where the coordinates move with the fluid), and the resolution of the method can easily be adjusted with respect to variables such as the density.

Method

The smoothed-particle hydrodynamics (SPH) method works by dividing the fluid into a set of discrete elements, referred to as particles. These particles have a spatial distance (known as the "smoothing length", typically represented in equations by $h$), over which their properties are "smoothed" by a kernel function. This means that the physical quantity of any particle can be obtained by summing the relevant properties of all the particles which lie within the range of the kernel. For example, using Monaghan's popular cubic spline kernel the temperature at position $\mathbf{r}$ depends on the temperatures of all the particles within a radial distance $2h$ of $\mathbf{r}$.

The contributions of each particle to a property are weighted according to their distance from the particle of interest, and their density. Mathematically, this is governed by the kernel function (symbol $W$). Kernel functions commonly used include the Gaussian function and the cubic spline. The latter function is exactly zero for particles further away than two smoothing lengths (unlike the Gaussian, where there is a small contribution at any finite distance away). This has the advantage of saving computational effort by not including the relatively minor contributions from distant particles.

The equation for any quantity $A$ at any point $\mathbf{r}$ is given by the equation

$A(\mathbf{r}) = \sum_j m_j \frac{A_j}{\rho_j} W(| \mathbf{r}-\mathbf{r}_{j} |,h),$

where $m_j$ is the mass of particle $j$, $A_j$ is the value of the quantity $A$ for particle $j$, $\rho_j$ is the density associated with particle $j$, $\mathbf{r}$ denotes position and $W$ is the kernel function mentioned above. For example, the density of particle $i$ ($\rho_i$) can be expressed as:

$\rho_i = \rho(\mathbf{r}_i) = \sum_j m_j \frac{\rho_j}{\rho_j} W(| \mathbf{r}_i-\mathbf{r}_j |,h) = \sum_j m_j W(| \mathbf{r}_i-\mathbf{r}_j |,h),$

where the summation over $j$ includes all particles in the simulation.

Similarly, the spatial derivative of a quantity can be obtained easily by virtue of the linearity of the derivative (del, $\nabla$).

$\nabla A(\mathbf{r}) = \sum_j m_j \frac{A_j}{\rho_j} \nabla W(| \mathbf{r}-\mathbf{r}_j |,h).$

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. Combined with an equation of state and an integrator, SPH can simulate hydrodynamic flows efficiently. However, the traditional artificial viscosity formulation used in SPH tends to smear out shocks and contact discontinuities to a much greater extent than state-of-the-art grid-based schemes.

The Lagrangian-based adaptivity of SPH is analogous to the adaptivity present in grid-based adaptive mesh refinement codes. In some ways it is actually simpler because SPH particles lack any explicit topology relating them, unlike the elements in FEM. Adaptivity in SPH can be introduced in two ways; either by changing the particle smoothing lengths or by splitting SPH particles into 'daughter' particles with smaller smoothing lengths. The first method is common in astrophysical simulations where the particles naturally evolve into states with large density differences.^[1] However, in hydrodynamics simulations where the density is often (approximately) constant this is not a suitable method for adaptivity. For this reason particle splitting can be employed, with various conditions for splitting ranging from distance to a free surface ^[2] through to material shear.^[3]

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,^[4] particle mesh, or particle-particle particle-mesh.

Uses in astrophysics

The adaptive resolution of smoothed-particle hydrodynamics, combined with its ability to simulate phenomena covering many orders of magnitude, make it ideal for computations in theoretical astrophysics.

Simulations of galaxy formation, star formation, stellar collisions, supernovae 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.^[5]

Uses in fluid simulation

Fig. SPH simulation of ocean waves using FLUIDS v.1 (Hoetzlein)

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-base technique 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 which 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^[6]
• M. Ihmsen et al., 2010, introduce boundary handling and adaptive time-stepping for PCISPH for accurate rigid body interactions^[7]
• K. Bodin et al., 2011, replace the standard equation of state pressure equation with a density constraint and applies a variational time integrator.^[8]
• R. Hoetzlein, 2012, develops efficient GPU-based SPH for large scenes in Fluids v.3^[9]
• 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 ^[10]
• M. Macklin et al., 2013 simulates incompressible flows inside the Position Based Dynamics framework, for bigger timesteps ^[11]
• 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.^[12]
• J. Kyle and E. Terrell, 2013, SPH Applied to Full-Film Lubrication^[13]

Uses in solid mechanics

In 1990, Libersky and Petschek ^[14][15] extended SPH to Solid Mechanics.

The main advantage of SPH 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 ^[16] and lack of consistency.^[17] Over the past years, different corrections have been introduced to improve the accuracy of the SPH solution. That is the case of Liu et al.,^[18] Randles and Libersky ^[19] and Johnson and Beissel,^[20] who tried to solve the consistency problem. Dyka et al.^[21][22] and Randles and Libersky ^[23] introduced the stress-point integration into SPH and Belytschko et al.^[24] showed later 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. Here are some recent examples of applications and developments of the method:

• Libersky and Petschek ^[14] modified the SPH method in order to solve Strength of Materials problems.
• Johnson and Beissel ^[20] and Randles and Libersky ^[19] applied SPH to impact phenomena.
• Bonet and Kulasegaram applied SPH to metal forming simulations.^[25]
• William G. Hoover developed a SPH-based method whose acronym is SPAM (smooth-particle applied mechanics) to study impact fracture in solids.^[26]
• Rabczuk and co-workers applied a modified SPH (SPH/MLSPH) to simulate fracture and fragmentation.^[27]
• Herreros and Mabssout developed the Taylor-SPH (TSPH) method to overcome the problem of shock wave propagation in solids.^[28]

Notes

1. http://arxiv.org/abs/astro-ph/9512078
2. http://dl.acm.org/citation.cfm?id=1568695
3. http://www.ase.uc.edu/~liugr/Storage/Journal%20Papers/2006/JA_2006_09.pdf
4. "The Parallel k-D Tree Gravity Code"; "PKDGRAV (Parallel K-D tree GRAVity code" use a kd-tree gravity simulation.
5. http://www.astro.ex.ac.uk/people/mbate/Cluster/clusterRT.html
6. Solenthaler (2009). "Predictive-Corrective Incompressible SPH".
7. Imhsen (2010). "Boundary handling and adaptive time-stepping for PCISPH". Workshop on Virtual Reality Interaction and Physical Simulation VRIPHYS.
8. Bodin (2011). "Constraint Fluids. http://www.physics.umu.se/english/research/statistical-physics-and-networks/complex-mechanical-systems/fluids-and-solids/". IEEE Transactions on Visualization and Computer Graphics. External link in |title= (help)
9. Hoetzlein (2012). "Fluids v.3, A Large scale, Open Source Fluid Simulator. http://fluids3.com". External link in |title= (help)
10. Akinci (2012). "Versatile Rigid-Fluid Coupling for Incompressible SPH http://www.nadir.tk/research". ACM TOG, SIGGRAPH proceedings. External link in |title= (help)
11. Macklin (2013). "Position Based Fluids http://blog.mmacklin.com/publications". ACM TOG, SIGGRAPH proceedings. External link in |title= (help)
12. Akinci (2013). "Versatile Surface Tension and Adhesion for SPH Fluids SPH http://www.nadir.tk/research". ACM TOG, SIGGRAPH proceedings. External link in |title= (help)
13. Journal of Tribology (2013). "Application of Smoothed Particle Hydrodynamics to Full-Film Lubrication".
14. Libersky, L.D.; Petschek, A.G. (1990). "Smooth Particle Hydrodynamics with Strength of Materials, Advances in the Free Lagrange Method". Lecture Notes in Physics 395: 248–257. doi:10.1007/3-540-54960-9_58.
15. L.D. Libersky; A.G. Petschek; A.G. Carney; T.C. Hipp; J.R. Allahdadi; F.A. High (1993). "Strain Lagrangian hydrodynamics: a three-dimensional SPH code for dynamic material response". J. Comput. Phys. 109: 67–75.
16. J.W. Swegle; D.A. Hicks; S.W. Attaway (1995). "Smooth particle hydrodynamics stability analysis". J. Comput. Phys. 116: 123–134. Bibcode:1995JCoPh.116..123S. doi:10.1006/jcph.1995.1010.
17. T. Belytschko; Y. Krongauz; J. Dolbow; C. Gerlach (1998). "On the completeness of meshfree particle methods". Int. J. Numer. Methods Eng. 43: 785–819. Bibcode:1998IJNME..43..785B. doi:10.1002/(sici)1097-0207(19981115)43:5<785::aid-nme420>3.0.co;2-9.
18. W.K. Liu; S. Jun; Y.F. Zhang (1995). "Reproducing kernel particle methods". Int. J. Numer. Methods Eng. 20 (8–9): 1081–1106. Bibcode:1995IJNMF..20.1081L. doi:10.1002/fld.1650200824.
19. P.W. Randles; L.D. Libersky (1997). "Recent improvements in SPH modelling of hypervelocity impact". Int. J. Impact Eng. 20: 525–532. doi:10.1016/s0734-743x(97)87441-6.
20. G.R. Johnson; S.R. Beissel (1996). "Normalized smoothing functions for SPH impact computations". Int. J. Numer. Methods Eng. 39: 2725–2741. Bibcode:1996IJNME..39.2725J. doi:10.1002/(sici)1097-0207(19960830)39:16<2725::aid-nme973>3.0.co;2-9.
21. C.T. Dyka; R.P. Ingel (1995). "An approach for tension instability in Smoothed Particle Hydrodynamics". Comput. Struct. 57: 573–580. doi:10.1016/0045-7949(95)00059-p.
22. C.T. Dyka; P.W. Randles; R.P. Ingel (1997). "Stress points for tension instability in SPH". Int. J. Numer. Methods Eng. 40: 2325–2341. Bibcode:1997IJNME..40.2325D. doi:10.1002/(sici)1097-0207(19970715)40:13<2325::aid-nme161>3.0.co;2-8.
23. P.W. Randles; L.D. Libersky (2000). "Normalized SPH with stress points". Int. J. Numer. Methods Eng. 48: 1445–1462. Bibcode:2000IJNME..48.1445R. doi:10.1002/1097-0207(20000810)48:10<1445::aid-nme831>3.0.co;2-9.
24. T. Belytschko; Y. Guo; W.K. Liu; S.P. Xiao (2000). "A unified stability analysis of meshless particle methods". Int. J. Numer. Methods Eng. 48: 1359–1400. Bibcode:2000IJNME..48.1359B. doi:10.1002/1097-0207(20000730)48:9<1359::aid-nme829>3.0.co;2-u.
25. J. Bonet; S. Kulasegaram (2000). "Correction and stabilization of smooth particle hydrodynamics methods with applications in metal forming simulations". Int. J. Numer. Methods Eng. 47: 1189–1214. Bibcode:2000IJNME..47.1189B. doi:10.1002/(sici)1097-0207(20000228)47:6<1189::aid-nme830>3.0.co;2-i.
26. W. G. Hoover; C. G. Hoover (2001). Computing in Science and Engineering 3 (2): 78–85. Missing or empty |title= (help)
27. T. Rabczuk; J. Eibl; L. Stempniewski (2003). "Simulation of high velocity concrete fragmentation using SPH/MLSPH". Int. J. Numer. Methods Eng. 56: 1421–1444. Bibcode:2003IJNME..56.1421R. doi:10.1002/nme.617.
28. M.I. Herreros; M. Mabssout (2011). "A two-steps time discretization scheme using the SPH method for shock wave propagation". Comput. Methods Appl. Mech. Engrg. 200: 1833–1845. Bibcode:2011CMAME.200.1833H. doi:10.1016/j.cma.2011.02.006.

References

• [1] R.A. Gingold and J.J. Monaghan, "Smoothed particle hydrodynamics: theory and application to non-spherical stars," Mon. Not. R. Astron. Soc., Vol 181, pp. 375–89, 1977.
• [2] L.B. Lucy, "A numerical approach to the testing of the fission hypothesis," Astron. J., Vol 82, pp. 1013–1024, 1977.
• [3] Hoover, W. G. (2006). Smooth Particle Applied Mechanics: The State of the Art, World Scientific.
• [4] Impact Modelling with SPH Stellingwerf, R. F., Wingate, C. A., Memorie della Societa Astronomia Italiana, Vol. 65, p. 1117 (1994).
• [5] 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).
• [6] 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).
• [7] Harada, T., Koshizuka, S. and Kawaguchi, Y. Smoothed Particle Hydrodynamics on GPUs. In Proceedings of Computer Graphics International (June 2007, Petropolis Brazil).
• [8] 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). Proceedings published as Lecture Notes in Computer Science v. 3994/2006 (Springer-Verlag).
• [9] M. Kelager. (2006) Lagrangian Fluid Dynamics Using Smoothed Particle Hydrodynamics, M. Kelagar (MS Thesis, Univ. Copenhagen).
• [10] Kolb, A. and Cuntz, N. (2005) ] Dynamic particle coupling for GPU-based fluid simulation, A. Kolb and N. Cuntz. In Proceedings of the 18th Symposium on Simulation Techniques (2005) pp. 722–727.
• [11] Liu, G.R. and Liu, M.B. Smoothed Particle Hydrodynamics: a meshfree particle method. Singapore: World Scientific (2003).
• [12] Monaghan, J.J. (1992). Smoothed Particle Hydrodynamics. Annu. Rev. Astron. Astrophys. (1992). 30 : 543-74.
• [13] Muller, M., Charypar, D. and Gross, M. ] Particle-based Fluid Simulation for Interactive Applications, In Proceedings of Eurographics/SIGGRAPH Symposium on Computer Animation (2003), eds. D. Breen and M. Lin.
• [14] Vesterlund, M. Simulation and Rendering of a Viscous Fluid Using Smoothed Particle Hydrodynamics, (MS Thesis, Umea University, Sweden).
• [15] Violeau, D., Fluid Mechanics and the SPH method. Oxford University Press (2012).

External links

• First large simulation of star formation using SPH
• SPHERIC (SPH European Research Interest Community)
• ITVO is the web-site of The Italian Theoretical Virtual Observatory created to query a database of numerical simulation archive.
• SPHC Image Gallery depicts a wide variety of test cases, experimental validations, and commercial applications of the SPH code SPHC.
• A derivation of the SPH model starting from Navier-Stokes equations

Software

• Algodoo is a 2D simulation framework for education using SPH
• pysph Open Source Framework for Smoothed Particle Hydrodynamics in Python (New BSD License)
• DualSPHysics is an open source SPH code based on SPHysics and using GPU computing
• Fluidix is a GPU-based particle simulation API available from OneZero Software
• FLUIDS v.1 is a simple, open source (Zlib), real-time 3D SPH implementation in C++ for liquids for CPU and GPU.
• GADGET is a freely available (GPL) code for cosmological N-body/SPH simulations
• GPUSPH SPH simulator with viscosity (GPLv3)
• SimPARTIX is a commercial simulation package for SPH and DEM simulations from Fraunhofer IWM
• SPLASH is an open source (GPL) visualisation tool for SPH simulations
• SPH-flow
• SPHysics is an open source SPH implementation in Fortran
• SYMPLER: A freeware SYMbolic ParticLE simulatoR from the University of Freiburg.
• Physics Abstraction Layer is an open source abstraction system that supports real time physics engines with SPH support
• Pasimodo is a program package for particle-based simulation methods, e.g. SPH
• Punto is a freely available visualisation tool for particle simulations
• AQUAgpusph is the free (GPLv3) SPH of the researchers, by the researchers, for the researchers