Finite Difference Evaluation of Moisture Profile in Boxed-heart Large-cross-section Square Timber of Pinus densiflora during High Temperature Drying

1. INTRODUCTION

Log-wood and green lumber contain a substantial amount of moisture. If they are used without drying, dry defects arise, such as surface checking, twisting, and bowing. Therefore, wood should be dried to the in-use moisture content. However, the wood drying process consumes considerable energy and time. Hence, it would benefit from optimization. Such optimization necessitates theoretically and experimentally explaining the amount of energy input by the size of the wood, moisture distribution and temperature distribution inside the wood.

There have been many studies on the moisture transfer in wood (Siau, 1995). In wood, the water exists in three phases (free/liquid water, bound water, and water vapor) and the moisture is moved mainly through the free water flow above the fiber saturation point (Tiemann, 1921; Salin, 2006). Due to the formation of the moisture gradient, the diffusion of bound water moves through the cell wall, and the movement of water vapor occurs due to the relative vapor pressure gradient (Siau, 1995). Studies have been carried out to predict the moisture transfer in the range of moisture content below the fiber saturation point and to analyze the moisture transfer regarding the permeability of the wood cell wall and the pits (Kang et al., 2008; Stamm, 1960). Analysis of moisture transfer in wood during high temperature drying has been studied by Soma et al. (2013) and Cai (2010). However, there are few studies on moisture transfer in wood at temperatures above 100°C for domestic woods.

The finite difference method (FDM) is a robust numerical approach to solving the heat, moisture and mass transfer problem. The moisture transfer can be expressed as a finite difference, using the mass balance in the wood. The FDM method has been used for one-dimensional (Avramidis et al., 1994), two-dimensional (Schnabel et al., 2017; da Silva et al., 2014) and three dimensional (Salin, 2006; da Silva et al., 2011; da Silva et al., 2014; Younsi et al., 2006) analysis of changes of temperature or moisture content in wood during drying. The disadvantage of the FDM is that it is difficult to formulate irregularly shaped objects. However, it is relatively easy to analyze moisture transfer in square timber using the FDM because it is dried to stereotyped form shape.

High temperature drying is a method for rapidly drying wood at 100°C or higher, which can shorten the drying period by about 4 to 10 times compared to conventional kiln drying. However, during high temperature drying of wood with a large cross-section, the moisture gradient between the center of the wood and the surface becomes very high, which can lead to defects due to dry stress. Therefore, evaluating the moisture transfer within the wood during the high temperature drying is essential for the analysis of the internal stress distribution of the wood. In this study, the moisture transfer in boxed heart square timber of red pine was evaluated, during high temperature drying, by the three-dimensional FDM.

2. MATERIALS and METHODS

2.1. Three-dimensional FDM

FDM was used to predict the internal moisture content of the wood. The governing equations for mass transfer in three dimensions can be expressed by Eq. (1), using the mass balance in the unsteady state.

$$D_x\frac{{\partial }^{2}C}{{\partial x}^2}+D_y\frac{{\partial }^{2}C}{{\partial y}^2}+D_z\frac{{\partial }^{2}C}{{\partial z}^2}=\frac{\partial C}{\partial \mathit{{\rm T}}}$$

················ (1)

where D_x = diffusion coefficient in x direction

(m²/s); D_y = diffusion coefficient in y direction

(m²/s); D_z = diffusion coefficient in z direction

(m²/s); C = mass concentration (kg/m³); x, y

and z = distance in the direction of flow (m),

and τ = time (s).

If the object is incremented in the x-, y-, and z-axis directions, the intersection can be designed as shown in Fig. 1. Furthermore, the derivative terms of Eq. (1) can be approximated as shown in Eq. (2)~(5).

Figure 1. Schematic of the nodes used in the three-dimensional finite difference method.

Download Original Figure

$$\frac{{\partial }^{2}C}{{\partial x}^2}\approx \frac{C_{m+1,n,l}^p+C_{m-1,n,l}^p-{2C}_{m,n,l}^p}{{\left(\mathrm{\Delta }x\right)}^2}$$

·············· (2)

$$\frac{{\partial }^{2}C}{{\partial y}^2}\approx \frac{C_{m,n+1,l}^p+C_{m,n-1,l}^p-{2C}_{m,n,l}^p}{{\left(\mathrm{\Delta }y\right)}^2}$$

·············· (3)

$$\frac{{\partial }^{2}C}{{\partial z}^2}\approx \frac{C_{m,n,l+1}^p+C_{m,n,l-1}^p-{2C}_{m,n,l}^p}{{\left(\mathrm{\Delta }z\right)}^2}$$

·············· (4)

$$\frac{\partial C}{\partial \mathit{{\rm T}}}\approx \frac{C_{m,n,l}^{p+1}-C_{m,n,l}^p}{\mathrm{\Delta }\mathit{{\rm T}}}$$

