Dear LB fellows,
I’m currently working on implementing LBM to solve the 2D lid driven cavity flow.
In the near future, I will be required to run simulations at relatively high Reynolds number and will be interested in computing the strain rate tensor S.
This guided me towards using regularized boundary conditions for the top moving lid. For consistency reasons, I also chose to implement the other boundary conditions (no-slip straight walls for south, west and east) with regularized BCs.
Even though I get a satisfying velocity field with regards to [1] for the stationary regime*, my simulation experiences a mass loss. At each time step, the total mass contained in the cavity decreases. For one time step, the mass loss is small, I’m talking of something around (1.10^-7)% of the theoritical mass for my cavity (Dimension_x * Dimension_Y). However, after 200000 timesteps, its effect can be clearly seen.
This post seem to indicate that the regularized scheme for BCs does not conserve mass, which, if true, is likely to be at the origin of my observations. Any thoughts or remarks on that matter ?
However I do not quite understand how can Reg. BCs not conserve mass.
Finally, let me describe the scheme I am implementing to apply these BCs, maybe I’m doing it wrong.
1- For every node on the lattice, including wet boundary nodes, collide and stream.
2 - For every wet bdry node - Except corners -
2.1 - Compute rho. Zou & He method : rho = 2.0*(known pops) + tangential pops
[code=“cpp”]
Example : North Wall
rho_ = 2.0*(f[x][y][5]+f[x][y][6]+f[x][y][2]) + (f[x][y][3]+f[x][y][1]+f[x][y][0]);
2.2 - Bounce-back of nonquilibrium parts : for unknown populations, fneq(k) = f(opp(k) - feq(opp(k))
2.3 - Compute the tensor Pi^(1)
[code="cpp"]
Pi_xx = 0.0; Pi_xy = 0.0; Pi_yy = 0.0;
for(int k=0;k<9;k++)
{
Pi_xx += fneq[k]*e[k][0]*e[k][0];
Pi_yy += fneq[k]*e[k][1]*e[k][1];
Pi_xy += fneq[k]*e[k][0]*e[k][1];
}
2.3 - Compute the tensor Q
[code=“cpp”]
for(int k=0;k<9;k++)
{
Q_xx = e[k][0]e[k][0] - 1.0/3.0;
Q_yy = e[k][1]e[k][1] - 1.0/3.0;
Q_xy = e[k][0]e[k][1];
QPi = 4.5w[k](Q_xxPi_xx+Q_yyPi_yy+2.0Q_xy*Pi_xy);
}
2.4 - Compute regularized populations
[code="cpp"]
for(int k=0;k<9;k++)
{
f[x][y][k] = eq[k] + QPi; // eq[k] is the equilibirum pop for direction k
}
3 - Deal with corners.
3.1 - Extrapolate rho from the nearest diagonal bulk nodes.
[code=“cpp”]
Example - Top Left corner. Looking for rho at # by extrapolate from *'s.
+ + +
-
-
+ +
-
-
-
* +
-
-
+ + *
3.2 - Bounce back of nonequilibrium. Same as before for conventional links.
For buried links (directions 5 and 7 for the top left corner), also apply the bonce back of nonequilibrium part and use the fact that the sum over directions of each non equilibrium part must be 0 by construction.
[code="cpp"]
Example - Top Left Corner
fneq[1] = fneq[3];
fneq[4] = fneq[2];
fneq[8] = fneq[6];
fneq[5] = -0.5*(fneq[0]+fneq[1]+fneq[2]+fneq[3]+fneq[4]+fneq[6]+fneq[8]);
fneq[7] = fneq[5];
3.3 - Compute Pi^(1), Q and regularized populations same as before.
4 - Compute macroscopic fields
---------> Back to step 1.
- I read somewhere that it might not be relevant to compare the results with reference papers in a point wise manner. If you know about this, feel free to tell me where I could read more about it, even though it’s not the purpose of this present post.
Useful material
6 2 5
3 0 1
7 4 8
[1] : Ghia, U. K. N. G., Kirti N. Ghia, and C. T. Shin. “High-Re solutions for incompressible flow using the Navier-Stokes equations and a multigrid method.” Journal of computational physics 48.3 (1982): 387-411.