Practical Fluid Dynamics: Part 1
September 20, 2014
September 20, 2014
Press Releases
September 20, 2014
PR Newswire
View All

If you enjoy reading this site, you might also want to check out these UBM Tech sites:

Practical Fluid Dynamics: Part 1

June 26, 2008 Page 1 of 3

[In this technical article originally printed in Game Developer magazine, Neversoft co-founder Mick West looks at how to efficiently implement fluid effects - from smoke to water and beyond - in video games, with example code.]

Fluid effects, such as rising smoke and turbulent water flow, are everywhere in nature but are seldom implemented convincingly in computer games. The simulation of fluids (which covers both liquids and gases) is computationally very expensive.

It's also mentally expensive, with even introductory papers on the subject relying on the reader to have math skills at least at the undergraduate calculus level.

In this two-part article, I will attempt to address both these problems from the perspective of a game programmer who's not necessarily conversant with vector calculus. I'll explain how certain fluid effects work without using advanced equations and without too much new terminology.

I'll also describe one way of implementing the simulation of fluids in an efficient manner without the expensive iterative diffusion and projection steps found in other implementations.

A working demonstration in ZIP form with source code accompanies this article; example output from this can be seen in Figure 1.

FIGURE 1 Smoke output is shown from the accompanying code.

# Grids and Particles

There are several ways of simulating the motion of fluids, but they all generally divide into two common styles: grid methods and particle methods. In a grid method, the fluid is represented by dividing up the space a fluid might occupy into individual cells and storing how much of the fluid is in each cell.

In a particle method, the fluid is modeled as a large number of particles that move around and react to collisions with the environment, interacting with nearby particles. Let's focus first on simulating fluids with grids.

The simplest way to discuss the grid method is in respect to a regular two-dimensional grid, although the techniques apply equally well in three dimensions. At the most basic level, to simulate fluid in the space covered by a grid you need two grids: one to store the density of liquid or gas at each point and another to store the velocity of the fluid.

Figure 2 shows a representation of this, with each point having a velocity vector and containing a density value (not shown). The actual implementation of these grids in C/C++ is most efficiently done as one-dimensional arrays. The amount of fluid in each cell is represented as a float.

The velocity grid (also referred to as a velocity field, or vector field) could be represented as an array of 2D vectors, but for coding simplicity it's best represented as two separate arrays of floats, one for X and one for Y.

FIGURE 2 The fluid density moves over a field of velocities, with a density stored at each point.

In addition to these two grids, we can have any number of other matching grids that store various attributes. Again, each will be stored as a matching array of floats, which can store factors such as the temperature of the fluid at each point or the color of the fluid (whereby you can mix multiple fluids together).

You can also store more esoteric quantities such as humidity, for example, if you were simulating steam or cloud formation.

Page 1 of 3

### Related Jobs

Infinity Ward / Activision — Woodland Hills, California, United States
[09.20.14]

Senior AI Engineer
Infinity Ward / Activision — Woodland Hills, California, United States
[09.20.14]

Lead Tools Engineer - Infinity Ward
Insomniac Games — Burbank , California, United States
[09.19.14]

Senior Engine Programmer
Nexon America, Inc. — El Segundo , California, United States
[09.19.14]

