# Lattice Boltzmann Method

(Redirected from Lattice Boltzmann Methods)
Post-publication activity

Curator: Sauro Succi

## Fluids, atoms and kinetic theory

The quantitative study of the flowing state of matter (fluids) relies upon two major levels of description: continuum and atomistic. In the continuum picture fluids are described in terms of space-filling fields, such as density, velocity, pressure, which vary smoothly in space and time. Such a description is fairly adequate for most applications, to the point that the physics of fluids is often implicitly identified with continuum fluid mechanics. Nevertheless, for over a century, it is known that fluids (gas and liquids) are ultimately composed by a collection of individual entities, atoms and molecules, whose discrete nature becomes apparent at scales around the nanometer and below. Much less known, at least to the general public, is the existence of a third, intermediate level of description, in which fluids are represented in terms of the probability (density) $$f(\vec{r},\vec{v};t)$$ of finding a given particle at a given position in space, $$r$$ and time $$t\ ,$$ with a given velocity, $$\vec{v}\ .$$ This intermediate (mesoscopic) level is described by kinetic theory, whose cornerstone is the celebrated Boltzmann's equation (Boltzmann, 1976, Cercignani, 1975):

$\tag{1} \frac{\partial f}{\partial t} + \vec{v} \cdot \nabla f = C(f,f)$

In the above equation, the left-hand side describes the free streaming of the particles (force-free for simplicity), while the right-hand sides describes the changes in the probability distribution function $$f(\vec{r},\vec{v};t)$$ induced by inter-particle collisions. Although much more economical than the atomistic description, the Boltzmann equation contains way more information than the continuum description, since the distribution function lives in a six-dimensional space (phase-space) spanned by the particle position $$\vec{r}$$ and velocity $$\vec{v}\ .$$ The Boltzmann equation is commonly regarded conceptually adequate only for rarefied gas dynamics, the reason being that Boltzmann originally derived it under the assumption of diluteness. Moreover, besides living in a six-dimensional space, the collision term involves a highly complicated and non-linear integral in velocity space which makes the equation very hard to treat even from the numerical point of view. Taken all together, it is no surprise to see kinetic theory still often regarded as a specialistic corner of fluid dynamics. Remarkably, modern developments in statistical mechanics and nonlinear dynamical systems, have conspired to change this picture. In particular, it is now accepted that kinetic theory, and Boltzmann-like equations, provide a fairly general mathematical framework to describe a wide host of non-equilibrium phenomena, straddling across the borders of different disciplines, such as physics, chemistry, biology and even the study of traffic flows. The key observation is that, even though the traditional kinetic theory of gases was originally derived for weakly-interacting (dilute) molecules, kinetic-like equations can be applied whenever strong interactions between elementary degrees of freedom (particles) can be cast in the form of weak interactions between collective degrees of freedom (quasi-particles). As we shall see, the lattice Boltzmann equation (LBE) falls precisely within this modern line of thinking.

## The Lattice Boltzmann equation: history

The lattice Boltzmann equation (LBE) is a minimal form of Boltzmann kinetic equation which is meant to simulate the dynamic behaviour of fluid flows without directly solving the equations of continuum fluid mechanics. Instead, macroscopic fluid dynamics emerges from the underlying dynamics of a fictitious ensemble of particles, whose motion and interactions are confined to a regular space-time lattice (Benzi and Succi, 1992; Chen and Doolen, 1998). Technically, the distinctive feature of LBE is a dramatic reduction of the degrees of freedom associated with the velocity space. In fact, particle velocities are restricted to a handful of discrete values $$\vec{v}=\vec{c}_i\ ,$$ $$i=0,b\ ,$$ by assuming that at each site the particles can only move along a finite number of directions. By endowing this set with sufficient symmetries to fulfill the basic conservation laws of mass, momentum and energy, the LBE can be shown to quantitatively reproduce the equations of motion (Chapman and Cowling, 1952) of continuum fluid mechanics, in the limit of long wavelengths as compared to the lattice scale. The idea of simplified Boltzmann equations with discrete speeds can be traced to the pioneering work of Broadwell, back in the 60's (Broadwell, 1964). However, to the best of our knowledge, these discrete-velocity Boltzmann equations were mostly intended to provide simpler, sometimes even analytically tractable, model equations for rarefied gas dynamics, but never meant as a computational alternative for the numerical solution of the Navier-Stokes equations of continuum fluid dynamics. The major conceptual twist in this direction had to wait another 20 years, with the advent of the celebrated Frisch-Hasslacher-Pomeau (FHP) lattice gas automaton (Frisch et al., 1986; Wolfram, 1986, Wolf-Gladrow, 2000). Originally developed in response to the major pitfalls of the Lattice Gas Cellular Automaton approach (Frisch et al., 1986; Wolf-Gladrow, 2000), the LBE rapidly developed into a vigorous self-standing research subject (Mc Namara and Zanetti, 1988; Higuera and Jimenez, 1989; Higuera et al., 1989; Chen et al., 1992). The distinctive features of LBE as a computational solver for fluid problems are its space-time locality, and the fact that information travels along the straight lines defined by the (constant) particle velocities associated with the lattice, rather than along the space-time dependent material lines defined by the flow speed. Due to these properties, the LB approach counts today an impressive array of applications across virtually all fields of fluid dynamics and allied disciplines, such as biology and material science (Succi, 2001, 2008).

