I’m currently trying to implement the Zou and He pressure boundary condition but the result that I get looks wrong. On the following pictures you see the z-velocity and the density in a channel that has uniform pressures defined at z=0 and z=z_max, the higher pressure beeing at z=0 (even though the z-velocity might suggest that it is the other way round):
The images correspond to the state after streaming. I’m trying to resolve this problem for three days now without success and hope that someone here can help me. The code I use looks like this (this code is for the boundary at z=z_max, D3Q19, numbering of f’s as in OpenLB):
//bounce back of nonequilibrium part:
//redistribution of Momentum:
feq=equilibrium(iPop, rhoIn, u0, t1 );
I don’t understand the very beginning of the code: you successively attribute two different values to u0. Which one is the right one? And, shouldn’t rhoIn contain the sum of the incoming f’s, instead of 1?
Oh, my mistake: I didn’t realize this was a pressure boundary condition. Then, my question is, are you using the Skordos trick to reduce round-off errors, or are you using plain BGK? In the latter case (plain BGK), I think there’s a 1 in excess in your formula. That means, the expression
should be replaced by
Also, I would have the tendency to say that this formula should be multiplied by a factor -1, but I am unsure of the orientation and sign of the wall in your problem.
Actually, we recently wrote a review paper on boundary conditions, in which all formulas are stated step-by-step. You may want to look at equation 27 in the preprint, and solve this equation for rho, to cross-verify the formula in your code.
I think that I have the same as you but with the D2Q9 model, in fact there a inreasing in pressure at vicinity of the outlet (boundary of low pressure).
have find a solution?
I noticed that the problem appears only when pressure gradient rises.