SPH Fluid Simulation

Recalled from the previous post, there are 2 major approaches used to describe the fluid motion, Eulerian and Lagrangian. Eulerian is grid-based, measuring properties at fixed points in space. Physical values are stored on grids. Lagrangian is particle-based, measuring properties as particles drift through the flow. Physical values are stored on each particle. Both methods have pros and cons. The prevalent methods used nowadays are mostly the hybrids. (e.g. FLIP is the hybrid between particle-based and volume-based).

Let's start from the particle-based approach and take SPH (Smoothed Particle Hydrodynamics) for example. Via discrete particles located arbitrary in space , SPH utilizes weighted summation from adjacent particles to interpolate values and derivative of continuous field.

SPH uses smoothing kernel (W) to determine which adjacent particles have how much contribution of values on one particular point in space. The selected smoothing kernel should be normalized, positive, even and last but not least, no kernel interaction outside the range of the radius "h". The derivatives of the field could be approximated with analytical differentiation:

Since the physical value at point j doesn't directly depend on other components in space, thus, the gradient of the field only has effect on the smoothing kernel. Same for Laplacian.

With above SPH equations, we can start to approximate each term in Navier-Stokes's Equation. Since the main idea of Lagrangian method is the particle itself moves with the fluid, thus the advection term can be ignored, which simplify the equation a lot. For other terms (e.g. density, pressure, viscosity, external forces, etc), they all can be approximated with SPH approach. I won't paste the equation for each term here since they all can be found from lots of material online, and the main purpose of this article is just for getting myself familiar with the big picture and techniques behind the scene. However, there is a thing worth mentioning,  if just using SPH equation above to directly approximate something like pressure and viscosity force, we will confront an asymmertric problem. Think about a situation, if there are only 2 particles with contact in space, the pressure of one particle only uses the pressure force of another, and since the pressure for each partilce may be different, the result pressure will not be symmetrical. These can be fixed with few ways to rewrite the equation.

OK. With the above knowledge, let's talk about the advantages and disadvantages of SPH. There are several benefits of using particle-based approach over pure grid-based method. First of all, the equation is more intuitive and easier to implement. Second, since SPH stores mass on particles, it guarantees the conservation of mass w/o extra calculation. Third, SPH doesn't have to track fluid boundaries; it supports free surface effect. But there still has to be a way for mesh generation for rendering since SPH doesn't contain topology info.

The obvious drawback of using particle methods is the need of large numbers of particles to get more accurate result, which enhances computational cost. The other thing is the density summation would cause problems when neighboring particles have different rest densities (i.e. different masses). This would lead to incorrect pressure and have weird gaps (break into pieces?). Some solutions have been introduced which involves alternative density computation, formulas rewriting and new surface tension model.

What if we need the boundary to limit the movement of particles? One way could be adding boundary particles to prevent penetration and inverse the velocity. Another way is to use Ghost Particles. There is a technical paper about Ghost SPH in the upcoming Siggraph 2012: http://www.cs.ubc.ca/~rbridson/docs/schechter-siggraph2012-ghostsph.pdf