COOLERG A STAQK 9F FRUST: FAME?) BEE} ANALYSSS Thesis for the Degree of Pix 0. mamas; STAEE HNWERSEW CfiRLGS EORRERO ANGEL 1%? Ith‘b 0-159 TE LIBRARIES Immunitfifif “um WUNIH\|\\Hl1‘|\1\\||\1\| 3 1293 01765 5451 '9' - N I" "’ L I D r; z . .i: 1 Miehign ;. $1.3m Univ crsity This is to certify that the thesis entitled Cooling a Stack of Fruit: Packed Bed Analysis presented by Carlos Borrero Angel has been accepted towards fulfillment of the requirements for A degree in_E.QQd_S£_LBnc e L‘.' n Afihy Date May 4, I967 I JUL" 30 2000 V‘J ‘12) COOLING A STACK OF FRUIT: PACKED BED ANALYSIS AN ABSTRACT Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOC TOR OF PHILOSOPHY Department of Food Science 1967 ABSTRACT COOLING A STACK OF FRUIT: PACKED BED ANALYSIS by Carlos Borrero Angel This project consists of two phases an experimental and an analytical study. In the experimental study we measure: the temperature history of the particles in the packed bed with respect to radial and longitudinal position, the pressure drop across the packed bed and the convective heat transfer coefficient for three air velocities. The analytical study consists of the formulation of the equations governing the system, namely the equation for heat transfer within the spherical particle and that of the flowing fluid, their solution, and obtaining a method for presenting the cooling rate data. The pressure data indicates that the pressure drop across the bed of fruit is a function of the height in the bed. The comparison of the friction factor calculated using the pressure dr0p with that obtained from correlations in the literature show the experimental values of the pressure drop to be much smaller than those predicted probably due to the lack of consolidation of the particles in the bed. The convective heat transfer coefficient was evaluated using a transient method. The values obtained for h, for each air velocity, are in the range of the values reported in the literature for similar cases. The experimentally determined h values were used in the computation of the temperature history of the packed bed by the finite differenc e model. CARLOS BORR ER O ANGEL The backward difference method is used for the simultaneous solution of the two equations governing the system. The total error for the numerical solution of a single sphere is obtained and the finite difference model is compared with the theoretical solution for a low NBi system with no heat generation. The experimental and calculated temperature-time data were plotted on semilogarithmic paper and the values of fa and ja, the apparent temperature response parameter and lag factor respectively, obtained. The first degree polynomial lines for the values of fa and ja versus height were obtained for both the experimental and calculated data. The results indicate that the values of fa and ja vary with height in the bed, the slope of the lines are similar for the calculated and experimental data for each air velocity. The experimental and calculated results suggest that the slope of the line of fa versus height in the bed is a function of the air velocity while the value of the intercept of the line is dependent on the particle and system properties. The heat generation contribution was evaluated using the finite difference model. At low air velocities the contribution due to heat generation becomes significant and its value is a function of the height of the particle in the bed. COOLING A STACK OF FRUIT: PACKED BED ANALYSIS A THESIS Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Food Science 1967 a mis padres AC KNOWLEDGEMENTS The author wishes to express his appreciation for the guidance and friendship of Dr. I. J. Pflug (Food Science), a truly inspiring teacher and counselor, during the graduate program. Grateful appreciation is also extended to Dr. C. M. Cooper (Chemical Engineering), Dr. L. R. Dugan (Food Science), Dr. C. W. Hall (Agricultural Engineering), and Dr. A. M. Pearson (Food Science) for serving on the author's guidance committee. The author is indebted to the American Society of Heating, Refrigerating and Air Conditioning Engineers (ASHRAE) for supporting project RP-42 of which this study is a part. I wish to express particular thanks to Mr. R. S. Buchanan, Director of Research for ASHRAE, and to Mr. E. J. Robertson, Chairman of the ASHRAE Task Group on Food Science and Refrigeration, for their interest in this project. Thanks go to Dr. B. S. Schweigert, Chairman of the Food Science Department, for making available the graduate assistantship. The author is particularly indebted to Dr. A. J. Kopelman (Food Science) for his friendship, encouragement and constructive criticism throughout the study and for his patience in editing the manuscript. Grateful acknowledgement is extended to Mr. Bruce D. Eder for reviewing the manuscript and to Mr. E. Gil and Mr, R. K. Fisher for their help in the laboratory. ii TABLE OF CONTENTS Page LIST OF TABLES v LIST OF FIGURES vi NOMENCLATURE viii 1. INTRODUCTION, LITERATURE REVIEW AND OBJECTIVES . . ........ . . . ........... l 2. THEORY O O 0 O O O ......... O O O O O O O O O O O 2.1. Differential equations for the packed bed system ................. . . . . Z. 2. Proposed method of solution . .......... . 2.3. Analysis of spheres . . ..... . . . ....... 21 2. 4. Generalized numerical solution for three one-dimension heat transfer geometries ...... 24 3. EXPERIMENTAL DESIGN AND PROCEDURE . . . . . . 29 3.1. a) Packed bed system considerations and deSign O O O O O O O O O O O O O O O O O O O O O O O C 29 . b) Description of the equipment ......... . 3O 3. 2. Spherical transducers . . . . . . . . ....... 33 3. 3. Measurement of pressure drop in the packed bed . 36 3. 4. Measurement of air velocity through the packed bed ............. . ........ 36 3. 5. Materials for the experiment ............ 37 3. 6. Temperature measurement , . . .......... 37 3. 7. Procedure for a single run ......... . . . . 38 4. BASIS OF ANALYSIS .................... 43 4.1 Physical and thermal properties of product and cooling medium ............ 43 Heat of respiration . . . . . . . . . . . . . . . 46 . Convective heat transfer coefficients. . ...... 47 Pressure drop calculations . . . . . . . . . . . . . 51 Calculation of the flow rate of air 159:“? «11an through the bed . . . . . . . . . . ....... . . 57 4. 6. Evaluation of the total error of the numerical solution . . . . . . ........... 57 5. RESULTS AND DISCUSSION ......... . ...... 6O 6. CONCLUSIONS . . . ............. . ...... 83 iii REFERENCES . ....... APPENDIX Al. A2. A3. EXPERIMENTAL AND CALCULATED DATA . . . PROGRAM PAC KB ED 3 PR OGRAM GEOMETR Y iv Page 84 88 92 95 LIST OF TAB LES Physical Properties of Metal Transducers . ..... Location of Center of Transducer in the Longitudinal Direction for the Three Heights Of PaCked Bed 0 9 o o o o o o o o o o o o o o o o o o 0 Heat Transfer Coefficients Obtained Experimentally and Calculated from Equation 4. 3. 2 for 2. 75 inch COpper and Aluminum Spheres . . . . . . . . . . . . Experimental Film Heat Transfer Coefficient, h, BTU/hr OF ft2 for the Air Velocities Tested . . . Coefficients for Equation 4. 4.12 versus Bed Height 0 o o o o o o 0 Q o o o o O o o o o 0 Q 9 o O Tabulation of the Values of f and j Obtained from the Numerical and the Exact Solution . . . . . Slope, A’, and Intercept, B', for the Lines of the Values of f versus Height for the Experimental and Calculated Data for the Three Air Velocities . . Values of f from the Experimental Data for Test 217; \felocity = 37 ft/min . . . . ........ Slope, A', and Intercept, B', for the Lines of the Values of the Experimental and Calculated Data ofj versus Height for the Three AirVelocities..................... Temperature of a Sphere Within the Packed Bed from Program PACKBED 3 and the Experimental Results; Values Calculated Using fa and ja arealsoShown.............. ....... Comparison of the Values of T -T/To- Tco Calculated Using the Finite Difference Model and those Obtained from the Charts of Schumann (1929) . Experimental Pressure Data . . . . . . . . . . . . . fa and ja Data from Experimental Model . . . . . . . fa and ja Data from Finite Difference Model . . . . . Page 33 35 48 50 53 59 63 73 79 88 89 91 Figure 3.1 3.2 3.3 3.4a 3.4c 4.1 4.2 5.1 5.2 5.3 LIST OF FIGUR ES Page Schematic diagram of the experimental system . . . 31 Test section. A.l. View of test crate. 2. Side view of test crate. B. 1. View of crate support. 2. Side view of crate support. . . . . . . . . . . . . 32 Copper and aluminum spherical transducers A. Copper sphere B. Aluminum sphere ..... 34 Apple (plain spheres) and transducer (textured spheres) location in the experimental three crate(45inch)packedbed .. . . . . . . . . . . . . 39 Apple (plain spheres) and transducer (textured spheres) location in the experimental two crate(30inch)packedbed . . . . . . . . . . . . . . 40 A. Apple (plain spheres) and transducer (textured spheres) location in the experimental one crate (15 inch) packed bed. B. Top view for all crates of apple (plain spheres) and transducer (textured spheres) location. . . . . . 41 Relationship between overall pressure drop and kinetic energy for the flow of air through the bed Of fruit 0 O O O O O O O O O O O O O O O O O O O O O 54 Relationship between friction factor and modified NRe for the packed bed of fruit ........... 56 Experimental (solid line) and calculated (broken line) apparent temperature response parameter versus height of the fruit in the bed for 245 ft/min air velocity ............ . ......... 64 Experimental (solid line) and calculated (broken line) apparent temperature response parameter versus height of the fruit in the bed for 173 ft/min air velocity . ..... . . . . . ......... . .. 65 Experimental (solid line) and calculated (broken line) apparent temperature response parameter versus height of the fruit in the bed for 37 ft/min airvelocity.. ........... 66 vi Figure Page 5. 4 Relationship between the slope of the line of fa versus height and the air velocity for the calculated (broken line) and experimental (solid line) systems . ....... . . . . . . . . . 6‘7 5. 5 Experimental (solid line) and calculated (broken line) apparent lag factor versus height of the fruit in the bed for 245 ft/min air velocity _. . . . . 70 5. 6 Experimental (solid line) and calculated (broken line) apparent lag factor versus height of the fruit in the bed for 173 ft/min air velocity . . . . . 71 5. 7 Experimental (solid line) and calculated (broken line) apparent lag factor versus height of the fruit in the bed for 37 ft/min air velocity . . . . . 72 5. 8 Relationship between the integrated heat of respiration and height of fruit in the bed ..... 74 vii NOMENCLATURE surface area of body; Ac , cross sectional area of packed bed; Ali’ A21, A4i surface area terms in the finite difference solution, ft‘Z constant in backward difference equation, dimensionless constant in equation Y = A' Z + B' interfacial area per unit volume of bed, l/ft total particle surface per volume of particles, l/ft constant in backward difference equation, dimensionless constant in equation Y = A' Z + B' empirical constant, dimensionless constant in backward difference equation, dimensionless heat capacity; Cpa , heat capacity of air, BTU/1b 0F particle diameter, ft. empirical constant, dimensionless constant in backward difference equation, 0F temperature response parameter; time required for the asymptote of the cooling (or heating) curve to cross one log cycle, the time required for a 90% reduction of temperature on the linear portion of the cooling or heating curve, fa, the apparent temperature response parameter, hr. constant in equation 4. 4. 9, dimensionless friction factor, dimensionless 2 gravitational constant, ft/ hr mass flow rate, lb/ftz hr viii x g l" Re Pr "U222 0.0 surface film coefficient, BTU/hr oF ft2 lag factor (Ta - T/To - Tl); jC , lag factor at the geometric center; jS , lag factor at the surface; jm' lag factor for the mass average; ja' apparent value of the lage factor for the packed bed, dimensionless drag force, 1b. Chilton-Colburn (1934) heat transfer factor, dimensionless thermal conductivity; ka thermal conductivity of air; kao thermal conductivity of air at O 0C, BTU/hr OF ft length of packed bed, ft molecular weight, lb/lb mole Biot Number, h R/k, dimensionless Reynolds Number, U D p/p. , dimensionless Prandtl Number (cpu/k)a , dimensionless pressure, lb/ftz heat flux, BTU/ftz heat generation term; Qi , heat generation term at node i, BTU/ft3 hr distance from center to point of measurement in a sphere, cylinder or slab, ft radius of sphere or cylinder; half the thickness of a slab, ft hydraulic radius, ft temperature; To' initial temperature; T,Jo , inlet air temperature; T , cooling medium temperature; Ta' air 1 temperature surrounding the sphere; Ti j' temperature of node i at time j, oF ix N><<3C3"’ 9 time, hr velocity; < U> , average velocity, ft/hr volume, ft constant in equation 2. 2.1, dimensionless height of bed, ft thermal diffusivity, k/p cp, ftZ/hr roots of equation 2. 3. lb, dimensionless fractional void space, dimensionless time displacement, Z/U, hr o At/(Ar)2, dimensionless angle in spherical coordinates,equation 2.1.1, radians constant in equation 2. 4. 6, dimensionless viscosity of air, lb/ft hr 3.14159 . . . , dimensionless density, 1b/ft3 angle in spherical coordi.nat68,equation 2.1.1, radians constant in equation 2. 3. 5, dimensionless 1. INTRODUCTION LITERATURE REVIEW AND OBJECTIVES The rapid reduction of the temperature of fruits and vegetables after harvest is an important step in maintaining the quality of these products. Rapid cooling before storage, so called precooling, has primarily a twofold function: (1) to slow down the metabolic process by retarding respiration; and (2) to arrest the development of decay- producing microorganisms. Various means are available to cool fresh produce, however, there is a large amount of interest in establishing a technique that will result in fast, effective and efficient heat removal from a variety of products. There are many factors influencing the rate of cooling. For example; size, shape, composition and surface heat transfer characteristics bear a significant relationship to the rate of heat removal. Surface heat transfer coefficient is influenced by the type and flow of the cooling medium as well as the geometry of the particle and the system used, while internal heat transfer is strictly a function of the physical and thermal properties of the particle and its geometry. The cooling system used is usually dependent on the physical and physiological properties of the product. Heat exchange from an object to a cooling fluid in motion is a combined conduction convection process where the controlling resistance to heat transfer may be the surface, low N system, Bi or the interior of the object, high N system. Bi Theoretical studies of cooling rates of single fruit have been made by Thevenot (1955), Hicks (1955) and Kopelman ital (1966). They have evaluated the effect of convective heat transfer coefficient, size and geometry of the object on the cooling rate. Cooling indices such as cooling rate, cooling coefficient, half cooling time as introduced and used by Sainsbury (1951), Guillou (1960) and Thevenot (1955) are based on the principle of Newtonian heating or cooling, low N system. This assumption requires Bi that the resistance to heat transfer be at the surface, which is not accurate for most fruit cooling systems. Since the value of the thermal conductivity of fruit or produce is small, k < 1, the NBi is in the moderate to high range. Pflug and Blaisdell (1963) and Kopelman e131 (1966) use the temperature response parameter, f, or time required for 90% reduction of temperature on the linear portion of the cooling curve, and j, the lag factor to define the cooling rate. This method of presentation of cooling data takes into account the internal resistance to heat flow and incorporates the lag in temperature response due to this resistance. Both Newtonian heating or cooling systems and high NBi systems may be presented in this form (section 2. 3). Kopelman and Pflug (1967) have shown the properties of f in the . ,' . . . > high and low NBi region to be. In the high NBi region, NBi 50, f was shown to be independent of N or h; proportional to o , Bi and proportional to R2 . In the low N region, N . < 0. 2, f was Bi B1 shown to be inversely proportional to N or h; proportional to pc Bi and proportional to R. These relationships are extremely important in the design of transient heat transfer apparatus and consequently are applicable to the design of fruit cooling systems. In this review we will consider only air cooling of produce. We can divide air cooling of produce into two parts, low velocity as in cold storage rooms and high velocity which is used in pre- cooling systems. In low velocity systems such as refrigerated storages the produce in lug boxes or crates are placed in a room and air is circulated in the room by fans, located at the ceiling or floor of the room. A major concern in a storage room of this type is to obtain uniform flow of air throughout the room. Unfortunately in large storage rooms uniform flow is very difficult to achieve because of gravity forces, temperature gradients and the deflection of the air stream by physical barriers. Many different designs of storage rooms have been tried and evaluated (Rostos, 1960) and there appears to be a need for a design which will insure the maximum air flow through the crate of fruit. Several different arrangements have been studied Rostos (1960) suggests a room with a grating in the floor where the air was circulated through the fruit to an evaporator and then returned to the floor by means of a plenum. Gurevitz (1966) suggested a similar design for apple storages where the air is also circulated through the fruit, the crates would be packed tightly against one another and the bottom of the crates would have a 50% porosity, which would approximate the porosity of a stack of apples. Gurevitz (1966) emphasizes the importance of having no channels for flow of the air to short circuit since the primary purpose of this design is to have the maximum air flow through the stack. High velocity air cooling systems are generally used for precooling the fruits before market or storage. Guillou (1960) outlined a method of forced air cooling whereby air was forced, under static head differential, across the fruit in containers of various types and varying degree of opening. Even though he does not list box dimension some interesting observations can be made from his data. From Figure 15 (Gillou, 1960) a plot of the static pressure head versus half-cooling time with the number of containers the air must pass through as a parameter, we can see that the half-cooling time, the time to achieve a 50% reduction in temperature, for grapes in lidded lugs for a static head of 4 inches of water using constant inlet temperature air is 0. 45, 0. 9, l. 5 and l. 9 hours for l, 2, 3 and 4 tiers, respectively. These data indicate an approximately linear relationship between the half-cooling time and the length of air travel. Graphs for other products suggest a Similar relationship. Some precooling studies have been made using a packed bed type of system. The packed bed arrangement allows the particles to have the maximum exposure to the air column passing through the box. Soule _e_t__a_l_ (1964) and Bennett 3211 (1966) have reported a bulk cooler for citrus fruits which was a lug box through which air was passed. This system is a packed bed, however, in their analysis they oversimplified their system by assmning the air temperature along the column of fruit to be the constant inlet air temperature. Soule 3311(1964) and Bennett e_t_al (1966) present their data in the form of the coefficients for the third degree poly- nomial fit of their cooling curve. Hood (1967) studied the precooling of pickling cucumbers using a system similar to that of Soule iii (1964). Hood presented his data in the form of arithmetic plots of the fruit center temperature versus time measured at three depths and three cucumber diameters. His data indicated that the cooling time increased with depth. Bakker-Arkema and Bickert (1966) studied cooling of sugar beets in a packed bed system using very low air velocities (< 10 ft/min) where the cooling system was assumed to be of a Newtonic nature and with no heat generation. The packed bed system has been widely studied for both heat and mass transfer applications. The majority of the papers available deal with the treatment of the packed bed for heat transfer in catalytic reactors. These systems are treated as being low NBi systems where particle resistance to heat transfer is negligible compared to the exterior film resistance. Analytical treatment of such systems were made by Schumann (1929), Furnas (1930, a, b), Wilhelm gt 2.1. (1943), Klinkenbeerg (1948, 1954), Coberly and Marshal (1951), Amundson (1956), and others. Schumann (1929) and Furnas (1930 b) present the solution in dimensionless form for the temperature of the particle as well as the fluid as a function of time and position. Beside the mathematical presentation, the solutions are presented in the form of curves covering a wide range of values. Schumann's (1929) classical paper has been the basis for a large number of analyses in the area of heat transfer in packed beds. The most exhaustive review to date of the literature concerning packed bed systems is that of Barker (1965). Barker also lists some of the many correlations for evaluating the surface film heat transfer coefficients in packed bed systems. Barker (1965) presents correlation equations from the literature for N versus Re J the Chilton-Colburn (1934) heat transfer factor, which for H, randomly packed spheres show a forty fold spread in values of J H' Some of the significant correlations in this area are those of Furnas (1930 a), Gamson _e_t_a_l_ (1943), Leva (1947), Bird SEQ (1960), De Acetis and Thodos (1960), Bradshaw and Myers (1963), and McConnachie and Thodos (1963). Convective heat transfer coefficients at the wall of the packed bed and effective thermal conductivities of the packing have been studied extensively for catalytic converters and regenerators. These systems are evaluated as low NBi systems because of the small particle size or high thermal conductivity of the packing. Some of this work is reported by Wilhelm _e__t_ 31(1943)’ Coberly and Marshall (1954), Amundson (1956), Hill and Wilhelm (1959), Kunii and Smith (1961), and Yagi $2.91 (1961). Flow through packed beds has also been widely studied. Generally there are two approaches used for studying pressure drops in packed beds. In one method the packed bed is considered as a bundle of tangled tubes of weird cross section; the second approach is to consider the bed as a collection of submerged objects. The tube and bundle concept has been more popular and we will only discuss this concept. Correlations for the friction factor for wide range of N particle dimensions and fluids are Re’ available from the literature. Some of these papers are: Leva (1947), Leva and Grummer (1947), Brownell and Katz (1947), Ergun (1952), McConnachie and Thodos (1963), and Wentz and Thodos (1963). Most of the papers available in the area of packed beds do not consider the interparticle conduction, nor the heat generation of the particles in the bed. We are interested in taking a column of fruit and evaluating, using a packed bed approach, the behavior of such a bed, measuring the variables of the system such as pressure drop, temperature gradients, and convective heat transfer coefficients. This study consists of two parts, an experimental and an analytical section. The objectives of the experimental section are: 1. To measure the pressure drop through the bed, specifically to see if existing correlations for packed beds can be used to predict and design such systems for fruit cooling. 2. To evaluate, experimentally, the convective heat transfer coefficient using a transient method for the various air velocities. 3. To measure the temperature history of the particles in the packed bed as a function of height in the bed, and air velocity during transient cooling. The objective of the second phase of this study is to formulate the equations governing the system, namely, the equation for heat transfer within the spherical particle and that of the flowing fluid and to solve these equations. The analytical model developed can then be used to predict the temperature history of the particles in the packed bed. Tfluarnodelshouldincorporatethose humorS'which had been ignored in previous studies such as interparticle conduction and heat generation, and should be of such a general nature that particle size and properties, air velocity, heat generation and heat transfer coefficient are variables in the solution. A method of expressmg the cooling rate of the particle as a function of height in the bed is to be determined which will allow for the calculation of particle temperatures as a function of time. The cooling rate should incorporate the particle and system properties. 2. THEOR Y 2.1. Differential equations for thepacked bed system Heat transfer studies of packed bed systems have been made primarily for low N systems, N . i 0. 2. In these systems, B1 B1 the major resistance to heat transfer is in the exterior film due to either the very small dimension of the particle or its high thermal conductivity. In some cases where the NBi of the system is moderate or high this type of analysis is not applicable and intra- particle conduction must be considered. In a packed bed model for cooling produce the system NBi is in the range of l to 10; suggesting that neglecting the temperature gradient within the particle is not a realistic assumption. Since we can expect the major portion of the heat exchange to be from the particle to the fluid flowing past it, we will assume no heat transfer in the direction perpendicular to flow of the fluid; that is to say we will have a one dimensional heat transfer system. Thus we can say that we have a radially adiabatic packed bed where the transfer of heat only takes place from the particle to the fluid. We will consider a packed bed whose particles are spheres having a radius R and where the fluid is introduced into the bed at Z = 0, at constant temperature, TOO . We will neglect particle to particle heat transfer and the temperature range of 40 - 100 0F is assumed to be low enough for the radiative heat transfer contribution to be negligible. Because the particles in the bed 10 are biological materials and are respiring there will be a heat generation contribution, Q . The differential equation for conduction heat transfer in a sphere is: 33.1".-k_1__3__ 211+...1 3’. 'nefl pcp at ‘ r2 at r at 2 as 51 as r sine 2 + 1 a T +0 (2.1.1) rZ sin2 0 8(1) 2 where Q is the heat generation term, BTU/ft3hr Since the particle dimension is small in comparison to the length of the bed in the direction of the flow, we are assuming constant fluid temperature surrounding the particle (sphere). This assumption will provide the particle radial symmetry. The equation for heat conduction in a sphere with radial symmetry is: —....-.—2- ar +Q :2 ‘3ch (2.1.2) The initial and boundary conditions for the sphere are: Initial condition T:T attzOf‘or OErER (2.1.3a) B ounda ry c onditions 321 z .13(T_T) for t> o (2.1.3b) 8r k a r:R g—E: 0 at r=0 for t_>_0 (2-1-3C) 11 where T is the temperature of the sphere and Ta that of the air. A heat balance of the fluid element surrounding one sphere will give the following differential equation. BTa 8Ta hA(T'Te) “U PecpaV‘s‘z“ = Vpacpast“ (“-4) where h = the convective heat transfer coefficient, BTU/F hr ftz U 2 the veloc1ty of the fluid past the particle, ft/hr Z = the height of the bed, ft A 3 . 0 pa, Cpa = the den51ty, 1b/ft , and spec1f1c heat, BTU/lb F, of air respectively Ta : air temperature, OF V = volume of air surrounding the particle, ft3 The initial conditions for the fluid can be specified as T :T at t=0 and 0: Z_<_L (2.1.5a) a 0 we can also specify the boundary condition that the inlet air temperature is constant and equal to a given value: TazTOO at 2:0 and t> o (2.1.5b) Is is now obvious that we have two one -dimensional heat transfer equations which must be solved simultaneously. For most cases the equation for conduction heat transfer in a sphere may be solved analytically, however, in our case we find that because of the presence of the heat generation term, Q, and the variable 12 surface temperature, equation 2.1.2 becomes a nonlinear differential equation and for this particular case there is no exact solution available. The general plan for the solution of equations 2.1. 2 and 2.1. 4 is to write them in finite difference form and to solve these two equations simultaneously. 2. 2. Proposed method of solution Based on the previous discussion (Section 2.1), the solution of equations 2.1. 2 and 2.1. 4 for the packed bed system will be made using the following as sumptions: 1. We have a radially adiabatic fixed bed; with no inter- particle heat transfer. 2. There is no radiative heat transfer. 3. The temperature surrounding each particle at any given time instant is constant. 4. Within the sphere there is radial symmetry, i. e., one dimensional heat flow. 5. The bed is ordered and the spherical particles rest one directly above the other. We can solve these equations for either a coupled or uncoupled system. In the coupled system we evaluate the temperature for each longitudinal position by freezing time and then iterate with respect to time. In the uncoupled system we calculate the temperature history of each sphere along the longitudinal axis and using the temperature history of the air surrounding the previous position to move up the column to the next sphere. 13 We first tried to solve the coupled system using an explicit numerical method. This method, the forward difference approximation, discussed in Schneider (1957), Dusinberre (1963), and others, was used to evaluate equations 2.1. 2 and 2.1. 4. It was found that the stability criterion for this solution required that the time increment be of the order of 10"5 hrs. The forward difference approximation, although adequate for the solution of the packed bed system is not efficient due to the magnitude of the time increments. To overcome this limitation an implicit finite difference method was used. (In an implicit method the time step may be much larger than in the explicit method because this method is unconditionally stable for any positive time increment.) An implicit method, the backward difference approximation, is similar to the forward difference approximation except that the temperatures are evaluated at the future time j + 1 . The backward difference method required the simultaneous solution of the equations for the temperature of all the nodes at time j + l . The Gauss elimination method described in Smith (1965) and Beck (1966) was applied to solve the tridiagonal equations obtained. The solution of these equations by the backward difference approxi- mation requires 40% more computer time (Beck, 1966) than the forward difference method, however, because of unconditional stability, the time interval can be larger and the total computing ~ time can be decreased. Generally the forward and backward difference methods have approximately the same accuracy (Beck, 1966). 14 The development of the finite difference equations will be presented in four parts. First, the evaluation of the fluid temperature surrounding the sphere; second, the convective boundary condition for the sphere; third, the calculation of the internal nodes of the sphere; and fourth, the calculation of the temperature at the geometric center, r = 0 . Before entering into the development of the equations for the various positions there are several points which should be taken into account. In the solution of the finite difference equations for the sphere where there is a volume and area change, the tridiagonal constants for the evaluation of each node vary. In order to simplify and to clarify the presentation of these equations the area constants are presented as Ali’ AZiand A41; the i subscript refers to the node. The computation formula of these constants for each node can be found in the program PACKBED 3 in Appendix (2). In the equations below all temperatures are designated as Ti,j where the i and j subscripts refer to position and time, respectively. The evaluation of the temperature 2f thigas surroundijithe ghere. The first longitudinal node was evaluated using the constant inlet air temperature, T . The succeeding nodes were evaluated co considering Ta to be the calculated temperature for the air in the previous element. As the height in the column is increased there is a time displacement in the calculation of the temperature history of each sphere and this displacement is n = Z/U where 15 Z is the height of the sphere being evaluated and AZ is the diameter of the sphere. —_— i U V T t pa cpa AZ i, j+1 Az —. —. , .Y. - 1 pa Cpa At(Ti, j+l Ti,j) 1 _Y_ U pana AZ Ta An energy balance over the longitudinal segment AZ gives: heat in by convection + heat in by flow 2 change in internal energy of segment XpacaUV P ( -T ) hA1i(Ti+l,j+l-Ti,j+l) + AZ Ta i,j+1 pac a V = +(T1JH' Ti’j) (2.2.1) where X = 1 or 2 and is dimensionless This equation holds for i = 1. For the entrance region it is important to obtain a realistic approximation for T. In 1,j+l' considering the energy output by flow we must consider that, in reality, for the first longitudinal segment AZ = D/Z and for all position other than the first AZ 2 D. We take care of this contingency by using the variable X whose value is equal to l for all spheres but the first where it is equal to 2. Taking equation 2. 2.1 and collecting terms we get: (pc UV) p c V _ _ a pa _ a a Ti,j+1 hAli X Az "KER" p c U V p __ a pa _ a a + Ti+l,j+1 ”1 Ali) ‘ ' X Ta ( 12 > K: Ti,j 16 dividing both sides of equation (2. 2. 2) by h Ali to make the equation dimensionless and separating terms we obtain the constants Ai' Bi’ C1 and Di in the equation Di : Ai T. 1-1,j+1 + Bi Ti,j+1 + Ci Ti+l,j+l Ai = o (2.2.3a) ’pacpa U V pac a V 13i = 1+x hAliAZ + ——R——AthA1i (2.2.3b) ci = -1 (2.2.3c) X Ta (pacpa U V) D1 = AZhAli (2.2.3d) The backward difference equation for the first node is then: Di : Bi Ti,j+1 + Ci Ti+1,j+1 (2°2'4) Convective boundary condition for the sphere air sphere \\ q 4—(9 o ®———j+1 lr-o— 2 I : l I | I O- — — — l- — —J «Ar—+4—Ar—p : : l I ,l . 1-1 1 1+1 In the cases where the surface resistance to heat transfer is finite and cannot be ignored, the condition at the boundary can be formulated in the following manner, when the surface of the sphere is cooled by a known heat flux. Equation 2.1. 2 may be written in the following form 17 ET _ 1 a ZBT -Qi+-é-£-—k :Z-E-r-(r 5?) (2.2.5) = h(T. 1,j+1 " Ti-1,j+1) (2'2'6) i, j+1 for the boundary condition Ti is the temperature of the air 1, j+l surrounding the sphere. Taking the right hand side of equation 2. Z. 5 and expressing it in finite difference form we obtain an equation where the heat flux removed by convection can be used: r- H 4n r2 g—T- - 4Tr r2 g—T-l it); _a__ r23: _ k ri+l,j+1 ri,j+1 411' 2 8r 8r _ 2 Ar r 41rri L. ._J (2.2.7a) where T. . - T. . 411' r2 35-:— = A1i IH’JHM 111“ (2.2.7b) 1+l,j+l hAZ. 2 8T 1 411' r -- = T. . - T. . (2.2. 7C) 8r 1’ j+1 k 1, J+1 1—1,J+1 and Ti .+1-Ti . 4Trr pc —- =A4ipcp )JAt Ll (2.2.7d) ','+1 18 Substituting equations (2. 2. 7 b, c, and (1) into equation 2. 2. 5 and multiplying through by r (Ar)Z/r (Ar)‘2 k we obtain A4101 A4i -_____— + —— T. . -T. . k 9(Ar)2 ( 1’3“ I’J) Ali A2i = (2:)? (Ti+1,j+1 ' Ti.j+1)- NBi TX: (Ti.j+1 - Ti'IJH) (2.2.8) where _ h r NBi- k and 9 _ (1 At ‘ ——2 (Ar) After separating terms in equations 2. 2. 8 we can obtain the constants A., B., C. and D. 1 1 1 1 NBl A21 A1 = - —;—-A—r- (2.2.93) A1i A41 NBi A21 Bl : (Ar)2 + m + W (2.2.913) Ali Ci : - (Ar)2 (2.2.9C) D 2 LiQi + —_A4i T <2 2 9d) 1 k 9 (Ar)2 i,j . . And using these constants the backward difference equation for node i = 2 is given by Di 2 A1 Ti-1,j+1 + Bi Ti,j+1 + Ci Ti+1,j+1 (2.2.10) 1‘? Calculation of the temperature of the internal nodes of the sfliere The finite difference equation for the inner nodes of the sphere (R> r> O) are obtained usmg equation (2.2. 5) for i > 2 . For the inner nodes we can say that (T - T. .) 4w r2 p c g3; = p c .A4. 1vl+1 1'1 (2.2.11a) p at . 1 At 1,J+1 . (T. . - T. . ) k 4w r2 g: = k A2i ”1'ng “J” (2. 2.11b) r 1+Lj+l r T. . - T. . ) k%ufi§£ = kAL “fil ““34 (221M) r . . 1 Ar 1,J+1 replacing these terms into equation 2. 2. 5 and multiplying by (Ar)2/(Ar)2k we get: A4i Qi A4i - ——-——-——- + —-—-— (T. . - T. .) k 8 (Ar)2 1,J+1 1,] A2. A1. = ——4L(T. . -11 . )- ——i-(T.. -'r. . ) (Ar)2 1+1,J+1 1,_]+1 (Ar)2 1,_]+1 1-1,J+1 (2.242) Separating terms in equation 2. 2.12 we obtain the constants for the backward difference equation for the internal nodes. Ali A1 : ' (Ar)2 A41 Ali B1 : 9(Ar)2 + (Ar)2 + c. = _ A2i (2. 2.13a) (2. 2.13b) (2. 2.13c) A4i Qi A4i ' D. = —----- T. .+ -——--- (2.2.13d) 1 6 (Ar)2 1,] k The backward difference equation for the inner nodes (i > 2) is then: Di = A1 Ti-1,j+1 + B1 Ti,j+1 + Ci Ti+1,j+l (2°2°14) Calculation of the center node for the sphere When r = O for one dimensional heat conduction in a sphere BT/ar = O , which occurs when there is symmetry around the origin, and where 2/r (GT/8r) assumes the indeterminate form 0/0 . By McLaurin's expansion T'(r) = T'(O) + r T"(O) +-;: r2 T"’(O) + where T' = aT/ar but T'(O) = 0, so the limiting value of Z/r (GT/8r) as r —’ O is the value of 2 EizT/ar2 . For r = 0 equation 2.1. 2 becomes 2 8T _ ar p cp }? - 3k arz + oi (2.2.15) Expressing equation 2. 2.15 in finite difference form and multiplying by At/p cp we obtain: At - . . - . . = 9 .. . - ‘ . . . . p cp Q1 + (Ti, J+1 T1,J) 3 (T1+1,J+1 2 T1,3+-1 + T1-1, 3+1) (2.2.16) At r = 0 we find that Ti+l,j+1 = T14“).+1 and equation 2.2.16 can be simplified to: 21 At - p Cp Qi + (Ti,j+1 - Ti,j) : 6 e(Ti-1,j+1- Ti,j+l) (2.2.17) Separating terms to obtain the coefficients A. = -6 9 1 B. = 1 + 6 9 1 (2.2.18) C. = O 1 D. = TF..-+Cl At 1 1,] 1 p cp We can now formulate the backward difference equation for the center point: Di:Ai Ti-1,j+1 +Bi T1,?r1 (2.2.19) 2.3. Analysis of spheres The analytical solution of equation 2.1. 2 for transient heat conduction in a sphere with no heat generation (Schneider, 1957) is -fifat . . r T-T1_ g 251nf3i-Bicosfii Sin i—R R2 T0 —T1 — 121 (Si - sin [3. B 1 e 1 1 R where Bi are the roots of the transcendental equation NBi : 1 - Bi cot Bi (2.3.1b) Equation 2.1. 2 can be presented in the following form (Ball, 1923) 00 ........_....—___ 2'2 je (2.3.2) 22 where j contains the location variable After sufficient time has elapsed, equation 2. 3. 2 converges rapidly and all terms but the first may be neglected. Taking the log of the first term form of equation 2. 3. 2 we obtain: 2 ‘31 ln(10) R2 at 10g(T-Tl):1ogj(TO-Tl)- (2.3.3) Equation 2. 3. 3 is the equation of the straight line plot of the temperature difference between the object and the fluid surrounding it versus time, where the slope of the line is {if a t/1n(lO) R2 and the intercept is j. Defining the unaccomplished temperature difference, T - Tl/TO — T to be 0.1 , the time required to accomplish this 1, temperature change is f, the temperature response parameter (Kopelman, 1966). Therefore 2 f _ 1n(10)R (2.3.4) Bf a For any material the functional solution can then be stated to be fa _ . .157: _ ”N131 (2.3.5) Pflug gt 31. (1965) plotted and tabulated values for jc’ jm and js (the value of j for the center of the mass average and surface, respectively) and f o/rZ were presented as a function of NBi . The solution presented above holds for all spheres, however, in systems where l/NBi > 6, a low NBi system, the internal 23 resistance to heat transfer is negligible and the object may be assumed to have a uniform temperature distribution. In low NBi systems we can make an energy balance over the sphere as follows: 8T Vp Cp .57:— Ah(T-Ta) (2.3.6) Integrating from T' to T" and t' to t" we obtain: TH -T fl = - -h—é—- (t” -t') (2.3.7) 1n _ 1 Vpcp Substituting in equation 2. 3. 7 for f and the unaccomplished temperature change, as defined above, we obtain: V p c R p c f :1n(10) A h =1n(10) 3 h (2.3.8) Multiplying and dividing equation 2. 3. 8 by RZ/k gives f l 10 __g_ : Tl”? (2.3.9) R Bi From equation 2. 3. 9 it is obvious that having the value of the temperature response parameter and the physical properties of the body we can readily evaluate the heat transfer coefficient of bodies whose N . is small (N B1 Kopelman ital. (1967) that low N Bi < O. 2). It has been shown by Bi systems can be used effectively to measure the film heat transfer coefficient using transient techniques. 24 2. 4. Generalized numerical solution for three one-dimensional heat transfer geometries The similarity between the equations for one dimensional transient heat transfer for the three geometries indicate that a general type of numerical solution may be written to cover all three cases. The equations for each of the three geometries are as follows: - Q + 21 - a _8__ QI- infinite slab (Z. 4.1a) p c at Br Sr P Q 8T _ 1_ §___ 11‘ . . . . - p c + at — a r [8r (r 8r) infinite cylinder (2. 4.1b) Q 8T _ _1_ L 2 8T ' 9 cp + 5%" ‘1 r2 8r (r 8r Sphere (2‘4'1C) where r is the characteristic dimension (half the thickness for a slab and the radius for the cylinder and sphere). The boundary conditions for equations 2. 4.1 are: tEO T=TO for O_<_ rER r=0 %§-= 0 for all t (2.4.2) r=R -k%rI:h(T-Tl) for t>0 At the center, r = O, we find, for the infinite cylinder and the sphere, that aT/Br = O which will give the indeterminate form for l__8_'_I_‘ and -Z— 2:— of 0/0. Using McLaurin's expansion as r 3r r 8r indicated above we find that, for r = O , equations 2. 4. lb and 2. 4.1c become Q. 2 _ 1192:“. 3T (2.4.3) p c at 2 p 81‘ for the infinite cylinder and Q 2 1 8T _ a T - C + W — 3 0. 2 (2. 4. 4) p 81' for the sphere. If we have symmetry about the center of the slab and radial symmetry in both the infinite cylinder and sphere, we can express the center node in the finite difference form as follows: at r z 0 Ti-l,j+l : Ti+1,j+1 for the slab + TiLj+l " T1,; _ 2 T1.-1,j+1 " T1,341 2 4 5 At a _ 2 ( ' ' ) (Ar) Q. ._i k for the infinite cylinder Qi Ti,j+1 " T1,; _ Ti- ,j+1 " Ti,j+l - — + -— Z X 2 (Z. 4o 5b) k At (1 (Ar)2 and for the sphere Q. T. . - T. . T. . - T. . _i 1, 1+1 1, 1 _ 1-1, 3+1 1, 3+1 (ArrZ 26 The general form of equations 2. 4. 5a, b, c is: T. . -T. . T. . -T. . 1,3+l 1,3 ___ X x 2 1-1.J+1 1:1“ (2.4.6) 01 - —- + k 0. At (Ar)2 where k = 1,2 and 3 for the slab, cylinder and sphere, respectively. At the boundary we can bring the convective contribution for the surface into equation 2. 4. la, b, c to get the finite difference forms: for the slab 9.. - 231 f; + Ti,j+1 'T1,j _ 5" i+1,3°+l 8" i,3+1 (2 4 7a) k c1 At _ Ar ° ' for the infinite cylinder P — 8T 8T 217 r —— - 211' r -— 3; + T1, 3+1 " Ti,3‘ _ 1 8r ,i+1.i+1 a i,j+1 - k 0. At _ 2w ri Ar (2. 4. 7b) and for the sphere 2 _ _.l + #3, 3+1 i,L : 1 i+l,j+1 i,3+1 (1 At 2 Ar 4w r 1 __ .1 (2. 4. 7c) ri is the variable heat transfer dimension. At the boundary we can define the following terms g3; _ T3+1,j+l ' 11,331 (2 4 8a) 8r - Ar ' ° i+l,j+1 27 ar g k __ _ (2. 4. 8b) i,j+1 (Ti,j+l ‘ Ti-l,j+l) By substituting equations 2. 4. 8a and b into equations 2. 4. 7, we obtain the finite difference equation for this node. For the internal nodes we can substitute the following terms into equations 2. 4. 7. (T. . -T.. ) .32 = A2. 1+1,1+15 “F'L (2.4.9a) 1' i+1,j+1 1 r (T.. -T. . ) .82 = Al. “1+1 4'14“ (2-4-9b) 8r . . 1 Ar 1,3+1 (T. . -T. .) 53—1; _ 141‘” 1,3 at _ A41 At (2.4.9c) where the terms Ali’ A21, and A4i are defined as follows: for the slab A1. = A2. = A4. 1 1 1 for the cylinder Ali = 21T riJrl A2. = 211' r. 1 1 A4. 1 (A1i + A21)/2 for the sphere 2 A2. = 441' r.2 1 1 A4. 1 (Ali + AZi)/Z 28 Program GEOMETRY which incorporates the solution of the three geometries is presented in Appendix 3. 3. EXPERIMENTAL DESIGN AND PROCEDURE 3.1a. Packed bed system considerations and desigg In the study of the packed bed approach for precooling fruits and vegetbables we could have used a model system. Using a model system has two distinct advantages; the individual particle can be constructed of stable materials of readily definable properties, of a simple shape for which an analytical solution may be available. By using the model system all of the particles used will be uniform which will avoid assumptions regarding shape and distribution otherwise needed. On the other hand, it was felt that the use of such a model would be oversimplification of the system mainly because the heat of respiration cannot be incorporated into the model. A fruit system is more difficult to work with since the system must be idealized by several assumptions. In working with apples we must assume; a simple geometry, sphere, with a given characteristic dimension and temperature independent thermal properties. Using the fruit system will enable us to evaluate the use of a model against the experimental data incorporating all of the differences caused by the deviation of the actual system from the numerical or finite difference model. Thus, after considering the advantages and disadvantages of a packed bed system made of fruits compared to a model system, we decided to use the fruit system. 29 3O 3.1b. Description of the equipment The experimental setup is shown in Figure 3.1, the apparatus was located in a 38 :1: 1 0F walk-in refrigerator. The apparatus was designed so that either 1, 2 or 3 boxes could be evaluated. A 96 :1: 2 0F walk-in warming room was located adjacent to the refrigerator. The height of the packed bed was determined by the number of crates stacked one above the other. Each crate (Figure 3. 2) was 15 inches high. The crates were raised into position by two hydraulic jacks attached to the entrance region (Figure 3.1). A centrifugal blower, downstream, was used to pull air through the bed with a damper before the blower, to regulate the air flow. The air was pulled through the bed because we wanted to eliminate gravity forces caused by convection currents and the energy contribution of the blower to the system. The crates were designed to, each, weigh approximately 100 pounds when full of fruit, a weight which could be readily carried by two men. The inside of the crate was lined with 7/8 inch expanded polystyrene insulation to validate our assumption of a radially adiabatic packed bed. A 1/2 inch layer of foam rubber was glued on the top of each crate so an air tight seal would be formed between crates when they were raised by the hydraulic jacks to joing the air tunnel. The bottom of the crate had three 1/2 inch rods spaced so as to give support to the wire mesh bottom (Figure 3. 2). 31 .839? Haunogiomxw 23 mo Edumdflu ofimgoflum .H.m oudwfim hoe/0H3 Hamzwfinufloo nomadp mocm> wfifiouswwmbm EL? Gowuomm HNMVll-DQ wouMHo mxomfl 3123?»: . uuommdm Baku 32 .unommdm macho mo 33> 33% .N .unommdm «5.an mo 33> .3 .m .ouduo «m3 .30 33> 33m .N .oudno 33 mo 33> .H .< .Gofloom amen. .N .m ondwfim " ‘4 (N \4 a as |_ L CI F‘l LJ —l .H E H 33 3 . 2. Metal rtrans duce rs Metal transducers were utilized to measure the convective film heat transfer coefficients in the packed bed system. The transducers were made from aluminum and copper, the physical properties of these metals are listed in Table 3. 2.1. Table 3. 2.1 Physical Properties of Metal Transducers ” 2 BTU 1b BTU ft Metal k (-————o—) p (7) c <—-;—) a (7;; hr ft F ft P 1b F Aluminum 70.0 169.0 0.214 1.936 Copper 222.0 559.0 0.0915 4.360 The transducers were 2. 75 :1: 1/64 inch diameter spheres with a 0. 073 inch hole drilled to the geometric center. A 30 gauge copper constantan thermocouple was placed in the hole and packed in place with filings of the same metal as the transducer by means of a hypodermic needle with a square cut end. To assure good contact between the thermocouple and the metal a drop of SAE 5 instrument oil was placed in the bottom of the hole before introducing the thermocouple. A 1/8 inch hole was drilled and tapped 3/16 inch from the thermocouple hole and a threaded rod was placed inside. This rod was used to support the transducer during the air tunnel tests and also to protect the thermocouple wire from mechanical failure. Figure 3. 3 is a photograph of the c0pper and aluminum transducers. 34 A B Figure 3. 3. Copper and aluminum spherical transducers. A. Copper sphere B. Aluminum sphere 35 Preliminary tests were made to obtain a relationship of NNu versus NRe for single spheres suspended in an air stream. The procedure for these tests was to warm the spheres to 96 :1: 1 c)F in an incubator and then to place them in a l x 1 ft. air tunnel equipped with straightening vanes to reduce free stream turbulence (Chessness, 1965) and to obtain the temperature history of these spheres at various air velocities. In the packed bed system only the aluminum transducers were used. The copper transducers were used only as a comparison to the aluminum in the study to obtain the NNu versus NRe relationship. The reason for not using the copper transducers was the 3. 3:1 ratio in weight between copper and aluminum. The heavy copper transducer would cause damage to the fruit with the possibility of the fruit yielding and conforming to the shape of the metal sphere causing additional errors. The aluminum transducers were spaced along the longitudinal axis of the packed bed at a fixed radial position (Figure 3. 4). Table 3. 2. 2 gives the longitudinal displacement of the transducers for the respective bed height. Table 3. 2. 2. Location of Center of Transducer in the Longitudinal Direction for the Three Heights of Packed Bed Height of Bed Transducer No. 15 inches 30 inches 45 inches 1 1.5in. 1.51n. 1.51n. 2 6. O in. 16. 5 in. 22. O in. 3 11. 0 in. 28. 5 in. 42. 5 in. 36 3.3 Measurement of pressure drop in theiacked bed The pressure drop through the bed was measured with static pressure taps mounted perpendicularly into the walls of the wind tunnel 5 inches upstream and 5 inches downstream of the packed bed test section. The upper tap contained a sliding probe which was actually a pitot tube measuring only the static pressure. The pressure taps were connected to an inclined manometer (i 0. 001 in.) with gauge liquid sp. gr. 0. 7970. The height of the bed for each case was the number of boxes in each test. The air velocity (for each of the tests) was varied by controlling the opening of the damper. Pressure drop measurements were made for several air velocities other than those of the heat transfer study. All pressure drop measurements were made after the bed had reached room temperature. 3. 4. Measurement of air velocity througli the packed bed The velocity of the air flowing through the packed bed was measured 8 ft. downstream from the test section using an Alnor thermoanenometer type 8500. The pipe was 12 inches in diameter at the point of measurement. This was the largest pipe diameter that could be used to obtain a velocity profile because the length of the anenometer probe was 6. 5 inches. In every test, velocity measurements were made every 0. 5 inch across the diameter of the 12 inch tube when the bed was at cold room temperature. Velocity measurements for the single sphere tests were made 1 ft. upstream from the sample, across the square 1 ft2 air tunnel. 37 3. 5. Materials for the experiment The product chosen for this experiment was Jonathan apples. This variety of apples was chosen because the fruits are relatively spherical in shape and they are stable enough to permit several heatings or in other words several tests for each quantity of fruit. The apples, grown in the Ionia district of Michigan, were stored under controlled atmosphere conditions in our laboratory until used in these tests. They were sorted for size and shape and placed randomly into the 18 x 18 x 15 inch crates. 3. 6. Temperature measurement Enameled copper constantan thermocouple wire glass wrapped with fiber glass overwrap was used throughout this study. Air temperatures were measured using 24 gauge and apple temperature using 30 gauge thermocouple wire. Temperatures were recorded at 5 or 15 minute intervals using a Bristol model 1P12G595-21-T7 potentiometer calibrated to give a 0-200 OF temperature output. The potentiometer was fitted with an analog to digital encoder which fed an input to a switching system and the data was recorded on eight channel punched paper tape and a printed tape. The punched paper tape was processed by the Michigan State University CDC 3600 computer and the data sorted by thermocouple number, time sequence and stored on a. magnetic tape in a form readily accessible for computation. Various programs were written for plotting and evaluating the data . 38 The apples for the temperature measurement were chosen for their spherical shape and weight(127 :t l g). The technique for placing the thermocouple in the apple was as follows: a hollow tube, the shaft of a 20 gauge hypodermic syringe, was inserted in the apple at the half way point between the stem and the calix. The 30 gauge wire was then inserted into the tube and the tube was pulled out from the apple leaving the thermocouple junction in the center of the apple. The 24 gauge thermocouple for measuring the air temperature was affixed to the 30 gauge wire near the entrance of this wire into the apple. In order to avoid probable temperature and drying errors, both lessions caused by the syringe were sealed with caulking compound. Apples with thermocouples were located in the packed bed as shown in Figure 3. 4 for each of the three packed bed heights. The textured spheres represent the transducers. There were two apples with thermocouples in each plane, in three planes there was an apple in the center of the crate one at the side of the crate and in one plane one of the apples was 4. 5 inches from the side and one at the center. 3.7. Procedure for a single run A detailed description of the procedure is as follows: The apples were loaded randomly into the crates except for the apples containing the thermocouples and the transducers, which were placed at the designed locations in the bed. The thermocouple wires were brought out of the crate, taped to the exterior wall to 39 20 45 42 6? s 27 ED 69 22 13 ® kg""“"“®‘ """r‘"‘F"““—“—-®"‘ GD Figure 3 . 4a. Apple (plain spheres) and transducer (textured spheres) location in the experimental three crate (45 inch) ‘ packed bed. 30 29 6? e 10 E91696} ® P———————————-———J————— M Figure 3. 4b. Apple (plain spheres) and transducer (textured spheres) location in the experimental two crate (30 inch) packed bed. 41 20 1‘19 9 i . 1 13 : a : 11 : l e e—@ : 15 6 : : | Q (B @ J A '20 .L—-10 —-—5.5 —55—-‘ ‘ 1mm mm 10 % 20 B Figure 3. 4c. A. Apple (plain spheres) and transducer (textured spheres) location in the experimental one crate (15 inch) packed bed. B. Top view for all crates of apple (plain spheres) and +1n'1v1nlq nnnnn l+n1r+111fintq nnnnnn \ 1nnn+1 r1“ 42 prevent loosening of the thermocouples during transport. The crate was then placed in the 96 i: 2 OF room; a fan was placed above the crate so the apples would reach uniform temperature faster. When the apples reached uniform temperature, 96 :1: 2 0F, six apples were weighed to :1: 0. 01 g; these apples were weighed before and after each test. The temperature recording system was turned on and the crate carried into the cold room. The crate was then placed on the mount and the jacks raised until the tunnel was sealed. The centrifugal fan was turned on and the thermocouple connectors coupled to the measuring system. The temperature history of the apples in the crate was measured until the temperature of the slowest cooling apple reached 43 OF, which is approximately 5 C)F above the cold room temperature. At the end of the test, the velocity and pressure drop were measured. For most cases the apples withstood 3 or 4 tests, however, after each test they were checked to eliminate the spoiled fruit. 4. BASIS OF ANALYSIS In this section we shall present the determination of some of the peripheral values necessary for the evaluation of the packed bed system. The following topics will be covered: 1) physical and thermal properties of product and cooling medium, 2) convective heat transfer film coefficients for both single spheres and the packed bed, 3) the equation for calculation of the heat of respiration, 4) pressure drop through the packed bed, 5) determination of the air flow rate through the packed bed, 6) determination of the total error of the numerical solution. Throughout this work many comparisons of experimental and calculated results are made on the basis of the equation for the least square lines with form Y = A‘Z + B' obtained from the data. The orthogonal coefficients for the least square lines were calculated using program UTEX-Z -LSC FWOP from the Michigan State University Computer Center Library. 4.1. Physical and thermal properties of product and cooling medium The temJLeratuIe In the numerical calculations the evaluation of temperature dependent properties was made for the temperature at the node at the time of the calculation. The temperature dependent properties for use in dimensionless numbers were evaluated at the average of the initial To, and final temperature, T1 . 43 44 Physicalproperties of the apples Total solids (by vacuum A.O.A.C., 1960) l6 % w/w Average density (by water displacement A. O. A. C. , 1960) . . 0. 82 gr/cc Specific heat of apples The specific heat for temperatures above freezing of most non fatty foods may be calculated from Siebels (1892) formula cp = 0. 008 x (percent moisture) + 0. 20 The value of the specific heat used throughout this study was 0.86 BTU/lb OF. Thermal conductivity of apple flesh The value of the thermal conductivity of an apple flesh used in this study was 0. 203 BTU/hr C)F ft; this value was determined by Kopelman (1966) for Jonathan apples grown in Michigan. Throughout this study we shall use this value and assume it to be temperature independent in the temperature range of our experiment, 38 - 100 OF. Physical properties of the cooling medium Density of air Air at atmospheric pressure approximates an ideal gas and its density can be expressed as follows (Zimmerman and Lavine, 1964): MP P : —— a R T 1545. 43 ft-lb/mole OR where R T=°R 45 Smcific heat of air The value of the specific heat of air from Zimmerman and Lavine (1964) for the temperature range of 0-100 OF is 0. 240 BTU/lb OF. Thermal conductivity of air The thermal conductivity of air, ka' was calculated as a function of temperature according to the equation (International Critical Tables, 1929) 3 2 k _ k (273+125)(Tc)/ a _ ao 'TC+125 273 T = OK C kao = thermal conductivity of air 0 0c = 0.0140 BTU/hr 0F ft Prgoerties of the packed bed Diameter of the wticles in the packed bed The average weight of the fruit in the bed was determined by dividing the weight of the fruit by the number of fruit in the crate. The average apple weighed 127 g and the radius was taken to be the radius of a sphere having the same volume. The average radius for the apples used throughout this study was 1. 32 inches. The porosity of the bed The porosity of the packed bed was evaluated as follows: volume of voids volume of bed 46 The value of the porosity used was obtained by filling a 1 ft3 box with apples and then pouring water in the box. The volume of water was then determined by weight difference. The value of the porosity used throughout this study was 6 = 0. 451 . 4. 2. Heat of respiration The evaluation of the cooling rates of apples must include the amount of heat liberated by the product as it continues to metabolize sugar 3. (’6HIZO6 + 6 O2 —‘ 6 H20 + 6 CO2 + Q The heat of respiration of apples is a function of the variety. Smock and Neubert (1950) report Q for apples in the range of 10 0c 2 to 3. 3. The data of heat of respiration for Jonathan apples is very limited and only three values are available (Rose 23%., 1963). The linear regression line of log of the heat of respiration versus temperature for Jonathan apples indicated the 010 0C to be 2. 46 for this variety. The following equation will be used to determine the heat of respiration: T. . -32 _}.z_.l.____ Qi,j = Q32 oF leO 0C 18 where Q10 0C : 2.46 ._ -4 . Q32 OF 2 2. 453 x 10 BTU/min lb Qi j : heat of respiration at node i at the 9 temperature Ti j, BTU/min lb 47 4. 3. Convective heat transfer c ocfficients The copper and aluminum transducers were used to determine the value of the convective heat transfer coefficient of a single sphere submerged in a moving fluid (section 3. 2) and also to determine the heat transfer coefficient in the packed bed. The temperature history of the single sphere was recorded every 15 seconds using a 12 point recording potentiometer. The data was plotted on semilog paper and the temperature response parameter required in the calculation of h (equation 2. 3. 8) was obtained. Most of the correlations available from the literature for single submerged objects are of the form: D! _ Y NNu . C NRe (4.3.1) where C' and D" are empirical constants. We selected equation 4. 3.2 (McAdams, 1954) for calculating h values which will be compared to the h values determined from our data. / . OC' I : I < NNu 0.37 NRe 17 \ NRe 70,000 (4.3.2) Table 4. 3.1 lists the experimental values of the film heat transfer coefficient obtained for each air velocity tested and also the values calculated from equation 4. 3. 2 for the corresponding NRe O 48 Table 4. 3.1. Heat Transfer Coefficients Obtained Experimentally and Calculated from Equation 4. 3. 2 for 2. 75 inch Copper and Aluminum Spheres h, BTUjhr OF ftZ Velocity Ex erimental Equation 4. 3. 2 ft/min Copper Aluminum (McAdams, 1954) 82 2. 42 2. 42 2. 60 135 2.90 2.92 3.73 260 4.10 4.20 5.05 430 5.20 5.20 6.78 The experimental values of the film coefficient are in very good agreement for the two metals. These are the average values from 2 and 5 tests for the copper and aluminum spheres, respectively. The values of h determined experimentally are in the range of the various correlations available. Since the experimental results were well replicated (less than 5% difference), we feel that the experimental values of the film coefficient obtained using the transducers give a good approximation of our system. Since we have a low NBi system, the convective heat transfer coefficients in the packed bed were evaluated by making an energy balance for the aluminum sphere (equation 2. 3. 6). The temperature of the transducers within the bed was measured every 5 minutes in each test. The average values obtained for the convective film coefficient for each air velocity using equation 2. 3. 6 were used in the numerical computation. The value of h for any one air velocity was observed to vary about 20%. The discrepancy in the 49 experimental value may be attributed to experimental error, but is primarily due to variations in the flow pattern around the transducers in the packed bed. The experimentally determined convective film coefficient was used to calculate the Chilton-Colburn (1934) factor, JH’ by the following equation: J = h (N )‘2/3 (4.3.4) The value of JH obtained when plotted versus the NRe , G pa/p'a , should correlate in the following form: 1 J : C'(N D (4.3.5) where C' and D' are empirical coefficients. Barker (1965) lists many correlations for JH versus NRe' The values of J obtained for a given NRe can vary as much as H 40 fold. Comparison of the experimental data with any one correlation from the literature was not found to be very meaningful because of this very wide range of possible values. We chose two correlations, those of Gamson 3.2.31.1. (1943) and Bird _e_t__a_1_. (1960), to compare with the experimental data and they are, respectively: _ 041 JH - 1.064 NRe (4.3.6) -0 41 : . o 07 .TH 0.61 NRe (4 3 ) We chose these two equations because from here it is obvious that for any given NRe the value of JH’ and consequently the value 50 of h, from these correlations will differ by a factor of 1. 75. Table 4. 3. 2 lists the experimental values of the heat transfer coefficient and also those of Gamson fifl‘ (1943). The values of h from equation 4. 3. 7 are exactly 0. 57 times those calculated from equation 4. 3. 6. From this table we can see that the experimental values lie between these two correlations. Table 4. 3. 2. Experimental Film Heat Transfer Coefficient, h, BTU/hr OF ft2 for the Air Velocities Tested o , Z J Velocity h, BTU/hr F ft H ft/min Experimental Gamson _e_t 31 Experimental Cams on ital 11943) @943) 37 2.9 (2.4-3.6) 2.9 0.0640 0.0645 173 4.9 (4.0-6. 0) 7.15 O. 0237 O. 0346 245 7.5 (6.0-8. 4) 8.8 0.0250 0. 0295 It can also be noted from Table 4. 3. 2 that the value of JH calculated using the experimental values of h for the 245 ft/min interstitial velocity is larger than the J for 173 ft/min. Since H theory predicts that JH will correlate with NRe in'the form of equation 4. 3. 5, we feel, that this value is larger than it must be in actuality. The spread of h values for this velocity is large (6. O to 8. 4 BTU/hr OF ftZ) and the value of h may actually lie in this range. We speculate that the value of h for the Z45 ft/min test is somewhat high and that a value of 6. O BTU/hr OF ft2 is probably 51 a better approximation. In the numerical calculation we shall try both of these values and observe the differences. 4. 4 Pressure drop calculations One of the approaches to evaluate the pressure drop through a packed bed is to assume that the packed bed is a tube of complicated cross section. The average velocity through this tube is given by (Bird 2311. , 1960): 2 AP Rh < > : ————-— U ZuL (4.4.1) where Rh is the hydraulic radius and is defined by volume of voids _ 6 volume of bed Rh _ a wetted surface (4' 4' 2) volume of bed and where 3”v a : ref- (4. 4. 3) av = total particle surface/volume of. particles 6 D : ;— for spheres (4. 4.4) V where D is the diameter of the sphere. The drag force of a fluid on the particles in the bed, F : ACAP can be related to the surface area of the particles and the kinetic energy of the fluid past the surface (Wentz and Thodos, 1963). 2 F = f A 3.1L (4.4.5) 1 S ch 52 constant which includes shear effects caused H where f1 by velocity gradients and surfaces 2 cross sectional area of packed bed, ft2 2 surface area of spheres in the bed, ft2 = volume of bed, ft3 AP : pressure drop through the bed lb/ft2 L = length of bed, ft A :aV 1: aA L (4-4-6) then AP : f aL—p-y-— (4.4.7) 1 ch Taking the logarithm of both sides we get: log AP : (3' log a; +D‘ (4.4.8) C where C' and D' are empirical constants. Substituting equations 4. 4. Z and 4. 4. 3 into 4. 4. 7 we get: :31: 6DL 11—6 I f1 ° “'4'” pU/Z gC Ergun (1952) has shown that for randomly packed spheres the value of the friction factor is given by the following modified form of equation 4. 4. 8. AP 10 U‘77-2 gC D c _ f 1-6 _ fk (4.4.10) Ergun (1952) determined the following relationship between U_I_>_2_ . , and the modified NRe’ ”(1%) . the friction factor, fk 53 fk = 1.75 + 1350 Re' (4.4.11) Replacing equation 4. 4.10 into equation 4. 4.11 we obtain the Ergun equation which is the sum of the Blake -Kozeny equation for laminar flow and the Burke-Plummer equation for turbulent flow. When the log of the pressure drop is plotted versus the log of the kinetic energy (equation 4. 4. 8), a straight line relationship is observed. The equation of this line is of the form U2 D' AP ___ CI (B—__) (4. 4.12) where C' and D' are empirical constants for the system. Figure 4.1 shows the three lines obtained for the three bed heights. The values of the constants C' and D' for the three lines are listed in Table 4. 4.1. Table 4. 4.1. Coefficients for Equation 4. 4.12 versus Bed Height Bed Height (inches) C' D' 15 0.65 0.8163 30 1.31 0.8133 45 2.10 0.8349 We can see from Table 4. 4.1 that the value of D', the exponent of the kinetic energy term, is very similar for the three packed beds. We would expect this since the velocity effects would be the same for 54 I l l I l I I I l l T l l I I l l '- -1 L : 45 inches 1.0 - 4 .. O -I. I ' I : . $L=30inchfs — . /_1 — o d p 9 L =15 inches NA — _ 11‘...” E 9-1 0 1 _ <1 . C .- I- —( l.- -1 O O. 01 1 1 1 1. 1 1 1 l‘ 1 1 1 1 1 1 1 1 1 0. 001 0. 0. 05 LL 2 gc Figure 4.1. Relationship between overall pressure dr0p and kinetic energy for the flow of air through the bed of fruit. 55 the three packed beds and we would expect the pressure drop to vary with the length of the bed. The values of C', change with the height of the bed and are approximately in the ratio 1:2:3, which is the corresponding height ratio. The relationship of the friction factor calculated using equation 4. 4.10 versus the N for the three bed heights is presented in Re' Figure 4. 2. The equation of the least square line for the experimental data is: -.447 rk _ 1.41 NR8, (4.4.13) The least square line obtained for the experimental data is shown in Figure 4. 2 as well as Equation 4. 4.11 of Ergun (1952). There is a large discrepancy between these two lines. The difference in the value of the friction factor as calculated from the two equations may in part be due to the loose packing of the fruit and to possible channeling through the bed causing significantly lower values for the friction factor. The experimental results suggest that the use of existing correlations for the friction factors may give values which are several times larger than the actual measured value and thus be the cause of over estimation of the pressure drop through the system. We cannot expect to have a correlation equation which will give a consistent value of the friction factor for produce systems because the bed is not shaken or consolidated. Consolidating the bed may cause damage to the produce, by bruising, consequently, in most 10 l-c Figure 4. 2. 56 [TIT TWFIIT JLlll Ergun (1952) [L141 1 _ UDE '— 11(1-6) Relationship between friction factor and modified N for the packed bed of fruit. Re- 57 systems the degree of consolidation or settling of the bed will vary from case to case depending on the method of filling the crates. 4. 5. Calculation of the flow rate of air throggh the bed The velocity distribution in the 12 inch cross section 8 feet downstream from the packed bed along the tunnel was measured as described in section 3. 4. The velocity profile of the duct was uniform and similar to that for turbulent flow. Since the velocity measurements were taken at even intervals, the flow rate through the bed was computed using the trapezoidal rule to sum up the fractional flow rates. The equation for the flow rate is: 12 2 2 Flow rate = 1r U (r - (r - Ar) ) n21 n n n The interstitial velocity or average velocity through the void spaces can then be obtained by dividing the flow rate by the product of the void fraction and the area of the bed. The air velocities used in this study were 245, 173 and 37 ft/min for our experimental packed bed. 4. 6. Evaluation of the total error of the numerical solution The amount by which the numerical solution differs from the exact solution of the partial differential equation is the total error. The total error is the sum of the discretization error and the stability error. The discretization error is the difference between the exact solution and the finite difference solution of the differential equation, this difference is due to neglecting terms in the Taylor expansion of the differential equation. This error does not lend itself easily 58 to computation since it is a function of the derivatives in the equation. The stability error is the difference between the finite difference solution of the differential equation and the actual numerical solution. This error is sometimes called the round off error and is caused by the rounding off of the values after each calculation. Numerical solutions are invariably more accurate than this error would indicate (Smith, 1965) since the stability analysis assumes all errors have the same sign, which is not necessarily true in practice. Smith (1965) states that for stable and convergent solutions the discretization error is the most influential of the two. We can decrease the discretization error by decreasing the size of Ar and At, however, a constant relationship between these two factors must be maintained. This method of improvement of the finite difference solution is limited by amount of time and the number of computations possible. In cases where the exact solution is available, the total error will be the difference between the exact and the numerical solution. We will, therefore, take the numerical solution for a sphere, apple (no heat generation), with a convective boundary condition and compare it to the exact solution (section 2. 3). We will take a range of values of Ar and At keeping 9 constant. Since we have chosen the temperature response parameter, f, and j to explain and correlate the cooling phenomenon, we will compare these values from the numerical solution with those of the exact solution. 59 Table 4. 6.1. Tabulation of the Values of f and j Obtained from the Numerical and the Exact Solution Vi Numerical Solution Exact Solution Parameter Ar = r/3 Ar = r/15 Ar = r/30 Equation 2.3.1 At = 0. 05 At = 0. 002 At = 0.0005 f 4.10 3.18 3.10 3.00 j 1.20 1.20 1.20 1.20 It becomes apparent that a judicious choice of values of Ar and At must be made, in order to decrease the magnitude of the total error, taking into account that after a given point the number of calculations to decrease this error below 5% becomes very large and the time required to calculate the temperature history of a packed bed becomes an important consideration. We will, in all the computations, use Ar : r/15 ft and At = O. 002 hr. It can be seen by comparing the value of f from the exact solution, 3. 0 hr. , against the value from the computed temperature data using the above mentioned values, 3.18 hr. , that the total error of the computation is approximately 6% for a single sphere. 5. RESULTS AND DISCUSSION In this section we shall give the results of the experimental study and those calculated by the finite difference model. It is our intent to show that the finite difference model simulates the experimental system, to point out how to use such a model for predicting the temperature of the particle in a packed bed fruit cooling system and to try to point out the limitations of the method. We feel that agreement between the calculated and experi- mental results will verify the basic assumptions made in the analysis. Differences between the experimental and calculated data are expected because of the nature of the problem, biological variability and number of parameters beyond our control. The experimental and the calculated values of fa and ja will be presented in the form of plots of apparent f and j versus the height of the sphere (apple) in the bed. The values of fa and ja determined from the experimental and calculated data are listed in Appendix 1. It is our opinion that comparisons between the calculated and experimental data for verifying the simulation can be made successfully by comparing the coefficients of the regression lines. Even though the calculated lines are higher degree than first order, but nearly straight lines, we decided to present the experimental and calculated data in the form of the first degree polynomial Y = A'Z + B' so that comparison can be made on the basis of the slope and the intercept. 60 61 Apparent tergperature response parameter The values of fa obtained from the experimental and calculated data are ploted in Figures 5.1, 5. 2 and 5. 3. From these graphs we can see that the slope of the line of the calculated values of fa versus height and that of the experimental values of fa for each air flow rate are in very good agreement, within 5%. It is interesting to note that for the three air flow rates the lines for the calculated data were below the lines of the experimental data, the value of the intercept of the line was smaller. We would have predicted that, on the basis of the total error of the numerical computation (section 4. 6), the opposite would occur, however, we are not greatly disturbed by this anomaly. We feel that the agreement in the slope of the lines of the calculated versus the experimental data, for any one air velocity, establishes the fact that the numerical model simulates the experimental system. The differences in the intercept of the lines of the calculated and experimental values of fa are due to experimental error and may, in part, be caused by the differences between the ideal system using the assumptions listed in section 2.1 and the experimental system. Table 5.1 lists the values of the intercept and the slope for the respective air velocities. The slope of the fa versus height lines appears to be a function of the air flow rate and the relationship seems to be of a logarithmic nature. We plotted the slope versus the air velocity on log -log paper (Figure 5. 4) for both the calculated and the experimentally determined values and obtained parallel straight lines for both cases. This limited amount of data does not give 62 conclusive evidence of the relationship, however, the possibility is exciting since it indicates that the relationship between the slope and the velocity may be expressed as: slope = 203.00 U-'O'85 where U is the interstitial velocity in ft/hr and slope = hr/ft. We can see from Table 5.1 that the spread of the fa versus height line varies with the air velocity. For the experimental data, the standard deviation is about 10 to 15% of the value of the intercept while for the calculated data it ranges 2-4% of the value of the intercept. In the experimental determination of the value of the heat transfer coefficient, we stated that the value h = 7. 5 BTU/hr OF ft‘2 was, apparently, too high for the 245 ft/min velocity. We speculated that a value of h = 6. 2 BTU/hr OF ftZ, still in the range of the experimental values, was more realistic. Both values of h were used to calculate the temperature history of the bed and the values of fa determined for each case are presented in Figure 5.1. The equations for these lines are: 2 h = 6.2 BTU/hr 0F ft ; f h = 7. 5 BTU/hr OF 11.2; f a 0.0287Z+l.20 a 0.0288Z+l.10 The two equations indicate that the slope of the two lines is the same and that a change in the value of h only affects the intercept of the line. The two lines having the same slope is further evidence that the slope is primarily a function of the velocity. We believe, that the particles in the packed bed behave as single spheres, how- 63 ever, in the bed their cooling rate is also a function of their position. If we calculate the value of f for a single sphere (apple), for h = 7. 5 and 6. 2 BTU/hr OF ftZ we obtain 0. 99 and 1. 08 hrs, respectively. Comparing the difference in f caused by the change in h, 0. 09 hrs, with the change in fa , 0.10 hrs, we can see that both systems appear to behave in a similar manner. Table 5.1. Slope, A', and Intercept, B', for the Lines of fa versus Height for the Experimental and Calculated Data for the Three Air Velocities Air Velocity A' B' :1: S. D. ft/min Experimental Calculated Experimental Calculated 37 0.1340 0.1385 2.88:1: 0.312 2.38:1: 0.090 173 0.0348 0.0359 1.56:1: 0.202 1.47:1: 0.038 245 0.0270 0.0287 a) 1.38:1: 0.084 1.20:1: 0.017 a) 0.0288 b) 1.10:1: 0.020 b) a) calculated using h : 6. 2 —_—%E_U—2— hr F ft b) calculated using h = 7. 5 ___EEOII_J__Z_ hr F ft A comparison of the fa for the apples, in the experimental study, located near the wall with that for the apples in the center of the crate, showed for the 173 and 245 ft/min cases no indication that the value of fa varied with radial position. However, in the 37 ft/min air velocity we did find some differences (Table 5. 2.). 64 H 0.0271 Z +1.38 0.0287 Z +1.21 0.0288 Z +1.10 A Experimental fa B Calculated h = 6.2 BTU/hr OF ft'Z fa C Calculated h 7.5 BTU/hr OF £12 1a II N H J 1 L 1 O 10 20 30 4O 50 Z (inches) Figure 5.1. Experimental (solid line) and calculated (broken line) apparent temperature response parameter versus height of the fruit in the bed for 245 ft/min air velocity. 65 0.0348 Z +1.56 0.0359 Z +1.47 A Experimental fa B Calculated h = 4. 9 BTU/hr OF {:2 fa l l L l 10 20 30 40 50 Z (inches) Figure 5. 2. Experimental (solid line) and calculated (broken line) apparent temperature response parameter versus height of the fruit in the bed for 173 ft/min air velocity. 66 Ha H A Experimental 3. 0.134 Z + 2. 89 B Calculated h = 2. 9 BTU/hr OE ftZ fa 0.138 Z + 2.38 A B 7 ”a? 5 H :5 O ,c: a: ‘4-4 3 l l J L 20 30 40 50 O 10 Z (inches) Experimental (solid line) and calculated (broken line) apparent temperature response parameter versus height of the fruit in the bed for 37 ft/min air velocity. Figure 5. 3. 6‘7 1'O__ l I I Tr1111 1 1 FIIIII. 1. .. 0.5- ... 1' '1 cu Q1 3 0.1- -1 U) " .1 0.05.- .. 0.01 1 1 1 1 11111 1 1 1 4 1111 10 50 100 500 air velocity (ft/min) Figure 5.4. Relationship between the slope of the line of f versus height and the air velocity for the calculated roken line) and experimental (solid line) systems. 68 Radial differences may be due to two major causes: radial heat flow through the sides of the crate and different local velocities. The effects of these two causes are self compensating because increasing of the cooling rate of the apple at the side of the crate due to heat flow to the side will have an opposing effect due to lower velocity at the side (wall effects). From Table 5. 2 where the fa for the apples at the various positions in the bed is tabulated, we can see the value at the sides to be slightly higher than the respective apples at the center. This suggests that at low flow rates, the effect of velocity differences between two apples in the same plane on their cooling rate is greater than the effect of heat flow to the side. Therefore, in high velocity, due to the high degree of turbulence, the effect of velocity differences will be smaller and the two effects will nullify each other; this was observed in the experimental study for the 173 and 245 ft/min cases. Table 5. 2. Values of fa from the Experimental Data for Test 217 Velocity : 37 ft/min Z (in) Radial. Position 1. 32 14. 0 27. 0 42. 0 Center of crate 2.7 4.3 6. 78 7.5 4. 5 in from side 6. 0 On side of crate 2. 76 5. 0 7. 62 69 The Lag Factor In a manner siznilar to the apparent temperature response parameter, the value of ja obtained from both the experimental and the calculated data were plotted versus the height of the particle in the bed. The least square lines for the three air velocities are shown in Figures 5. 5, 5. 6 and 5. 7. Table 5. 3 lists the values of the slope, intercept and the standard deviation for the lines of ja versus height for both the expe rimental and the calculated data. The slope of the least square lines for the 173 and 245 ft/min air velocity (Table 5. 3) is small, showing a very small change in ja with height in the bed. The 3 '7 ft/min case on the other hand, shows that there is a significant change in ja with height in the bed. There is agreement in the values of the slopes of the experimental and calculated lines, however, we cannot predict, as in the case of fa , a relationship between the lepe and the air velocity. We. can see, Table 5. .3, that the fluctuation of ja values about the experimental line is 10-15% while that for the calculated data may be as high as 9% which is suggestive that the possible error in the determination of ja from the calculated data is quite high. The lines of ja versus height in the bed for the 245 ft/min air velocity case where. h : 7. 5 and 6. 2 BTU/hr OF ftZ were used to calculate the temperature history of the bed are shown in Figure 5. 5. The equations of the lines are, respectively 70 .00868 Z +1.25 .00784 Z +1.53 .00669 Z +1.59 A Experimental ja B Calculated h z 6. 2 BTU/hr OF 112 ja c Calculated h 7. 5 BTU/hr OF {1:2 ja 3 )— a: 2 _ C "-) ”g. #" ..__——-—-:‘:—""———--—""' 4—A ——_——-—'—-' B 1 II— J l l l 10 20 30 40 50 Z (inches) Figure 5. 5. Experimental (solid line) and calculated (broken line) apparent lag factor versus height of the fruit in the bed for 245 ft/min air velocity. 71 A Experimental ja B Calculated h : 4. 9 BTU/hr OF ftz ja 0.00729 Z +1.19 0.00735 Z +1.53 H 1 1 l I 10 20 30 50 50 Z (inches) Figure 5. 6. Experimental (solid line) and calculated (broken line) apparent lag factor versus height of the fruit in the bed for 173 ft/min air velocity. 72 0. 0522 Z + 0. 98 0.0686 Z +1.06 A Experimental ja B Calculated h = 2. 9 BTU/hr OF ftz ja 1 1 1 1 10 20 3O 4O 50 Z (inches) Figure 5. 7. Experimental (solid line) and calculated (broken line) apparent lag factor versus height of the fruit in the bed for 37 ft/min air velocity. 73 ja 0. 00699 2 +1.59 ja 0.00814Z+1.53 The difference in the slope is unexpected and may be caused by the possible error in the determination of ja . Table 5. 3. Slope A', and Intercept B', for the Lines of the Values of the Experimental and Calculated Data of ja versus Height for the Three Velocities Air Velocity A' B' :1: S. D. ft/min Experimental Calculated Experimental Calculated 37 0.0522 0.0686 0.98i 0.190 1.06:1: 0.085 173 0.00729 0.00715 1.19i0.120 1.53:0.010 245 0.00868 0.00814 a) 1.25:1: 0.110 1.53:1: 0.030 a) 0. 00669 b) 1. 59 :1: 0. 098 b) 2 a) calculated using h : 6. 2 BTU/hr OF ft b) calculated using h = 7. 5 BTU/hr OF ft2 Heat generation If we consider the integrated values of the heat gene rated by each fruit at the separate positions in the bed during the time required to reach 53 OF (Figure 5. 8) we can observed that for both the 245 and 173 ft/min air velocities the heat generation term along the column is rather small, the 510pe of the line of BTU produced versus height in the bed is 0. 00142 and 0. 00211 respectively. The 37 ft/min test, on the other hand, shows a much larger slope, 0. 0106, and a 74 A h = 6.2 BTU/hr OF ftZ Q = 0.00142 Z + 0.025139 B h = 4.9 BTU/hr OF ftz Q = 0.00211 Z + 0.02777 C h = 2.9 BTU/hr OF ftZ Q = 0.0106 Z + 1.03176 ;/ / 0.3 )- / // :5 / [-4 0.2 r E // O // B ’ I I 0.1 - ’— I I ,q—I" ’ ‘I I ’ I / I I I, ’— I, "' J 1 1 l 0 10 20 30 40 Z (inches) Figure 5. 8. Relationship between the integrated heat of respiration and height of fruit in the bed. 75 significant quantity of heat produced by each fruit as the height in the column is increased. When we consider commercial fruit cooling system where the height of the fruit stack is large, and air velocities usually very low (< 10 ft/min), we can see that as the height of the stack increases the contribution of heat by respiration for each fruit becomes quite significant. The overall effect of this contribution on the cooling rate of fruit is difficult to evaluate since it is not the heat of respiration which will affect the cooling rate of the fruit, but the net rate of heat production which is the heat of respiration minus heat of evaporation. We observed experimentally that in our lowest air flow rate we had an average weight loss of approximately 0.1 gr. from an apple at Z : 42 in during cooling of the fruit. Considering the weight loss to be moisture this would be a loss of 0. 2 BTU. The energy produced by the fruit at this height was 0. 52 BTU (Figure 5. 8). It appears that the likelihood of a cancellation of the effect of the heat generation contribution would be greater at the higher air flow rates, however, as the cooling time increases the heat generation would become more significant. The magnitude of the net rate of heat production would, of course, be dependent on the actual value of the rate of heat generation. f andj a a f It would be quite helpful to have a single value, even if it is of an arbitrary nature, obtained from experimental temperature time data, which would represent the cooling rate of the particle. 76 This is why we plot the logarithm of the temperature difference between the center of the fruit and the constant inlet air temperature (rather than the variable local temperature) versus time. This method will provide a straight line curve, the slope of which we will designate fa' the apparent temperature response parameter. We make the distinction between the temperature response parameter, f, and the apparent temperature response parameter, fa , since the latter would be a function of longitudinal position. This method of presenting the data will incorporate all of the parameters affecting the temperature time relationship (k, p Cp' the effect of geometry, the thermal properties of the system, h), and for the packed bed fa will vary with the height of the particle in the bed since the temperature of the cooling medium is a function of the height of the particle in the bed. Throughout this study all references to 3' will be for the j at the center of the particle. The value of ja is similar to j and is calculated using the constant inlet temperature as the cooling medium temperature. Experimental determination of fa and ja is not of equal accuracy. We can, by analogy, using the properties of f and j speak of the reliability of these values. The value of fa is more reliable than ja since ja is subject to more experimental errors. For example: ja is quite sensitive to the location of the temperature sensing element and to the initial temperature distribution of the 77 particle, whereas fa , while not location independent due to the heat generation term is for all practical purposes, in our experimental system, independent of both of these parameters. Conduction errors along the thermocouple wires will have large effects on ja’ however, will influence fa only slightly. There is a great deal of subjectivity involved in the determination of the value of ja since a slight difference in drawing the straight line asymptote will cause a large difference in ja’ this is true for both the ja determined experimentally and that obtained from the data of the finite difference model. However, fa will be affected only slightly. The fa and ja data can be used, for a given system, to evaluate the temperature of the particles within the bed, at a given time. This type of data is useful since it allows the processor or designer to predict the behavior of the system. In a manner similar to that for determining the temperature of a single object (equation 2. 3. 3) we can establish the temperature of a particle in the packed bed by substituting the corresponding value of fa and ja for f and j. It must be kept in mind, however, that we can only calculate the temperature of the particle when the temperature lies on the straight line asymptote. If we attempt to calculate the temperature of a particle deep within the bed after a short time from the initiation of the test we will find that since the temperature of the particle is in the lag phase we will get an erroneous value, which may even be higher than the initial tempe ratur e . 78 In Table 5. 4 we present the temperature of a sphere (apple) at three locations in the bed for the three air velocities. Temperatures from the experimental and calculated data are presented along with those calculated using the following form of equation 2. 3. 3 log (Tl - T) : -t/fa + log 33 (T1 - To) The experimental value of the temperature is the average of several readings and the fa and ja values used for each case were obtained from the equation of the least square lines for the respective value. We can see from Table 5. 4, that the use of the values of fa and ja from the least square lines in the modified form of equation 2. 3. 3, that there is good agreement for both the calculated and experimental data. It is apparent, however, that there is a greater discrepancy between the average and the actual temperature for the experimental system than for the finite difference model. We can attribute this difference to the small values of ja obtained from the experimental data. In the case of the 37 ft/min for the finite difference model it appears that the slope of the line of ja versus height is too large since the calculated temperature differs from that obtained from the output of program PACKBED 3 by as much as 9%. 79 Table 5. 4. Temperature of a Sphere Within the Packed Bed from the Output of Program PACKBED 3 and the Experimental Results; Values Calculated Using fa and 5a are also Shown CALCULATED EXPERIMENTAL calculated from calculated Time Velocity Z from from Z temp. from (hr) ft/min (in) PACKBED 3 fa ja data (in) data fa ja data 4.0 37 14.52 54.8 53.9 14.0 56.0 53.9 4.0 37 19.80 61.7 62.0 20.0 60.8 57.1 4.0 37 27.72 70.6 77.0 27.0 73.3 74.0 2.0 173 14.52 49.0 49.1 14.0 49.7 47.3 2.0 173 19.80 51.7 51.1 20.0 53.0 49.6 2.0 173 27.72 56.0 56.0 27.0 58.3 50.5 1.0 245 14.52 a) 62.1 a) 61.0 14.0 63.0 56.0 1.0 245 19.80 a) 65.4a) 5.9 20.0 68.4 64.4 1.0 245 27.72 a) 70.0a) 71.0 27.0 72.2 68.2 1.0 245 14.52 b) 60.7 b) 60.7 1.0 245 19.80 b) 64.2 b) 64.3 1.0 245 27.72 b) 69.2 b) 69.5 a) calculated data h b) calculated data h H 11 6.2 BTU/hr OF ft‘2 7.5 BTU/hr OF ft2 80 The reliability of the numerical model There is no theoretical solution available for evaluating the longitudinal temperature profile of a packed bed of solids of low thermal conductivity with heat generation, therefore, we do not have a comparison for our finite difference model. A solution for a packed bed of materials of high thermal conductivity, low NBi , was presented by Schumann (1929). Schumann presented the results of this solution as a plot of the dimensionless temperature versus time with the height in the bed as a parameter. Furnas (193 0), expanded the range of values of Schumann's graphs. Since we did not have a way of evaluating our experimental model, we decided to compare the numerical solution for a packed bed of aluminum spheres whose diameter was 2. 64 in and whose properties were those given in Table 3. 2.1, a low N system, Bi (NBi < 0. 01) with the solution of Schumann (1929) for such a packed bed. In Table 5. 5 we list some of the values for the unaccomplished temperature To - T/TO - TCD , versus the dimensionless time and height in the bed for both the numerical system and the theoretical curves. It may be seen that the values from the numerical model approximate those of the theoretical case with a difference of less than 10%. It is also encouraging to see that for most of the cases the value of the unaccomplished temperature obtained from the numerical solution was smaller than that obtained from Schumann's graphs. We would have predicted, in our system, that this would occur on the basis of the total error (section 4. 6) of the numerical 81 solution, since the total error would cause the temperature of the sphere to be slightly larger than theory would predict at any given time. On the basis of this agreement, we feel that the finite difference model estimates the conditions within the packed bed system and presents a method for calculating the temperature history of the spheres within the packed bed. 82 Table 5. 5. Comparison of Value of TO-T/TO - T Calculated Using the Finite Difference Model and those Obtained from the Charts of Schumann (1929) TO-T TO - T % Y 4) Height Time T -T T - T difference (1n) (hrs) 0 °°calc. 0 0° Schumann 1.90 1.02 9.24 0.45 0.1952 0.20 2.40 1.90 2.04 9.24 0.90 0.3929 0.40 1.78 1.90 3.06 9.24 1.35 0.5623 0.59 4.70 1.90 4.08 9.24 1.80 0.6949 0.72 3.98 1.90 5.10 9.24 2.25 0.7928 0.82 3.32 1.90 6.12 9.24 2.70 0.8622 0.90 4.20 2.96 1.02 14.52 0.45 0.0979 0.09 8.50 2.96 2.04 14.52 0.90 0.2352 0.25 5.80 2.96 3.06 14. 52 1.35 0.3825 0.41 6.70 2.96 4.08 14.52 1.80 0.5207 0.57 9.80 2.96 5.10 14.52 2.25 0.6397 0.69 7.20 Where Y = —A}3— g- , dimensionless V p c U P and <19 ill-— t, dimensionless V p c P T0 - T and -——-—— , dimensionless; T, particle To - Tco temperature at time t; To' initial particle temperature, and T the inlet air temperature all in OF . m ’ 6. C ONC LUSIONS l. The finite difference model presented in this study can be used for predicting the cooling of fruit in a packed bed arrange- ment. 2. The apples or fruit in the stack cool indpendently of one another, however, their cooling rates are a function of their position in the stack. 3. The slope of the fa versus height line is a function of the interstitial velocity and their relationship appears to be: slope = 203. 00 U'O' 85. 4. The heat generation contribution in high air velocity, high h, is negligible, however, depending on the product, may contribute a substantial quantity of energy to the system in low air velocity systems. 5. fa and ja may be used to express the cooling curve of a particle in a packed bed and to calculate the temperature of this particle during the cooling cycle. 6. Pressure data indicated that the pressure drop across the bed of fruit is a function of the height of the fruit, as predicted by theory. However, the friction factor correlations for packed beds give values of the pressure drop which is much too large. The lack of consolidation of the stack of fruit is probably the cause. 83 REFER ENC ES A.O.A.C. 1960. Official Methods of Analysis. 9th Ed. Assoc. Offic. Agr. Chemists, Washington, D.C. Amundson, N. R. 1956. Fixed beds with large particles. Ind. Eng. Chem. 48, 35. Bakker-Arkema, F. N. and W. G. Bickert. 1966. Deep bed sugarbeet cooling. Paper No. 66-350 Annual Meeting ASAE. University of Massachusetts. Ball, C. O. 1923. Thermal process time for canned food. National Research Council Bull. 37. Barker, J. J. 1965. Heat transfer in packed beds. 13E. Eng. Chem. 57, 443. Beck, J. V., 1966. Personal communication. Bennett, A. H., J. Soule and G. E. Yost. 1966. Temperature response of florida citrus to forced air precooling. ASHRAE Journal 14, 1966. Bird, R. B., W. E. Stewart and E. N. Lightfoot. 1960. Transport Phenomena. John Wiley Co., New York. Bradshaw, R. D. and J. E. Myers, 1963. Heat and mass transfer in fixed and fluidized beds of large particles. A. I. Ch. E. Journal 9, 590. Brownell, L. E. and D. L. Katz, 1947. Flow of fluids through porous media. Chem. Eng. Prog. Trans. 43, 10. Coberly, C. A. and W. R. Marshall, 1951. Tenlperature gradients in granular beds. Chem. Eng. Prog. 47, 141. Chesness, J. L. 1965. A Study of a Leaf Model used in the Sprinkling Method of Freeze Protection of Plants, Ph. D. Thesis. Michigan State University, East Lansing, Michigan. Chilton, T. H. and A. P. Colburn, 1934. Mass transfer (absorption) coefficients. Ind. Eng. Chem. 26, 1183. Ergun, S., 1952. Fluid flow through packed columns. Chem. Eng. Prog. 48, 89. Furnas, C. C. 1930. Heat transfer from a gas stream to a bed of broken solids, I. Ind. Eng. Chem. 22, 26. 84 85 Furnas, C. C. 1930. Heat transfer from a gas stream to a bed of broken solids, 11. Ibid, 721. Gamson, B. W., G. Thodos and O. A. Hougen, 1943. Heat mass and momentum transfer in the flow of gasses through granular solids. A. I. Ch. E. Journal 39, l. Guillou, R. 1960. Coolers for fruits and vegetables. California Agr. Exp. Sta. Bull. 773. Guillou, R. 1963. Pressure cooling for fruits and vegetables. ASHRAE Journal 5, (11), 45. Gurevitz, D. 1966. Engineering analysis of the controlled atmosphere storage. Ph. D. Thesis, Michigan State University, East Lansing, Michigan, 1966. International Critical Tables, 1929. McGraw Hill Book Co. , Inc. 5:214. Hill, F. B. and R. H. Wilhelm, 1959. Radiative and convective heat transfer in a quiescent gas-solid bed of particles. A.I. Ch. E. Journal 5, (4), 486. Hicks, E. W. 1955. Precooling of fruits and vegetables: some theoretical consideration. Proceedings Ninth International Congress of Refrigeration, Paris, France. Hood, C. E. 1967. A forced air precooler for pickling cucumbers. Food Technol. 21 (2), 190. Klinkenbeerg, A. 1948. Numerical evaluation of equations describing transient heat and mass transfer in packed solids. Ind. Eng. Chem. 40, 1992. Klinkenbeerg, A. 1954. Heat transfer in cross flow heat exchangers and packed beds. Ind. Egg. Chem. 46, 2286. Kopelman, I. J., J. L. Blaisdell, and I. J. Pflug. 1966. Influence of fruit size and coolant velocity on the cooling of jonathan apples in water and air. ASHRAE Trans. 72, 1. Kopelman, I. J. and I. J. Pflug, 1967. Characteristics of transient heat conduction systems in the high and low biot number region. ASHRAE Trans. 73, 1. Kopelman, I. J., C. Borrero, and I. J. Pflug, 1967. Evaluation of surface film heat transfer coefficients using transient method. Paper to be presented at the Institute International du Froid Meeting, Madrid, Spain. 86 Kunii, D. and J. M. Smith, 1961. Thermal conductivity of porous rocks filled with stagnant fluid. Soc. Petroleum Eng. Journal, March 1961. Leva, M. 1947. Pressure drop through packed tubes: part I a general correlation. Chem. Engr. Prog. 43 (10), 549. Leva, M. 1947. Heat transfer to gases through packed tubes. Ind. Eng. Chem. 39, 857. Leva, M. and M. Grummer, 1948. Heat transfer to gases through packed tubes. Ind. Eng. Chem. 39, 857. McAdams, W. H., 1954. Heat Transmission, 3rd Ed. McGraw Hill Book Co. , New York. McConnachie, J. T. L. and G. Thodos, 1963. Transfer processes in the flow of gases through packed and distended beds of spheres. A. I. Ch. E. Journal 9, 1. Pflug, I. J. and J. L. Blaisdell, 1963. Methods of analysis of precooling data. ASHRAE Journal 5 (11), 33. Pflug, I. J., Blaisdell, J. L., and I. J. Kopelman, 1965. Developing temperature-time curves for objects that can be approximated by a sphere, infinite plate and infinite cylinder. ASHRAE Trans. 71, 1. Rose, D. H., Wright, R. C., T. M. Whiteman, 1963. The commercial storage of fruits, vegetables and florist and nursery stocks. U. S. Dept. Agric. Handbook 66. Raney, J. L. 1962. E-2 UTEX LSCFWOP. Michigan State University Computer Center Library, East Lansing, Michigan. Rostos, G. M. 1960. Survey of fruit cool stores. Bulletin 282 Commonwealth Scientific and Industrial Research Organization, Australia. Sainsbury, G. F. 1951. Improved fruit cooling methods. Refrig. Eng. 59, 5. Sainsbury, G. F. 1955. Commercial experience with high speed precooling of soft fruits. Proceedings of Ninth International Congress of Refrigeration, Paris, France. Schneider, P. J. 1957. Conduction Heat Transfer. Addison-Wesley Publishing Co., Inc. 2nd Ed. 87 Siebel, E. 1892. Specific heats of various products. Ice and Refrige_ration 2, 256. Smith, G. D. 1965. Numerical Solution of Partial Differential Equations. Oxford University Press, London. Smock, R. H. and A. M. Neubert, 1950. Apples and Apple Products. Interscience Publishers, New York. Schumann, T. E. W., 1929. A liquid flowing through a porous prism. l. Franklin Inst. 208, 404. Soule, J., G. E. Jost, and A. H. Bennett, 1964. Certain heat characteristics of organges, grapefruit and tangelos during forced air precooling. ASAE Paper No. 64-311. Thevenot, R., 1955. Precooling. Proceedings of Ninth International Congress of Refrigeration, Paris, France. Wentz, C. A. and G. Thodos, 1963. Pressure drOps in the flow of gases through packed and distended beds of spherical particles. A.I. Ch. E. Journal 9, 1. Wilhelm, R. H., W. C. Johnson, and F. B. Acton, 1943. Conduction convection and heat release in catalytic converters. Ind. Eng. Chem. 35, 562. Yagi, S., D. Kunii and N. Wakao, 1961. Radially effective thermal conductivities in packed beds. Intl. Heat Transfer Conf. ASME, New York. Zimmerman, O. T., and I. Lavine, 1964. Psychometric Tables and Charts. Industrial Research Service, Inc. Dover, New Hampshire. APPENDD€1 Table A1.1 Experimental Pres sure Data Height Velocity KEZ NRe' AP fk (in) (ft/hr) 923;) ( U D 9 ) (lb/fiz) gC H(1-€) 15 2268 0.0005 1751 0.0166 0.4991 15 3492 0.0012 2965 0.0249 0.3158 15 4170 0.0017 3219 0.0415 0.3691 15 5160 0.0025 3893 0.0498 0.2893 15 7740 0.0057 5974 0.0995 0.2571 15 9450 0.0086 7363 0.1369 0.2327 15 12500 0.0149 9679 0.2405 0.2367 15 15000 0.0214 11577 0.2779 0.1911 30 2268 0.0005 1751 0.0166 0.4495 30 7260 0.0050 5604 0.1617 0.2374 30 9540 0.0086 7363 0.2696 0.2292 30 11820 0.0133 9123 0.3733 0.2067 30 13020 0.0161 10049 0.5143 0.2347 30 13620 0.0176 10512 0.2571 0.2145 30 15000 0.0214 11578 0.5309 0.1826 45 2184 0.0005 1686 0.0373 0.4037 45 7260 0.0050 5604 0.2903 0.2841 45 7740 0.0057 5974 0.3028 0.2607 45 10020 0.0095 7734 0.4272 0.2195 45 12000 0.0137 9262 0.6345 0.2273 45 12720 0.0154 9818 0.6968 0.2221 45 13620 0.0176 10512 0.7465 0.2076 45 15000 0.0214 11578 0.9124 0.2092 88 89 Table A1. 2 fa and ja Data from Experimental Model Air Velocity = 245 ftlmin Air Velocity = 173 ft/min Air Velocith 37 ft min Z fa ja Z fa ja Z fad ja 1.32 1.51 1.35 1.32 1.76 1.10 1.32 2.50 1.30 1.32 1.51 1.35 1.32 1.70 1.00 1.32 2.60 1.30 1.32 1.49 1.25 1.32 1.63 1.16 1.32 2.90 1.25 1.32 1.49 1.00 1.32 1.66 1.16 1.32 3.00 1.00 1.32 1.50 1.34 1.32 1.63 1.00 1.32 2.70 1.25 1.32 1.50 1.34 1.32 1.47 1.25 1.32 3.10 1.25 1.32 1.27 1.25 1.32 1.46 1.00 1.32 3.00 1.00 1.32 1.27 1.25 1.32 1.36 1.30 1.32 2.76 1.05 1.32 1.32 1.30 1.32 1.33 1.16 1.32 2.70 1.13 1.32 1.16 1.25 6.00 1.60 1.36 6.00 4.00 1.27 1.32 1.35 1.44 6.00 1.56 1.44 6.00 3.70 1.27 1.32 1.27 1.30 6.00 1.90 1.34 6.00 3.60 1.15 6.00 1.51 1.35 6.00 1.80 1.34 6.00 4.40 1.25 6.00 1.52 1.35 6.00 1.73 1.37 10.00 4.30 1.70 6.00 1.42 1.25 10.00 2.00 1.25 10.00 4.30 1.43 6.00 1.34 1.45 10.00 1.90 1.15 10.00 4.80 1.25 6.00 1.54 1.52 13.00 2.15 1.34 10.00 4.60 1.35 10.00 1.75 1.35 13.00 2.15 1.40 13.00 5.20 1.44 10.00 1.75 1.25 13.00 2.10 1.41 13.00 5.20 1.30 13.00 1.75 1.40 13.00 2.00 1.40 13.00 5.00 1.45 13.00 1.75 1.40 14.00 2.19 1.25 13.00 5.00 1.40 13.00 1.82 1.30 14.00 2.26 1.25 14.00 4.96 1.70 13.00 1.82 1.30 14.00 2.29 1.16 14.00 4.80 1.34 14.00 1.89 1.25 14.00 2.26 1.34 14.00 4.30 1.25 14.00 1.88 1.43 14.00 1.73 1.25 14.00 5.00 1.61 14.00 1.72 1.33 14.00 1.89 1.10 20.00 5.18 1.84 90 Table A1. 2 continued Air Velocity = 245 ft/min Air Velocity = 173 ft/min Air Velocity = 37 ft/min fl z fa ja z fa ja z fa ja 14700 1.80 1.33 20.00 2.58 1.15 20.00 5.04 1.61 20.00 2.00 1.34 20.00 2.46 1.35 27.00 6.73 2.70 20.00 2.00 1.25 27.00 2.69 1.42 27.00 6.48 2.60 27.00 2.24 1.43 27.00 2.59 1.34 27.00 6.78 1.88 27.00 2.14 1.63 27.00 2.95 1.44 27.00 6.00 1.84 27.00 2.14 1.35 27.00 2.82 1.43 29.00 6.30 1.84 27.00 2.34 1.35 27.00 2.09 1.52 29.00 6.60 1.65 27.00 2.04 1.65 27.00 2.02 1.34 42.00 7.36 3.56 27.00 1.89 1.70 29.00 2.77 1.25 42.00 7.50 2.50 29.00 2.44 1.34 29.00 2.86 1.50 42.00 7.62 2.70 29.00 2.60 1.25 42.00 2.95 1.61 42.00 2.44 1.79 42.00 2.86 1.52 42.00 2.60 1.61 42.00 3.49 1.43 42.00 2.42 1.61 42.00 3.32 1.50 42.00 2.42 1.61 42.00 2.42 1.79 42.00 2.37 1.70 91 Table Al.3 fa and ja Data from Finite Difference Model Air Velocity z 745 f't/min Air Velocity = 173 ft/min Air Velocity = 37 it/mi- A B c D Z fa Ja fa Ja fa J8. fa Ja 1.32 1.20 1.59 1.10 1.61 1-40 1.56 2.30 1.29 3.96 1.30 1.60 1.19 1.61 1.52 1.55 2.82 1.29 6.60 1.40 1.59 1.28 1.63 1.75 1.55 3.40 1.43 9.24 1.49 1.61 1.38 1.63 1.80 1.56 3.90 1.61 11.90 1.50 1.61 1.46 1.61 1.92 1.56 4.25 1.79 14.50 1.64 1.61 1.54 1.63 2.05 1.59 4.50 2.06 17.20 1.71 1.64 1.62 1.74 2.11 1.72 4.75 2.32 19.80 1.79 1.64 1.70 1.76 2.21 1.79 5.00 2.52 22.40 1.88 1.64 1.76 1.81 2.31 1.80 5.40 2.65 25.10 1.94 1.72 1. 2 1.81 2.40 1.82 5.80 2.68 27.70 2.00 1.79 1.88 1.81 2.50 1.83 6.20 2.75 30.40 2.05 1.79 1.94 1.82 2.60 1.82 6.52 2.94 33.00 2.14 1.84 2.00 1.83 2.66 1.84 6.90 3.15 35.70 2.20 1.84 2.06 1.83 2.71 1.83 7.25 3.50 38.30 2.30 1.84 2.11 1.85 2.80 1.84 7.65 3.70 40.90 2.40 1.84 2.17 1.85 2.91 1.86 8.00 3.85 43.70 2.46 1.89 2.22 1.86 3.00 1.86 8.40 4.00 46.30 2.53 1.93 2.35 1.86 3.10 1.86 8.65 4.15 o 2 A h = 6.2 BTU/hr 1“ ft 13 h = 7. 5 BTU/hr 017 it2 c h = 4.9 BTU/hr OF f’tZ , D h = 2.9 BTU/hr OF it“ 0(W()O APPENDIX 2 PROGRAM PACKDEDJ DIMENSION A(35)Ib<35)9C(35)oA1(35)9A2(35)9A3(35)9A4(35)9D(35)9 10(35)0T(35)9CC(3b)QDDC3519TA(750019TAB(7SOO) COMMON AoboCquTQQRESRQNNQDELTAU OIZ COMMON COEHTQRADQDELTARQRHOQSRHTQX .ZZ QTAOTAHQKKK THIS PROGRAM CALCULATES THE TEMPERATURE HISTORY OF A PACKED DED OF SPHERES ISSIE=NUMDER OF LONGITUUINAL STAGESONN=NUMDER OF NUUquNN*Z=THE NUMBER OF DIVIEIONb FOR THE prEREOUELTAuzTIME INTERVAL IN HOUR: ISSIE=18 NN=17 DELTAU=OOOE BEGINNING OF PROGRAM CALCULATIONS DO 7 I2=IOISSIE ORESP=OO KKK=O INITIAL TEMPERATURE DISTRIDUTION DO 10 N=10NN IO T(N)=960 RRINT 300 OIZ 300 FORMATIIH10*CALCULATION OF PACKED DEU FOR STAGE NO: *012) ZZ=IZ~1 TIME = 00 TD$=TIME DISPLACEMENT FOR SPHERESQZZ=NUMUER OF STAGEb ADOVE THE FIRST SPHERE TIME=TIME+TD5*ZZ DO 7 KK=1025 DO 9 K=195O KKK=KKK+I CALL SPHERE CALL TRIDI TIME =TIME+DcLTAU TABIKKK1=T(I) 9 CONTINUE PRINT IOOOTIME9T(I)0(T(I)9I=2917Ib) CORESP 100 FORMAT(1H 9F100696(E120503X)) 7 CONTINUE END SUBROUTINE SPHERE DIMENSION A(35)od(35)9C(35)1A1(35)OAZ(35)QA3(35)OA4(35)9D(35)9 10(35)OT(35)9CC(35)IDD(35)9TA(7SOO)9TAd(7bOO) COMMON AvbvaDoTQORESPQNNQUELTAU OIZ COMMON COEHTORADQDELTARqRHOQSRHT9X 911 QTAOTAUQKKK A d C AND D ARE THE CONbTANTS FOR THE UACKWARU DIFFERENCE EUUATIUN FOR EACH NODE 92 12 13 18 IF(IZOEO¢I)1COIJ TA(KKK)=40. X320 GO TO 18 TA(KKK)=TAU(KKK) x=10 CONTINUE RHO=DEN$ITY OF SPHERE LDS/CU FT 9 THERMC=THERMAL CONDUCTIVITY BTU/FT HR F o SPHT= SPECIFIC HEAT OF SPHERE UTU/Lb F0 RAD: CHAPACTEPISTIC DIMENSION IN FT RHO=5102 THERMC=o203 SPHT=086 PAD=1.32/lco DELTAP=PAD/(NN-C) ALPHA=THERMC/(RHO*SPHT) CP=SPECIFIC HEAT OF AIR bTU/LQ F QRHUA3UCNSITY OF AIR Ld/CU FT. SPVV SPECIFIC VOLUME OF AIR AT TEMPERATURE T CU FT/Ldo VEL=VELOCITY IN FPMO COEHTZCUNVECTIVE FILM COEFFICIENT IN oTU/HR SQ FT F. CP=024 COEHT=602 VEL=245 $QIO=Zo46 BIOT=COEHT*RAD/THEPMC THETA=DELTAU*ALPHA/(DELTAR**2.) SPVOL=1545o4J*(4éOo+T(1))/(25o967*14o7*144o) RHOA=lo/SPVOL PI=301416 K=l AI(1)=4*PI*(RAD**2) VOL=(180*180*2054/(I7280*b40))“(40/Jo*pI*(RAD**3)) VEL=VEL*60.*RHOA*CP/b4o A2(1)=Oo A4(1):Oo A3(1)=Oo A‘I)=Oo 8(1): 1. +X *VEL/(COEHT*A1(1))+CP*RHOA*VOL/(DELTAU*COEH T* 1A1(1)) C(I)='10 Q(1)=Oo D<1)=X *TA(KKK)* ( VEL)/(CUEHT*A1(1))+ CP*HHUA*VUL/(DELTAU*CUEH 1T*A1(I))*T(1) A1(2)=4.*PI*(HAu-(2*K-1)*UELTAR/2.)**2 A212)=4.*PI*RAD**2 A4(2)=(A1(2)+A2(2))/2o A(2)=BIOT*(-A2(2))/(RAD*DELTAR) 8(2): (A4(2)/THtTA+ A1(2))/(DELTAR**Zo)+BIOT*A2(2)/(RAD*DELTAR) C(2)=’ A1(a)/(DELTAP**2) 0(2)=A4(2)*20455*(100**(“4))*U10 **((T(2)—32o)/18o)*60o*HHO 75 94 C(2)=O(2)/THERMC 0(2): A4I2)/(THETA*(DELTAQ**2))*T(2) +O(2) NIP=NN~1 DO 5 N=3yNIP K=N~2 Al(N)=4.*PI*(RAD-(2*K-1)*DELTAR/2.)**2 A2(N)=4.*PI*(RAD-(2*K+l)*DELTAR/2.)**a A4(N)=(A1(N)+A2(N))/2. A(N)=-A1(N)/DELTAR**2 B(N)=(A4(N)/THETA+A1(N)+A2(N))/DELTAR**2 C(N)=-A2(N)/(DELTAR**d) 0(N)=A4(N)*2o455*10.**(-4)* U10 **((T(N)-32.)/18o)*60.*RHO Q(N)=O(N)/THERMC D(N)= A4(N)*T(N)/(THETA*(DELTAR**2)) + U(N) CONTINUE A(NN)=-6 *THETA B(NN)=1.+6.*THETA C(NN):OO O(NN)= 2.453*10.**(~4)* ulO **((T(N)-Ja.)/IB.)*60.*RHU l*DELTAU/(RHO*SPHT) D(NN)=T(NN) +Q(NN) DO 75 NzloNIP O(N)=0(N)*THERMC*UELTAU*DELTAH ORESR=D(N) + URLSP Q(NN)=Q(NN)*RHO*5RHT*a./3.*(DELTAR**3) ORESR=ORESR+U ORESRzHEAT OF RtsRIRATION IN bTU RETURN END SUBROUTINE TRIDI DIMENSION A(33).d(35)aC(35).A1(35)oAa(35)oA3(35)sA4(35)90(35). 10(35)9T(35).LC(JD)oDD(Jb) COMMON AQUQCOUQT9ORESP9NN9UELTAU THIS SUUROUTINE CALCULATES TH: TEMPERATURES OF THE NOUES USING THE CONSTANTS FROM SUUROUTINE SPHERE bY THE GAUSS ELIMINATION TECHNIUUE CC(1)=C(1)/B(1) DD(1)=D(1)/B(1) DO IK=29NN M3K-1 CCIK):C(K)/IU(K)“A(K)*CC(M)) DD(K):(D(K)"A(K)*DD(M))/(D(K)-A(K)*CC(M)) CONTINUE T(NN)=DD(NN) DO 2 KK=29NN KKK=NN~KK+1 KKKK =KKK+1 T(KKK):DD(KKK)“CC(KKK)%T(KKKK) RETURN END 0000 10 100 APPENDIX 3 PROGRAM GEOMETRY DIMENSION A(35)9D(35)9C(35)9A1(35)9A8(3b)9A3(35)9A4(35)9D‘35)9 IO(35)9T(3O)9CC(JS)9OD(J5) COMMON A9D9C9Q9T9UHESP9NN9UELTAU COMMON ISHAPE THIS PROGRAM CALCULATES THE TEMPERATURE HISTORY OF A SLAB9CYLINUEH OR A SPHERE NN=NUMBEP OF NOUES9NN-IZNUMDER OF DIVISIONS IN THE ODJECT DELTAU=TIME OIVISION9 NN=IO TIME = O. DELTAU=900E INITIAL TEMPERATURE UISTRIUUTION DO IO N=I9NN T(N)=960 BEGINNING OF PROGRAM CALCULATIONS DO 7 KK=I9ZO CALL SHAPE CALL TRIOI TIME =TIME+DELTAU PRINT 1009 TIME9(T(I)9I:19NN94) FOPMAT(IH 9FIO.O94(EIZ.593X)) CONTINUE END SUBPOUTINE SHAPE DIMENSION A(35)9b(35)9C(35)9A1(35)9A¢I3b)9A3(35)9A4(35)9DI3S)9 IQ(35)9T(35)9CC(JD)9DO(JD) . COMMON A 9 E59C9 I.) 9 T 9 QHESPo IVN9LJELTAU COMMON ISHAPE THIS SUDROUTINE CALCULATES THE COEFFICIENTS FOR THE UACKWAHU UIFF EPENCE EQUATION TA=4OQ PHO=5102 THERMC=.ZO3 CP3024 ' PHO=DENSITY LUS/CU FT 9 THCRMC2THEHMAL CONQUCTIVITY BTU/FT HP F 9 SPHT= SPECIFIC HEAT LTU/Lb F9 PAD: CHAPACTEPISTIC UIMENSION IN FT SPHT3086 ALPHA=THERMC/(RHO*5PHT) PAO=I.3E/Ic. DELTAP=PAO/(NN—I) COEHT=299 CLle 95 96 BIOT=COEHT*RAD/THERMC THETA=DELTAU*ALPHA/(UELTAR**Z¢) PI=391416 ' VOL=(180*18.*2964/(I7289*b40))‘(49/3.*PI*(RAD**3)) NIP=NN~1 ISHAPE=1 ISHAPE =1 FOR SLAB92 FOR CYLINDER9 3 FOR SPHERE IF(ISHAPE9EO¢I)J9II 3 DO 1 N=I9NIP Al‘N)=lO A2IN)=IO A4(N)=Io I CONTINUE A(1)=Oo 8(1)=Io/THETA+I9+UIOT*UELTAR/RAD C(I)=‘Io D‘l)=Io/THETA*T(I)+bIOT*DELTAR/RAD*TA GO TO 2 II CONTINUE IF (ISHAPEoEO02)17912 l7 CONTINUE CL=LENGTH OF CYLINDER AI(1)=29*PI*(HAU“(2*K-l)*OELTAH/29)*CL A2(1)=29*PI*HAU*CL A3(1)=Oo A4(1)=(A1(1)+A2(1))/£o A(1)=O. 8(1): (A4(1)/THETA+ A1(1))/(DELTAH**2.)+BIUT*A2(1)/(RAD*DELTAH) C(l)=- A1(1)/(DELTAR**2) D(1)= A4(1)/(THETA*(DELTAR**2))*T(1)+dIOT*A2(1)*TA/(DELTAR*RAD I) DO 101 N329NIP K=N-I AI(N):20*PI*(HAD*(E*N‘II*OELTAR/Eo)*CL A2(N)=20*PI*(HAd-(2%K+l)*UELTAH/ao)*CL A4IN)=(AI(N)+A2(N))/Eo 101 CONTINUE GO TO 2 I2 IF‘ISHAPEOEOOJ)I39I4 13 CONTINUE K31 AI(I)=4Q*PI*(RAD‘(2*K-I)*DELTAR/20)**E A2(1)=40*PI*RAD**2 A4(I)=(AI(I)+A2(I))/do A‘I)=Oo 8(1): (A4(1)/THETA+ A1(1))/(DELTAH**2.)+BIOT*A2(1)/(RAD*DELTAH) C(l)=— Al