FLIP fluid simulation with sample code

I've been continuing to play around with incompressible fluid simulation and have implemented the FLIP method (FLuid Implicit Particle). The FLIP algorithm represents the fluid concentration and velocity field using a large set of particles.  At every iteration these particles are transferred to an auxiliary grid using kernel density estimation which defines a grid-based velocity field which is made incompressible by pressure projection and then used to advect the particles. Finally the velocity of each particle is corrected using the difference between the projected velocity field and the initial velocity field obtained by kernel density estimation.

The primary benefit of this approach is that it has extremely low numerical viscosity, making simulating fluids such as smoke and water possible at sane grid resolutions.  The downside is considerably more complexity than a purely grid-based solver; you need to track particles, perform the mapping between particles and grid (and vice versa) and have some sensible scheme to reseed particles when either there are too few or too many in a region.

To play with the method I implemented it in C++, largely following Robert Bridson's Fluid Simulation for Computer Graphics.  It is a canned example demonstrating variable density flow under the influence of gravity using an body force proportional to the fluid 'concentration', a sort of poor-man's Boussinesq approximation.  The simulation below was performed on a 100x100 auxiliary grid with 16 particles per grid-cell. This is more than the 2x2 typically used, but helps get around reseeding.  As this is not free-surface flow, there are particles everywhere which allows causes the smoke to form nice instabilities and finally turbulent mixing (although I did not run it particularly far).

In order to have an unmangled version of the code available to myself, I am posting it here.  It is uncommented and unstructured, but (hopefully) fairly understandable once you understand the basic algorithm.  The code is bundled with the GMM++ sparse matrix library for the pressure solves and required VTK and CMake for output and as a build-system respectively.

You can download the code here: FLIP2.zip