- Email: [email protected]
Predicting heat conduction during solidification of a food inside a freezer due to natural convection
Journal of Food Engineering 56 (2002) 17–26 www.elsevier.com/locate/jfoodeng
Predicting heat conduction during solidification of a food inside a freezer due to natural convection Nelson O. Moraga *, Hern an G. Barraza Departamento de Ingenierıa Mec anica, Universidad de Santiago de Chile, Av. Lib. Bdo. OÕHiggins 3363, Casilla 10233, Correo 2., Santiago, Chile Received 30 June 1999; received in revised form 14 December 2000; accepted 13 March 2002
Abstract A new mathematical model for numerical simulation of two dimensional food freezing due to natural convection is presented. Fluid mechanics and heat transfer by natural convection between air and a solid food in a freezer are predicted along with the heat conduction inside a plate shaped food. The mathematical model used includes continuity, linear momentum and energy partial differential equations for air and the non-linear heat diffusion equation for the solidification of the water content in the food. Unsteady two dimensional results of the velocity and temperatures distributions in the air around the food obtained by using the finite volume method are presented. Experimental data of temperature variation in time for a slab shaped food inside a freezer are used to assess the accuracy of the mathematical model and the solution procedure. It is found that the use of this model, that does not require heat transfer coefficients, allows to simulate two dimensional freezing experiments in a slab shaped food by natural convection with maximum and mean deviations of 5.2% and 1.2%, respectively. Ó 2002 Elsevier Science Ltd. All rights reserved.
1. Introduction An accurate knowledge of freezing processes is needed in the design of refrigeration systems to assure high quality of foods at reduced cost. From the heat transfer point of view the task to be performed is to determine the transient temperature distribution in the food during the cooling process and to predict the freezing time. Such a simple task is by far not easy since there are at least four main source of uncertainties: (a) thermal properties of foods; (b) freezing heat transfer coefficients; (c) approximations made in the prediction method and (d) the mathematical model used. Freezing processes are difficult to predict because the non-linear nature of the diffusion equation, in which the food properties: density, specific heat, enthalpy and thermal conductivity vary with temperature. Experimental methods to measure and model thermal properties of frozen foods have been reviewed by Singh (1995). Estimated properties with the Box–Kanemasu method and experimental thermal conductivity and apparent volumetric specific heat of carbohydrate and protein rich foods have been found that do not agree *
Corresponding author. Tel./fax: +56-2-682-5498. E-mail address: [email protected] (N.O. Moraga).
well for methylcellulose and wheat gluten, Saad and Scott (1996). Predicted values from food composition data by using the SchwartzbergÕs bound water model for freezing point and bound water for meats, surimi, fruits, cheese and Tylose gels have been found to be within 10% of actual values, Pham (1996). Sensitivity of freezing time calculation for Alaska pollock surimi performed with a commercial finite element program have shown that: apparent specific heat, geometric dimension, heat transfer coefficient, fluid ambient temperature, thermal conductivity and density influence the calculations, De Qian and Kolbe (1994). The influence of change of volume expansion, thermal conductivity, specific heat, unfrozen water, initial freezing point on the predicted freezing curve with an analytical model has been reported by Hung and Thompson (1983) to account for differences of freezing times between 2% and 13%. The imprecise knowledge of the heat transfer coefficient h has been found to be the major error source in freezing time prediction methods, Cleland and Earle (1984). Small deviations in the heat convective coefficient were found, by the use of a finite element method, to result in large deviations in the core temperature of food chilling using air, Nicolai and Baerdemacker (1996). Spatial variations of h and the measurement techniques used in the experiments to determine h are
0260-8774/02/$ - see front matter Ó 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 2 6 0 - 8 7 7 4 ( 0 2 ) 0 0 1 3 5 - 8
18
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Nomenclature a b c Cp d fi g k n p t T u v x y a b q / U
length of the freezer chamber (m) height of the freezer chamber (m) length of the plate shaped food (m) specific heat (J/kg °C) height of the plate shaped food (m) functions to define time variations of freezer wall temperatures (°C) gravitational acceleration (m2 /s) thermal conductivity (W/m °C) normal direction to food external area pressure of air (N/m2 ) time (s) temperature (°C) air velocity in axial x-direction (m/s) air velocity y-direction (m/s) horizontal cartesian coordinate (m) vertical cartesian coordinate (m) under-relaxation factor coefficient of thermal expansion (1/°C) density (kg/m3 ) transport dependent variable water content in food on wet basis
two issues that affect the accuracy of the freezing prediction methods, Cleland, Cleland, and Jones (1994). The approximations commonly used in the predictions methods are: (a) that no changes of volume occurs in the food during freezing, even though density of ground beef changes about 4,8%; (b) weight losses in the range of 2–3% caused by water evaporation at the food surface are considered to be negligible; (c) the food is assumed to be an homogeneous isotropic material; (d) food texture is considered to be similar to the materials used to measure the heat transfer coefficient by the lumped system method. Physical, mathematical and computational aspects of freezing and melting processes that have been examined by Alexiades and Solomon (1993), indicated that the common practice used to construct a mathematical model has been based on solving the heat diffusion equation with boundary conditions of the third kind. The mathematical model used to find the time evolution of temperature distributions in the food has been based on the non-linear, second order, partial differential heat diffusion equation with heat convective boundary conditions on the food surfaces. Approximate analytical method to predict freezing times have been derived using shape factors, Salvadori, Michelis, and Mascheroni (1997), Cleland, Cleland, and Earle (1987), Cleland and Earle (1982) with percentage differences with experimental data in the range between 5% and 15%. Finite
l
dynamic viscosity (N s/m2 )
Subscripts A geometrical center of food a air F geometrical center of air in the freezer f food oa initial state in air of initial state in food p pressure ref reference state T temperature u air velocity in horizontal direction v air velocity in vertical direction / transport dependent variable Superscripts n number of iteration n 1 number of previous iteration new new calculated value guessed values for dependent variables 0 correction values for dependent variables
difference methods have been more often used to predict freezing times with errors of 10% for elliptical foods, Ilicali, Cetin, and Cetin (1996) and for rectangular meat patties, Rubiolo de Reinick (1996). A two-time-step finite difference method has been shown to have smaller maximum errors (7–19%) than the Crank–Nicolson method (16–37%) in the prediction of freezing of codfish, Saad and Scott (1997). Freezing of irregular shaped foods has been predicted by using finite element method by Comini, Del Guidice, Lewis, and Zienkiewiez (1974), Abdalla and Singh (1985), Caroll, Mohtar, and Segerlind (1996), Puri and Anantheswaran (1993), among others. This paper presents a new mathematical model that uses a numerical procedure for predicting heat conduction with solidification of the water content for solid foods, being cooled by natural convection in air inside a freezer. Instead of using a standard model of heat diffusion for the food with heat convection boundary conditions the mathematical model proposed replaces the external boundary condition for a set of four partial differential for the simulation of the natural convection in the surrounding air. The finite volume method is used to solve the system of discretized equations. Experimental measurement with thermocouples, scanner and a PC of time histories of temperature inside a standard package, used to test freezers, that has properties very close to the beef meat is used to validate the numerical solution procedure.
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
This original procedure has been used to simulate two dimensional unsteady heat conduction during freezing of a slab shaped piece of salmon meat under a forced convective heat transfer cooling, Moraga and Salinas (1999). Unsteady local heat transfer coefficients by forced convection have been successfully calculated at the four walls of a slab shaped food in forced air around it after the temperature fields in the food and in the air were determined, Moraga and Medina (2000). In these two previous papers the cooling mechanism was forced convection from the food to the surrounding air and the fluid motion of air is independent of the temperature field in the cooling fluid. In the present work the momentum equation in the vertical direction is coupled to the energy equation by the buoyancy forces that are calculated in terms of the temperature difference and the prediction of the solidification due to natural convection is more complicated. By the other hand, the heat transfer coefficients for natural convection are generally known with lower accuracy than in forced convection. Therefore, a procedure such as the one proposed that does not require the knowledge of heat transfer coefficient should be more interesting to use in food freezing processes by natural convection.
2. Physical situation and mathematical model A slab shaped food of height equal to d and length equal to c is inside a freezer of dimensions b by a initially at temperature Tof in contact with air at rest that has a uniform initial temperature Toa , as it is shown in Fig. 1. The refrigeration system is activated, and the freezer walls start to cool, along with, the air inside the rectangular cavity and later the food located on top of the bottom wall. In order to define the mathematical model the following assumptions are made. The air is assumed to be Newtonian fluid, that flows in a laminar flow. Unsteady natural convection in the air and heat conduction in the
19
food are assumed to be two dimensional. Thermal properties of the food change with temperature and a temperature dependent apparent specific heat is used to predict the phase change from liquid water to ice. Mass transfer of water inside the food and the change of volume have been neglected. Continuity, linear momentum in x and y spatial directions and energy equations describe the natural convection in air, ou ov þ ¼0 ð1Þ ox oy ou ou ou op o ou o ou l l q þu þv ¼ þ þ ot ox oy ox ox ox oy oy ð2Þ ov ov ov op o ov o ov q l l þu þv ¼ þ þ ot ox oy oy ox ox oy oy þ gbðT Tref Þ o qCp T o qCp T oT o oT þv ¼ þu k qCp ox oy ot ox ox o oT k þ oy oy
ð3Þ
ð4Þ
where according to BoussinesqÕs approximation air density has been assumed to change linearly with temperature. The values of the constant properties used for air were: q0 ¼ 1:287 kg/m3 , Cp ¼ 1:006 kJ/kg °C, k ¼ 0:024 W/m °C, l ¼ 1:71 105 kg/m s, b ¼ 3:5 103 1/K. The non-linear heat diffusion equation inside the food assuming temperature changes in the thermal properties is oT o oT o oT ðqCp Þf ¼ kf kf þ ð5Þ ot ox ox oy oy where a temperature dependent apparent specific heat Cpf has been used to capture the phase change of liquid water to ice inside the food, Hsiao (1985). Initially, air at rest inside the freezer and food are at uniform temperature: for
t ¼ 0;
u ¼ v ¼ 0;
Tair ðx; y; 0Þ ¼ Toa ;
Tfood ðx; y; 0Þ ¼ Tof ;
ð6Þ
Boundary conditions in the freezer include the no slip condition and temperature variations in time and space in all the internal walls of the freezer, u¼v¼0
Fig. 1. Physical situation of food inside freezer.
at x ¼ 0;
x ¼ a;
y ¼ 0;
y¼b
T ðx ¼ 0; y; tÞ ¼ f1 ðtÞ;
T ðx ¼ a; y; tÞ ¼ f2 ðtÞ;
T ðx; y ¼ 0; tÞ ¼ f3 ðtÞ;
T ðx; y ¼ b; tÞ ¼ f4 ðtÞ
ð7Þ ð8Þ
where the boundary condition of the fourth kind can be determined by experimental measurements of the wall temperature variations with time.
20
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Table 1 Thermal properties of food used in calculations 8 3 > < q ¼ 1053 kg=m T P Tcr : Cp ¼ 1:448ð1 UÞ þ 4:187U J=kg °C > : K ¼ 0:09705 þ 0:501U þ 5:052 104 UT W=m °C 8 > > >
1053 kg=m3 0:98221 þ 0:1131U þ 0:25746ð1 UÞ=T T < Tcr : > Cp ¼ 3:874 2:534U þ 902:893ð1 UÞ=T 2 J=kg °C > > : K ¼ 0:388 þ 1:412U þ 0:9575=T W=m °C
Conjugate boundary conditions are imposed between the air and the food. oTa oTf ðTa Þw ¼ ðTf Þw ; ka ¼ kf ð9Þ on w on w Table 1 shows the thermal properties of the material, assumed to be those of the beef meat, that were calculated in terms of the initial water content on wet basis U as reported by Salvadori and Mascheroni (1991).
3. Method of numerical solution The physical domain of air and food inside the freezer was discretized with different non-uniform grids of 27 19, 39 31, 51 36, 69 51 and 81 78 nodes. Temperature, pressure, specific heat and density were calculated at the grid points. A staggered grid was used to evaluate the velocity components and thermal conductivity on the faces of the control volumes (Harlow & Welch, 1965). Control volume faces were located between the air the food and the size of the grid was refined near the food walls. The governing Eqs. (1)–(5) were integrated in time and space according to the practice of the finite volume method. The diffusion terms of the viscous forces in the linear momentum equations and the heat conduction in energy equations were calculated by using linear interpolation functions between the dependent variables (/ ¼ u; v; T Þ and the independent spatial variables ðx; yÞ. Convective terms in linear momentum and energy equations were evaluated by the power-law differencing scheme (Patankar, 1980). The semi-implicit method for pressure-linked equations, SIMPLE algorithm was used to solve the governing equations (Patankar, 1980). The guess and correct procedure was initiated with a guessed value for pressure p , that was used to solve the discretized linear momentum equations to obtain the velocity components u and v . These values were introduced into the conti-
nuity equation to calculate the correction term for pressure p0 and for the velocity components u0 , v0 . The corrected values for the dependent variables were evaluated as (/ ¼ / þ /0 ), where /0 were the corrections for / ¼ u, v and p. An iterative procedure was used to find a convergent solution for pressure, velocity and temperature. The iteratively improved pressures, velocities and temperatures were calculated using under-relaxation: pnew ¼ p þ ap p0 ; / ¼ u; v; T
/n ¼ a/ / þ ð1 a/ Þ/n1 ; ð10Þ
where ap and all a/ ð/ ¼ u, v, T Þ were under-relaxation factors with the following values between 0 and 1, au ¼ av ¼ 0:03;
ap ¼ 0:3;
aT ¼ 0:9
ð11Þ
A line by line method, that results as a combination of the direct tri-diagonal matrix algorithm (TDMA) for one dimensional situations and the iterative Gauss– Seidel method with successive under relaxation was the numerical procedure used to solve the system of discretized governing equations (Versteeg & Malalasekera, 1996). The iterative procedure was ended when the discretized continuity equation was satisfied with an error smaller than 106 in every one of the control volumes. The numerical procedure used was based on the fact that the heat diffusion equation, Eq. (5), is an asymptotic case of the energy equation, Eq. (4). Therefore, the space occupied by the solid food was assumed to be filled with a fluid with an extremely large value of the dynamic viscosity, l ¼ 1020 , Pa s, and hence that on that region the numerical solution for the fluid mechanics problem gave as the result that the velocity components were zero. Temperature dependent values for density, apparent specific heat and thermal conductivity of the food were assigned to the nodes on that part of the computational domain and the energy equation was solved with the zero values for the two velocity components.
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
4. Experimental procedure A standard package of 1.0 kg, with dimensions equal to 0.05 m by 0.10 m in length, to test domestic freezers was used in the experiments. The package composition was 467.2 g of water, 230 g of oxietilmetil cellulose, 5 g of NaCl and 0.6 g of 6-Cl-m cresol. An extremely thin plastic membrane was used to seal the surface in order to avoid a loss of weight by mass transfer in its surface. Thermal properties of such a substance are very close to the beef meat. The package was put inside a freezer chamber of 0.49 m in length by 0.35 m in height, as shown in Fig. 1. Initially, the still air and the freezer walls were equal to 28 °C, while the package was at 26 °C. The refrigeration cycle was started and temperatures in the walls began to decrease in while the temperature data was recorded every 10 s until the solid reached a final state with a uniform temperature of 33 °C, after 18,000 s. Temperature data from 18T-type thermocouples (0.3 mm in diameter) located on the internal surface of the freezer walls as shown in Fig. 2, were recorded by a personal computer. Fig. 3 depicts the wall temperature histories
21
recorded during the first 8000 s, where the large temperature changes are seen to occur during the first 2000 s. The non-symmetric nature of the wall temperatures was caused by the direction of the refrigerant fluid flow inside the evaporator located on the walls of the freezer. The unsteady wall temperature data measured were used in the numerical simulation as the known boundary conditions on the walls of the freezer.
5. Results The program written to solve the system of equation according to the SIMPLE algorithm in finite volume was tested to produce numerical results that were independent of the grid used and of the time step. A total of five different meshes with different number of nodes were used. The physical domain that included the air in the freezer and the solid food was discretized by using five alternative non-uniform grids of 27 19, 39 31, 51 36, 69 51 and 81 78 nodes. These five grids were selected based on previous experiences solving conjugate a forced convection and food freezing problem, Moraga & Salinas (1999). Fig. 4 shows that the grid used with 81 78 nodes was more refined in the vicinity of the food walls where temperature and velocity gradients were larger. A total of three control points were chosen to test that the numerical solution was grid independent. Fig. 4 depicts the location of the two regions in the food that coincide with the position where the time history of food temperature was measured, TA and TB , and in the space occupied by the air, TF . The
Fig. 2. Position of thermocouples on freezer walls.
Fig. 3. Temperature history of freezer walls.
Fig. 4. Mesh used and thermocouples locations for air and food.
22
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
thermocouples TA and TF were located at the geometrical center of the food and of the freezer, respectively. The results of the temperature histories predicted in the geometrical center of the food and of the air in the center of the freezer for the five grids used in the numerical solution procedure are shown in Fig. 5. A time step of 1 s was used in these calculations. Temperature histories could be predicted in the food with grids of at least 51 36 nodes while 81 78 nodes were needed to calculate the temperature variations with time for the air independently of the grid size. Therefore it was concluded that a grid of 81 78 was needed to solve the discretized natural convection of air with freezing of the water content of the food located inside the freezer under time changing wall temperature boundary conditions. The effect of the time step used in the numerical results predicted for velocities and temperatures was also studied. The mathematical model defined by Eqs. (1)–(9) was discretized with the selected grid of 81 78 nodes. Four different time steps of 0.1, 1, 10 and 100 s were used to integrate the basic discretized equations in time. Fig. 6 shows the time evolution of temperatures calculated with these four different time steps at the geometrical centers in the food (TA ) and in the air (TF ). The calculations with Dt larger or equal to 10 s originates temperatures differences of 5 °C at the end of the freezing process in the food geometrical center and dif-
Fig. 5. Effect of grid size on temperatures measured at geometrical center of food (TA ) and of air in freezer (TF ), (Dt ¼ 1 s).
Fig. 6. Effect of time step on temperatures measured at geometrical center of food (TA ) and of air in freezer (TF ), (grid size: 81 78 nodes).
ferences of 4 °C in the prediction of the air temperature in the center of the freezing chamber. The results for the temperature histories calculated with different time steps were seen to be time step independent when Dt was equal to 1 s. The evolution in time of the air velocity vectors in the freezer is shown in Fig. 7 for 300, 1200, 2700, 5400, 9000 and 14,400 s. At time equal to 1200 s the center of the two main vortexes is seen to be in the lower position and two secondary non-symmetric flows are seen to appear on top of the food. Tertiary vortexes start to develop in the two lower corners for time equals to 2700 s. A more complicated fluid flow pattern with six vortexes is seen for times of 9000 and 14,400 s. The velocity at the center of center of the freezer reaches a maximum value of 0.5 m/s after 6500 s and it remains in this value for the next 6000 s until it begins to decrease. The velocity near the surface of the food changes with time and in space. The maximum values of the velocities occur at the two upper corners of the food and they are reached after 2000 s. Temperature distributions in the air and in the food are shown for six time instants, in Fig. 8. A uniform temperature equal to 26 °C is observed for almost all the lower half part of the freezer at time equal to 300 s. After 2700 s temperatures lower than 0 °C are noticed in the air that is in the half upper part of the freezer. Air temperature is seen to reach temperatures below 0 °C in all the freezer area for a time of 9000 s, while the liquid
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
23
Fig. 7. Time evolution of air velocity in freezer.
water had changed to ice in the lower half portion of the food. Finally at 14,400 s temperature in the food is seen to be below 20 °C. The non-symmetrical initial cooling process caused by the direction of the circulation of the refrigerant liquid in the external walls of the freezer last until 5400 s and from then on the isotherms tends to be almost symmetrical in y-direction. The validation procedure used include a comparison of temperature histories predicted in two locations in the food (TA , TB ) and one in the air at the center of the freezer (TF ) against measurement performed in the same locations with thermocouples. A comparison between predicted values and data measured with thermocouples for the temperature history at two locations in the food was performed and the results are shown in Fig. 9. The predicted freezing curve at two locations in the food are seen to describe the experimental results during all the freezing process. The predicted curves for TA and TB are seen to have an unusual increment around 14,000 s. At this instant of time, the temperature distribution in the food and in the surrounding air shown in Fig. 7 becomes
more uniform. The sudden changes in the temperature at TA and TB after 14,000 s seems to be caused by the sudden changes of the fluid temperature at the same time (see Fig. 6). Maximum deviation for the temperature at the geometrical center of the food TA and at the lower positions located at one fourth of the food height TB occurs near the beginning of the phase change process. The maximum and mean deviations of the predicted temperatures and the measured temperatures at the positions TA and TB inside the food were calculated as percentages of the total experimental temperature range (56 °C). The maximum and mean deviations for TA were 4.6% and 1.2%, respectively and 5.4% and 0.9% for TB . The accuracy in the numerical simulation for onedimensional food freezing of an infinite slab has been discussed by Saad & Scott (1997). They have found that finite differences with a two time step method has maximum percent errors for simulations of the temperature data measured in codfish during freezing in the range between 7.2% and 19% and mean percent
24
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Fig. 8. Unsteady temperature distributions in air and in food.
errors in the range from 1.8% to 6.4%. The temperature simulated using the finite difference method along with the Crank–Nicolson scheme showed larger discrepencies, with maximum percent errors between 16.5% and 37.4% and mean percent errors in the range from 7.8% to 13.5%. It is evident that the full two dimensional conjugate natural convection–conduction model solved by the finite volume method exhibits accurate results for temperature predictions of food freezing inside a freezer. Another advantage of the proposed method it is that the knowledge of the spatial and time variations of the heat transfer coefficient around the food surface are not needed in the calculations.
6. Conclusions A new mathematical model to simulate the unsteady two dimensional natural convection in air and heat conduction with solidification of water content in a slab shaped solid food located inside a freezer has been presented. This model eliminates the need of specifying a convective heat transfer coefficient and replaces the convective boundary condition by four partial differential equations for the non-isothermal laminar flow of air around the food. The numerical predictions of cooling curves in a slab of a piece of material with the properties of beef meat have been found to have mean and maximum percent deviations of 5.4% (3 °C) and 2.1% (1.2
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Fig. 9. Experimental and predicted cooling curves for at the food geometrical center (TA ) and at 1/4 of food height (TB ).
°C) for simulations of the experimental data measured in the food with thermocouples. The finite volume method with the SIMPLE algorithm has been implemented to solve the discretized system of equations. A sensitivity analysis of the discretization process allows to recommend the use of a non-uniform grid with 81 78 nodes and time steps of 1 s to calculate the unsteady distributions of velocity and temperature in the air inside the freezer and the temperature variations in time of a slab shaped food.
Acknowledgements This work was supported by CONICYT (Chile) through FONDECYT 1000207 project and by Universidad de Santiago de Chile, project DICYT 019316 MB.
References Abdalla, H., & Singh, R. (1985). Simulation of thawing of foods using finite element method. Journal of Foods Process Engineering, 7, 273–286.
25
Alexiades, V., & Solomon, A. D. (1993). Mathematical modeling of melting and freezing processes. Washington: Hemisphere Publishing. Caroll, N., Mohtar, R., & Segerlind, L. J. (1996). Predicting the cooling time for irregular shaped food products. Journal of Food Process Engineering, 19, 385–402. Cleland, A. C., & Earle, R. L. (1982). Freezing time predictions for foods––a simplified procedure. International Journal of Refrigeration, 5, 134. Cleland, A. C., & Earle, R. L. (1984). Assessment of freezing time prediction methods. Journal of Food Science, 49, 1034– 1042. Cleland, D. J., Cleland, A. C., & Earle, R. L. (1987). Prediction of freezing and thawing times for multidimensional shapes by simple formulae, Part 2: irregular shapes. International Journal of Refrigeration, 10, 234–240. Cleland, D. J., Cleland, A. C., & Jones, R. S. (1994). Collection of accurate experimental data for testing the performance of simple methods for food time prediction. Journal of Food Process Engineering, 17, 93–119. Comini, G., Del Guidice, S., Lewis, R. W., & Zienkiewiez, O. C. (1974). Finite element solution of non linear heat conduction problems with special reference to phase change. International Journal Numerical Methods in Engineering, 8, 613–624. De Qian, W., & Kolbe, E. (1994). Analysis of food block freezing using a P-C based finite element package. Journal of Food Engineering, 21, 521–530. Harlow, F. H., & Welch, J. E. (1965). Numerical calculation of timedependent viscous incompressible flow of fluid with free surface. Physics of Fluids, 8, 2182–2189. Hsiao, J. S. (1985). An efficient algorithm for finite difference analysis of heat transfer with melting and solidification. Numerical Heat Transfer, 8, 653–666. Hung, Y. C., & Thompson, D. R. (1983). Freezing time prediction for slab foodstuffs by an improved analytical method. Journal of Food Science, 48, 555–560. Ilicali, C., Cetin, M., & Cetin, S. (1996). Methods for the freezing time of ellipses. Journal of Food Engineering, 28, 361–372. Moraga, N. O., & Medina, E. E. (2000). Conjugate forced convection and heat conduction with freezing of water content in a plate shaped food. International Journal of Heat and Mass Transfer, 43, 495–517. Moraga, N. O., & Salinas, C. H. (1999). Numerical model for heat and fluid flow in food freezing. Numerical Heat Transfer, part A, 35, 495–517. Nicolai, B. M., & Baerdemacker, J. de. (1996). Sensitivity analysis with respect to the surface heat transfer coefficient as applied to thermal process calculations. Journal of Food Engineering, 28, 23– 33. Patankar, S. V. (1980). Numerical Heat Transfer and Fluid Flow. Washington: Hemisphere. Pham, Q. T. (1996). Prediction of calorimetric properties and freezing time of foods from composition data. Journal of Food Engineering, 30, 95–107. Puri, V. M., & Anantheswaran, R. C. (1993). The finite element method in food processing: a review. Journal of Food Engineering, 19, 242–274. Rubiolo de Reinick, A. C. (1996). Average and center temperature vs time evaluation for freezing and thawing rectangular foods. Journal of Food Engineering, 30, 299–311. Saad, Z., & Scott, E. P. (1996). Estimation of temperature dependent thermal properties of basic food solutions during freezing. Journal of Food Engineering, 28, 1–19. Saad, Z., & Scott, E. P. (1997). Analysis of accuracy in the numerical simulation of the freezing process in food materials. Journal of Food Engineering, 31, 95–111. Salvadori, V. O., & Mascheroni, R. H. (1991). Prediction of freezing and thawing times of foods by means of a
26
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
simplified analytical method. Journal of Food Engineering, 13, 67– 78. Salvadori, V. O., Michelis, A. de, & Mascheroni, R. H. (1997). Prediction of freezing times for regular multidimensional foods using simple formulae. Lebensmittel–Wissenchaft und–Technologie, 30, 30–35.
Singh, R. P. (1995). Thermal properties of frozen foods. In M. A. Rao & S. S. H. Rizvi (Eds.), Engineering Properties of Foods (2nd ed.). New York: Marcel Dekker (chapter 4). Versteeg, H. K., & Malalasekera, W. (1996). An introduction to computational fluid dynamics. The finite volume method. Essex: Longman Scientific and Technical.
Predicting heat conduction during solidification of a food inside a freezer due to natural convection Nelson O. Moraga *, Hern an G. Barraza Departamento de Ingenierıa Mec anica, Universidad de Santiago de Chile, Av. Lib. Bdo. OÕHiggins 3363, Casilla 10233, Correo 2., Santiago, Chile Received 30 June 1999; received in revised form 14 December 2000; accepted 13 March 2002
Abstract A new mathematical model for numerical simulation of two dimensional food freezing due to natural convection is presented. Fluid mechanics and heat transfer by natural convection between air and a solid food in a freezer are predicted along with the heat conduction inside a plate shaped food. The mathematical model used includes continuity, linear momentum and energy partial differential equations for air and the non-linear heat diffusion equation for the solidification of the water content in the food. Unsteady two dimensional results of the velocity and temperatures distributions in the air around the food obtained by using the finite volume method are presented. Experimental data of temperature variation in time for a slab shaped food inside a freezer are used to assess the accuracy of the mathematical model and the solution procedure. It is found that the use of this model, that does not require heat transfer coefficients, allows to simulate two dimensional freezing experiments in a slab shaped food by natural convection with maximum and mean deviations of 5.2% and 1.2%, respectively. Ó 2002 Elsevier Science Ltd. All rights reserved.
1. Introduction An accurate knowledge of freezing processes is needed in the design of refrigeration systems to assure high quality of foods at reduced cost. From the heat transfer point of view the task to be performed is to determine the transient temperature distribution in the food during the cooling process and to predict the freezing time. Such a simple task is by far not easy since there are at least four main source of uncertainties: (a) thermal properties of foods; (b) freezing heat transfer coefficients; (c) approximations made in the prediction method and (d) the mathematical model used. Freezing processes are difficult to predict because the non-linear nature of the diffusion equation, in which the food properties: density, specific heat, enthalpy and thermal conductivity vary with temperature. Experimental methods to measure and model thermal properties of frozen foods have been reviewed by Singh (1995). Estimated properties with the Box–Kanemasu method and experimental thermal conductivity and apparent volumetric specific heat of carbohydrate and protein rich foods have been found that do not agree *
Corresponding author. Tel./fax: +56-2-682-5498. E-mail address: [email protected] (N.O. Moraga).
well for methylcellulose and wheat gluten, Saad and Scott (1996). Predicted values from food composition data by using the SchwartzbergÕs bound water model for freezing point and bound water for meats, surimi, fruits, cheese and Tylose gels have been found to be within 10% of actual values, Pham (1996). Sensitivity of freezing time calculation for Alaska pollock surimi performed with a commercial finite element program have shown that: apparent specific heat, geometric dimension, heat transfer coefficient, fluid ambient temperature, thermal conductivity and density influence the calculations, De Qian and Kolbe (1994). The influence of change of volume expansion, thermal conductivity, specific heat, unfrozen water, initial freezing point on the predicted freezing curve with an analytical model has been reported by Hung and Thompson (1983) to account for differences of freezing times between 2% and 13%. The imprecise knowledge of the heat transfer coefficient h has been found to be the major error source in freezing time prediction methods, Cleland and Earle (1984). Small deviations in the heat convective coefficient were found, by the use of a finite element method, to result in large deviations in the core temperature of food chilling using air, Nicolai and Baerdemacker (1996). Spatial variations of h and the measurement techniques used in the experiments to determine h are
0260-8774/02/$ - see front matter Ó 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 2 6 0 - 8 7 7 4 ( 0 2 ) 0 0 1 3 5 - 8
18
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Nomenclature a b c Cp d fi g k n p t T u v x y a b q / U
length of the freezer chamber (m) height of the freezer chamber (m) length of the plate shaped food (m) specific heat (J/kg °C) height of the plate shaped food (m) functions to define time variations of freezer wall temperatures (°C) gravitational acceleration (m2 /s) thermal conductivity (W/m °C) normal direction to food external area pressure of air (N/m2 ) time (s) temperature (°C) air velocity in axial x-direction (m/s) air velocity y-direction (m/s) horizontal cartesian coordinate (m) vertical cartesian coordinate (m) under-relaxation factor coefficient of thermal expansion (1/°C) density (kg/m3 ) transport dependent variable water content in food on wet basis
two issues that affect the accuracy of the freezing prediction methods, Cleland, Cleland, and Jones (1994). The approximations commonly used in the predictions methods are: (a) that no changes of volume occurs in the food during freezing, even though density of ground beef changes about 4,8%; (b) weight losses in the range of 2–3% caused by water evaporation at the food surface are considered to be negligible; (c) the food is assumed to be an homogeneous isotropic material; (d) food texture is considered to be similar to the materials used to measure the heat transfer coefficient by the lumped system method. Physical, mathematical and computational aspects of freezing and melting processes that have been examined by Alexiades and Solomon (1993), indicated that the common practice used to construct a mathematical model has been based on solving the heat diffusion equation with boundary conditions of the third kind. The mathematical model used to find the time evolution of temperature distributions in the food has been based on the non-linear, second order, partial differential heat diffusion equation with heat convective boundary conditions on the food surfaces. Approximate analytical method to predict freezing times have been derived using shape factors, Salvadori, Michelis, and Mascheroni (1997), Cleland, Cleland, and Earle (1987), Cleland and Earle (1982) with percentage differences with experimental data in the range between 5% and 15%. Finite
l
dynamic viscosity (N s/m2 )
Subscripts A geometrical center of food a air F geometrical center of air in the freezer f food oa initial state in air of initial state in food p pressure ref reference state T temperature u air velocity in horizontal direction v air velocity in vertical direction / transport dependent variable Superscripts n number of iteration n 1 number of previous iteration new new calculated value guessed values for dependent variables 0 correction values for dependent variables
difference methods have been more often used to predict freezing times with errors of 10% for elliptical foods, Ilicali, Cetin, and Cetin (1996) and for rectangular meat patties, Rubiolo de Reinick (1996). A two-time-step finite difference method has been shown to have smaller maximum errors (7–19%) than the Crank–Nicolson method (16–37%) in the prediction of freezing of codfish, Saad and Scott (1997). Freezing of irregular shaped foods has been predicted by using finite element method by Comini, Del Guidice, Lewis, and Zienkiewiez (1974), Abdalla and Singh (1985), Caroll, Mohtar, and Segerlind (1996), Puri and Anantheswaran (1993), among others. This paper presents a new mathematical model that uses a numerical procedure for predicting heat conduction with solidification of the water content for solid foods, being cooled by natural convection in air inside a freezer. Instead of using a standard model of heat diffusion for the food with heat convection boundary conditions the mathematical model proposed replaces the external boundary condition for a set of four partial differential for the simulation of the natural convection in the surrounding air. The finite volume method is used to solve the system of discretized equations. Experimental measurement with thermocouples, scanner and a PC of time histories of temperature inside a standard package, used to test freezers, that has properties very close to the beef meat is used to validate the numerical solution procedure.
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
This original procedure has been used to simulate two dimensional unsteady heat conduction during freezing of a slab shaped piece of salmon meat under a forced convective heat transfer cooling, Moraga and Salinas (1999). Unsteady local heat transfer coefficients by forced convection have been successfully calculated at the four walls of a slab shaped food in forced air around it after the temperature fields in the food and in the air were determined, Moraga and Medina (2000). In these two previous papers the cooling mechanism was forced convection from the food to the surrounding air and the fluid motion of air is independent of the temperature field in the cooling fluid. In the present work the momentum equation in the vertical direction is coupled to the energy equation by the buoyancy forces that are calculated in terms of the temperature difference and the prediction of the solidification due to natural convection is more complicated. By the other hand, the heat transfer coefficients for natural convection are generally known with lower accuracy than in forced convection. Therefore, a procedure such as the one proposed that does not require the knowledge of heat transfer coefficient should be more interesting to use in food freezing processes by natural convection.
2. Physical situation and mathematical model A slab shaped food of height equal to d and length equal to c is inside a freezer of dimensions b by a initially at temperature Tof in contact with air at rest that has a uniform initial temperature Toa , as it is shown in Fig. 1. The refrigeration system is activated, and the freezer walls start to cool, along with, the air inside the rectangular cavity and later the food located on top of the bottom wall. In order to define the mathematical model the following assumptions are made. The air is assumed to be Newtonian fluid, that flows in a laminar flow. Unsteady natural convection in the air and heat conduction in the
19
food are assumed to be two dimensional. Thermal properties of the food change with temperature and a temperature dependent apparent specific heat is used to predict the phase change from liquid water to ice. Mass transfer of water inside the food and the change of volume have been neglected. Continuity, linear momentum in x and y spatial directions and energy equations describe the natural convection in air, ou ov þ ¼0 ð1Þ ox oy ou ou ou op o ou o ou l l q þu þv ¼ þ þ ot ox oy ox ox ox oy oy ð2Þ ov ov ov op o ov o ov q l l þu þv ¼ þ þ ot ox oy oy ox ox oy oy þ gbðT Tref Þ o qCp T o qCp T oT o oT þv ¼ þu k qCp ox oy ot ox ox o oT k þ oy oy
ð3Þ
ð4Þ
where according to BoussinesqÕs approximation air density has been assumed to change linearly with temperature. The values of the constant properties used for air were: q0 ¼ 1:287 kg/m3 , Cp ¼ 1:006 kJ/kg °C, k ¼ 0:024 W/m °C, l ¼ 1:71 105 kg/m s, b ¼ 3:5 103 1/K. The non-linear heat diffusion equation inside the food assuming temperature changes in the thermal properties is oT o oT o oT ðqCp Þf ¼ kf kf þ ð5Þ ot ox ox oy oy where a temperature dependent apparent specific heat Cpf has been used to capture the phase change of liquid water to ice inside the food, Hsiao (1985). Initially, air at rest inside the freezer and food are at uniform temperature: for
t ¼ 0;
u ¼ v ¼ 0;
Tair ðx; y; 0Þ ¼ Toa ;
Tfood ðx; y; 0Þ ¼ Tof ;
ð6Þ
Boundary conditions in the freezer include the no slip condition and temperature variations in time and space in all the internal walls of the freezer, u¼v¼0
Fig. 1. Physical situation of food inside freezer.
at x ¼ 0;
x ¼ a;
y ¼ 0;
y¼b
T ðx ¼ 0; y; tÞ ¼ f1 ðtÞ;
T ðx ¼ a; y; tÞ ¼ f2 ðtÞ;
T ðx; y ¼ 0; tÞ ¼ f3 ðtÞ;
T ðx; y ¼ b; tÞ ¼ f4 ðtÞ
ð7Þ ð8Þ
where the boundary condition of the fourth kind can be determined by experimental measurements of the wall temperature variations with time.
20
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Table 1 Thermal properties of food used in calculations 8 3 > < q ¼ 1053 kg=m T P Tcr : Cp ¼ 1:448ð1 UÞ þ 4:187U J=kg °C > : K ¼ 0:09705 þ 0:501U þ 5:052 104 UT W=m °C 8 > > >
1053 kg=m3 0:98221 þ 0:1131U þ 0:25746ð1 UÞ=T T < Tcr : > Cp ¼ 3:874 2:534U þ 902:893ð1 UÞ=T 2 J=kg °C > > : K ¼ 0:388 þ 1:412U þ 0:9575=T W=m °C
Conjugate boundary conditions are imposed between the air and the food. oTa oTf ðTa Þw ¼ ðTf Þw ; ka ¼ kf ð9Þ on w on w Table 1 shows the thermal properties of the material, assumed to be those of the beef meat, that were calculated in terms of the initial water content on wet basis U as reported by Salvadori and Mascheroni (1991).
3. Method of numerical solution The physical domain of air and food inside the freezer was discretized with different non-uniform grids of 27 19, 39 31, 51 36, 69 51 and 81 78 nodes. Temperature, pressure, specific heat and density were calculated at the grid points. A staggered grid was used to evaluate the velocity components and thermal conductivity on the faces of the control volumes (Harlow & Welch, 1965). Control volume faces were located between the air the food and the size of the grid was refined near the food walls. The governing Eqs. (1)–(5) were integrated in time and space according to the practice of the finite volume method. The diffusion terms of the viscous forces in the linear momentum equations and the heat conduction in energy equations were calculated by using linear interpolation functions between the dependent variables (/ ¼ u; v; T Þ and the independent spatial variables ðx; yÞ. Convective terms in linear momentum and energy equations were evaluated by the power-law differencing scheme (Patankar, 1980). The semi-implicit method for pressure-linked equations, SIMPLE algorithm was used to solve the governing equations (Patankar, 1980). The guess and correct procedure was initiated with a guessed value for pressure p , that was used to solve the discretized linear momentum equations to obtain the velocity components u and v . These values were introduced into the conti-
nuity equation to calculate the correction term for pressure p0 and for the velocity components u0 , v0 . The corrected values for the dependent variables were evaluated as (/ ¼ / þ /0 ), where /0 were the corrections for / ¼ u, v and p. An iterative procedure was used to find a convergent solution for pressure, velocity and temperature. The iteratively improved pressures, velocities and temperatures were calculated using under-relaxation: pnew ¼ p þ ap p0 ; / ¼ u; v; T
/n ¼ a/ / þ ð1 a/ Þ/n1 ; ð10Þ
where ap and all a/ ð/ ¼ u, v, T Þ were under-relaxation factors with the following values between 0 and 1, au ¼ av ¼ 0:03;
ap ¼ 0:3;
aT ¼ 0:9
ð11Þ
A line by line method, that results as a combination of the direct tri-diagonal matrix algorithm (TDMA) for one dimensional situations and the iterative Gauss– Seidel method with successive under relaxation was the numerical procedure used to solve the system of discretized governing equations (Versteeg & Malalasekera, 1996). The iterative procedure was ended when the discretized continuity equation was satisfied with an error smaller than 106 in every one of the control volumes. The numerical procedure used was based on the fact that the heat diffusion equation, Eq. (5), is an asymptotic case of the energy equation, Eq. (4). Therefore, the space occupied by the solid food was assumed to be filled with a fluid with an extremely large value of the dynamic viscosity, l ¼ 1020 , Pa s, and hence that on that region the numerical solution for the fluid mechanics problem gave as the result that the velocity components were zero. Temperature dependent values for density, apparent specific heat and thermal conductivity of the food were assigned to the nodes on that part of the computational domain and the energy equation was solved with the zero values for the two velocity components.
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
4. Experimental procedure A standard package of 1.0 kg, with dimensions equal to 0.05 m by 0.10 m in length, to test domestic freezers was used in the experiments. The package composition was 467.2 g of water, 230 g of oxietilmetil cellulose, 5 g of NaCl and 0.6 g of 6-Cl-m cresol. An extremely thin plastic membrane was used to seal the surface in order to avoid a loss of weight by mass transfer in its surface. Thermal properties of such a substance are very close to the beef meat. The package was put inside a freezer chamber of 0.49 m in length by 0.35 m in height, as shown in Fig. 1. Initially, the still air and the freezer walls were equal to 28 °C, while the package was at 26 °C. The refrigeration cycle was started and temperatures in the walls began to decrease in while the temperature data was recorded every 10 s until the solid reached a final state with a uniform temperature of 33 °C, after 18,000 s. Temperature data from 18T-type thermocouples (0.3 mm in diameter) located on the internal surface of the freezer walls as shown in Fig. 2, were recorded by a personal computer. Fig. 3 depicts the wall temperature histories
21
recorded during the first 8000 s, where the large temperature changes are seen to occur during the first 2000 s. The non-symmetric nature of the wall temperatures was caused by the direction of the refrigerant fluid flow inside the evaporator located on the walls of the freezer. The unsteady wall temperature data measured were used in the numerical simulation as the known boundary conditions on the walls of the freezer.
5. Results The program written to solve the system of equation according to the SIMPLE algorithm in finite volume was tested to produce numerical results that were independent of the grid used and of the time step. A total of five different meshes with different number of nodes were used. The physical domain that included the air in the freezer and the solid food was discretized by using five alternative non-uniform grids of 27 19, 39 31, 51 36, 69 51 and 81 78 nodes. These five grids were selected based on previous experiences solving conjugate a forced convection and food freezing problem, Moraga & Salinas (1999). Fig. 4 shows that the grid used with 81 78 nodes was more refined in the vicinity of the food walls where temperature and velocity gradients were larger. A total of three control points were chosen to test that the numerical solution was grid independent. Fig. 4 depicts the location of the two regions in the food that coincide with the position where the time history of food temperature was measured, TA and TB , and in the space occupied by the air, TF . The
Fig. 2. Position of thermocouples on freezer walls.
Fig. 3. Temperature history of freezer walls.
Fig. 4. Mesh used and thermocouples locations for air and food.
22
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
thermocouples TA and TF were located at the geometrical center of the food and of the freezer, respectively. The results of the temperature histories predicted in the geometrical center of the food and of the air in the center of the freezer for the five grids used in the numerical solution procedure are shown in Fig. 5. A time step of 1 s was used in these calculations. Temperature histories could be predicted in the food with grids of at least 51 36 nodes while 81 78 nodes were needed to calculate the temperature variations with time for the air independently of the grid size. Therefore it was concluded that a grid of 81 78 was needed to solve the discretized natural convection of air with freezing of the water content of the food located inside the freezer under time changing wall temperature boundary conditions. The effect of the time step used in the numerical results predicted for velocities and temperatures was also studied. The mathematical model defined by Eqs. (1)–(9) was discretized with the selected grid of 81 78 nodes. Four different time steps of 0.1, 1, 10 and 100 s were used to integrate the basic discretized equations in time. Fig. 6 shows the time evolution of temperatures calculated with these four different time steps at the geometrical centers in the food (TA ) and in the air (TF ). The calculations with Dt larger or equal to 10 s originates temperatures differences of 5 °C at the end of the freezing process in the food geometrical center and dif-
Fig. 5. Effect of grid size on temperatures measured at geometrical center of food (TA ) and of air in freezer (TF ), (Dt ¼ 1 s).
Fig. 6. Effect of time step on temperatures measured at geometrical center of food (TA ) and of air in freezer (TF ), (grid size: 81 78 nodes).
ferences of 4 °C in the prediction of the air temperature in the center of the freezing chamber. The results for the temperature histories calculated with different time steps were seen to be time step independent when Dt was equal to 1 s. The evolution in time of the air velocity vectors in the freezer is shown in Fig. 7 for 300, 1200, 2700, 5400, 9000 and 14,400 s. At time equal to 1200 s the center of the two main vortexes is seen to be in the lower position and two secondary non-symmetric flows are seen to appear on top of the food. Tertiary vortexes start to develop in the two lower corners for time equals to 2700 s. A more complicated fluid flow pattern with six vortexes is seen for times of 9000 and 14,400 s. The velocity at the center of center of the freezer reaches a maximum value of 0.5 m/s after 6500 s and it remains in this value for the next 6000 s until it begins to decrease. The velocity near the surface of the food changes with time and in space. The maximum values of the velocities occur at the two upper corners of the food and they are reached after 2000 s. Temperature distributions in the air and in the food are shown for six time instants, in Fig. 8. A uniform temperature equal to 26 °C is observed for almost all the lower half part of the freezer at time equal to 300 s. After 2700 s temperatures lower than 0 °C are noticed in the air that is in the half upper part of the freezer. Air temperature is seen to reach temperatures below 0 °C in all the freezer area for a time of 9000 s, while the liquid
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
23
Fig. 7. Time evolution of air velocity in freezer.
water had changed to ice in the lower half portion of the food. Finally at 14,400 s temperature in the food is seen to be below 20 °C. The non-symmetrical initial cooling process caused by the direction of the circulation of the refrigerant liquid in the external walls of the freezer last until 5400 s and from then on the isotherms tends to be almost symmetrical in y-direction. The validation procedure used include a comparison of temperature histories predicted in two locations in the food (TA , TB ) and one in the air at the center of the freezer (TF ) against measurement performed in the same locations with thermocouples. A comparison between predicted values and data measured with thermocouples for the temperature history at two locations in the food was performed and the results are shown in Fig. 9. The predicted freezing curve at two locations in the food are seen to describe the experimental results during all the freezing process. The predicted curves for TA and TB are seen to have an unusual increment around 14,000 s. At this instant of time, the temperature distribution in the food and in the surrounding air shown in Fig. 7 becomes
more uniform. The sudden changes in the temperature at TA and TB after 14,000 s seems to be caused by the sudden changes of the fluid temperature at the same time (see Fig. 6). Maximum deviation for the temperature at the geometrical center of the food TA and at the lower positions located at one fourth of the food height TB occurs near the beginning of the phase change process. The maximum and mean deviations of the predicted temperatures and the measured temperatures at the positions TA and TB inside the food were calculated as percentages of the total experimental temperature range (56 °C). The maximum and mean deviations for TA were 4.6% and 1.2%, respectively and 5.4% and 0.9% for TB . The accuracy in the numerical simulation for onedimensional food freezing of an infinite slab has been discussed by Saad & Scott (1997). They have found that finite differences with a two time step method has maximum percent errors for simulations of the temperature data measured in codfish during freezing in the range between 7.2% and 19% and mean percent
24
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Fig. 8. Unsteady temperature distributions in air and in food.
errors in the range from 1.8% to 6.4%. The temperature simulated using the finite difference method along with the Crank–Nicolson scheme showed larger discrepencies, with maximum percent errors between 16.5% and 37.4% and mean percent errors in the range from 7.8% to 13.5%. It is evident that the full two dimensional conjugate natural convection–conduction model solved by the finite volume method exhibits accurate results for temperature predictions of food freezing inside a freezer. Another advantage of the proposed method it is that the knowledge of the spatial and time variations of the heat transfer coefficient around the food surface are not needed in the calculations.
6. Conclusions A new mathematical model to simulate the unsteady two dimensional natural convection in air and heat conduction with solidification of water content in a slab shaped solid food located inside a freezer has been presented. This model eliminates the need of specifying a convective heat transfer coefficient and replaces the convective boundary condition by four partial differential equations for the non-isothermal laminar flow of air around the food. The numerical predictions of cooling curves in a slab of a piece of material with the properties of beef meat have been found to have mean and maximum percent deviations of 5.4% (3 °C) and 2.1% (1.2
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
Fig. 9. Experimental and predicted cooling curves for at the food geometrical center (TA ) and at 1/4 of food height (TB ).
°C) for simulations of the experimental data measured in the food with thermocouples. The finite volume method with the SIMPLE algorithm has been implemented to solve the discretized system of equations. A sensitivity analysis of the discretization process allows to recommend the use of a non-uniform grid with 81 78 nodes and time steps of 1 s to calculate the unsteady distributions of velocity and temperature in the air inside the freezer and the temperature variations in time of a slab shaped food.
Acknowledgements This work was supported by CONICYT (Chile) through FONDECYT 1000207 project and by Universidad de Santiago de Chile, project DICYT 019316 MB.
References Abdalla, H., & Singh, R. (1985). Simulation of thawing of foods using finite element method. Journal of Foods Process Engineering, 7, 273–286.
25
Alexiades, V., & Solomon, A. D. (1993). Mathematical modeling of melting and freezing processes. Washington: Hemisphere Publishing. Caroll, N., Mohtar, R., & Segerlind, L. J. (1996). Predicting the cooling time for irregular shaped food products. Journal of Food Process Engineering, 19, 385–402. Cleland, A. C., & Earle, R. L. (1982). Freezing time predictions for foods––a simplified procedure. International Journal of Refrigeration, 5, 134. Cleland, A. C., & Earle, R. L. (1984). Assessment of freezing time prediction methods. Journal of Food Science, 49, 1034– 1042. Cleland, D. J., Cleland, A. C., & Earle, R. L. (1987). Prediction of freezing and thawing times for multidimensional shapes by simple formulae, Part 2: irregular shapes. International Journal of Refrigeration, 10, 234–240. Cleland, D. J., Cleland, A. C., & Jones, R. S. (1994). Collection of accurate experimental data for testing the performance of simple methods for food time prediction. Journal of Food Process Engineering, 17, 93–119. Comini, G., Del Guidice, S., Lewis, R. W., & Zienkiewiez, O. C. (1974). Finite element solution of non linear heat conduction problems with special reference to phase change. International Journal Numerical Methods in Engineering, 8, 613–624. De Qian, W., & Kolbe, E. (1994). Analysis of food block freezing using a P-C based finite element package. Journal of Food Engineering, 21, 521–530. Harlow, F. H., & Welch, J. E. (1965). Numerical calculation of timedependent viscous incompressible flow of fluid with free surface. Physics of Fluids, 8, 2182–2189. Hsiao, J. S. (1985). An efficient algorithm for finite difference analysis of heat transfer with melting and solidification. Numerical Heat Transfer, 8, 653–666. Hung, Y. C., & Thompson, D. R. (1983). Freezing time prediction for slab foodstuffs by an improved analytical method. Journal of Food Science, 48, 555–560. Ilicali, C., Cetin, M., & Cetin, S. (1996). Methods for the freezing time of ellipses. Journal of Food Engineering, 28, 361–372. Moraga, N. O., & Medina, E. E. (2000). Conjugate forced convection and heat conduction with freezing of water content in a plate shaped food. International Journal of Heat and Mass Transfer, 43, 495–517. Moraga, N. O., & Salinas, C. H. (1999). Numerical model for heat and fluid flow in food freezing. Numerical Heat Transfer, part A, 35, 495–517. Nicolai, B. M., & Baerdemacker, J. de. (1996). Sensitivity analysis with respect to the surface heat transfer coefficient as applied to thermal process calculations. Journal of Food Engineering, 28, 23– 33. Patankar, S. V. (1980). Numerical Heat Transfer and Fluid Flow. Washington: Hemisphere. Pham, Q. T. (1996). Prediction of calorimetric properties and freezing time of foods from composition data. Journal of Food Engineering, 30, 95–107. Puri, V. M., & Anantheswaran, R. C. (1993). The finite element method in food processing: a review. Journal of Food Engineering, 19, 242–274. Rubiolo de Reinick, A. C. (1996). Average and center temperature vs time evaluation for freezing and thawing rectangular foods. Journal of Food Engineering, 30, 299–311. Saad, Z., & Scott, E. P. (1996). Estimation of temperature dependent thermal properties of basic food solutions during freezing. Journal of Food Engineering, 28, 1–19. Saad, Z., & Scott, E. P. (1997). Analysis of accuracy in the numerical simulation of the freezing process in food materials. Journal of Food Engineering, 31, 95–111. Salvadori, V. O., & Mascheroni, R. H. (1991). Prediction of freezing and thawing times of foods by means of a
26
N.O. Moraga, H.G. Barraza / Journal of Food Engineering 56 (2002) 17–26
simplified analytical method. Journal of Food Engineering, 13, 67– 78. Salvadori, V. O., Michelis, A. de, & Mascheroni, R. H. (1997). Prediction of freezing times for regular multidimensional foods using simple formulae. Lebensmittel–Wissenchaft und–Technologie, 30, 30–35.
Singh, R. P. (1995). Thermal properties of frozen foods. In M. A. Rao & S. S. H. Rizvi (Eds.), Engineering Properties of Foods (2nd ed.). New York: Marcel Dekker (chapter 4). Versteeg, H. K., & Malalasekera, W. (1996). An introduction to computational fluid dynamics. The finite volume method. Essex: Longman Scientific and Technical.