CS 348b Final Project - Where There's Smoke...
My project focused on rendering images of physically based images of smoke. Originally I had proposed using similar techniques based on simulating the fluid flow of gases to render fire, but getting the fluid simulator proved more challenging that I estimated, so in the end I decided to render the simpler smoke model. The basic method used is to model the flow of the gas in which the smoke is suspended using the Navier-Stokes equations, generating a density field for the smoke, and then rendering the density field using standard volume rendering techniques.
The primary task of the implementation was developing the fluid simulator to simulate the flow of the gas transporting the smoke. In this I followed the method described in Fedkiw et al.. Specifically, the gas is modeled as an incompressible, inviscid fluid which transports a density field used to represent the smoke. To do this I implemented a solver using the semi-Lagrangian method, which was first described in the computer graphics context by Stam .
Conceptually, the structure of the solver is fairly straightforward and can be broken down into four main components, based on accounting for different terms in the equations of flow. The space to be simulated is discretized into a regular rectangular grid, with each cell recording the velocity and temperature of the gas and density of the smoke. At each time step, the flow of the fluid is then simulated by 1) Adding external forces to the flow (e.g. buoyant forces from temperature differences and gravity forces acting on the smoke density) 2) Advecting (transporting) the velocity fields 3) enforcing conservation of mass by solving a discrete Poisson equation for the pressure 4) Advecting the temperture field and smoke density.
Each of these components is implemented as follows:
External forces: The velocity field is simply updated by adding the acceleration over the given time step: v += F dt
Advection: This is the key step of the semi-Lagrangian method. One can show that for equations of the form dS/dt = v . grad S the value of S at t+1 and a point x can be computed by simply taking the value of S at the point which will flow to x at time t. Thus to implement this, we can just trace the particle backwards in time from xusing the velocity field to find the next value.  Describes this in some detail.
Conservation of mass: To ensure conservation of mass, we need to enforce the constraint div v = 0. This can be shown to be equivalent to solving for the pressure, using Poisson's equation. This is a standard problem in numerical analysis, and there are several techniques available. Here I used the successive overrelaxation method.
- Advection of scalar fields: This part is essentially the same as for updating the velocity.
In addition to the standard semi-Lagrangian method described above, Fedkiw et al. introduced a "vorticity confinement" term, which I also implemented. The basic idea is that the semi-Lagrangian method tends to excessively damp the flow of the fluid, resulting in overly smooth looking flows. This small scale turbulence can be recovered, however, by adding appropriate external forces to the flow, specifically, forces related to the curl (vorticity) of the velocity field, since this is an indicator of the amount of turbulence locally at the flow.
Once I had the fluid solver, I used it to generate density fields of smoke. To render the smoke, I used PBRT's single scattering volume integrator. I tried several different values for the scattering, absorption and Henyey-Greenstein parameters. As expected, high sigma_a values generated a nice black sooty looking smoke, while high sigma_s values generated a whiter looking smoke. In general the single scattering approximation appeared to perform quite well for smoke, though for other smoke like phenomena (e.g. steam) multiple scatter is likely necessary to achieve the best results.
Below are some images and animations I obtained using these techniques:
Some images comparing different amounts of vorticity confinement, showing the smooth flow and the increased turbulence, with a black smoke. These are simulated using a 40x40x40 grid, which is probably too coarse for the resolution rendered here, and leads to a loss of some of the fine scale detail.
An animation of a similar smoke plume with a simulated "wind" blowing and different light scattering properties, yielding a lighter looking smoke.
In image and animation of smoke rising from a container. Rather than model the smoke inside the container, I using a decreasing source term of the smoke density from the edge of the opening, to simulate the smoke "bunching up" toward the edges. This was simulated with a 100x100x100 grid and higher vorticity confinment, and appears with significantly more fine scale detail than the previous images.
While conceptually straightforward, implementing the fluid simulator proved to be somewhat more challenging than I had anticipated, though mostly I think these difficulties were more from lack of experience with numerical methods. For example, neither of my main references on the method specified very explicitly how they were performing the "particle tracing" for the advection step of the process. Initially I used a simple linear step (v * dt); however, this ended up causing problems in regions where the velocity was changing rapidly (smoke density would disappear). I achieved better results, though, after implementing a crude second order Runge-Kutta integrator with an adaptive time step.
Another difficulty was implementing the Poisson solver for the pressure simulation. Foster and Fedkiw  recommend using a conjugate gradient solver to solve the sparse linear system of equations formed by discretizing Poisson's equation. There are some quite good implementations of this freely available; however, it was not clear to me exactly how properly include the boundary conditions in this system, and fearing debugging that might be difficult, I decided to implement the slower, but more straightforward successive over relaxation method described in Foster and Metaxis . Unfortunately, I had problems with the form of the equation they have in their paper, which was resulting in some very odd behavior for the smoke. Strauss , though, has a nice elementary description of the method, and after modifying my system based on that my results seemed reasonable. Still, it required many more iterations to achieve convergence than  reported.
Overall, I was quite happy with the results, though I was I had been able to complete the fluid solver more quickly and render the more visually interesting fire simulations. Probably, the most interesting feature of this method, I think, is its ability to generate smooth realistic looking animations: all we need to do is set up the conditions we would like to see, and let the physics take their course.
 Fedkiw, R., Stam, J. and Jensen, H.W., "Visual Simulation of Smoke", SIGGRAPH 2001, 23-30 (2001)
 J. Stam. Stable Fluids. In SIGGRAPH 99 121–128 (1999)
 N. Foster and D. Metaxas. Modeling the Motion of a Hot, Turbulent Gas. SIGGRAPH 97 181–188 (1997)
 N. Foster and R. Fedkiw. "Practical Animation of Liquids", SIGGRAPH 2001, 15-22 (2001)
 W. Strauss, Partial Differential Equations. John Wiley & Sons (1992)