## The Lattice Boltzmann equation: mathematical formulation

The lattice Boltzmann equation reads as follows (Wolf-Gladrow, 2000; Succi, 2001)$\tag{2} f_i(\vec{r} +{\vec{c}}_i \Delta t, t+\Delta t) = f_i(\vec{r};t) + \Omega_{ij} (f^{eq}_j(\vec{r};t)-f_j(\vec{r};t))$

In the above $$f_i(\vec{r};t)$$ represents the probability of finding a particle at position $$\vec{r}$$ and time $$t$$ with velocity $$\vec{v}=\vec{c}_i$$ (see Figure 1).

Here, the subscript $$i$$ labels a set of discrete speeds ($$i=0...b$$) connecting the nodes of a regular lattice. A typical choice in two and three dimensions are the $$9$$-speed and $$19$$-speed lattices depicted in Figure 2 and Figure 2.

The left hand side of (2) corresponds to an exact discrete-velocity representation of the Boltzmann streaming operator $$\partial_t f + {\vec{v}} \cdot \nabla f\ ,$$ while the right-hand side stands for particle-collisions, represented through a relaxation to a local equilibrium. Such relaxation is controlled by the scattering matrix $$\Omega_{ij}\ ,$$ whose eigenvalues define the time-scale for local equilibration of the various kinetic moments. The local equilibrium $$f_i^{eq}$$ is the lattice analogue of a local Maxwellian with density $$\rho=\sum_i f_i$$ and flow velocity $${\vec{u}} =\sum_i f_i {\vec{c}}_i/\rho$$ (see below). In practical applications, the scattering matrix is often taken in its simplest diagonal form, resulting in the so called Lattice BGK (LBGK) scheme (Qian et al., 1992; Chen et al., 1991), after the continuum model Boltzmann equation introduced by Bhatnagar-Gross-Krook, back in 1954 (Bhatnagar et al., 1954). The discrete local equilibria are typically given by a second-order expansion in the Mach-number of a local Maxwellian, $\tag{3} f_i^{eq}(\vec{r};t) = w_i \rho \left( 1 + \frac{u_a c_{ia}}{c_s^2} + \frac{u_a u_b Q_{iab}}{2 c_s^4} \right)$

