I’ve been wondering about the equilibrium distribution function f[sup]eq[/sup][sub]i[/sub] used in the BGK collision step. As I understand the evolution of f ist determined by
f[sub]i[/sub](x+c[sub]i[/sub], t+1) = (1-1/tau)f[sub]i[/sub](x, t) + 1/tauf[sup]eq[/sup][sub]i[/sub](x, t),
f[sup]eq[/sup][sub]i[/sub] = rhot[sub]i[/sub](1 + 3*c[sub]i[/sub]u + 9/2(e[sub]i[/sub]u)[sup]2[/sup] - 3/2u[sup]2[/sup]).
The last equation can be found in bgk_ma2_equilibrium, but if I look at bgk_ma2_collision_base something different seems to be used (
Can anybody explain me what's going on? Thanks!
Actually, the code you are referring to implements the equilibrium, using exactly the formula you are mentioning.
It is possible that you were unable to match the code with the formula because you were unaware of the fact that the cs2 (“square-speed-of-sound”) constant in the code has value of 1/3 (and corresponds to all “3” terms in your equation).
Or, more likely, you couldn’t do the match because the code has been optimized for speed and is therefore difficult to read.
All collision terms and similar things in Palabos have a generic and a few specialized versions. The generic version is usually nice to read, and works with all lattice structures (D3Q15, D3Q19, etc.) . The specialized versions are written for a specific lattice structure (e.g. D3Q15). They are algorithmically equivalent to the generic version, but usually much faster. The drawback is that they look somewhat cryptic. From a didactic point of view, it is therefore better to look at the generic version. In your case, you can find it at the URL
Note that there is a variable called “rhoBar”. Simply think of this variable as being the density rho, and you’ll find an exact match with the theory. In practice, the variable “rhoBar” can be equal to either rho or rho-1, depending on the round-off policy chosen in a given Palabos program. But that’s another story, and it’s explained on the following page: