Binary mixture flow with free energy lattice Boltzmann methods

We use free energy lattice Boltzmann methods (FRE LBM) to simulate shear and extensional flow of a binary mixture in two and three dimensions. To this end, two classical configurations are digitally twinned, namely a parallel-band device for binary shear flow and a four-roller apparatus for binary extensional flow. The FRE LBM and the test cases are implemented in the open-source C++ framework OpenLB and evaluated for several non-dimensional numbers. Characteristic deformations are captured, where breakup mechanisms occur for critical capillary regimes. Though the known mass leakage for small droplet-domain ratios is observed, suitable mesh sizes show good agreement to analytical predictions and reference results.


I
Mixture flows are omnipresent in nature and essential to many industrial processes. Taylor [37] proposed machinery to examine the deformations of droplets induced by shear and extensional flow of binary multicomponent mixtures. The deformation is governed by the balance of outer forces and surface tension. Once this force balance is in favor of deformation, the droplet will break. Modifying the properties of the system, the breakup process can be adjusted in terms of number and size of resulting droplets. These phenomena are essential in manufacturing processes, for example in order to maximize the efficiency of creating emulsions [5,39]. For the computer simulation of binary mixture flow several methods have been employed in the past. Due to the intrinsic parallelizability which enables the outsourcing of high performance computing (HPC) machinery, the lattice Boltzmann method (LBM) emerged as an unconventional alternative for multicomponent computational fluid dynamics (CFD). The popularity of LBM for CFD and beyond has increased significantly [22]. Several data structures are available commercially and opensource. Exemplarily for the latter, the highly parallel C++ framework OpenLB [17] has been successfully used for simulations of various transport processes on Top500 HPC machines (e.g. [17,9,28,25,6,34,30,8,33,29,31,4]). Simulating multiphase and multicomponent flows in LBM is mostly based on a phase field model with diffuse interfaces. The interfacial zone brings forth additional physics captured by the Cahn-Hilliard equation (CHE), though in turn upholds the high parallelizability of LBM. Several approaches for the underlying mixture dynamics exist [13], for example the free energy model (FRE) [36,26]. Tunable physical effects and top-down configuration of thermodynamics, are advantages of models akin to FRE. Albeit a high potential for numerical simulations, applications with FRE LBM for flows relevant to industrial processes are still rare.
The dynamic effects on an immiscible binary component mixture can be abstracted into shear-and extension-dominated flows. For these two types of dynamic mixture flows, the present work implements and tests the FRE LBM with a simplistic binary fluid composition (equal density and viscosity). In particular, deformation as well as breakup phenomena are distinctively assessed to determine the models applicability for more complex applications. As such, we approve the suitability of the presented FRE LBM for numerically simulating these types of binary mixture flows via digitally twinning classical devices and comparing the results to references.
This manuscript is structured as follows. Section 2 summarizes the methods, the numerical results are described in Section 3 and Section 4 draws conclusions and closes the paper.
where : Ω × → R denotes the density, : Ω × → R is the fluid velocity, and > 0 is a kinematic viscosity, and Ω ⊆ R and ⊆ R + denote space and time, respectively. Here, P th = P chem + I : Ω × → R × is the thermodynamic pressure tensor. For single phase and single component flow, P th reduces to the isotropic pressure I . In case of a multicomponent flow the corresponding thermodynamics are introduced by the anisotropic chemical pressure tensor P chem [14].
To model a multicomponent system, capturing additional physics of the diffuse interface, the Cahn-Hilliard equation (CHE) is introduced where is the order parameter, denotes the chemical potential and is related to the mobility of the interface. Complemented with respective initial and boundary conditions, (1) Discrete velocity sets.

