Inclusion of both a forcing term and pressure BC

Lattice Boltzmann Method LBM: Theory
#1

Dear all,

Based on the paper of Guo et al. 2002 named “Discrete lattice effects on the forcing term in the lattice Boltzmann method” I implemented a simulation of a steady Poiseuille flow driven simultaneously by a body force term (rho*g) and a pressure gradient term placed between the inlet and the outlet boundaries ((Pin-Pout)/Lx).
My main objective with this work was to verify that the model proposed by Guo was in fact the most accurate one since as demonstrated in that paper, lattice effects did not appear when recovering the incompressible Navier-Stokes equation.

To do so I followed this procedure:

(1) I fixed the Reynolds number, which provides me the values for rhog and Pin-Pout)/Lx (which I considered to have equal magnitude), and changed the Mach number, i.e I increased/decreased the number of nodes in my domain.
(2) I did this for the Guo model and for a much simpler one that considers the force to be Fi=wiei*F/(cs^2).
(3) Finally, I computed the deviations of the numerical solutions (of both models) from the expected analytical result.

What is amazing is that for the two models and using tau=1 I always achieved the analytical solution to machine accuracy regardless the Mach number used.
Have anyone got the same result? Or am I doing something wrong? If not, what is the explanation for that?

Thank you in advance for the help.

Goncalo Silva

P.S. At inlet and oultet boundaries I used the Zou and He boundary conditions. For the walls I used the traditional bounce-back approach placing the wall haf way between the solid and fluid nodes. My lb model is the D2Q9.

#2

Hello Goncalo,

it would be interesting to see the curve “error as function of tau” for fixed Re and spatial resolution. Have you checked whether you are second order accurate for fixed tau (other than 1)?
How do you compute the numerical error? Do you take the entire computational grid into account? Of course tau = 1 is a special choice, but I find it strange that you recover machine accuracy.

I am curious what your final result will be.
Timm

#3

Hello Tim

First of all thank you for your interest in my problem.

The way I have defined my numerical error was as the sum (over all computational nodes) of the absolute value of the difference between the numerical solution and the analytical one over the absolute value of the sum of the analytical solution, i.e. I used a relative error with the conventional L1 norm definition, as it was used for instance for He and Luo in the paper: “Lattice Boltzmann Model for the Incompressible Navier-Stokes Equation” (1997).

The results that I am reporting were obtained with the model that considers the body force to be Fi=wieiF/cs^2. My implementation with the Guo model for the forcing term yields equivalent results…

So as you suggested me, here it is the relation between tau and the error for Re=10, LY=30 (channel width) and LX=15 (channel length) at steady state is as follows

tau Error Relat
0.51 0.0019225
0.55 0.0019304
0.60000 0.0018953
0.70000 0.0016857
0.80000 1.29450E-03
0.90000 7.29980E-04
0.99000 8.02250E-05
0.99900 8.09370E-06
0.99990 8.10080E-07
0.99999 8.10E-08
1.00000 5.32550E-14
1.00001 8.10170E-08
1.00010 8.10240E-07
1.00100 8.10950E-06
1.01000 8.18040E-05
1.10000 8.87810E-04
1.50000 5.86270E-03
1.70000 9.08980E-03
2.00000 1.46570E-02
2.50000 2.52710E-02
3.00000 3.66480E-02

The error increases as tau value departs from 1 (with the exception of tau=0.51).

Keeping tau=1 and still with LY=30 and LX=15 the steady state solutions yield the following error values as function of Reynolds:

Re Error Relat
1 9.0978E-13
5 2.3235E-13
10 5.32550E-14
15 2.4712E-14
20 3.8433E-14
25 1.1738E-14
30 4.5003E-14
35 2.8743E-14
40 1.286E-14
45 1.2921E-14
50 2.352E-14
60 2.931E-14
70 0.061668
80 0.36923
90 NaN
100 NaN

Note that for Re=70 the Mach is above 0.58. However for Re=60 Mach is 0.53 and the numerical solution is equal to the analytical one to machine accuracy.

Using LY=30, LX=15 but now tau=0.9 the Mach vs Error evolution is:

0.5 0.020547
0.4 0.017771
0.3 0.0099619
0.2 0.0050707
0.1 0.0013894
0.09 0.0011588
0.08 0.00091038
0.07 0.00068389
0.06 0.0005326
0.05 0.00036605

