Entropy-Assisted Computing of Low-Dissipative Systems

Entropy feedback is reviewed and highlighted as the guiding principle to reach extremely low dissipation. This principle is illustrated through turbulent flow simulations using the entropic lattice Boltzmann scheme.


Introduction
The operation of many natural and engineered systems depends crucially on their ability to function at very low dissipation rates, the lower most often the better [1].Zero-dissipation, however, is an ideal limit which could only be reached if these systems could operate at virtually infinite processing speed.Hence, a very general question arises: how low can one keep dissipation in a given thermodynamic system?
Here we show that the ability of a specific class of fluid-kinetic systems [2][3][4] to function at a very low dissipation is dramatically enhanced by enforcing the second principle of thermodynamics in the form of an entropic feedback [5].Through concrete examples of turbulent flows, we highlight how entropy-assisted simulation maintains the system at low viscosity, through a highly orchestrated and self-consistent interplay between local enhancement and reduction of the dissipation.Balancing of these dissipation fluctuations leads to a spatial distribution of the average effective viscosity which keeps the simulation "alive and well".We envisage the entropy-assisted computing procedure to offer a general paradigm for the computer simulation of a wide class of low-dissipative complex phenomena, such as classical and quantum turbulence and wave propagation in active media.

Survival below Minimum Dissipation Threshold
The second principle of thermodynamics stands out as one of the most general and inescapable laws of physics, with profound bearings on the time evolution of virtually all natural systems [1].In its essence, it states that any natural system is driven towards a state of maximum entropy (equilibrium), characterized by a maximum number of microscopic configurations.However, it says little about a most relevant question: how long does it take for a given system to reach its equilibrium state?This question, the heart of non-equilibrium thermodynamics, is all but academic, since most natural phenomena, life in the first place, depend on the time the system is able to borrow from temporary elusion of the second principle [1].The rate of decay to equilibrium is measured by transport coefficients, such as kinematic viscosity, and can change widely from system to system, from seconds in an ordinary gas, to years and centuries in glassy materials.
The typical form of kinematic viscosity is given by ν ∼ v 2 T τ, where v T = √ k B T/m is the thermal speed, τ a typical relaxation time.The kinematic viscosity measures the diffusivity of momentum across the system: such diffusivity results from the competition between kinetic energy, which sustains the free motion of molecules, and potential energy, which controls their interactions (collisions).Kinetic energy drives the system out of equilibrium, while molecular collisions pull it back to a local equilibrium, in which entropy is locally maximized (Boltzmann's H-theorem) [6][7][8].Thus, a low-viscous fluid is not one with nearly no collisions, but one where collisions are so frequent and effective that they inhibit any migration of momentum from place to place, which is the source of macroscopic dissipation.From the above argument, it is seen that zero-viscosity is a mirage because it would imply instantaneous relaxation, i.e., τ → 0. Given that strictly zero-viscosity is a chimera in a real (finite-speed) world, a natural question arises: what is the minimum dissipation which can be sustained by a given physical system?While the answer depends on the specific system in mind, here we shall focus on discrete dynamical systems, i.e., featuring a fundamental minimal length scale a and minimal time scale h.
For the case of simulated fluids, for instance, the condition is that the smallest coherent structures (eddies) capable of surviving dissipation be resolved by the discrete grid, i.e., l K > a, where l K ∼ L/Re 3/4 is the so-called Kolmogorov's length [9], the smallest active scale in the game and Re = UL/ν is the Reynolds number, i.e., the ratio of nonlinear energy transfer to dissipation, for a fluid moving at a macroscopic speed U on a domain of macroscopic size L. As a result, the minimal viscosity is given by ν min = (u/N 1/3 )ν l , where N = L/a is the grid size, u ≡ U/U l , U l = a/h and ν l = a 2 /h being the natural lattice speed and lattice viscosity, respectively.Given that u < 1 for reasons of numerical stability, we see that the minimal viscosity is always smaller than the lattice viscosity, the ratio of the two decreasing like 1/N 1/3 , so that the minimum viscosity can be brought to zero only in the continuum limit N → ∞.Another face of the same chimera.The message is that fluids cannot support viscosity below their minimum bound ν min .Breaking such constraint leads to two basic scenarios: a mild reaction (loss of accuracy), whereby the resolved eddies, l > a, still survive, although with a corrupted dynamics, the degree of corruption increasing as they approach a.The second, more dramatic, possibility is loss of realizability: the system develops disruptive instabilities, typically in the form of an uncontrolled growth of the smallest eddies.This is nothing short of a survival problem, except that it concerns a discrete dynamical system.From the practical point of view, the art of keeping the system alive and well in the forbidden regime ν < ν min is known as turbulence modeling, a topic of utmost practical and conceptual importance.Essentially, the idea is to replace the nominal viscosity with an effective one, representing the effects of unresolved eddies as "random" collisions on the resolved ones.This picture explicitly draws upon an analogy with kinetic theory, where there is a clear scale separation between molecular and hydrodynamic degrees of freedom.Turbulence, on the contrary, features a continuum spectrum of scales, hence the notion of eddy viscosity, although very useful, still resists a rigorous justification.

Minimum Viscosity in Discrete Phase-Space-Time
However, a modern formulation of continuum fluid mechanics in a form which explicitly ingrains the discreteness of space-time is known as the lattice Boltzmann equation [4] In the above f i is the probability of finding a "particle" at position x in the lattice at time t, moving with discrete velocity c i along b lattice links; f mirr i is the so-called mirror state.In the simplest case [4], it is taken as f mirr i = 2 f eq i − f i , with f eq i the local equilibrium, which is a universal non-linear function of the local order parameters.For standard fluids, these are the fluid density n(x, t) = ∑ i f i (x, t) and velocity u(x, t) = n −1 ∑ i c i f i (x, t).The left-hand side of Equation (1) represents the free-streaming step, while the right-hand side describes the interactions among the discrete populations f i at each given lattice site.On condition that the lattice obeys proper symmetries, the lattice Boltzmann Equation (1) reproduces fluid dynamics with the viscosity The negative shift, −1/2 in Equation ( 2) is crucial; indeed, if only in principle, it permits to achieve zero viscosity in the limit τ → h/2, i.e., without sending τ → 0, unlike in the continuum.This negative shift is the result of the broken time-symmetry, which contributes a negative viscosity (sometimes called propagation viscosity) to the overall momentum diffusivity, besides the conventional contribution due to the collisional relaxation.Thus, in a discrete world, the viscosity receives contributions from both dynamical steps of the kinetic description: free streaming and collisions.They carry opposite signs, hence, if only in principle, they can cancel each other, leaving the time step and relaxation time both finite.This property, typically regarded as a very useful numerical artifact, has played a major role in the lattice Boltzmann simulation of a variety of complex flows, and most notably turbulent ones [10].Amazingly, suitably designed discrete kinetic systems keep describing correct fluid behavior several orders of magnitude below the hydrodynamic minimum viscosity bounds mentioned earlier in this paper.How come the minimal viscosity can be eluded by several orders of magnitude?
The key is the second principle in fully discrete setting.Indeed, the lattice Boltzmann Equation ( 1) is compatible with a discrete-time H-theorem, based on the H-function (negative of the entropy), , where w i are suitable positive-definite weights.Lattice Boltzmann systems equipped with the H-theorem are known as entropic [5] and function on a feedback mechanism, whereby the local relaxation τ is adjusted in space and time, so as to secure the entropic bound, where f and f are the pre-and post-collisional states, respectively.
The working principle is explained in Figure 1 and amounts to using in Equation ( 1) the entropy-supervised mirror state f mirr i = (1 − α) f i + α f eq , where the stretch α is found from the isentropic constraint, H[ f mirr ] = H[ f ].This can be interpreted as the effective viscosity, with the effective relaxation time, τ eff = 2τ/α.The entropy-assisted computation thus informs the pre-collision state f about its isentropic mirror f mirr and stipulates the single condition that the second law is respected by the post-collision state f .Whenever non-equilibrium effects become strong enough to endanger realizability, the entropic constraint adjusts the relaxation time so as to secure compliance with the second principle.This feedback is self-activated "on demand", i.e., only whenever and wherever the need arises.And when the danger is gone, most elegantly, the entropic feedback, leaves the stage unsolicited.The second principle decides by itself: sometimes viscosity is increased (ν eff > ν) to smooth out sharp features, sometimes it is reduced (ν eff < ν) to sharpen the dying ones.In the most demanding cases, the effective viscosity may even drop negative (τ eff < h/2) to promote local instabilities and sustain the system against dissipative death.The effective viscosity self-adapts to the actual state of turbulence to literally protect it against defective evolution and disruptive instabilities.
In Figure 2 we illustrate the above by the vorticity field of a flow past a circular cylinder at Re ∼ 3300, in which many active scales of motion are visible.The Reynolds number, Re = UD/ν, is based on the diameter of the cylinder, which is here taken as D = 30 a, while the mean flow velocity is U = 0.03 (a/h), corresponding to a viscosity ν = UD/Re ∼ 2.7 × 10 −4 (a 2 /h).With these parameters, the minimum viscosity is ν min = 0.03/30 1/3 ∼ 0.01 (a 2 /h), so that ν/ν min ∼ 1/40.The flow structures in Figure 2a are colored with the effective viscosity, normalized as R = (ν eff − ν)/ν.The high quality of resolution of the flow structures (vortex tubes, tangles etc.) is maintained by a concerted action of dampers (R > 0) and promoters (R < 0).The tiniest structures would not be able to survive unless the effective viscosity is enabled to go negative from time to time, in order to compensate for over-dissipation and "regenerate" small scale structures otherwise doomed by over-damping.Also to be noted (Figure 2b) is the spottiness of the effective viscosity, with a highly fine-grained mixing of dampers and promoters.All of the above configures a very elegant preemptive scenario which we can take as the hallmark of entropic computing: very attentive "guardian angels".Amazingly, the spatial pattern of the time-averaged effective viscosity shown in Figure 3 resembles indeed a "guardian angel", protecting the system against numerical crisis!   2 demonstrates a larger variation of the effective viscosity, the time-average picture is much milder: most of the activity (strongest damping neighboring the strongest promotion) is concentrated at and around the twin shear layers, just behind the cylinder.In the rest of the domain, the deviation of the effective viscosity from its nominal value is less than a fraction of a percent.

Entropic Lattice Boltzmann Algorithm
Here we describe the essentials of the simulation method used above.In the entropic lattice Boltzmann scheme, populations associated with the discrete velocities c i evolve according to (1).The local equilibrium f eq i was found by minimizing the entropy function H[ f ].The entropic mirror state is specified by the stretch α, which is computed as the positive root of the entropy condition: . Whenever the simulation is resolved at a particular grid node x, the stretch α becomes fixed automatically to the value α = 2 at that node, and the effective viscosity ν eff (3) reduces to the nominal viscosity ν (2).The stretch was obtained numerically at each grid point using Newton-Raphson method.For the simulation presented above, we used the lattice with b = 15 discrete velocities [11].Apart from an entropy-supported kinetic equation, we require augmenting boundary conditions that are capable of simulating both resolved and under-resolved flows.Existing boundary conditions such as the bounce-back scheme [4], provide reliable results for resolved simulations but with reducing grid sizes and increasing the Reynolds number, the quality of the simulations is lost due to shock-like instabilities generated at the walls.Hence, we used here the recently proposed Tamm-Mott-Smith boundary condition [12] for circumventing these instabilities.Entropy-supported kinetic Equation (1) provides reliable simulations for any choice of lattice and flow velocity as long as the Mach number remains small.Armed with stable boundary condition, the present scheme was extensively tested for various flow setups such as decaying turbulence (ν/ν min ∼ 10 −4 ) [13], turbulent channel (ν/ν min ∼ 10 −3 ) [12], grid generated turbulence (ν/ν min ∼ 10 −3 ), flow past an airfoil (ν/ν min ∼ 10 −4 ) and others.For the particular simulation of the flow past a circular cylinder presented in Figures 2 and 3, we used a computational domain that is 9D long in the span-wise direction, 35D along in the stream-wise direction with 10D upstream of the cylinder and 25D downstream of it; along vertical direction the domain was 21D long with cylinder axis in the mid-plane.The cylinder was resolved with the diameter D = 30 grid points.Apart from the simulation of turbulent flows, the entropy feedback has significantly improved stability of thermal flows with temperature gradients [14], multiphase flows [15] and other fluid dynamics problems.This gives us strong confidence that entropy-guidance can be extended to other low-dissipative physical systems.

ELBM: Questions and Answers
Here we answer some typical questions to the entropic lattice Boltzmann method (ELBM) reviewed above.For the sake of convenience, we rewrite the ELBM equation setting the time step h = 1: where β ∈ [0, 1], while the entropic mirror state f mirr is The stretch α is defined by the entropy balance between the pre-collision state f and the mirror state f mirr , A discrete-time H-theorem states: If the non-trivial solution α exists for the entropy balance (6), then the total entropy H(t) = ∑ x H( f (x, t)) is not increasing, H(t + 1) ≤ H(t).Note that the validity of the H-theorem requires not just the equilibrium to be evaluated through the minimization of H but also, and most importantly, the fulfillment of the entropy balance condition (6).

Is the entropic feedback in ELBM a stabilizing technique or a physically sound subgrid-scale model for turbulence?
A: The ELBM should be viewed as a built-in subgrid model rather than a mere stabilization technique.Stabilization in ELBM is a by-product of the discrete-time H-theorem.Instead of a mere addition of artificial viscosity, the ELBM allows the effective viscosity to fluctuate around the target value ν.In order to clarify this point, note a few general features of the entropic stretch α.
• Hydrodynamic limit: whenever the simulation is resolved (populations stay close to the local equilibrium), the stretch α tends to the fixed value α = 2 (and so does also the dual stretch, α * = 2, according to (7)).Then ELBM self-consistently becomes equivalent to the lattice Bhatnagar-Gross-Krook (LBGK) equation (α = 2) and recovers the Navier-Stokes equations with the kinematic viscosity, where c s is speed of sound (a O(1) lattice-dependent constant).
Note that the above is a direct implication of the built-in H-theorem.Indeed, the resolved simulation, at the kinetic level, is characterized by the fact that all populations are asymptotically close to the local equilibrium.Then, the entropy function becomes well represented by its second-order approximation: at fixed locally conserved fields (density and momentum here), if The levels of the entropy are then asymptotically close to the levels of the above quadratic form.It is under such condition that the entropy estimate ( 6) results in α = 2.Note that the standard Chapman-Enskog approximation is valid under precisely the same condition of closeness to the local equilibrium, thereby the viscosity ν is the same for both ELBM and LBGK.
• Effective viscosity and self-averaging: The effective viscosity in the above notation reads, Depending on the outcome for stretch α, the effective viscosity ν eff (αβ) is larger than the viscosity ν ≡ ν eff (2β) if α < 2, and it is smaller than ν if α > 2. In the first case, the (larger) effective viscosity leads to smoothing the velocity gradient at the given node, while in the second case, the smaller viscosity leads to a sharpening of the velocity gradient.Note that, when β → 1 (vanishing viscosity ν → 0), the effective viscosity (9) can drop to even negative values if α > 2. This asymmetry between the over-relaxation being "shorter" (α < 2) or "longer" (α > 2) than the LBGK over-relaxation α = 2 is the crucial implication of the compliance with the H-theorem: even if the effective viscosity becomes negative at some lattice nodes, this does not lead to numerical instability because even in that case the H-theorem (and the proper behavior of the total entropy) remains valid.
Parameterization with the effective viscosity ν eff (αβ) can be seen as an alternative to the parameterization with the over-relaxation α.Let us note that, if a pair {α, α * } is connected by the duality relation (7), then the mean value of the corresponding effective viscosity is equal to the viscosity (8), The relation ( 10) is termed self-averaging, and provides important albeit heuristic argument that the averaged-in-time effective viscosity in ELBM simulation is close to the viscosity ν.In other words, we expect that it is only the matter of resolution that the Entropy 2015, 17, 8099-8110 average effective viscosity deviates from ν.This assertion, while not rigorous, is supported by simulation (see Figure 3).The rapid fluctuations of the stretch α around α = 2 at a given monitoring point chosen at random in the simulation domain are clearly seen in Figure 4.In summary, the ELBGK exploits the self-adaptive mechanism of effective viscosity by choosing automatically the over-relaxation α at each node to guarantee the H-theorem at all sites and all discrete time-steps.When the grid is coarsened, over-relaxation α becomes "smeared" in an interval, [α min , α max ], with 1 < α min < 2, and α max > 2. The self-adapted over-relaxation set up by (6), results in two oppositely directed effects: if α < 2, the effective viscosity is larger than ν, and the ELBM will tend to smoothen any flow perturbation.On the other hand, if α > 2, the flow perturbation is enhanced (effective viscosity is smaller than ν).In ELBM simulations, these two effects act simultaneously on various nodes, with the net effect combining stabilization (through smoothing, α < 2) with the preservation of the resolution (through sharpening, α > 2).Note that, as β → 1, the effective viscosity can even drop to negative values when α > 2. This, however, does not lead to instabilities as the total entropy balance remains under control by the discrete-time H-theorem.This all is very different from a conventional perspective on "eddy viscosity" turbulence modeling, and it is not surprising that ELBM does not reduce to familiar large eddy simulation (LES) models [16].

What is the relation of ELBM to the entropic stabilizing techniques proposed in CFD?
A: During the last four decades numerous entropic stabilizing techniques have been proposed in computational fluid dynamics (CFD) (see, e.g., Refs [17][18][19][20][21][22][23][24][25] and references therein).The idea behind is, roughly speaking, to maintain an appropriate amount of artificial viscosity through the analysis of discretization of the entropy balance (physical or artificial).In this regards, ELBM is based on a different premise: it applies to strictly discrete systems (in velocity-space-time), and the discrete-time H-theorem does not reduce to the estimate of the entropy production (cf., e.g., [26]).

How ELBM performs in comparison to other stabilizing techniques proposed for LBM?
A: The closest analog of the conventional stabilization techniques in the LBM setting is perhaps the method of entropic limiters [27][28][29].The idea behind is to measure the closeness of the pre-collision state to the corresponding local equilibrium (in the sense of the entropy difference), and to apply equilibration instead of over-relaxation if the difference exceeds a user-defined threshold.This is similar to conventional artificial viscosity stabilization techniques in CFD.
4. What is the main numerical mechanism promoting stability in ELBM?
A: Stability is promoted by the discrete-time H-theorem.Note that the implication of the H-theorem in the presence of the over-relaxation allows post-collision distributions to be both closer to the equilibrium than the LBGK outcome (α < 2) or further away from the equilibrium (α > 2).
It must be noted that, in principle, for some pre-collision states, the corresponding entropic mirror state may not exist (and hence no entropy balance is possible).However, this happens beyond the domain of validity of the lattice Boltzmann models, and is of no concern in practice.
In particular, pathological cases (no solution for α) occur in none of the simulations referred to in this paper.A: ELBM is based on the discrete-time H-theorem which is imposed in a rather "orthodox" manner through the entropy balance condition (6) for the over-relaxation.A different realization of the entropic control was introduced recently by three of the present authors in [30] (we refer to this as KBC model).The idea is to replace the entropic over-relaxation on all the non-conserved moments as it is done in ELBM by a combination of the standard (unsupervised) over-relaxation of the stresses with the proper equilibration of the rest of the non-conserved moments.More specifically, if we write a moment representation of the populations, where k i is the contribution of locally conserved fields, s i are the stresses and h i are the remaining higher-order moments, then the mirror state for KBC models reads, where γ is the entropic stabilizer which is found by minimizing the entropy in the post-collision state (4) with the mirror state (11): The rationale behind is this: The over-relaxation of the stresses in the mirror state is the only formal condition to recover the viscosity ν (8); hence, an optimal post-collision state should minimize the entropy under this constraint.Thus, the KBC post-collision state is a quasi-equilibrium which corresponds to the minimum of the entropy function once all the relevant constraints are applied.Moreover, Equation ( 12) admits the following approximate solution, with X|Y = ∑ i,j X i [∂ 2 H/∂ f i ∂ f j ] eq Y j the entropic scalar product, and δs i = s i − s eq i , δh i = h i − h eq i .While (11) lumps together all the higher-order moments in the h-part of the populations, a generalization which makes a distinction within these moments is Entropy 2015, 17, 8099-8110 straightforward: For h i = ∑ m h mi with m labeling the different higher-order moments (or groups of such moments), we have instead of (11), while the formula (13) generalizes to with C −1 the inverse of the correlation matrix C mn = δh m |δh n .While the H-theorem is not directly imposed in the KBC models (unlike the ELBM), simulations of various setups demonstrated they are 'virtually indestructible' (Ref.[31]).
Note that in both ELBM and KBC models (and eventually in any lattice Boltzmann model) a statement that "it recovers the viscosity ν" refers only to a fully resolved simulation.Validity of the Navier-Stokes equation at small scales for a given simulation is checked independently, for example, by measuring the viscosity in the energy and enstrophy balance equations.For a detailed analysis of these aspects for the KBC models we refer to recent papers [32,33].

Conclusion
The second law of thermodynamics provides a parameter-free solution to the problem of controlling the effective turbulent viscosity, so as to tame numerical disruption.Besides turbulence, to which the above findings have an immediate impact on, the general notion of entropy-assisted computing, likely with different realizations of the entropy feedback, is expected to apply to other states of matter characterized by extremely low dissipation, such as superfluids [34] and cosmological fluids [35] near black-hole horizons.It is also of interest to explore whether a similar paradigm might inform the behavior of active matter systems [36].Finally, one may extrapolate even further and conjecture that entropy-assisted feedback systems, functioning according to the feedback loop discussed above, may be engineered outside the realm of fluid mechanics, typically at the intersection of information, biology and statistical physics [37][38][39].In an even broader perspective, we surmise that entropy-assisted procedures might also inspire the design of novel active feedback systems in natural, biological and possibly also medical sciences.

Figure 1 .Figure 2 .
Figure 1.Entropy-assisted computing.The initial state f is over-relaxed to the state f with the entropy function H value strictly below the value at the entropy mirror state f mirr .The zigzag trajectory of over-relaxations eventually ends up at the bottom of the well-at the equilibrium f eq .

Figure 3 .
Figure 3. Distribution of the time-averaged normalized effective viscosity R = (ν eff − ν)/ν at the mid-section of the flow past a round cylinder.Red/Yellow: Promoters (R < 0); Blue: Dampers (R > 0); Green: Nominal (R = 0).While the snapshot in Figure2demonstrates a larger variation of the effective viscosity, the time-average picture is much milder: most of the activity (strongest damping neighboring the strongest promotion) is concentrated at and around the twin shear layers, just behind the cylinder.In the rest of the domain, the deviation of the effective viscosity from its nominal value is less than a fraction of a percent.

Figure 4 .
Figure 4. Top: History of the entropic strech α at a monitoring point; Middle: Histogram of α; Bottom: Close-up of the hystogram around the dominant value α = 2.

5 .
Very recently, in Ref. [30], Karlin et al. presented a new entropic stabilizer for LB schemes.How is it different from the ELBM?