where $$w_i$$ is a set of weights normalized to unity and $$Q_{iab}=c_{ia}c_{ib}-c_s^2 \delta_{ab}$$ is the quadrupole projector along the $$i$$-th direction and $$c_s^2 = \sum_i w_i c_i^2$$ is the lattice sound speed. In the above, subscripts a,b denote cartesian components, and repeated indices are summed upon. The reason for expanding to second order in the flow field, more precisely in the local Mach-number $$Ma=u/c_s\ ,$$ is related to the lack of symmetry of a given lattice to ensure the isotropy of higher tensors coupling to higher powers of the Mach-number. Fortunately, isotropy of fourth-order tensors $$T_{abcd} = \sum_i w_i c_{ia} c_{ib} c_{ic} c_{id}\ ,$$ is sufficient to recover the correct behaviour of low-Mach, quasi incompressible flows. More demanding applications, with major thermal/compressibility effects, call for higher-order lattices, securing the correct symmetry of higher-rank tensors. The analysis and the identification of suitable lattices fulfilling the correct symmetry relations is facilitated by systematic expansions of the generating function of Hermite polynomials: $e^{-(v-u)^2/2v_T^2} = e^{-v^2/2 v_T^2} \sum_{n=0}^{\infty} \left(\frac{u}{v_T} \right)^n H_n \left(\frac{v}{v_T} \right)$ where $$v_T = \sqrt{kT/m}$$ is thermal speed and $$H_n$$ is the $$n$$-th order tensorial Hermite polynomial (the scalar notation is just for simplicity).

### Hydrodynamics from LBE

Once the discrete populations $$f_i$$ are known, the main fluid quantities are obtained by simple (linear) summation upon the discrete speeds, namely $\rho(\vec{r};t)= m \sum_{i=0}^b f_i(\vec{r};t)$

$\rho(\vec{r};t) \vec{u}(\vec{r};t) = m \sum_{i=0}^b f_i(\vec{r};t) \vec{c_i}$

${\vec{\vec{P}}}(\vec{r};t) = m \sum_{i=0}^b f_i(\vec{r};t) \vec{c_i} \vec{c_i}$

which define the mass density, current density and momentum-flux tensor of the fluid, respectively. It should be appreciated that, while $$\rho$$ and $$\rho \vec{u}$$ are conserved quantities, only the trace of the momentum-flux tensor is conserved, which corresponds to the fluid energy (pressure). The non-diagonal components of $$P_{ab}$$ are zero at equilibrium, and represent the fluid shear. Chapman-Enskog analysis (Chapman and Cowling, 1992) of the LBE shows that in the limit of long-wavelengths, low-frequency, the eq. (2) reproduces exactly the Navier-Stokes equation for quasi-incompressible flows, with an ideal equation of state $\tag{4} P(\rho) =\rho c_s^2,$

and a kinematic viscosity $\tag{5} \nu = c_s^2 (\tau-\Delta t/2),$

$$\tau$$ being the leading non-zero (inverse) eigenvalue of the matrix $$\Omega_{ij}\ .$$

The conceptual and practical simplicity of this computational scheme lies at the heart of the computational efficiency of the LBE, and form a solid basis for its success as a Navier-Stokes simulator 'in kinetic disguise'. The main computational highlights of LB are as follows:

i) Non-locality (streaming) is linear and non-linearity (collision) is local. This should be contrasted with Navier-Stokes equations, in which the transport term $$({\vec{u}}\cdot {\vec{\nabla}}) {\vec{u}}$$ is non-local and non-linear at a time. In particular, this term shows that the fluid momentum is transported along the material streamlines $$d {\vec{r}} = {\vec{u}} dt\ ,$$ defined by the fluid velocity itself. For complex flows, such streamlined can become very "wild" and demanding in terms of numerical stability.

ii) Streaming is exact, since it proceeds along constant streamlines $$d {\vec{r}}_i={\vec{c}}_i dt$$ (lightcones);

iii) Complex boundary conditions are readily formulated in terms of elementary mechanical rules (bounce-back, reflections) of interactions of the LB 'molecules' with solid walls;

iv) Fluid pressure and the strain tensor are available locally, on-the-fly, as linear combinations of the equilibrium and non-equilibrium components of the discrete distribution function, respectively. This means that pressure can be obtained from the flow field configuration with no need of solving a (usually very expensive) Poisson problem.

v) Nearly ideal amenability to parallel computing (low communication/computation ratio)

These aspects configure LBE as a very special finite-difference scheme for hydrodynamics and lie at the heart of its impressive growth as an alternative numerical technique for complex fluid dynamic problems.

## Boundary conditions

One of the major strengths of the LB technique is its capability of handling boundary conditions associated with highly irregular geometries by means of elementary mechanical operations on the discrete distribution function. This capability is best illustrated by the important case of non-slip boundary conditions, whereby the fluid velocity at solid (non-moving) walls is set to zero: $\tag{6} \vec{u}(\vec{r}_W) = 0$

