Transient temperature profiles of frozen logs subjected to axisymmetric thawing and heating were calculated by a temperature method and an enthalpy method. The present paper discusses only the temperature method, which uses the conventional (temperature) formulation of the nonlinear heat conduction equation. This approach required the specification of a thawing temperature interval over which the latent heal was incorporated in the specific heat. Thermal properties were varied with position and temperature, and changed discontinuously with the phase. The log surface temperature was specified. The computerized finite-difference program HEAT was used in conjunction with this method. Computed temperature profiles were in overall agreement with experimental data obtained from heating logs in agitated water

Steinhagen, H. Peter

English

Wood and Fiber Science (E-Journal)

COMPUTERIZED FINITE-DIFFERENCE METHOD TO CALCULATE TRANSIENT HEAT CONDUCTION WITH THAWING H. Peter Steinhagen Associate Professor Department of Forest Products University of Idaho, Moscow, ID 83843 (Received July 1985) ABSTRACT Transient temperature profiles of frozen logs subjected to axisymmetric thawing and heating were calculated bv a temoerature method and an enthalov method. The oresent oaoer discusses onlv the . . . . temperature method, which uses the conventional (temperature) fornulation of the nonlinear heat conduction equation. This approach required the specification of a thawing temperature interval over which the latent heat was incorporated in the specific heat. Thermal properties were varied with position and temperature, and changed discontinuously with the phase. The log surface temperature was specified. The computerized finite-difference program HEAT was used in conjunction with this method. Computed temperature profiles were in overall agreement with experimental data obtained from heating logs in agitated water. Keywords: Heat transfer computer program, thawing model, phase change, nonlinear heat conduc- tion. INTRODUCTON Computing transient heat conduction with thawing is important in a number of wood processing operations. As an example, northern veneer and flakeboard mill operators request estimates of log thawing times for wood frozen during the winter months. The additional time required for raising veneer block temperatures to the levels recommended for veneer cutting must also be computed, as tem- perature measurements inside the logs are not feasible from an economic stand- point. Transient temperature profiles were calculated for frozen logs subjected to axi- symmetric thawing and heating. The two computational methods found suitable were a temperature method after Bonacina et al. (1973) and an enthalpy method derived from Voller and Cross (1981). The temperature method uses the con- ventional (temperature) formulation of the nonlinear heat conduction equation, with temperature as the only dependent variable. In the enthalpy method, the dependent variables are enthalpy and temperature. This paper discusses the tem- perature method. The method was used in conjunction with the computer program HEAT, which was prepared by Beckman (1972) in FORTRAN V language, for use on a UNIVAC 1108 computer. HEAT uses a finite-difference technique with predictor-corrector integrations (Hamming 1962) and can solve multidimensional heat conduction problems with various boundary and initial conditions. Thermal property data may vary with position and time or temperature. These data may also be dis- continuous, which is important for phase change calculations. Computations have been compared with data from experiments in which logs were thawed and heated in agitated water (Steinhagen 1977a). A more explicit Wwdond Fiber .Science. 18(3). 1986. DO. 460467  0 1986 by tho Soeicly of Wood Scicncc and Tcchnolwy S ~ ~ ~ ~ ~ ~ ~ - C A L C U L A T I N G  TRANSIENT HEAT CONDUCTION 46 1 t -03 FIG. 1. Log transient with imponant energy terms. description of the computational procedure is given in Steinhagen (1978), from which this paper was abstracted. THEORY Consider a long, cylindrical log with radius re (Fig. 1). The length of the log allows one to ignore any heating end effects; practically, this requires a length-to- diameter ratio in excess of four, according to data published by MacLean (1946). The log's initial temperature (To) is below 0 C, the freezing point of free water. At time zero, the log's surface is exposed to the heating bath temperature (T,) so that the surface temperature (T,) is immediately raised to the level of the bath temperature, involving surface thawing. The question, then, is how much log heating time (t) is required to reach the temperature (T) at the radial depth (r) along the log's center cross section. Note that the concept of two phases simultaneously existing during thawing is eliminated from consideration, given the assumption that the log's thermal con- ductivity (k), specific heat (c), and density ( p )  are functions of position and tem- perature and will discontinuously change with the phase. The latent heat is in- corporated in the specific heat over a finite thawing temperature interval (concept 462 WOOD AND FlBER SCIENCE. 1UI.Y 1986. V. 18(3) FIG. 2. Control systems with finite-difference terms of "equivalent specific heat," Bonacina et al. 1973). In this way, the thawing problem becomes equivalent to a problem of nonlinear heat conduction without phase change. An energy balance on a control system (Fig. 1, broken lines) gives where q 1, is the rate of heat flow at position r, q I ,+, is the rate of heat flow at position r + Ar, and E is the rate of energy storage. Substituting the conventional rate equations into the energy balance yields where A is the area through which heat flows, and V is the volume. Further development will yield the nonlinear equation for axisymmetric heat conduction: This equation cannot be solved because of the assumed temperature dependence of k, p, and c. However, an approximate solution is available via the finite- difference approach described below. FINITE-DIFFERENCE APPROACH A number of equidistant nodes are devised along the log radius (Fig. 2). Each node is successively considered the focal node i; nodes j l ,  j2 are the adjacent nodes. The heat flow from j into i is qji (=-q,), and the energy stored in the system surrounding node i is E,. One may then write an energy balance on each system and substitute the rate equations into the energy balances. The negative expression of the derivative in the conduction rate equations, -aT/ar, is approx- imated by the difference notation AT/Ar, for which one may write (Ti - Ti)/dj,, where d,, is the distance between j and i. This leads to a set of solvable equations analogous to Eq. (2): S I P ~ ~ ~ ~ ~ ~ ~ - C A L C U L A T I N C  TRANSIENT HEAT CONDUCTION 463 TO Te FIG. 3. Typical profile of C;, over temperature. dT, 4iAjii ,,,V.c. -= - kj2iAj2i (T,, - Ti) + -' d t  djii (Tj2 - Ti) dj2i Equation (4) may be rewritten as d ~ ,  A-jn pYfi - = 2 C..(T JI J - T,) dt j, where Cji is the conductance (kjiAj,/dji), and n is the number of nodes immediately surrounding focal node i (i.e., 2 in our one-dimensional system, and 4, 6 in two- and three-dimensional systems, respectively). Equation (5) is of the type that HEAT (Beckman 1972) will solve simultaneously for all nodes in a multidimen- sional system. As HEAT allows only the conductance (C,,) but not the capacitance (piVici) to be temperature-dependent, one must circumvent this shortcoming via dividing both sides of Eq. (5) by p,Vici, which yields .- ,=in where C,: is the conductance per unit capacitance (kj,Aji/dji)/(piVic,). The quantity C,;, in addition to the initial temperature (To) and surface tem- perature (T.), is computer input data and must be evaluated for the temperature points shown in Fig. 3 by the dots. The thawing temperature interval (AT,) shown in Fig. 3 is a computational parameter that must be adjusted so that agreement is obtained with experimental results. This thawing interval may or may not be fictitious. Within the thawing interval, the effective specific heat of a node (c,) is Li c. = c. + - ' "T, 464 WOOD AND FIBER SCIENCE. JULY 1986. V. 18(3) Tnem 1. Thermalproperties oflog No. 10. lnclvdcl latent heat. 'Surfam node. where c,,, is the true specific heat at the thawing temperature, and Li is the latent heat. For wood, we have where L, is the latent heat of water fusion (334,000 J/kg), and MC, is the node's moisture content expressed in percent of the dry mass of wood. This concept considers only the free water (MC, > 30%) as frozen, the bound water phase change being negligible (Kubler 1962). TABLE 2. Computer input data C,:, in seconds-', for log No. 10. j 12c - 1 . 5 - c  I . ~ ( + ) c  1 .5 ( - )C 1 . 5 ( + ) c  54 c 2 Ib  - - - - - - I 2 0.0018 0.0016 0.000075 0.000054 0.00088 0.00098 3 2 0.0016 0.0014 0.000064 0.000048 0.00077 0.00086 2 3 0.0018 0.0015 0.000082 0.000062 0.00089 0.00099 etc. etc. etc. etc. etc. etc. etc. etc. 19 18 0.0010 0.0009 0.000034 0.000023 0.00043 0.00048 18 19 0.0064 0.0057 0.000200 0.000140 0.00270 0.00310 Include8 l a t ~ n t  heal. Surface node I - I with rpcifird ternpraturc T. - 54 C needs no C;, vvalvo specification. S ~ ~ ~ ~ ~ ~ ~ ~ - C A L C U L A T I N O  TR NS ENT HEAT CONDUCTION 465 HEATING TIME ( h l  FIG. 4. Transient temperature profile of log No. 10. COMPUTATIONS Transient temperature profiles were computed for four logs identified in the above-mentioned experimental paper (Steinhagen 1977a) as Nos. 1 (aspen), 4 (black cherry), 7 (red oak), and 10 (eastern white pine). In accordance with the known log specifications, thermal conductive properties as a function of temper- ature, moisture content, and specific gravity were determined from the literature (Steinhagen 1977b). The thermal property values of log No. 10 are given in Table 1, for 19 nodes and an assumed thawing interval of 3 degrees (- 1.5 to 1.5 C),  with To = -22 C and T, = 54 C. These values were used to prepare the data of C;, (Table 2). These data were submitted to the computer in tabulated form (the format is given in the user manual for program HEAT, Beckman 1972). HEAT performed linear interpolation of the table values. As the smallest value of the Euler stability limit criterion (I/? C,:) was 1/0.0064 seconds (Table 2) which must J not he exceeded by the integration time step, a value of 0.040 hours was chosen as a time step in the computations for log No. 10. Note that the computer input for capacitance had to be unity, due to the transformation of Eq. (5) into Eq. (6). 466 WOOD AND FIBER SCIENCE, JULY 1986, V. 18(3) RESULTS AND DISCUSSION The transient temperature profile of log No. 10 is given in Fig. 4. The solid lines were established in the earlier experiments (Steinhagen 1977a). The symbols represent the comparable computed data. For a given position and temperature, the discrepancy between measured and computed time was usually within the 10% target, for all four logs studied in the project. The computational error seemed to increase when the latent heat entered the calculations. This is apparently due to the discontinuity and use of a fictitious thawing temperature interval to accommodate the latent heat effect. From a prac- tical point of view, it appears, however, that the method has provided acceptable temperature-time solutions to the problem of axisymmetric thawing and heating of logs with position- and temperature-dependent thermal properties. Increasing the number of nodes had no significant effect on the quality of results. Oscillations among the computed values were not observed in these computations. Calculations with constant, though distinct, thermal properties in each phase were considerably less cumbersome to prepare than with variable properties. The resulting error was usually again within target. Solutions have been generalized into charts (Steinhagen et al. 1980) covering a large range of values that the nondimensionalized parameters can assume. CONCLUSIONS Utilization of the temperature method in conjunction with computer program HEAT has provided temperature-time solutions to axisymmetric log thawing and heating eases. The agreement between calculated and measured temperature profiles appears satisfactory. Specification ofa thawing temperature interval for wood appears to be a problem with this method. ACKNOWLEDGMENTS The study was financially supported by the University of Wisconsin-Madison (UW) School of Natural Resources, College of Agriculture and Life Sciences, by funds provided under McIntire-Stennis Public Law No. 87-788 and by a grant from the University of Idaho Research Council. Special thanks go to John F. Lutz, U.S. Forest Products Laboratory, Madison, Wisconsin (now retired), and Professor Hans Kuhler, UW Dept. of Forestry, for instigating and supporting this research in all phases, as well as to Professor Glen E. Myers, UW Dept. of Me- chanical Engineering, for advice in the area of heat transfer analysis. REFERENCES BECKMAN, W. A. 1972. HEAT (Heat transfer computer program in FORTRAN V language). EES Repon No. 37, Engineering Experiment Station, University of Wisconsin, Madison, W1. BONACINA, C., G. COMINI, A. FASANO, AND M. PRIMICERIO. 1973. Numerical solution of phase- change problems. Int. 1. Heat Mass Transfer 16:1825-1832. HAMMING, R. W. 1962. Numerical methods for scientists and engineen. McGraw-Hill, New York, w ..., KUBLER, H. 1962. Shrinkage and swelling of wood by coldness. Holz Roh-Werkst. 20(9):364-368. MACLEAN, 1. D. 1946. Rate of temperature change in shon-length round timbers. ASME Transactions 68(1):1-16. Slcinhapm-CALCULATING TRANSIENT HEAT CONDUl3ION 467 STEINHAGEN, H. P. 1977a. Heatinp: times for frozen veneer loas-New experimental data. For. Prod. . J. 27(6):24-28. . 1977b. Thermal conductive ~roaerties of wood. ereen or dm. from -409 to +IOODC: A . . . - . . literature review. USDA For. Serv. Gen. Tech. Report FPL9, Forest Products Laboratory, Mad- ison, WI. . 1978. Heating times for frozen and nonfrozen veneer logs. Ph.D. thesis, University of Wisconsin-Madison. MJN78-04884, University Microfilms, Ann Arbor, MI. -, G. E. MYERS, AND H. KURLER. 1980. Heating time charts for frozen and nonfrozen veneer logs. For. Prod. J. 30(4):27-37. VOLLER, V., AND M. CROSS. 1981. Accurate solutions of moving boundary problems using the enthalpy method. Int. J.  Heat Mass Transfer 24545-556. 

Computerized Finite-Difference Method to Calculate Transient Heat Conduction with Thawing

https://wfs.swst.org/index.php/wfs/article/download/1260/1260

Computerized Finite-Difference Method to Calculate Transient Heat Conduction with Thawing

Abstract

Similar works

Full text

Available Versions

Wood and Fiber Science (E-Journal)