## Services on Demand

## Journal

## Article

## Indicators

- Cited by SciELO
- Access statistics

## Related links

- Cited by Google
- Similars in SciELO
- Similars in Google

## Share

## Brazilian Journal of Chemical Engineering

##
*Print version* ISSN 0104-6632

### Braz. J. Chem. Eng. vol.30 no.4 São Paulo Oct./Dec. 2013

#### http://dx.doi.org/10.1590/S0104-66322013000400013

**FLUID DYNAMICS; HEAT AND MASS TRANSFER; AND OTHER**

**Heat transfer analyses using computational fluid dynamics in the air blast freezing of guava pulp in large containers**

**W. M. Okita; M. J. Reno; A. P. Peres; J. V. Resende ^{*}**

Federal University of Lavras, Department of Food Science, Phone: + (55) (35) 3829-1659, Fax: + (55) (35) 3829-1401, Laboratory of Food Refrigeration, P.O. Box 3037, Zip Code: 37200-000, Lavras - MG, Brazil. E-mail: jvresende@dca.ufla.br

**ABSTRACT**

Heat transfer during the freezing of guava pulp conditioned in large containers such as in stacked boxes (34 L) and buckets (20 L) and unstacked drums (200 L) is discussed. The air velocities across the cross-section of the tunnel were measured, and the values in the outlet of the evaporator were used as the initial conditions in computational fluid dynamics (CFD) simulations. The model tested was turbulent standard k-ε. The CFD-generated convective heat transfer coefficients were mapped on the surfaces for each configuration and used in procedures for the calculation of freezing-time estimates. These estimates were compared with the experimental results for validation. The results showed that CFD determined representative coefficients and produced good correlations between the predicted and experimental values when applied to the freezing-time estimates for the box and drum configurations. The errors depended on the configuration and the adopted mesh (3-D grid) construction.

**Keywords:** CFD; Airflow; Heat transfer; Fruit pulp; Freezing time.

**INTRODUCTION**

In Brazil, fruit pulp falls into the important category of frozen foodstuffs that are used as raw materials for ice cream, yogurt and jelly production. The freezing method used in the processing of guava is typically implemented in small agro industries. In these installations, the power consumption related to the freezing time is an important factor because it increases the overall production costs. Fruit pulps are packaged in a variety of package types and sizes, and a large increase in the production and consumption of fruit pulps has been recorded over the last decade. In industrial situations, the process of freezing fruit pulp inside boxes, buckets and drums that contain large amounts of product can be very complex. The homogeneity of the product is not guaranteed; voids or air bubbles can occur inside the product and packages. The product can have a standard form and size, but this rarely occurs. The external heat transfer may be constant across the entire surface of the product; however, wide variations in the heat transfer coefficients, mainly the convective heat transfer coefficient, may occur in different positions (Pham and Willix, 1987; Reno *et al*., 2011; Resende *et al*., 2013).

The rate of heat transfer depends on airflow conditions, and local variations in the heat transfer coefficients are expected along the surface; these variations produce local differences in temperature (Pham *et al*., 2009). In food freezing, the air velocity distribution determines the efficiency and the homogeneity of the treatments to which the product is being submitted. In equipment used in food processing, the airflow is generally turbulent and transient. Free-stream turbulence influences the transfer phenomena in porous media. Heat transfer heterogeneity is related to the local velocity and turbulence intensity, which influences the heat transfer coefficients in food processing (Verboven *et al*., 1997; Verboven *et al*., 2006). The measurement of local heat transfer coefficients for a complex geometry is difficult, and values associated with the different geometric shapes, surfaces and package space arrays have not been determined. For the freezing of foodstuffs in boxes, buckets and drums, the literature values show that the surface heat transfer coefficients vary when the measurements are made in different locations in the stacking. Thus, the values of the coefficients are different between the layers of product in the stack; therefore, studies that ignore these variations should be used only with due care (Pham and Willix, 1987; Mannaperuma *et al*., 1994a, 1994b; Verboven *et al*., 2006; Kondjoyan, 2006).

Proper methods of cooling may ensure the temperature uniformity of products stacked in a cold store and improve the heat transfer efficiency between the products and air (Tanaka *et al*., 2012). To improve the cooling performance of cold stores, number of studies modeling the airflow pattern and temperature distribution have demonstrated the effectiveness of various methods (Nahor *et al*., 2005; Chourasia and Goswami, 2007; Konishi *et al*., 2009). Tanaka *et al*. (2012) developed a transient three-dimensional computational fluid dynamics model to investigate the cooling performance of a partially loaded cold store in the cooling process. The model accounted for turbulence by means of the standard k-ε model with standard wall profiles (Ghasemi and Hosseini, 2012). The model was validated successfully in cooling experiments with a mean error of 0.36 m/s and 1.4 °C for the air velocity in a refrigerated room and the temperature of the product, respectively. Based on the model, the effects of different loading patterns of corrugated fiberboard containers filled with products in a cold store on cooling effectiveness were studied. As a result, flat loading with air gaps was the optimal configuration to achieve high uniformity of temperature and rapid cooling. Tapsoba *et al*., (2006) investigated experimentally and numerically the airflow patterns within a refrigerated truck loaded with slotted but empty pallets using CFD. Air velocity measurements were carried out on a reduced scale model with a Laser Doppler Anemometer above and inside the pallets. The numerical predictions obtained with a Computational Fluid Dynamics package using the RSM turbulence model showed a satisfactory agreement with experimental data in high velocity zones.

Computational fluid dynamics (CFD) can be used to simulate the local surface heat transfer coefficients on the surfaces. On the packages' surfaces, the heat transfer coefficients vary from one cross-section to another and within each cross-section (Pham *et al*., 2009). The use of simulated CFD data allows the packages to be divided into subsections and the heat transfer to be averaged over each subsection. The methods of predicting freezing time could be used to validate the average heat transfer coefficients. According to Cleland and Özilgen (1998), these methods have improved significantly over the last few decades, and they are already not the limiting factor in the accurately design of freezing systems.

The modeling of rapid cooling equipment is an obvious application of CFD; determining optimized airflow patterns can improve equipment effectiveness and energy efficiency (Smale *et al*., 2006). The objectives of this study were as follows: (i) to investigate, both experimentally and numerically, the heat transfer process during the freezing of guava pulp packed in a stack of large containers (boxes, buckets and drums) in small agro industries; (ii) to simulate the airflow pattern and air velocity using CFD and experimental measurements; (iii) to simulate the local convective heat transfer coefficients using CFD; and (iv) to validate the results of the convective heat transfer coefficients used in the freezing-time estimation model and to compare these results with the experimental data.

**MATERIAL AND METHODS**

**CFD Simulation**

**Mathematical Modeling**

The conservation equations can be written in a generalized form. The time-averaged mathematical models, together with the Reynolds decomposition model, can be written as follows:

where S_{M} represents possible momentum sources and the stress tensor, τ, is related to the shear rate by:

The conservation equation for energy that takes into account density changes that result from temperature variations can be written in a generalized form:

in which S_{E} represents possible energy sources and H_{tot} is the total enthalpy, related to the static enthalpy H (T,p) by:

The term ∇.(u.t) represents the work due to viscous stresses and is called the viscous work term. The term u.S_{M} represents the work due to external *momentum* sources and is currently neglected.

**Turbulence Model**

Numerical methods have been developed to solve Equations (1), (2) and (3), in which the complexity increases significantly because of the Reynolds stress terms on the right-hand side of Equation (2). The turbulence model used in this work is the standard k-ε model, which is frequently used in industrial applications because of its numerical robustness and low computing times. The Reynolds stress is obtained from:

Here, µ_{T} is the turbulent viscosity, δ is the identity tensor and the superscript T is the transposition operation.

The turbulent kinetic energy (k) and the dissipation rate of the turbulent kinetic energy (ε) are obtained with:

where λ is the production term by the shear effect,

and the viscosity is calculated from:

C_{1} (1.44), C_{2} (1.92), C_{µ}(0.09), σ_{k} (1.00) and σ_{ε} (given by Equation (10)) are the model constants.

**Numerical Methods**

The equations were solved numerically by the finite volume method using the commercial CFD code CFX (v. 11), in which the control volume method is used to discretize the transport equations. The pressure-velocity coupling algorithm SIMPLEC (SIMPLE Consistent) and the higher-order upwind interpolation scheme were used in all numerical experiments. More details on these schemes can be found in Patankar (1980).

The numerical experiments were performed with an accuracy of 10^{-5} for the Euclidean norm of the source mass in the pressure-velocity coupling. A transient run was performed using the steady-state results as the initial conditions with five iterations for each time step. Time steps of 0.001 seconds and a total simulation time of 40 seconds were used for the configurations of boxes, buckets and drums.

