Note, there is already a question here, however the answers don't answer the original question, let alone specific considerations when dealing with nearly in-compressible fluids (water). Another similar question here, but does not speak to water, and really only talks about performance considerations. I'm not too concerned about performance right now.

I've been trying to figure out what method to use for CFD with water. I've tried simulating with LBM, and no matter what I do, I get massive instability with infinities and NaN's everywhere after some time, it seems inevitable (using float32's, as non scientific Nvidia Gpus can have 1/32 double precision throughput). The thing is LBM is extremely easy to implement, and Navier Stokes methods are not.

So what I'm trying to figure out is LBM even the right tool for the job? I see people talk about LBM in terms of compressible fluid, but for water and the scale I'm dealing with this would appear moot, it isn't necessary to deal with. So if I'm dealing with water, and I don't want to deal with compressibility, should I just stick with Navier Stokes? Or is there some way to get LBM to behave for water?

  1. Note, I'm primarily concerned about the force of water on macroscopic objects, and what effect the water has on that single object's movement.
  2. I'll be using a GPU, and I will be implementing everything (and have already done so for LBM on the GPU), my simulations are specific enough that other peoples frameworks will not work (though I also don't need the precision of simulation that they need)
  • $\begingroup$ Hi there. When you say you want to simulate water, what kind of properties do you want to include in your model. You already named compressibility. If you want to study the interaction of water flow, then most likely the effects of compressibility will not be very important (except if you are in some high velocity setting). I would start out with dismissing compressibility! $\endgroup$
    – MPIchael
    Commented Oct 4, 2019 at 9:30
  • $\begingroup$ @MPIchael Yes compressibility is not important ( I mention this in my question ), in fact the only thing I'm actually concerned about is the effect of water on the movement of a moving rigid body object within water. The object will move which will also displace water, which will also have an effect of the movement of said object. These kinds of effects are the ones I'm hoping to model, ie like fish movement through water. $\endgroup$
    – Krupip
    Commented Oct 4, 2019 at 12:09
  • $\begingroup$ I use LBM for almost three years in my PhD research and before that for many years I used FEM and FVM for CFD. In my opinion: LBM is only suitable for creepy flows in structured grid geometries. Otherwise, even in incompressible regime, if you have intermediate Reynolds number ~ 100-200 you will be in trouble... Unless you have a really powerful cluster with many core hours queue reserved specially for you to do simulation. I strongly recommend to use FEM or FVM for any kind of CFD simulations (incompressible, compressible, etc.). $\endgroup$ Commented Oct 4, 2019 at 13:00
  • $\begingroup$ After three years days and nights using LBM for more than 1000 really large-scale CFD simulations (my number of lattice points are 100 Millions to 200 Millions and number of timesteps ~ 3000000, so imagine the nightmare that I have right now just for simulating a flow in intermediate laminar regime of Re ~ 500 to 1000): Don't waste your time with LBM... Yes it's really easier to implement and parallaelize it (even it's not true really, I can refer you to my code which has ~ 2 millions C++ lines parallelized with MPI, which only god knows what's going on there...) but it doesn't worth it... $\endgroup$ Commented Oct 4, 2019 at 13:03
  • $\begingroup$ And about GPU: LBM and GPU is just waste of time... Unless your model is small enough to be able to fit it in a single GPU without any need to have hybrid MPI+GPU architecture. MPI+GPU that we recently implemented in our code with hope of getting better performance: It really sucks... why? Because of really expensive handshaking of communications between GPUs. In fact, recently we realized that communication between GPUs through MPI is about 1000 times slower... So, don't waste your time and go ahead with many powerful FEM or FVM based CFD solvers like OpenFOAM and Nek5000. Believe me bro... $\endgroup$ Commented Oct 4, 2019 at 13:05

3 Answers 3


First, I'm not sure why you emphasize on water? I mean I understand that you are looking for a CFD scheme that works for your special case, but you need to know that water fluid is not a special fluid at all. Water is a incompressible fluid and you can simulate its movement by using incompressible Navier-Stokes equation as long as your Mach number is not high, typically smaller than 0.3. So, by knowing your Mach number and the Reynolds number for your specific case, you can pick incompressible Navier-Stokes equation to simulate your water flow.

You have two options here:

  1. Solve Incompressible Navier-Stokes equation: $$\nabla \cdot \mathbf{u} = 0$$

$$\rho\frac{\partial \mathbf{u}}{\partial t} + \rho \mathbf{u} \cdot \nabla \mathbf{u} = \nabla \cdot \tau + \mathbf{f}$$

For Newtonian fluid: $$\tau = -P\mathrm{I} + \mu(\nabla \otimes \mathbf{u} + (\nabla \otimes \mathbf{u})^{T})$$

Finally: $$\rho \frac{\partial \mathbf{u}}{\partial t} + \rho \mathbf{u} \cdot \nabla \mathbf{u} = -\nabla P + \mu \nabla^{2} \mathbf{u} + \mathbf{f}$$

You can solve it by using FEM or FVM methods, which is successfully implemented and verified in OpenFOAM and Nek5000 packages. I recommend Nek5000 which even has industrial scale accuracy.

  1. Use lattice Boltzmann method (LBM): $$f_{i}(\mathbf{r}+\Delta t \mathbf{c}_{i},t+\Delta t) = f_{i}(\mathbf{r},t) - \frac{f_{i}^{eq}-f_{i}}{\tau}$$

For ensuring conditional stability you need to have relaxation time $\tau > 0.5$ and $\mathrm{Mach} < 0.02$. Note that LBM $\mathrm{Mach}$ is not equal to real Mach number. These are the conditions which may work for your case or for countless cases will not work at all.

If I was in your shoe and I had the knowledge about LBM that I have right now three years ago, I would definitely choose the first option. Any sane person will go for the first option I believe. The most famous argument that we LBM people tell other people including traditional CFD people is: LBM is easy to implement cause it just contains two operations of collision and streaming and in collision step everything could be calculated locally and in streaming step you need to just copy everything in certain directions in comparison to traditional CFD techniques that you have nonlinear momentum acceleration term ($\mathbf{u} \cdot \nabla \mathbf{u}$) which is a pain to discretize and also for the pressure you need to solve time-consuming Poisson pressure equation: $\nabla^{2} P = - (\nabla \otimes \mathbf{u})^{T} : (\nabla \otimes \mathbf{u})$. Also nowadays, we argue that meshing step is more time-consuming in comparison to LBM, cause we call LBM a mesh free method. All of them are true for creepy flows ($Re << 1$) in structured grids like cube, which of course it just limits us to toy model and is not suitable for real world engineering applications. In fact, one iteration of LBM is much much faster than one iteration of FEM or FVM but the problem is that in order to satisfy those stability conditions for LBM you need so many iterations which is at least three order of magnitude higher than conventional FEM or FVM flow solvers. If your adviser is a LBM expert and you are forced to use LBM and of course if you have lots of computational resources, go for it, tough luck! But, if you don't, I recommend to pick a famous open source or commerical CFD solver and more focus on engineering aspects instead of thinking why your simulation become unstable. Already, you see that it is really difficult to get a stable LBM simulation, so just write down the cons pros of these two approaches and take the one that is more inclined to your knowledge, resources, etc.


Let me put some calculation here to show that indeed LBM needs a huge number of meshes and timesteps to work properly at mid range Re flow regime. Let's say we have a water pipe ($\nu = 10^{-6}$ $\frac{\mathrm{m}^{2}}{\mathrm{s}}$) with 1.0 mm diameter and 1 $\frac{\mathrm{m}}{\mathrm{s}}$ velocity in it, which is not really a big deal to solve with FEM or FVM. The Reynolds number is: $Re = \frac{D u}{\nu} = \frac{10^{-3} \times 1}{10^{6}} = 1000$. For LBM stability, a rule of thumb is: $Mach < 0.02$ and $\tau = 0.54$ for Mach number and relaxation time respectively. So:

$$\tau - \frac{1}{2} = \frac{3 \nu \Delta t}{\Delta x^{2}}$$

$$Mach = \frac{\sqrt{3} u \Delta t}{\Delta x}$$

If $Mach < 0.02$:

$$\frac{\Delta t}{\Delta x} < \frac{0.02}{\sqrt{3} u} = 0.012$$


$$0.04 = \frac{3 \times 10^{-6} \Delta t}{\Delta x^{2}}$$


$$\frac{\Delta t}{\Delta x^{2}} = 13333.33$$


$$\frac{\Delta t}{\Delta x} = 13333.33 \Delta x$$


$$13333.33 \Delta x < 0.012$$