·············· (5)

where, m, n and l = number representing locations in each direction; p = number representing pth time, p+1 = (p+1)th time.

Therefore, when Eq. (1) is expressed in the finite difference form, the differential equation at the internal node inside the wood is as follows [Eq. (6)]:

$$\begin{array}{c}{D}_x\frac{C_{m+1,n,l}^p+C_{m-1,n,l}^p-2C_{m,n,l}^p}{{\left(\mathrm{\Delta }x\right)}^2}\\ +D_y\frac{C_{m,n+1,l}^p+C_{m,n-1,l}^p-2C_{m,n,l}^p}{{\left(\mathrm{\Delta }y\right)}^2}\\ +D_z\frac{C_{m,n,l+1}^p+C_{m,n,l-1}^p-2C_{m,n,l}^p}{{\left(\mathrm{\Delta }z\right)}^2}\\ =\frac{C_{m,n,l}^{p+1}-C_{m,n,l}^p}{\mathrm{\Delta }\mathit{{\rm T}}}\end{array}$$

·············· (6)

The node of the surface, the edge and vertex of the wood is different from that described in Eq. (6), and the moisture is emitted from the surface of the wood to the ambient air. Thus, it can be added to the mass balance to design the equation. At the node, representing the outer surface boundary, where surface moisture emission occurs in one direction, we can formulate the following [Eq. (7)]:

$$\begin{array}{c}{D}_x\mathrm{\Delta }y\mathrm{\Delta }z\frac{C_{m-1,n,l}^p-C_{m,n,l}^p}{\mathrm{\Delta }x}\\ +D_y\frac{\mathrm{\Delta }x}{2}\mathrm{\Delta }z\frac{C_{m,n+1,l}^p+C_{m,n-1,l}^p-{2C}_{m,n,l}^p}{\mathrm{\Delta }y}\\ +D_z\frac{\mathrm{\Delta }x}{2}\mathrm{\Delta }y\frac{C_{m,n,l+1}^p+C_{m,n,l-1}^p-{2C}_{m,n,l}^p}{\mathrm{\Delta }z}\\ +S\mathrm{\Delta }y\mathrm{\Delta }z\left(C_e-C_{m,n,l}^p\right)\\ =\frac{\mathrm{\Delta }x}{2}\mathrm{\Delta }y\mathrm{\Delta }z\frac{C_{m,n,l}^{p+1}-C_{m,n,l}^p}{\mathrm{\Delta }\mathit{{\rm T}}}\end{array}$$

·············· (7)

where S = surface emission coefficient (m²/s); C_e = moisture concentration in equilibrium (kg/m³).

The mass transfer at the edge, where the moisture is emitted in two different direction can be expressed by Eq. (8).

$$\begin{array}{c}{D}_x\frac{\mathrm{\Delta }y}{2}\mathrm{\Delta }z\frac{C_{m+1,n,l}^p-C_{m,n,l}^p}{\mathrm{\Delta }x}\\ +D_y\frac{\mathrm{\Delta }x}{2}\mathrm{\Delta }z\frac{C_{m,n-1,l}^p+C_{m,n,l}^p}{\mathrm{\Delta }y}\\ +D_z\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }y}{2}\frac{C_{m,n,l+1}^p+C_{m,n,l-1}^p-{2C}_{m,n,l}^p}{\mathrm{\Delta }z}\\ +S\frac{\mathrm{\Delta }y}{2}\mathrm{\Delta }z\left(C_e-C_{m,n,l}^p\right)\\ +S\frac{\mathrm{\Delta }x}{2}\mathrm{\Delta }z\left(C_e-C_{m,n,l}^p\right)\\ =\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }y}{2}\mathrm{\Delta }z\frac{C_{m,n,l}^{p+1}-C_{m,n,l}^p}{\mathrm{\Delta }\mathit{{\rm T}}}\end{array}$$

·············· (8)

The differential equation at the vertex node, at which water is emitted from three surfaces of an object, is shown in Eq. (9).

$$\begin{array}{c}{D}_x\frac{\mathrm{\Delta }y}{2}\frac{\mathrm{\Delta }z}{2}\frac{C_{m+1,n,l}^p-C_{m,n,l}^p}{\mathrm{\Delta }x}\\ +D_y\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }z}{2}\frac{C_{m,n-1,l}^p+C_{m,n,l}^p}{\mathrm{\Delta }y}\\ +D_z\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }y}{2}\frac{C_{m,n,l+1}^p-C_{m,n,l}^p}{\mathrm{\Delta }z}\\ +S\frac{\mathrm{\Delta }y}{2}\frac{\mathrm{\Delta }z}{2}\left(C_e-C_{m,n,l}^p\right)\\ +S\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }z}{2}\left(C_e-C_{m,n,l}^p\right)\\ +S\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }y}{2}\left(C_e-C_{m,n,l}^p\right)\\ =\frac{\mathrm{\Delta }x}{2}\frac{\mathrm{\Delta }y}{2}\frac{\mathrm{\Delta }z}{2}\frac{C_{m,n,l}^{p+1}-C_{m,n,l}^p}{\mathrm{\Delta }\mathit{{\rm T}}}\end{array}$$

