Porous medium vacuum drying is a complicated heat and mass transfer process. Based on the theory of heat and mass transfer, a coupled model for the porous medium vacuum drying process is constructed. The model is implemented and solved using COMSOL software. The water evaporation rate is determined using a nonequilibrium method with the rate constant parameter values of 1, 10, 1000, and 10000 are simulated. The effects of vapor pressures of 1000, 5000, and 9000 Pa; initial moistures of 0.6, 0.5, and 0.4 water saturation; heat temperatures of 323, 333, and 343 K; and intrinsic permeability of 10−13, 10−14, and 10−15 m2 are studied. The results facilitate a better understanding of the porous medium vacuum drying process.

1. Introduction

Vacuum drying is an excellent drying method for vegetables and fruits, among others. This approach is a low-temperature, nonpolluting method that produces good results. However, vacuum drying requires a complicated device and entails high costs. Scientists and engineers are currently studying vacuum drying equipment that could be used in corn drying [13]. However, the corn vacuum drying theory remains unclear. Given that corn is a porous medium, the vacuum drying of corn is a complicated heat and mass transfer process that has been the subject of intensive research [47]. All vacuum drying models have to address the water phase change during numerical solving. In one method, the vapor pressure is equal to its equilibrium value [811]. Another method is nonequilibrium method [1216].

Erriguible studied convective and vacuum drying and identified the couple problem of heat and mass transfer between the inner and the outer porous medium [8, 9]. The numerical code applied to the porous medium and the computational fluid dynamics (CFD) Software (FLUENT) was obtained. Murugesan presented the same problem [10].

Perré and Turner used a dual-scale modeling approach to describe the coupling of the drier (large-scale) and the porous medium (macroscale) throughout the drying process [11]. The model was used to investigate the vacuum drying of a softwood board placed in an experimental vacuum chamber heated by two infrared emitters.

Torres et al. proposed a coupled model to describe the vacuum drying of oak wood at the laboratory scale [12, 13]. This model describes the physics of wood-water relations and interactions with a vacuum dryer. The results provided important information on liquid and gas phase transport in wood.

Warning implemented a multiphase porous media model involving heat and mass transfer within a potato chip in a commercial CFD program [14]. The simulations were run at different frying pressures of 1.33, 9.89, 16.7, and 101 kPa. Good agreement between the predicted and literature experimental moisture, oil, and acrylamide content was achieved. Regardless of fryer pressure, the model showed that the core pressure increased and became approximately 40 kPa higher than the surface. The model modified Darcy’s law to account for the Klinkenberg effect.

Halder et al. developed an improved multiphase porous media model involving heat and mass transfer with careful consideration given to the selection of input parameters [15, 16]. The nonequilibrium formulation for evaporation, which provides a better description of the physics and is easier to implement in a typical CFD software, is used because it can explicitly express the evaporation rate in terms of the concentration of vapor and temperature. External heat transfer and mass transfer coefficients are estimated to reflect the different frying phases accurately, that is, the nonboiling phase and surface boiling and falling rate stages in the boiling phase.

As noted by Halder et al., water evaporation during frying or during other drying-like processes implemented using an equilibrium formulation may not always occur [14, 15]. The equations resulting from an equilibrium formulation cannot be implemented in any direct manner in the framework of most commercial software. A nonequilibrium formulation provides a better description of the physics and is also easier to implement in software, thus appearing to be the obvious alternative. Additionally, details on the nonequilibrium model have not been discussed thoroughly.

In this paper, heat and mass transfer of porous medium in the vacuum drying process is implemented by using a nonequilibrium method. First, different phase change rates were studied to understand their effect on the drying process. The effects of vapor pressure, initial moisture, heat temperature, and intrinsic permeability on the drying process were then examined.

2. Physical Model

A physical one-dimensional (1D) model that explains the drying process is shown in Figure 1. The bottom of the porous medium is heated by a hot plate. The top of the porous medium is subjected to gas pressure. The height of the porous medium is 1 cm.

3. Mathematical Model

The porous medium consists of a continuous rigid solid phase, an incompressible liquid phase (free water), and a continuous gas phase that is assumed to be a perfect mixture of vapor and dry air, considered as ideal gases. For a mathematical description of the transport phenomenon in a porous medium, we adopt a continuum approach, wherein macroscopic partial differential equations are achieved through the volume averaging of the microscopic conservation laws. 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.

The moisture movement of the inner porous medium is liquid water and vapor movement; that is, the liquid water could become vapor, and the vapor and liquid water are moved by the pressure gradient.

The compressibility effects of the liquid phase are negligible, and the phase is homogeneous:

The solid phase is rigid and homogeneous:

The gaseous phase is considered an ideal gas. This phase ensures that

The assumption of the local thermal equilibrium between the solid, gas, and liquid phases involves

Mass conservation equations are written for each component in each phase. Given that the solid phase is rigid, the following is given:

The averaged mass conservation of the dry air yields

For vapor,

For free water,

For water, the general equation of mass conservation is obtained from the sum of the conservation equations of vapor () and free water (). The general equation is written as follows:

For the Darcy flow of vapor,

For the Darcy flow of air, where the gas and free water velocity is given by The effective diffusion coefficient [8] is given by

The vapor fraction in mixed gas is given by The pressure moving the free water is given by For capillary pressure, The saturation of free water and gas is Free water relative permeability is given by Gas relative permeability is given by The water phase change rate is expressed as Water saturation vapor pressure is given by

By considering the hypothesis of the local thermal equilibrium, the energy conservation is reduced to a unique equation:

4. Boundary Condition and Parameters

The air pressure on the external surface at the top of the porous medium is fixed, and the boundary condition for air is given by

The boundary condition for vapor at the top of the porous medium is given by

To simulate the vapor pressure of the vacuum drying chamber effect on the drying process, four different vapor pressure boundary values are used.

The boundary condition for free water at the top of the porous medium is

The boundary condition for energy at the top of the porous medium is

The boundary condition at the bottom of the porous medium is

Three different values are used in the simulation.

The initial moisture of the porous medium is represented by the liquid water saturation; different initial water saturation values are used. To compare the effects, drying base moisture content was also used, as shown in (3.9). The water phase change rate is studied using four different rate constant parameter values. The modeling parameters are shown in Table 1.

5. Numerical Solution

COMSOL Multiphysics 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 relative initial and boundary conditions of each. COMSOL offers three possibilities for writing the equations: (1) using a template (Fick’s Law, 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 an unsymmetric-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) grid was used to solve the equations using COMSOL Multiphysics 3.5a. Given the symmetry condition setting at the left and the right sides, the 1D model shown in Figure 1 was, in fact, the model that was applied. The mesh consists of 4 × 100 elements (2D), and time stepping is 1 (0 s to 200 s of solution), 10 (200 s to 100000 s of solution), and 100 (100000 s to 500000 s of solution). Several grid sensitivity tests were conducted to determine the sufficiency of the mesh scheme and to ensure that the results are grid-independent. The maximum element size was established as 1e−4. A backward differentiation formula was used to solve time-dependent variables. Relative tolerance was set to 1e−3, whereas absolute tolerance was set to 1e−4. The simulations were performed using a Lenovo Thinkpad X200 with Intel Core 2 Duo processor with 2.4 GHz processing speed, and 2048 MB of RAM running Windows XP.

6. Results and Discussion

6.1. Effect of Phase Change Rate

The phase change rate of water could not be determined using any method for porous medium drying [14, 15]. The rate constant parameter has a dimension of reciprocal time in which phase change occurs. A large value signifies that phase change occurs within a small time frame. For assumption of equilibrium, either is infinitely large or phase change occurs instantaneously. A very high value, however, makes the convergence of the numerical solution difficult.

In the simulations, is set as 1, 10, 100, 1000, and 10000. The other parameters are = 0.6, = 10−13, = 1000 Pa, and = 323 K. However, when is 100, the numerical solution is not convergent even the time step is reduced and the grid is refined, the reason for which is unknown. The results of the other values are shown in Figures 2 to 8. Figure 2 shows the moisture curves of different values. Moisture is obviously affected by the value. When is set as 1, the drying process is longer (approximately 9 h). However, when ≥ 10, the total drying time remains almost the same (approximately 6 h). Under the vacuum conditions, the free water evaporated easily because it was boiling, which resulted in a faster drying rate. For the quick drying process, a higher value ( > 100) is typically adopted [13].

Figures 3 and 4 show the temperature curves at = 5, 7.5, and 10 mm at different values. The temperature is increased rapidly at the start of drying and is then lowered gradually. As the drying process continues, the temperature is increased slowly. The temperature is rapidly increased when the free water has evaporated. When the drying is nearly finished, the temperature remains unchanged. The end temperatures at the same final position are the same for all values. At = 5 mm, the end temperature is approximately 312 K; at = 7.5 mm, the end temperature is approximately 315 K; and at = 10 mm, the end temperature is approximately 317K. The simulation results do not coincide with those of the equilibrium method, wherein temperature is increased at the beginning, remained unchanged throughout the process, and then increased slowly until a steady value was obtained [8].

The temperature change is near to same as the value increases. That is especially noticeable when comparing = 1 with = 1000 and 10000. The temperature curve is almost the same for = 1000 and 10000.

Figures 5 to 8 show the moisture change curves of different values along the direction at 0.5, 1, 2, 3, 4, 5, 6, 7, 8, 9, and 10 h. The moisture curve of = 1 obviously differs from the curve of = 10, 1000, and 10000. Throughout most of the drying time, the moisture of = 1 is lower near the heat position ( = 0), and the moisture is higher farther away from the heat position ( = 0.01 m). Until the latter part of the drying process (7 h), the middle section near the top has higher moisture because of the relatively low free water evaporation rate and the relatively large free water movement. Pressure gradient moves the moisture from the bottom to the top. However, the moisture values at , 1000, and 10000 are higher in the middle near the top throughout most of the drying time; the moisture change appears in Figure 6 as a “Ω” symbol.

The moisture curves of = 1000 and 10000 are almost the same. However, the curves are not as smooth as that of = 1 because of the higher free water evaporation rate. This issue could be resolved by reducing the time step in the resolving process. The moisture is greater than the initial moisture at the initial drying time because of the free water movement. As shown in Figures 2 to 8, the simulation results of = 1000 and 10000 are almost the same. = 1000 is adequate for the quick drying process [14, 15].

6.2. Effect of Vapor Pressure in Vacuum Drying Chamber

The pressure of a vacuum drying chamber, especially vapor pressure, plays an important role in the vacuum drying process and is also linked to the drying cost. The moisture curves of = 1000, 5000, and 9000 Pa are shown in Figure 9. The other simulation parameters are = 0.6, = 10−13, = 1000, and = 323 K. The vapor pressure has a greater effect on the drying process; a lower vapor pressure results in greater pressure degradation. The drying times are approximately 6, 8, and 12 h. The movements of free water and vapor, as well as the free water evaporation rate, are quicker, as given by (3.13), and (3.21), respectively.

The different vapor pressure effects on the temperature curve are shown in Figure 10. Compared with the temperature curve at = 1000 Pa, no significant change was observed in the temperature curve at = 5000 and 9000 Pa. The temperature is increased at the start of drying, then it is lowered which is insignificant, especially at = 9000 Pa. The moisture changes at = 5000 and 9000 Pa are shown in Figures 11 and 12, respectively, as compared with the moisture change at = 1000 Pa (Figure 7).

The curve was smoother because the increasing vapor pressure lowers the water evaporation rate, as shown in (3.21). For = 9000 Pa in particular, the maximum moisture value appears at = 0.01m in the drying process. The moisture curve does not appear as the “Ω” symbol.

6.3. Effect of Initial Moisture Content

The effect of initial moisture content on the moisture curve is shown in Figure 13. To compare the results, moisture is represented by the moisture ratio . The drying times are 4.5, 5.5, and 6 h, respectively, for the initial = 0.6, 0.5, and 0.4. The other parameters are = 1000, = 10−13, = 1000 Pa, and = 323 K.

The temperature curve is shown in Figure 14. The temperature was obviously affected by the initial moisture, especially on the surface of the porous medium, at = 0.01 m. Higher initial moisture resulted in lower temperature during the drying process.

The moisture change is shown in Figures15 and 16 for initial moisture = 0.5 and 0.4, respectively, as compared with = 0.6 (in Figure 7). Except for the value, the changes are similar.

6.4. Effect of Heat Temperature

The effect of heat temperature on moisture is shown in Figure17. The drying end times are 6, 5, and 4.5 h, respectively, for the heat temperature = 323, 333, and 343 K. The other parameters are = 0.6, = 10−13, = 1000, and = 1000 Pa. The temperature changes at different points were similar except for the value, as shown in Figure18. A similar result was found for the moisture changes at different times in Figures 19 and 20 compared with Figure 7.

6.5. Effect of Intrinsic Permeability

Intrinsic permeability of porous medium is an inherent property and cannot be changed, and measuring it is difficult. Figures 21 and 22 show the moisture and temperature curves at different permeabilities = 10−13, 10−14, and 10−15. Intrinsic permeability has a greater effect because the transfer of free water and vapor is affected by (3.13). The drying time became evidently longer as the intrinsic permeability was reduced because the moisture movement velocity was lowered at the same pressure gradient. The temperature is increased because the mass transfer rate is lower.

The effect of intrinsic permeability on moisture change at different time is shown in Figures 23 and 24. The moisture in the middle of the porous medium was not lowered quickly, even at the drying end time for = 10−14 (Figure 23) compared with = 10−13 (Figure 7).

7. Conclusion

A coupled model of porous medium vacuum drying based on the theory of heat and mass transfer was implemented in this paper. The model was implemented and solved using COMSOL. The temperature increased quickly at the start of drying and then lowered gradually. As the drying process continued, the temperature increased slowly. In the absence of free water, temperature increased rapidly. As the drying process concluded, the temperature remained unchanged. The water evaporation rate could not be obtained during the porous medium vacuum drying process. The rate constant parameter is essential to the nonequilibrium method. When ≥ 1000, the simulation of the drying process was not evidently affected. Vapor pressure and heat transfer affected the transfer of mass. A similar effect was found in the initial moisture and the heat temperature. Intrinsic permeability had a greater effect on the drying process.


Diagonal tensor
:Diffusivity ()
:Diffusion tensor ()
:Gravity vector (m)
:Intrinsic averaged enthalpy (J)
:Water phase rate (kg )
:Intrinsic permeability ()
:Relative permeability
:Mass (kg)
:Outer unit normal to the product
:Pressure (Pa)
:Capillary pressure (Pa)
:Universal Gas constant (J )
:Time (s)
:Temperature (K)
:Moisture content (in dry basis).
Greek Letters
:Latent of phase change (J)
:Effective thermal conductivitytensor (W)
:Viscosity (kg)
:Density (kg)
:Vapor fraction.
:Dry air
:Vapor saturation.
Mathematical Operators
:Gradient operator
:Divergence operator.


This research was supported by the National Natural Science Foundation of China (Grants no. 31000665, no. 51176027).