where $$\vec{r}_W$$ denotes a generic location on the solid wall. The simplest way to implement such boundary conditions in a LB framework is to reflect (bounce-back) the populations living on the wall-site $$\vec{r}_W\ ,$$ namely $\tag{7} f_{in}(\vec{r}_W) = f_{out}(\vec{r}_W)$

where subscript "in" and "out" denote the set of incoming and outgoing populations at the wall, respectively. Figure 3: The bounce-back implementation in the lattice Boltzmann models. The discrete nature in the boundary condition is illustrated. The incoming populations at a wall are equal to their mirror partners, coming out from a buffer node at the boundaries. For the sake of concreteness the case of D2Q9 is reported. The wall location can be thought as being located halfway between the fluid and buffer nodes (He et al., 1997).

The rationale behind the bounce-back rule is quite simple (see also Figure 3); by reflecting the populations at each given wall-site, the net current on that site is forced to zero by construction, regardless of the geometrical complexity of the boundary. This rule illustrates well the spirit of simplicity which inspires the set-up of boundary conditions in the LB framework. However, it also hides the subtle complexities which arise whenever the boundary is not aligned with the lattice (misaligned boundaries) (Ginzbourg and d’Humieres, 1996), which is the rule rather than the exception in most engineering applications. Dealing with misaligned boundaries calls for various types of interpolation techniques, too technical to be described in this introductory work. The crucial point, though, is that these interpolations techniques are greatly aided by the fact that the populations/particles keep moving along well-prescribed straight lines, whose intersections with arbitrarily complex boundaries are conceptually straightforward to be computed. In actual practice, however, this may become quite labour-intensive for real-life, complex three-dimensional geometries. A general formulation of such boundary conditions can be found in (Chen et al., 1998). For other formulations, merging LB with finite-volume techniques, see (Peng and Duncan, 2003; Ubertini et al., 2003).

Besides sheer geometrical complexity, another double-edged aspect of LB boundary conditions, is their exposure to potential mathematical ambiguities, whenever the conditions to be imposed on physical quantities are outnumbered by the LB populations to be specified at the boundary. Indeed, most common hydrodynamic boundary conditions involve the fluid density, pressure, velocity and its gradients, at a given boundary. There is no general guarantee that the number of specified macroscopic conditions be in exact match with the number of incoming LB populations. The excess of LB variables must then be resolved through educated guesses on the non-equilibrium distribution function (the equilibrium component is uniquely specified by the hydrodynamic variables). The literature on LB boundary conditions has grown considerably over the last decade, and understandably so, with an increasing penchant towards technical details beyond the scope of this introductory work. For an informed, although necessarily non-exhaustive, overview see Boghosian in (Ubertini et al., 2009).

## Post-classical LBE's

The foundations of the LB method were laid down in a short series of trailblazing papers between 1988-1992, culminating with the all-popular LBGK scheme. A number of important conceptual developments followed in the subsequent years. Among them, a fairly prominent one is the so-called entropic LB (ELB), i.e. Lattice Boltzmann schemes equipped with a H-theorem.

### Restoring the second principle: Entropic Models

One the most profound results of Boltzmann's kinetic theory was the establishment of a direct link between mechanics and thermodynamics, in the form of his celebrated H-theorem. This theorem states that the functional, $$H(t)= \int f(\vec{r},\vec{v};t) \log f(\vec{r},\vec{v};t) d \vec{v} d \vec{r}\ ,$$ cannot increase in the course of the evolution of $$f$$ under the Boltzmann dynamics. It can only decrease, until a steady-value is reached, which signals the attainment of thermodynamic equilibrium in which no change can take place any longer. The identification of $$-H$$ with the system entropy is immediate. The profound implications of the $$H$$-theorem have raised long-lasting controversies. To the purpose of LB theory, however, the interest towards a lattice H-theorem is motivated mostly by reasons of physical realizability and ensuing numerical stability issues. Indeed, by specifying the local equilibria at the outset, instead of realizing them through a sequence of collisions, it is clear that compliance with the original Boltzmann's H-theorem is no longer guaranteed. The earliest computationally competitive LBE was originally designed top-down, in compliance with the basic conservation laws of the physics of fluids, with somewhat less attention to the second principle of thermodynamics. A major result in LB theory was to realize that some LBE models do indeed support a discrete H-theorem. This means that the lattice equilibria, (3), can be shown to represent second-order expansions in the Mach-number of the minimizers of the following lattice $$H$$-function (Karlin et al., 1999)$H=\sum_i f_i \log (f_i/w_i)\ .$ This is an important result, since it puts the LBE back on a sound physical basis. More precisely, it shows that the top-down relaxation dynamics represented by the equation (2) is realizable, i.e. not only it fulfills local conservation laws, but also the evolutionary constraints imposed by the 2nd principle.

### Lattice Boltzmann from Boltzmann: the Hermite-Gauss connection

Years later, it was also realized that LB can be formally derived by a Gauss-Hermite projection of the corresponding continuum BGK equation (Abe, 1997; He and Luo, 1997; Shan and He, 1998): This realization elucidates the nature of the discrete velocities and place LBE on a more systematic mathematical basis (Grad, 1949). For a comprehensive treatment of higher-order LB schemes based on Hermite projections, see (Shan et al., 2006; Philippi et al., 2006).

### Multi-time relaxation LBE

As previously mentioned, practical applications based on the LBE are largely dominated by the single-time relaxation, LBGK form. This is mostly a matter of simplicity, for it is known that LBGK, besides suffering of some obvious physical limitations (with just a single relaxation parameter, the transport of mass, momentum and energy are forced to proceed at the same pace), is also exposed to a number of numerical inaccuracies, (Lallemand and Luo, 2000). To cope with these problems, a number of authors have developed optimized versions of the original collision matrix in scattering form (Higuera and Jimenez, 1989; Higuera et a., 1989) now known as MRT (Multi-Relaxation Time) and demonstrated enhanced stability. Essentially each hydrodynamic mode is made to relax towards its equilibrium at a different rate. The separate relaxation was also demonstrated in a more general lattice-independent fashion (Shan and Chen, 2007). By and large, the above limitations and inaccuracies are not perceived to overshadow its simplicity, which is why LBGK still holds the lion's share in the field. Nevertheless, se MRT approach is gaining increasing popularity for engineering applications.

## Idiosyncrasies, lucky strikes and some misbeliefs

### Weak compressibility and pressure drops

An often asked question is how can LB possibly describe liquids, given the fact that it obeys an ideal equation of state, eq. (4). The answer is that, for incompressible flows (incompressibility here refers to the flow, not to the fluid!), the equation of state does not play any role, as long as the incompressibility condition $${\vec{\nabla}} \cdot {\vec{u}} = 0$$ is fulfilled. As a matter of fact, LB is not exactly incompressible, for it supports density changes of the order of $$\delta \rho/\rho \sim O(Ma^2)\ .$$ Recall that lattice equilibria are a second-order expansion of a local Maxwellian and consequently, at variance with the latter, they are positive-definite only for sufficiently low Mach numbers, typically $$Ma < 0.2\ .$$ As a result, standard LB can only support pressure changes of the order of $$\delta P/P \sim Ma^2\ .$$ Whenever larger pressure gradients are required, the common practice is to represent them through an external force $$F^{ext} = -\nabla P\ .$$

### The mirage of zero viscosity

The formula (5) of the LB viscosity contains a slight, and yet crucial, departure from the continuum analogue, namely the negative contribution $$-c_s^2 \Delta t/2\ .$$ This negative viscosity, which stems from Taylor expansion of the discrete free-streaming operator to second order in the lattice time-step, is a golden nugget for LB applications to low-viscous and turbulent flows. The reason is that it permits to achieve vanishingly low viscosities, $$\nu \sim \epsilon <<1\ ,$$ with time steps $$\Delta t \sim O(1)\ .$$ This magic is realized by choosing $$\tau = \Delta t/2 + \epsilon\ ,$$ which delivers $$\nu = c_s^2 \epsilon\ .$$ This expression naturally invites a greedy question: how low can one go with $$\epsilon\ ?$$ Numerical practice shows that, with ordinary lattice resolutions ($$100-1000$$ grid points per linear dimension), the viscosity can be brought down to $$0.01-0.001$$ lattice units. Below such values, non-linear instabilities are usually triggered, due to the fact that the lattice viscosity is too low to dissipate the shortest wavelengths. This figure can however be lowered by several orders of magnitude by the Entropic LB previously mentioned, thus taking LB much closer to the "mirage" of zero viscosity.

