# Heat and Mass Transfer of Droplet Vacuum Freezing Process Based on Dynamic Mesh

Authors: Lili Zhao; Yuekai Zhang; Zhijun Zhang; Xun Li; Wenhui Zhang
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

A numerical simulation using dynamic mesh method by COMSOL has been developed to model heat and mass transfer during vacuum freezing by evaporation of a single droplet. The initial droplet diameter, initial droplet temperature, and vacuum chamber pressure effect are studied. The surface and center temperature curve was predicted to show the effect. The mass transfer rate and radius displacement were also calculated. The results show the dynamic mesh shows well the freezing process with the radius reduction of droplet. The initial droplet diameter, initial droplet temperature, and vacuum pressure have obvious effect on freezing process. The total freezing time is about 200 s, 300 s, and 400 s for droplet diameter 7.5 mm, 10.5 mm, and 12.5 mm, respectively. The vacuum pressure less than 200 Pa is enough for the less time to freezing the droplet, that is, the key point in freezing time. The initial droplet temperature has obvious effect on freezing but little effect on freezing temperature.


---

## Body

## 1. Introduction

Droplet vacuum freezing is fast temperature decrease process with the phase change from liquid to vapor and from solid to vapor. It belongs to the evaporative supercooling method for ice production, which is a research hot spot in this field [1–3]. That is the basic principle of the evaporative freezing method: the water evaporates under vacuum conditions. If the vapor pressure of the atmosphere is below 611 Pa (the vapor pressure is around 611 Pa at 0°C), the water will keep on evaporating.The heat and mass transfer is a complicated process, which has been studied by numerical and experiment method. Kim and Satoh et al. have studied the water spray evaporation method for ice particle production by theoretical and experimental [4, 5]. In their study, the conditions for the formation of ice particles were investigated theoretically by the diffusion-controlled evaporation model. The experiment was with pure water and 7% ethylene glycol. Li et al. have studied the cooling/freezing phenomena of a water droplet due to evaporation in an evacuated chamber by experiment and the heat transfer dominating the evaporation-freezing phenomena in order to estimate the pressure in the evaporator [6]. Asaoka et al. [7, 8] also proposed an alternative vacuum freezing method to make ice slurry by evaporating ethanol solution instead of water. This system improves the ice forming stability but still needs large amount of energy for vacuum. The model and simulation above used the analytical model; the finite element method is also used to examine the process [9]. It is based on the diffusion model with the porous medium; radius decrease by the water evaporating is not considered.In this paper, we developed a model using dynamic mesh method to follow the tracks of droplet boundary. The diameter, initial temperature, and vacuum chamber pressure effect are studied.


## 2. Model and Validation

### 2.1. Material

A physical two-dimensional axis symmetry model that explains the vacuum freezing process of droplet is shown in Figure1. The property of the water is shown in Table 1.Table 1
Parameters used in the simulation process.
Parameter
Value
Unit
ρ
998
Kg m−3
c
l
3420
J kg−1 K−1
c
i
1720
J kg−1 K−1
k
i
0.6
W m−1 K−1
k
l
2
W m−1 K−1
λ
i
2.8 × 106
J kg−1
λ
l
2.5 × 106
J kg−1Figure 1
2D axis symmetry model of droplet.


### 2.2. Heat Transfer Model

The value of any physical quantity at a point in space is given by its average value on the averaging volume centered at this point.By considering the hypothesis of the local thermal equilibrium, the energy conservation is reduced to a unique equation:(1)
ρ
c
∂
T
∂
t
-
∇
(
k
·
∇
T
)
=
0
.Boundary condition for (1) on the symmetric surface is
(2)
k
∇
T
=
0
.Boundary condition for (1) on the outer surface is
(3)
k
∂
T
∂
n
=
-
λ
·
m
˙
s
.Initial condition for (1) is
(4)
T
-
=
T
0
.


### 2.3. Mass Transfer Coefficient

The phase change rate of water is difficultly decided by experiment method. The rate constant parameter has the dimension of reciprocal time in which phase change occurs. A large value of signifies that phase change occurs in a small time. The mass transfer could be as,(5)
m
˙
s
=
K
(
P
sat
⁡
(
T
s
)
-
P
v
)
.


### 2.4. Resolution and Dynamic Mesh

COMSOL Multiphysics 3.5a was used to solve the set of equations. COMSOL is 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 relative initial and boundary conditions of each. COMSOL offers three possibilities for writing the equations:(
1
) using a template (Fick Law, Fourier 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 an unsymmetrical-pattern multifrontal method. We used a direct solver for sparse matrices (UMFPACK), which involves significantly more complicated algorithms than solvers used for dense matrices. The main complication is the need to handle the fill-in in factors L and U efficiently.A two-dimensional (2D) axis symmetry grid was used to solve the equations using COMSOL Multiphysics 3.5a. The mesh consists of 1068 elements (2D), and time stepping is freely taken by solver. Several grid sensitivity tests were conducted to determine the sufficiency of the mesh scheme and to ensure that the results are grid independent. A backward differentiation formula was used to solve time-dependent variables. Relative tolerance was set to1
×
10
-
4, whereas absolute tolerance was set to 1
×
10
-
6. The simulations were performed using a Tongfang PC with Intel Core 2 Duo processor with 3.0 GHz processing speed and 4096 MB of RAM running Windows 7.A deformed mesh can be useful if the boundaries of your computational domain are moving in time or are a function of a parameter. The point is that a new mesh needs not be generated for each configuration of the boundaries—instead, the software simply perturbs the mesh nodes so they conform with the moved boundaries.In COMSOL Multiphysics, the movement of the interior nodes is in three ways:(
1
) by propagating the moving boundary displacement throughout the domain to obtain a smooth mesh deformation everywhere—this is done by solving PDEs for the mesh displacements (a Laplace or Winslow smoothing PDE) with boundary conditions given by the movement of the boundaries; (
2
) by specifying an explicit formula for the mesh deformation; (
3
) by letting the mesh movement be determined by some physical deformation variables, such as the displacement components of structural mechanics.


### 2.5. Validation

Figures2, 3, and 4 are the simulation results by the parameters of Table 1. It is compared with the experiment result from [9]. The simulation results are in good agreement with the experiment results that proved the model and method explain well the heat and mass transfer process. The details of discussion will be given later.Figure 2
The simulaiton results and experiment results (from [9]) of droplet diameter 7.5 mm (P
v
=
200 Pa, T
0
=
297.15).Figure 3
The simulaiton results and experiment results (from [9]) of droplet diameter 10.5 mm (P
v
=
200 Pa, T
0
=
297.15).Figure 4
The simulaiton results and experiment results (from [9]) of droplet diameter 12.5 mm (P
v
=
200 Pa, T
0
=
297.15).


## 2.1. Material

A physical two-dimensional axis symmetry model that explains the vacuum freezing process of droplet is shown in Figure1. The property of the water is shown in Table 1.Table 1
Parameters used in the simulation process.
Parameter
Value
Unit
ρ
998
Kg m−3
c
l
3420
J kg−1 K−1
c
i
1720
J kg−1 K−1
k
i
0.6
W m−1 K−1
k
l
2
W m−1 K−1
λ
i
2.8 × 106
J kg−1
λ
l
2.5 × 106
J kg−1Figure 1
2D axis symmetry model of droplet.


## 2.2. Heat Transfer Model

The value of any physical quantity at a point in space is given by its average value on the averaging volume centered at this point.By considering the hypothesis of the local thermal equilibrium, the energy conservation is reduced to a unique equation:(1)
ρ
c
∂
T
∂
t
-
∇
(
k
·
∇
T
)
=
0
.Boundary condition for (1) on the symmetric surface is
(2)
k
∇
T
=
0
.Boundary condition for (1) on the outer surface is
(3)
k
∂
T
∂
n
=
-
λ
·
m
˙
s
.Initial condition for (1) is
(4)
T
-
=
T
0
.


## 2.3. Mass Transfer Coefficient

The phase change rate of water is difficultly decided by experiment method. The rate constant parameter has the dimension of reciprocal time in which phase change occurs. A large value of signifies that phase change occurs in a small time. The mass transfer could be as,(5)
m
˙
s
=
K
(
P
sat
⁡
(
T
s
)
-
P
v
)
.


## 2.4. Resolution and Dynamic Mesh

COMSOL Multiphysics 3.5a was used to solve the set of equations. COMSOL is 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 relative initial and boundary conditions of each. COMSOL offers three possibilities for writing the equations:(
1
) using a template (Fick Law, Fourier 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 an unsymmetrical-pattern multifrontal method. We used a direct solver for sparse matrices (UMFPACK), which involves significantly more complicated algorithms than solvers used for dense matrices. The main complication is the need to handle the fill-in in factors L and U efficiently.A two-dimensional (2D) axis symmetry grid was used to solve the equations using COMSOL Multiphysics 3.5a. The mesh consists of 1068 elements (2D), and time stepping is freely taken by solver. Several grid sensitivity tests were conducted to determine the sufficiency of the mesh scheme and to ensure that the results are grid independent. A backward differentiation formula was used to solve time-dependent variables. Relative tolerance was set to1
×
10
-
4, whereas absolute tolerance was set to 1
×
10
-
6. The simulations were performed using a Tongfang PC with Intel Core 2 Duo processor with 3.0 GHz processing speed and 4096 MB of RAM running Windows 7.A deformed mesh can be useful if the boundaries of your computational domain are moving in time or are a function of a parameter. The point is that a new mesh needs not be generated for each configuration of the boundaries—instead, the software simply perturbs the mesh nodes so they conform with the moved boundaries.In COMSOL Multiphysics, the movement of the interior nodes is in three ways:(
1
) by propagating the moving boundary displacement throughout the domain to obtain a smooth mesh deformation everywhere—this is done by solving PDEs for the mesh displacements (a Laplace or Winslow smoothing PDE) with boundary conditions given by the movement of the boundaries; (
2
) by specifying an explicit formula for the mesh deformation; (
3
) by letting the mesh movement be determined by some physical deformation variables, such as the displacement components of structural mechanics.


## 2.5. Validation

Figures2, 3, and 4 are the simulation results by the parameters of Table 1. It is compared with the experiment result from [9]. The simulation results are in good agreement with the experiment results that proved the model and method explain well the heat and mass transfer process. The details of discussion will be given later.Figure 2
The simulaiton results and experiment results (from [9]) of droplet diameter 7.5 mm (P
v
=
200 Pa, T
0
=
297.15).Figure 3
The simulaiton results and experiment results (from [9]) of droplet diameter 10.5 mm (P
v
=
200 Pa, T
0
=
297.15).Figure 4
The simulaiton results and experiment results (from [9]) of droplet diameter 12.5 mm (P
v
=
200 Pa, T
0
=
297.15).


## 3. Results and Discussion

At first, the diameter effect was studied in Figures2, 3, and 4. The surface temperature and center temperature of the droplet were given. Both have the quick temperature decrease at the initial stage, called first cooling stage. And then it will get the second stable freezing stage, because the phase change and latent of phase changed. And then in the third frozen stage, the temperature will be decreased quickly. The experiment results have some differences that are not discussed here because it is a large title.The diameter of dropelt increases from 7.5 mm to 12.5 mm, and the freezing time is about 200 s, 300 s, and 400 s. It is shown that the diameter has larger effect on freezing time. Because, most of the time, we want the less time to freeze the droplet, the less diameter droplet is needed. The surface time temperature of droplet was lower than the center temperature of it most of the time because the freezing power is the surface vapor and sublimation. Both temperatures are same in the end.Figure5 is the surface total mass transfer rate of droplet that is decided by (5). The parameters include mass transfer coefficent, saturated vapor pressure of droplet surface, and vacuum chamber pressure. The rate is larger at the initial freezing stage, slowly decreased, and then reached zero. It is corresponding to the temperature curve of the droplet surface and center. The radius displacement is calculated with the moving mesh method in Figure 6. The more the diameter droplet is, the more the radius displacement is.Figure 5
Mass transfer rate of different diameters (d
=
7.5 mm, 10.5 mm, and 12.5 mm, P
v
=
200 Pa, and T
0
=
297.15).Figure 6
Radius displacement of different diameters (d
=
7.5 mm, 10.5 mm, and 12.5 mm, P
v
=
200 Pa, and T
0
=
297.15).As the important effect factor, the pressure of vacuum chamber is changed from 100 to 500 Pa. The results were shown in Figure7. The pressure is lower than 300 Pa; the temperature curve has the obvious three stages discussed above. But when pressure is higher than 400 Pa, the temperature just has the two stages. And compared with the temperature curve of 100 Pa and 200 Pa, the freezing time to the end time is almost the same. But compared with the temperature curve of 200 Pa and 300 Pa, the freezing time is less than 200 Pa. This is very meaningful for freezing water into ice in order for energy storage.Figure 7
Center temperature of different vacuum chamber pressure (T
0
=
297.15 K, d
=
7.5 mm).Figures8 and 9 are radius displacement and total mass transfer rate with different vacuum chamber pressure. The mass transfer rate has some different characteristic than normal. The higher the vacuum chamber pressure is, the higher the initial rate is. But the later stage rate of mass transfer rate is contrary. The detailed reason should be as follows: because the surface freezing is earlier at lower pressure.Figure 8
Radius displacement of different vacuum chamber pressure (T
0
=
297.15 K, d
=
7.5 mm).Figure 9
Mass transfer rate of different vacuum chamber pressure (T
0
=
297.15 K, d
=
7.5 mm).Figures10, 11, and 12 are the different initial temperature results. From Figure 10, the initial temperature will obviously affect the first cooling stage. But little is the third frozen stage. The radius displacement and total mass transfer rate are affected with the temperature increased.Figure 10
Center temperature of different initial temperature (P
v
=
200 Pa, d
=
7.5 mm).Figure 11
Radius displacement of different initial temperature (P
v
=
200 Pa, d
=
7.5 mm).Figure 12
Mass transfer rate of different initial temperature (P
v
=
200 Pa, d
=
7.5 mm).


## 4. Conclusion

The initial droplet diameter, initial droplet temperature, and vacuum chamber pressure effect are studied. The surface and center temperature curve were predicted to show the effect. The mass transfer rate and radius displacement were also calculated. The results showing the dynamic mesh well show the freezing process with the radius reduction of droplet. The initial droplet diameter, initial droplet temperature, and vacuum pressure have the obviously effect on freezing process. Usually droplet vacuum freezing has the quick temperature decrease at the initial stage, called first cooling stage. And then it will get the second stable freezing stage, because the phase change and latent of phase changed. And then in the third frozen stage, the temperature will be decreased quickly. The total freezing time is about 200 s, 300 s, and 400 s for droplet diameters 7.5 mm, 10.5 mm, and 12.5 mm, respectively. The vacuum pressure less than 200 Pa is enough for the less time to freeze the droplet that is the key point in freezing time. The initial droplet temperature has the obvious effect on freezing but little effect on freezing temperature.



---

*Source: 798040-2014-07-20.xml*