Physics of puffing and microexplosion of emulsion fuel droplets

The physics of water-in-oil emulsion droplet microexplosion/puffing has been investigated using high-fidelity interface-capturing simulation. Varying the dispersed-phase (water) sub-droplet size/location and the initiation location of explosive boiling (bubble formation), the droplet breakup processes have been well revealed. The bubble growth leads to local and partial breakup of the parent oil droplet, i.e., puffing. The water sub-droplet size and location determine the after-puffing dynamics. The boiling surface of the water sub-droplet is unstable and evolves further. Finally, the sub-droplet is wrapped by boiled water vapor and detaches itself from the parent oil droplet. When the water sub-droplet is small, the detachment is quick, and the oil droplet breakup is limited. When it is large and initially located toward the parent droplet center, the droplet breakup is more extensive. For microexplosion triggered by the simultaneous growth of multiple separate bubbles, each explosion is local and independent initially, but their mutual interactions occur at a later stage. The degree of breakup can be larger due to interactions among multiple explosions. These findings suggest that controlling microexplosion/puffing is possible in a fuel spray, if the emulsion-fuel blend and the ambient flow conditions such as heating are properly designed. The current study also gives us an insight into modeling the puffing and microexplosion of emulsion droplets and sprays.

As can be seen, explosive boiling is a key phenomenon in microexplosion of emulsion fuel droplets.Explosive boiling of water has been of great interest since it occurs in many engineering devices.Several interesting phenomena have been revealed by boiling-droplet experiments. 15,16 iling vapor bubbles are initiated by nucleation.Bubble nucleation was observed near the hightemperature boundary of a droplet prior to explosive boiling.Violent boiling started at the edge of the droplet and the vapor-liquid interface was wrinkled by the Landau(-Darrieus) instability mechanism. 15,16 ubble oscillations and droplet fragmentation were also observed, as well as the Rayleigh-Taylor (RT) instability following boiling. 15,16 t is expected that similar phenomena will occur in microexplosion/puffing of emulsion droplets.
Microexplosion and puffing are beneficial to achieving enhanced secondary atomization.In fuel sprays, these eruptive secondary atomization mechanisms will help to meet conflicting spray requirements of longer penetration (achieved by large droplets) and better evaporation/mixing (achieved by small droplets).Another important benefit is that generation of evaporated water vapor reduces the local flame temperature, thereby reduces the emissions of NO x . 1 Addition of OH radicals due to evaporated water vapor also helps soot reduction. 1espite its potential benefit to fuel spray atomization, the dynamics of microexplosion/puffing is, however, not fully understood yet.Accurate prediction of emulsion-fuel sprays requires more knowledge on physical processes such as droplet heating, boiling, breakup, and also on their complex interactions with the ambient turbulent gas flow.[6][7][8][9] But its relevance to droplets of O (10 μm) in a realistic fuel spray has not been well addressed. 1Currently, a few experimental photos of microexplosion/puffing of emulsion fuel droplets are available under realistic fuel spray conditions, 4,5 but the detailed processes are not resolved and physical mechanisms remain unknown.
8][19][20][21] One approach used Rayleigh's growth model, which solves the growth of one single quasi-steady symmetric vapor bubble at the center of the droplet. 18Although marginally acceptable for miscible multicomponent droplets, 2,14 the assumptions made in such an approach may largely deviate from what has been discovered for emulation droplets in fuel sprays.For example in an emulsion droplet, multiple bubble formation can occur at multiple water sub-droplets. 2,5  recent experiment 4 implies that the dominant mode of eruptive secondary breakup in a fuel spray may be puffing.Some other models mainly predict the final outcome of droplet breakup by microexplosion/puffing, such as the secondary droplet size distribution. 20,21 t is clear that in order to predict microexplosion, more knowledge on the physics of microexplosion of emulsion fuel droplets is needed.
In this study, high-fidelity numerical simulation is used to investigate the microexplosion phenomena from first principles.The objective is to identify dominant mechanisms of microexplosion/puffing and obtain modeling insights for Eulerian-Lagrangian simulation of practical-scale fuel sprays under microexplosion conditions.Based on our multiphase flow code, [22][23][24] the gas-liquid interfaces (between the oil droplet and the ambient gas, between the water sub-droplet and evaporated water vapor, and between the oil droplet and the water vapor) are directly solved.To the best of our knowledge, such a simulation study is the first of its kind.To obtain detailed information, a single emulsion droplet at a characteristic scale in a fuel spray is considered and the evolution of microexplosion/puffing is investigated.Comparison with theoretical predictions is performed at each stage of the results discussion to obtain an accurate, quantified understanding of the physics of puffing/microexplosion phenomena.
Since it is difficult and expensive to directly solve the physical phenomena that initiate explosive boiling such as emulsion formation and droplet heating, in this study we instead use physically reasonable initial conditions, for which past findings [1][2][3][4][5][6][7][8][9][10][11][12][13] have been referred to, to perform parametric studies to investigate the effects of these initial conditions on emulsion droplet  and puffing.First, the distribution of dispersed water sub-droplets plays an important role in microexplosion. 1,6,8 B changing the initial dispersed-phase sub-droplet size, Suzuki et al. 9 showed that the occurrence probability of microexplosion could be varied.2][13] Second, droplet heating causes thermocapillary (Marangoni) migration of inner sub-droplets toward the hot side, 1,7,25,26 which may lead to sub-droplet coalescence.In experiments with a longer time scale (∼O(1 s) for heating) using larger droplets, complete phase separation of the oil and the water was observed before microexplosion. 7However, in a realistic fuel spray where the injected fuel does not have a sufficient time for complete phase separation, sub-droplet coalescence occurs only partially. 4,5 2][13] In such a case, the observed explosion-type breakup (∼10 μs from the beginning of droplet deformation to explosion) of small emulsion droplets (∼20-50 μm) is considered due to explosive boiling of multiple water sub-droplets. 5Under these considerations, sub-droplet coalescence is not directly considered in this study, and the number and location of water sub-droplets are varied to investigate their effects on microexplosion and puffing.
Nucleation is the initial formation of a small vapor bubble nucleus, which finally leads to explosive boiling.The initiation of nucleation can be modeled by the homogeneous nucleation theory. 27,28 he nucleation probability is sensitive to the local temperature, which has been confirmed by experiments. 27,28 t was also observed that the nucleation location was always close to the wateroil interface. 27Typically for miscible multicomponent droplets, one large bubble is formed at the center of the droplet. 2,14 ut for immiscible emulsion droplets to be studied in this paper, multiple bubbles can be formed simultaneously. 2,27 ollowing these findings, initial bubbles are placed at the surface of water sub-droplets.Predicting the nucleation site and resolving tiny nuclei are difficult and not attempted in this study.
The rest of this paper is organized as follows.In Sec.II, the computational setup is explained.In Sec.III, numerical results are presented and physical analyses are made for each stage of puffing and microexplosion.The effect of initial conditions on microexplosion and puffing is also discussed.Finally in Sec.IV, concluding remarks are given.

