Parameter Sensitivity of the Microdroplet Vacuum Freezing Process

The vacuum freezing process of microdroplets (<100μm in diameter) is studied by dynamic mesh method. The mass transfer coefficient was studied using the results of related papers that considered droplet diameters exceeding 1mm. The diameter, initial temperature, and vacuum chamber pressure effects are also discussed. To estimate parameter sensitivity, the effects of material density, specific heat, and thermal conductivity in 20% scope, as well as latent evaporation/sublimation in 5%, were simulated. The results show that the mass transfer coefficient K is essentially different between microdroplets (<100μm) and macrodroplet (>1mm). Pressure and droplet diameter have an effect on cooling and freezing stages, but initial temperature only affects the cooling stage. The thermal conductivity coefficient kl affected the cooling stage, whereas ki affected the freezing stage. Heat capacity Cl affected the cooling stage, but Ci has virtually no effect on all stages. The actual latent heat of freezing ΔH was also affected. Higher density corresponds to lower cooling rate in the cooling stage.


Introduction
The use of evaporative supercooling method for ice production presents an immensely popular topic for research [1,2].Droplet vacuum freezing involves a rapid decrease in temperature in which the phase changes from liquid to solid, from liquid to vapor, and from solid to vapor [3].The basic principle is as follows.Water evaporates under vacuum conditions.If the vapor pressure of the atmosphere is below 611 Pa (vapor pressure at 0 ∘ C), then the water continues to evaporate.Heat must be provided during evaporation to decrease droplet temperature.Thus, the temperature reaches below 0 ∘ C and then ice is produced.
The heat and mass transfer of droplet vacuum freezing are complicated processes and have been studied by numerical and experiment methods.Kim et al. [4,5] conducted theoretical and experimental studies on the use of water spray evaporation method to produce ice particles.In their study, the conditions for the formation of ice particles were theoretically investigated using a diffusion-controlled evaporation model.The experiment included pure water and 7% ethylene glycol.Li et al. [6] experimentally studied the cooling/freezing phenomena of a water droplet due to evaporation in an evacuated chamber, as well as the heat transfer dominating the evaporation/freezing phenomena, to estimate the pressure in the evaporator.Asaoka et al. [7,8] also proposed an alternative vacuum freezing method to produce ice slurry by evaporating ethanol solution instead of water.This system improves ice forming stability but requires a large amount of energy for vacuum production.These models and simulations have been used as analytical model, and finite element method has been used to examine the process [9].This method is based on the diffusion model through a porous medium, and radius reduction attributed to water evaporation is disregarded.We have studied the droplet vacuum freezing process using the dynamic mesh method [10], in which the droplet diameter was reduced during the freezing process.The diameter of droplet was set as 7.5, 10.5, and 12.5 mm for comparison with experiment results [9].To continuously produce ice from water droplets, the diameter of the droplet should be less than 100 m because of the height limitation of the spray vacuum chamber [3].
In the present study, a microdroplet (<100 m diameter, compared with the conventional >1 mm) vacuum freezing process with dynamic mesh method is studied.The mass transfer coefficient was studied along with the results of related studies [3,10].The diameter, initial temperature, and vacuum chamber pressure effect were also studied.To estimate parameter sensitivity, the effect of material density, specific heat, thermal conductivity, and latent heat of evaporation was simulated in 20% or 5% scope.

Material and Method
2.1.Material.The microdroplet is spherically shaped and its diameter is reduced during vacuum freezing process because of surface evaporation.To simplify computations, the two-dimensional axis symmetry model (Figure 1) was adopted.The solid and liquid properties of the microdroplet are provided in Table 1.To estimate sensitivity, parameters including ,   ,   ,   , and   were changed by ±20%.By contrast, the parameters   and   were just changed by ±5%.These adjustments are adopted because   must be larger than   and these two parameters are highly sensitive.These reasons will be discussed later.

Heat Transfer Model.
By considering a hypothesis of local thermal equilibrium, energy conservation is reduced to a unique equation: The boundary conditions for (1) on the axis symmetric surface are as follows: The boundary conditions for (1) on the outer surface are defined as The initial condition for (1) is Outside surface Axis symmetric

