*3.2. Computational Domain of Droplet Breakup*

The schematic diagram of computational domain and a droplet passing through metal foam is shown in Figure 3, which corresponds to the experimental measurement. A metal foam of 60 PPI was inserted into a rectangular channel of 50 (width) × 50 (height) lattice cells. The physical unit of one lattice is 73 μm. The length of the channel (*L*) is 500 lattice cells, and the distance between the droplet and metal foam (*L1*) is 175 lattice cells. The 3D matrix of the metal foam was restructured by stacking the 2D cross-section images obtained by micro-CT. It should be noted that this reconstructed method perfectly matches the quadrilateral grids of the LBM and has a high computational efficiency, although its accuracy is less than the curve surface mesh. The flow is assumed to be steady, incompressible, and laminar. The rectangular channel is filled with silicon oil, and the water droplet is initially placed at the channel inlet. Constant velocity is imposed on the silicon oil at the inlet boundary of the channel to drive the droplet forward. The half-way bounce back rule is utilized for the solid wall [39]. The velocity inlet boundary condition with uniform velocity and the convective boundary condition are applied at the inlet and outlet of the channel, respectively [40]. It is noted that the conventional SC model suffers several limitations such as the numerical instability at a high-density ratio [41]. In this study, a density ratio of 1.0 is taken for the silicon oil and water system. Such a density ratio of 1.0 is present for the multi-phase related to the micro-reactors and microfluidic devices and is easy to simulate using the SC model. This study does not touch the cases with high-density ratios.

**Figure 3.** Schematic diagram of a droplet passing through metal foam in a computational domain of 500 (length) × 50 (width) × 50 (height) lattice cell and L1 of 175 lattice cell.

#### **4. Discussion**

#### *4.1. Dimensionless Parameters*

The physical values are scaled using several dimensionless parameters. One parameter is the capillary number (*Ca*) which describes the relative importance of viscosity and the interfacial tension. Here the *Ca* number is defined as Equation (14) [21]

$$\mathbf{C}\_{a} = \frac{\boldsymbol{\mu}}{\sigma} \tag{14}$$

where *u*, η, and σ are the average inlet velocity of the silicon oil phase, the dynamic viscosity of the silicon oil, and the interfacial tension between the silicon oil and the water.

The second dimensionless parameter *Dd (Dd* = *dh*/*dp*) describes the relative size between the pores of the metal foam and the droplet where *dp* is the average pore diameter of the metal foam, and *dh* is the droplet diameter at the channel inlet. The third parameter, *Ld*, is defined as the ratio of the thickness of the metal foam (*Lh*) to *dp*. The fourth parameter *Db* is defined as *Db* = *D*/*dh* where *D* is the distance between the two successive droplets when the first one reaches the position of *L*1 (Figure 3).

In this study, the change of droplet superficial area (Δ*S*) is defined as

$$
\Delta S = \frac{S\_t}{S\_i} \tag{15}
$$

where *si* is the droplet's original superficial area before the droplet gets into the metal foam; *st* is the superficial area of the droplet which has entered the metal foam. *st* is obtained by averaging the corresponding values of three simulations with different metal foams (the PPI and porosity are the same). In this way, the generality of the simulation is promoted. Theoretically, when Δ*S* exceeds 1.0, the droplet definitely deforms or breaks up.

#### *4.2. Modeling Coe*ffi*cients and Validation of LBM Simulation*

The LBM with the SC model requires identifying a correct *Gw*,*<sup>k</sup>* before simulating the flow pattern of two immiscible fluids. The water contact angle surrounded by the silicon oil on the surface of a pure nickel material surface was measured using a Drop Shape Analyzer (KRUSS, Hamburg, Germany; see Figure S1a). The relation between the static contact angle and *Gw*,*<sup>k</sup>* is investigated in independent LBM simulations by equilibrating a water droplet surrounded by a flat solid surface. The value of the interaction strength between the silicon oil and water (*G*o,w) is fixed and *Gw,k* is altered to examine the wettability of two fluids. As shown in Figure 4b, when *Gw,k* is ±0.06, the simulated contact angle is 135◦, which is in good agreement with the experimental value of 132. Therefore, in the following simulation, *Gw,k* is chosen as ±0.06. It is known that *G*w,k is related to the grid size in the SC model. The physical unit size of one lattice in this study was set as 73 μm, which balanced the simulation accuracy and the calculation time, considering that a high spatial resolution is favorable for improving the simulation accuracy but meanwhile significantly increases the calculation time especially for a 3D simulation. In this study, the accuracy of the grid size value is verified by comparing the simulated and experimental results regarding the contact angle (see Figure S1) and the droplet shape when the droplet leaves the metal foam (see Figure 4).

**Figure 4.** Comparison of the lattice Boltzmann method (LBM) simulation and the experimental results. (**a**) The experimental and (**b**) the simulation results for *Ca* = 0.063, *Dd* = 1.06, and *Ld* = 4.88. (**c**) The experimental and (**d**) the simulation results for *Ca* = 0.084, *Dd* = 1.31, and *Ld* = 2.44.

The experimental results using a Plexiglass tube (Figure 1) were used for the validation of the LBM simulated ones. The droplet deformation and breakup after the droplet left the metal foam were recorded by high-speed video. Correspondingly, the LBM code was used to simulate the hydrodynamics of the droplet. Quantitative comparisons between the simulation and the experiments were performed. When *Ca* = 0.063, *Dd* = 1.06, and *Ld* = 4.88, as shown in Figure 4a, the experimental result showed that the droplet broke into two droplets with the diameters of 0.7 and 1.13 mm, respectively, after it left the foam. The same phenomenon can be observed in the simulation result (see Figure 4b) with the two droplets' diameters of 0.7 and 1.04 mm, respectively. When *Ca* = 0.084, *Dd* = 1.31, and *Ld* = 2.44, only one droplet with the diameter of 1.35 mm was detected after the droplet passed through the metal foam (see Figure 4c). The corresponding simulation result exhibits the same tendency with the droplet's diameter of 1.30 mm (see Figure 4d). Based on the validations above, conclusions can be drawn that the LBM model and code are adequate to simulate the hydrodynamics of droplets passing through metal foam.

In the numerical study, we considered the case of equiviscous droplets with medium viscosity, i.e., ν*<sup>w</sup>* = ν*<sup>s</sup>* (dynamic viscosity of water and silicon oil), although the experimental setup shows ν*w*/ν<sup>s</sup> = 1/48. It is known that a large viscosity ratio likely results in an instable LBM simulation. To solve the instability, the multi relaxation time (MRT) approach has been reported to be feasible [42]. However, the MRT method requires a large computer resource that is most likely inapplicable for the 3D LBM simulation on the hydrodynamics of one droplet passing a porous structure like metal foam. Here, despite that a viscosity ratio should be the real value for a more reliable comparison, the droplet dynamics are somehow not so significantly different. From the numerical viewpoint, furthermore, the equiviscous case has a great numerical stability. That is, the efficient calculation compensates for any loss of accuracy.

#### *4.3. The E*ff*ect of Ca*

For the metal foam reactor with immiscible fluids, smaller droplets produced in the reactor exhibits higher activity. For example, in our previous studies on biodiesel synthesis by transesterification with methanol, methyl ester yield is strongly dependent on droplet size [4]. This is because the mass transfer of triglycerides (TG) from the oil phase to the methanol/oil interface limits the rate of methanolysis reaction and controls the kinetics at the beginning of the reaction [43,44]. The overall volumetric TG mass transfer coefficient increase is attributed to an increase in the specific interfacial area by the decrease in droplet size, which leads to an increase in the TG reaction rate. For the metal foam with two immiscible fluids, when the continuous phase drives the droplets to pass through the porous material, the solid walls stretch the droplets, which deform into an elongated shape and break up. The deformation and breakup increase the superficial area of the droplets. When the distance between two droplets is smaller than one critical valve, they coalesce [45], resulting in a decrease in the superficial area.

Figure 5 shows the effect of *Ca* on droplet dynamics when the droplet passes through a metal foam. For the *Ca* values of 0.031 and 0.046, as shown in Figure 5b,c, the droplet is stretched when it collides with the walls of the metal foam. Then, the droplet shape restores at the pores of the metal foam. When the droplet leaves the metal foam, an obvious elongation can be observed. It should be noted that no droplet breakup is detectable for the *Ca* values of 0.031 and 0.046. Several previous studies [46–49] have indicated, when the viscosity ratio of the droplet to the matrix fluid λ is less than 4, there is a "critical capillary number" *Cac*, above which the droplet continues to deform and finally breaks in the creeping flow regime. The critical capillary number for droplet breakup in shear flow is the lowest for λ, roughly around 0.647, and its value (Ca*<sup>c</sup>* ≈ 0.4) is slightly less than the case for λ = 1, where *Cac* ≈ 0.41 [50]. The droplet hydrodynamics passing through obstructions in confined microchannels were explored by Chung et al. both numerically and experimentally [8]. They found that for the cylinder obstruction, the *Cac* was around 0.1. In our case, as shown in Figure 5d, for a *Ca* of 0.061, the droplet breaks up and forms one daughter droplet. The *Cac* in our case is around 0.061, which is lower than that of 0.1 for the case of one droplet passing through one cylinder obstruction. It was found that the stretching rate in shear flows may be increased by incorporating periodic reorientations

in the flow [51]. Wen et al. [4] also reported that droplet size decreases with increasing reorientations by turns in a zigzag micro-channel for the mixing of ethanol and soybean oil. The smaller *Ca*c for the case of metal foam than that for the cylinder obstruction in a micro-channel is mostly likely attributed to the frequent periodic reorientation in the flow. When *Ca* further increases up to 0.092, three daughter droplets can be observed, as shown in Figure 5f.

**Figure 5.** Hydrodynamics of one droplet passing a metal foam with different *Ca* values (*Ld* = 7.32, *Dd* = 1.82): (**a**) schematic diagram, (**b**) *Ca* = 0.031, (**c**) *Ca* = 0.046, (**d**) *Ca* = 0.061, (**e**) *Ca* = 0.078, and (**f**) *Ca* = 0.092.

As shown in Figure 6, when *Ca* increases from 0.031 to 0.092, Δ*S* rises nonlinearly from 1.0 to 1.2 at the channel length *L* = 18 mm where the droplet has left the metal foam, and its shape evolves to be stable. The rise in Δ*S* at *L* = 18 mm is obvious, with the increase in *Ca* from 0.031 to 0.061. For *Ca* = 0.031, the Δ*S* of 1.0 at *L* = 18 mm suggests that no daughter droplet is formed when the droplet leaves the porous structure. The high *Ca* value strengths the droplet breakup and restrains the coalescence.

**Figure 6.** Change of the superficial area of the droplet passing through metal foam with different *Ca* values (*Ld* = 7.32; *Dd* = 1.82).

#### *4.4. The Impact of Dd*

PPI and porosity determine the cell diameter of metal foam. The cell average diameter of the metal foam with 60 PPI and the porosity of 95% was calculated to be 0.817 mm. The effect of *Dd* on droplet dynamics when the droplet passes through a metal foam is show in Figure 7. For the *Dd* values of 0.89 and 1.07, as shown in Figure 7b,c,the droplet passes though the metal foam without splitting into daughter droplets. For *Dd* = 0.89, no obvious droplet deformation can be observed. This is because the droplet is small enough for it to pass through the metal foam without a strong interaction with the three-dimensional porous structure of the metal foam. For *Dd* = 1.07, the droplet is stretched into a linear shape. For *Dd* ≥ 1.25(see Figure 7d–g), a polliwog-shaped droplet is formed, which is the characteristic of the case that an obstacle intrudes into a moving droplet followed by the breakup of droplet. The formation of the tail of the polliwog-shape is attributed to the interfacial area between the droplet and the obstacle. The polliwog-shaped droplets are only present at the outlet of the metal foam for the *Dd* values of 1.82 and 2.24. In contrast, for 1.25 ≤ *Dd* ≤ 1.33, the polliwog-shaped droplets are formed inside the metal foam. This difference is because the polliwog-shaped droplet requires enough space to be stretched until it breaks up. For the big droplet of *Dd* = 1.82 or 2.24, the tail is so long that the required space for stretching the droplet to break up exceeds the pore size of metal foam, resulting in the formation of polliwog-shaped droplets at the metal foam outlet.

**Figure 7.** Hydrodynamics of one droplet passing a metal foam with different *Dd* values (*Ld* = 7.32, *Ca* = 0.061): (**a**) schematic diagram, (**b**) *Dd* = 0.89, (**c**) *Dd* = 1.07, (**d**) *Dd* = 1.25, (**e**) *Dd* = 1.33, (**f**) *Dd* = 1.82, and (**g**) *Dd* = 2.24.

.

To further elucidate the details of droplet breakup evolution in the metal foam, the case of *Dd* = 1.25 is taken as a representative one (see Figure S2). Droplet splitting in a microfluidic channel by the use of an obstacle was investigated numerically by Lee and Son [9]. The droplet is elongated around a cuboid obstacle and its portion becomes narrow near the front of the obstacle. The elongated portion breaks off at the front corner of the obstacle and then at the rear corner of the obstacle. Chung et al. [8] found that for cylinder obstruction, the thread becomes uniformly thinner around the cylinder, and finally breaks up at the front parts, independent of Ca. They concluded that the thread breakup was attributed to the velocity gradients induced by the geometric effect of the obstructions. Different from the hydrodynamics of droplets passing an obstacle in a micro-channel, for the droplet passing through metal foam, a waist is formed between the head and tail (see Supplementary Figure S2a). The droplet

continually lengthens and breaks into two droplets, which are then driven into a spherical shape by surface tension (see Figure S2b,c). That is, the breakup of the droplet occurs at the waist not on the surfaces of the porous obstacles. In this regard, this droplet breakup is to some extent similar to that in a simple shear flow [8].

Figure 8 shows the effect of *Dd* on Δ*S*. When *Dd* increases from 0.89 to 1.33, the Δ*S* at *L* = 18 mm rises from 1.0 to 1.2. With a further increase in *Dd* from 1.33 to 1.82, *Dd* declines to 1.15 at *L* = 18 mm. For *Dd* = 0.89, the Δ*S* of 1.0 at *L* = 18 mm suggests that no daughter droplet is formed when the droplet leaves the porous structure.

**Figure 8.** Change of the superficial area of the droplet passing through metal foam with different *Dd* values (*Ld* = 7.32; *Ca* = 0.061).

#### *4.5. The Impact of Ld*

The effect of *Ld* on droplet hydrodynamics is shown in Figure 9. For the case of *Ca* = 0.061, *Dd* = 1.25, and *Ld* = 2.44 (see Figure 9a), the polliwog–shaped droplet is formed at the metal foam outlet and then it breaks up into two droplets. With the increase in *Ld* from 2.44 to 9.76, the droplet breaks up inside the metal foam and the number of formed droplets is fixed at two. Therefore, no difference in the number of the generated daughter droplets is observed regardless of various *Ld*. As mentioned in Section 4.4, when *Dd* < 1.25, the droplet passes though the metal foam without splitting into daughter droplets. For *Dd* = 1.25, the *Dd* values of the two broken droplets are both less than 1.25, thus no splitting of the two droplets takes place in the following metal foam part. For the cases of *Ca* = 0.061, *Dd* = 1.82, and *Ld* in the range of 2.44 to 9.76, as shown in Figure 9e–h, the polliwog–shaped droplet is always formed at the metal foam outlet and then it breaks up into two droplets. *Ld* shows a negligible effect on the number of the generated daughter droplets in this case. In contrast, for a microchannel with obstructions, Chung et al. reported that the number of satellite droplets equals to the number of cylinder obstructions when the penetration of droplet fluid into the cylinder interval occurs [8]. The number of obstructions plays an important role in the number of generated daughter droplets. The difference in the droplet hydrodynamics is mainly because for metal foam, the breakup of the droplet occurs at the waist of the deformed droplet, whereas the breakup takes place on the surfaces of the obstructions in the case of a microchannel with obstructions.

**Figure 9.** Hydrodynamics of one droplet passing a metal foam with different *Ld* values (*Ca* = 0.061): (**a**) *Ld* = 2.44, *Dd* = 1.25; (**b**) *Ld* = 4.88, *Dd* = 1.25; (**c**) *Ld* = 7.32, *Dd* = 1.25; (**d**) *Ld* = 9.76, *Dd* = 1.25; (**e**) *Ld* = 2.44, *Dd* = 1.82; (**f**) *Ld* = 4.88, *Dd* = 1.82; (**g**) *Ld* = 7.32, *Dd* = 1.82; and (**h**) *Ld* = 9.76, *Dd* = 1.82.

#### *4.6. The Impact of Db*

The above discussions are related to the hydrodynamics of a single droplet passing through metal foam. Actually, metal foam is usually used for the mixing of two immiscible liquids where a series of droplets are formed [4]. Therefore, in this section, the hydrodynamics of two successive droplets passing through metal foam is investigated. The effect of *Db* on droplet hydrodynamics is illustrated in Figure 10. As shown in Figure 10, when *Db* is below 3.00, the Δ*S* less than 1 can be found at the inlet of the metal foam. A *Db* of 3.00 is the critical value to avoid the coalescence of two droplets at the inlet zone of the metal foam in the case of *Ca* of 0.061, *Dd* of 1.33, and *Ld* of 7.32. To further elucidate the details of two successive droplets hydrodynamics at the inlet zone of the metal foam, the case of *Db* of 2.00 is taken as the representative one. A *Db* of 2.00 at the inlet of the rectangular channel decreases to 0.35 at the position of 16 mm (1 mm before the inlet of the metal foam). Then, the two droplets merge into a big droplet at the inlet zone of the metal foam (see Figure S3).

**Figure 10.** Change of the superficial area of two successive droplets when they pass through a metal foam with different *Db* values (*Ld* = 7.32, *Ca* = 0.061, and *Dd* = 1.33).

#### **5. Conclusions**

We developed a LBM of droplet deformation and breakup and performed validation experiments. The simulated results are in good agreement with the experimental ones. We defined the four dimensionless parameters: capillary number (*Ca*), the relative size between pores of the metal foam and the droplet (*Dd*), and the metal foam thickness (*Ld*). The results show that for *Ld* of 7.32 and *Dd* of 1.82, the *Cac* for droplet passing metal foam of 60 PPI is around 0.061. This *Cac* is lower than that of 0.1 for the case of one droplet passing through one cylinder obstruction because of the frequent periodic reorientation in the flow of the metal foam. For *Ld* of 7.32 and *Ca* of 0.061, when *Dd* = 0.89, no obvious droplet deform can be observed. For *Dd* = 1.07, the droplet is stretched into a linear shape. When *Dd* exceeds 1.25, a polliwog-shaped droplet is formed. The breakup of the droplet occurs at the waist of the deformed droplet and not on the surfaces of the porous obstacles. For *Ca* of 0.061 and *Dd* of 1.25, *Ld* in the range of 2.44 to 9.76 shows a negligible effect on the number of the generated daughter droplets. A *Db* of 3.00 is the critical value to avoid the coalescence of two droplets at the inlet zone of the metal foam in the case of *Ca* of 0.061, *Dd* of 1.33, and *Ld* of 7.32.

**Supplementary Materials:** The following are available online at http://www.mdpi.com/2227-9717/7/12/877/s1, Figure S1: (a) Image of a water droplet surround by silicon oil on pure nickel surface and (b) a droplet simulated by LBM with Gw,k of ± 0.06, Figure S2: Evolution of one droplet breakup: (a) t = 0.09 s, (b) t = 0.096 s, and (c) t =0.102 s, Figure S3: Hydrodynamics of two successive droplets passing a metal foam with Db of 2.00 (Ld = 7.32, Ca = 0.061, Dd = 1.33): (a) t = 0.12 s, (b) t = 0.144 s, and (c) t = 0.168 s.

**Author Contributions:** Conceptualization, J.Z. and X.Y.; methodology, J.Z.; software, J.Z.; validation, J.Z. and X.Y.; formal analysis, J.Z. and X.Y.; investigation, J.Z.; resources, J.Z.; data curation, J.Z.; writing—original draft preparation, J.Z. and X.Y.; writing—review and editing, J.Z. and X.Y.; visualization, J.Z.; supervision, S.-T.T.; project administration, X.Y.; funding acquisition, X.Y.

**Funding:** This research was funded by the China Natural Science Foundation, grant number 21476073 and 21176069; and the Sub Topic of Major Science and Technology Project, grant number 2017ZX06002019-003.

**Conflicts of Interest:** The authors declare no conflict of interest.
