The fastest and most memory efficient lattice Boltzmann CFD software, running on any GPU via OpenCL.

• CFD model: lattice Boltzmann method (LBM)
  • streaming (part 2/2):

    f₀^temp(x,t) = f₀(x, t)
    f_i^temp(x,t) = f_{(t%2 ? i : (i%2 ? i+1 : i-1))}(i%2 ? x : x-e_i, t)   for   i ∈ [1, q-1]

  • collision:

    ρ(x,t) = (Σ_i f_i^temp(x,t)) + 1
    
    u(x,t) = ¹∕_ρ(x,t) Σ_i c_i f_i^temp(x,t)
    
    f_i^eq-shifted(x,t) = w_i ρ · (^(u_°c_i)2∕_(2c⁴) - ^(u_°u)∕_(2c²) + ^(u_°c_i)∕_c²) + w_i (ρ-1)
    
    f_i^temp(x, t+Δt) = f_i^temp(x,t) + Ω_i(f_i^temp(x,t), f_i^eq-shifted(x,t), τ)

  • streaming (part 1/2):

    f₀(x, t+Δt) = f₀^temp(x, t+Δt)
    f_{(t%2 ? (i%2 ? i+1 : i-1) : i)}(i%2 ? x+e_i : x, t+Δt) = f_i^temp(x, t+Δt)   for   i ∈ [1, q-1]

• peak performance on most GPUs (datacenter/gaming/professional/laptop), validated with roofline model
• up to 4.29 billion (2³²) lattice points, or 1624³ resolution, on a single GPU (if it has 225 GB memory)
• optimized to minimize memory demand to 55 Bytes/node (~⅙ (~⅓) of conventional FP64 (FP32) LBM solvers)
  • in-place streaming with Esoteric-Pull: almost cuts memory demand in half and slightly increases performance due to implicit bounce-back boundaries; offers optimal memory access patterns for single-node in-place streaming
  • decoupled arithmetic precision (FP32) and memory precision (FP32 or FP16S or FP16C): all arithmetic is done in FP32 for compatibility on all hardware, but LBM density distribution functions in memory can be compressed to FP16S or FP16C: almost cuts memory demand in half again and almost doubles performance, without impacting overall accuracy for most setups
• DDF-shifting and other algebraic optimization to minimize round-off error
• updating density and velocity in the stream_collide() kernel is optional (higher performance if disabled)
• velocity sets:
  • D2Q9
  • D3Q15
  • D3Q19 (default)
  • D3Q27
• collision operators
  • single-relaxation-time (SRT/BGK) (default)
  • two-relaxation-time (TRT)
• only 8 flag bits per lattice point (can be used independently / at the same time):
  • TYPE_S (stationary or moving) solid boundaries
  • TYPE_E equilibrium boundaries (inflow/outflow)
  • TYPE_T temperature boundaries
  • TYPE_F free surface (fluid)
  • TYPE_I free surface (interface)
  • TYPE_G free surface (gas)
  • TYPE_X remaining for custom use or further extensions
  • TYPE_Y remaining for custom use or further extensions

• boundary types
  • stationary mid-grid bounce-back boundaries (stationary solid boundaries)
  • moving mid-grid bounce-back boundaries (moving solid boundaries)
  • equilibrium boundaries (non-reflective inflow/outflow)
  • temperature boundaries (fixed temperature)
• global force per volume (Guo forcing), can be modified on-the-fly
• local force per volume (force field)
  • optional computation of forces from the fluid on solid boundaries
• state-of-the-art free surface LBM (FSLBM) implementation:
  • volume-of-fluid model
  • fully analytic PLIC for efficient curvature calculation
  • improved mass conservation
  • ultra efficient implementation with only 4 kernels additionally to stream_collide() kernel
• thermal LBM to simulate thermal convection
  • D3Q7 subgrid for thermal DDFs
  • in-place streaming with Esoteric-Pull for thermal DDFs
  • optional FP16S or FP16C compression for thermal DDFs with DDF-shifting