Mass Transfer Model.
Mass transfer occurs through a phase change of water from liquid or solid to gas on the droplet surface.However, mass transfer is not considered inside the droplet.A general expression of nonequilibrium evaporation rate used for phase change modeling in porous medium is given by [11,12] İ =  ( V,eq −  V ) .
This equation could also be expressed as [13] In vacuum cooling, the equation could be expressed as [14][15][16][17][18] Here, phase change occurs based on the movement of the pressure difference between the saturation vapor pressure in the microdroplet surface and vacuum chamber pressure [11], which is also considered the mass transfer boundary condition: The mass transfer coefficient, , is a rate constant parameter with the dimension of reciprocal time in which phase change occurs.A large value of  signifies that phase change occurs within a short amount of time.The  value is conventionally obtained from experiments.In the present paper, we decided to use experimental data from [3].Detailed information will be provided later.
Water partial pressure,  V , is assumed to be equal to the total vacuum chamber pressure.The saturation pressure,  sat , is related to the surface temperature of the microdroplet.This value is estimated using the equation when the surface temperature is higher than freezing temperature [9] but uses another equation for estimation when the surface temperature is lower than freezing temperature [9].The details could be found in related papers [9,19].
The initial condition for ( 5) is set using the initial diameter of the droplet: Thus, as in the dynamic mesh method, the boundary conditions of ( 8) could be rewritten as 2.4.Resolution and Validation.COMSOL Metaphysics 3.5a was used to solve the set of equations.COMSOL is an advanced software used for modeling and simulating any physical process described by partial derivative equations.The set of equations introduced above was solved using the corresponding relative initial and boundary conditions.COMSOL offers three possibilities for writing the equations: (1) using a template (Fick's law or Fourier's law), (2) using the coefficient form (for mildly nonlinear problems), and (3) using the general form (for most nonlinear problems).Differential equations in the coefficient form were written using a nonsymmetrical pattern multifrontal method.We used a direct solver for sparse matrices (UMFPACK), which involves significantly more complicated algorithms than those solvers used for dense matrices.The need for efficient handling of the fill-in factors  and  is considered the method's major complication.
A two-dimensional (2D) axis symmetry grid was used to solve the equations using COMSOL Metaphysics 3.5a.The mesh consists of 1068 elements (2D), and time stepping is freely taken by the solver.A backward differentiation formula was used to solve time-dependent variables.Relative tolerance was set at 1 × 10 −4 , whereas absolute tolerance was set at 1 × 10 −6 .The simulations were performed using a Tongfang PC running on Windows 7, with an Intel Core 2 Duo processor clocked at 3.0 GHz processing speed and 4096 MB of RAM.
Several grid sensitivity tests were conducted to determine the sufficiency of the mesh scheme and to ensure that the results are grid independent.Figures 2 and 3 present the grid test results.The vacuum chamber pressure would be higher or lower than 600 Pa, depending on the different physical process.In the transition involved in Figure 2, liquid water is changed to vapor.In Figure 3, liquid water changes into solid ice, and then the conversion from ice to vapor is attributable to the physical property difference.The test results in Figures  2 and 3 show that all of the three elements yielded accurate results, although certain numerical results have step-changed.Simulation times of 0.5, 1, and 3 h presented 267, 1068, and 4602 elements in Figure 3.In the following simulation, 1068 elements were used.
As presented above, the mass transfer coefficient or rate constant  should be obtained in the experiment.In the present simulation, the experimental data from [3] was used for decision-making.In Figure 4, different  values were used, but the results do not show the difference when  exceeds 0.1 because, at  ≥ 0.1, mass and heat transfer is limited by inner microdroplet process.In addition, the Mathematical Problems in Engineering  obtained simulation results were used to approximate the experimental data.The calculation results in the reference are also shown in Figure 4.The bias between the experimental present simulation results and reference results [3] is attributed to the fact that the microdroplet is not spherical during vacuum freezing.The droplet surface is in the boiling state, and thus the increased surface area results in larger mass and heat transfer than in the simulation.Certain information on the inner controlled results could be obtained later.
In our previous study [10], the rate constant  is 4 × 10 −6 , which is considerably lower than 0.1 in the present study.The former diameter of droplet is 7.5 mm or larger.However, the later diameter is less than 100 m.The droplet dimension was reduced, thereby altering the mass and heat transfer process.