$$\Delta x < 10^{-6}$$

I pick the upper bound as $\Delta x = 10^{-6}$ m. So:

$$\Delta t = 1.3 \times 10^{-8}$$

First, look at the insane timestep that we calculated. It means if you want to simulate water flow in a pipe for 1 second in reality it needs almost 100 millions timesteps! Wow!

Let's say your 1 mm diameter pipe has 10 mm or 1 cm length, which is not that uncommon. So, the total volume of your pipe is:

$$V = \frac{\pi}{4}D^{2} L = \frac{\pi}{4} (10^{-3})^{2} \times 10 \times 10^{-3} = 8 \times 10^{-9}$$

Remember LBM works with voxels and it means your mesh is just a tiny cube and has same size everywhere in your domain (ridiculous right?! forget about adaptive mesh refinement in LBM...). So, each voxel volume is:

$$V_{voxel} = \Delta x^{3} = 10^{-18}$$

Finally, you would find the number of voxels that you need to fill this 1 mm diameter and 10 mm length pipe as:

$$N = \frac{V}{V_{voxel}} = \frac{8 \times 10^{-9}}{10^{-18}} = 8 \times 10^{9}$$

Another insane finding... You need 8 billions voxels to do this simulation and still you can't be 100% sure that your simulation will remain stable or not... Come on... FEM or FVM people do this simulation as a piece of cake less than a couple of seconds on their crappy old laptops by using ANSYS Fluent or any other open source code if you want to argue that ANSYS is commercial. Want to say something about MRT? That's garbage... It's just good for publishing nonsense articles in simple 2D rectangular geometries to show that they were able to simulate lid driven cavity problem with higher Reynolds number. Just look at the CFD market and count besides XFlow and Exa who produce LBM software for serious industrial usage? For your information Exa is acquired by Dassault Systemes and I'm pretty sure they just bought Exa to close it.

  • $\begingroup$ This in combination with everything else I've read, tells me that LBM is a no go for the kind of simulation I'm doing. I think I'm going to use PIC/FLIP methods, as I need to be fast (not necessarily 100% accurate), though all my simulations are submerged. Also those resources you've linked to will be a good resource for me to verify the accuracy (or lack of) of my solutions! $\endgroup$
    – Krupip
    Commented Oct 4, 2019 at 20:03
  • $\begingroup$ Since you mention the very heart of the LBM (collision & streaming) I think it would be helpful if you present the universal streaming step and a general form of the collision step - i.e., not limited to the LBGK collision operator. $\endgroup$
    – Dan Doe
    Commented Jan 6, 2022 at 7:19

I disagree with the answer given by Alone Programmer. His reasoning is not objective and seems to be based on a very primitive lattice-Boltzmann model with BGK collision operator. LBM has significantly matured over the years and is in my opinion a very attractive numerical solver. Nonetheless it depends on the application if it really fits your needs. Irrespectively of the numerical method chosen you have to think about whether you really want to maintain your own code base or if it is not better to use and extend existing software tools.