This results in a 1.79 accuracy. If I stick to Mach numbers below 0.1 the slope is 1.94 which is much closer to 2. In my opinion this result is acceptable.

First I thought that the problem was related to the forcing term, so I removed it and checked the results.

So, setting the body force to zero and using just the pressure boundary conditions to drive the flow (for Re=10, LY=30 and LX=15) I will have a similar behavior, i.e. the minimum error is not for the minimum tau (where I have the minimum Mach number) but now for tau=0.9!:

tau Error Relat
0.51 0.0019097
0.55 0.0018667
0.60000 0.0017678
0.70000 0.0014308
0.80000 0.00091209
0.90000 0.0002201
1.00000 0.00063735
1.10000 0.0016526
1.50000 0.0071377
1.70000 0.01062
2.00000 0.016571
2.50000 0.027825
3.00000 0.039844

Doing exactly the same thing, i.e. using only pressure gradients to move the flow, but setting constant tau=1, Re=10 I have:

Ma Error Relat
0.5 0.019608
0.4 0.013699
0.3 0.0077519
0.2 0.0029499
0.1 0.00079936
0.09 0.00063735
0.08 0.00052002
0.07 0.00038565
0.06 0.00028337
0.05 0.00019996

This gives a slop of 2.0155 which agrees very well with the expected 2nd order accuracy.

So my code produces logical results for tau=1 if no body force is used or if a body force is used but tau is different from 1.

If you want I can send you my codes (they are simple and short MATLAB codes)…

Again thank you for your help.

Regards

Goncalo Silva

P.S. I don’t know if this is relevant for the understanding of my problem but, to have a poiseuille flow without a developing region, i.e. to have a linear pressure drop from the inlet to the outlet I defined my body force as:

g=0.5*12.Rekvisc.^2./((LY-2).^3) (0.5 factor is just to account that this body force will only contribute for the flow movement in 50%)

and my outlet pressure drop value as:

rho_outlet=(rho_inlet-3.12.Rekvisc.^2.(rho0).*(LX-1)./((LY-2).^3)) + 3.*g(1).rho0.(LX-1) (in my case rho_inlet was simply set to 1).

So my body force/pressure difference between inlet/outlet changes whenever I change Re, tau, Ma, etc…

#4

Hello Goncalo,

the error behavior for flow without body force is reasonable. I agree. That you have a minimum of your error at tau = 0.9 is not surprising. We have just published a paper, where this is also discussed. Have a look here. Actually, tau around 0.9 is usually the best choice for a LB simulation. There is also a paper dealing with this issue by Holdych et al. (see ref. 16 in our paper).
But I really do not understand why you have machine accuracy for body force driven flow at tau = 1. Right now I do not have so much time. I will have a closer look at your final statements in your last post (concerning choice of body force and pressure gradient) in the next days.

Have a nice time,
Timm

PS: Just a little try: What happens if you compute the L2 error instead of the L1 error?

#5

Hello Tim

Thank you for the paper suggestion. I will read it thoroughly.
Actually, the fact that tau=0.9 is the optimal tau choice is a surprise for me. When I say that I am thinking about Chen, Martinez and Mei paper (On boundary conditions in lattice Boltzmann methods) where they say something like: tau=1 works very well when bounce-back wall boundaries are used but other values of tau may work not so well.
I think I will have to study more on this subject…

Regarding the use of an L2 norm for the error evaluation the results are more or less the same.
Here they are for the same parameters as before, i.e. tau=1, LY=30 and LX=15:

Re Error Relat L2 norm
1 1.69E-13
5 9.64E-14
10 5.31820E-14
15 1.77E-14
20 3.19E-14
25 1.13E-14
30 4.64E-14
35 3.14E-14
40 1.54E-14
45 1.69E-14
50 3.11E-14
60 4.25E-14
70 0.049575
80 0.13625
90 NaN
100 NaN

Machine accuracy is still recover regardless the Re, and consequently the Ma, I use (under certain ranges as Re=90 testify).

I will also spend some more time trying to understand what is happening…I think I will spend some more time looking to the equations and trying to figure out why when the end result of the collision process only considers the equilibrium distribution function plus the forcing term I get machnine accuracy. But when the forcing term is not present and so the distribution functions resulting from the collision process are simply equal to the equilibrium distribution functions the result that I get follow the expected second order accuracy rule…

Once more thanks for your help.

Goncalo Silva

#6