• Smagorinsky-Lilly subgrid turbulence LES model to keep simulations with very large Reynolds number stable
  equations

  Π_αβ = Σ_i e_iα e_iβ (f_i - f_i^eq-shifted)
  
  Q = Σ_αβ Π_αβ²
                       ______________________
  τ = ½ (τ₀ + √ τ₀² + ^(16√2)∕_(3π²) ^√Q∕_ρ )

• on Windows: real time interactive rasterization and raytracing graphics
• on Windows and Linux (even in WSL and/or in a remote console through SSH): real time interactive ASCII console graphics
• with interactive graphics mode disabled, image resolution can be as large as VRAM allows for (132 Megapixel (16K) and above)
• (interacitive) visualization modes:
  • flags (and force vectors on solid boundary nodes if the extension is used)
  • velocity field
  • streamlines
  • velocity-colored Q-criterion isosurface
  • rasterized free surface with marching-cubes
  • raytraced free surface with fast ray-grid traversal and marching-cubes, either 1-4 rays/pixel or 1-10 rays/pixel

1. Check the settings and extensions in src/defines.hpp by uncommenting corresponding lines.

2. Write a C++ setup skript as main_setup() function in src/setup.cpp (get inspiration from existing setups):

   • For unit conversion, use the units struct.
   • For initializing the box, use call LBM lbm(...); constructor.
   • Set the initial condition in a loop that iterates over the entire lattice by writing to lbm.rho[n]/lbm.u.x[n]/lbm.u.y[n]/lbm.u.z[n]/lbm.flags[n].
   • Call lbm.run() to initialize and execute the setup (infinite time steps) or lbm.run(time_steps) to execute only a specific number of time steps.
   • As long as the lbm object is in scope, you can access the memory. As soon as it goes out of scope, all memory associated to the current simulation is freed again.

3. When done with the setup, hit compile+run on Windows in Visual Studio Community or execute ./make.sh on Linux; this will automatically select the fastest installed GPU. Alternatively, you can add the device ID as the first command-line argument, for example ./make.sh 2 to compile+run on device 2, or bin/FluidX3D.exe 2 to run the executable on device 2. Compile time for the entire code is about 10 seconds.

4. On Windows with WINDOWS_GRAPHICS/CONSOLE_GRAPHICS or on Linux with CONSOLE_GRAPHICS, press the P key to start/pause the simulation.

   Toggle rendering modes with the keyboard:

   • 1: flags (and force vectors on solid boundary nodes if the extension is used)
   • 2: velocity field
   • 3: streamlines
   • 4: vorticity / velocity-colored Q-criterion isosurface
   • 5: rasterized free surface
   • 6: raytraced free surface

   Camera movement:

   • Mouse or I/J/K/L: rotate camera
   • Scrollwheel or ,/.: zoom (centered camera mode) or adjust camera movement speed (free camera mode)
   • U: toggle rotation with Mouse and angle snap rotation with I/J/K/L
   • Y/X: adjust camera field of view
   • R: toggle camera autorotation
   • F: toggle centered/free camera mode
   • W/A/S/D/Space/C: move free camera
   • V: toggle stereoscopic rendering for VR
   • B: toggle VR-goggles/3D-TV mode for stereoscopic rendering
   • N/M: adjust eye distance for stereoscopic rendering
   • H: print current camera position/rotation in console as copy/paste command
   • Esc/Alt+F4: quit

• works in Windows, Linux and Android with C++17
• runs on any hardware that supports OpenCL 1.2, from any vendor (Nvidia, AMD, Intel, ...):
  • world's fastest datacenter GPUs like H100, A100, MI250(X), MI210, MI100, V100(S), P100, ...
  • gaming GPUs (desktop or laptop)
  • "professional"/workstation GPUs
  • integrated GPUs
  • Xeon Phi
  • CPUs
  • even smartphone ARM GPUs
• supports importing and voxelizing triangle meshes from binary .stl files
• supports exporting volumetric data as binary .vtk files
• supports exporting rendered frames as .png/.qoi/.bmp files; time-consuming image encoding is handled in parallel on the CPU while the simulation on GPU can continue without delay

Here are performance benchmarks on various hardware in MLUPs/s, or how many million lattice points are updated per second. The settings used for the benchmark are D3Q19 SRT with no extensions enabled (only LBM with implicit mid-grid bounce-back boundaries) and the setup consists of an empty cubic box with sufficient size (typically 256³). Without extensions, a single lattice point requires:

• a memory capacity of 93 (FP32/FP32) or 55 (FP32/FP16) Bytes
• a memory bandwidth of 153 (FP32/FP32) or 77 (FP32/FP16) Bytes per time step
• 326 (FP32/FP32) or 364 (FP32/FP16S) or 1241 (FP32/FP16C) FLOPs per time step (FP32+INT32 operations counted combined)

In consequence, the arithmetic intensity of this implementation is 2.13 (FP32/FP32) or 4.73 (FP32/FP16S) or 16.12 (FP32/FP16C) FLOPs/Byte. So performance is only limited by memory bandwidth.

• What physical model does FluidX3D use?

  
  FluidX3D implements the lattice Boltzmann method, a type of direct numerical simulation (DNS), the most accurate type of fluid simulation, but also the most computationally challenging. Optional extension models include volume force (Guo forcing), free surface (volume-of-fluid and PLIC), a temperature model and Smagorinsky-Lilly subgrid turbulence model.
  
  
• FluidX3D only uses FP32 or even FP32/FP16, in contrast to FP64. Are simulation results physically accurate?

  
  Yes, in all but extreme edge cases. The code has been specially optimized to minimize arithmetic round-off errors and make the most out of lower precision. With these optimizations, accuracy in most cases is indistinguishable from FP64 double-precision, even with FP32/FP16 mixed-precision. Details can be found in this paper.
  
  
• Why is the simulation box size limited to 2³² grid points?

  
  The 32-bit unsigned integer grid index will overflow above this number. Using 64-bit index calculation would slow the simulation down by ~20%, as 64-bit uint is calculated on special function units and not the regular GPU cores. 2³² grid points with FP32/FP16 mixed-precision is equivalent to 225GB memory and single GPUs currently are only at 80GB, so it should be fine for a while to come.
  
  
• Comparted to the benchmark numbers stated here, efficiency seems much lower but performance is slightly better for most devices. How can this be?

  
  In that paper, the One-Step-Pull swap algorithm is implemented, using only misaligned reads and coalesced writes. On almost all GPUs, the performance penalty for misaligned writes is much larger than for misaligned reads, and sometimes there is almost no penalty for misaligned reads at all. Because of this, One-Step-Pull runs at peak bandwidth and thus peak efficiency.
  Here, a different swap algorithm termed Esoteric-Pull is used, a type of in-place streaming. This makes the LBM require much less memory (93 vs. 169 (FP32/FP32) or 55 vs. 93 (FP32/FP16) Bytes/node for D3Q19), and also less memory bandwidth (153 vs. 171 (FP32/FP32) or 77 vs. 95 (FP32/FP16) Bytes/node per time step for D3Q19) due to so-called implicit bounce-back boundaries. However memory access now is half coalesced and half misaligned for both reads and writes, so memory access efficiency is lower. For overall performance, these two effects approximately cancel out. The benefit of Esoteric-Pull - being able to simulate domains twice as large with the same amount of memory - clearly outweights the cost of slightly lower memory access efficiency, especially since performance is not reduced overall.
  
  
• Why don't you use CUDA? Wouldn't that be more efficient?

  
  No, that is a wrong myth. OpenCL is exactly as efficient as CUDA on Nvidia GPUs if optimized properly. Here I did roofline model and analyzed OpenCL performance on various hardware. OpenCL efficiency on modern Nvidia GPUs can be 100% with the right memory access pattern, so CUDA can't possibly be any more efficient. Without any performance advantage, there is no reason to use proprietary CUDA over OpenCL, since OpenCL is compatible with a lot more hardware.
  
  
• Why no multi-relaxation-time (MRT) collision operator?

  
  The idea of MRT is to linearly transform the DDFs into "moment space" by matrix multiplication and relax these moments individually, promising better stability and accuracy. In practice, in the vast majority of cases, it has zero or even negative effects on stability and accuracy, and simple SRT is much superior. Apart from the kinematic shear viscosity and conserved terms, the remaining moments are non-physical quantities and their tuning is a blackbox. Although MRT can be implemented in an efficient manner with only a single matrix-vector multiplication in registers, leading to identical performance compared to SRT by remaining bandwidth-bound, storing the matrices vastly elongates and over-complicates the code for no real benefit.
  
  
• How about multi-GPU?

  
  Multi-GPU would allow for larger simulation domains, but is currently not yet supported in FluidX3D. It makes code development and maintaining much more effort and many extensions close to impossible, for example ray-tracing graphics. With the in-place streaming implementation and FP16 memory compression combined, FluidX3D already only requires 1/3 the memory of existing LBM GPU implementations, so large domain size is already possible on a single GPU. Nevertheless, I will look into multi-GPU support in the future, so stay tuned.
  
  

• I'm on a budget and have only a cheap computer. Can I run FluidX3D on my toaster PC/laptop?

  
  Absolutely. Today even the most inexpensive hardware, like integrated GPUs or entry-level gaming GPUs, support OpenCL. You might be a bit more limited on memory capacity and grid resolution, but you should be good to go. I've tested FluidX3D on very old and inexpensive hardware and even on my Samsung S9+ smartphone, and it runs just fine, although admittedly a bit slower.
  
  
• I don't have an expensive workstation GPU, but only a gaming GPU. Will performance suffer?

  
  No. Efficiency on gaming GPUs is exactly as good as on their "professional"/workstation counterparts. Performance often is even better as gaming GPUs have higher boost clocks.
  
  
• Do I need a GPU with ECC memory?

  
  No. Gaming GPUs work just fine. Some Nvidia GPUs automatically reduce memory clocks for compute applications to almost entirely eliminate memory errors.
  
  
• My GPU does not support CUDA. Can I still use FluidX3D?

  
  Yes. FluidX3D uses OpenCL 1.2 and not CUDA, so it runs on any GPU from any vendor since around 2012.
  
  
• I don't have a dedicated graphics card at all. Can I still run FluidX3D on my PC/laptop?

  
  Yes. FluidX3D also runs on all integrated GPUs since around 2012, and also on CPUs.
  
  
• I need more memory than my GPU can offer. Can I run FluidX3D on my CPU as well?

  
  Yes. You only need to install the OpenCL Runtime for Intel CPUs.
  
  
• In the benchmarks you list some very expensive hardware. How do you get access to that?

  
  I'm a scientist (PhD candidate in computational physics) and I use FluidX3D for my research, so I have access to BZHPC, SuperMUC-NG and JURECA-DC supercomputers.
  
  

• I don't have an RTX/DXR GPU that supports raytracing. Can I still use raytracing graphics in FluidX3D?

  
  Yes, and at full performance. FluidX3D does not use a bounding volume hierarchy (BVH) to accelerate raytracing, but fast ray-grid traversal instead, implemented directly in OpenCL C. This is much faster than BVH for moving isosurfaces in the LBM grid (~N vs. ~N²+log(N) runtime; LBM itself is ~N³), and it does not require any dedicated raytracing hardware. Raytracing in FluidX3D runs on any GPU that supports OpenCL 1.2.
  
  
• I have a datacenter/mining GPU without any video output or graphics hardware. Can FluidX3D still render simulation results?

  
  Yes. FluidX3D does all rendering (rasterization and raytracing) in OpenCL C, so no display output and no graphics features like OpenGL/Vulkan/DirectX are required. Rendering is just another form of compute after all. Rendered frames are passed to the CPU over PCIe and then the CPU can either draw them on screen through dedicated/integrated graphics or write them to the hard drive.
  
  
• I'm running FluidX3D on a remote (super-)computer and only have an SSH terminal. Can I still use graphics somehow?

  
  Yes, either directly as interactive ASCII graphics in the terminal or by storing rendered frames on the hard drive and then copying them over via `scp -r [email protected]:"~/path/to/images/folder" .`.
  
  

• I want to learn about programming/software/physics/engineering. Can I use FluidX3D for free?

  
  Yes. Anyone can use FluidX3D for free for public research, education or personal use. Use by scientists, students and hobbyists is free of charge and well encouraged.
  
  
• I am a scientist/teacher with a paid position at a public institution. Can I use FluidX3D for my research/teaching?

  
  Yes, you can use FluidX3D free of charge. This is considered research/education, not commercial use. To give credit, the references listed below should be cited. If you publish data/results generated by altered source versions, the altered source code must be published as well.
  
  
• I work at a company in CFD/consulting/R&D or related fields. Can I use FluidX3D commercially?

  
  No. Commercial use is not allowed with the current license.
  
  
• Is FluidX3D open-source?

  
  No. "Open-source" as a technical term is defined as freely available without any restriction on use, but I am not comfortable with that. I have written FluidX3D in my spare time and no one should milk it for profits while I remain uncompensated, especially considering what other CFD software sells for. The technical term for the type of license I choose is "source-available no-cost non-commercial". The source code is freely available, and you are free to use, to alter and to redistribute it, as long as you do not sell it or make a profit from derived products/services, and as long as you do not use it for any military purposes (see the license for details).
  
  
• Will FluidX3D at some point be available with a commercial license?

  
  Maybe I will add the option for a second, commercial license later on. If you are interested in commercial use, let me know. For non-commercial use in science and education, FluidX3D is and will always be free.
  
  

• OpenCL-Headers for GPU parallelization (Khronos Group)
• marching-cubes tables for isosurface generation on GPU (Paul Bourke)
• src/lodepng.cpp and src/lodepng.hpp for .png encoding and decoding (Lode Vandevenne)
• SimplexNoise class in src/utilities.hpp for generating continuous noise in 2D/3D/4D space (Stefan Gustavson)
• skybox/skybox8k.png for free surface raytracing (HDRI Hub)

• Lehmann, M.: Esoteric Pull and Esoteric Push: Two Simple In-Place Streaming Schemes for the Lattice Boltzmann Method on GPUs. Computation, 10, 92, (2022)
• Lehmann, M., Krause, M., Amati, G., Sega, M., Harting, J. and Gekle, S.: Accuracy and performance of the lattice Boltzmann method with 64-bit, 32-bit, and customized 16-bit number formats. Phys. Rev. E 106, 015308, (2022)
• Lehmann, M.: Combined scientific CFD simulation and interactive raytracing with OpenCL. IWOCL'22: International Workshop on OpenCL, 3, 1-2, (2022)
• Lehmann, M., Oehlschlägel, L.M., Häusl, F., Held, A. and Gekle, S.: Ejection of marine microplastics by raindrops: a computational and experimental study. Micropl.&Nanopl. 1, 18, (2021)
• Lehmann, M.: High Performance Free Surface LBM on GPUs. Master's thesis, (2019)
• Lehmann, M. and Gekle, S.: Analytic Solution to the Piecewise Linear Interface Construction Problem and Its Application in Curvature Calculation for Volume-of-Fluid Simulation Codes. Computation, 10, 21, (2022)

• For any questions, feedback or other inquiries, don't hesitate to contact me at [email protected].
• Updates will be posted on Twitter via @FluidX3D and @ProjectPhysX, under the hashtag #FluidX3D or on my YouTube channel.