I have been programming LBM multi-threaded OpenMP algorithms now for 2 years (update 2022: now over 4 years off and on) and I am still convinced that the method or a modification of it will be a very attractive tool for transient incompressible flows and in the future also for compressible simulations (update 2022: this has proven true, a lot of LBM start-ups and companies have emerged in the meanwhile and XFlow has proven to be very competitive, e.g. see their simulation of the Bugatti Chiron). The following thoughts consider positive and negative aspects of LBM:

  • LBM is based on the kinetic theory of gases. In the kinetic theory of gases one creates fictive dilute gas models of discrete particles or particle clusters interacting with each other and observes their dynamics. Boltzmann introduced a probability distribution over phase space, an extended concept of the traditional density that is only strictly valid for dilute gases. Surprisingly the continuous Boltzmann equation even with simplified collision models such as the BGK approximation preserves the complete set of compressible Navier-Stokes equations (Chapman-Enskog expansion). This means the theoretical framework can be used to make predictions about compressible and incompressible flows. One can make use of the law of similarity, meaning it does not even matter if you simulate a liquid or a fluid: As long as the assumptions about the material parameters such as density and the dimensionless parameters of interest in the differential equations such as the Reynolds and Mach-number are identical, the behaviour is similar. Traditional LBM methods apply a discretisation of the space with distinct single-speed (they only connect neighbouring cells) velocity subsets- so called lattices. These have a certain number of speeds and thus a number of "symmetries". These symmetries are connected to tensors and these again to conserved quantities and conservation equations. The common discretisations have too little speeds in order to preserve the full Navier-Stokes equations and you are left with the conservation equations for incompressible flow with an error term that rises with the model Mach number. The model is therefore used to preserve the incompressible conservation equations asymptotically. You can simulate substances in the limit of a low-Mach number incompressible approximation wherever that is an incompressible fluid like water or the flow around a car. That being said there is room for improvement: There has been continuous effort in particular from the research group of Prof. Ilya Karlin in Zürich to make modifications to the basic idea and apply it to (highly) compressible flows. The more basic nature of LBM - namely the distribution functions- as supposed to the direct application of the continuum hypothesis and the Navier-Stokes equations - might allow for a simulation of highly compressible fluids with shock waves where the continuum hypothesis clearly breaks down. If such an algorithm will really be faster than finite-volume-based methods has still to be seen but I have the confidence it can be.

  • Incompressible LBM struggle with medium Reynolds numbers already when used with the standard formulation with BGK operator and basic boundaries such as the 2D Zou/He boundaries. But there are significantly more stable formulations such as MRT, TRT, cascaded and entropic collision operators and furthermore several modifications of the Smagorinsky turbulence model and wall functions that can be implemented efficiently. These formulations also require for more sophisticated boundary conditions and strategies like block-structures mesh refinement are available. This means while the basic algorithm is simple incorporating a stable model is significantly more difficult, doing so effectively even more so. Nonetheless the limitations discussed by Alone Programmer are not there. You can simulate highly turbulent flows even on desktop setups. The limitations arise from the nature of the models (DNS, LES, DES vs RANS) and not from LBM as a CFD solver!

  • LBM traditionally works on a regular structured grid. This means meshing is easy but boundaries are stair-cased if one uses the traditional bounce-back algorithm. Implementing more precise boundaries like interpolated bounce-back violates though exact mass-conservation.

  • In LBM you have limited options for turbulent inlet conditions. You have to use experimentally determined velocity profiles, trip wires to trigger turbulence, pressure-periodic boundary conditions or a periodic box.

  • The LBM is comparably simple and predictable. Thus the algorithm is very fast and its parallel scalability is excellent, unrivaled by finite-volume based methods if you apply it to transient problems. You can find implementations in the literature on CUDA GPUs and Xeon Phi coprocessors that update up to 1 billion lattice nodes a second. On a 12-core processor I have archieved 150 million lattice node updates with a D3Q19 lattice. (Update 2022: With increasing processing power of CPUs and GPUs modern implementations are able to update between 300 million cells on a modern CPU and 3.5 billion cells on a modern accelerator card, see e.g. the benchmark of the University of Geneva) Do not compare apples to oranges: lattice-Boltzmann is inherently transient and requires a lot of times steps you can control the accuracy and temporal resolution of the algorithm by setting the macroscopic velocity and therefore a particular model Mach number. The model Mach number should be lower than $0.3$ across the majority of the flow field to obtain usable results (see e.g. Aslan et al.). You can obtain faster results with finite-volume method if you calculate steady-state processes and you can use RANS models but let's face it: Eddy-viscosity RANS models such as $k$-$\epsilon$ or $k$-$\omega$, have severe drawbacks. They are inherently isotropic and thus struggle with anisotropy such as in swirling flows. $k$-$\epsilon$ further overestimates the turbulent intensity rendering it pretty much useless for most tasks. On the other hand you have very expensive RSM models that have 6 distinct transport equations for the turbulent stresses. They handle anisotropy way better but they have stability issues and so only few modern publications on this topic can be found. URANS formulations altogether offer very limited insight into transient turbulent processes. DNS is almost impossible for technical applications and LES is still restricted to large-scale clusters. DES + LES is computationally very expensive but delivers good results. Now where does LBM go? It is a solver like FVM and FDM. The traditional formulation is basically DNS but there exist RANS and as mentioned LES formulations. Without a doubt LBM is way faster than FVM when performing DNS or LES but the question is would a RANS model suffice? There is though a reason truck aerodynamics where separation prevails are dominated by LBM: If one needs a transient resolution of turbulent structures like vorticity LBM is clearly to be preferred! Just look at what XFlow can do and look at how much of that traditional commercial software such as Ansys Fluent simply can't because of the different modelling paradigm.

  • Nonetheless LBM is currently actively researched, it lacks a wide acceptance and readily available code. There are only few optimised options available such as the open-source projects Palabos, OpenLB and the commercial softwares XFlow and ProLB that allow somebody to use LBM without programming it by yourself from scratch. LBM codes are often used in research, custom made, tailored to a particular platform by skilled personnel. *(update 2022: This has gotten better, in particular Palabos and OpenLB have worked on their GPU support! XFlow and SimScale have gained traction.).

