Equilibrium functions are interpreted as particle directional probability and therefore are allowed to take positive values; however we all know that during an LBM solution, sometimes they take negative values. Dellar (2002) has a paper that discusses this issues and believes that it is ok. However I still have problem understanding that. Can anyone help me?

For BGK the values of the probability distribution function in a bulk (don’t examine boundary conditions) for the advection-diffusion and simple hydrodynamic flow have to be positive for stability issues. It’s pure mathematical von Neumann stability analysis.

In terms of physics it has to be positive as well. It’s exactly where the entropic LBM (I. Karlin) emerges or limiters of the LBM (A. Gorban).

However, for TRT/MRT the values can be negative giving stable and physical stability picture.

However, in principle to be on a safe side it’s always good to insure them to be positive or in other words obtain such conditions that they are positive.