Lattice Boltzmann Modeling-based Design of a Membrane-free Liquid-liquid Microseparator

The benefits of continuous processing and the challenges related to the integration with efficient downstream units for end-to-end manufacturing have spurred the development of efficient miniaturized continuously-operated separators. Membrane-free microseparators with specifically positioned internal structures subjecting fluids to a capillary pressure gradient have been previously shown to enable efficient gas-liquid separation. Here we present initial studies on the model-based design of a liquid-liquid microseparator with pillars of various diameters between two plates. For the optimization of in silico separator performance, mesoscopic lattice-Boltzmann modeling was used. Simulation results at various conditions revealed the possibility to improve the separation of two liquids by changing the geometrical characteristics of the microseparator.


Introduction
Introduction of microflow devices in chemical processing offers not only intensification of transport phenomena and more controlled process conditions, but also more efficient approaches for a transfer from lab-scale to industrial production, thus resulting in paradigm shift in chemical engineering. The quest for continuous and lately even end-to-end production in pharma industry has revealed new challenges as well as the interest in implementing micro-and meso-flow systems, particularly regarding the reactor integration with the downstream processing units. 1,2 In the (bio)chemical processes, removal of gaseous (by)products from the liquid phase or the separation of two-liquid phase systems after extraction is often required. In the past decade, microflow-based phase separators have increasingly gained attention due to several advantages comprising high surface-to-volume ratio, very efficient mass and heat transport, and the possibility of integrating micro-and even nano-structures. 3 Membrane-based devices are typically a matter of choice to achieve separation of phases with different hydrophobicity. [4][5][6][7] Another approach is to have a capillaries-perforated-wall inside a microchannel. 8 However, both solutions are prone to clogging, which might be prevented by the use of microsettlers. 9 One way to achieve continuous phase separation without using either a membrane or external forces is by subjecting the fluids flowing through the separator to a capillary pressure gradient. 10 Accordingly, both wetting and non-wetting phases are forced to different parts of the separator so that end-users merely have to introduce two separate outlets in those parts. The original design of such a device was used for separating a gas-liquid mixture, 10 although a similar setup has been studied for liquid-liquid separation. 11 In its essence, the design is fairly simple as a "forest" of different-sized pillars is introduced into a channel "between-twoplates". These pillars are equally spaced with respect to their centers while their diameters decrease in the direction perpendicular to the flow thus creating a capillary pressure gradient.
The lattice Boltzmann (LB) method is a mesoscopic approach to solving computational fluid dynamics (CFD) problems. It has numerous additional applications such as modeling heat and mass transfer as well as multiphase systems. [12][13][14] Lattice-Boltzmann (LB) simulation is different from the conventional Navier-Stokes-equations-solving procedure, as it does not treat the fluids as a continuum. Conversely, LB treats any medium (whether fluid or solid) as comprised by discrete particle populations of the same substance, occupying a fictitious discrete space (namely, a lattice structure) with statistically the same momentum. At the micro scale, it is important to consider the use of LB since the continuum assumption of Navier-Stokes equations may not be justified at this scale.
The multiphase LB algorithms have been around for almost as long as the method itself, and the first such algorithm was the so-called "color-gradient" (CG) approach. 15 Although newer algorithms have been developed (with prospective advantages), 13 it is still quite favorable to use CG as it is fairly simple to implement while it offers a rather seamless control over the simulation parameters, especially in liquid-liquid applications.
This communication presents results from ongoing research, which uses CFD to improve the performance of a membrane-free liquid-liquid microseparator. This work aims at showing how modeling-based design can, in principle, be applied to test and improve existing equipment as an attempt to reduce development costs by cutting off material expenses.