Front-End Developer

 Anonymous
 Semi-lagrangian advection has been used since a long time in the field of numerical weather forecast. Semi-Lagrangian schemes presented in the literature have used either two or three time levels for time integration (Staniforth & Cot e 1991); the latter are more expensive but provide greater accuracy. In graphic application one can use a faster (second order in time) trapezoidal rule for the integration : x_i - x_i* = (u(x_i, t+dt) + u(x_i*, t)) / 2 This is an implicit equation for x_i* (the position back in time) which could be solved by estimating u using the 2 first term of the Taylor expansion around x_i x_i* = x_i - dt ( 1 + dt/2 Grad u(x_i,t)^-1 * ((u(x_i,t) + u(x_i, t+dt)/2) where all terms on the RHS are computed at grid point x_i. The biggest problem with Stam algoritm is the use of linear interpolation, which is not appropriated for rendering of fluid. One early drawback of semi-Lagrangian schemes was that, since the basis functions of common high-order interpolants such as cubic Lagrange or Hermite are not positive definite, a scalar that is initially bounded between 0 and 1 will not necessarily remain so. However, this problem can be easily fixed by a local clipping, as proposed by Bermejo and Staniforth (1992). This local clipping produces a quasi-monotone scheme. The term quasimonotone means that the scalar field can develop wiggles, but it cannot exceed its range at the previous time step. This is different from traditional clipping, where under-shootsand over-shoots are reset to minimum and maximum background values chosen, based on physical considerations, at the beginning of the computation. As shown by Bermejo (2001) this clipping is equivalent to a linear combination of a higher order interpolant and a low-order (monotone) Lagrange interpolant. To guarantee (discrete) conservation of scalar, two approaches are proposed by Priestly (1993) and Bermejo and Conde (2002). Personnally I prefer the one from Priestly (1993), since it can be used with the quasi-monotone interpolation scheme previously described. The complete reference is Priestly, A. 1993 A Quasi-Conservative version of the Semi-Lagrangian Advection Scheme. Mon. Wea. Rev., 121

 Stéphane Gaudreault
 According to google, patent number: 6266071 is issue since July 24, 2001. http://www.google.com/patents?id=ySgIAAAAEBAJ&printsec=abstract&zoom=4

 Anonymous
 Normally, you should always avoid the case where two trajectories arrive at the same point back in time. The semi-lagrangian method is stable even for a CFL number larger than 1, but this does not mean that it is stable for any timestep. In fact, you are limited by your framerate, but there is also a theorical limit. When the advection velocity varies, the time-integration interval should be bounded for both numerical accuracy and computational stability. A stability condition for semi-Lagrangian schemes in terms of the so-called Lipschitz number or deformational CFL (DeCFL) number has been discussed by Smolarkiewicz and Pudykiewicz (1992) and Lin and Rood (1996). see http://ams.allenpress.com/perlserv/?request=get-abstract&doi=10.1175%2F1520-0493(1996)124%3C2046%3AMFFSLT%3E2.0.CO%3B2 Using this stability condition for the choice of the timestep, any two trajectories cannot cross each other. In this case, there is no need to keep a list of points of origin and mass fraction.

 Anonymous
 Could you post a sample that shows it as a liquid?

 Lou Hayt
 During my education at FullSail I made a liquid dynamics project that was based on Stam's articlesource code. I used the densities to form a dynamic height map (this was mentioned in this article), added a floating boat (not polished though), waves, rain effects etc you can download a zip file with my source code together with Stam's article and source code here: http://louhayt.com/Documents/5.%20Fluid%20Dynamics%20Project%20(One%20month%20class%20project).zip from my understanding this code is light enough to perform well on mobile devices, I was happy to read the possible improvements suggested in this article! :)

 Richard Mathie
 Hey why is the pressure term completely neglected in the naiver stokes equation, it should read Du/Dt =du/dt + (u.grad)u= + nu grad^2 u it dosent matter in the end as the "vorticity confinement, and projection" steps are effectively taking care of it, but the physical reason they are needed is the pressure!

 Richard Mathie
 odd that didnt come out right try Du/Dt=du/dt + (u.grad)u= ## -1/rho grad P## + nu grad^2 u

 Bryan Rayner
 Using this article as a guideline, I've made a simple fluid simulation for a physics assignment. However, there are a few problems which arise later - particularly dealing with the decompression step. A special case exists wherein a thick wave gradually gains speed - the cells before it cause the wave to move faster towards them - since the wave is more than one cell thick, it isn't balanced off by the rear. As it moves to a faster and faster velocity, it begins to move across multiple cells at the same time - resulting in an ever increasing velocity which is not slowed in the slightest. This causes multiple glitches in the system - as soon as a thick wave starts, it will grow and fill the entire fluid (as the decompression algorithm will continue to spread it out in a given direction.) My quick fix - dissalow cells that have more than a given amount of fluid to decompress - the backward advection will pull them out, and they will eventually slow down enough to spread out to a normal state.

 none Comment: