Accounting Advection

If a velocity field is not mass conserving, then some points will have multiple velocity vectors from other points pointing toward them. If we simply move our scalar quantities (like density) along these vectors, there will be multiple quantities going to (or coming from) the same point, and the result will be a net loss or gain of the scalar quantity. So, the total amount of something such as the density would either fade to zero or gradually (or perhaps explosively) increase.

The usual solution to this problem is to make sure the vector field is incompressible and mass conserving. But as mentioned before, it's computationally expensive to implement. One partial solution is to make the advection step mass conserving, regardless of whether the velocity field is actually mass conserving. The basis of this solution is to always account for any movement of a quantity by subtracting in one place what is added in another.

Advection uses a source and destination buffer to keep it independent of update order. In Stam's implementation, the destination buffer is simply filled one cell at a time by combining a value from four cells in the source buffer and placing this value into the destination buffer.

To properly account for compressible motion, we need to change this copying to accumulating, and initially make the destination buffer a copy of the source buffer. As we move quantities from one place to another, we can subtract them in the source and add them in the destination.

With the forward advection in Figure 3, we're moving a quantity from point P to points A, B, C, and D. To account for this, we simply subtract the original source value in P from the destination value in P, and then add it (interpolated appropriately) to A, B, C, and D. The net change on the destination buffer is zero.

With the reverse advection in Figure 4, as used by Stam, the solution would initially seem to be symmetrically the same: Subtract the interpolated source values in E, F, G, and H from the destination buffer, and add them to P.

While this method works fine for signed quantities such as velocity, the problem here is that quantities such as density are positive values. They cannot drop below zero, as you cannot have a negative quantity of liquid.

Suppose that point E was one source point for two destinations P1 and P2, both of which wanted 0.8 of E. If we follow our initial plan and subtract 0.8*E from E and add 0.8*E to both P1 and P2, the net effect is zero, but now the value at E is negative. If we clamp E to zero, then there is a net gain of 0.6*E. If we subtract 0.8*E from the source value of E after updating P1, then when we update P2, it will only get 0.8*0.2*E, when clearly both P1 and P2 should get equal amounts. Intuitively, it seems they should both get 0.5*E, and the resulting value in E should be zero, leading to a net zero change.

To achieve this result, we create a list that for each point, noting the four points of origin for each and the fraction of each point they want. Simultaneously, we can accumulate the fractions asked of each source point. In an ideal world, this would add up to one, as the entire value is being moved somewhere (including partially back to where it started).

But with our compressible field, the amount of the value in each point that is being moved could be greater than or less than one. If the total fraction required is greater than one, then we can simply scale all the requested fraction by this value, and the total will be one. If less than one, the requesting points can have the full amount requested. We should not scale in this case, as it will lead to significant errors.

With the mass conservation of advection fully accounted for in both directions, it turns out that neither forward nor backward linear advection alone will produce smooth results. After some experimentation, I determined that applying forward advection followed by backward advection worked very well and gives a smooth and artifact-free flow of fluid over a compressible velocity field.

Now What?

We can now perform both forward and reverse advection in a mass-conserving manner, meaning we can move fluid around its own velocity field. But even though our velocity field does not need to be mass-conserving, we actually still want it to be, since the velocity fields of real world fluids generally are incompressible.

Stam solves this problem by expensively forcing the field to be fully mass conserving after every change. It's necessary, since the reverse advection requires it. The key difference is that since our advection step does not require the field to be mass-conserving, we're really only doing it for cosmetic purposes.

To that end, any method that rapidly approaches that state over several time steps will suit our purpose. That method, and the method of diffusion, can be found in the accompanying code, and I will discuss how they work in the next instalment of the article.

[EDITOR'S NOTE: This article was independently written and commissioned by Gamasutra's editors, since it was deemed of value to the community. Its publishing has been made possible by Intel, as a platform and vendor-agnostic part of Intel's Visual Computing microsite.]