Color-gradient lattice Boltzmann method
The lattice Boltzmann (LB) method follows the movement of discrete particle populations, where those particles refer to the same substance and statistically share the same momenta. These groups of particles are mathematically represented by distribution functions modeling, at time t and about position x  . In computations, each substance (or phase) s is represented by a separate set of functions f. In doing simulations with the LB method in 2D cases the computational domain is usually broken down into a square lattice, where at each node f it is discretized into a finite amount of lattice directions a. In this study, a 3D LB model was used with a cubic lattice with 19 lattice directions. The multisubstance LB equation at spatial coordinate x  can thus be written as follows: 15 where e  is the basic lattice velocity, t is the discrete time, and Ω is the collision operator. The left-hand side of Eq. 1 represents the streaming, in which after each time step Δt particle populations represented by f move to the neighboring lattice nodes according to e  . The 19-velocity LB model is commonly referred to as D3Q19 (Fig. 1), where ( ) , , x y z e e e e =     has the following values: 0,1, 1,0,0,0,0,1, 1,1, 1,0,0,1, 1,1, 1,0,0 0,0,0,1, 1,0,0,1, 1,0,0,1, 1, 1,1,0,0,1, 1 . 0,0,0,0,0,1, 1,0,0,1, 1,1, 1,0,0, 1,1, 1,1 In the color-gradient multiphase LB model, Ω actually combines three suboperators: 17 is the perturbation operator and is the source of interfacial tension, whereas ( ) ( ) 3 Ω represents the "recoloring" step, which keeps the phases pure, and the interfaces sharp. The BGK collision operator is sufficiently accurate for this study as microfluidic devices operate at low Reynolds numbers. 16 It takes the following form: where τ * is the effective relaxation parameter, and f eq is the equilibrium distribution function. In the color-gradient model, f eq is defined as: 17 ρ s is the local macroscopic density of the substance s, and T ( , , ) is the local macroscopic velocity, which are respectively computed as: where ρ=ρ 1 +ρ 2 . f in Eq. 5 is related to an arbitrarily chosen parameter α according to the streaming direction in the lattice (with a=0 corresponding to null stream-  To prevent unphysical behavior, 0<α s <1 should be satisfied. Through α one can control the density ratio of two fluids in binary systems: 18 w in Eq. 5 is lattice-based weighting factor, which in case of D3Q19 takes the following values: τ * from Eq. 4 controls the kinematic viscosity v of the fluid through: By depending itself on both species densities (via phase field ρ s ), τ * is obtained by combining the relaxation parameters of the two species: 18 where: Assuming the same A for both species, this parameter controls the interfacial tension σ in the binary system through the following expression: 17 4 f * is the distribution function post BGK collision and perturbation, 0≤β≤1 controls the interface thickness, whereas f eq = f eq,1 + f eq,2 . θ is the angle between and e  , whose cosine can be computed with the help of the dot product:

Boundary conditions
At solid walls, the bounce-back rule was used to give the no-slip boundary condition. 19 The wetting of the two phases was controlled by the fictitious-density boundary, 20 and was set to completely favor the wetting of only one of the phases, i.e., one phase was completely wetting, while the other was completely non-wetting. The inlet and outlet boundary conditions used in the computations were the non-equilibrium extrapolation density and velocity boundaries. 21 To create segregated flow at the inlet, a constant density boundary was used there. The density would be varied with a step function controlling which substance was let in at the time, thus creating a segregated flow. At the outlet, a constant velocity boundary was imposed, which controlled the total flow rate in the system. To simplify the computations, especially in terms of pressure regulation, only a single outlet was present. In the space leading to this outlet, a triangular obstacle was placed, thus creating two paths for fluids to flow. Either north or south side of this triangle represented the two outlets that would have been present in a physical microseparator. This is depicted in Fig. 2.

Simulation setup
A sketch of the geometry was drawn in Inkscape, and exported to a PNG file. The resulting PNG file was then processed with a Python script to return a 2D wall-mapping array, which could then be loaded into our in-house developed LB solver. The latter is written in CUDA C++ programming language, and it incorporates the model described above. The solver takes an input geometry file (created by the Python script), and produces VTK output files for later analysis in ParaView. Parameters, other than the geometry, are set in the code itself. The simulation parameters were the following:

Results
In this study, a system similar to the one proposed by Wheeler et al. 10 was studied by utilizing the LB method. However, initially, a physical mi-croseparator was laser-cut into a polycarbonate plate and tested (Fig. 3 should serve as a visual aid). The microseparator comprised an inlet at one side, an outlet at the opposite side, and the pillar "forest" in between. Preliminary testing of the prototype microseparator with selected biphasic system at different process conditions regarding total fluid flow rate (0.36 -3.6 mL min -1 ) and its ratio (2:1, 1:1, 1:2) did not reveal satisfactory results as the highest phase purity achieved at either exit was ~90 %. The modeling-based optimization of the device was therefore selected as the next step.
By consoling with the experimentalists, a somewhat changed geometry was drawn-up and tested using CFD. Corresponding numerical results showed that the proposed separator was yet unable to separate the liquids (Fig. 4).
These results showed that there were two major issues with this geometry. The first issue was that the non-wetting fluid had not been forced far enough to its outlet, and consequently, it exited through both outlets. The second issue was that, before the two outlets, the pillars were discontinued, which actually gave favorable conditions for either fluid to exit at either outlet. The former of these issues could be solved by moving the inlet further to the side, where the non-wetting fluid is expected to exit the separator (where pillar diameters are smaller), while the latter problem could be solved by expanding the pillar "forest" further down the channel.
This newly suggested geometry was drawn up and further CFD computations were performed, also at higher Re (Fig. 5, Fig. 6, and Fig. 7).
Indeed, this geometry performed better in CFD testing than the previous concept. It was able to successfully separate the two phases for three different flow ratios (wetting : non-wetting = 2 : 1 ; 1 : 1 ; 1 : 2). However, at the flow ratio 2:1, the non-wetting phase formed long slugs in the separator, which would not always coalesce, thus leading the wetting phase into the wrong outlet. The flow rate ratio has this impact because smaller pillars have greater gaps between them than the larger ones. While this has the desired effect of creating a capillary pres-  sure gradient, it also results in greater void volume around the smaller pillars, which requires more fluid to fill the area. This contributes to the flow rate ratio limitations of the device. As mentioned previously, the issue needs to be studied in more detail, further CFD simulations at broader range of conditions will be conducted in the ongoing research. This newly proposed geometry is currently awaiting laboratory testing.

Conclusion
The LB method was successfully applied to computationally model a multiphase-flow separation problem. Through LB computations, a new geometry for an existing membrane-free separator was proposed. The new geometry proved to be able to separate two immiscible fluids better than the original geometry. An issue that remains unsolved is that separation apparatus seems to be effective only at certain flow ratios of the two phases. As an attempt to solve the flow-ratio issue, future work will include experimental testing of the newly proposed geometry (and possibly other geometries).

ACKnowlEDGmEnTS
The financial support of the Slovenian Research Agency through PhD Grant MR-39080 (FS) and Grant P2-0191 as well as to projects J7-1816 and N2-0067 is acknowledged. The authors also acknowledge the support of EC H2020 project COM-PETE, Grant 811040. The authors also thank Mark Selan for his contributions to the study.