I have some questions regarding your definition of symbols.
What is kvisc? Is it the numerical viscosity kvisc = (tau - 0.5) / 3?
How is your Reynolds number defined?
What is g(1)?
What is rho0?
Please write down the equations with fully defined symbols.

Now, what happens if you only use a body force? Does your code have periodic boundary conditions?
Or to put it into other words: What happens if you change the weight of body force and pressure boundary condition? You have only tested the cases 50% / 50% and 0% / 100%. Just use 30% / 70% and so on and find out in which cases you see machine accuracy.
I also wonder whether this would work in 3D.

Another point: The use of the simple equation Fi=wieiF/cs^2 and Guo’s model has to yield the same result in your case, since the body force is both homogeneous and constant. In this special case, both models are equivalent (if I am not mistaken).

Have a nice Christmas,
Timm

#7

Hello Tim

You are right, kvisc is the kinematic viscosity, defined as you wrote.
I do define my Reynolds number as the product of the mean velocity by the channel width over the kinematic viscosity.
The mean velocity is defined as 2/3 of the maximum velocity being this last parameter defined by me, when I set the Mach number magnitude.
The channel width is given by the number of computational nodes throughout the transverse flow direction. Since Re, Ma and kvisc (tau) are my simulation input parameters the channel width is simply defined as a function of them, i.e LY= round(3/2Rekvisc/Umax)+2 (the plus 2 is just to account for the fact that the wall nodes are not part of my fluid domain and are located halfway).

g(1) is the streamwise component of my acceleration vector, i.e. g(1)=12.Rekvisc.^2./((LY-2).^3) and g(2)=0

rho0 is the constant part of the fluid density since I am using the incompressible model of He and Luo (1997). Therefore, rho=sum(fi) and rho0u=sum(eifi).

I hope now all symbols meaning is clear :slight_smile:

Using only a body force with periodic b.c. at inlet/outlet results in logical L2 errors, i.e. the order of convergence of the model is around 2.
Now it comes the odd part (at least for me). For fixed Re, Ma and tau, but changing the relative contribution of the body force (g) to the overall flow driving forces, expressed by the body force and the pressure gradient (dp/dx), my error changes!
Here it is an example for Re=10, Ma=0.09, tau=0.1 (and LY=30 and LX=15):
g/(dp/dx+g) L2
0 0,00034197
0,1 0,00027358
0,25 0,00017099
0,4 0,000068395
0,49 6,8395E-06
0,499 6,8395E-07
0,4999 6,8395E-08
0,5 3,23E-11
0,5001 6,84E-08
0,501 6,8395E-07
0,51 6,84E-06
0,55 0,000034197
0,6 0,000068395
0,75 0,00017099
0,9 0,00027358
1 0,00034197

