Once, before I knew what a partial differential equation is, I tried to program a fluid simulation. The resulting fluid was highly compressible, as for want of necessary knowledge no Poisson equation was employed.
Several years later, I once more took in hands the little old program and made something more extensive out of it. Now one can choose among various predefined setups. The most interesting extension was that the fluid now could interact with itself via potential fields. With this, I could realize things like star clusters and crystal lattices.
A couple of years later still, I ported the program from Java to C++ and added incompressible fluids and generation of videos, furthermore I added modes that had little to do with fluids, like SmoothLife, reaction–diffusion–advection systems and animations based on Perlin noise.
You can download the C++ source code for the program here. The program uses Qt and in most modes runs smoothly on my 4 core 8 hyperthread 2.2GHz 64bit processor.
You can influence the fluid with the mouse and modifier keys. The exact meaning of the mouse buttons and combinations may vary with the selected mode, but for most of the fluid modes the assignment is as follows:
The button Clear
resets everything to zero in most modes, however in some modes it just stops the motion.
In the following I present some of the predefined modes of the program. Some of the example videos I did not embed from YouTube, as YouTube's video quality is too low.
This mode essentially corresponds to the original program, except for coloring and size.
Example video:
Here I have added two centers of gravity that let the fluid whirl around. Additionally, the fluid is accelerated in the direction of its motion, so that after a while it becomes fast enough to leave the vortex.
Example video:
If the fluid interacts with itself via a regularized Coulomb potential, it forms clusters that attract each other like stars according to the Newtonian law of gravity.
Example video:
That even works to a degree if certain simulation parameters are set to extreme values. Now, the stars sometimes explode.
Example video:
Similar to the previous, but with more and smaller stars.
Example video:
Here the stars are tiny and sometimes explode into gas. The resulting blast leaves an empty space, which leads to interesting patterns.
Example video:
The same again with a bump mapping shader. The empty spaces from the explosions now look like little craters.
Example video:
Here is the result of a big simulation with many tiny stars. The simulation itself did not run as smoothly, however.
A Coulomb potential pulls the fluid into the center of an agglomeration. But if one uses a potential that pulls the fluid into a ring shaped area around the center, the result is hexagonal crystal lattices. These crystals can collapse under their own gravity. What remains is a degenerate lump. If two such lumps meet, one or both are torn apart into a gas cloud from which new crystals form.
Example video:
Here is the whole thing similarly in larger size.
If one uses a potential that has several ring shaped minima with the right radii, other crystal structures can arise. Here are some with pentagonal lattice defects. Additionally, you can see the overall potential (green) and the strength of the force field (blue).
Example video:
In this mode there are various kinds of objects: Gas, small suns
, crystals, seemingly rigidly rotating lumps and pulsating mega-suns
in which the lumps are melted down.
Example video:
This mode bears a superficial similarity to starburst
, but does not feature lumps or mega-suns.
Example video:
Here I used repulsive instead of attracting ring potentials for a change.
Example video:
If one additionally uses a small attracting ring, one gets crystal structures as soon as the density of the fluid is high enough. If the density decreases, the crystals dissolve again.
Example video:
If you take two attracting ring potentials with the golden section as their radius ratio, not only hexagonal crystals are possible, but also pentagonal quasicrystals. In this mode, the hexagonal crystals are stable at high densities whereas the quasicrystals dominate at low densities. If the density decreases, the hexagonal phase morphs subtly into the quasicrystalline phase.
Example video:
This complex self-interaction force field with repelling and attracting components leads to a multitude of stable forms that can bump into each other. Some of them are asymmetrical in such a way that they accelerate in a certain direction.
Example video:
In this mode the fluid has an interesting structure, and there are many crystalline phases. Some of these are very rare and hard to observe. Here's an interesting thing you can do: Add fluid until it crystallizes, then generate two horizontal streams in opposite directions. The crystal liquefies because of the motion, but after a while the middle regions of the two streams turn crystalline and rigid again, because there the flow has almost no shear. A little later an instability arises: Both crystalline beams bend and finally liquefy completely. The flow takes the form of two counter-rotating eddies.
In the video you can also observe the behavior of quasicrystal domain walls, starting from second 1:34.
Example video:
This is the mode I just described, but with bump map visualization. With this, the crystalline phases don't look too great, however the fluid does.
Example video:
... or is it coffee with milk? But I do not drink coffee. So fire it is.
In this mode, the fluid is only almost incompressible, there is a Coulomb-like self attraction of bright areas. Therefore dark areas grow while bright areas shrink. The advection step does not conserve the amount of bright substance. If you use the mouse to paint a bright trail, it contracts while whirling interestingly.
Example video:
This mode is defined almost the same way as waft / oversaturated / many phases
, but with incompressible, substance conserving advection.
The simulation runs rather slowly in size 1024×1024, so better watch the video:
This mode implements the Gray-Scott model with parameter values that lead to replicating dots. (k = 0.0570 und f = 0.0220).
Example video:
I didn't come up with a better name for that. This reaction–diffusion equation does not go along well with advection: If you move the fluid, everything goes dark. After you stop the motion, the system recovers. In the process, the streamlines of the motion become visible.
Example video:
To this reaction–diffusion equation I have given the official name The striving for harmony bears in itself the seed of discord
.
Why? There are blue and orange regions that constantly change their form.
Whenever too big a homogeneous region forms, a part of it spontaneously
changes to the other color.
Example video:
Technical details: The reaction part of the equation is, in every point, the differential equation for the Lorenz attractor. This attractor consists of two loops (represented in the program by the two colors). The state vector of a system ruled by this differential equation moves around one of the loops and from time to time changes to the other loop if it gets to a critical location. The diffusion part of the RD equation couples each Lorenz system belonging to a point to the neighbouring points. When the system gets to the critical part of the loop, but its neighbour systems are at a non-critical part of the same loop, the influence of the neighbour systems keeps it on the loop. However if the systems in a region have sufficiently similar phase, they arrive together at the critical location and it is possible that they all change to the other loop together.
Under lightning and shock waves, a black hole grows bigger and bigger.
Video:
SmoothLife is a generalization of
Conway's Game of Life
from discrete grids to continuous functions.
How the aliveness
of a point changes with time depends on the total aliveness inside one disk shaped and one surrounding ring shaped neighbourhood.
You can zoom the simulation data using the mouse wheel. When zooming, the sizes of the neighbourhoods are adjusted, so that the patterns behave similarly to before.
This mode implements SmoothLife as it is described in this paper, but with a bump map shader used for rendering.
Video:
This mode implements SmoothLife, as it was demonstrated in this YouTube video. The formula for the transition function has a subtle difference to the one from the paper. Also, the time evolution is theoretically continuous (Practically, it still boils down to discrete time steps of Euler integration).
Video:
This SmoothLife transition function produces comically wobbling gliders.
Video:
These modi are described on their own page: Wave equation with fiber potential.