So in my opinion the reason why LBM is becoming increasingly more popular is because of its capabilities to resolve transient flows accurately and effectively. Writing a well performing LBM algorithm from scratch requires though time and I would not advise anybody to do so for their PhD thesis. Most supervisors will expect you to model physics and likely won't give you the time you would have to put into a solver itself. If you already have a well-performing framework then it is an excellent choice for transient flows in particular in complex geometries. But you should ask yourself: Does a cheaper model suffice? If you do not need the transient information of LBM simulations then every time step is a burden! If the smeared time-averaged results of cheap RANS models are sufficient turn to OpenFOAM.

  • $\begingroup$ Unfortunately, I should disagree. You pretty much described LBM and its various applications/implementations, etc., which are not necessarily interesting to OP or even answering this question. The main point in your answer is that you provided inherent transient nature of LBM as some sort of a unique feature, which I think you just exaggerating it... There are a couple of well defined, validated, verified, transient CFD frameworks based on FEM or FVM, so LBM does not add anything here or does not have any unique feature. $\endgroup$ Commented Oct 30, 2019 at 1:50
  • $\begingroup$ You talked about more complicated collision operators such as MRT, TRT, etc. which mathematicians can understand what's going on with them, but in my experience as someone use all of these collision operators not just simple LBGK, in my opinion these are just nonsense for engineering applications. If you are a mathematicians or a physicist that study a fundamental problem of fluid dynamics, probably LBM with its all LBGK, MRT, TRT, filtering, BFL boundary conditions, etc might be a good option for you, but for engineering application, never LBM... $\endgroup$ Commented Oct 30, 2019 at 1:53
  • $\begingroup$ I agree about kinetic theory, Chapman-Enskog expansion, etc that might give a better understanding of fluid dynamics fundamental aspects or probably from a non-equilibrium thermodynamics point of view, but I can't understand why an engineer should use LBM when there are accurate, fast, scalable, multiphysics, parallel implementation of CFD based on FEM or FVM such as OpenFOAM or Nek5000. This debate will never get stopped cause some people just love LBM... But my opinion is based on my three years of experience working days and nights with LBM. Nothing is in LBM right now that I did not try... $\endgroup$ Commented Oct 30, 2019 at 1:58
  • $\begingroup$ So again: if your flow is creepy ($Re << 1$) and steady state and you don't want to spend money to buy license of Fluent or CFX, the best option is to pick an open source LBM implementation such as OpenLB or Palabos. Other than that using LBM is just torturing yourself... I have a friend that her research is about creepy flows in porous media and she's pretty OK with using LBM cause she could fit her simulations even in her desktop computer without the need for any cluster, but for me to deal with 100 to 200 millions lattice points with 3000000 timesteps, it's just an endless nightmare... $\endgroup$ Commented Oct 30, 2019 at 2:13
  • $\begingroup$ @Alone Programmer I have listed advantages and disadvantages of the method and left it to the questioner to decide. The transient nature and the small required time steps are unique features of LBM which help it for transient applications. You are comparing a framework to a method. LBM as a method is excellent for transient applications, just the easily accessible frameworks are limited. And for certain applications there are no alternatives. E.g. find me a solver that is able to mesh complex geometries such as realistic porous media. Nobody will use FVM for that for sure. $\endgroup$
    – 2b-t
    Commented Oct 30, 2019 at 8:20

Just happened to have a look at this old post and happened to notice that the given calculation in the accepted answer is scientifically completely (!) wrong such that I had to outline it in this post:

  • The given computational burden is wrong by four orders of magnitude: A calculation which actually takes less than 2 hours and a half on a modern GPU would take close to three years and close to 40 years on a modern CPU with his assumptions! Cringe!
  • The answer does not differentiate between a model, a numerical method and an implementation giving a skewed view and rendering the given answer virtually worthless.

