A 3D multi-scale hygro-mechanical model of oak wood

A multi-scale framework is proposed for the prediction of the macroscopic hygro-elastic properties of oak wood. The distinctive features of the current multi-scale approach are that: (i) Four different scales of observation are considered, which enables the inclusion of heterogeneous effects from the nano-, micro-, and meso-scales in the effective constitutive behavior of oak at the macro-scale, (ii) the model relies on three-dimensional material descriptions at each considered length scale, and (iii) a moisture-dependent constitutive assumption is adopted at the nano-scale, which allows for recovering the moisture dependency of the material response at higher scales of observation. In the modeling approach, oak wood is assumed as homogeneous at the macro-scale. The meso-scale description considers the cellular structure of individual growth rings with three different densities. At the micro-scale, the heterogeneous nature of cell walls is described by the characteristics of the primary and secondary cell wall layers. Finally, the nano-scale response is determined by cellulose micro-fibrils embedded in a matrix of hemicellulose and lignin. The oak properties at the four length scales are connected via a three-level homogenization procedure, for which, depending on the geometry of the fine-scale configuration, an asymptotic homogenization procedure or Voigt averaging procedure is applied at each level to determine the effective hygro-elastic properties at the corresponding coarse scale. In addition, the moisture adsorption isotherms at each scale are constructed from a volume-weighted averaging of the moisture adsorption characteristics at the scale below. The computational results demonstrate that the macro-scale moisture-dependent, hygro-elastic behavior of oak wood is predicted realistically, thereby revealing the influence of the material density, the micro-fibril orientation, and the hygro-elastic properties from the underlying scales. The computed macro-scale properties of oak are in good agreement with experimental data reported in the literature.


Introduction
Oak is one of the most used hardwoods in structural and building applications (Foster and Ramage 2017; Brischke et al. 2012), due to its relatively high strength and hardness and its high resistance to biological attacks (Carmona Uzcategui et al. 2020).Moreover, its aesthetic qualities make oak an appealing material for art and decoration purposes, including furniture, veneer decorations, and panel paintings (Ross 2010;van Duin and Kos 2014).Oak wood has a heterogeneous, hierarchical material structure that extends over multiple length scales, from the macroscopic, structural scale down to the nanoscopic, molecular scale (Blanchette et al. 1994;Han et al. 2019;Gibson 2012;Saavedra Flores et al. 2016).At the macro-scale, oak wood is characterized by concentric rings that are periodic in the radial material direction, known as annual growth rings (Stokes and Smiley 1996).Each growth ring consists of earlywood and (denser) latewood, which can be respectively distinguished by light and dark colors that reflect the different meso-structural morphologies.At the meso-scale, oak is characterized by a cellular structure consisting of various cell types, which include libriform fibers, fiber tracheids, axial parenchyma, vessels, and rays.Libriform fibers and fiber tracheids are long and slender hollow tubes oriented mainly along the axial (or longitudinal) direction of the stem of a tree, and serve to provide mechanical resistance.Libriform fibers typically have a thicker cell wall than fiber tracheids.Axial parenchyma cells have thinner cell walls and a shorter cell length than the fibers in longitudinal direction, and mainly lie in the regions between vessels.Vessels are tubular structures oriented along the longitudinal direction of wood, and have a relatively wide diameter and a shorter length than libriform fibers and fiber tracheids.Ray parenchyma, or rays, are prismatic cells oriented along the radial direction of the growth rings.In Pedunculate oak (Quercus robur), fiber tracheids, axial parenchyma, and large vessels can be found in earlywood regions.Conversely, latewood contains libriform fibers, axial parenchyma, and small vessels, while rays can be observed in either a uni-seriate or multi-seriate form in both earlywood and latewood (Elaieb et al. 2019;Kim and Daniel 2016;Wheeler et al. 1989).At the micro-scale, cell walls are distinguished by primary and secondary walls.The secondary wall comprises three layers, referred to as the S1, S2, and S3 layers.The S2 layer is the thickest of these three layers, and contributes the most to the mechanical resistance of oak.The primary wall is connected to the middle lamella layer, which bonds different adjacent wood cells.At the nanoscale, the cell wall layers consist of cellulose micro-fibrils bundled in fibril aggregates, which are embedded in a hemicellulose and lignin matrix (Rowell 1984).The micro-fibrils are oriented along single or multiple directions that are specific for a cell wall layer (Bergander and Salmén 2002).Cellulose, hemicellulose, and lignin are organic polymers.Cellulose is a polysaccharide consisting of a long linear chain of several glucose monomer units linked together by glycosidic bonds.Due to the high chain organization of these molecules, cellulose is hardly accessible to water, which makes it hydrophobic (Bergenstråhle et al. 2008;Matthews et al. 2006).Hemicellulose is a polysaccharide consisting of a short chain of different carbohydrate monomers.Hemicellulose is typically characterized by a considerable degree 1 3 Wood Science and Technology (2023) 57:1215-1256 of chain branching that is responsible for its non-crystalline structure, which results in a highly hydrophilic behavior (Siau 1983;Berry and Roderick 2005).Lignin is a complex phenolic polymer consisting of phenylpropane units commonly linked by ether bonds.Its molecular arrangement is strongly cross-linked, resulting in an amorphous structure.Lignin is hydrophilic, with a moderate capacity for moisture uptake (Cousins 1976).
The mechanical and hygro-expansive behavior of cellulose, hemicellulose, and lignin, together with the complex morphologies of oak across the different length scales, determine its macroscopic hygro-mechanical response under relative humidity variations.Relative humidity variations may result in dimensional changes of oak wood (Derome et al. 2012;Luimes et al. 2018a) that can promote fracture at the macroscopic level (Ross 2010).For example, recent studies on oak museum objects have shown that moisture-induced dimensional changes in cabinet door panels and panel paintings can potentially lead to cracking mechanisms in the panel itself, and in the decorative layers and paint layers (Luimes et al. 2018a, b;Bosco et al. 2020bBosco et al. , 2021;;Luimes and Suiker 2021).These studies also suggest that the understanding of the mechanical and hygro-expansive behavior of oak museum objects can be further improved via multi-scale approaches that incorporate the characteristics and diversity of the underlying material structures.
In the literature, several multi-scale models for wood have been proposed.For the analysis of softwood (e.g., spruce), which is characterized by a relatively regular cellular configuration, the proposed modeling approaches consider both twodimensional (Persson 2000;Rafsanjani et al. 2012Rafsanjani et al. , 2013;;Scheperboer et al. 2019;Livani et al. 2021) and three-dimensional (Qing and Mishnaevsky Jr 2009;Saavedra Flores et al. 2015, 2016;Pina et al. 2019;Rojas Vega et al. 2022) material structures that include various scales of observation, ranging from the scale of cellulose micro-fibrils to that of annual growth rings.Analytical and/or computational homogenization techniques have been applied for obtaining the effective material properties, occasionally incorporating the dependency of the mechanical properties of the constituent phases on the moisture content.Conversely, multi-scale models of hardwood account for the influence of more complex morphologies on the effective material properties, by including the contributions of ray cells (Perré and Badel 2003;de Borst and Bader 2014;Malek andGibson 2015, 2017) and/or vessels (Hofstetter et al. 2005;Livani et al. 2021).Nevertheless, most models assume somewhat idealized meso-structural geometries (Hofstetter et al. 2005;de Borst and Bader 2014;Malek andGibson 2015, 2017), such as those mimicking the complex behavior of balsa wood (Malek andGibson 2015, 2017;Livani et al. 2021).In Badel and Perré (2002), Perré and Badel (2003), and Badel and Perré (2006), more detailed meso-structural geometries of a single oak wood growth ring are considered.These models, however, rely on a two-dimensional description that provides the effective material properties with respect to the transverse anatomical plane of oak, thereby neglecting the material behavior along the radial and tangential planes.Furthermore, the homogenization procedure is applied at a single level and only considers the coupling between the cellular scale and the macroscopic, structural scale.The dependence of the constitutive properties on the moisture content is also ignored.
In the current paper, a novel multi-scale framework for oak wood is presented, which overcomes the limitations addressed above by including the following aspects: (1) Four different scales of observation are considered, which enables the inclusion of heterogeneous effects from the nano-, micro-, and meso-scales in the effective constitutive behavior of oak at the macro-scale, (2) the model relies on three-dimensional material descriptions at each considered length scale, and (3) a moisture-dependent constitutive assumption is adopted at the nano-scale, which allows for recovering the moisture dependency of the material response at higher scales of observation.The oak material properties at the four length scales are connected via a three-level homogenization procedure, for which, depending on the geometry of the fine-scale configuration at the specific level, an asymptotic homogenization procedure or Voigt averaging procedure is applied to determine the effective hygro-elastic properties at the corresponding coarse scale.Asymptotic homogenization determines the effective response of a heterogeneous solid from the response of the underlying, fine-scale unit cell by subjecting it to periodic boundary conditions (Bakhvalov and Panasenko 2012).Conversely, Voigt averaging straightforwardly estimates the effective response of a heterogeneous solid based on the assumption of a uniform strain field in the underlying material structure (Hill 1963).In addition, the moisture adsorption isotherms at a scale are constructed from a volume-weighted averaging of the moisture adsorption characteristics at the scale below.
The three-level homogenization procedure consists of three subsequent steps.In the first homogenization step, asymptotic homogenization is applied to obtain the hygro-mechanical properties of individual oak cell wall layers from the response of corresponding, three-dimensional periodic nano-scale unit cells.These nano-scale unit cells consist of a cellulose micro-fibril embedded in a matrix of hemicellulose and lignin.For the hydrophilic hemicellulose and lignin, the elastic properties are assumed to be moisture-dependent, and moisture adsorption isotherms are specified to relate the moisture content to the ambient relative humidity.These nano-scale moisture adsorption isotherms are further used to construct the effective moisture adsorption isotherms of the micro-scale cell wall layers by using a weighted volume averaging scheme.In the second homogenization step, a micro-scale representation of a cell wall is considered, which is modeled as a three-dimensional composite plate consisting of a stack of cell wall layers.Each cell wall layer is characterized by the hygro-elastic properties and the characteristic moisture adsorption isotherm computed from the nano-scale model.With this input, the effective hygro-mechanical properties of the cell walls are determined through a Voigt averaging procedure.Further, the effective moisture adsorption isotherm is constructed from a weighted volume average of the moisture adsorption characteristics of the individual cell wall layers.In the third homogenization step, accurate three-dimensional meso-scale unit cells describing growth rings of different densities are generated from detailed microscopy images of oak samples.The application of asymptotic homogenization to these meso-scale configurations ultimately determines the hygro-mechanical properties of oak at the macroscopic scale.Additionally, a weighted volume average of the meso-scale moisture adsorption behavior provides the macro-scale moisture adsorption isotherm.

3
Wood Science and Technology (2023) 57:1215-1256 The paper is organized as follows.The proposed multi-scale modeling approach is described in detail in Sect."Multi-scale modeling of oak wood".The nano-scale unit cell and the (moisture-dependent) material properties of the nano-scale constituents are presented in Sect."Nano-scale model".The micro-scale model and the moisturedependent hygro-elastic properties of individual cell wall layers obtained from the first homogenization step are discussed in Sect."Micro-scale model".The meso-scale unit cell and the moisture-dependent hygro-mechanical characteristics of cell walls deduced from the second homogenization step are presented in Sect."Meso-scale model".The moisture-dependent hygro-elastic properties of oak at the macro-scale, as determined from the third homogenization step, are treated in Sect."Macro-scale results".Finally, Sect."Conclusion" summarizes the main conclusions of the study.
Throughout the paper the following notation for tensors and tensor products is employed: a, a, A and n A denote, respectively, a scalar, a vector, a second-order tensor, and an nth-order tensor.By making use of Einstein's summation convention, vector and tensor operations are defined as follows: A dyadic product is represented by a ⊗ b = a i b j e i ⊗ e j , and inner products are defined as , with e i ( i = x, y, z ) the unit vectors of a Cartesian basis.The symbols ∇ and ∇⋅ respectively denote the gradient and the divergence operators, as defined by ∇f

Multi-scale geometry
For the analysis of the mechanical and hygro-expansive properties of oak wood within a multi-scale framework, four characteristic length scales are distinguished, see Fig. 1.At the macro-scale, which is defined by the characteristic length scale L , oak wood is assumed to be a three-dimensional homoge- neous medium characterized by effective hygro-elastic properties, see Fig. 1a.The meso-scale unit cell, which is defined by the characteristic length scale L, reflects the cellular structure of an entire growth ring in the transverse plane, which includes earlywood and latewood, see Fig. 1b.The meso-scale geometry is generated from detailed microscopy images taken from oak samples, which include different types of cells (i.e., libriform fibers, fiber tracheids, axial parenchyma, and rays) and vessels.The two-dimensional domain constructed from the microscopy image is extruded along the third, longitudinal (axial) direction, and is subsequently made periodic in all three directions.At the micro-scale, which is defined by the length scale l, the heterogeneous structure of cell walls is considered, see Fig. 1c, consisting of primary and secondary walls.The primary wall is connected to the middle lamella, which bonds different adjacent wood cells.Since it is practically cumbersome to distinguish the primary walls of two adjacent wood cells from the middle lamella in between, similar to other works in the literature, e.g., Malek and Gibson (2017), these three cell wall layers are considered as a single object, referred to as the compound middle lamella (CML) layer.The secondary wall is composed of S1, S2, and S3 cell wall layers that have different properties.The micro-scale unit cell is idealized as a composite laminate consisting of two stacks of secondary walls with a compound middle lamella layer in between.At the nano-scale, which is characterized by the length scale , cell wall layers are composed of periodic bundles of cellulose micro-fibrils that are embedded in a matrix of hemicellulose and lignin, see Fig. 1d.The composite fiber illustrated in Fig. 1d is analogous to the configurations considered in Terashima et al. (2009) and Malek and Gibson (2017).

Review of the applied homogenization methods
The effective hygro-mechanical properties of oak wood are determined through a multi-level homogenization procedure, where, depending on the geometry of the fine-scale configuration at a specific level, an asymptotic homogenization procedure (Palencia 1980;Bakhvalov and Panasenko 2012;Peerlings and Fleck 2004;Bacigalupo et al. 2016;Fantoni et al. 2017) or Voigt averaging procedure (Jones 1998;Gao and Oterkus 2019) is applied to determine the effective hygro-elastic properties at the corresponding coarse scale.In the following, a review of the two homogenization methods is presented.Asymptotic homogenization applies at the scales at which the underlying fine scale can be considered periodic, i.e., at the nano-and mesoscales.Conversely, Voigt averaging needs to be applied at the micro-scale, since the composite cell wall is non-periodic in the direction normal to the cell wall.

Asymptotic homogenization procedure
Consider a composite material with a heterogeneous fine-scale structure, see Fig. 2b.The coarse-scale three-dimensional domain Ω is constructed as a peri- odic repetition of the fine-scale unit cell Q depicted in Fig. 2a.Denoting L and L = L as the characteristic lengths of the coarse-and fine-scale domains, the requirement of a strong separation of scales is represented by the condition  ≪ 1 .The field quantities (displacement, strain, and stress) governing the hygro- mechanical behavior vary smoothly at the coarse scale, while they are periodic at the fine scale.Consequently, the field quantities may be assumed to explicitly depend on two variables: A slow variable X that varies at the coarse-scale level Fig. 2 Reference domain for the asymptotic homogenization procedure: a coarse-scale domain Ω with a characteristic length scale L and b underlying periodic fine-scale unit cell Q with a characteristic length scale L. The coordinate systems used at the coarse and fine scales are X = (X, Y, Z) and x = (x, y, z) , respectively and a fast variable x = X∕ that varies at the fine-scale level.In the asymptotic homogenization procedure, which is based on the framework developed in Bosco et al. (2017aBosco et al. ( , 2017bBosco et al. ( , 2020a)), the heterogeneous material at the fine scale is modeled as an equivalent homogeneous solid at the coarse scale.The procedure starts from the mechanical equilibrium equation, which, in the absence of body forces, is expressed by where is the Cauchy stress tensor that satisfies the hygro-elastic constitutive relation: Here, u(X) is the displacement field and ΔRH(X) is the ambient, coarse-scale relative humidity variation, which is equal to the difference between the current relative humidity RH and the initial reference relative humidity RH 0 , i.e., ΔRH = RH − RH 0 .The material properties at the fine-scale level determine the coarse-scale response and are represented by local constitutive tensors that are periodic functions of the fine-scale variable x , i.e., the elastic stiffness tensor 4 C = 4 C(x) and the hygroscopic expansion tensor = (x) .Although not explicitly indicated in Eq. ( 2), the local stiffness tensor 4 C and the hygroscopic expansion tensor can also depend on the current value of the relative humidity RH .Note that the hygroscopic strain is commonly defined in terms of a moisture content variation, instead of a relative humidity variation as done in Eq. ( 2).The formulation in Eq. ( 2) is imposed by the homogenization procedure, which requires defining the effective properties as a function of a uniform field over the fine-scale domain.The moisture content depends on the water adsorption capacity of each of the constituent phases, and thus is non-uniform across the fine-scale configuration.Conversely, the relative humidity is uniform across the scales, and therefore is used as the reference field for expressing the hygroscopic strain in Eq. (2).
The displacement field is formulated as an asymptotic expansion in terms of : Expression (3) is inserted into the constitutive relation (2), after which the result is substituted into the equilibrium equation ( 1), leading to a set of boundary value problems-the so-called "cell problems"-defined at the fine scale, i.e., where 4 I S is the fourth-order symmetric identity tensor, given by I S ijkl = ( il jk + ik jl )∕2 , with ij the Kronecker delta symbol.The solution of these cell problems provides the fine-scale influence functions 3 N 1 and b 1 , which represent the periodic fluctuations of the displacement field within the fine-scale domain, as (1) Wood Science and Technology (2023) 57:1215-1256 caused by respective variations of the coarse-scale deformation and the relative humidity.
The cell problems in Eqs. ( 4) and ( 5) are completed by periodic boundary conditions expressed in terms of the unknown influence functions 3 N 1 and b 1 .Peri- odic boundary conditions require that the values assumed by each influence function are equal at corresponding points on opposite boundaries (left and right, top and bottom).In addition, for warranting the uniqueness of the solution, the influence functions are prescribed to vanish on average across the dimensions of the unit cell.Finally, the cell problems are solved numerically by the finite element method, thereby considering a specific material configuration that represents the fine-scale domain.Once the influence functions 3 N 1 and b 1 have been obtained from Eqs. ( 4) and ( 5), the effective elastic tensor 4 C and hygroscopic expansion tensor ̄ at the coarse scale can be computed as Further details on the derivation of the effective properties, the solution of the cell problems, and the finite element implementation of the method can be found in Bosco et al. (2017a) and Yuan and Fish (2008).

Voigt averaging procedure
Consider a three-dimensional laminate Ω composed of a finite number of layers.Each layer k is characterized by uniform material properties, defined through the elastic stiffness tensor 4 C k and the hygroscopic expansion tensor k .The effective properties of the laminate are obtained from a Voigt averaging procedure, which assumes full kinematic compatibility between the composing layers.Hence, the strains in the individual layers are supposed to be equal, and thereby are equal to the effective strain in the laminate.Furthermore, the effective stress in the laminate is consistently determined by the average of the stresses in the individual layers weighted by the layer thicknesses.The effective elastic stiffness tensor 4 C and the effective hygroscopic expansion tensor β of the multi-layer laminate then follow as (Gao and Oterkus 2019) with t k the relative thickness of layer k, as determined with respect to the total thick- ness of the laminate.

Three-level homogenization approach
The proposed multi-scale model of oak wood includes four different scales of observation, see Fig. 1, which are connected via a three-level homogenization approach.Accordingly, in the first homogenization step the effective micro-scale hygro-elastic properties of individual cell wall layers are computed by applying asymptotic homogenization to the response of the nano-scale unit cell composed of a cellulose fibril embedded in a matrix of hemicellulose and lignin.Following the terminology introduced in Sect."Asymptotic homogenization procedure", the micro-scale thus corresponds to the coarse scale and the nano-scale to the fine scale.The nano-scale unit cell has been depicted in Fig. 1d and has periodic boundary conditions in three directions.The asymptotic homogenization procedure departs from generating a nano-structural fiber geometry for each cell wall layer, which is characterized by specific volume fractions of cellulose, hemicellulose, and lignin.The hygro-elastic properties of each constituent phase are reflected by the tensors 4 C and , where the elastic parameters of the hydrophilic hemicellulose and lignin phases are dependent on the relative humidity.Next, the unit cell problems given by Eqs. ( 4) and ( 5) are solved for these nano-scale structures, resulting in the values of the influence functions.These values are subsequently substituted into relations ( 6) and ( 7) to calculate the effective hygroelastic properties 4 Clay and ̄ lay of each cell wall layer.In the second homogenization step, the effective properties of oak cell walls are computed at the meso-scale by applying a Voigt averaging of the micro-scale properties of the cell wall layers.The meso-scale thus corresponds to the coarse scale and the micro-scale to the fine scale.Note that asymptotic homogenization cannot be applied here, since the micro-structure is non-periodic in the direction normal to the cell wall.The assumed micro-scale domain is a heterogeneous plate composed of two stacks of secondary layers with a compound middle lamella layer in between, as shown in Fig. 1c.Two different micro-structural geometries are considered, which refer to the micro-fibril orientations and layer thicknesses characteristic of i) axial cells (i.e., libriform fibers, fiber tracheids, axial parenchyma) and ii) ray cells.The hygro-elastic properties 4 Clay and ̄ lay of each cell wall layer depend on its specific micro-fibril orientation and wall thickness.To obtain the effective properties of the cell walls, the elastic tensor 4 Clay and hygro- scopic expansion tensor ̄ lay computed for the cell wall layers in the first homogenization step are converted to a global reference system associated with the cell wall, providing 4 Cx lay and ̄ x lay .The Voigt averages defined by Eqs. ( 8) and ( 9) subsequently determine the effective elastic tensor 4 Ccw and hygroscopic expan- sion tensor ̄ cw of the cell wall.In the third homogenization step, the effective properties of oak are computed at the macro-scale from asymptotic homogenization of the meso-scale response.The macro-scale thus is interpreted as the coarse scale and the meso-scale as the fine scale.As illustrated in Fig. 1b, the periodic meso-scale unit cell reflects the cellular structure of an oak growth ring, for which three representative mesostructural geometries of different densities are generated based on microscopy 1 3 Wood Science and Technology (2023) 57:1215-1256 images of oak.The cell walls are characterized by their effective properties 4 Ccw and ̄ cw computed in the second homogenization step, which refer to the local reference system of the cell wall.These local constitutive tensors are transformed to a global reference system associated with the meso-scale unit cell, and, as such, are used as input to solve the unit cell problems represented by Eqs. ( 4) and ( 5).The obtained influence functions are finally employed to compute the effective material properties 4 C and ̄ at the macro-scale via Eqs.( 6) and ( 7).Accordingly, with the current three-level homogenization procedure, the effective constitutive properties calculated at the macro-scale include meso-structural, micro-structural, and nano-structural information of oak wood.

Geometry
As illustrated in Fig. 3, the nano-scale unit cell is represented by a prismatic composite fiber with a length and a square cross section with edges w, for which the composite fiber consists of a cellulose micro-fibril core embedded in a matrix of hemicellulose and lignin, see also Malek and Gibson (2017).The nano-scale unit cell has periodic boundary conditions in the three directions of the Cartesian (1, 2, 3) coordinate system, which correspond to the principal material orientations Fig. 3 Nano-scale model: nano-scale unit cell of a prismatic composite fiber, consisting of a cellulose micro-fibril embedded in a matrix of hemicellulose and lignin.The cellulose micro-fibril comprises a crystalline core embedded in non-crystalline regions of the cellulose micro-fibril.The cellulose micro-fibril is composed of an inner crystalline region that is covered by a thin non-crystalline layer; the precise structure and distribution of these two constituent phases, however, are still a topic of debate in the wood science community (Hearle 1963;Roberts 1991;Terashima et al. 2009;Nishiyama 2009;Fernandes et al. 2011;Kulasinski et al. 2014).The crystalline regions are characterized by cellulose molecules with highly organized glucose chains linked together by glycosidic bonds.The non-crystalline cellulose is characterized by a lower chain ordering and a lower number of hydrogen bonds than crystalline cellulose, thereby reflecting an intermediate state between crystalline and amorphous cellulose (Khodayari et al. 2020a;Ioelovich et al. 2010).The crystalline and non-crystalline cellulose are indicated in Fig. 3.The inner crystalline cellulose core is defined by a length cc and a width w cc .The length of the surrounding, non- crystalline cellulose cover is nc = − cc , and its thickness is w nc .The outer het- erogeneous matrix, in which the hemicellulose and lignin are randomly distributed (Xu et al. 2007), has a length and a thickness (w − w cc − 2w nc )∕2 .Compared to crystalline cellulose, hemicellulose is based on carbohydrates with shorter branched chains and a lower degree of polymerization (Rowell 1984).Lignin is a complex phenolic polymer consisting of phenylpropane units linked by ether bonds (Mendez et al. 2019;Persson 2000;Qing and Mishnaevsky 2009).
The nano-scale unit cell is used to compute the effective hygro-elastic tensors 4 Clay and ̄ lay of the composite fibers of the cell wall layers (S1, S2, S3, and CML) at the micro-scale.Each cell wall layer is characterized by different volume fractions of the cellulose, hemicellulose, and lignin phases, which are respectively denoted as v c , v h and v lg .Accordingly, a different nano-scale geometry is generated for each specific cell wall layer, with their composition and geometrical parameters summarized in Table 1.The nano-scale geometry is constructed by first generating a unit cell consisting of a cellulose core that is embedded in a uniform matrix.The geometrical parameters of the crystalline cellulose region are kept fixed, with the width and length of the domain respectively taken as w cc = 3 nm and cc = 25 nm (Sved- ström et al. 2012).For the non-crystalline cellulose cover, a thickness w nc = 2 Å is selected (Andersson 2007).The length nc of the non-crystalline region is computed based on the crystallinity index (CI), which is the volume fraction of the crystalline region with respect to the total cellulose domain.In this work, the crystallinity index is taken as CI = 0.73 (Kránitz 2014), which results in a length of the non-crystalline region of nc = 1.66 nm.Correspondingly, the length of the nano-scale unit cell equals = cc + nc = 26.55nm.Next, for each nano-scale geometry associated with a specific cell wall layer, the volume fraction of the matrix is obtained as the sum of the volume fractions of the hemicellulose and lignin phases, according to the compositions listed in Table 1.From the specific geometry of the micro-fibril and the value of the matrix volume fraction, the cell width w for the S3, S2, S1, and CML layers is computed as w = 5.07, 4.81, 5.75, 19.73 nm, respectively.Based on the above dimensions, a nano-scale unit cell is generated and discretized by 17,496 finite elements, using linear eight-node brick elements equipped with a 2 × 2 × 2 Gauss integration scheme.Additional simulations not presented here have confirmed that this finite element discretization is sufficiently fine for providing meshindependent results.As mentioned, the matrix is initially considered homogeneous and composed of hemicellulose only, see also Malek and Gibson (2017).Subsequently, hemicellulose finite elements are selected randomly and replaced by lignin elements, where the total number of lignin elements is determined by the prescribed volume fraction v lg .

Hygro-elastic properties of the nano-scale constituents
The effective moisture-dependent elastic response of oak is governed by the capability of moisture uptake of the different constituent phases at the nano-scale, which is strongly influenced by their molecular structure (Korkmaz and Büyüksarı 2019).Consider first the crystalline and non-crystalline cellulose phases.Due to the densely packed, highly ordered glucose chains in crystalline cellulose, this constituent phase is poorly accessible to moisture, as a result of which moisture penetration remains limited to a small region near the outer surface of the cellulose phase (Bergenstråhle et al. 2008;Matthews et al. 2006).Molecular dynamics simulations have shown that the non-crystalline cellulose phase also has this feature (Khodayari et al. 2020b).Accordingly, in this work the crystalline and non-crystalline cellulose phases are considered to be hydrophobic (see also Cave 1978;Persson 2000;Mendez et al. 2019), by which their properties are independent of the moisture content.The crystalline cellulose is further assumed as a transversely isotropic material, with its plane of isotropy perpendicular to the micro-fibril axis (Persson 2000).
The elastic behavior of crystalline cellulose is thus characterized by five independent parameters, which are the Young's moduli of the cellulose micro-fibril in the transverse plane, E 1 (= E 2 ) , and in the longitudinal direction, E 3 , the Poisson's ratios 12 and 13 (= 23 ) , and the shear modulus G 13 (= G 23 ) .The values selected for the five independent elastic parameters are taken from Persson (2000), in correspondence with E 1 = 17.5 GPa, E 3 = 150.0GPa, 12 = 0.50 , 13 = 0.01 and G 13 = 4.5 GPa, see Table 2.Note that value listed for the shear modulus G 12 fol- lows from the relation G 12 = E 1 ∕(2(1 + 12 )) .The non-crystalline cellulose region commonly exhibits a somewhat stiffer response along the longitudinal micro-fibril  direction than in the transverse plane (Khodayari et al. 2020a, b).However, due to a lack of experimental data, non-crystalline cellulose is assumed to be isotropic, with the elastic properties equal to those of amorphous cellulose, i.e., a Young's modulus E = 6.9 GPa (Sundberg et al. 2013) and a Poisson's ratio = 0.15 (Kulasinski et al. 2015), see Table 2.
Hemicellulose is a hydrophilic phase with a relatively high capacity for moisture uptake.Its moisture content m h (defined by the mass of moisture per unit volume divided by the mass of dry material per unit volume) reaches a maximum percentage value of approximately 52% at the saturation point.Lignin has a lower capacity for moisture uptake, with its moisture content m lg reaching a maximum percentage value of 20% .Figure 4a shows the moisture adsorption isotherms of hemicellulose and lignin, which provide the relation between the equilibrium moisture content of each constituent phase, m ∈ {m h , m lg } , and the ambient relative humidity, RH .The adsorp- tion isotherms were measured over a relative humidity range between 0 and 97.4% and refer to native Pinus radiata hemicellulose (glucomannan) (Cousins 1978) and periodate lignin (Cousins 1976).As argued in Cousins (1976Cousins ( , 1977)), periodate lignin well approximates the behavior of native lignin.Although moisture equilibrium states generally are defined along both an adsorption isotherm and a desorption isotherm, the modeling of hysteresis effects, which requires a distinction into these two isotherms (Luimes and Suiker 2021), for simplicity is left out of consideration.Cousins (1978) and Cousins (1976), respectively; b experimentally measured elastic modulus of regenerated hemicellulose (relative to its value in the dry state) as a function of its moisture content m h (Cousins 1978); c experimentally measured elastic modulus of regenerated periodate lignin (relative to its value in the dry state) as a function of its moisture content m lg (Cousins 1977) The hygro-mechanical properties of the hydrophilic hemicellulose and lignin phases are quantified by considering a reference moisture content m = 12 % , with m ∈ {m h , m lg } .Note that this specific moisture content value is widely adopted in the literature for representing the macroscopic hygro-elastic properties of wood.For the sake of consistency, m = 12 % will be used as the reference moisture content for all four considered length scales shown in Fig. 1.The lignin phase has an amorphous structure for which an isotropic elastic behavior is assumed, defined by the Young's modulus E and Poisson's ratio .In Cousins (1976), these material parameters were measured as E = 3.1 GPa and = 0.26 at the reference moisture content m = m lg = 12% , and are adopted accordingly, see Table 2. Hemicellulose molecules tend to be aligned with the cellulose chains, resulting in a transversely isotropic material with the isotropic plane oriented perpendicular to the cellulose micro-fibril axis (Cave 1978).At the reference moisture content m = m h = 12% , the five inde- pendent elastic parameters are E 1 = 3.5 GPa, E 3 = 16.0GPa, 12 = 0.40 , 13 = 0.10 , and G 13 = 1.5 GPa (Persson 2000), as listed in Table 2.
From the elastic parameters presented in Table 2, for each constituent phase the elastic stiffness tensor 4 C can be constructed.For the hydrophilic hemicellulose and lignin phases, the elastic stiffness tensor relates to the reference moisture content m = 12% , and can be generalized toward arbitrary values of the moisture content by applying the experimental trends reported in Cousins (1978) for hemicellulose and in Cousins (1977) for lignin, as depicted in Fig. 4b, c, respectively.The figures illustrate the measured elastic moduli of hemicellulose and lignin-normalized by the corresponding values in the dry state-as a function of the moisture content.Except for the initial, increasing trend measured for lignin, the elastic moduli of the constituent phases decrease under an increasing moisture content.This decrease can be mainly attributed to hydrogen bond breakage generated during the hydration process, resulting in an increase in porosity and a mechanical weakening (Kulasinski et al. 2015).It is noted that these experimental results were obtained for regenerated samples with isotropic elastic properties; however, in the absence of more specific experimental data, these results will be used for determining the elastic stiffness tensors 4 C of both hemicellulose (which is considered transversely isotropic) and lignin.The elastic stiffness tensor can be determined at an arbitrary value of the moisture content, by multiplying the stiffness tensor corresponding to reference moisture content m = 12% by the ratio between the value of the curve in Fig. 4b or c at an arbitrary, alternative moisture content value and the value at m = 12%.
The hygroscopic properties of the hydrophilic constituents of oak wood at the nano-scale are considered to be independent of the moisture content.The hygroexpansion coefficients summarized in Table 2 refer to a percentage change in strain per percentage change in moisture content, and are based on the values reported in Cave (1978), which are regularly used in the literature (Persson 2000;Mendez et al. 2019).For hemicellulose, a transversely isotropic hygro-expansive response is considered, where the minor hygro-expansive deformation along the direction of the micro-fibril axis typically is ignored (Cave 1978;Persson 2000).In the other two directions, the values of the hygro-expansion coefficients are taken equal to 1 = 2 = 1.37 (Persson 2000).Finally, the response of lignin is considered Wood Science and Technology (2023) 57:1215-1256 isotropic, with the hygro-expansive coefficients taken as 1 = 2 = 3 = 0.35 (Pers- son 2000).
As mentioned in Sect."Asymptotic homogenization procedure", in the asymptotic homogenization procedure the nano-scale elastic tensor 4 C of the hydrophilic hemicellulose and lignin phases needs to be specified in terms of the relative humidity.This is, since the relative humidity is uniform at the nano-scale (and the scales above), while the moisture content in the fine-scale configuration depends on the moisture adsorption capacity of the constituent phases, and thus is non-uniform at the nano-scale.The elastic tensors of the hemicellulose and lignin phases can be calculated in terms of the relative humidity as follows.For a specific value of the relative humidity, the moisture contents of the hemicellulose and lignin phases follow from the moisture adsorption isotherms in Fig. 4a.The elastic stiffness tensors 4 C corresponding to these moisture contents, which thus also apply to the specific relative humidity value, can be subsequently determined by using the curves in Fig. 4b, c, in accordance with the procedure explained earlier in this chapter.Note that for the nano-scale hygroscopic expansion tensor this conversion is not necessary, since the hygro-expansion coefficients listed in Table 2 are independent of the moisture content.

Geometry
At the micro-scale level of oak wood, a cell wall is considered as a composite laminate consisting of two secondary walls (which include the S1, S2, and S3 layers) and a compound middle lamella (CML) layer in between.As mentioned in Sect."Multi-scale geometry", the CML layer represents the effective behavior of a middle lamella layer connected to a primary layer at each side.Correspondingly, the cell wall is represented by the seven-layer laminate depicted in Fig. 5. Generally, cell walls composed of libriform fibers, fiber tracheids, parenchyma, and rays may have different characteristics.In the present work, the so-called axial cells, which refer to the libriform fibers, fiber tracheids and parenchyma that are (mainly) oriented along the axial (or longitudinal) material direction of oak wood, for simplicity are considered to have the same cell wall structure, while for the ray cells, or rays, which are oriented along the radial material direction, a different cell wall structure is assumed.The geometrical parameters that distinguish the axial cells from the ray cells are the thicknesses of the secondary and CML cell wall layers, and the orientations of the corresponding micro-fibril angles (MFAs).The MFA is defined as the angle between the micro-fibrils' longitudinal axis and the (vertical) z-axis of the global (right-handed) x-y-z coordinate system depicted in Fig. 5, as measured about the (out-of-plane) y-axis.The micro-fibrils' longitudinal axis is indicated by the 3-axis of the local 1-2-3 coordinate system, while the 1-axis denotes its transverse axis in the plane of the cell wall, and the 2-axis denotes its transverse axis in the out-of-plane direction.Note that the local 1-2-3 coordinate system corresponds to that of the nano-scale unit cell shown in Fig. 3. Table 3 summarizes the values of the relative layer thickness, t lay (defined relative to the total thickness of the composite cell wall), and the orientations of the MFAs adopted in the axial and ray cell wall models.The thickness of the cell wall layers is based on the average thickness of earlywood and latewood cell walls measured from the SEM images of axial and ray cells in Pedunculate oak, as presented in Kim and Daniel (2016) and Harada and Wardrop (1960), respectively.As mentioned above, the cell wall layers of axial cells and rays are characterized by different values of the MFA, where for both cell types the MFAs of two opposing secondary layers of the same type are oriented mirror-symmetrically with respect to the vertical y-z plane, see Fig. 5.For the axial cells, the micro-fibril angles in the S3 and S1 layers are selected as 75 • and −60 • , respectively (Mendez et al. 2019).The MFA in the S2 layer of the axial cells may vary considerably, from a value of approximately 0 •  for libriform fibers (Lichtenegger et al. 1999), to 20 • for fiber tracheids (Lichteneg- ger et al. 1999), and 40 • for parenchyma cells (Chafe and Chauret 1974).Since the MFA in the S2 layer of axial cells strongly influences the effective response of wood (Cave 1978), two different MFAs will be considered in the numerical analyses in Sect."Macro-scale results", namely 10 • and 20 • .Note that these values fall within the range of values mentioned above.For ray cells, in the S3 layer the micro-fibril angle has been observed to vary within a range of 30 • − 60 • (Harada and Wardrop 1960;Chafe and Chauret 1974).The MFA in the S1 layer of a ray cell typically has the opposite value of the MFA in the S3 layer.In the S2 layer, the MFA value varies between 0 • (Harada and Wardrop 1960) and 39 • (Chafe and Chauret 1974).As a working assumption, the MFA in each secondary wall layer of ray cells is taken as the average value of the corresponding range of literature values mentioned above.Accordingly, the MFAs for the S1, S2, and S3 layers are equal to 45 • , 20 • , and −45 • , respectively.Finally, the MFA of the CML layer is assumed to be randomly distributed for both axial cells and ray cells (Chafe and Chauret 1974).

Hygro-elastic properties of micro-scale cell wall layers
The effective elastic stiffness tensor 4 Clay and hygro-expansion tensor ̄ lay of each micro-scale cell wall layer are obtained from asymptotic homogenization of the material response at the nano-scale using Eqs.( 6) and ( 7).As discussed in Sect."Hygro-elastic properties of the nano-scale constituents" above, the calculated effective properties depend on the relative humidity, which can be assumed to be a uniform field.At the micro-scale, each cell wall layer is considered homogeneous and thus is characterized by a specific, constant moisture content value that differs per cell wall layer.Therefore, in the following, the micro-scale elastic and hygroscopic properties of the cell wall layers will be presented as a function of their specific moisture content.For this purpose, the effective moisture content mlay in a micro-scale cell wall layer is computed as the sum of the moisture contents of the nano-scale hydrophilic phases weighted by their volume fraction: Here, v h and v lg , respectively, are the volume fractions of hemicellulose and lignin in each specific cell wall layer, as listed in Table 1.Further, m h and m lg represent the moisture contents in hemicellulose and lignin, respectively, which can be obtained from the moisture adsorption isotherms depicted in Fig. 4a using the value of the relative humidity as input.With this procedure, from Eq. ( 10) the effective moisture adsorption isotherm of each micro-scale cell wall layer can be constructed over the complete range of relative humidities, see Fig. 6.Observe from Fig. 6 that the adsorption characteristics are different for the individual cell wall layers, as a result of their typical dependency on the underlying nano-scale properties.In specific, the CML layer has the largest capacity for moisture uptake, which is caused by the large volume fractions of the hydrophilic hemicellulose and lignin phases, see Table 1.
Conversely, the S2 layer has the smallest moisture uptake capacity, due to the relatively large volume fraction of the hydrophobic cellulose phase.Table 4 summarizes the effective micro-scale hygro-elastic parameters of the composite fibers-as represented by the nano-scale unit cell in Fig. 3-in each of the cell wall layers at the effective reference moisture content mlay = 12% .The effective fiber properties for each layer relate to the nano-scale composition and cell width listed in Table 1, and are expressed in the local Cartesian coordinate system defined along the micro-fibril principal material directions (1, 2, 3), as illustrated in Fig. 3.Note that the effective fiber properties for the secondary cell wall layers refer to a specific orientation of the micro-fibril angle, see Table 4, while for the CML layer, they will be used to compute the cell wall layer properties considering randomly oriented micro-fibrils, as explained below in more detail.In correspondence with the nano-scale constitutive behavior and geometries adopted for the crystalline cellulose and hemicellulose phases, the effective fiber properties for the secondary cell wall layers appear to be approximately transversely isotropic, as also reported in Mendez et al. (2019).In the CML layer, lignin is the predominant component, which causes the elastic modulus Ē3 in the longitudinal (3-) direction of the micro-fibril to be much lower than for the secondary layers.Further, for all layers the effective hygro-expansion coefficient β3 in the longitudinal (3-) direction of the micro-fibril is considerably lower than the coefficient β1 (= β2 ) in the transverse (1,2) plane.
Figure 7 illustrates the dependence of the micro-scale hygro-elastic properties (normalized with respect to their values in the dry state) of the composite fibers in the S2 layer on the micro-scale moisture content mS2 .For the other cell wall lay- ers similar trends are obtained, so these results have been omitted here for brevity.The material response is approximately transversely isotropic and thus is characterized by the values of the five independent elastic parameters Ē1 , Ē3 , ν13 , ν12 and Ḡ13 .For completeness, the shear modulus Ḡ12 = Ē1 ∕(2(1 + ν12 )) has also been depicted.Observe that for all hygro-elastic properties an increase in moisture content globally corresponds to a decrease in the hygro-elastic response of the cell wall layers, which results from the moisture-dependent constitutive behavior at the nano-scale level, see Fig. 4b, c.
The hygro-elastic properties of the cell wall layers will be used as input for the homogenization procedure between the micro-and meso-scales, see Sect."Meso-scale model" below.For this purpose, the properties first need to be expressed in the global reference system x = (x, y, z) of the micro-scale unit cell, depicted in Fig. 5.For the S1, S2, and S3 layers, the constitutive tensors in the local (1-2-3) coordinate system are transformed to the global (x-y-z) reference system by applying the usual operations for a coordinate transformation of second-order and fourth-order tensors.Hence, in the global reference system the stiffness tensor 4 Cx lay and the hygroscopic expansion tensor ̄ x lay of the cell wall layers are obtained as Here, the second-order tensor R represents the standard rotation tensor.Its elements are defined by the cosines and sines of the orientation angles MFA of the microfibrils, as listed in Table 3.The rotation tensor meets the requirements R T = R -1 and det R = 1 , where the superscripts "T" and " −1 " denote the transpose and inverse of a tensor, respectively, and " det " indicates the determinant.In the CML layer, the orientations of the micro-fibrils are assumed to be randomly distributed.Following a similar approach as in Derome et al. (2013) and Ogierman and Kokot (2017), the elastic stiffness tensor 4 Cx lay and hygro-expansion tensor ̄ x lay in the global coordinate system are calculated by averaging the hygro-elastic tensors along N discrete orientations of the MFA.This is done in accordance with the expressions where the rotation matrix R MFA i is computed for each specific angle of the micro-fibril, MFA i = 2i ∕N , with i ∈ {1, 2, ..., N} .Similar as in Derome et al.  (2013) and Ogierman and Kokot (2017), the total number of micro-fibril orientations N is assumed to be equal to 36.

Image preparation
The meso-scale model is based on high-resolution microscopy images of an oak wood sample.For obtaining the images, a cubic oak sample was cut from a larger specimen along the three anatomical (L, T, R) planes, which refer to the longitudinal (L), transverse (T), and radial (R) material directions of wood.Using a disposable blade and a sliding Reichert microtome, three sections (i.e., one per principal material direction) were cut and subsequently prepared into high-quality permanent slides containing several growth rings.The surface of the slides was equal to 14 × 8.5 mm 2 , and the thickness ranged between 15 and 25 μ m.A solution of 1 g of Safranin powder in 100 ml of distilled water was prepared for staining the thin slides.After staining, the slides were washed with water to remove the residual stain, until the running water was clear.Next, the slides were dehydrated with an ethanol concentration that was gradually increased in a stepwise fashion from 50% to 75% , and finally to 96% , in order to prevent damage of cell walls under rapid shrinkage (Yeung et al. 2015).The slides were embedded in Canada balsam under a cover glass and were dried for a period of 12 to 24 hours at a temperature of 60 • C (Yeung et al. 2015).Subsequently, the permanent slides were scanned using an automated digital microscope 100× to create high-quality images of the complete slides.
Figure 8a illustrates the transverse (R-T) section of the oak sample.The image reveals the complex cellular structure of oak, in which libriform fibers, fiber tracheids, axial parenchyma, vessels, and rays are present.Further, three different growth rings A, B, and C with width W and height H have been selected for the construction of meso-scale unit cells.The variations in morphological features within the growth ring are clearly visible.Specifically, a pattern of large earlywood vessels can be identified across the width W e of the earlywood region of growth ring B. In contrast, much smaller vessels are present across the width W − W e of the latewood region, resulting in a relatively dense material structure.Moreover, parallel regions of multi-seriate and uni-seriate rays between the axial cells can be detected.respectively.Here, the large voids with an elongated oval shape correspond to vessel sections.Additionally, the relatively dark regions indicate fiber sections, while the light regions correspond to parenchyma cells.It can be further observed that the cross section of rays in the R-T plane is substantially different from that in the T-L plane, in correspondence with respective cuts made along the longitudinal and transverse axes of the ray cells.The resolution of the original images is very high so that cell walls can be accurately identified when zooming in at higher magnification.The high resolution enables an accurate image processing required for the construction of three-dimensional finite element configurations of oak growth rings at the mesoscale, as described in Sect."Image-based meso-scale unit cell" below.

Image-based meso-scale unit cell
The meso-scale unit cells are generated based on the microscopy images presented in Fig. 8. Accordingly, three regions of dimensions W × H are selected in the R-T plane of the meso-scale material structure, as indicated by the black rectangular boxes A, B, and C in Fig. 8a.Note that the three regions correspond to three distinct growth rings, where each growth ring is characterized by a specific proportion and distribution of the various cell types, and thus has a specific local density, denoted by A , B , and C .Under the assumption of a periodic meso-structure, see Fig. 2, the three selected regions may be considered as representative of three distinct oak samples that are characterized by the overall densities A , B , and C .Each of the three two-dimensional regions extracted from the R-T plane image is next converted into a binary image with an appropriate pixel intensity threshold.The threshold is fine-tuned such that the black pixels of the binary image represent cell wall material, as shown in Fig. 9a.An extra homogeneous layer with a thickness of three black pixels is added at the rectangular outer boundary of each of the three regions to impose geometrical periodicity on the domain.Subsequently, the cell wall edges are identified using the image processing toolbox for boundary detection of MATLAB, which contains an edge detection scheme based on localizing the variations of graylevel images.Since the resolution of the obtained images is insufficient to provide smooth edges for the cell walls, jagged-shaped edges for the black pixel regions are retrieved instead, as illustrated in Fig. 9b.From the jagged geometry, the smooth edges of the cell walls are reconstructed using a second-order B-spline interpolation function, for which the corner points along the interface between the black and white pixels are selected as the control points of the B-spline function, as indicated by the solid circles in Fig. 9b.The smooth cell edges identified with this procedure are provided as input for the generation of a finite element mesh.Figure 9c shows a detail of the two-dimensional finite element mesh based on triangular elements.
From the discretized two-dimensional domain, the three-dimensional FEM mesh can be subsequently constructed.In the literature, three-dimensional FEM meshes for wood are generally established through the extrusion of the two-dimensional mesh in the out-of-plane direction (Malek and Gibson 2017;Mendez et al. 2019).However, when extruding the current 2D image in the (longitudinal) out-of-plane direction, the section of rays becomes fully "open" along the longitudinal direction, see the local geometry within the dashed rectangle shown in Fig. 9e, which may lead to an underprediction of the effective, macro-scale Young's moduli, especially the effective Young's modulus ĒT in the tangential direction (Malek and Gibson 2017).
To minimize this discrepancy, the following strategy is proposed for generating the three-dimensional meso-scale domain.First, the two-dimensional image is extruded along the longitudinal (L) direction, see Fig. 9e, where the thickness of the extruded volume is set equal to the average thickness of rays in the longitudinal direction, as determined from the image in the T-L plane presented in Fig. 8b.Next, in the regions containing rays, the domain is closed in the longitudinal direction by adding elements with ray cell wall properties in the two outer element layers, see Fig. 9f.The three-dimensional unit cell obtained is depicted in Fig. 9g, and thus reflects the presence of single ray cells in the longitudinal direction.Upon extrusion, the triangular finite elements initially defining the 2D mesh in the R-T plane are replaced by linear, six-node wedge elements equipped with a two-point Gauss integration scheme to construct the 3D mesh.Depending on the material density, the 3D finite element configurations contain between 3.2 × 10 6 and 4.0 × 10 6 elements.A mesh- refinement study not presented here has shown that the constructed meshes are sufficiently fine for providing mesh-independent simulation results.Table 5 summarizes the characteristics of the generated three-dimensional meso-scale unit cells.The width W and height H of the unit cell and the width W e of the earlywood region are measured from the microscopy image in the transverse plane, as indicated in Fig. 8a for sample B. Further, the porosity , the ray volume fraction v ray and the vessel volume fraction v vess,EW in the earlywood (EW) of samples A, B, and C are determined by measuring the corresponding regions in the discretized three-dimensional meso-scale unit cell.The porosity values range between [0.45 − 0.60] ; assuming a cell wall density in the dry state of cw = 1440 kg/m 3 (Kellogg and Wangaard 1969), the densities = (1 − ) cw of the generated meso-scale unit cells follow as A = 582 kg/m 3 , B = 685 kg/m 3 , and C = 798 kg/m 3 .These densities are in good agreement with the den- sities reported in the literature for Pedunculate oak, which vary in the range [461 − 874] kg/m 3 (Badel and Perré 2006).The solid volume fraction of ray cells is defined as v ray = V ray ∕(V ax + V ray ) , with V ray and V ax being the volumes of axial cells and ray cells.From the value of v ray , the solid volume fraction of axial cells,

Hygro-elastic properties of meso-scale cell walls
The effective mesoscopic hygro-elastic properties of the cell walls follow from a Voigt average of the material characteristics at the micro-scale.Specifically, the cell wall layer properties 4 Cx lay and ̄ x lay , which are expressed in the global reference system of the cell wall via Eqs.( 11)-( 14), are used as input to compute the overall cell wall stiffness tensor 4 Ccw and hygro-expansion tensor ̄ cw by using Eqs.( 8)-( 9).In addition, the effective mesoscopic moisture content mcw in cell walls composed of axial cells ( mcw = max ) and ray cells ( mcw = mray ) is calculated as the weighted vol- ume average of the moisture contents in the micro-scale cell wall layers, i.e., Here, t lay = t S1 , t S2 , t S3 , t CML are the relative thicknesses of the cell wall layers, as listed in Table 3 for the axial cells and ray cells, and mlay = mS1 , mS2 , mS3 , mCML are the moisture contents in the cell wall layers, which, for an arbitrary value of the relative humidity, can be determined from the micro-scale moisture adsorption isotherms depicted in Fig. 6.By performing this procedure over the full relative humidity range, the effective meso-scale adsorption isotherms of axial cell walls and ray cell walls can be constructed, see Fig. 10.Observe that the axial cell walls and ray (15 cell walls are defined by comparable moisture adsorption isotherms, with the ray cell walls having a slightly higher capacity for moisture uptake. Table 6 summarizes the effective hygro-elastic properties of the cell walls composed of axial cells and ray cells, referring to an effective reference moisture content mcw = 12% .The results correspond to an MFA of 10 • in the S2 layer of axial cells; the specific influence of this MFA on the macro-scale properties of oak will be investigated in more detail in Sect."Macro-scale results".In correspondence with the geometry selected for the micro-scale unit cells, the properties of the cell wall layers are virtually orthotropic, with the stiffer direction oriented along the (vertical) z-axis of the cell wall, see Fig. 5.This behavior is consistent with the elastic response of cell walls determined from nano-indentation experiments (Konnerth et al. 2010).Further, for the axial cell wall the Young's modulus in z-direction is almost two times larger than that of the ray cell wall, which is mainly due to a thicker (relatively stiff) S2 layer in the axial cell wall, see Table 3.Additionally, the hygro-expansion coefficient is the largest in the direction perpendicular to the cell wall, i.e., the y-direction, and is negligible in the z-direction.The latter result has also been reported by other investigators (Rafsanjani et al. 2014).
Figure 11 illustrates the meso-scale hygro-elastic properties (normalized with respect to the value in the dry state) as a function of the moisture content.Since for the axial cell wall and ray cell wall similar trends are found, only the results for the axial cell walls are depicted here.Apart from a slight initial increase in some of the curves that is caused by the moisture-dependent behavior of lignin, see Fig. 4c, the coefficients typically decrease under an increasing moisture content.
The meso-scale hygro-elastic properties of the cell walls serve as input for the asymptotic homogenization procedure that provides the macro-scale properties.Accordingly, the meso-scale properties listed in Table 6, which refer to the principal material directions (x, y, z) of a cell wall, see Fig. 1c, need to be converted to the global reference system (R, T, L) used at the macro-scale, see Fig. 1a.Since for the axial cells the vertical z-axis of the cell wall coincides with the longitudinal (L) direction at the macro-scale, only the x-and y-directions in the local material points of the cell walls need to be identified in more detail, after which they can be expressed in terms of the R-and T-directions of the global reference system by applying a standard coordinate transformation.The identification of the cell wall x-and y-directions is performed for all triangular elements used for constructing the two-dimensional mesh for the cell wall configuration in the R-T plane, see Fig. 9d.For each triangular element, the cell wall edge point located closest to the element center is identified, from which the x-direction of the cell wall follows from the corresponding derivative of the B-spline function parameterizing the cell wall edge in this edge point.Trivially, the y-direction is oriented perpendicular to the x-direction.For the ray cells, a similar procedure is followed.

Moisture adsorption isotherm
At the macro-scale, the effective moisture content m is computed as the weighted volume average of the meso-scale moisture contents max and mray of the cell walls composed of axial and radial cells: Here, max and mray follow from Eq. ( 15), and v ray and v ax = (1 − v ray ) are the solid volume fractions of ray cells and axial cells in the considered growth rings A, B, or C displayed in Fig. 8a, with the corresponding values of v ray listed in Table 5.For an arbitrary value of the relative humidity, from the meso-scale moisture adsorption isotherms depicted in Fig. 10 the meso-scale moisture contents max and mray can be determined, which thus result in the macro-scale moisture content m via Eq.( 16).Performing this procedure over the full range of relative humidities leads to the macro-scale moisture adsorption isotherm depicted in Fig. 12.The macro-scale moisture adsorption isotherm is presented together with experimental data for oak reported in Luimes and Suiker (2021), showing that the results obtained by the model prediction and the experimental measurements are in reasonably good agreement.The discrepancy in the 50-90% relative humidity range is probably due to the fact that the oak wood tested in Luimes and Suiker 2021) may have somewhat different characteristics at the small scale than the oak sample analyzed in the present communication.