### LB is -not- a dilute-gas approximation

Another common catch, sometimes fueled even by published work, is that LBE should apply only to dilute gases, the reason being that it derives from an approximation of the true Boltzmann equation, which is itself valid only for dilute gases. This line of thinking is overly restrictive. First, it fails to recognize that, although kinetic theory is originally meant to describe weakly-interacting (dilute) systems, it can also be applied whenever strong interactions between elementary degrees of freedom can be cast in the form of weak interactions between appropriate collective degrees of freedom (quasi-particles). LB falls within this category, since each LB population is hydrodynamically representative of a large collection of real molecules. The second source of confusion is to assume that, since BGK is a model approximation to the true Boltzmann equation, it should inherit all of its limitations, and more. This is not necessarily the case, because BGK is based on two fairly universal assumptions, namely:

i) There exist a local attractor to the collisional dynamics (local equilibrium)

ii) In the vicinity of this local equilibrium, the collisional dynamics can be described in terms of a simple relaxation process, the relaxation time-scale $$\tau\ ,$$ fixing the viscosity of the fluid.

Assumption i) is certainly true for any statistical-dynamics system supporting microscopic invariants, while assumption ii) is also fairly plausible, possibly even rigorous, for generic systems with short-range interactions. For sure, neither of the two is restricted to dilute gases.

## LBE for non-ideal fluids

One of the most brilliant outgrowths of the basic LB theory is the extension to the case of interacting (non-ideal) fluids. Many approaches have been devised to this purpose, (Gustensen, 1991; Shan and Chen, 1993; Shan and Chen, 1994; Swift et al., 1996), but for the sake of coinciseness, here we shall refer to the pseudo-potential approach developed by Shan and Chen (Shan and Chen, 1993; Shan and Chen, 1994). The Shan-Chen (SC) approach augments the LB equation with the following force term

$$\tag{8} F_a(\vec{r};t) = -G \Psi[\rho(\vec{r},t) ] \sum_i w_i \Psi[\rho(\vec{r}+\vec{c}_i \Delta t,t)] c_{ia}$$

where $$G$$ is a coupling strength and $$\Psi$$ is a density-dependent pseudo-potential, incorporating the effects of potential energy interactions through a simple nearest-shell coupling. Positive/negative values of $$G$$ code for repulsion/attraction respectively. Figure 5: Sketch of the isothermal Shan-Chen bulk equation of state $$P_b(\rho)$$ as a function of the inverse density, i.e. the specific volume of the fluid. For sufficiently large (attractive) interaction strength the model is able to describe a phase-transition with liquid-vapor coexistence. The constant pressure in the region of the coexistence is provided by the equivalent of the Maxwell rule for the lattice Boltzmann model (Shan and Chen, 1993; Shan and Chen, 1994).

In spite of its simplicity, the SC model features two essential ingredients of non-ideal flows, namely i) a non-ideal bulk equation of state, $$P_b(\rho)=c_s^2 \rho + c_s^2 G \frac{\Psi^2}{2}\ ,$$ and a non-zero surface tension $$\gamma \propto -G \int (\nabla \Psi)^2 dl\ ,$$ where the integral runs across the fluid interface (Shan and Chen, 1993; Shan and Chen, 1994). Being purely attractive ($$G<0$$), the SC model is exposed to density collapse instabilities. The natural antidote in a real fluid is hard-core repulsion, but this is known to hinder the numerical stability of the scheme. Therefore, in the SC approach, the attractive instability is tamed through a high-density saturation of the pseudo-potential, which is typically taken in the form: $\tag{9} \Psi[\rho]=\rho_0 (1-e^{-\rho/\rho_0})$