·············· (9)

2.2. Moisture Diffusion Coefficient

The moisture diffusion coefficient of the softwood was calculated by the geometric moisture diffusion model. This model is a modification of the electric circuit model of Stamm (1960), in consideration of the water vapor permeability in the pits and the rays (Kang et al., 2008). Fig. 2 shows the structure of water diffusion in the geometric moisture diffusion model. The moisture diffusion coefficient for each direction can be expressed as shown in Table 1. In the case of wood with a boxed heart, the radial direction and tangential direction are not fixed according to each node because the pith position is not constant. Consequently, the transverse moisture diffusion coefficient was determined, using the mean moisture diffusion coefficient of the radial and tangential directions. Assuming that the diameter of the tracheid (35 μm), the length of the rays (250 μm), and the size of the pits (1 μm) were constant, the ratio of pits to wood cells (p) and the ratio of rays to wood cells (r) were calculated, using the values measured by Eom (2015). The transverse bound-water diffusion coefficient of the cell wall (D_BT), the diffusion coefficient of the water vapor in the air within the lumens (D_v), and the longitudinal bound-water diffusion coefficient of the cell wall (D_BL) were determined by Eq. (10)~(13) (Siau, 1995).

Figure 2. Geometric moisture diffusion model (Kang et al., 2008).

Download Original Figure

Table 1. Three-directional moisture diffusion coefficient (Kang et al., 2008) Where a² = volume fraction of air in wood; p = ratio of pits to wood cells; r = ratio of rays to wood cells; D_BT = transverse bound-water diffusion coefficient of the cell wall (m²/s); D_v = diffusion coefficient (m²/s) of water vapor in the air within in the lumens; D_BL = longitudinal bound-water diffusion coefficient of the cell wall (m²/s)

Download Original Figure

$$D_{\mathit{BT}}=7\times {10}^{-6}exp\left[-\frac{38,500-290\mathit{MC}}{\mathit{RT}}\right]$$

·············· (10)

$$D_{\mathit{BL}}=2.5D_{\mathit{BT}}$$

·············· (11)

$$D_{\mathit{{\rm Y}}}=\frac{0.018D_a{P}_0}{G_M^w{\mathit{{\rm P}}}_w\mathit{RT}}\frac{\partial H}{\partial \mathit{MC}}$$

·············· (12)

$$D_a=\frac{2.2}{P}{\left(\frac{T}{273}\right)}^{1.75}$$

·············· (13)

where, D_BT = transverse diffusion coefficient of bound water in the cell wall (m²/s); D_BL = longitudinal diffusion coefficient of bound water in the cell wall (m²/s); D_v = diffusion coefficient of water vapor in the air within the lumen (m²/s); D_a = diffusion coefficient of water vapor in the air (m²/sec); R = universal gas constant (8314 m³Pa/(kgmol K)); T = temperature (K); MC = moisture content (%); G^w_M = specific gravity of cell-wall at moisture content MC; ρ_w = water density (kg/m³); P₀ = saturated vapor pressure (Pa); H = relative humidity (%), and P = air pressure (Pa).

2.3. Surface Emission Coefficient

The surface emission coefficient was determined using the convective mass transfer coefficients. The surface emission can be determined by applying the surface moisture content and the internal moisture content change to Eq. (14) (Yeo and Smith, 2005). At this time, it was assumed that the water vapor pressure of the outer surface of the wood is very similar to the saturated water vapor pressure at the wet bulb temperature.

$$S=h_{H_{2}0,\mathit{air}}\cdot \frac{\frac{M_{H_{2}0}{p}_s}{{\mathit{RT}}_s}-\frac{M_{H_{2}0}{p}_e}{{\mathit{RT}}_{\mathit{dry}}}}{G_{\mathit{SMC}}{\mathit{{\rm P}}}_w\frac{\mathit{SMC}}{100}-G_{\mathit{EMC}}{\mathit{{\rm P}}}_w\frac{\mathit{EMC}}{100}}$$

·············· (14)