Effective hygro-elastic properties of oak
The macro-scale hygro-elastic response of oak is obtained from asymptotic homogenization of the material response at meso-scale level.The meso-scale properties presented in Sect."Hygro-elastic properties of meso-scale cell walls" depend on the relative humidity, and are used as input for computing the macro-scale stiffness tensor 4 C and the hygroscopic expansion tensor ̄ via Eqs.( 6) and ( 7), respectively.The macro-scale properties are thus obtained as a function of the relative humidity and are subsequently expressed in terms of the moisture content by using the macroscopic moisture adsorption isotherm depicted in Fig. 12.

Effective elastic properties of oak
The macroscopic elastic response of oak following from the three-level homogenization procedure is approximately orthotropic and thus is represented by the nine (16) m = v ax max + v ray mray .16), using the meso-scale moisture adsorption isotherms in Fig. 10 as input.The experimental data have been taken from Luimes and Suiker (2021) independent elastic parameters shown in Figs. 13, 14, and 15.In these figures, solid, dashed, and dotted lines refer to the meso-scale unit cells A, B, and C, which are respectively characterized by a relatively low density ( A = 582 kg/m 3 ), a medium density ( B = 685 kg/m 3 ), and a high density ( C = 798 kg/m 3 ), see also Table 5.Further, the results are depicted for MFAs of 10 • and 20 • in the S2 layer of axial cells.
Figure 13 illustrates the effective Young's moduli of oak as a function of the macro-scale moisture content m , together with the experimental data reported in Ozyhar et al. (2016), Güntekın et al. (2016a), Luimes et al. (2018b) and Güntekın et al. (2016b).The Young's moduli in the radial direction ĒR , tangential direc- tion ĒT , and longitudinal direction ĒL are presented in Fig. 13(a), Fig. 13(b), and Fig. 13(c), respectively, and clearly show that a higher material density generally leads to a higher stiffness value.Specifically, the Young's modulus ĒL increases almost proportionally with the material density , as also reported in other modeling works on hardwood (Malek and Gibson 2017).Additionally, all three Young's moduli show a slightly increasing trend at relatively small moisture content, followed by a monotonic decrease for higher moisture content values.The initial increase in the curves is caused by the moisture-dependent behavior of lignin at the nano-scale, see Fig. 4(c).Despite that the experimental data in Fig. 13 represent some scatter, they lie close to the model predictions and show a comparable decreasing trend under a growing moisture content.It can be further observed that the Young's modulus ĒL in the longitudinal direction has the largest value, and is less influenced by moisture content variations than the Young's moduli ĒR and ĒT in the radial and tangential directions.This feature can be ascribed to the predominantly longitudinal orientation of the stiff, hydrophobic cellulose micro-fibrils in the S2 layer of axial cells, i.e., MFA = 10 • , 20 • , see Table 3.Note hereby that an increase in the MFA orientation from 10 • to 20 • in the S2 layer of axial cells indeed leads to a substantial drop in the Young's modulus ĒL in the longitudinal direction, and has a negligible influ- ence on the Young's moduli ĒT and ĒR in the transverse and radial directions.Fig- ures 13(a) and 13(b) also illustrate that the radial Young's modulus ĒR generally has a higher value than the tangential Young's modulus ĒT .This feature can be ascribed to the prevalent alignment and the denser configuration of fiber cells in the radial direction compared to the transverse direction and the stiffening effect of the rays in the radial direction (Perré and Badel 2003).At the highest moisture content considered, m = 21% , the ratio ĒR ∕ ĒT ranges between 2.16 and 2.32, with the exact value depending on the material density.This result is comparable to the experimental value ĒR ∕ ĒT ≈ 2 reported in Burgert et al. (2001) for a fully saturated Pedunculate oak wood sample with a dry density of 530 kg/m 3 .
Figure 14 shows the effective Poisson's ratios of oak as a function of the macroscale moisture content m , together with the experimental data reported in Ozy- har et al. (2016), Güntekın et al. (2016a), Luimes et al. (2018b) and Güntekın et al. (2016b).The influence of the material density on the Poisson's ratios νLT (Fig. 14(a), tangential plane), and νRT (Fig. 14(c), transverse plane) is mild, and on the Poisson's ratio νLR (Fig. 14(b), radial plane) is negligible.Further, in correspond- ence with the experimental results presented in Mendez et al. (2019), the Poisson's ratios only slightly increase under an increasing moisture content, with the exception being the somewhat more substantial increase for νLT in case the S2 layer of axial cells is characterized by the larger micro-fiber angle of 20 • , see Fig. 14(a).The effect on the Poisson's ratio by a change in MFA from 10 • to 20 • is the largest for νLT , fol- lowed by νLR and νRT .Despite their spread, the experimental data for the Poisson's ratios show a consistent, approximately constant trend under a varying moisture content, for which the values in most cases are comparable to the modeling results obtained for an MFA of 20 • in S2 layer of axial cells.
Figure 15 finally shows the effective shear moduli of oak as a function of the macro-scale moisture content m , together with experimental data taken from the lit- erature (Ozyhar et al. 2016;Güntekın et al. 2016b;Luimes et al. 2018b).The shear moduli associated with the tangential plane ḠLT , the radial plane ḠLR , and the trans- verse plane ḠRT are presented in Fig. 15(a), Fig. 15(b), and Fig. 15(c), respectively.Similar to the Young's moduli, the shear moduli increase with increasing material density .Further, all shear moduli show a comparable trend under a growing moisture content, which is characterized by a small, initial increase, followed by a substantial, monotonic decrease at larger moisture content.Additionally, the MFA orientation in the S2 layer of axial cells does not influence the effective shear modulus ḠRT in the transverse plane, while the moduli ḠLT and ḠLR in the tangential and radial planes become larger for a larger micro-fiber angle.This result has also been reported in Astley et al. (1998).Observe further that the (limited) experimental data points lie relatively close to the modeling results, especially to those related to the larger micro-fiber angle of 20 • .

Effective hygro-expansion coefficients
The effective hygroscopic coefficients of oak following from the three-level homogenization procedure are depicted in Fig. 16 as a function of the macro-scale moisture content m , where the solid, dashed, and dotted lines refer to the meso-scale unit cells A, B, and C with densities A = 582 kg/m 3 , B = 685 kg/m 3 , and C = 798 kg/m 3 , respectively.The results correspond to an MFA of 10 • and 20 • in the S2 layer of axial cells.The hygroscopic behavior of oak is approximately orthotropic and thus is essentially defined by the three coefficients βR , βT and βL depicted in Fig. 16(a), Fig. 16(b)  and Fig. 16(c), respectively.These coefficients represent the diagonal elements of the hygroscopic expansion tensor; it has been indeed confirmed that the computed offdiagonal elements are negligible.The experimental results presented in Fig. 16 are based on the test data presented in Ozyhar et al. (2016), Badel and Perré (2006) and Kránitz (2014); since the values of the moisture content at which the hygro-expansion coefficients were measured are not reported in these references, only the spread in experimental data has been displayed by means of a gray region of constant bandwidth.Note that the experimental bandwidths for the three hygroscopic coefficients indeed mostly fall within the ranges of the modeling results.It is further seen that the hygroexpansion coefficient is hardly affected by both the material density and the MFA, with the latter result being confirmed by the experimental-numerical study presented in Yamamoto et al. (2001).In addition, the hygro-expansion coefficients typically decrease as a function of the moisture content; only for the longitudinal component βL a slight initial increase caused by the typical moisture-dependent behavior of lignin can be observed.Note further that the hygro-expansion coefficient in the transverse direction, βT , is larger than the hygro-expansion coefficient in the radial direction, βR , for which, depending on the material density, the ratio βT / βR lies between 1.69 to 1.90.This result corresponds well with experimental results presented in Badel and Perré (2006) for oak samples with an average density of 700 kg/m 3 , for which βT / βR ≈ 1.78 .Finally, the hygro-expansion coefficient in the longitudinal direction is much smaller than in the two other directions, which is a typical characteristic of (oak) wood (Ozyhar et al. 2016;Kránitz 2014).The spread of the experimental data, which have been taken from Ozyhar et al. (2016), Badel and Perré (2006) and Kránitz (2014), is indicated by the gray region 1 3 Wood Science and Technology (2023) 57:1215-1256

