choice of units in thermal lattice Boltzmann

Dear all,
I’m doing simulations of natural convection in a square cavity, using the model from [1]. The velocity field and temperature field are solved using two independents BGK equations, coupled through a forcing term (boussinesq term).
My goal is to compare the results from my Lattice Boltzmann code to those of commercial CFD solvers (like Fluent for example).
To do these comparisons, a square cavity with a hot wall on one side and a cold wall on the other side is simulated for different Rayleigh number and Prandtl number.
So far, the simulations between the Lattice Boltzmann method and the others method don’t agree. But I think the only problem is that I’m not computing the Rayleigh number properly, as the results of my LBM simulation match those of Fluent but for different Rayleigh number.
Can someone help me in computing the correct Rayleigh number and the corresponding LB parameters?
This is how I compute the Rayleigh number at the moment:
The article [1] states that the kinematic viscosity nu is linked to the relaxation time tau through

nu = (tau-0.5)/3 * dx^2 / dt

and the thermal diffusivity is determined by

D = (tauT-0.5)/2 * dx^2 /dt

(where tauT is the relaxation time for the thermal lattice BGK)
The Rayleigh number Ra is defined by:

Ra = g*betta*deltaT*L^3 / (nu*D)

where g is the gravity force, betta the thermal expansion, deltaT the temperature difference between the hot and cold wall, and L the size of the cavity.
So in term of lattice units, if we suppose dx/dt = 1 and dx = 1/(Nx-1) (where Nx is the number of node along x axis) it should be:

Ra = g*betta*deltaT*(Nx-1)^2 * 24 / ( (2tau-1)*(2tauT-1))

and the Prandtl number Pr should be

Pr = nu/D = 2/3 * (2tau-1)/(2tauT-1)

What is wrong with this analysis, and how the Rayleigh number should be computed?
Thank you

Note: I’ve already read the article “choice of units in lattice Boltzmann simulation” from Jonas latt, but it didn’t help.
[1] A coupled lattice BGK model for the Boussinesq equations. Guo et al.

Does anyone have any idea?
I studied the evolution of the Rayleigh number when a parameters varies (by matching Fluent simulations).
varying g*betta -> Ra varies linearly
varying DeltaT -> Ra varies linearly
varying resolution Nx -> Ra variation is cubic
The last variation contradicts my formula, should the formula be the following?

Ra = g*betta*deltaT*(Nx-1)^3 * 24 / ( (2tau-1)*(2tauT-1))

This doesn’t seems right either…



nu = (tau-0.5)/3 * dx^2 / dt


D = (tauT-0.5)/2 * dx^2 /dt

are lattice model dependant
make sure that you use the lattice model corespond to those relations
you can use d2q9 for fluid d2q4 for energy.

Second, there is a relation to insure that you will get stable and correct simulation.

u_th=sqrt(g*betta*deltaT*(Nx-1)) = sqrt(Ra*D*nu/(Nx-1)^2)

where u_th be less than 0.1

always I follow the following steps:
1- we can assume Nx. Ra and Pr should be known

D=(Nx-1)*sqrt( * 0.01/(Ra*Pr))

2-once you know D you can get tauT and nu

3-once you know nu you can get tau

4-check for the taus if they are in the valid range. if they are not repeat the previus steps with another Nx


5-deltaT is scalar assume it is unity


where gbeta is the forcing term
it can be appliey in different ways.

With those steps your units should be ok

another issue if you want to compare to commercial pakage you need to convert LBM results
to nondimentional or dimensional results. if this is your problem I can not help you I am not very familair
with this stuff.