where $$\rho_0$$ is a reference density marking the border between the dense and light phases. The role of this pseudo-potential is to produce an exponentially vanishing force as the density increases. This saturation is key to the numerical stability of the method. The SC method is simple, elegant and efficient, but in its original form, it also comes with a number of limitations. First, both the equation of state and surface tension are governed by the same parameter, $$G\ ,$$ which means that they cannot be changed independently. Second, the method is affected by the so-called spurious currents, i.e. artificial structures around the interface, which owe their existence solely to the lack of (higher order) isotropy of the lattice. These spurious currents prevent stable operation with density jumps between the liquid and vapour phases above $$50-100\ .$$ Recently the origin of numerical instability has been attributed to the high sound speed dictated by the non-realistic EoS (9), which results in a greater-than-one acoustic Courant number (Yuan and Schaefer, 2006, Kupershtokh 2010). With more realistic EoS, density-ratios as high as realistic values of normal fluids have been reportedly achieved (Yu and Fan, 2009). Recent extensions of the SC model, accounting for multi-range interactions extending beyond the first Brillouin region, have recently shown promising potential for a variety of complex flows, such foams, emulsions and ferro-fluids (Falcucci et al, 2007, Sbragaglia et al., 2007, Shan, 2007, Falcucci et al, 2009, Benzi et al., 2009). Figure 6: Rupture of a ferrofluid droplet in external magnetic fields using a single-component lattice Boltzmann model for nonideal fluids. From (Falcucci and Succi, 2009), and Phys. Rev. E Kaleidoscope, May 2009, pre.aps.org/node/2731.

## LB beyond Navier-Stokes: meso and nano-hydrodynamics

The LB method was born as a computational alternative to the discretization of the partial differential equations of continuum mechanics (Navier-Stokes equations), hence in connection with the numerical simulation of macroscopic fluids. In the recent years, however, LB has made proof of a remarkable versatility in describing a wide host of micro and nano-fluidic situations for which the continuum picture may become inaccurate or even physically untenable, and yet sufficiently universal to dispense with atomistic details. For instance, by suitable adaptation of the boundary conditions, LB has proven capable of describing typical beyond-hydro phenomena, such as slip-flows and other features of microfluidics (Ansumali, 2002), with quantitative accuracy.

Another important extension of the basic LB method concerns fluids with suspended objects, both rigid and flexible ones, typically colloids and polymers respectively (Ladd, 1994). The suspended body is usually treated as a collection of grid-free points (markers), linked together by either rigid constraints or some form of coarse-grained potentials. In the former case, the fluid-body interaction is treated as a collision of the LB particles with a moving solid surface, while in the latter the markers are evolved according to effective Newton's equation of motion, where the force terms account for the various interactions with the fluid (viscous drag). The major advantage of the LB approach for this type of flows is the fact that the fluid interactions with the suspended bodies are treated explicitly, with no resort to long-range, action-at-distance, hydrodynamic interactions. The second advantage is the usual freedom from Poisson solvers for the pressure. The price is small time-steps. Once again, capitalizing on these assets in a quantitative fashion requires a great deal of attention to technical details (Ahlrichs and Duenweg, 1998; Duenweg and Ladd, 2008). In general, a comprehensive LB theory of beyond-hydrodynamic regimes remains to be developed. Figure 7: Two-strand biopolymer translocation through a membrane. The calculation couples LB for the hydrodynamic solvent with Molecular Dynamics for the biopolymer (Fyta et al., 2006).

## Whither LB?

LB methods can describe virtually all transport phenomena controlled by advection-diffusion-reaction mechanisms, as well as non-ideal gas effects driven by potential energy interactions. Whence their impressively broad range of applications. Amazingly, LB can also describe quantum fluids (Succi and Benzi, 1993), although, at the time of this writing, the overwhelming majority of LB applications is on classical fluids, especially multiphase flows and flows with suspended bodies, for which the solution of the field equations of continuum fluid mechanics is either physically inadequate or conceptually opaque and computationally over-demanding. Besides the natural realm of macroscopic fluid dynamics, such as turbulence, aerodynamics, multiphase flows, hemodynamics (Chopard 2005, Hoekstra 2006, Melchionna 2009) the LB portfolio of applications keeps expanding across scales of motion, particularly towards micro and nano-fluidics (Harting et al. 2006, Hyvaluoma et al. 2009). Realizing the full potential of the method, up to engineering standards of accuracy, however, still raises many challenges ahead.