Conclusion
A three-dimensional, multi-scale model is presented that can be used for the analysis of the moisture-dependent hygro-elastic response of oak wood.The model considers four different length scales, which correspond to material structures at the nano-, micro-, meso-, and macro-scales.From these scales, the effective, hygro-mechanical behavior of oak at the macro-scale is determined using a threelevel homogenization procedure, which is applied in three subsequent steps.In the first homogenization step, for each oak cell wall layer a nano-scale unit cell is defined that describes the hygro-mechanical behavior of a three-dimensional cellulose micro-fibril embedded in a matrix of hemicellulose and lignin.The unit cell is subjected to an asymptotic homogenization procedure to compute the hygromechanical properties of individual cell wall layers.The cell wall layers appear to be characterized by a transversely isotropic hygro-elastic constitutive behavior, for which the effective material properties generally decrease in value under an increasing moisture content.The effective moisture adsorption isotherms of the micro-scale cell wall layers are constructed from the moisture absorption characteristics of the nano-scale constituents by using a weighted volume averaging scheme.The computational results of this first homogenization step are used as input for the second homogenization step, in which oak cell walls are idealized at the micro-scale as three-dimensional composite plates consisting of a stack of cell wall layers.The effective meso-scale hygro-elastic properties are calculated by applying a Voigt averaging procedure, and the moisture adsorption isotherms of the cell walls are computed from a weighted volume average of the microscale moisture adsorption characteristics.The resulting cell wall properties are approximately orthotropic, which is consistent with literature results.In the third homogenization step, accurate three-dimensional meso-scale unit cells describing growth rings of different densities are generated from detailed microscopy images of oak samples.With the results of the second homogenization step used as input, the application of asymptotic homogenization to the meso-scale models ultimately provides the hygro-elastic material properties of oak at the macroscale.Further, a weighted volume average of the meso-scale moisture adsorption behavior results in the macro-scale moisture adsorption isotherm.
The numerical results illustrate that the proposed multi-scale model realistically predicts the effective moisture-dependent hygro-elastic behavior of oak wood at the macro-scale.The constitutive behavior appears to be approximately orthotropic and corresponds well with experimental data taken from the literature.The elastic and hygro-expansion coefficients globally decrease for an increasing moisture content, and the elastic behavior becomes stiffer for a higher material density.Further, the micro-fibril orientation in the S2 layer of axial cells appears to only noticeably influence the effective properties associated with the longitudinal material direction, i.e., the Young's modulus in the longitudinal direction and the Poisson's ratios and shear moduli in the radial and tangential planes.
The effective mechanical and hygroscopic properties calculated with the current model can be used as input for numerical analyses of engineering and cultural heritage applications aimed at predicting the moisture-dependent mechanical response of oak structures under complex hygro-mechanical loading paths.To further validate the computational results presented in this paper, more systematic experimental studies are needed at the different length scales considered in the homogenization approach.In addition, by using micro-tomography scanning techniques, the three-dimensional meso-scale unit cell configurations defined in this study can be further refined in order to accurately represent the complex meso-scale morphology, especially along the longitudinal material direction.

Fig. 1
Fig. 1 Model representation of oak wood involving four different length scales: a macroscopic scale with a characteristic length L , corresponding to a homogeneous oak wood sample; b mesoscopic scale with a characteristic length L, reflecting the cellular structure of an entire growth ring; c microscopic scale with a characteristic length l, representing a heterogeneous cell wall layer; d nanoscopic scale with a characteristic length , reflecting a cellulose micro-fibril embedded in a hemicellulose and lignin matrix Hygro-elastic properties of the nano-scale phases, expressed in the local 1-2-3 coordinate system.The elastic parameters of the hydrophilic hemicellulose and lignin phases correspond to a reference moisture content of m = 12% Constituent phase Young's modulus (GPa) Poisson's ratio (-) Shear modulus (GPa) Wood Science and Technology (2023) 57:1215-1256

Fig. 4
Fig. 4 Nano-scale model: a moisture adsorption isotherm of hemicellulose ( m = m h ) and lignin ( m = m lg ), taken fromCousins (1978) andCousins (1976), respectively; b experimentally measured elastic modulus of regenerated hemicellulose (relative to its value in the dry state) as a function of its moisture content m h(Cousins 1978); c experimentally measured elastic modulus of regenerated periodate lignin (relative to its value in the dry state) as a function of its moisture content m lg(Cousins 1977)

Fig. 5
Fig. 5 Micro-scale model: The micro-scale unit cell consists of two stacks of secondary cell wall layers (S1, S2, and S3) and a compound middle lamella (CML) layer in between.The micro-fibril angles (MFAs) of the individual cell wall layers are designated by black solid lines and are representative of axial cells, in agreement with the specific values listed in Table3for the front and back layers of the depicted cell wall.The MFAs in two opposing (front and back) secondary layers of the same type are oriented mirror-symmetrically with respect to the vertical y-z plane

Fig. 7
Fig. 7 Micro-scale model: Effective hygro-elastic properties of the S2 layer (normalized with respect to the value in the dry state) as a function of the moisture content.a Young's moduli; b Poisson's ratios; c shear moduli; d hygro-expansion coefficients Figure 8b, c show the radial (R-L) and tangential (T-L) sections of the oak sample,

Fig. 8
Fig. 8 Meso-scale model: microscopy images of oak sections taken along the three principal material directions.a Section along the transverse (R-T) plane.The cellular structure of the growth rings reveals a pattern of large earlywood vessels and much smaller latewood vessels.The three rectangular regions of width W and height H indicated by the black boxes A, B, and C refer to different growth rings and are used for the construction of the meso-scale unit cells.For growth ring B, the transition between the earlywood and latewood regions is indicated by a dashed line, with their widths equal to W e and W − W e , respectively.b Section along the radial (R-L) plane.The dark and light regions indicate fibers and light parenchyma cells, respectively.The elongated oval-shaped voids represent the vessels.c Section along the tangential (T-L) plane.The specific cross-sections of multi-seriate and uni-seriate rays can be identified

Fig. 9
Fig. 9 Meso-scale model: procedure for the generation of the meso-scale unit cell.a Binary image extracted from growth ring B in the microscopy image of the meso-structure in the R-T plane, shown in Fig. 8a.The cell wall material is represented by the black pixels.b Smoothing of the jagged-shaped edges of the cell walls using a second-order B-spline interpolation function.The solid circles represent the control points of the B-spline function.c Detail of the two-dimensional cell wall domain discretized by triangular finite elements.d Detection of the in-plane principal material directions of the axial cells in each finite element.e Detail of the three-dimensional model after extrusion, in which the ray cells are fully "open" along the longitudinal, out-of-plane direction.f Closure of the open sections through the addition of elements with ray cell wall properties in the two outer element layers in the longitudinal direction.g The final three-dimensional meso-scale unit cell of a growth ring, which is used for computing the effective properties at the macro-scale

Fig. 11
Fig. 11 Meso-scale model.Effective hygro-elastic properties of cell walls of axial cells (normalized with respect to their value in the dry state) as a function of the moisture content: a Young's moduli; b Poisson's ratios; c shear moduli; d hygro-expansion coefficients

Fig. 12
Fig.12Macro-scale results: effective moisture adsorption isotherm of oak wood as obtained from Eq. (16), using the meso-scale moisture adsorption isotherms in Fig.10as input.The experimental data have been taken fromLuimes and Suiker (2021)

Fig. 13
Fig. 13 Macro-scale results: macro-scale Young's moduli of oak as a function of the moisture content.a Effective Young's modulus ĒR in the radial direction; b effective Young's modulus ĒT in the tangential direction; c effective Young's modulus ĒL in the longitudinal direction.The results relate to an MFA of 10 • and 20 • in the S2 layer of the axial cells.The experimental data have been taken from Ozyhar et al. (2016), Güntekın et al. (2016a), Luimes et al. (2018b) andGüntekın et al. (2016b)

Fig. 14
Fig. 14 Macro-scale results: macro-scale Poisson's ratios of oak as a function of the moisture content.a Effective Poisson's ratio νLT in the tangential plane; b effective Poisson's ratio νLR in the radial plane; c effective Poisson's ratio νRT in the transverse plane.The results correspond to an MFA of 10 • and 20 • in the S2 layer of axial cells.The experimental data have been taken from Ozyhar et al. (2016); Güntekın et al. (2016a); Luimes et al. (2018b); Güntekın et al. (2016b)

Fig. 15
Fig. 15 Macro-scale results: macro-scale shear moduli of oak as a function of the moisture content.a Effective shear modulus ḠLT in the tangential plane; b effective shear modulus ḠLR in the radial plane; c effective shear modulus ḠRT in the transverse plane.The results correspond to an MFA of 10 • and 20 • in the S2 layer of axial cells.The experimental data have been taken from Ozyhar et al. (2016), Güntekın et al. (2016b) and Luimes et al. (2018b)

Fig. 16
Fig. 16 Macro-scale results: effective macro-scale hygro-expansive properties of oak as a function of the moisture content.a Effective hygro-expansion coefficient βR in the radial direction; b effective hygroexpansion coefficient βT in the tangential direction; c effective hygro-expansion coefficient βL in the longitudinal direction.The results correspond to an MFA of 10 • and 20 • in the S2 layer of axial cells.The spread of the experimental data, which have been taken fromOzyhar et al. (2016),Badel and Perré (2006) andKránitz (2014), is indicated by the gray region

Table 1
Malek and Gibson (2017)tion and cell width of the nano-scale unit cells characterizing the composite fibers in the S1, S2, S3, and CML layers.The compositions have been taken fromMalek and Gibson (2017)

Table 2
Nano-scale model

Table 3
Micro-scale modelRelative layer thickness, t lay and micro-fibril angle, MFA, of cell wall layers containing axial cells and ray cells

Table 5
Meso-scale modelCharacteristics of the meso-scale unit cells.The unit cells are generated from microscopy images of the three different oak growth rings A, B, and C displayed in Fig.8a.The table lists the porosity , density , unit cell width W, ray volume fraction v ray , width of the earlywood (EW) region W e and vessel volume fraction v vess,EW in EW