**Computational Grid and Boundary Conditions**

The Courant number is of fundamental importance for transient flows. For a one-dimensional grid, it is defined by:

in which u is the local air velocity, Δt is the time step and Δx is the mesh size. The calculated Courant number is a multidimensional generalization of this expression where the velocity and length scale are based on the mass flow into the control volume and the dimension of the control volume. As an implicit code, ANSYS CFX does not require the Courant number to be small for stability. However, for transient runs, as used in this work, the maximum and RSM Courant number are written to the output file at every time step. After the grid tests were performed, meshes were constructed with adequate refining and the time step was chosen based on a Courant number below 3 in order to ensure mesh independency and optimization of the computational time.

The freezing tunnel numerical grids were constructed using the mesh building code ICEM CFD. The final grids each had approximately 122,086 cells for the box configuration, 607,511 cells for the bucket configuration and 158,964 cells for the drum configuration. The test cases were simulated using hexahedral grids for boxes and drums and using non-structured tetrahedral grids for bucket configurations. Fig. 1 shows the 3-D surface meshes used for the configurations in a freezing tunnel for boxes, buckets and drums.

The model tested was turbulent standard k-ε. The airflow was assumed to be at 253 K (-20 °C), and the boundary conditions were an inlet air velocity of 3.5 ms-1, an outlet relative pressure of 0 atm and adiabatic conditions at the walls for the simulation of the heat transfer process.

**Experiments**

**Guava Pulp**

The guava pulp used in this study was acquired from a small agribusiness located in the municipal district of Lavras/MG, Brazil. To preserve the microbiological quality of the pulp from the time of experimental assembly to the moment of processing, 0.5 g of potassium sorbate was added per 100 g of pulp.

**Freezing of the Guava Pulp**

The pulp was batch-frozen in a freezing tunnel with forced air directed onto the product. Three types of containers were used for the freezing of the guava pulp: high-density polyethylene (HDPE) plastic boxes, HDPE buckets and metal drums. In each experimental setup, 600 kg of pulp was distributed according to the following configurations and then frozen: 1) 40 boxes, each containing 15 kg of guava pulp; 2) 40 buckets, each containing 15 kg of guava pulp; and, 3) three drums, each containing 200 kg of guava pulp. In the drum and box configurations, the pulp was wrapped in polyethylene bags prior to being introduced into the container. The containers were then stacked in the center of the freezing tunnel with equal lateral, front and bottom spacing for all configurations. The dimensions of the HDPE boxes were 0.60 m (length), 0.40 m (width) and 0.17 m (height). Spaces were left at each edge to allow the passage of air around the box. The dimensions of the polyethylene bags that contained the pulp were 0.47 m (length), 0.36 m (width) and 0.086 m (height).

The dimensions of the HDPE buckets were 0.32 m (diameter) and 0.37 m (height). The dimensions of the metallic drums were 0.82 m (height) and 0.57 m (diameter). Fig. 2 shows the configurations inside the tunnel, the dimensions, and the system of rectangular coordinates (X, Y, Z) that was used as the reference coordinate system throughout this work.

**Monitoring of Sample and Air Temperatures**

To obtain the temperature history, type-T (copper-constantan AWG-24) thermocouples and resistive thermal devices (RTDs) were installed throughout the experimental configuration at specific locations. The locations were chosen among the samples to register the transient behavior of the heat transfer during the process in the most comprehensive manner possible (Reno *et al*., 2011). The thermocouples and RTDs were connected to a signal conditioning system

(National Instruments - Model SCXI - Hungary). The temperature measurements were obtained at 2-minute intervals using the LabVIEW 8.5 software package. The average temperature of the cooling air was monitored using a parallel arrangement of type-T thermocouples distributed evenly in front of the evaporator. The registration of the temperatures was performed using a data acquisition system.

**Monitoring of Sample and Air Temperatures**

To obtain the temperature history, type-T (copper-constantan AWG-24) thermocouples and resistive thermal devices (RTDs) were installed throughout the experimental configuration at specific locations. The locations were chosen among the samples to register the transient behavior of the heat transfer during the process in the most comprehensive manner possible (Reno *et al*., 2011). The thermocouples and RTDs were connected to a signal conditioning system (National Instruments - Model SCXI - Hungary). The temperature measurements were obtained at 2-minute intervals using the LabVIEW 8.5 software package. The average temperature of the cooling air was monitored using a parallel arrangement of type-T thermocouples distributed evenly in front of the evaporator. The registration of the temperatures was performed using a data acquisition system.

**Air Velocity Measurement**

The air velocity was determined at various loca-tions along the freezing tunnel using a hot-wire anemometer (TSI - model 9535 - USA). These devices were established at different heights (coordinate Y(m) = 0.420; 0.845; 1.270; 1.695; and 2.120), lengths (coordinate X(m) = 0.000; 0.090; 0.180; 0.270; 0.360; 0.450; 0.540; 0.630; 0.720; 0.810; 0.900; 0.990; 1.080; 1.170; 1.260; 1.350; 1.440; 1.530; 1.610; 1.700; 1.790; 1.880; 1.970; 2.060; 2.150; 2.240; and 2.330) and widths (coordinate Z(m) = 0.135 and 1.195) in the freezing tunnel along the horizontal and vertical directions of the air circulation.

The vertical (u_{y}) and horizontal (u_{z}) velocities were measured and the average air velocity was determined using the equation:

Validation of Heat Transfer Coefficients from CFD Through a Freezing Time Prediction

The values for the convective heat transfer coefficients used in the freezing-time calculations were those obtained from the CFD simulation for the different configurations. For the freezing-time prediction of the samples conditioned in boxes, buckets and drums, a geometric factor (E) was used (Equation (13)):

The geometric factor (E) can be defined as the ratio between the freezing time of one test body and the real freezing time for the product with similar dimensions, with both being exposed under the same conditions.

One of the most recommended equations (Cleland, 1990; Delgado and Sun, 2001; Resende *et al*., 2007) for predicting the freezing time of one-dimensional foodstuffs (t_{f,slab}) is the method of Pham (1986).

**Boxes**

For the boxes, the geometric factor (E) was calculated using two methods: i) method EHTD (equivalent heat transfer dimensions), considered by Cleland and Earle (1982) and given by the equation:

The functions W1 and W2 depend on the Biot number (Bi_{s}) and on the dimensional ratios β_{1} and β_{2} in accordance with the equations:

and, ii) the geometric factor (E_{AN}) derived from the analytical solution of heat transfer problems that involve a phase change (Hossain *et al*., 1992):

where the values of Z_{n} and Z_{m} are the roots of the equations:

The values of Z_{nm} in Equation (20) are given by:

The six first roots of Equations (18) and (19) are tabulated in Carslaw and Jaeger (1959) as a function of the values of Bi_{s} and Bi_{s}×β_{1}.

**Buckets and Drums**

The E_{AS} form factor was calculated using Equation (21) for cylinders with heights greater than their diameters (Hossain, 1995):

where β_{1} = (height /diameter). The expressions in Equation (21) are valid for low Biot values (Bi<3) and have the advantage of simplicity (Hossain, 1995).

**Thermophysical Properties**

The thermophysical properties for the guava pulps required for calculation of the freezing times were obtained from Reno *et al*. (2011) and are shown in Table 1.

Comparison Between the Experimental and Calculated Results

A comparison of the difference between the predicted and experimental values was made according to the following equation:

**RESULTS AND DISCUSSION**

**Airflow Pathways and Air Velocities Between the Containers**

Fig. 3 illustrates the comparison of the actual measured air velocity profile with the airflow pathways simulated using CFD (standard k-ε) in the freezing tunnel for the configurations of guava pulp conditioned inside the stacked boxes, the stacked buckets and the metallic drums.

Figs. 3(a2), (b2) and (c2) show the airflow pathways directed by the baffle plates installed in the tunnel ceiling on the inlet and outlet of the evaporator. These figures also reveal that the magnitude of the air velocities and the pathways obtained from the simulations are compatible with those from the experimental measurements for the boxes, and drum configurations. Figs. 3(a1)-(a2) and 3(c1)-(c2) also show that the higher air velocity is dependent on the Y coordinate, which can be explained by the evaporator that was located in the superior position of the tunnel wall. Figures 3(a1)-(a2) and 3(c1)-(c2) further indicate that higher velocities occurred at positions with higher Y values (above the stacking and drums), where a larger free area was available for the return of air to the evaporator. The obstacles formed by the containers induced preferential airflow pathways above the stacked boxes and unstacked drums in the air return to the evaporator. The simulated profile for the buckets showed uniform air distribution along the Y coordinate that was different from the experimental measurements. In this case, lower values for the air velocity were found close to the floor and higher values at the top of the stacking.

The graph in Fig. 4 shows a comparison between the simulated and experimental results for the air velocity in the coordinates related to the tunnel center (Y = 1.27 m and Z = 0.135 m) in the box, bucket and drum configurations. Figs. 4 (a) and 4 (c) show that the simulated values scatter around the experimental results, which indicates that the simulated average values approach the experimental average values across the cross-section of the freezing tunnel. Fig. 4 (b) indicates that, in the bucket configurations, the air velocity values simulated using CFD are underestimated compared to the experimental data. In these configurations, the scattered simulated values of velocity levels are lower than the values determined by experimental measurements.

In these figures the high values of the standard deviation indicate turbulence and the peaks and oscillations observed are characteristics of streams that follow preferential channels in the gaps between the containers. These peaks were not detected in the experimental measures that were made in this coordinate. The oscillation frequency is related to the number of gaps. The frequency is lower and the length is higher in the configuration of drums (Fig. 4 (c)) where the size of the containers is larger and the number of gaps is smaller. Figure 5 shows the comparison between experimental and simulated results.

The average errors along the Coordinate (X) were 18.22%, -20.70% and -13.26% for the box, bucket and drum configurations, respectively.

**Convective Heat Transfer Coefficients Simulated Using CFD**

Figure 6 shows the maps of convective heat transfer coefficients on the containers' surfaces. These coefficients were simulated using CFD for the configurations of guava pulp processed in boxes (Fig. 6 (a)), buckets (Fig. 6 (b)) and drums (Fig. 6 (c)) as a function of the location in the stacking. Higher air velocities result in higher h values and, consequently, a higher heat transfer rate between the product and air cooling. Figs. 6 (a)-(b) show coefficient values distributed on the box and bucket surfaces which have a uniform air velocity profile along the Y coordinate (Fig. 3 (a)-(b)). In the drum configuration the highest values are placed at the top of the containers (Fig. 6 (c)). The validation is possible through an analysis and comparison of the freezing-time estimates calculated using the average values of the coefficients shown in Fig. 6 with the experimental results during cooling obtained from the readings of the temperature sensors installed in all the configurations.

**Validation of the Convective Heat Transfer Coefficients Through Freezing-Time Estimates**

The predicted freezing time is an important tool for industrial installations because it allows the residence time of processes to be calculated without extensive experimental efforts that demand extra time. One of the most significant challenges in the estimation of freezing times continues to be knowledge of the convective heat transfer coefficient, which plays a decisive role in processes that involve convection. In this work, the simulated average values of h using CFD were used for freezing-times estimates. For the box configuration, each box was considered to be a compact rectangular brick. The simulated average values of h were obtained from 2000 points extracted on the surface (Fig. 7) for each layer of the box configuration. For the bucket configuration, 5000 points in each layer were evaluated, and 1000 points on the surface were evaluated for metallic drums.

Table 2 shows the simulated average values of h obtained from CFD, the freezing-time estimates using these coefficients, and the experimental freezing times measured for the box, bucket and drum configurations. Table 2 also shows a comparison between the predicted and simulated values as function of the location in the configuration. The thermal center was defined as the last point in the configuration to undergo temperature changes or the last point to freeze. In all the experiments, the final temperature of the thermal center was -18 °C and experimental freezing times were taken when all thermocouples installed in the configuration reported temperatures immediately lower than -18 °C. This temperature was chosen because it is the most commonly used commercially and because is a parameter that represents minimal enzymatic and microbiological activities in food.

Table 2 shows that, for the box configuration, where the layer was considered to be a compact rectangular brick, the freezing times predicted by models using the average h values simulated and extracted from the whole surface are underestimated. The estimate results for the freezing times that were calculated using simulated data were shown to be accurate according to the results of the experimental measurements. The mean errors were less than 15.8% when the geometric factor based on the EHTD method (Cleland and Earle, 1982) was tested. These errors were reduced when the semi-analytical method (Hossain *et al*., 1992) was used to calculate the geometric factor. With this procedure, the error was less than 10.5%. The reduced error values can be explained. For the box configuration, the dimensional ratios β_{1} and β_{2} are greater than 3 (β_{1} = 4.19 and β_{2} = 5.47). According to Hossain (1995), when the dimension ratios are high or when the Biot number is high (Bi>3), various methods became inexact; in these cases, the full analytical method should be used to provide the best estimates of freezing times. Another observation that can explain the underestimated values is that, in this configuration, the empty spaces between the packages were not taken into account. Therefore, the heat transfer surface areas applied to the calculations were reduced and did not accurately describe the real behavior of the air passing between the product configurations.

The errors for the bucket configuration between the experimental and predicted freezing times obtained using Equations (13) and (21) and the simulated h values from CFD were less than 26.3% (Table 2) and were also underestimated. The freezing-time estimates calculated with Equations (13) and (21) for buckets were deemed insufficient for commercial acceptance (greater than 20%). These simulations were not made with individual buckets, and the runs were performed with consideration of the layer of buckets that included the front and rear lines in the stacking. In this case, local flow characteristics were not taken into account. Table 2 shows that the errors increase as a function of the stacking height. In these coordinates, the air velocity is higher and, with the explicit presence of the airflow, preferential pathways were not considered in the evaluation of local h values.

Equations (13) and (21) were also tested for the freezing-time estimates of metallic drums that contained guava pulp, and the average h values were obtained considering all details of the airflow around the cylindrical geometry. For the metallic drums, the use of the simulated average h values yielded accurate estimates with mean errors of 2.3% for the drum located closer of the door, 6.7% for the central drum, and 11.3% for the drum located at a greater distance from the door. In these simulations, the drums were considered individually, and the conditions were established to express the airflow pathways.

**CONCLUSION**

In this work, the predicted values provided satisfactory results to represent the airflow pathways and heat transfer coefficients that occur in the air blast freezing of guava fruit pulps in large containers. The airflow pathways and the magnitude of the air velocities from simulations were similar to those experimentally measured for the box, bucket and drum configurations. The air velocities were highest as a function of the height increase of the tunnel. The accuracy was due to the evaporator location in the superior position of the tunnel wall. The containers with the guava pulp form obstacles to the air flow and induce preferential pathways above the stacked product in the air return to the evaporator.

The surface heat transfer coefficient values were mapped, and the results showed that the highest values corresponded to those locations where the air velocities were greatest. The validations of these results obtained from CFD were made to compare the estimated values for freezing times with experimental measurements. The results indicated good correlations between the predicted and experimental values when applied for the freezing-time estimates for the box and drum configurations. The mean errors were less than 15.8% using the EHTD method and were less than 10.5% using a semi-analytical method applied to compact rectangular bricks in the box configuration. For the bucket configuration, the mean errors were obtained and were compared to the predicted values using the simulated average h values; the errors of the experimental values were approximately 26.5% and were underestimated. In this case, the poor results were attributed to the use of a 3-D grid that was built without consideration of each bucket as an individual object. For the drums, the use of simulated h values obtained from CFD produced mean errors of 3.3% for the drum located near the door, 6.7% for the central position, and 11.3% for the drum positioned in the bottom of the tunnel.

**ACKNOWLEDGMENTS**

The authors wish to thank the Fundação de Amparo à Pesquisa do Estado de Minas Gerais (FAPEMIG- Brazil), Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq - Brazil) and Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES - Brazil) for financial support for this research.

**REFERENCES**

Carslaw, H. S. and Jaeger, J. C., Conduction of Heat in Solids. 2th Ed., Oxford University Press, London (1959). [ Links ]

Chourasia, M. K. and Goswami, T. K., CFD simulation of effects of operating parameters and product on heat transfer and mass loss in the stack of bagged potatoes. Journal of Food Engineering, 80, 947-960 (2007). [ Links ]

Cleland, A. C., Food Refrigeration Process**:** Analysis, Design and Simulation. Elsevier Applied Science, London/New York (1990). [ Links ]

Cleland, A. C. and Earle, R. L., Freezing time prediction for foods - a simplified procedure. International Journal of Refrigeration, 5, 134-140 (1982). [ Links ]

Cleland, A. C. and Özilgen, S., Thermal design calculations for food freezing equipment- past, present and future. International Journal of Refrigeration, 21, 359-371 (1998). [ Links ]

Delgado, A. E. and Sun, D-W., Heat and mass transfer models for predicting freezing process - a review. Journal of Food Engineering, 47, 157-174 (2001). [ Links ]

Ghasemi, H. and Hosseini, S. H., Investigation of hydrodynamics and transition regime in an internal loop airlift reactor using CFD. Brazilian Journal of Chemical Engineering, 29, 821-833 (2012). [ Links ]

Hossain, M. M., A simple method of freezing time calculation for foodstuffs of various shapes. Food Australia, 47, 109-112 (1995). [ Links ]

Hossain, M. M., Cleland, D. J. and Cleland, A. C., Prediction of freezing and thawing times for foods of two-dimensional irregular shape by using a semi-analytical geometric factor. International Journal of Refrigeration, 15, 235-240 (1992). [ Links ]

Kondjoyan, A., A review on surface heat and mass transfer coefﬁcients during air chilling and storage of food products. International Journal of Refrigeration, 29, 863-875 (2006). [ Links ]

Konishi, Y., Tanaka, F., Uchino, T. and Hamanaka, D., CFD prediction of the temperature distribution within a refrigerated truck filled with fruit and vegetables during transport. Transactions of the Japan Society of Refrigerating and Air Conditioning Engineers, 26, 159-165 (2009). [ Links ]

Mannaperuma, J. D., Singh, R. P. and Reid, D. S., Effective surface heat transfer coefficients encountered in air blast freezing of whole chicken and chicken parts, individually and in packages. International Journal of Refrigeration, 17, 263-272 (1994a). [ Links ]

Mannaperuma, J. D., Singh, R. P. and Reid, D. S., Effective surface heat transfer coefficients encountered in air blast freezing of single plastic wrapped whole turkey. International Journal of Refrigeration, 17, 273-280 (1994b). [ Links ]

Nahor, H. B., Hoang, M. L., Verboven, P., Baelmans, M. and Nicolaï, B. M., CFD model of airflow, heat and mass transfer in cooling stores. International Journal of Refrigeration, 28, 368-380 (2005).

Patankar, S., Numerical Heat Transfer and Fluid Flow. Ed. Hemisphere Publishing Corporation, New York (1980). [ Links ]

Pham, Q. T., Simplified equation for predicting the freezing time of foodstuffs. Journal of Food Technology, 21, 209-219 (1986). [ Links ]

Pham, Q. T. and Willix, J., Heat transfer coefficients in the air blast freezing of rows of cartons. In: Proceedings 17th Int. Congress of Refrigeration, 350-357 (1987). [ Links ]

Pham, Q. T., Trujillo, F. J. and McPhail, N., Finite element model for beef chilling using CFD-generated heat transfer coefficients. International Journal of Refrigeration, 32, 102-113 (2009). [ Links ]

Smale, N. J., Moureh, J. and Cortella, G., A review of numerical models of airflow in refrigerated food applications. International Journal of Refrigeration, 29, 911-930 (2006). [ Links ]

Reno, J., Resende, J. V., Peres, A. P., Giarola, T. M. O. and Prado, M. E. T., Heat transfer and energy consuption in the freezing of guava pulp in large containers. Applied Thermal Engineering, 31, 545-555 (2011). [ Links ]

Resende, E. D., Kieckbusch, T. G., Toledo, E. C. V. and Maciel, M. R. W., Discretisation of the non-linear heat transfer equation for food freezing processes using orthogonal collocation on finite elements. Brazilian Journal of Chemical Engineering, 24, 399-409 (2007). [ Links ]

Resende, J. V., Prado, M. E. T. and Silveira Junior, V., Non-uniform heat transfer during air-blast freezing of a fruit pulp model in multilayer boxes. Food and Bioprocess Technology, 6, 146-159 (2013). [ Links ]

Tanaka, F., Konish, Y., Kuroki, Y., Hamanaka, D. and Uchino, T., The use of CFD to improve the performance of a partially loaded cold store. Journal of Food Process Engineering, 35, 874-880 (2012). [ Links ]

Tapsoba, F., Moureh, J. and Flick, D., Airflow patterns in an enclosure loaded with slotted pallets. Int. J. Refrig., 29, 899-910 (2006). [ Links ]

Verboven, P., Nicolai, B., Scheerlinck, N. and Baerdemaeker, J., The local surface heat transfer coefficient in thermal food process calculations: A CFD approach. Journal of Food Engineering, 33, 15-35 (1997). [ Links ]

Verboven, P., Flick, D., Nicolaï, B. M. and Alvarez, G., Modelling transport phenomena in refrigerated food bulks, packaged and stacks: Basics and advances. International Journal of Refrigeration, 29, 985-997 (2006). [ Links ]

*(Submitted: August 7, 2012 ; Revised: January 11, 2013 ; Accepted: January 18, 2013)*

* To whom correspondence should be addressed