How to explain this? Note that for example if I fixe g/(dp(dx+g) equal to 0.75 or to 0.25 I have 2.27 order of convergence, which is acceptable, I say.
So, based on your last suggestion Tim I found out that everything is logical except when tau=1 and g/(dp/dx+g)=0.5 are set.
Furthermore, and regarding your last comment, I think I will have to agree with you. Analyzing eqs. (21a) and (21b) of Guo’s paper and since g (or F) is constant (spatially and temporally) the Navier-Stokes equations are recovered without any additional terms when using Fi=wieiF/cs^2.
However, Guo et al. used this particular flow model to demonstrate that their body force definition was the most accurate one. Theoretically, this body force definition should be as accurate as any other present in their paper, shouldn’t it?
Finally I would like to thank you once for suggesting me your paper. Now that I read it I must say that, in contrast with some LB papers, it is really easy to follow. I learned a lot from reading it :slight_smile:
Regards,
Goncalo

#8

Hi Goncalo,

it is really interesting what you have found out. I am not aware of any publication reporting similar results. I have never seen anybody using both body force and pressure BC either. However, there is one paper I have read which is about analytic solutions of the LB equation. Maybe it might be useful. Unfortunately, I do not know the author nor the title. And for the next ten days I won’t be in my office. But when I am back I can tell you.
However, I believe that your results can only be recovered using an incompressible model. I am pretty sure that the standard LB model will not show a machine accuracy under similar circumstances. Moreover I believe that in a more complex geometry the error might be larger. And it is an open question whether your observation could be recovered in 3D.
The reason for the machine accuracy is unclear at this point. Maybe some simple thinking could explain it? Also note that tau = 1 is a special case in the LBM. You should look for some more theoretical articles analyzing the background of the LBM.

Timm

#9

Hi Tim

I downloaded a paper where analytical solutions of the LBM for some simple flows are studied, the paper name is: “Analytical Solutions of Simple Flows and Analysis of Nonslip Boundary Conditions for the Lattice Boltzmann BGK Model” by He, Zou, Luo and Dembo. I don´t know if is this the paper that you were refereeing to.

You are right Tim, the machine accuracy situation only takes place when the incompressible model is used. If the standard equilibrium distribution function is used the model accuracy behaves as it should, i.e. second order regardless the tau and g/(g+dp/dx) value. This raises me a somehow “philosophical” question: Considering that the general purpose of the LB model is to simulate the incompressible Navier-Stokes equations and that the incompressible model is the one that generally yields more accurate results (at least it is that what my experience suggests) why people still use the standard model? For example, why in your paper did you use the standard model? Wouldn´t the results with the incompressible model be much closer to the analytical ones? Does that have anything to do with the fact that you are simulating in 3D? (this explanation does not make much sense to me but since I am still sticking my studies to 2D problems maybe there is something with the extra dimension that I am missing.)

Finally, and for those that have already read the paper “Discrete lattice effects on the forcing term in the lattice Boltzmann method” by Guo et al. (2002), is there someone that could explain me why have these authors used a Poseuille flow driven simultaneously by a pressure BC and a constant body force as a benchmark test to demonstrate that their model was superior? I do not know if you ever thought about it but I think that the pressure BC is not doing anything, periodic boundary conditions at inlet/outlet would be enough…

Thank you for the help.

Goncalo

#10

I downloaded a paper where analytical solutions of
the LBM for some simple flows are studied, the
paper name is: “Analytical Solutions of Simple
Flows and Analysis of Nonslip Boundary Conditions
for the Lattice Boltzmann BGK Model” by He, Zou,
Luo and Dembo. I don´t know if is this the paper
that you were refereeing to.

No, this is not the paper I meant, but I think that the paper I want to show you is not very important for your needs anyway.

You are right Tim, the machine accuracy situation
only takes place when the incompressible model is
used. If the standard equilibrium distribution
function is used the model accuracy behaves as it
should, i.e. second order regardless the tau and
g/(g+dp/dx) value. This raises me a somehow
“philosophical” question: Considering that the
general purpose of the LB model is to simulate the
incompressible Navier-Stokes equations and that
the incompressible model is the one that generally
yields more accurate results (at least it is that
what my experience suggests) why people still use
the standard model? For example, why in your paper
did you use the standard model? Wouldn´t the
results with the incompressible model be much
closer to the analytical ones? Does that have
anything to do with the fact that you are
simulating in 3D? (this explanation does not make
much sense to me but since I am still sticking my
studies to 2D problems maybe there is something
with the extra dimension that I am missing.)

Incompressible models only work for steady state simulations, if I remember correctly. I am not completely sure, and I am also very interested why people do not use the incompressible models more often. Could someone please enlighten us?

Finally, and for those that have already read the
paper “Discrete lattice effects on the forcing
term in the lattice Boltzmann method” by Guo et
al. (2002), is there someone that could explain me
why have these authors used a Poseuille flow
driven simultaneously by a pressure BC and a
constant body force as a benchmark test to
demonstrate that their model was superior? I do
not know if you ever thought about it but I think
that the pressure BC is not doing anything,
periodic boundary conditions at inlet/outlet would
be enough…

I must have missed this point when I (diagonally) read the paper. Pressure boundary conditions of course have an effect on the flow, so it is not correct to say that they do nothing.
I see no reason why they use both methods simultaneously. I have to read the paper again (more thoroughly) when I am back in my office.
Some funny thing happened to me when I first implemented Jonas Latt’s velocity boundary conditions: I forgot to adjust the density at the inlet and outlet. As a consequence, the inlet and outlet densities have been the same and the resulting flow velocity was exactly 50 % of the expected value. This means that 50 % of the flow is driven by the pressure gradient and 50 % by the velocity of the incoming and outgoing fluid portions.
I think one could spend years with understanding the subtleties of LBM.

Timm

#11

Hello Tim

Even if you think that the paper that you first wanted to suggest me may not help… please, tell me its name when you come to back to your office and have time to… I always like to give papers a quick read… sometimes I find some interesting things… :slight_smile:

Regarding the incompressible model I believe that the one that I am using (He and Luo) is valid also for unsteady state flow regime…
He and Luo (1997) besides testing their incompressible for a steady Poisueille flow (where, by the way, using pressure bc and tau =0.75 they also achieve machine accuracy and with tau=1 or 2 they claim to have the same results) they also test their model for a Womersley flow achieving better accuracy than with the standard one. The authors say that the applicability of their model to unsteady flows is always valid provided that the time scale L/cs (i.e. ratio of the characteristic macroscopic length to the sound speed), which represents the time needed for a sound signal to travel a distance L, is much smaller than the time in which temporal macroscopic changes are expected to occur.

However when we talk about incompressible models we may be talking about different ones. At least I know three different ones:

  • Zou, Hou, Chen, Doolen, J. Stat. Phys. 81:35 (1995) - The model proposed in this paper condensates the parameter rho*u into a single parameter J (named frequently as current density, see for instance Succi’s book). The model only recovers steady state Navier-Stokes equations. (Maybe it was this one that you were talking about…)

  • He and Luo, J. Stat. Phys. 88: 927 (1997) - This is the model that I use. Since I am just interested in steaty state solutions I could also use the first model of Zou et al. However I tested both and found out that He and Luo model besides being numerically more robust it is also more accurate. (This conclusion was solely based on a simple Poiseuille flow test…)

-Guo Z et al. J. Comput. Phys. 165: 288 (2000) - According to the authors this model is theoretically an improvement of He and Luo incompressible model as it does not have the constraint of the macroscopic time scale being much higher than L/cs. I am not 100% sure but I think that this is the model major advantage… In this paper the authors perform a lot of unsteady state flow tests with excellent results. Nevertheless, I have never seen any other article, besides this one, where this model is employed…I don’t know the reason for that… I also have never implemented this model so I don’t have much of an opinion about it… (I found it recently…)

I understand that in situations where there are no density (pressure) gradients inside the flow domain, as it happens when using a constant body force and periodic b.c. in a straight geometry, the standard lb model and all incompressible models should yield exactly the same results with the same accuracy. However, my experience says that when pressure gradients inside the domain occur all incompressible models, for equivalent Ma numbers, yield more accurate results (being ones more accurate than others).

Regarding my last comment where I say that the pressure BCs are not doing anything on Guo et al. (2002) test. What I meant to say was that for the purpose of their study the use of both a body force and pressure boundary conditions is, in my opinion, unjustified. If I understood it well, the fact that pressure b.c. are being used do not add or change anything to their conclusions…

Thank you for your suggestions and comments Tim. Enjoy these week holidays.

Regards

Goncalo

#12

Hello

I understood why Guo et al. used the pressure boundary condition.
Since F=rho*g and g is constant (in their test) they introduce the pressure b.c. so that rho changes throughout the flow field which will make F changing in the same proportion…

Goncalo

#13

Hi Goncalo,

the paper I was talking about:
Gabor Hazi: Accuracy of the lattice Boltzmann method based on analytical solutions, Phys. Rev. E 67, 056705 (2003)

Timm

#14

Hello

Thank you for the paper suggestion Tim.

I will spend some more time trying to figure out the subtleties of Guo’s model to account for the body force presence and maybe extend my study to 3D cases.

I will keep in touch whenever I find some more interesting things.
Moreover, if I understand (or find some article speaking of) the advantages/disadvantages of the incompressible model I will also say something to the forum… I think this is really an interesting open question, at least for me it is. I mean, if the incompressible model is so much better (or at least slightly better) than the standard one why people still worry to use the standard model in their studies of incompressible fluid flow?

Have a nice 2009 :slight_smile:

Goncalo

#15

hi GoncaloSilva
in this text at 24-12-2008 you said that you used tau=0.1 ?!!!
I read that tau must be greater than 0.5
if I am wrong please tell me
thanks

#16

Dear Seyedmehdi

You are absolutly right! I errouneously wrote tau=0.1. In fact what I meant to write was tau=0.9. For the kind of flow in which I was performing tests, i.e. a Poiseuille flow with both Pressure BC and a body force driving it, I found that this tau value was the most accurate one.
This conclusion is supported by R. Noble, J. G. Georgiadis, R. O. Buckius, J. Comput. Phys. 193, 595 (2004) and the recent work of T. Krüger, F. Varnik, and D. Raabe, Phys. Rev. E (accepted for publication) (you can find the link for this last paper somewhere else in this forum as it is Timm paper)

I hope my post helps and sorry about the confusion that I have probably created.

Regards

GoncaloSilva