Even though the author claims to have extended experience with the method, I have the impression that his/her frustration (which is visible throughout the answer resulting in a highly subjective answer) is caused by a lack of profound knowledge of the theoretical foundations of the method as well as its efficient implementation. I will outline this in the post below.

  • The post lacks a clear differentiation between numerical methods and implementations: Finite element method (FEM), finite volume method (FVM) and the lattice-Boltzmann method (LBM) are methods and should be clearly separated from their implementations which are software tools/libraries such as OpenFOAM, Nek5000 or Ansys Fluent. There are also open-source (e.g. Palabos, OpenLB) and commercial solvers (e.g. Simulia XFlow by Dassault Systemes or Pacefish used by SimScale) based on the lattice-Boltzmann method which offer also advanced features like grid-refinement, large-eddy and RANS turbulence models. It is sufficient to have a look what modern LBM codes with large-eddy turbulence models can do, such as this aerodynamics simulation of the flow around a Bugatti by XFlow, to see that the comparison is completely lacking context: You can't compare your own buggy, in-house code to a commercially developed and fully tested professional framework in order to decide which method is more suitable than the other. This is in particular true for LBM where a main source of instability are low-order boundary conditions that fail to model higher-order contributions that still have an effect on $f^{(neq)}$ and therefore on the evolution of the Boltzmann equation.

  • Finite element methods are not popular in computational fluid dynamics due to conservativity and are generally significantly slower. Neither Nek5000 (finite volume, finite differences, spectral elements) nor OpenFOAM (finite volume) are based on finite element.

  • The answer mixes methods (see point above) and levels of modelling and turbulence models (such as Direct Numerical Simulation (DNS), large-eddy turbulence models and (unsteady) Reynolds-averaging (RANS)). Directly solving the Navier-Stokes equations (DNS) for a transient flow is quite challenging in particular with traditional numerical methods such as finite element methods. A commercial solution based on FVM will only be able to solve a pipe flow at Reynolds $1000$ reasonably fast when assuming steady flow (so not transient, which is reasonably for the given Reynolds number but would not be reasonable for higher Reynolds numbers as the flow will transition to turbulent flow for $Re_D \geq 2300$) and deploying an unstructured grid. Even then lattice-Boltzmann will be about as fast when calculating transiently, in particular if taking into consideration grid generation. Said scenario is actually the worst case for a comparison and not representative: LBM can't accelerate towards steady-state as it is inherently transient - even though the solution is steady you have to perform a transient LBM simulation.

  • The local Mach number required for stability given by $0.02$ is a magnitude smaller than you would actually need (see e.g. Aslan et al. - "Investigation of the Lattice Boltzmann SRT and MRT Stability for Lid Driven Cavity Flow") making the required number of time steps a magnitude (and more) larger.

  • The Mach number in lattice-Boltzmann is defined as

    $$Ma := \frac{u_{lb}^{max}}{c_s}$$

    where $c_s$ is the lattice speed of sound which for the most common single-speed lattices is given by $c_s = \frac{1}{\sqrt{3}} \frac{\Delta x}{\Delta t}$. But $\Delta x$ and $\Delta t$ are lattice units, the dimensionless system in LBM which chooses $\Delta x = 1$ and $\Delta t = 1$ (You are never ever calculating any $\Delta x$ or $\Delta t$ with the standard incompressible LBM! Just have a look at Latt or Krueger et al.! You are working on the basis of a dimensionless system with convenient dimensionless units that is related to the dimensioned system through the law of similarity. What is then done in the derivation above just hurts if you have seen any introductory video to the lattice-Boltzmann equation.). This means your highest possible over-velocity in the fluid domain is

    $$|u_{lb}^{max}| = c_s \, Ma \approx \frac{1}{3} \frac{\Delta x}{\Delta t} \, 0.3 = 0.1$$

    If we assume that there are no significant contractions inside the domain we can assume that the for a pipe flow the characteristic velocity can be chosen to the maximal admittable velocity $U_{lb} = |u_{lb}^{max}|$.

    Even keeping the very conservative stability limiy $\tau \approx 0.54$ as a stability limit (which is only required for simple collision operators like BGK without any turbulence models, while other collision operators are closer to $0.51$ and turbulences models allow arbitrary close relaxation times to $0.5$ that are only limited by the floating-point accuracy so something like $10^{-15}$ for the double data-type) this gives with the formula for the model viscosity

    $$ \nu = \left( \frac{\tau}{\Delta t} - \frac{1}{2} \right) c_s^2 \Delta t $$

    a limiting model viscosity $\nu_{lb} = 0.01\overline{3}$ which means that with you have to choose your characteristic length according to the Reynolds number $Re_D = \frac{U_{lb} \, D_{lb}}{\nu_{lb}}$ resulting in $D_{lb} = \frac{Re \, \nu_{lb}}{U_{lb}}$ which in our case gives $$D_{lb} = \frac{Re \, \nu_{lb}}{U_{lb}} = \frac{1000 \cdot 0.01\overline{3}}{0.1} = 133.\overline{3}$$ for $D_{lb}$ which is the diameter of the pipe $D$ in lattice units. Assuming a rectangular computational domain that is $C$ times longer than its diameter this gives a number of cells $$N = NX \cdot NY \cdot NZ \approx C L_{lb}^3$$ which for a mentioned $C = 10$ gives $2.37 \cdot 10^7$ nodes, so 23 million nodes which are ehhhhm, a bit less than the mentioned 8 billion. And as mentioned these stability limits are not relevant for more advanced collision operators, turbulence models and in any case resolution requirements might be reduced with grid-refinement approaches.

  • Same goes for the time step: For unit conversion we can introduce the conversion factor between the dimensionless $lb$ system and the dimensioned system $d$ in SI-units: $c_L = \frac{L_d}{L_{lb}}$ for length, $c_U = \frac{U_d}{U_{lb}}$ for velocity resulting in a conversion factor for time $c_T = \frac{c_L}{c_U}$ and therefore in $T_{lb} = \frac{T_d}{c_T} = \frac{c_U \, T_d}{c_L}$. For the given numbers $c_L = \frac{10^{-3}}{133} \approx 7.5 \cdot 10^{-6}$, $c_U = \frac{1}{0.1} = 10$ and therefore $T_{lb} \approx \frac{10 \cdot 1}{7.5 \cdot 10^{-6}} = 1.\overline{3} \cdot 10^6$. So roughly 1.3 million time steps are required for a stable simulation of 1 second for the given system, sounds like a lot but it is not that much: A well-programmed three-dimensional lattice-Boltzmann implementations can update around 200-300 million nodes within one second on a single CPU (see e.g. Wittmann et al. or Bauer et al.), implementations on the old Xeon Phis can reach around 1 billion node updates per second (see Robertsén et al.) and on a graphics/accelerator card well over 3.5 billion nodes per second (see Latt et al.). This means even for what seems to be an unreasonable large number of time steps the method will be very fast!

  • Therefore the simulation of one second real-time with a primitive BGK collision operator for the given Reynolds number of $1000$ resulting in a grid with $23$ million computational nodes and $10$ million time steps on a single GPU will remain stable and complete in approximately $\frac{23 \cdot 10^6 \cdot 1.3 \cdot 10^6}{3500 \cdot 10^6} = 8542.\overline{857142}$ seconds, so around 2 hours and 20 minutes which is reasonable for a transient simulation.

  • If you are though implement a more stable collision operator such as multiple-relaxation times (MRT, see d'Humieres) or cumulant (see Geier), or if you add a turbulence model such as the Smagorinsky large-eddy model (see e.g. Sagaut) your stability requires significantly less computational nodes (something like $67$ or $33$ cells or even lower - likely anything less than $67$ is too little for resolving the physics with a uniform grid, you would have to calculate your dimensionless wall distance $y+$ and see, probably one would reduce the simulation domain to $5 \, D$ in favour of having a larger $D_{lb}$ anyways). Assuming a simulation domain of $67 \cdot 67 \cdot 667$ we would have $8$ times less cells, with just under $3$ million nodes, $c_L$ would drop to the half $1.49 \cdot 10^{-5}$ and so would the required time steps $666666.\overline{6}$. Our computational burden decreases by a factor of $16$ resulting in $571$ seconds or less than 10 minutes simulation for 1 second of real-time (on a single A100 GPU or less than 2 hours on a powerful desktop). This means the LBM simulation on a GPU will voxelise the geometry (sub-seconds) and complete the simulation before you have even created your computational mesh for a finite volume simulation.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.