 | Journal of Renewable Materials |  |
DOI: 10.32604/jrm.2023.021808
ARTICLE
Numerical Modelling of Drying Induced Cracks in Wood Discs Using the Extended Finite Element Method
1Research Institute of Wood Industry, Chinese Academy of Forestry, Key Lab of Wood Science and Technology of Status Forestry and Grassland Administration, Beijing, 100091, China
2School of Materials Science and Engineering, Harbin Institute of Technology, Harbin, 150001, China
*Corresponding Author: Yun Lu. Email: y.lu@caf.ac.cn
Received: 06 February 2022; Accepted: 30 March 2022
Abstract: Drying crack is a common phenomenon occurring during moisture discharge from wood, reducing efficient wood utilization. Drying crack is primarily caused by drying stress, and the reasonable methods for determining drying stress are sparse. In this study, the initiation and propagation of cracks during wood discs drying were simulated using the extended finite element method (XFEM). The distribution of drying stress and displacement was analyzed at different crack conditions based on the simulation results. This study aimed to solve the problem of the limitation of drying stress testing methods and provide a new idea for the study of wood drying stress. The numerical simulation results are found in good agreement with the experimental results, thus corroborating the feasibility of XFEM in modeling drying crack of wood discs. The stress concentration was observed at the crack tip region, while a minor stress was presented in the region of crack passing through, indicating that the crack formation process was also a process of releasing drying stress. Further, more energy was required to form double cracks in comparison with the single crack mode.
Keywords: Wood drying; queensland peppermint; drying cracks; numerical simulation; experimental validation; XFEM
Abbreviations
| XFEM | Extend finite element method |
| MPS | Maximum principal stress |
| MC | Moisture content |
| FSP | Fiber saturation point |
Drying crack is a common phenomenon for wood, which usually occurs during both the manufacturing process and use. The microcracks on wood surfaces affect the aesthetic appearance, while larger cracks can decrease wood strength, thus bringing about potential safety risks while in use. Due to the cracks on the surface, the interior of the timber structure is vulnerable to water and fungal attacks [1,2]. Therefore, the durability of wood in use is largely dependent on the cracks that develop during wood processing and due to surroundings. With the rise of wooden construction, wood logs are used as beams, columns, and other load-bearing components. For safety, better drying quality needs to be matched with this application. With wood as an anisotropic material, drying cracks emerge more easily in wood logs than in lumber, mainly presented as the radial cracks in the cross-section of wood logs [3,4]. The cracks are attributed to drying stress, generated by the uneven moisture distribution during water removal from wood and the anisotropic shrinkage between tangential and radial [5,6]. The effective means for detecting drying stress during the occurrence of cracks are scarce. However, the numerical simulations can essentially make up for the experimental limitations.
The finite element method (FEM) is a powerful tool for stress field analysis. However, there are some deficiencies while dealing with the fracture mechanics problems using the traditional FEM. The element boundary must coincide with the geometry edge of the cracks, but it is challenging to perform refinement of the mesh around a crack tip. Furthermore, it requires a remeshing technique to update the mesh to match the geometric discontinuities as the crack progresses, reducing the computational efficiency. To address these limitations, the extended finite element method (XFEM) was first applied to the crack growth problems by Belytschko et al. [7] which is considered a breakthrough in modeling crack propagation. The basic idea of XFEM is to add degrees of freedom in elements with special displacement functions, and the description for the discontinuous field is entirely mesh-independent [8,9]. Hence, XFEM can be used to simulate the propagation of a discrete crack along an arbitrary and solution-dependent path without requiring remeshing [10].
XEFM has been a convenient numerical tool to model fracture problems in wood science. Based on XFEM, Qiu et al. simulated the crack propagation in Glulam beams, and it was shown that the numerical results correlated well with experiments in terms of crack propagation direction and the applied forces for crack growth [11]. To estimate the load-bearing capacity of glued laminated timber, Berg et al. modeled the cracking behavior of beams by XFEM [12]. Hu et al. simulated the crack propagation in cross-sections of wood logs caused by relaxed growth stress and accumulated drying stress [13]. It was demonstrated that the crack behavior is closely related to the maximum principal stress and critical fracture energy of wood logs. The aforementioned research studies fully confirm that the XFEM is advantageous in modeling crack behavior. However, the research findings on cracks during wood drying are limited.
In this study, the XFEM was used to simulate the drying crack of Queensland peppermint wood discs. Based on the simulation results, the distribution of stress and displacement was analyzed at different drying crack conditions. This study aims to solve the problem of the limited drying stress testing methods and provide a new idea for the study of wood drying stress.
The Queensland peppermint (Eucalyptus exserta) wood discs were used to investigate crack behavior during the drying process. The diameter and thickness of the wood discs were about 200 and 50 mm, respectively. The initial moisture content (MC) of the wood was about 60%. In this experiment, eighteen wood discs were chosen (with great care to avoid defects and reaction wood) as test specimens, out of which three were used to determine the distribution of initial MC, whereas the fifteen were used for observing the initiation and propagation of drying cracks. The drying experiment was conducted in a 101–2AB type drying oven kept at the constant temperature of 80°C. The test specimens were taken out from the drying oven and observed cracking every 30 min.
The micro-mechanical testing machine (Instron 5848 Crop, Instron, USA) equipped with a load cell with a capacity of 2000 N was used to determine the tensile modulus of elasticity and modulus of rupture of wood in the tangential direction. Fifteen faultless specimens with the dimension of 76 mm (Tangential) × 15 mm (Radial) × 1.5 mm (Longitudinal) were selected for the tensile test. As shown in Fig. 1, the neck areas of the specimens were 2 mm × 30 mm (R × T), which ensured that there was a growth ring in the tensile area. The loading rate was set as 0.15 mm/min during the test. Before the tensile test, the specimens were put in a conditioning chamber for at least two weeks until the weight reached equilibrium. To obtain the test specimens with an MC of 12%, the temperature and relative humidity of the conditioning chamber were set as 20°C and 65%, respectively.
Figure 1: Specimen preparation (a) and micro-mechanical testing machine (b)
XFEM is based on the partition of unity concept, which is an extension of the traditional FEM [14]. XFEM enables local enrichment functions to be easily incorporated into a finite element approximation. For cracks, it consists of near-tip asymptotic functions that capture the singularity around the crack tip and a discontinuous function that represents the jump in displacement across the crack surfaces. The approximation for a displacement vector function with the partition of unity enrichment is shown in Eq. (1) [15].
where Ni(x) is the common nodal shape function; ui is the usual nodal displacement vector associated with the continuous part of the finite element solution; ai is the nodal enriched degree of freedom vector; H(x) is the associated discontinuous jump function across the crack surface; Fα(x) is the associated elastic asymptotic crack-tip functions; bai is the nodal enriched degree of freedom vector; N is the number nodes. In Eq. (1), ui applies to all the nodes in the model, H(x)ai is valid only for the nodes that are in the crack support domain, and Fα(x)bai is used only for nodes whose shape function support is cut by the crack tip.
Fig. 2 illustrates the normal and tangential coordinates for a smooth crack. The discontinuous jump function across the crack surface and the asymptotic crack tip function is given by Eqs. (2) and (3) [13,15].
Figure 2: Illustration of normal and tangential coordinates for a smooth crack
Cracks are generally categorized into three major fracture modes. Mode I fracture means the opening mode, where the tensile stress is normal (perpendicular) to the plane of the crack [16]. Mode II and Mode III represent a crack sheared due to in-plane shear stress and transverse shear stress respectively. Owing to the three orthogonal principal directions defined as longitudinal (L), radial (R), and tangential (T) axes, wood has six crack propagation systems. A pair of letters in the combination of L, R, and T is used to identify a crack propagation system, i.e., RL, RT, LR, LT, TL, and TR, respectively [17]. The first letter indicates the direction normal to the crack surface, while the second letter specifies the direction of crack propagation [11,18]. The radial crack along the wood ray is the major fracture of wood discs during drying. The crack is produced when the tangential tensile stress exceeds the tensile strength, where the stress is perpendicular to the plane of the cracks. Thus, the TR type crack is the primary fracture propagation of wood discs, which belongs to the opening mode (Mode I).
This study employs the theory of damage mechanics based on traction separation laws to simulate crack initiation and propagation. The crack initiation is determined by the principle of maximum principal stress (MPS). As shown in Eq. (4), when the MPS in the structure exceeds the principal stress value of the crack initiation, the crack is initiated [2]. The direction of the crack surface is perpendicular to the direction of MPS. As TR-type crack is the primary fracture propagation of wood disc, the cracks always propagate radially. Thus, the MPS value for crack initiation is equivalent to the modulus of rupture in the tangential direction.
where δtT and ftT are the tensile stress and tensile strength of wood, respectively, in the tangential direction.
The damage evolution law based on energy is used to describe the crack propagation. The crack propagation is determined by the energy release rate G, which is defined as the energy required to form a crack of unit area. GIC is considered as the resistibility to crack propagation, which is an inherent characteristic of materials. According to the energy balance relation in crack propagation, certain elastic strain energy is stored in the material when the object is deformed by force. This part of the energy will be released, accompanied by the propagation of cracks. The crack will grow only when the released elastic strain energy is greater than the surface energy required to form a new surface [2,15]. The precondition for crack propagation is illustrated as Eq. (5) [15].
where G is strain energy release rate, and GIC is critical strain energy release rate.
In this study, a two-dimensional plane geometric model was employed. According to the practical dimension of wood discs used in this experiment, the diameter of the circular region was set as 200 mm. Due to the distinctiveness of the pith position, a round area with a 20 mm diameter was removed from the whole circular region (Fig. 3). The material parameters required in this model are given in Table 1. These parameters used in this model are the standard values at 12% MC since some materials parameters, such as elastic modulus and tensile strength, are closely related to moisture content. The elastic modulus for tangential tension and MPS were obtained from the extension test. According to the previously published literature [19], the Poisson’s ratio of common hardwood species ranges from 0.231 to 0.496, and the average value is 0.332. Thus, the μTR was taken as 0.332 in this model.
Figure 3: Gridding partition method (a) and applied loading condition (b)
According to the previous literature, the density of Queensland peppermint wood at 12% MC is 773 kg/m3 [20]. As known, the GIC was determined by the density and MC of wood. The relationship is expressed in Eq. (6) [2,21]. Thus, the GIC of Queensland peppermint wood was determined as 0.669 at 12% MC.
where ρ12 is the wood density at 12% MC, and μ is wood MC.
Fig. 3a illustrates the mesh division of the assembly unit in this model, which includes a total of 484 elements. An element CPE4R (a 4-node bilinear plane strain quadrilateral, reduced integration element with hourglass control) in ABAQUS is used for element partition of the log cross-section. The load and boundary condition in this simulation model is presented in Fig. 3b. The boundary condition is set as a gradually increased radial displacement at the inner boundary. That is so because the loading method controlled by displacement is more stabilized than the force-controlled loading method in XFEM analysis [13].
4.1 Analysis of Drying Cracks in Wood Discs
For Queensland peppermint wood, there is much heartwood and only a small proportion of sapwood. The initial MC of Queensland peppermint wood discs is depicted in Fig. 4a. The region near the pith showed the highest initial MC, followed by the heartwood and sapwood area. The transition region between heartwood and sapwood came last. The location of crack initiation is shown in Fig. 4b. The results showed that the crack initiation was focused in the range of 10–30 mm from the pith, and then gradually expanded to the outer region as the drying proceeded. This phenomenon can be explained by the shrinkage difference between heartwood and sapwood. As Jiang et al. reported, the tangential air-drying shrinkage ratio at the region near the pith and bark was 7.44% and 4.84%, respectively [22]. The value of the shrinkage ratio also decreased from the pith to bark. Though the MC at the region near the bark dropped to fiber saturation point (FSP) ahead of the region near the pith, only a few minor cracks appeared here due to the smaller shrinkage ratio in this region. For the region near the pith, it is easy to produce one or more major cracks because of the unique material properties and greater shrinkage ratio in this region once the MC decreases to FSP. Further, these major cracks will propagate from the location of crack initiation to the outer region.
Figure 4: Initial MC of wood discs (a) and location distribution of crack initiation (b)
4.2 Simulation Results of a Single Crack
A single crack is the most common crack pattern in the drying process of wood discs, which is characterized by one main crack originating from the region near the pith and propagating to the outer region along the radial direction. The simulation results of stress and displacement during the formation process of single crack are shown in Fig. 5. The distribution of MPS at the three-stage, i.e., crack initiation, crack extending to the outer boundary, and after fracture conditions, is illustrated in Figs. 5a, 5b, and 5c, respectively. The corresponding displacement distribution of the above three stages is shown in Fig. 5d through Fig. 5f.
Figure 5: Simulation results of stress (a–c) distribution and displacement (d–f) distribution under single crack mode
As shown in Fig. 5, the values of MPS and displacement were smaller at the crack initiation stage, which were slightly larger at the region near the pith than in other regions. When the crack extended to the boundary, the stress value in the crack tip region was large, but the displacement value was small. This is because of the stress concentration occurring in the crack tip region before cracking, resulting in a considerably large stress value in this region. However, there was no cracking at this moment, so the displacement value was small. For the region of crack passing through, the stress value was small. However, the corresponding displacement was large, which is mainly attributed to the partial stress release after crack formation. When the fracture occurred in the wood discs, the stress was released entirely. At this stage, the stress value in the cracked area was small, and the displacement increased from the internal region to the outer region. Fig. 6 shows the experimental and simulation results under single crack mode. As can be seen, the numerical simulation results agreed well with the experimental results, demonstrating the applicability of XFEM in modeling drying cracks and drying stress.
Figure 6: Failure pattern with single crack: experimental results (a) and simulation results (b)
4.3 Simulation Results of Double Cracks
The diameter splitting crack pattern is another crack mode in wood discs, characterized by two opposite main cracks originating from the pith and propagating along the radial direction to the outer region. The final result is that the wood cross-section splits into two parts. Fig. 7a through Fig. 7c show the simulation results of MPS distribution at the crack initiation stage, cracks extended to the radius center, and cracks extended to the outer region, respectively. Fig. 7d through Fig. 7f show the corresponding displacement distribution for each stage. Similar to the single crack mode, the values of stress and displacement were small at the crack initiation stage. As the cracks extended to the radius center, the phenomenon of stress concentration was observed in the crack tip region, and minor stress was presented in the region of crack passing through. There is a larger displacement value in the cracked region for the displacement distribution at this stage but a smaller displacement value in the uncracked region. When the crack extended to the outer region, the cracking region was enlarged. Compared to the previous stage, the stress value on both sides of the cracks further decreased; however, the displacement value increased. The simulation results under double cracks mode were consistent with the experimental results, as shown in Fig. 8. When the same displacement load with single crack mode was applied to wood discs under the double-crack mode, the wood cross-section was not splitting into two parts. In effect, the double cracks were less likely to appear than the single crack in the drying process.
Figure 7: Simulation results of stress (a–c) distribution and displacement (d–f) distribution under double cracks mode
Figure 8: Failure pattern with double cracks: experimental results (a) and simulation results (b)
There are various kinds of crack patterns during wood discs drying, among which the single crack and double cracks are the two main cracking modes. In Queensland peppermint wood discs, the cracks quickly occurred near the pith region and then gradually expanded to the outer region with the drying proceeded. In this study, the initiation and propagation of single crack and double cracks in wood discs during drying were simulated by the XFEM. The numerical simulation results agreed well with the experimental results, demonstrating the useful applicability of XFEM in modeling drying cracks and drying stress. According to the simulation results of drying cracks, the stress concentration was observed at the crack tip region, while a minor stress was presented in the region of crack passing through. This indicated that the crack formation process is also a process of releasing drying stress. In addition, the wood cross-section was not splitting into two parts in double cracks simulation, if the applied load as much as single crack simulation.
Funding Statement: This study was financed by a Grant-in-Aid for Scientific Research from the Youth Program of National Natural Science Foundation of China (No. 31800478) and the National Natural Science Foundation of China (Grant Nos. 31870535 and 32122058)
Conflicts of Interest: The authors declare that they have no conflicts of interest to report regarding the present study.
References
1. Sandberg, D., Söderström, O. (2006). Crack formation due to weathering of radial and tangential sections of pine and spruce. Wood Material Science and Engineering, 1(1), 12–20. DOI 10.1080/17480270600644407. [Google Scholar] [CrossRef]
2. Chen, K. (2019). The study on decay and shrinkage cracks of timber components (Ph.D. Thesis). School of Civil Engineering Southeast UniversityChina. [Google Scholar]
3. Kang, W., Lee, N. H. (2002). Mathematical modeling to predict drying deformation and stress due to the differential shrinkage within a tree disk. Wood Science and Technology, 36(6), 463–476. DOI 10.1007/s00226-002-0153-5. [Google Scholar] [CrossRef]
4. Fu, Z., Zhao, J., Lv, Y., Huan, S., Cai, Y. (2016a). Stress characteristics and stress reversal mechanism of white birch (Betula platyphylla) disks under different drying conditions. Maderas Ciencia y Tecnología, 18(2), 361–372. DOI 10.4067/S0718-221X2016005000033. [Google Scholar] [CrossRef]
5. Hasalani, M., Itaya, Y. (1996). Drying-induced strain and stress: A review. Drying Technology, 14(5), 1011–1040. DOI 10.1080/07373939608917138. [Google Scholar] [CrossRef]
6. Fu, Z., Zhao, J., Yang, Y., Cai, Y. (2016b). Variation of drying strains between tangential and radial directions in Asian white birch. Forests, 7(3), 59. DOI 10.3390/f7030059. [Google Scholar] [CrossRef]
7. Belytschko, T., Black, T. (1999). Elastic crack growth in finite elements with minimal remeshing. International Journal for Numerical Methods in Engineering, 45(5), 601–620. DOI 10.1002/(SICI)1097-0207. [Google Scholar] [CrossRef]
8. Moës, N., Dolbow, J., Belytschko, T. (1999). A finite element method for crack growth without remeshing. International Journal for Numerical Methods in Engineering, 46(1), 131–150. DOI 10.1007/s00226-002-0153-5. [Google Scholar] [CrossRef]
9. Shi, J., Chopp, D., Lua, J., Sukumar, N., Belytschko, T. (2010). Abaqus implementation of extended finite element method using a level set representation for three-dimensional fatigue crack growth and life predictions. Engineering Fracture Mechanics, 77(14), 2840–2863. DOI 10.1002/(SICI)1097-0207. [Google Scholar] [CrossRef]
10. Zhuang, Z., Liu, Z. L., Cheng, B. B., Liao, J. H. (2014). Extended finite element method. Beijing: Tsinghua University Press. [Google Scholar]
11. Qiu, L. P., Zhu, E. C., Van de Kuilen, J. W. G. (2014). Modeling crack propagation in wood by extended finite element method. European Journal of Wood and Wood Products, 72(2), 273–283. DOI 10.1007/s00107-013-0773-5. [Google Scholar] [CrossRef]
12. Berg, S., Sandberg, D., Ekevad, M., Vaziri, M. (2015). Crack influence on load-bearing capacity of glued laminated timber using extended finite element modelling. Wood Material Science & Engineering, 10(4), 335–343. DOI 10.1080/17480272.2015.1027268. [Google Scholar] [CrossRef]
13. Hu, X. Y., Liu, Z. L., Zhuang, Z. (2017). XFEM study of crack propagation in logs after growth stress relaxation and drying stress accumulation. Wood Science and Technology, 51(6), 1447–1468. DOI 10.1007/s00226-017-0943-4. [Google Scholar] [CrossRef]
14. Melenk, J. M., Babuška, I. (1996). The partition of unity finite element method: Basic theory and applications. Computer Methods in Applied Mechanics and Engineering, 139(1–4), 289–314. DOI 10.1016/S0045-7825(96)01087-0. [Google Scholar] [CrossRef]
15. A.S.U. Manual, Abaqus 6.14 (2014). [Google Scholar]
16. Larsen, F., Ormarsson, S. (2014). Experimental and finite element study of the effect of temperature and moisture on the tangential tensile strength and fracture behavior in timber logs. Holzforschung, 68(1), 133–140. DOI 10.1515/hf-2012-0149. [Google Scholar] [CrossRef]
17. Larsen, F., Ormarsson, S., Olesen, J. F. (2011). Moisture-driven fracture in solid wood. Wood Material Science and Engineering, 6(1–2), 49–57. DOI 10.1080/17480272.2010.532234. [Google Scholar] [CrossRef]
18. Smith, I., Vasic, S. (2003). Fracture behaviour of softwood. Mechanics of Materials, 35(8), 803–815. DOI 10.1016/S0167-6636(02)00208-9. [Google Scholar] [CrossRef]
19. Kretschmann, D. (2010). Mechanical properties of wood. In: Wood handbook: Wood as an engineering material. Madison: US Dept. of Agriculture, Forest Service, Forest Products Laboratory. [Google Scholar]
20. Research Institute of Wood Industry, Chinese Academy of Forestry (1982). Wood physical and mechanical properties of major tree species in China. Beijing: China Forestry Press. [Google Scholar]
21. Danielsson, H. (2013). Perpendicular to grain fracture analysis of wooden structural elements-models and applications (Ph.D. Thesis). Lund UniversitySweden. [Google Scholar]
22. Jiang, X., Ye, K., Lv, J., Zhao, Y., Yin, Y. (2007). Wood properties and processing of eucalyptus and acacia plantation in China. Beijing: Science Press. [Google Scholar]
 | This work is licensed under a Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. |