where, h_H20,air = convective mass transfer coefficient; G_SMC and G_EMC = specific gravities of wood based on oven dry weight and volume at SMC and EMC; SMC = surface moisture content (%); EMC = equilibrium moisture content (%); M_H20 = water molecular weight (18 kg/kgmol); p_s = water vapor pressure in air adjacent to surface of wood (Pa); T_s = surface temperature (K); p_e = vapor pressure equilibrated to chamber condition (Pa); T_dry = dry bulb temperature (K). Fig. 3.

Figure 3. Illustration of nodes in wood with a large cross section.

Download Original Figure

3. RESULTS and DISCUSSION

During high temperature drying, the temperatures of red pine at half the point between the surface of the wood, the center, and the surface are shown in Fig 4. The air temperature between the woods was measured to be about 5 to 10°C higher than the target dry bulb temperature suggested by drying schedule. As a result of measuring the temperature of the center of the wood and temperature at the half point between the surface and the center, the temperatures inside the wood did not reach the dry bulb temperature until the end of drying, unlike the wood surface temperature. This result is because the heat of evaporation influences the heat transfer inside the wood, up to the later stage of high temperature drying. The temperature of the half point gradually increased after 24 h and remained almost the same as the surface temperatures after 196 h.

Figure 4. Temperature of red pine timber according to drying time.

Download Original Figure

Based on the measured temperature and coefficients relevant to mass transfer, the moisture distribution in the wood was predicted. The moisture change according to drying time is expressed in Fig. 5, using the fractional moisture content (E bar, E = (AMC-EMC) / (IMC-EMC), where, IMC = initial moisture content (%); EMC = equilibrium moisture content (%); AMC = average moisture content (%)). Fractional moisture content is defined as the ratio of the average moisture content to initial moisture content.

Figure 5. Fractional moisture content of red pine timber during drying.

Download Original Figure

The measured fractional moisture content at 24 h was higher than the predicted value, and the predicted 72 h result was similar to the measured value. However, the measured fractional moisture content at 196 and 240 h were lower than the predicted results.

Fig. 6 and Fig. 7 show the measured temperature profile, predicted and measured moisture distributions in longitudinal direction and transverse direction of red pine timber. Fig. 6(a) and Fig. 7(a) show temperature profile in longitudinal direction at each time. The temperature of each node was used to determine moisture diffusion coefficient and surface emission coefficient. Fig. 6(b) and Fig 7(b) show the predicted moisture distribution in longitudinal direction. The moisture content difference between the end part and the center part of the wood in longitudinal direction was not large. This is because the end coating was applied to the cross section of the wood.

Figure 6. Measured temperature profile (a), predicted moisture distribution (b), and predicted (c) and actual (d) transversal moisture distributions in red pine timber at 72 hours.

Download Original Figure

Fig. 6(c) and 6(d) show the predicted and measured moisture distribution of the red pine timber at 72 h of drying, respectively. When the 72 h moisture content distributions were compared, the predicted results showed an average moisture content of about 13% and an actual average moisture content of about 10%, generating a 3% error. The simulation predicted that the moisture content gradient at each site would be comparatively gentle.

Fig. 7(c) and 7(d) display the predicted and measured moisture content distribution in the wood at 240 h of drying. After drying, the final moisture content was 7.22%, and the experimental result was 2.67%. The moisture content was very low, being less than 2% on the surface, but it was predicted that a high moisture content of 10% or more would still exist inside the wood because of the surface checking, as at 72 and 196 h of drying, several surface checks were observed. Drying checks might be occurred by great moisture content gradient formed in the wood dried through high temperature drying schedule suggested in this study. Consequently, surface checks occurred, and surface emission was also generated inside the wood, which rapidly decreased the total moisture of the wood. This results show that the prediction of moisture content profile in wood calculated by using the moisture diffusion coefficient and the surface emission coefficient is possible only for the wood without checking.

Figure 7. Measured temperature profile (a), predicted moisture distribution (b), and predicted (c) and actual (d) transversal moisture distributions in red pine timber at 240 hours.

Download Original Figure

4. CONCLUSION

This study predicted the moisture movement in wood by the 3-D FDM during at high temperature drying. Based on the 3-D FDM, the moisture movement in the transverse and longitudinal direction of the wood was applied to the analysis of the moisture profile in the wood. The moisture diffusion coefficient was determined, by using the geometric moisture diffusion coefficient. The surface emission coefficient was calculated by considering the convective mass transfer coefficient. In order to compare the predicted results with the actual measurement results, we dried red pine timber at a high temperature and measured the average moisture content and the sectional moisture content of the wood at the time of the scheduled temperature change. Consequently, the prediction of the internal water movement of the wood up to 24 h, with a small wet bulb depression, showed a small error compared to the actual value, but thereafter, the actual fractional moisture content result was lower than the predicted result. This behavior is thought to be due to the checks on the surface of the wood. Therefore, the theoretical moisture diffusion coefficient and the surface emission coefficient used in this study are applicable for wood without checking.