Results and Discussion
The process of microdroplet vacuum freezing includes three stages: cooling, freezing with phase change, and frozen stages.In the cooling stage, the temperature of the microdroplet is lowered to the phase change temperature from the initial temperature.In the freezing stage, the liquid water is changed to solid ice, and the latent heat of the phase change is released.The temperature was maintained at a stable value.When all liquid water is turned into solid ice, the temperature of the microdroplet is rapidly decreased to the saturation temperature corresponding to the vacuum chamber pressure.This state is the frozen stage.
When the external conditions, such as pressure, initial temperature, and microdroplet diameter, or internal conditions, such as the property of material, were changed, the three stages of microdroplet vacuum freezing are also changed.The pressure is the only condition that could be controlled in the vacuum freezing process.The saturation temperature corresponds to the pressure, with lower pressure corresponding to lower temperature and vice versa.However, the pressure not only affects the end temperature but also affects the process.When the pressure is lowered, the freezing time is evidently shortened (Figure 5).Except at 600 and 500 Pa, the three stages presented above were all altered.The phase change is a process.The change from liquid water to solid ice occurs at a range of temperature.Pressures at 600 and 500 Pa are located near the triple-phase temperature, and thus the process occurs differently for other pressures.In the present simulation, the temperature scope of the phase change is 2 K from triple-phase temperature point.Results show that the freezing time is evidently different when the pressure is lower than 500 Pa.The three stages are affected by the pressure in the process.However, the initial temperature effect is different from the pressure effect in Figure 6.The initial temperature mainly affects the cooling stage, almost exerting no effect on the freezing and frozen stages.Based on the results, if pressure is not changed, then the end temperature of the microdroplet remains the same regardless of the initial temperature of the material or microdroplet.
Similar to pressure and initial temperature differences, the initial diameter of microdroplet has a considerable effect on the freezing process.The diameter of microdroplet ranges from 10 m to 100 m and is a one-order dimension.However, the freezing time exceeds two orders in Figure 7.In reality, a lower initial diameter results in flash evaporation [5] and a decreasing diameter increases this probability.
Unlike rate constant , the process is externally controlled.When  > 0.1, the process is internally controlled.These changes indicate process complexity.
Unlike external conditons, material property is related to the temperature.However, the accuracy of the property   parameters is the basis for the validity of simulation and models.Thermal conductivity is related to temperature and varies with the process.The manner by which this parameter affects the process should be presented.Thermal conductivity is related to both temperature and physical state, varying in the liquid and the solid states.To estimate parameter sensitivity, the value of parameter is changed in the 20% scope, except for latent heat of phase change.The results show the evident effect of thermal conductivity   and   , which are different between Figures 8 and 9. Larger   results in faster cooling rate during the cooling stage.By contrast, the frozen   stage is less affected.The value of   is less affected in the cooling stage and is evidently affected in the freezing stage.Larger   results in decreased freezing time.
As seen in Figures 10 and 11, the heat capacities   and   have different effects on the freezing process.  almost has no effect on the freezing process.By contrast,   evidently affects the first stage (or the cooling stage).Larger   results in a higher cooling rate.However, the end time is nearly similar.
As the main parameter, the latent heats   and   present an evident effect on all processes, as seen in Figures 12 and 13.However, their effects are different on the freezing process.In fact, in the second stage, namely, the freezing stage, liquid  water is changed into solid microdroplets.The released latent heat of freezing, Δ, is calculated as By increasing   , the Δ value decreases.Given that the latent heat of freezing Δ is decreased, cooling and freezing time decreased (Figure 12) and vice versa.The latent heat of sublimation   also varies.Increasing   results in the increase in latent heat of freezing Δ, and then the cooling rate in cooling stage is decreased.Thus, freezing time is also increased because of the increased Δ and vice versa.In the heat loss in (3), given the increase in latent heat values   and   , the cooling rates in the cooling stage shown in Figures 12 and 13 are both increased.As discussed above, the rate constant  is fixed, and the process is internally controlled.Here, (3) is unaffected by the heat transfer process by the changed   and   , obtained as the highest value as in  > 0.1.These heat values represent the latent heat of freezing Δ.The latent heat values   and   changed in the 5% scope because of higher sensitivity.
Finally, density is discussed.In the phase change from liquid water to solid ice, density changes.As shown  in Figure 14, higher density is correlated with lower cooling rate in the cooling stage, and freezing time is consequently shortened.This trend is true vice versa.

Conclusion
The mass transfer on the droplet surface varies between microdroplets and macrodroplets, as represented by the mass transfer coefficient .The pressure of the vacuum chamber and microdroplet diameter affect the cooling and freezing stages.However, the initial temperature of the microdroplet only affects the cooling stage.The thermal conductivity coefficient   affected the cooling stage, whereas   affected the freezing stage.Heat capacity   affected the cooling stage, but   almost has no effect in all stages.Latent heat of   and   affected cooling and freezing stages, but their effects are different.Both latent heat values affected the latent heat of freezing Δ.Finally, higher density results in decreased cooling rate during the cooling stage.
results with K = 10 Simulation results with K = 1 Simulation results with K = 0.1 Simulation results with K = 0.01 Simulation results with K = 0.001 Dropt center temperature (K)Experiment result from[3] Simulation results from[3]

Figure 8 :
Figure 8: Effect of thermal conductivity   at an initial temperature of 286 K, vacuum chamber pressure 300 Pa, and microdroplet initial diameter of 50 m.

Figure 9 :
Figure 9: Effect of thermal conductivity   at initial temperature of 286 K, vacuum chamber pressure of 300 Pa, and microdroplet initial diameter of 50 m.

Figure 10 :Figure 11 :
Figure 10: Effect of heat capacity   at initial temperature of 286 K, vacuum chamber pressure of 300 Pa, and microdroplet initial diameter of 50 m.

Figure 13 :
Figure 13: Effect of latent heat   at an initial temperature of 286 K, vacuum chamber pressure 300 Pa, and microdroplet initial diameter of 50 m.

Figure 14 :
Figure 14: Effect of density  at an initial temperature of 286 K, vacuum chamber pressure of 300 Pa, and microdroplet initial diameter of 50 m.

Table 1 :
Parameters used in the simulation process.
Figure 2: Grid sensitivity tests with 267, 1068, and 4272 elements at 600 Pa vacuum chamber pressure.
Figure 12: Effect of latent heat   at an initial temperature of 286 K, vacuum chamber pressure of 300 Pa, and microdroplet initial diameter of 50 m.