A. Mathematical formulations and numerical procedures
The governing equations for density ρ, velocity u, temperature T, and species mass fraction Y i (oil, H 2 O, O 2 , and N 2 ) are solved. 24They are where with P TH = T(∂p/∂T) ρ .For an ideal gas, P TH = p.Q u includes the viscous and surface tension terms.Q T includes the work by viscous forces and heat conduction modeled by Fourier's law.Q Y i is the mass diffusion term modeled by Fick's law.S * are the source terms due to boiling and will be defined later.
The equations for level-set functions are added to the above system to capture interfaces. 29,30 o level-set functions F i ( = F 1 , F 2 ) are used and follow: where S F,i denotes the surface regression due to boiling.Each F i is a signed distance function and F i = 0 represents an interface.The level-set method is combined with the MARS (Multi-interface Advection and Reconstruction Solver) method, 31 a kind of VOF (Volume of Fluid) methods, to improve the mass conservation. 22,32 rface tension is modeled by the CSF (Continuum Surface Force) method.33 The surface tension force is where σ is the surface tension coefficient, κ is the local surface curvature, n is the surface normal unit vector, and δ is the delta function to identify the surface.This force is included in Q u as F s /ρ.Boiling can be formulated as jump conditions at the interface.In the general form of evaporation, the jump conditions of heat and mass fractions are 34 where h l is the latent heat of evaporation, ω is evaporation rate, λ is thermal conductivity, and D dif is diffusion coefficient.The subscripts L and G denote liquid and gas, respectively.The square brackets denote the difference of a variable f between the liquid and gas phases at the interface, 34 i.e., [f] The surface velocity u S is defined as the sum of the liquid velocity and the surface regression velocity, i.e., u S = u L + s L and s L = s L n = ( ω/ρ L )n.Therefore, the velocity jump is 34 For boiling, the species mass fraction does not change across the interface (Y H 2 O,G = Y H 2 O,L = 1), so Eq. ( 6) is automatically satisfied.The evaporation rate is thus determined by the heat condition in Eq. (5). 35,36 erived from the above jump conditions, the boiling source terms are Solving Eq. ( 1) is split into three stages: the advective phase, the non-advective phase, and the acoustic phase.In the advective phase, is solved by the CIP (Cubic Interpolated Pseudo-particle or Constrained Interpolation Profile) method. 37,38 he third-order scheme uses a third-order polynomial curve fitting for flow variables and their derivatives to solve advection. 37,38 n the non-advective and acoustic phases, is solved by the CUP (Combined and Unified Procedure) method. 38The CIP-CUP method can be applied to both incompressible and compressible flows, using a unified treatment of the speed of sound. 38In this study, this method is used for simulating microexplosion, where explosive boiling occurs and compressibility needs to be considered.
In the acoustic phase, the pressure is determined by solving a Poisson equation.Omitting S * and Q * for simplicity and discretizing Eq. ( 11) in time, one obtains where the superscripts n denotes the time step n and * the advanced time step in the acoustic phase.
The ρ and T equations in Eq. ( 12) can be rewritten as and the divergence of the u equation is Using the thermodynamic relation one finally obtains the Poisson equation for p * , where c T = (∂ p/∂ρ) is the isothermal speed of sound and c s = (c 2 s is the adiabatic speed of sound.If c s is infinite, Eq. ( 16) reduces to the incompressible formulation.If c s is finite, it includes the compressible effect.
The equation of state (EOS) for gas is given by the ideal-gas law, and the EOS for liquid water and oil is given by Tait's empirical formulation 39 where the subscript 0 denotes the reference state, and A w = 296.3MPa and γ w = 7.415.Depending on the flow conditions of interest, the above numerical method can solve both incompressible and compressible flows. 38For the turbulent atomization cases in our previous research, 22,23 the same code was used with parameters fit for incompressible flows.For the evaporating spray case, 24 the spray was solved as a compressible flow, although the compressibility effect was minor.In this study of microexplosion, compressibility is taken into account by solving Eq. ( 16).

B. Code validation
For microexplosion/puffing, there is currently no suitable experimental or numerical data available for direct comparison.Since boiling and gas-liquid interface oscillation are two additional key phenomena for microexplosion, we studied the code performance on reproducing these key phenomena by comparing simulation results in simple configurations against analytical solutions.

Boiling
A sucking interface problem is solved. 40A boiling interface is located in a domain that separates superheated water from vapor, as shown in Fig. 1(a).The left boundary is a wall.The interface moves toward the right direction as vapor is generated from superheated water.The vapor temperature is

Linear and nonlinear droplet oscillation
Droplet oscillation affects the water sub-droplet dynamics (see Sec. III C).The analytical angular frequency for linear oscillation is 41 in two-dimensional (2D) and three-dimensional (3D) configurations, respectively, where R 0 is the equilibrium droplet radius and n is the oscillation mode (the lowest mode n = 2 is used).Figure 2 exemplifies the flow fields at two time instants in a 3D configuration.The initial velocity is u = 0 and the domain size is 5R 0 × 5R 0 × 5R 0 .All the boundaries are open boundaries.For R 0 = 10.6 × 10 −6 m and σ = 50 × 10 −3 N/m, the theoretical oscillation period (T osc = 2π /ω) is T osc = 10.0 × 10 −6 s for 3D cases.The simulation result is T osc = 9.8 × 10 −6 s.Again, excellent agreement is obtained.
Figure 3 presents the comparison of nonlinear droplet oscillations in large amplitude. 42For this example, the initial droplet shape is a prolate spheroid of a/b = 3, where a is the length of its semi-major axis and b is the length of its semi-minor axis.The volume of this prolate spheroid is  42 The time shown is normalized by (ρ l r 3 /σ ) 1 / 2 .V = 4πab 2 /3 and the equivalent radius of a sphere of the same volume is r = 3 1/3 b.The parameters are ρ l = 700 kg/m 3 , σ = 20 × 10 −3 N/m, r = 2.9 × 10 −4 m, μ l = 635 × 10 −6 Pa s and the Reynolds number Re = (σ r/ρ l ) 1/2 /ν l is 100.A 3D grid system of 135 × 135 × 135 points is used.Figure 3 shows the simulated droplet shapes.The results predicted by the current code show excellent agreement with the published results 42 both in droplet shapes and in the oscillation period.It should be stressed that even if the shape deformation is considerable, the current code can capture the oscillation dynamics well.
It should also be mentioned that in our previous paper, 22 our code has been validated for dynamic capillary wave propagation and droplet pinch-off from a ligament.By these and the comparisons performed in this study, it can be said that the current code accurately predicts explosive boiling and interface dynamics.

C. Simulation parameters and initial conditions
In order to investigate the detailed dynamics of microexplosion, parametric studies are conducted.The initial conditions are listed in Table I. Figure 4 illustrates droplet configurations.Cases A1-A5 have only one water sub-droplet inside the parent oil droplet and the size of the water sub-droplet is varied.The water sub-droplet surface depth d/R 0 and the bubble inclination angle θ b are identical for A1-A5.To have a different characteristic boiling time scale, the latent heat of evaporation h l is doubled in Case A6.In A7, the bubble location is deeper than in A5.In Cases B1-B4, the size of the water sub-droplet R 2 /R 0 is identical, but the locations of the water sub-droplet and the bubble are varied.Cases C1-C4 have multiple water sub-droplets.Each of C1-C4 has a corresponding case in series A. The multiple water sub-droplets and initial bubbles are placed symmetrically.Each 3D case also has a corresponding 2D case.Case A1-3D is a 3D version of Case A1.The water sub-droplet radii of A1 and A1-3D are identical, but the volume fraction of the water sub-droplet is reduced by the three-dimensionality in A1-3D.Similarly, Cases B3-3D and C3-3D correspond to Cases B3 and C3, respectively.In Case C3-3D, the number of water sub-droplets is increased to six to keep the symmetry in the 3D configuration.
Considerable computational cost is required to numerically solve a 3D case, especially when the sub-droplet shape, explosive boiling interface, and nonlinear droplet shape oscillation need to be resolved in the present study.On the other hand, 2D cases can help clearly reveal the underlying physics with a simple configuration and thus reduced computational cost.Therefore, in the present study, we have taken the strategy of using 2D cases as the main configuration for parametric study with a 3D counterpart in each of the case categories to take advantage of what 2D and 3D simulations can offer.It will be demonstrated that 2D simulations can largely determine the physics of puffing and microexplosion.By investigating the relationship between 2D and 3D cases, the complete physical mechanisms of microexplosion can be obtained.For all the cases, the physical properties are given as follows.The liquid properties are set to be similar to those of hexadecane and water, for which published data on vapor nucleation are available in Ref. 27.The oil properties used are within the property ranges of general liquid hydrocarbon fuels.The oil density is ρ oil = 770 kg/m 3 , water density ρ water = 850 kg/m 3 , the ambient air pressure 3 MPa, and the water boiling temperature 503 K.The ambient air temperature and the liquid water/oil temperature are set to be 553 K. Therefore, the superheat degree is 50 K. 27The latent For the liquids, constant properties are used due to the short time scale and limited temperature range of a microexplosion process.The time scale of microexplosion t∼O (1 μs) is much shorter than that of heat transfer within the droplet t∼O (100 μs).Therefore, temperature variation in water sub-droplets due to heat transfer by boiling is confined to a very thin layer and can be neglected.The temperature of the rest of the droplet volume remains at the initial temperature, and their viscosity remains at the initial value.Additionally, the temperature range in the present study is 503-553 K with no external heat source.In this range, the viscosity variation is smaller than that within a lower temperature range such as 300-350 K.Even in the thin layer, the temperature variation is limited.Therefore, the microexplosion dynamics of an emulsion droplet is not significantly affected by the variation of fluid properties of the liquids, such as viscosity.For the gas phase, the physical properties are derived by the NIST (National Institute of Standards and Technology) database. 43The initial velocity is zero.Evaporation of the oil droplet is neglected because its time scale is not on the same order of magnitude as that of microexplosion, and thus its effects are secondary in the present study.
The initial oil droplet diameter is set at D = 30 μm, which is a typical scale in a real-scale fuel spray.The computational domain size is 4.2D × 4.2D in 2D cases and 4.0D × 4.0D × 4.0D in 3D ones.All the boundaries are open boundaries.The total number of grid points is 381 × 381 in 2D cases and 381 × 381 × 381 in 3D ones.The minimum grid spacing is 0.26 μm.The grid spacing near the boundaries is gradually stretched towards the domain boundaries.Grid resolution tests have been conducted and the current resolution has been found to be sufficient for the phenomena in interest.Figure 5 shows the profiles of the displacement of the centers of gravity (CG) of the parent oil droplet and the water sub-droplet of Case A1 with different gird resolutions.This statistics will be discussed in detail in Sec.III C 2. It can be seen that the results obtained by the 85 × 85 grid cannot reproduce the droplet motion well.With the finer resolution of the 381 × 381 grid, the droplet motion is well reproduced with deviations around 4%-7% from the results of the 789 × 789 grid.The grid number 381 × 381 is thus chosen.It should be finally pointed out that grid resolution tests have been also performed for the other statistics in this paper and the same conclusion is attained.
As discussed in the literature, 2,27,28 vapor bubbles are generated by nucleation.Initially, the nucleation critical radius is quite small (<0.01 μm) compared to the water sub-droplet size, and the inner vapor pressure p v is high (∼100 atm) due to surface tension.The first-stage growth from a nucleus to a small bubble is "inertia controlled." 2 But very quickly (typically ∼0.1 μs compared to t∼O (1-10 μs) for microexplosion), this high pressure is lowered by the bubble growth.Then, the next growth stage is nearly isobaric and "diffusion controlled." 2 As it is impossible to resolve the initial nuclei and as this initial stage quickly ends, the simulations in this study start from a finite-size bubble.The bubble size is initialized to 0.07R 0 .The effect of the initial bubble size has been checked.If a larger bubble is placed, the basic dynamics is identical, but the early development phase of the bubble is more skipped and shorter.To include a sufficient length of the bubble growth, this value was chosen.

A. Overall dynamics
The results of Cases A1 and A3 are used to explain the overall dynamics.Figure 6 shows the temporal sequence of the initial puffing and subsequent after-puffing dynamics.As a single water sub-droplet is placed near the parent droplet surface, the primary mode is puffing.A solid line is used to show an interface either between a liquid species (the oil or water) and a gas species (water vapor or air) or between the two liquids (the oil and water).The two liquids are distinguished by darker gray (red) for water and lighter gray (orange) for the oil.In the gas phase, the contours show the mass fraction of oxygen, so lightest gray (yellow) illustrates air (O 2 /N 2 ) and gray (blue) water vapor.
Initially, the bubble near the water-oil interface grows by boiling.As the liquid oil "wall" is thin near the parent droplet surface, the bubble will push this side strongly toward outside.Consequently, the liquid oil wall soon ruptures and liquid fragments and water vapor are ejected (Figs.6(a) and 6(b): t = 0.4 μs).This is puffing.After puffing, the boiling surface keeps ejecting vapor, which pushes the water sub-droplet itself.At the edge of the water/vapor liquid-gas interface near the water/oil boundary (marked by gray (blue) arrows in Fig. 6), the edge becomes unstable by the Landau-Darrieus (LD) instability and the troughs therefore deepen.The edge regression proceeds to wrap up the water sub-droplet and finally detaches it from the parent droplet.During this process, the RT instability development is also observed at the oil-vapor inert (non-boiling) interface (marked by solid black arrows in Fig. 6(b)).After the detachment, the water sub-droplet stops pushing the oil droplet and the breakup process mostly ends.In Case A1, the boiling surface of the water subdroplet oscillates strongly.By this oscillation, the water sub-droplet changes the shape (indicated by dashed arrows in Fig. 6(b) at t = 2.2 and t = 3.2 μs).In such a case, the shape of the boiling surface is flattened and some detached part of the water sub-droplet re-merges with the oil droplet.The interface edge regression is thus suspended and then restarts.It can be expected that the time for the detachment of the whole water sub-droplet becomes longer.In 3D cases, the basic dynamics is the same, as illustrated in Fig. 7 for Case A1-3D.Similarly, puffing occurs on the right side and vapor is ejected from the hole created by puffing.The water sub-droplet penetrates into the oil droplet, along with its shape oscillation.The edge regression can be also seen, which finally leads to the water sub-droplet detachment.The three-dimensional LD and RT instability development creates several bubble cells around the water-oil edge (in 2D cases the edge bubble is a donut-like tube).
In the following, each stage of the evolution is discussed in quantified detail.While the qualitative similarity of the puffing dynamics between 2D and 3D cases has been shown, the quantitative difference between 2D and 3D cases are analyzed in detail.

B. Vapor bubble growth and oil droplet breakup: Puffing
Based on Rayleigh's general formulation of bubble growth, one obtains 2, 44-46 where 2 , and p v are the bubble pressure.The solution to Eq. ( 19) is R(t) ∼ At in the inertia-controlled stage and R(t) ∼ Bt 1 / 2 in the diffusioncontrolled stage. 2,44,46 I the formulation of Mikic et al., 46 A = (2h l ρ v T /3ρ L T sat ) 1/2 and B = (12α L /π ) 1 / 2 • Ja, where Ja is the Jacob number J a = T c L ρ L / h l ρ v , ρ v is the vapor density, T sat is the saturation temperature, and α L is the liquid thermal diffusivity.Figure 8 shows the time history of the bubble size (in equivalent radius) and the bubble pressure for several cases.If the bubble is placed close to the parent droplet surface, it quickly bursts (Cases A1, A3, A4, and A1-3D).If the initial bubble is placed near the parent droplet center (Cases B2, B3, and B3-3D), the bubble growth takes a longer time.In the early stage, there is pressure effect.Especially in longer-growth cases (Cases B2, B3, and B3-3D), the later bubble growth is close to the theoretical curve predicted by R(t) ∼ Bt 1/2 for the diffusion-controlled stage.
During the bubble growth, the bubble vapor pressure oscillates.The natural oscillation frequency of an isothermal bubble surrounded by a uniform liquid can be predicted by 45 where R b0 is the equilibrium bubble radius.Equation ( 20) has been derived for a non-growing bubble, but can give an estimate of the oscillation frequency of a growing bubble.For Cases B3 and B3-3D, the oscillation period T osc predicted by Eq. ( 20) is T osc = 2π/ω N ∼ 0.66 μs, which is close to the observed periods (Fig. 8(b)).The bubble therefore oscillates at about the natural frequency during its growth.The slightly longer period of Case B2 than B3 is due to the larger density of water on the droplet surface side.The size of the secondary droplets produced by initial puffing is determined by the bubble growth.Because the bubble pushes out adjacent liquid, its depth has a strong impact on the secondary droplet size.Figure 9 shows the breakup time and the secondary droplet size against the initial depth of the bubble center d bub .The open triangles denote Cases A7 and B2, where in the outer region of the bubble is water, and the filled squares denote Cases A2-A5, B1, B3, and C1, where in the outer region of the bubble is the oil.The general trends are the same.As the initial bubble depth is larger, the breakup time is longer and the secondary droplet size is larger.The breakup characteristics are therefore determined by which liquid (the oil or water) is pushed by the bubble.

C. Water sub-droplet detachment: After-puffing
After the initial puffing, a large portion of the water sub-droplet still remains attached to the parent droplet.Continuing vapor ejection generates "thrust" and affects the droplet dynamics.During this phase, the boiling surface exhibits dynamic motion and its area increases.Finally, the water sub-droplet is detached from the parent droplet.
Figure 10 shows the time history of the boiling surface area and the amplitude of shape oscillation of the water sub-droplet.The surface area is non-dimensionalized by the equilibrium value of Case A1 in Fig. 10(a) and Case A1-3D in Fig. 10(c).In 2D cases, the increase rates of the boiling surface area S boil in the first edge-regression phase are similar among A1, A2, A3, and A4.The evaporation rate is halved in Case A6, so the increase rate of S boil is smaller.In Cases A3 and A4, the water sub-droplet is detached shortly after the initial puffing.In Cases A1 and A2 (and A6), the boiling surface area increases in two stages.After the initial rise, S boil maintains or even decreases.Then it increases again.After the first cycle of edge regression, the water sub-droplet shape is largely deformed by oscillation (Fig. 10(b)), which makes the water sub-droplet re-merge with the parent droplet (see the shape change in Fig. 6(b) between t = 2.2 and t = 3.2 μs indicated by dashed arrows).The re-merging is especially evident in Case A1, where the boiling surface area decreases considerably.The already detached edge is re-attached, and the edge regression process starts again.Finally, the water sub-droplet is wrapped by the exposed boiling surface and detached from the parent oil droplet.The only difference among Cases A1-A4 is the water sub-droplet size.In A6 where the evaporation rate is reduced due to the increase of the latent heat of evaporation h l , the initial re-merging similarly occurs during t = 1.2-2.4μs.However, the extent of the boiling surface area cancellation is not as considerable as that in Case A1, because the initial rise of the boiling surface area decreases due to the reduced evaporation rate.The 3D case A1-3D shows the same puffing and after-puffing dynamics.The boiling surface area shows a two-stage increase (Fig. 10(c)) and the shape oscillation can be also seen (Fig. 10(d)).The difference in the interaction between the regression and oscillation is due to the difference between a 2D case and its 3D counterpart, which will be discussed later.Next, the observed droplet dynamics due to puffing is studied in detail.First, the edge regression is investigated.Then, the one-stage detachment and two-stage re-merging/detachment of the water sub-droplet are examined.

Edge regression
It is known that a boiling interface exhibits the Landau-Darrieus (LD) instability due to evaporating mass flux at the interface, 15,16,47 as illustrated in Fig. 11.Due to the evaporating mass flux, the streamlines are bent at the liquid troughs and crests.This motion works to increase the shape deformation amplitude.Even if the liquid motion is irrotational, vorticity is generated in the gas phase as illustrated by the dashed arrows. 47This vorticity generation plays a role in the instability development.
Figure 12 shows the vorticity component normal to the paper before and after puffing for Case A1.Positive and negative vorticities mean counterclockwise and clockwise rotations, respectively.At t = 0.2 μs when the initial puffing will soon occur, the vapor bubble is pushing the thin liquid "wall" toward the outside air in the right direction.The observed strong symmetric vorticity is due to this outward jet motion and not due to the LD instability mechanism.This motion is strong enough to pull the boiling surface of the water sub-droplet toward the same direction.As a result, the concave shape becomes convex after puffing (marked by bold gray (red) arrows in Fig. 12(a)).This induces oscillation of the water sub-droplet as a rebound motion (see also Fig. 6(b) at t = 2.2 μs and t = 3.2 μs).At t = 0.8 μs and t = 1.6 μs, the vorticity structures near the edge locations in trough shapes are marked by small solid black arrows.Commonly, in these concave regions, similar positive and negative vorticity pairs are observed.The vorticity generation is due to the local evaporative mass flux from the boiling surface.Note in the immediate vicinity the oil/vapor interface is nonboiling and inert, and therefore does not contribute to vorticity generation.The vorticity pairs work to deepen the trough as shown in Fig. 11.In this sense, the mechanism of the regression initiation is the LD instability.As a result, the boiling surface of the water sub-droplet gradually extends and regresses toward inside.During this regression process, the edge regions exhibit growing bubble-like structures and leave a wavy pattern (a chain of bubbles and necks) on the detached inert interface between the water vapor and the liquid oil (Fig. 6(b) at t = 5.5 μs).This is induced by the RT instability, whose mechanism is due to density difference and surface acceleration. 48In the neck regions created by RT plumes (marked by bold arrows in Fig. 12(a)), the locally outgoing vapor flow is accelerated by the nozzle effect and the vorticity magnitude is large, as shown schematically in Fig. 12(b).This motion eventually pushes the developed neck plumes outward.By repeating these cycles, the edge regression proceeds and the water sub-droplet will eventually detach itself from the parent oil droplet.
Although it is difficult to quantitatively estimate the speed of edge regression, which is initiated by an edge bubble, the bubble growth rate predicted by Mikic et al. 46 suggests that the growth rate variation is not large, if the time scale of puffing is considered.The growth rate is given by 46 where , and B = (12α L /π ) 1/2 • Ja.This equation was originally derived for a single bubble in a non-uniform temperature field.In the current cases, bubbles are formed one by one in sequence (Fig. 12(b)).With the time scale of 1-2 μs to generate each bubble, the speed variation is slow, so the average growth speed can be estimated by Eq. ( 21) as a constant, i.e., d R/dt ∼ 5m/s.Applying this value to 2D cases (A1-A4), it roughly gives d S boil /dt ∼ 2 • d R/dt = 2.6 × 10 −12 m 2 /μs.The actual observation gives 3.9 × 10 −12 m 2 /μs, which is on the same order of magnitude as the estimation.Therefore, repeated bubble generation is a reasonable physical mechanism, which will be used in Subsection III C 2. In A1-3D, the edge bubble generation evolves in three dimensions.Meanwhile, the RT instability also grows in three dimensions.Therefore, multiple bubble cells are created, as shown in Fig. 7.The overall increase rate of the boiling surface area does not vary rapidly and is almost constant, as shown in Fig. 10(c).

One-stage detachment and two-stage penetration/detachment of a water sub-droplet
Figure 13 shows the temporal history of the distance between the CGs of the parent droplet x CG, o and the water sub-droplet x CG,w , i.e., d ow = x CG,o − x CG,w .Here, x CG,w = 0 at t = 0 and therefore x CG,w > 0 means that the water sub-droplet is moving toward the parent droplet (Fig. 14).To set the relative distance to be zero at t = 0, d ow − d ow,0 = d ow − d ow (t = 0) is normalized by R 0 and plotted in Fig. 13.If the relative distance d ow − d ow,0 is negative, the water sub-droplet is penetrating into the oil droplet.If d ow − d ow,0 is positive, the water sub-droplet is departing from the oil droplet.Among Cases A1-A4, the only difference of initial parameters is the water sub-droplet size.As shown here, if the water sub-droplet size is large (A1, A2, and A6), the water sub-droplet tends to penetrate into the parent droplet.For the same diameter, the water sub-droplet of Case A1-3D penetrates faster than that of its 2D counterpart case A1.Factors affecting the water sub-droplet motion include the inertia of the initial bubble burst, the thrust force generated by boiling, and the shape oscillation.To understand the competing effect of these phenomena, a simplified setting is considered, as shown in Fig. 14.Note in the following analysis, the effect of the shape oscillation is neglected.
In both 2D and 3D cases, the water sub-droplet motion can be modeled as and Here, any effect on or due to the oil droplet is neglected for simplicity (x CG, o = const.and d ow − d ow,0 = −x CG,w ).L is the hypothetical depth in the third direction in 2D cases.f th is the thrust force due to boiling, and f p is the net force in the direction opposite to the thrust after the bubble burst.For simplicity, f th and f p are defined as the force per unit area and assumed to be constant.S eff is the effective (net) projected area of the boiling surface or the bubble surface, namely, S eff, boil = 2h boil L and S eff, bub = 2h bub L in 2D cases and S e f f,boil = π h 2 boil and S e f f,bub = π h 2 bub in 3D cases.h bub is assumed to be a constant and h boil = R sin θ , where θ in Fig. 14(a) is measured in the same way as θ b in Fig. 4(a).Therefore, S eff, boil increases with the increase of θ for 0 < θ < π/2, and decreases for θ > π/2.From the observation in Fig. 10, the rate of increase of the boiling surface area by edge regression is set to be a constant for simplicity, i. with in 2D cases, and in 3D cases.The time range is t < t w , where t w is the instant when the edge regression completes and the entire water sub-droplet detaches from the parent oil droplet.Figure 14(b) schematically shows the solutions of Eq. ( 24).When R is small, the first term −At 2 is larger, i.e., the bubble pulling effect is dominant, and thus the water sub-droplet tends to move away from the parent droplet.When R is large, the second term B(t − C) is larger, i.e., the boiling vapor thrust effect is dominant, and thus the water sub-droplet tends to penetrate toward the parent droplet.Cases A3 and A4 in Fig. 13 correspond to small-R cases.The water sub-droplet thus departs from the parent droplet easily.On the other hand, Cases A1, A2, A6, and A1-3D in Fig. 13 are large-R cases.The water sub-droplet thus penetrates into the parent droplet.Faster penetration of the 3D water sub-droplet (Case A1-3D) can be understood as follows.For a large R, the bubble effect is relatively small (h bub ∼0), and then the ratio of acceleration is ẍCG,w(3D) / ẍCG,w(2D) = 3π h boil /8R.This gives ẍCG,w(3D) / ẍCG,w(2D) > 1 if h boil /R > 8/3π = 0.85.Therefore, if the boiling surface area is sufficiently large, the water sub-droplet penetration is faster in 3D cases.
It is clear now that the direction of the water sub-droplet motion determines whether the water sub-droplet detachment is a one-stage or two-stage process.In each case, the water sub-droplet starts shape oscillation (the rebounding motion) after the bubble bursts.If the penetration is negative, this shape oscillation effect cannot lead to re-merging of the water sub-droplet with the parent droplet.In fact, in Cases A3 and A4, the water sub-droplet has already been detached largely before the rebounding shape oscillation occurs.
If the penetration is positive (Cases A1, A2, A6, and A1-3D), the water sub-droplet is flattened by the rebounding shape oscillation and re-merges with the parent droplet (Fig. 6(b)).The re-merging cancels out part of the already generated boiling surface and leads to the decreasing boiling surface area (Fig. 10).The degree of the re-merging depends on the interaction between the edge regression and shape oscillation.Particularly in Case A1, the rebounding motion pushes back the water subdroplet before it largely detaches (Fig. 6(b)), which results in considerable cancellation of the boiling surface in the middle of the process (Figs.10(a) and 10(b)).In Case A6, this cancellation process also occurs during t = 1.2-2.4μs, but it is not as extensive as in Case A1, because the initial increase of the boiling surface area is reduced due to the reduced evaporation rate.In Case A2, the re-merging starts around t = 1.4 μs, when the water sub-droplet has almost completely detached (Figs.10(a) and 10(b)).Because the penetration is positive, this interaction results in a subsequent fluctuation of the boiling surface area around the equilibrium value before the final detachment (Fig. 10(a)).
Cases A1 and A1-3D have the same water sub-droplet radius, but the re-merging is weaker in A1-3D due to the difference in the timing of the edge regression and shape oscillation, as shown in Fig. 10.The shape oscillation period is shorter in 3D than in 2D; it is about 0.8 times of that of A1 (see also Eq. ( 18)).And the water sub-droplet penetration is faster in the 3D case than in the 2D case, as discussed above.These factors lead to the weaker interaction between the edge regression and the shape oscillation in A1-3D.
In the long run, the edge instability mechanism exists as long as boiling continues, and the surface regression restarts after the re-merging.Finally, the water sub-droplet detaches completely from the parent oil droplet and the breakup due to puffing ends.

D. Effect of initial locations of the water sub-droplet(s) and bubble formation
As seen in Subsection III C, the water sub-droplet dynamics, mainly the penetration direction of the sub-droplet, is a key to determining the outcome of the parent droplet breakup.In this subsection, the initial locations and number of the water sub-droplet(s) and vapor bubble(s) are varied to investigate their effects on the puffing and microexplosion dynamics.In series B, the water sub-droplet size is the same as that of Case A1.Its position and the initial bubble location are varied.
If the location of the bubble formation is inclined (Case B1; Fig. 15(a)), the bubble development is asymmetric.The upper right region of the oil droplet, where the liquid "wall" is the thinnest, ruptures first.From where the puffing occurs, the vapor leaks outside.This slowly bends the trajectory of the water sub-droplet.The water sub-droplet gradually rotates counterclockwise and eventually breaks up a portion of the oil droplet.In Case B2 (Fig. 15(b)), the bubble is formed at a deeper location inside the oil droplet.The water sub-droplet moves toward the oil droplet surface on the right side.The dynamics develops more slowly than in Case A1 due to ρ water > ρ oil .Once the puffing occurs, the water sub-droplet itself also goes out of the parent droplet.In Cases B3 and B4 (Figs. 15(c) and 15(d)), the water sub-droplet moves leftward.The difference of the water sub-droplet depth d/R 0 between the two cases determines which side of the oil droplet ruptures first.
So far, single-bubble cases have been considered.In a realistic spray, the bubble formation and explosive boiling in multiple water sub-droplets in a single oil parent droplet is also likely.If bubbles form simultaneously at multiple locations, the subsequent combined effects may break up the entire parent droplet to a much greater extent compared to the breakup by single puffing.The droplet breakup thus appears as complete microexplosion. 5In fact, Lasheras et al. 2 pointed out that simultaneous nucleation of a large number of water sub-droplets is possible.Therefore it is worth examining the effect of multiple explosions on an oil droplet.
Figure 16 shows the 2D C-series cases with multiple explosions and Fig. 17 a corresponding 3D case.Initially, when the influence of puffing is confined to its local region, each explosion proceeds almost independently as "local puffing."As the dynamics develops, mutual interactions between puffing start, and the whole parent droplet may break up, whereas single puffing in general does not break up the whole parent oil droplet.
Figure 17 presents two instantaneous flow fields of the 3D case C3-3D.When each puffing is independent (Fig. 17(a)), the oil droplet shape is similar to that in Fig. 16(c) at t = 0.5 μs.After interactions among the six water sub-droplets start at t = 3.0 μs in Fig. 17(b), extensive fragmentation of the oil and water liquids will occur.The droplet shape at a later time t = 3.0 μs (Fig. 17 somewhat different from that of the 2D Case C3-2 (Fig. 16(c) at t = 3.0 μs) due to mainly two reasons.First, the water sub-droplets penetrate faster in the 3D case (see Sec. III C 2).Therefore, in the 2D Case C3-2 the water sub-droplet detachment occurs before the collision of sub-droplets at t = 3.0 μs (Fig. 16(c)).In contrast in the 3D Case C3-3D, the water sub-droplets collide with each other inside the parent droplet before they detach from the parent droplet (Fig. 17(b)).Second, the breakup is easier in the 2D case.As shown in Subsection III E, water vapor in the 3D case need to push relatively more oil during the breakup compared to the 2D case.By these reasons, the 3D breakup is slightly different from the 2D breakup.

E. Degree of breakup
As already seen, the degree of the parent droplet breakup depends on the water sub-droplet dynamics.Here, the degree of breakup is quantified and defined as where V oil indicates the volume of oil whose distance away from the CG of the initial oil droplet is within the initial droplet radius R 0 .t end is taken as the time instant when the water sub-droplets are detached and the major deformation of the oil parent droplet ends.This is a global value for estimation only and does not contain detailed shape deformation information.A large D b means that the parent droplet is largely deformed and blown far away by the water vapor.
Figure 18 shows D b for both 2D and 3D cases with a single or multiple water sub-droplets.The solid lines connect cases with the same water sub-droplet radius.For the A series (A1-A5 with a single water sub-droplet), the trend is clear.The larger D b in A1 and A2, compared to those in A3-A5, is due to a longer penetration time of the water sub-droplet.The series C cases have multiple water sub-droplets.With more water sub-droplets, D b increases, but not linearly.For the series B cases (B1-B4), D b is not the same even if the water sub-droplet size is the same.As shown in Sec.III D, the initiation location of explosive boiling and the subsequent water sub-droplet motion determine the degree of breakup.D b of A1-3D is smaller than that of A1.This is because the ratio of the mass of the surrounding oil to that of explosive vapor r = m o /m v is different in the 2D and 3D cases.The parent oil droplet radius is R 0 and the water sub-droplet radius R. The vapor mass m v is m v,2D = ω t • 2π RL in 2D cases and m v,3D = ω t • 4π R 2 in 3D cases.t is a time period.Similarly, the mass of the surrounding oil m o is m o,2D = ρ oil • π (R 2 o − R 2 )L in 2D cases and m o,3D = ρ oil • 4π (R 3 o − R 3 )/3 in 3D cases.Denoting ξ = R/R 0 , the ratio F = r 2D /r 3D = 3ξ (1 + ξ )[2(1 + ξ + ξ 2 )] −1 is always F(ξ ) < 1 for 0 < ξ < 1.Therefore, a 2D droplet in general breaks up faster and more easily than a 3D droplet, if the water sub-droplet diameter is the same.Note this estimation does not consider the asymmetric deformation, dynamic deformation, and vapor leakage.In the actual dynamics with multiple water sub-droplets, it is more complicated to estimate the relation of D b between a 2D case and its counterpart 3D case.In appearance, Case C3-3D seems to have a larger D b than C3, but the number of water sub-droplets is different between the two cases.So, there is no direct relation on D b between the two cases.
Even though the increase of D b does not linearly relate to the number of active water sub-droplets, the above trend indicates that multiple explosions can break up the entire droplet if the water subdroplet size is sufficiently large and the bubble growth occurs simultaneously.This suggests that control of secondary atomization by microexplosion/puffing is possible by properly designing the emulsion fuel mixture, the combustor configuration and its operating conditions.Consequently, the combustion performance and emission characteristics may be also optimized.Experimentally, such an attempt has been made by varying the fuel/water mixing time and thus the size of dispersed water sub-droplets at a certain volume fraction of water. 9is article is copyrighted as indicated in the article.Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP:

103302- 6
FIG. 1.A sucking interface due to boiling.(a) Schematic of the problem and (b) temperature profiles at different times.

FIG. 7 .
FIG. 7. Puffing and after-puffing dynamics of Case A1-3D.(a) t = 1.5 μs and (b) t = 2.5 μs.The gray surface indicates the liquid oil droplet surface.The darker gray and lighter gray (red and pink, respectively) surfaces show the water sub-droplet's boiling surface and inert non-boiling interface attached to the liquid oil, respectively.The contours on the plane that is perpendicular to the z-axis and contains the CG of the oil droplet show the oxygen mass fraction.The darkest gray and lightest gray (blue and yellow, respectively) correspond to water vapor and air, respectively.To clearly show the droplet surface, the contours on the plane are not shown in the regions inside the oil droplet or very close to the boiling surface.

FIG. 8 .
FIG. 8. Temporal traces of (a) the bubble size and (b) the normalized pressure difference between the bubble and ambient air before puffing.

48 TimeFIG. 13 .
FIG. 13.Distance between the CGs of the parent oil droplet and the water sub-droplet.t w corresponds to the time when the water sub-droplet fully detaches in each case.
FIG. 14.Schematic of water sub-droplet dynamics.(a) Model configuration for water sub-droplet dynamics and (b) schematic of solutions.In (b), t w corresponds to the maximal time of each case in Fig. 13.

TABLE I .
Simulation cases and parameters.
a An "active" water sub-droplet means the sub-droplet contains a bubble initiated at the sub-droplet surface (or equivalently the water/oil interface).