Free energy lattice Boltzmann method.
We assume a classical discretization of the phase space and the time domain [19] such that the following derivation is done completely in lattice units The present LBM is based on the lattice Boltzmann equation (LBE) where = 0, 1, . . . , −1 counts the discrete velocities in and ( , ) = ( ( , )) T denotes the populations in discrete space-time ( , ) ∈ (Ω ℎ , ℎ ) and maps to R . The here used velocity stencils are given in Figure 1. The classical Bhatnagar-Gross-Krook (BGK) collision model [2] reads where > 0.5 denotes the relaxation time at which relaxes towards the equilibrium with denoting the lattice weights and s is the lattice speed of sound. The term ( , ) obeys Guo's forcing scheme [7] ( , where is a force field. The macroscopic flow variables are recovered via the discrete moments of the populations respectively. As such, we use (2.4) to approximate the conservative variables in the NSE (2.1) and (2.2). Formal Chapman-Enskog expansions are given in [19] and limit consistency is proven force-free in [32]. Both establish the relation The appropriate body force to recover the chemical pressure is defined as the residual is specified below. For coupling the approximations of (2.1), (2.2) and (2.3) through the force (2.11), the thermodynamics can be consistently derived for a binary fluid mixture via the free energy functional [15,26] where 1 , 2 , and are tunable parameters for the interface tension, and arguments are neglected where unambiguous. The auxiliary variables and are defined as respectively, where 1 and 2 are the concentration fractions of the respective components. The chemical potentials and are defined as functional derivatives of the free energy 16) respectively. The free energy of the system is minimized at equilibrium through the thermodynamic force induced by the chemical pressure. In case of a planar interface, the minimization of the free energy yields a simplified interface solution where is the interface width and the bulk components are identified by = ±1 at = ±∞ [21]. For approximating the CHE (2.3), a second population ( , ) is required such that its zeroth moment yields the order parameter This population evolves according to a second LBE where > 0.5. To recover the coupled CHE (2.3) in the continuum limit [19,26], the corresponding equilibrium reads ( + , + 1) = ★ ( , ) (2.24) through the space-time cylinder. The computing steps implemented in one collision are summarized in Algorithm 1 and the streaming is realized as a mere pointer shift [20]. The present work illustrates only the bulk solver.
Standard LBM boundary methods are applied to impose velocity boundary conditions for the binary mixture. The macroscopic initial conditions are implemented via initialization of the populations to the corresponding equilibrium and alignment of kinetic moments through collisions [24] preceding the actual simulation time horizon.

N
To assess the capability of the FRE LBM for recovering shear and extensional binary flows, we emulate Taylor's parallel-band and four-roller devices [37] by means of numerical simulations in two dimensions (2D). The former is also tested for three dimensions (3D). The geometric setup of the 2D simulations is sketched in Figure 2.
All computations are done with OpenLB release 1.4 [18] on several HPC machines, either using up to 16 nodes with five quad-core Intel Xeon E5-2609v2 cores each, or a maximum of 75 nodes with respectively two deca-core Intel Xeon E5-2660v3.
The deformation of the 1 droplet where is the longer axis halved and the shorter one, is measured via intrinsic functors of OpenLB [17]. In case of a horizontally measured inclination angle = 0 • , an interpolation along the space directions recovers the location of the interface. If a simultaneous deformation and inclination of the droplet occurs, and are approximated through concentric circles and is computed at the intersection point.
Though essential differences between 2D and 3D deformations are known [35], certain non-dimensional regimes still allow a side-by-side comparison. Based on that we compare the FRE LBM solutions to 3D reference computations [16,23] and analytical predictions [27,37] (2) Geometric setup of numerical test cases for binary flow in two dimensions. Scales differ for the purpose of representation. ℎ = , (3.5) where , , , , , and are shear rate, droplet radius, viscosity, surface tension, interface thickness and a mobility parameter, respectively. The ratios of viscosity and density of the components are unity. The bounds between these regimes however are not sharp, such that the droplet may pass through multiple types during the breakup process. Due to differences between 2D and 3D droplet deformations, c in 2D is significantly larger than in 3D. The latter agrees well with the literature [16], and so does the breakup scenario (see Figure 5). In 2D for = 5 c at = 1, = 0.2 and = 40 we observe end-pinching as well as a capillary wave breakup (see Figure 4).

Steady state validation.
The droplet radius is set to = 20, which corresponds to a ratio of 40 between domain length and radius. Figure 6 summarizes the deformation in the subcritical capillary regime. For = 0.01, 0.02, 0.04 the droplet shows little to no deformation. Beginning at = 0.05 the deformation becomes significant and increases rapidly with increasing and with a considerably faster rate than in the shear flow. Our simulation results and the reference data from [11,12] agree from the perspective of an overall trend but differ at individual values. Based on the same reasoning as above, we conclude that a 3D extensional flow produces a higher deformation at lower than in 2D. Components 1 (red) and 2 (blue) are plotted at normalized time steps.

Breakup.
reproduction of these phenomena, the results show heavy satellite shrinkage due to the high domain-droplet ratio [41].

C
We set up a FRE LBM algorithm implementation in OpenLB and test it for shear and extensional binary mixture flow in two and three dimensions. Taylor's parallel-band and four-roller devices are digitally twinned in a simplified manner and used for validation of the numerical scheme. To the authors' knowledge, the present work is the first application of LBM for simulating a four-roller apparatus. Characteristic deformations for steady states and breakup scenarios in critical capillary regimes are captured. Though the known satellite loss for very small droplet-domain ratios is observed, with suitably fine meshes we find good agreement to references.