high prandtl number simulation!

Dear friends

I work with lattice Boltzmann methods, but have been only using the
single relaxation time BGK approximation for my collision operator. As
a result, I am stuck with Prandtl no of 0.5. Now, I would like to
simulate flows with very large Prandtl number, like oils, which have a
Pr=1000, and also liquid metals which have a Pr of 0.01.

Is there a suitable LB scheme which allows me to do so in a simple


what kind of Lattice Boltzmann algorithm you use?  My understanding is that you are using the single distribution function approach .. and then you only have one relaxation time that you use to calculate both thermal diffusivity and kinematic viscosity... and then no way to  change your Prandtl number.

Would you like to take a look at the 2 distribution approach? In this forum there should be a matlab code that uses this method to deal with Rayleigh Benard Convection.

Using this method you should be able to change your Pr number in a fair range of values.


hi tovarish
no, i use two disterbution and use tow relaxation, i use viscosity to calculate the relaxation time for flow field,and for thermal,use thermal diffusivity,so when the prandtl number is high,the solution converge and unstable!!bequase one of the relaxaion time get too close to 0.5!!


I am using two disterbution model . I want to simulate heated square cavity at different Rayleigh number.

the error in Xvelocity is proportional to Rayleigh number.

Ra My resule Papers result
1e3 3.643 3.644
1e4 16.98 16.178
1e5 43.21 34.73

I think two disterbution model has a limitation.
In the papers the best way to couple LBM and Themal application by using Hybird LBM

you will find a good information in those papers:

1- Viscous flow computations with the method of lattice Boltzmann equation
2-A Hybrid Lattice Boltzmann Finite Dierence Scheme for the Diusion Equation
4-Theory of the lattice Boltzmann method: Acoustic and thermal properties in two and three dimensions

Hi eshanf…

  can  ask you what's the range of Rayleigh number that you would like to play with?


Hi Andrea

In low Ra,1e3,1e4,1e5!
I use this scheme!
Tau(V)=3viscosity+0.5 (for flow field)
diffusivity+0.5 (thermal field)

I set the viscosity and Pr number and then calculate diffusivity and two relaxation time whit them!
and use these equilibrium distribution function:


whats your idea about that?

Hi SaS
I’ll try those paper
thank u a lot

Hi Ehsan

   I think you are doing fine.   I played a lot with the model that you use to simulate Rayleigh Benard Convection and I could get good results for a faire range of Pr and Ra values  (up to Pr=10^3 Ra=10^7 without too much numerical effort)  ..    To increase the resolution of the grid and use smaller time steps helps to push more your simulations but then you have to ask a lot more to your computer. This is the case when you want play with high Ra number. The size of your thermal boundary  layer, in fact, will decrease with increasing Ra and you definitely want to be able to catch some of the features in the boundary layer (thermal instabilities will arise from there).....

BUT in the range of Ra number that you are playing everything should be relatively smooth… For the problem of Rayleigh Benard Convection I didn’t have really any stability problem when I was playing with the relaxation time for the numerical scheme of the temperature of the order of 0.501 and smaller.

What I was doing actually was to fix the relaxation time of the fluid = 1 and then, using the definition of Pr to get the value of the relaxation time for the other scheme.

Anyway …
I remember that once Mr Latt posted an inspiring .pdf where he was showing how the choice of Rayleigh and Prandtl and the relaxation times where related to deltaX and delta time. I think that this .pdf could be helpful… just to be sure that setup of your simulation looks fine. Is this .pdf still around Mr Latt? how could not find it in the lbmethod pages.

hi all

tovarish can I ask you one question?

can you tell me the normalized x and y velocity at pr=0.71 and ra=1e5 ?

I want to compare with my results.

best regards

Hi SaS…

         you mean the average velocity of the system during Rayleigh-Benard convection right?  umm ..    yes I could but why you do not compare your results with the Nusselt number instead of my numerical results? 

the Nusselt number shoud be:
N u = 1.46(Ra/Ra) ^0.281 see 5. F. M. Richter, H. Nataf, and S. F. Daly, J. Fluid Mech. vol.129, (1983) pp.173 for example.
(you can also find in books the exponent= 0.3 obtained from linear stability analysis of the problem)

Given that the Nusselt number is defined as the convective divided by the conductive heat transfer, you could know if the average vertical velocity in your system is right.
Using the matlab language, you could calculate the Nusselt number as:

Nu = 1. + sum(sum(uy_Nu.T))/(lxk*(Thot-Tcold));

where uy_Nu is the vertical velocity of your system (same direction of the applied gradient of temperature),
k is the thermal conductivity and Thot and Tcold are the temperature applied to the lower and upper wall respectively.

let me know if it helps

I just want to add something SaS …

I would not be surprised if with Pr=0.71 and Ra10^5 you would already get a periodic solution with time (You would not get a constant Nu number then ) … but if you average over a large enough window of time then you would get results which should be close to the one given by the relation that I wrote before.

In Raleigh Benard Convection there are some different flow regimes corresponding to different bifurcation points.
look at this paper, I like it a lot. It is about lab. experiment with Rayleigh Benard Convection (they are not in Boussinesq approximation ) … but it shows pretty well the transition in convection behavior that you should expect increasing the Ra number (The also play with high Pr number) and gives many ideas for the investigation of the Rayleigh Benard convection problem

here a link of the paper from Pr. Manga home page: http://seismo.berkeley.edu/~manga/paper34.pdf

Hi Andrea
thank you very very much!its really helpful for me!
i use dlx=1,dly=1,dlt=1 in my simulation, when i run my simulation in low tau(near to 0.5,for example 0.52!), one of node in the domain get high value!and get the solution unstable!!(for example i use 100*100 grid size for cavity simulation)

i don’t know the following way get me a wrong value for Ra,Pr,… and other non dimensional parameter?
i use the arbitrary value for V,D,v(viscosity) for calculating these, that D is the number of grid(for example the height of the cavity),V<0.05 and use v(viscosity) ~0.01 for reynolds number, or for natural convection, first i calculate the gbeta(=viscosity^2Ra/(D^3Pr)), and then calculate the force term Force=gbetawi3*DeltaT, and add it to collision step!
i don’t know these work is true or i should choice them in other way!?

I should say that i simulate the natural convection in cavity and i could get a good Nu number in Ra=1e3,1e4,1e5,1e6!!
i don’t know why when i use high Pr number, my simulation get unstable!!!

well …

In some way what you see is normal.  If you play with too small relaxation time for the Advection Diffusion (AD) scheme.. the solution will start to get unstable. When with thermal force convection problem this is often the case... and the value of relaxation time where the solution tends to get unstable is around 0.51.

But now when I played with Rayleigh Benard Convection. Indeed the stability of the advection diffusion scheme is problem dependent… and for the case of Rayleigh Benard convection I was able to decrease a bit more the value of relaxation time and play with Pr number =1000. But I don’t know for the case of the cavity problem. I would not be surprised if the limit in time for this problem is close to 0.52. Anyway. What you can try is … first to increase the size of the grid and see if, with a better resolution you can decrease the value of tau for the AD scheme. Second try to change boundary conditions.

I underline again that I played with high Pr number fluid for the case of Rayleigh Benard convection only.
Aynway … the AD scheme that we are using is very very simple but not so good maybe. You would like to take a look at other AD Lattice Boltzmann schemes?


hi again
thank u, i’ll try!

Hi Andrea
can you tell me whats the restriction in choosing viscosity value in LBM?
or what is the main cause of unstability in LBM?

Hi Dears
I want to calculate Nu Number over a circular cylinder with constant temperature.
but I don’t know how use interpolation for that.
best regard