Index transformation between bounding boxes and uniform grids

I end up rewriting this code all the time. It transforms from points in space to grid coordinates world_to_grid() and back grid_to_world(), such that p = grid_to_world( aabb, dim, world_to_grid( aabb, dim, p ) ).  It's simple, but bugs here can mess up lots of things in ways that are hard to detect.  No range checking is performed in order to make it easy to apply custom boundary conditions.

    template< typename real, typename index, typename real3 >
    inline real3 _world_to_grid( const real *aabb, const index *dim, const real3 &pos ){
        return real3( real(dim[0]-1)*(pos[0]-aabb[0])/(aabb[1]-aabb[0]), real(dim[1]-1)*(pos[1]-aabb[2])/(aabb[3]-aabb[2]), real(dim[2]-1)*(pos[2]-aabb[4])/(aabb[5]-aabb[4]) );
    }
    
    template< typename real, typename index, typename real3 >
    inline real3 _grid_to_world( const real *aabb, const index *dim, const real3 &pos ){
        return real3( aabb[0]+real(pos[0])*(aabb[1]-aabb[0])/real(dim[0]-1), aabb[2]+real(pos[1])*(aabb[3]-aabb[2])/real(dim[1]-1), aabb[4]+real(pos[2])*(aabb[5]-aabb[4])/real(dim[2]-1) );
    }  

The code assumes that the dim variable holds the number of points in each direction (e.g. int dim[] = { nx, ny, nz };) and that the first point is coincident with the minimum value of the axis-aligned bounding box with the last point coincident with the maximum value of the bounding box. The aabb parameter is the axis-aligned bounding box defining the grid in world space, e.g. int aabb[] = { xmin, xmax, ymin, ymax, zmin, zmax };

A follow-up to fluid simulation on non-uniform grids

In my last post, I discussed preliminary results for fluid simulation on non-uniform Cartesian grids.  In that post I showed some preliminary results, but there were some bugs that added disturbing artifacts.

I have fixed a number of bugs and now have a solver based on BFECC advection using min-max limited cubic interpolation for non-uniform and often highly anisotropic meshes for the velocity and pressure solve, with high-resolution uniform density fields for the density field.  The results are fairly impressive:

This image shows a volume rendering (in Paraview) of a simulation computed using the grid in the background.  Near the area of interest, the grid is uniform, but it grows very quickly (geometrically, with growth rate ~1.5) outside this region.  The velocity/pressure grid is 133x127x134, but covers nearly a 10x6x10 m cubic volume, with 2cm cells in the fine region.  The density field is 1x6x1 m with 1cm uniform resolution.

Being able to run different resolutions and gradings for the velocity and density fields is extremely helpful: fine fluid details distract from a poor velocity solution, and high-resolution densities help avoid diffusion in the physics.  The image above shows the density as resolved by the fluid grid. It is terrible.  However the density as resolved by the density grid is -way- better:

It's still not perfect, but given the cell size and anisotropy, I think it does extremely well.  Although there are definite artifacts, the payoff is in the memory usage and runtime. The whole simulation is 5 seconds in real time and takes approximately 22 seconds per output frame, meaning my 5 second simulation completes in under an hour.  These simulations used to take on the order of 4-5 hours.

The results look pretty good.  I think the grading is too steep to get really nice results, but it's an excellent proof-of-concept.

Fluid Simulation for Graphics on Non-Uniform Structured Grids

I've been playing around with fluid simulation on non-uniform structured grids.  They have some charms in that it is easy to have very large domains with isolated regions of interest; e.g. near the camera or salient fluid features.  One of the big advantages is that it makes far-field boundaries easy, you simply extend the mesh far away from the domain.

My solver pretty run-of-the-mill; only the pressure solver was updated in order to handle anisotropic cells.  I've used a finite-volume formulation to derive the pressure-correction equation for this case, but it suffices to say that only the weights used in forming the Poisson equation change.

Here is an example, a grid that is roughly 10x10x5 meters, with 2cm fine cells for roughly 100x100x100 cells:

The aspect ratios get quite high, 100:1 is not uncommon.   You can see a buoyant flow simulation that I'm running as the contour values.  The simulation takes about 20 seconds per frame, of which 7 seconds is writing the VTK output files that I use to post-process.  In spite of that, the visual detail in the region of interest is excellent (given the resolution):

I find that it helps to advect a uniform resolution density field that is roughly 2X the fine cell resolution.  Writing the code with general-purpose fields that can be point-sampled arbitrarily makes this a trivial feature to implement, simply allocate different grids for density and velocity.

Finally here is a rendering of the final simulation in Paraview, total simulation time ~35 minutes. There are still some bugs to track down, but the results are pretty promising: