Tutorial · CFD · Algorithms · JavaScript
📅 July 2026 ⏱ ≈ 30 min 🎯 Intermediate – Advanced

Lattice-Boltzmann Method in About 200 Lines of JavaScript

You don't need a Navier-Stokes solver, a Poisson pressure equation, or a mesh generator to simulate real fluid flow. The Lattice-Boltzmann Method gets you flow around a cylinder — with a proper Kármán vortex street — from a handful of local, embarrassingly parallel array operations.

1. The D2Q9 Lattice

"D2Q9" means 2 spatial dimensions, 9 discrete velocities. A particle population can either sit still, or move to one of its 8 neighbours (4 axis-aligned + 4 diagonal) each timestep. Every direction has a weight w[i] used in the equilibrium formula: 

// D2Q9 velocity set: index 0 = rest, 1-4 = axis, 5-8 = diagonal
const ex = [0, 1, 0, -1, 0, 1, -1, -1, 1];
const ey = [0, 0, 1, 0, -1, 1, 1, -1, -1];
const w  = [4/9, 1/9, 1/9, 1/9, 1/9, 1/36, 1/36, 1/36, 1/36];
const opp = [0, 3, 4, 1, 2, 7, 8, 5, 6]; // opposite direction, for bounce-back
Lattice units: everything here is in "lattice units": Δx = Δt = 1, lattice sound speed cs = 1/√3. Physical viscosity maps to the relaxation time τ via ν = cs²(τ − 0.5).

2. Data Layout

Nine flat Float64Arrays (one per direction) hold the whole grid's distributions. A boolean mask marks obstacle cells — here, a circle placed a third of the way into the domain: 

const NX = 200, NY = 80;
const N  = NX * NY;
const f  = Array.from({length: 9}, () => new Float64Array(N));
const feq = Array.from({length: 9}, () => new Float64Array(N));
const rho = new Float64Array(N).fill(1);
const ux  = new Float64Array(N);
const uy  = new Float64Array(N);
const obstacle = new Uint8Array(N);

const cx = NX / 5, cy = NY / 2, R = NY / 9;
for (let y = 0; y < NY; y++)
  for (let x = 0; x < NX; x++)
    if ((x - cx)**2 + (y - cy)**2 < R*R)
      obstacle[y * NX + x] = 1;

3. Macroscopic Variables

Density and velocity fall directly out of the distributions as the zeroth and first velocity moments — no extra equation to solve:

ρ = Σᵢ fᵢ u = (1/ρ)·Σᵢ fᵢ·eᵢ
function computeMacroscopic() {
  for (let n = 0; n < N; n++) {
    let r = 0, vx = 0, vy = 0;
    for (let i = 0; i < 9; i++) {
      const fi = f[i][n];
      r += fi; vx += fi * ex[i]; vy += fi * ey[i];
    }
    rho[n] = r;
    ux[n]  = vx / r;
    uy[n]  = vy / r;
  }
}

4. Equilibrium Distribution

The BGK equilibrium is a second-order truncated expansion of the Maxwell-Boltzmann distribution in local velocity u:

fᵢ^eq = wᵢ·ρ·[1 + 3(eᵢ·u) + 4.5(eᵢ·u)² − 1.5|u|²]
function computeEquilibrium() {
  for (let n = 0; n < N; n++) {
    const r = rho[n], vx = ux[n], vy = uy[n];
    const usq = vx * vx + vy * vy;
    for (let i = 0; i < 9; i++) {
      const eu = ex[i] * vx + ey[i] * vy;
      feq[i][n] = w[i] * r * (1 + 3 * eu + 4.5 * eu * eu - 1.5 * usq);
    }
  }
}

5. Collision (BGK Relaxation)

Every distribution relaxes toward its equilibrium at rate 1/τ. This single subtraction is where viscosity lives — larger τ means slower relaxation and a more viscous fluid:

fᵢ ← fᵢ − (1/τ)·(fᵢ − fᵢ^eq)
const tau = 0.6;  // τ > 0.5 required for stability; ν = (τ-0.5)/3
function collide() {
  for (let i = 0; i < 9; i++)
    for (let n = 0; n < N; n++)
      if (!obstacle[n])
        f[i][n] -= (f[i][n] - feq[i][n]) / tau;
}

6. Streaming

Streaming shifts every distribution one lattice step along its own direction — a pure array-copy with wraparound skipped at the domain edges (handled separately by boundary conditions): 

function stream() {
  for (let i = 0; i < 9; i++) {
    const src = f[i], dst = new Float64Array(N);
    for (let y = 0; y < NY; y++)
      for (let x = 0; x < NX; x++) {
        const xs = x - ex[i], ys = y - ey[i]; // pull from upstream neighbour
        if (xs < 0 || xs >= NX || ys < 0 || ys >= NY) continue;
        dst[y * NX + x] = src[ys * NX + xs];
      }
    f[i] = dst;
  }
}
Performance note: production LBM codes stream in place with a clever index trick (or two buffers swapped per direction) to avoid the full array allocation shown here — but the pull-based formulation above is the clearest to read and to get correct first.

7. Boundary Conditions

Three different rules close the domain: a fixed-velocity inlet (simplified Zou-He), a zero-gradient outlet, and bounce-back for solid walls and the cylinder — by far the simplest no-slip condition in LBM, just reverse the incoming populations in place: 

function applyBoundaries() {
  // Bounce-back: obstacle nodes reflect every incoming direction
  for (let n = 0; n < N; n++) {
    if (!obstacle[n]) continue;
    const tmp = new Float64Array(9);
    for (let i = 0; i < 9; i++) tmp[i] = f[i][n];
    for (let i = 0; i < 9; i++) f[i][n] = tmp[opp[i]];
  }

  // Left edge (inlet): force a uniform horizontal velocity u0
  const u0 = 0.08; // keep well below c_s/√3 ≈ 0.577 for stability
  for (let y = 0; y < NY; y++) {
    const n = y * NX;
    ux[n] = u0; uy[n] = 0;
    let s = 0;
    for (let i = 0; i < 9; i++) if (ex[i] <= 0) s += f[i][n];
    rho[n] = s / (1 - u0);
    for (let i = 0; i < 9; i++) {
      const eu = ex[i] * u0;
      f[i][n] = w[i] * rho[n] * (1 + 3 * eu + 4.5 * eu * eu - 1.5 * u0 * u0);
    }
  }

  // Right edge (outlet): zero-gradient — copy the column just upstream
  for (let y = 0; y < NY; y++) {
    const n = y * NX + (NX - 1), nPrev = n - 1;
    for (let i = 0; i < 9; i++) f[i][n] = f[i][nPrev];
  }

  // Top / bottom: bounce-back (solid channel walls)
  for (let x = 0; x < NX; x++) {
    for (const y of [0, NY - 1]) {
      const n = y * NX + x;
      const tmp = new Float64Array(9);
      for (let i = 0; i < 9; i++) tmp[i] = f[i][n];
      for (let i = 0; i < 9; i++) f[i][n] = tmp[opp[i]];
    }
  }
}

8. Rendering Vorticity

Raw velocity is hard to read visually; vorticity (curl of velocity) makes the alternating Kármán vortices pop as red/blue bands. A simple central-difference curl and a diverging red-white-blue colour map do the job: 

function renderVorticity(ctx, img) {
  const data = img.data;
  for (let y = 1; y < NY - 1; y++)
    for (let x = 1; x < NX - 1; x++) {
      const n = y * NX + x;
      const dvdx = uy[n + 1] - uy[n - 1];
      const dudy = ux[n + NX] - ux[n - NX];
      const curl = (dvdx - dudy) * 40; // scale for visibility
      const p = (y * NX + x) * 4;
      if (obstacle[n]) { data[p] = data[p+1] = data[p+2] = 40; }
      else if (curl > 0) { data[p] = clamp(curl * 255); data[p+1] = 20; data[p+2] = 20; }
      else                { data[p] = 20; data[p+1] = 20; data[p+2] = clamp(-curl * 255); }
      data[p+3] = 255;
    }
  ctx.putImageData(img, 0, 0);
}
function clamp(v) { return Math.max(0, Math.min(255, v)); }

9. Wiring It Together — the Main Loop

Each animation frame runs the four steps in order — macroscopic, collide, stream, boundaries — then renders. That's the entire LBM algorithm, under 200 lines including setup and rendering: 

function step() {
  computeMacroscopic();
  computeEquilibrium();
  collide();
  stream();
  applyBoundaries();
}

const canvas = document.getElementById('lbm');
canvas.width = NX; canvas.height = NY;
const ctx = canvas.getContext('2d');
const img = ctx.createImageData(NX, NY);

// Initialize at rest with a small rightward velocity everywhere
computeEquilibriumAt(0.05);
function loop() {
  for (let s = 0; s < 4; s++) step(); // a few LBM steps per rendered frame
  renderVorticity(ctx, img);
  requestAnimationFrame(loop);
}
loop();
Reading the result: at Re ≈ UD/ν ≈ 100–200 (tune τ and u0), the wake behind the cylinder should start shedding alternating red/blue vortex pairs — the Kármán vortex street described in Steady and Unsteady Flow Around a Cylinder. For the full theory behind why FDM, FVM, and LBM give the same physics through very different discretizations, see Solving Navier-Stokes Numerically: FDM vs FVM vs LBM.
▶ Live Demo

Frequently Asked Questions

What will I learn in this tutorial?

Build a working D2Q9 Lattice-Boltzmann CFD solver from scratch in about 200 lines of JavaScript: streaming, BGK collision, bounce-back walls, and a vorticity-colored canvas render of flow around a cylinder.

What topics are covered in this tutorial?

This tutorial covers: The D2Q9 Lattice, Data Layout, Macroscopic Variables, Equilibrium Distribution, Collision, Streaming, Boundary Conditions, Rendering Vorticity.

How long does this tutorial take?

This tutorial takes approximately 30 minutes to complete.

What prerequisites do I need before starting?

This is a Intermediate – Advanced-level tutorial — no special preparation beyond basic JavaScript is assumed.