Evaluating Timoshenko Method for Analyzing CLT under Out-of-Plane Loading

Cross-laminated timber (CLT) is an engineered wood product made up of layers of structurally graded timber, where subsequent layers are oriented orthogonally to each other. In CLT, the layers oriented in transverse direction, generally termed as cross-layer, are subjected to shear in radial–tangential plane, which is commonly known as rolling shear. As the shear modulus of cross-layers is significantly lower than that in other planes, CLT exhibits higher shear deformation under out-of-plane loading in contrast to other engineered wood products such as laminated veneer lumber (LVL) and glue laminated timber (GLT). Several analytical methods such as Timoshenko, modified gamma and shear analogy methods were proposed to account for this excessive shear deformation in CLT. This paper focuses on the effectiveness of Timoshenko method in hybrid CLT, in which hardwood cross-layers are used due to their higher rolling shear modulus. A comprehensive numerical study was conducted and obtained results were carefully analyzed for a range of hybrid combinations. It was observed that Timoshenko method could not accurately predict the shear response of CLTs with hardwood cross layers. Comprehensive parametric analysis was conducted to generate reliable numerical results, which were subsequently used to propose modified design equations for hybrid CLTs.


Introduction
Structural timber products can be classified into structurally graded timber boards and engineered wood products. Products such as laminated lumber veneer (LVL), glue laminated timber (GLT), and plywood have been used for decades [1]. Cross-laminated timber (CLT), on the other hand, is relatively new to the construction industry. The concept of CLT was originated in Austria in the early 1980s. Unlike GLT in CLT subsequent lamina are orthogonal to each other along thickness direction. Most cross-laminated timber is made up of 3 to 5 layers although 7 and 9 layered CLTs are also used in construction [2]. Due to its orthogonal make up, CLT can sustain bi-directional bending in contrast to other structural timber products which typically carry the load in one direction [3,4]. The makeup of CLT could be either homogeneous or heterogeneous. In homogeneous CLT, timber boards used in longitudinal layers and cross-layers are of same structural grades, whereas in heterogeneous CLT, the cross-layers are made of timber with inferior mechanical properties [5].
CLT offers relatively lower bending stiffness when compared to an equivalent GLT section. This can be attributed to the orthogonal layout of cross-layers and the lower rolling shear modulus in timber [6]. Timber boards are orthotropic in nature, where the mechanical properties vary along the longitudinal, tangential, and radial direction [7]. To improve performance of CLT under out-of-plane Timoshenko beam theory is an extension of Bernoulli-Euler beam theory to account for shear deformation in thick composite beams. Geometric properties of a typical three-layered CLT is illustrated in Figure 1. Unlike Euler-Bernoulli theory, the cross-section is not assumed to remain straight and perpendicular to the neutral axis in Timoshenko theory which uses a shear form factor k to correct shear stiffness of the cross-section. The bending and shear stiffness of laminated beams are calculated using Equations (1) and (2), respectively [19,20].
where Ei and Gi are the elastic modulus and the shear modulus in MPa of the i-th layer. bi, hi, zi and a are width of the panel, height of the i-th layer, distance of the neutral axis of the i-th layer from the neutral axis of the panel in mm, and the distance between the neutral axis of the outer layers, respectively. Shear form factor, k is determined by using Equation (3) [21].
S (z) = First moment of area (mm 3 ). Deflection of the simply supported panel under uniformly distributed loading (UDL) condition can be determined using Equation (4).
where q is the uniformly distributed loading in N/mm and L is the span length of the panel in mm [22].

Shear Form Factor
As mentioned before, Timoshenko beam theory is an extension of Euler-Bernoulli theory, where the deflection due to shear deformation is also considered for calculation of total deflection. The kinematics assumptions used by Timoshenko beam theory results in presence of constant shear deformation through thickness of the section under out-of-plane bending. To compensate for the excess resulting shear deformation against results obtained from elasticity solutions, shear form factor has been introduced by Timoshenko. A range of shear form factor values has been proposed literature. The shear form factor suggested by Timoshenko originally was 3/2 for isotropic homogeneous solid section [23]. The most commonly used shear form factor value is 6/5 proposed by Reissner [23] based on First order shear deformation theory (FSDT). Shear form factor value of 12/π 2 was proposed by Mindlin [24] and shear form factor was proposed to be a function of Poisson's ratio by Cowper [25].
The estimation of shear form factor was extended for laminated composite beams by Barbero et al. [26], where he proposed a mechanics for thin-walled laminated beams. A theoretical approach Unlike Euler-Bernoulli theory, the cross-section is not assumed to remain straight and perpendicular to the neutral axis in Timoshenko theory which uses a shear form factor k to correct shear stiffness of the cross-section. The bending and shear stiffness of laminated beams are calculated using Equations (1) and (2), respectively [19,20].
where A i = bh i .
where E i and G i are the elastic modulus and the shear modulus in MPa of the i-th layer. b i , h i , z i and a are width of the panel, height of the i-th layer, distance of the neutral axis of the i-th layer from the neutral axis of the panel in mm, and the distance between the neutral axis of the outer layers, respectively. Shear form factor, k is determined by using Equation (3) [21].
S (z) = First moment of area (mm 3 ). Deflection of the simply supported panel under uniformly distributed loading (UDL) condition can be determined using Equation (4).
where q is the uniformly distributed loading in N/mm and L is the span length of the panel in mm [22].

Shear Form Factor
As mentioned before, Timoshenko beam theory is an extension of Euler-Bernoulli theory, where the deflection due to shear deformation is also considered for calculation of total deflection. The kinematics assumptions used by Timoshenko beam theory results in presence of constant shear deformation through thickness of the section under out-of-plane bending. To compensate for the excess resulting shear deformation against results obtained from elasticity solutions, shear form factor has been introduced by Timoshenko. A range of shear form factor values has been proposed literature. The shear form factor suggested by Timoshenko originally was 3/2 for isotropic homogeneous solid section [23]. The most commonly used shear form factor value is 6/5 proposed by Reissner [23] based on First order shear deformation theory (FSDT). Shear form factor value of 12/π 2 was proposed by Mindlin [24] and shear form factor was proposed to be a function of Poisson's ratio by Cowper [25]. The estimation of shear form factor was extended for laminated composite beams by Barbero et al. [26], where he proposed a mechanics for thin-walled laminated beams. A theoretical approach was proposed by Davalos et al. [27] to evaluate the shear form factor for thick composite laminate with arbitrary lay-up configurations of layers. The evaluation of shear form factor in accordance to Equation (3) was proposed by Augustine et al. [28].

Finite Element Modelling
Three dimensional models of CLT panels with varying span-to-depth ratio and rolling shear modulus for the cross-layer were developed using commercial FE (Finite Element) package Abaqus 2016. All CLT panels were subjected to uniformly distributed loading and resulting deflections were captured from FE models to evaluate the performance of the considered analytical technique, i.e., Timoshenko method in predicting those deflections.

Material Model
CLT is a composite material and exhibits linear elastic behavior until some form of mechanical failure takes place. In practice, non-linearity in structural response is typically initiated due to either tensile failure in the bottom layer or rolling shear failure in the cross-layer. This study intends to investigate deflection behavior of CLT panels within linear elastic range, and hence linear elastic material model was used.
Timber is an orthotropic material, whose elastic behavior is different in three orthogonal directions [29]. The three orthogonal directions in timber are longitudinal, tangential, and radial direction, as illustrated in Figure 2.
Buildings 2020, 10, x FOR PEER REVIEW 4 of 21 was proposed by Davalos et al. [27] to evaluate the shear form factor for thick composite laminate with arbitrary lay-up configurations of layers. The evaluation of shear form factor in accordance to Equation (3) was proposed by Augustine et al. [28].

Finite Element Modelling
Three dimensional models of CLT panels with varying span-to-depth ratio and rolling shear modulus for the cross-layer were developed using commercial FE (Finite Element) package Abaqus 2016. All CLT panels were subjected to uniformly distributed loading and resulting deflections were captured from FE models to evaluate the performance of the considered analytical technique, i.e., Timoshenko method in predicting those deflections.

Material Model
CLT is a composite material and exhibits linear elastic behavior until some form of mechanical failure takes place. In practice, non-linearity in structural response is typically initiated due to either tensile failure in the bottom layer or rolling shear failure in the cross-layer. This study intends to investigate deflection behavior of CLT panels within linear elastic range, and hence linear elastic material model was used.
Timber is an orthotropic material, whose elastic behavior is different in three orthogonal directions [29]. The three orthogonal directions in timber are longitudinal, tangential, and radial direction, as illustrated in Figure 2. Nine independent elastic constants are required to determine the mechanical response of orthotropic materials. These constants are three elastic moduli, three Poisson's ratios, and three shear moduli. Orthotropic material property for timber can be reduced to transversely isotropic material Nine independent elastic constants are required to determine the mechanical response of orthotropic materials. These constants are three elastic moduli, three Poisson's ratios, and three shear moduli. Orthotropic material property for timber can be reduced to transversely isotropic material by considering the radial-tangential plane to be isotropic, as the material properties are reported to be of the same order [30]. This assumption means that elastic moduli (E R and E T ) and Poisson's ratio (ν RT and ν TR ) normal to each other in the radial-tangential plane becomes equal. Further, the shear moduli in the other two perpendicular planes become equal. The additional constraint can be expressed as where L, R, and T represent longitudinal, tangential, and radial material direction, respectively. Due to the constraint imposed by Equation (5), two independent Poisson's ratios are reduced to one as two elastic modulus values are equal, resulting in reduction of nine-independent constants to 5. The consequent constitutive relationship is shown in Equation (6).
where R plane is considered to be isotropic and subscripts T are replaced by subscript R.
Elastic constants for timber material with reference to Figures 2b and 3 are listed in Table 1.
Buildings 2020, 10, x FOR PEER REVIEW 6 of 21 The parametric study consisted of 52 models of varying combination of span and thickness. Span lengths of the CLT panels were varied between 2 and 5 m considering current design practice for engineering wood products. Span lengths were changed by 250 mm increment between 2 and 5 m. Three-layered CLT panels with total thickness of 90, 120, 150, and 180 mm were used for each span length. These variations produced a wide range of L/D ratio ranging between 11.11 and 58.33. The uniform loading was considered constant at 1.5 N/mm so that the maximum deflection was not excessive.

Material Properties
Mechanical properties used for longitudinal layers both in FE models and in analytical calculations were mostly taken from standard values recommended for C24 grade structural timber [31]. The longitudinal elastic modulus (EL) used in the study is slightly higher than that of C24 grade, which is representative of other species including hardwoods. CLT was modelled assuming timber being transversely isotropic, and the corresponding material properties are shown in Table 2. The rolling shear modulus (GRT) of the cross-layers in each model was varied between 50 MPa and 250 MPa with an increment of 50 MPa. The range of rolling shear modulus was chosen based on the values reported for Norway spruce (a softwood species) [32] and those for hardwood species [33]. The radial and tangential elastic moduli in the cross-layers of models were varied according to their relationship with rolling shear modulus and Poisson's ratio as shown in Table 1. The Poisson's ratio values were obtained from experimental research conducted on Norway spruce [34].

Shear modulus
Longitudinal-tangential plane The shear modulus in radial-tangential plane (G RT ), commonly referred to as the rolling shear modulus, is dependent on the Poisson's ratio and the modulus of elasticity of that plane. It is worth noting, that determining rolling shear modulus using physical tests is easier than determining Poisson's ratio.

Geometry and Boundary Conditions
FE models were developed using fully integrated solid elements (C3D20) available in Abaqus 2016. Cubic elements of 8 mm were chosen based on conducted mesh convergence study. The width of all the models were kept constant at 200 mm. Span and thickness of the developed models were varied to investigate their effects on the behavior of panels under out-of-plane loading. The models were subjected to uniformly distributed load for simply supported boundary condition. Figure 3 shows the model geometry and applied loading.
The parametric study consisted of 52 models of varying combination of span and thickness. Span lengths of the CLT panels were varied between 2 and 5 m considering current design practice for engineering wood products. Span lengths were changed by 250 mm increment between 2 and 5 m. Three-layered CLT panels with total thickness of 90, 120, 150, and 180 mm were used for each span length. These variations produced a wide range of L/D ratio ranging between 11.11 and 58.33. The uniform loading was considered constant at 1.5 N/mm so that the maximum deflection was not excessive.

Material Properties
Mechanical properties used for longitudinal layers both in FE models and in analytical calculations were mostly taken from standard values recommended for C24 grade structural timber [31]. The longitudinal elastic modulus (E L ) used in the study is slightly higher than that of C24 grade, which is representative of other species including hardwoods. CLT was modelled assuming timber being transversely isotropic, and the corresponding material properties are shown in Table 2. The rolling shear modulus (G RT ) of the cross-layers in each model was varied between 50 MPa and 250 MPa with an increment of 50 MPa. The range of rolling shear modulus was chosen based on the values reported for Norway spruce (a softwood species) [32] and those for hardwood species [33]. The radial and tangential elastic moduli in the cross-layers of models were varied according to their relationship with rolling shear modulus and Poisson's ratio as shown in Table 1. The Poisson's ratio values were obtained from experimental research conducted on Norway spruce [34].

Verification of the Numerical Modelling Technique
FE models were verified using experimental results obtained from four point bending tests conducted in Deakin University. The equipment used is Universal open structure flexural frame with 300 kN maximum capacity of loading (UNIFLEX 300). It was sourced from Italy (Milano) through GEO-CON Products Pty. Ltd. The actual test setup and geometric details are shown in Figure 4.  The CLT used in the test was sourced from Europe and was manufactured from Norway spruce representing structural grade of C24. In the developed FE model, material properties were taken from Table 3 and the rolling shear modulus of cross-layer was taken as 50 MPa as this is representative of Norway spruce [35]. Tests were conducted under displacement control and the rate of displacement was kept at 1 mm/minute. Figure 5 shows the load-displacement behavior of 5 test specimens along with the response obtained from the developed elastic FE model. As previously mentioned, the FE model developed as part of the current study was elastic, and hence produced a straight line response as shown in Figure 5. It should also be noted that slip between the interfaces is considered to be negligible and consequently tie constrained is used between the interfaces. Figure 5a shows developed shear stresses in the radial-tangential plane of the cross layer due to applied loading, and it is obvious that significant rolling shear stresses developed at the end spans due to shear deformation of the beam. Figure 5b shows good agreement between the FE model and test results within the elastic range, and hence the adopted modelling technique may safely be used to predict deflection of CLT panels as long as they remain linearly elastic, which represents the typical design criteria in timber structures. The CLT used in the test was sourced from Europe and was manufactured from Norway spruce representing structural grade of C24. In the developed FE model, material properties were taken from Table 3 and the rolling shear modulus of cross-layer was taken as 50 MPa as this is representative of Norway spruce [35]. Tests were conducted under displacement control and the rate of displacement was kept at 1 mm/min. Figure 5 shows the load-displacement behavior of 5 test specimens along with the response obtained from the developed elastic FE model.

Numerical Modelling of Glue Laminated Timber (GLT)
Timoshenko method in the current study is reported to produce accurate predictions for glue laminated timber (GLT) beams. To assess the accuracy of the adopted FE modelling technique, numerical models were developed for GLT beams and obtained results were compared against those calculated using the analytical techniques. The mechanical properties used for both theoretical calculations and FE modelling are shown in Table 4 with respect to those used in CLT. GLT is also made of structural timber boards but contrary to CLT, the boards aligned in one direction. This entails that all the modulus values used in the analytical methods are the same, as listed in Table 4. The load used for calculation and modelling was 1.5 N/mm. As previously mentioned, the FE model developed as part of the current study was elastic, and hence produced a straight line response as shown in Figure 5. It should also be noted that slip between the interfaces is considered to be negligible and consequently tie constrained is used between the interfaces. Figure 5a shows developed shear stresses in the radial-tangential plane of the cross layer due to applied loading, and it is obvious that significant rolling shear stresses developed at the end Buildings 2020, 10, 184 8 of 21 spans due to shear deformation of the beam. Figure 5b shows good agreement between the FE model and test results within the elastic range, and hence the adopted modelling technique may safely be used to predict deflection of CLT panels as long as they remain linearly elastic, which represents the typical design criteria in timber structures.

Numerical Modelling of Glue Laminated Timber (GLT)
Timoshenko method in the current study is reported to produce accurate predictions for glue laminated timber (GLT) beams. To assess the accuracy of the adopted FE modelling technique, numerical models were developed for GLT beams and obtained results were compared against those calculated using the analytical techniques. The mechanical properties used for both theoretical calculations and FE modelling are shown in Table 4 with respect to those used in CLT. GLT is also made of structural timber boards but contrary to CLT, the boards aligned in one direction. This entails that all the modulus values used in the analytical methods are the same, as listed in Table 4. The load used for calculation and modelling was 1.5 N/mm. Deflections obtained for a range of span-to-depth ratios is presented in Table 5. Which shows excellent agreement between deflections obtained using FE models and those predicted analytically. In GLT sections, due to higher shear modulus values of cross-layer the shear deformation is very much limited. The results according to Table 5, indicates that there is good agreement between the analytical and Numerical results where shear deformation is insignificant.

Reassessing the Suitability of Timoshenko Method for CLT Panels
In this section, results obtained through Timoshenko method and FE models are compared. Then the different parameters that affects the deviation between Timoshenko method and Numerical models are identified and their impact on the results are discussed.

Comparison between Timoshenko Method and CLT Panels
FE models of CLT panels were used to estimate maximum deflections within the considered span, and obtained FE predictions were compared against those determined using Timoshenko method. Observed deviations in deflection predictions for a range of span-to-depth ratio are presented in Figure 6. The range of span length and panel depth used in the study was varied between 1026 to 2500 mm and 60 to 120 mm, respectively. The mechanical properties used for the modelling of the panels Buildings 2020, 10, 184 9 of 21 are listed in Table 2 and the schematic of the model is in illustrated in Figure 3. It is worth noting that Figure 6 presents observed trends in terms of span-to-depth (L/D) ratio.

Comparison between Timoshenko Method and CLT Panels
FE models of CLT panels were used to estimate maximum deflections within the considered span, and obtained FE predictions were compared against those determined using Timoshenko method. Observed deviations in deflection predictions for a range of span-to-depth ratio are presented in Figure 6. The range of span length and panel depth used in the study was varied between 1026 to 2500 mm and 60 to 120 mm, respectively. The mechanical properties used for the modelling of the panels are listed in Table 2 and the schematic of the model is in illustrated in Figure 3. It is worth noting that Figure 6 presents observed trends in terms of span-to-depth (L/D) ratio. Figure 6. Deviations in deflection predictions obtained using FE models vs. analytical methods.
In Figure 6, the percentage error was calculated based on numerical predictions as shown in Equation (7).
In Figure 6, the deviations in deflection obtained as percentage for FE models against those calculated using Timoshenko were plotted vs. rolling shear modulus of cross-layer as this has been reported to be one of the key parameters in dictating CLT deflection. Timoshenko method tends to underestimate the deflections varying within 0 percent to 15 percent. Timoshenko method in accordance to Figure 6, is found to be inaccurate as the rolling shear modulus of the cross-layer increases. It is also observed that, with increasing span-to-depth ratio, the deviation between Timoshenko method and FE models decrease considerably. From this study, it is concluded that span-to-depth ratio and rolling shear modulus are the parameters that affect the accuracy of Timoshenko method.

Parameters Affecting the Performance of Timoshenko Method
Numerical models were carefully used to conduct a comprehensive parametric study to investigate the suitability of Timoshenko method in predicting deflections of CLT panels. In Figure 6, the percentage error was calculated based on numerical predictions as shown in Equation (7).
Percentage error = de f lection f rom FE model − de f lection f rom analytical method de f lection f rom FE model × 100 In Figure 6, the deviations in deflection obtained as percentage for FE models against those calculated using Timoshenko were plotted vs. rolling shear modulus of cross-layer as this has been reported to be one of the key parameters in dictating CLT deflection. Timoshenko method tends to underestimate the deflections varying within 0 percent to 15 percent. Timoshenko method in accordance to Figure 6, is found to be inaccurate as the rolling shear modulus of the cross-layer increases. It is also observed that, with increasing span-to-depth ratio, the deviation between Timoshenko method and FE models decrease considerably. From this study, it is concluded that span-to-depth ratio and rolling shear modulus are the parameters that affect the accuracy of Timoshenko method.

Parameters Affecting the Performance of Timoshenko Method
Numerical models were carefully used to conduct a comprehensive parametric study to investigate the suitability of Timoshenko method in predicting deflections of CLT panels. Figure 7 shows observed deviations in maximum deflection predictions against varying length and for 5 different rolling shear modulus considering four panel thickness (90, 120, 150, and 180 mm). It is observed that Timoshenko method always under-predicts deflections in CLT panels, and predictions become more inaccurate as the panels become thick and span becomes shorter for a given panel thickness producing small values for L/D ratio. different rolling shear modulus considering four panel thickness (90, 120, 150, and 180 mm). It is observed that Timoshenko method always under-predicts deflections in CLT panels, and predictions become more inaccurate as the panels become thick and span becomes shorter for a given panel thickness producing small values for L/D ratio. indicates that the observed deviation in deflection predictions is somewhat inversely proportional to the L/D ratio for a given rolling shear. Lower L/D ratio results in significant contribution from shear deformation towards total deflection, and Timoshenko method cannot capture this contribution producing conservative results. It is also observed that for a given span, the deviation in analytical predictions are higher for panels with higher rolling shear modulus. These observations clearly indicate that the Timoshenko method should be modified before it can be used to predict deflections in CLT, especially in hybrid CLT with cross-layers of higher characteristic rolling shear modulus.
Since the observed deviation in deflection prediction is contributed primarily by shear deformation, a back calculation approach was adopted to calibrate Timoshenko method. Shear stiffness values calculated from Timoshenko method were compared to those estimated using Finite element model results for deflection. The effective shear stiffness of the section was calculated using the Equations (2) and (3) for Timoshenko method, whilst Equation (4) was used to determine apparent shear stiffness according to the deflection obtained from FE models. Figure 8 presents deviations between analytical vs. FE shear stiffness of CLT sections considered in the current study.  Figure 7 indicates that the observed deviation in deflection predictions is somewhat inversely proportional to the L/D ratio for a given rolling shear. Lower L/D ratio results in significant contribution from shear deformation towards total deflection, and Timoshenko method cannot capture this contribution producing conservative results. It is also observed that for a given span, the deviation in analytical predictions are higher for panels with higher rolling shear modulus. These observations clearly indicate that the Timoshenko method should be modified before it can be used to predict deflections in CLT, especially in hybrid CLT with cross-layers of higher characteristic rolling shear modulus.
Since the observed deviation in deflection prediction is contributed primarily by shear deformation, a back calculation approach was adopted to calibrate Timoshenko method. Shear stiffness values calculated from Timoshenko method were compared to those estimated using Finite element model results for deflection. The effective shear stiffness of the section was calculated using the Equations (2) and (3) for Timoshenko method, whilst Equation (4) was used to determine apparent shear stiffness according to the deflection obtained from FE models. Linear trends observed in Figure 8 resemble Figure 7, and clearly justifies the hypothesis that the deviations observed between analytically obtained deflections vs. FE results are primarily due to the inaccurate estimation of shear stiffness using the current form of Timoshenko method.

Effect of Rolling Shear Modulus
To further substantiate the above hypothesis, the results for CLT panel of span-to-depth ratio of 11.11 is considered, since lower L/D ratio was found to be less accurate for deflection calculation. Figure 9 compares the deflection obtained from Timoshenko method and the FE models. From Figure  9, it is observed that the difference between deflection estimated from Timoshenko method and FEM increases as the rolling shear modulus increases. For shear stiffness, the same trend is observed. Linear trends observed in Figure 8 resemble Figure 7, and clearly justifies the hypothesis that the deviations observed between analytically obtained deflections vs. FE results are primarily due to the inaccurate estimation of shear stiffness using the current form of Timoshenko method.

Effect of Rolling Shear Modulus
To further substantiate the above hypothesis, the results for CLT panel of span-to-depth ratio of 11.11 is considered, since lower L/D ratio was found to be less accurate for deflection calculation. Figure 9 compares the deflection obtained from Timoshenko method and the FE models. From Figure 9, it is observed that the difference between deflection estimated from Timoshenko method and FEM increases as the rolling shear modulus increases. For shear stiffness, the same trend is observed. Buildings 2020, 10, x FOR PEER REVIEW 12 of 21 (a) (b) Figure 9. Variation of (a) deflection and (b) shear stiffness of the section of panel with length to depth ratio of 11.11. Figure 9a demonstrates that with increasing rolling shear modulus of the cross-layer the deviation between deflections obtained from Timoshenko method and FE models for panels of spanto-depth ratio of 11.11, increases. The increase in deviation of deflection however has decreasing gradient with increasing rolling shear modulus, suggesting that increasing rolling shear modulus has prominent effect on the percentage deviation of deflection for lower values. In Figure 9b, the shear stiffness of the CLT panel according to Timoshenko method and FE models are compared for changing rolling shear modulus. With increasing rolling shear modulus the shear stiffness values are found to increasing at different rate resulting in increasing deviation between shear stiffness obtained from Timoshenko method and FE models. The percentage deviation against rolling shear modulus for both deflection and shear stiffness has been found to show similar trend.

Effect of Length and Depth
The effect of span-to-depth ratio on the deviations of deflection and shear stiffness between Timoshenko method and FEM are shown in Figure 10 for the rolling shear modulus of 250 MPa along with the deviation in estimation of shear stiffness of the sections.
(a) (b) Figure 10. Deviation of (a) deflection and (b) shear stiffness against changing span-to-depth ratio of Timoshenko method from FE models.
In Figure 10, the deviation as percentage for deflection and respective shear stiffness obtained between FEM results and Timoshenko method expressed against the whole range of span-to-depth ratio used in the study. It is observed in Figure 10a that the higher the span-to-depth ratio, the lower the error in deflection obtained from Timoshenko method with respect to the FE models. This reduction of error with increasing span-to-depth ratio is due to relatively smaller contribution of  Figure 9a demonstrates that with increasing rolling shear modulus of the cross-layer the deviation between deflections obtained from Timoshenko method and FE models for panels of span-to-depth ratio of 11.11, increases. The increase in deviation of deflection however has decreasing gradient with increasing rolling shear modulus, suggesting that increasing rolling shear modulus has prominent effect on the percentage deviation of deflection for lower values. In Figure 9b, the shear stiffness of the CLT panel according to Timoshenko method and FE models are compared for changing rolling shear modulus. With increasing rolling shear modulus the shear stiffness values are found to increasing at different rate resulting in increasing deviation between shear stiffness obtained from Timoshenko method and FE models. The percentage deviation against rolling shear modulus for both deflection and shear stiffness has been found to show similar trend.

Effect of Length and Depth
The effect of span-to-depth ratio on the deviations of deflection and shear stiffness between Timoshenko method and FEM are shown in Figure 10 for the rolling shear modulus of 250 MPa along with the deviation in estimation of shear stiffness of the sections. Figure 9a demonstrates that with increasing rolling shear modulus of the cross-layer the deviation between deflections obtained from Timoshenko method and FE models for panels of spanto-depth ratio of 11.11, increases. The increase in deviation of deflection however has decreasing gradient with increasing rolling shear modulus, suggesting that increasing rolling shear modulus has prominent effect on the percentage deviation of deflection for lower values. In Figure 9b, the shear stiffness of the CLT panel according to Timoshenko method and FE models are compared for changing rolling shear modulus. With increasing rolling shear modulus the shear stiffness values are found to increasing at different rate resulting in increasing deviation between shear stiffness obtained from Timoshenko method and FE models. The percentage deviation against rolling shear modulus for both deflection and shear stiffness has been found to show similar trend.

Effect of Length and Depth
The effect of span-to-depth ratio on the deviations of deflection and shear stiffness between Timoshenko method and FEM are shown in Figure 10 for the rolling shear modulus of 250 MPa along with the deviation in estimation of shear stiffness of the sections.
(a) (b) Figure 10. Deviation of (a) deflection and (b) shear stiffness against changing span-to-depth ratio of Timoshenko method from FE models.
In Figure 10, the deviation as percentage for deflection and respective shear stiffness obtained between FEM results and Timoshenko method expressed against the whole range of span-to-depth ratio used in the study. It is observed in Figure 10a that the higher the span-to-depth ratio, the lower the error in deflection obtained from Timoshenko method with respect to the FE models. This reduction of error with increasing span-to-depth ratio is due to relatively smaller contribution of Deviation of deflection (%) Figure 10. Deviation of (a) deflection and (b) shear stiffness against changing span-to-depth ratio of Timoshenko method from FE models.
In Figure 10, the deviation as percentage for deflection and respective shear stiffness obtained between FEM results and Timoshenko method expressed against the whole range of span-to-depth ratio used in the study. It is observed in Figure 10a that the higher the span-to-depth ratio, the lower the error in deflection obtained from Timoshenko method with respect to the FE models. This reduction of error with increasing span-to-depth ratio is due to relatively smaller contribution of deflection due to shear deformation of cross layer to the total deflection as span-to-depth ratio increases. The deviation of shear stiffness against span-to-depth ratio illustrated in Figure 10b shows pronounced noise. To explain this noise, the effect of length and depth on the deviation of shear stiffness is further investigated separately. The shear stiffness obtained from Timoshenko method and FE models against the whole range of length from 2 to 5 m for the panels of thickness 90 mm and 180 mm with rolling shear modulus of cross-layers of 250 MPa is shown in Figure 11.
Buildings 2020, 10, x FOR PEER REVIEW 13 of 21 deflection due to shear deformation of cross layer to the total deflection as span-to-depth ratio increases. The deviation of shear stiffness against span-to-depth ratio illustrated in Figure 10b shows pronounced noise. To explain this noise, the effect of length and depth on the deviation of shear stiffness is further investigated separately. The shear stiffness obtained from Timoshenko method and FE models against the whole range of length from 2 to 5 m for the panels of thickness 90 mm and 180 mm with rolling shear modulus of cross-layers of 250 MPa is shown in Figure 11. The shear stiffness obtained from Timoshenko method is independent of length of the panel as can be seen from Figure 11, whereas the shear stiffness estimated from the FE models shows increasing trend with increasing length of the panels. Figure 12 also suggests that the FE model converge towards the Timoshenko method as length of the panel increases. The shear stiffness obtained from Timoshenko method is independent of length of the panel as can be seen from Figure 11, whereas the shear stiffness estimated from the FE models shows increasing trend with increasing length of the panels. Figure 12 also suggests that the FE model converge towards the Timoshenko method as length of the panel increases.
Buildings 2020, 10, x FOR PEER REVIEW 13 of 21 deflection due to shear deformation of cross layer to the total deflection as span-to-depth ratio increases. The deviation of shear stiffness against span-to-depth ratio illustrated in Figure 10b shows pronounced noise. To explain this noise, the effect of length and depth on the deviation of shear stiffness is further investigated separately. The shear stiffness obtained from Timoshenko method and FE models against the whole range of length from 2 to 5 m for the panels of thickness 90 mm and 180 mm with rolling shear modulus of cross-layers of 250 MPa is shown in Figure 11. The shear stiffness obtained from Timoshenko method is independent of length of the panel as can be seen from Figure 11, whereas the shear stiffness estimated from the FE models shows increasing trend with increasing length of the panels. Figure 12 also suggests that the FE model converge towards the Timoshenko method as length of the panel increases. In Figure 12 the deviation between shear stiffness from Timoshenko method using Equation (2) to that obtained from FE models is presented. It is observed that the deviation ranges between 38 percent and 22 percent for panels with depth of 90 mm, whereas for panels with 180 mm thickness, the deviation ranges between 58 percent and 35 percent. This suggests that as the depth increases the deviation between the shear stiffness obtained from Timoshenko method and FE models increases. From Figures 11 and 12 it is apparent that both length and depth of the panels have effect on the shear stiffness obtained from FE models. In contrast, the length has no effect on the shear stiffness obtained from the original Timoshenko method according to Equation (2). Further, for higher length, decrease of deviation between results obtained from Timoshenko method and FE models is observed, in contrast higher depth results higher deviation between the analytical and numerical values. Hence, it is concluded that the effect of length and depth needs to be considered independently.
From the above discussions in Sections 3 and 4, it can be summarized that rolling shear modulus, length and depth are mostly responsible for the inaccurate prediction of shear stiffness proposed in the Timoshenko method which resulted in significant error for the deflection calculation.

Modified Shear Form Factor
In this section, a new coefficient is introduced to accommodate the effect of length and depth of a CLT panel to determine the effective shear stiffness. Effective shear stiffness is calculated from Equation (2) in the Timoshenko method which contains a shear form factor, k determined using Equation (3). The methodology used for the determination of the modification coefficient is illustrated in Figure 13.
Buildings 2020, 10, x FOR PEER REVIEW 14 of 21 In Figure 12 the deviation between shear stiffness from Timoshenko method using Equation (2) to that obtained from FE models is presented. It is observed that the deviation ranges between 38 percent and 22 percent for panels with depth of 90 mm, whereas for panels with 180 mm thickness, the deviation ranges between 58 percent and 35 percent. This suggests that as the depth increases the deviation between the shear stiffness obtained from Timoshenko method and FE models increases. From Figures 11 and 12 it is apparent that both length and depth of the panels have effect on the shear stiffness obtained from FE models. In contrast, the length has no effect on the shear stiffness obtained from the original Timoshenko method according to Equation (2). Further, for higher length, decrease of deviation between results obtained from Timoshenko method and FE models is observed, in contrast higher depth results higher deviation between the analytical and numerical values. Hence, it is concluded that the effect of length and depth needs to be considered independently.
From the above discussions in Sections 3 and 4, it can be summarized that rolling shear modulus, length and depth are mostly responsible for the inaccurate prediction of shear stiffness proposed in the Timoshenko method which resulted in significant error for the deflection calculation.

Modified Shear Form Factor
In this section, a new coefficient is introduced to accommodate the effect of length and depth of a CLT panel to determine the effective shear stiffness. Effective shear stiffness is calculated from Equation (2) in the Timoshenko method which contains a shear form factor, determined using Equation (3). The methodology used for the determination of the modification coefficient is illustrated in Figure 13.   (2), the form factor km for the model is estimated

Modification Co-Efficient for Shear Form Factor
A modification co-efficient (c), function of shear modulus, span length, and depth of a panel, is introduced to modify the shear form factor proposed in the Timoshenko method. The modification co-efficient is defined as the ratio of shear form factors calculated from FE models (km) to the shear correction factor obtained from the Timoshenko method (k), as shown below.

=
(8) Figure 13. Methodology to determine the shear correcting factor. Figure 13 illustrates that the shear stiffness of the panel for the FE models are estimated using the deflection value obtained from the FE models and the corresponding geometric properties and mechanical properties. For estimations of the effective shear stiffness (GA) eff the corresponding bending stiffness calculated using Equation (1) along with respective length and uniformly distributed load (q) is used in Equation (4). Using the estimated shear stiffness (GA) eff in Equation (2), the form factor k m for the model is estimated

Modification Co-Efficient for Shear Form Factor
A modification co-efficient (c), function of shear modulus, span length, and depth of a panel, is introduced to modify the shear form factor proposed in the Timoshenko method. The modification co-efficient is defined as the ratio of shear form factors calculated from FE models (k m ) to the shear correction factor obtained from the Timoshenko method (k), as shown below.
The dependency of this modification co-efficient, c on the length, depth, and rolling shear modulus is illustrated in Figure 14. It can be observed that, rolling shear modulus and depth have positive correlation with the modification co-efficient, whereas length is found to have an inverse effect on the modification factor.
Buildings 2020, 10, x FOR PEER REVIEW 15 of 21 The dependency of this modification co-efficient, c on the length, depth, and rolling shear modulus is illustrated in Figure 14. It can be observed that, rolling shear modulus and depth have positive correlation with the modification co-efficient, whereas length is found to have an inverse effect on the modification factor.

Formulating Equation for Modification Co-Efficient
The FEM results obtained over the whole range of geometry and rolling shear modulus as mention in Sections 3.2 and 3.3 are used to establish the relationship between the identified parameters from Section 5.1 and the modification factor in the following section.

Relationship between Modification Factor, Length and Rolling Shear Modulus
Modification co-efficient for each of the rolling shear moduli is plotted against length for a specific depth of the panel and is depicted in Figure 14. From Figure 14a,c it is observed that, modification co-efficient demonstrates somewhat linear relationship with increasing rolling shear modulus and depth. In contrast, with increasing length the modification co-efficient approaches 1, as can be observed from Figure 14b. This observations are used to develop an equation for modification co-efficient.

Formulating Equation for Modification Co-Efficient
The FEM results obtained over the whole range of geometry and rolling shear modulus as mention in Sections 3.2 and 3.3 are used to establish the relationship between the identified parameters from Section 5.1 and the modification factor in the following section.

Relationship between Modification Factor, Length and Rolling Shear Modulus
Modification co-efficient for each of the rolling shear moduli is plotted against length for a specific depth of the panel and is depicted in Figure 14. From Figure 14a,c it is observed that, modification co-efficient demonstrates somewhat linear relationship with increasing rolling shear modulus and depth. In contrast, with increasing length the modification co-efficient approaches 1, as can be observed from Figure 14b. This observations are used to develop an equation for modification co-efficient.

The form of the equation was chosen to be
where d is the depth of the panel, G RT is the rolling shear modulus, and L is the span length. The rationale for this form of equation for Equation (9) are the conclusions that were explained in Section 5.1, which are that the magnitude of c with increasing length approaches 1 and is positively correlated to depth of the panel and the rolling shear modulus. The summary of the results obtained from the parametric study is illustrated in Figure 15. The form of the equation was chosen to be where d is the depth of the panel, GRT is the rolling shear modulus, and L is the span length. The rationale for this form of equation for Equation (9) are the conclusions that were explained in Section 5.1, which are that the magnitude of c with increasing length approaches 1 and is positively correlated to depth of the panel and the rolling shear modulus. The summary of the results obtained from the parametric study is illustrated in Figure 15. The relationship between the modification co-efficient with length, depth and rolling shear modulus was formulated as shown in Equation (10).
where L is length of the panels, d is the depth of the panels, and is the rolling shear modulus of the cross-layer. It is should be noted that the proposed method is limited to three layered CLT where layers of same thicknesses. Finally, the final form of the modified equation for Timoshenko method for a simply supported CLT panel under uniformly distributed loading is proposed as shown in Equation (11  The relationship between the modification co-efficient with length, depth and rolling shear modulus was formulated as shown in Equation (10).
where L is length of the panels, d is the depth of the panels, and G RT is the rolling shear modulus of the cross-layer. It is should be noted that the proposed method is limited to three layered CLT where layers of same thicknesses. Finally, the final form of the modified equation for Timoshenko method for a simply supported CLT panel under uniformly distributed loading is proposed as shown in Equation (11).
where w mod , is the modified equation for mid-span deflection of three-layered CLT panel; w s = 5 384 qL 4 (EI) e f f , is the shear free mid-span deflection; c, is the modification co-efficient factor evaluated by Equation (10); L, is the span length of the panel; (EI) eff , is the effective bending stiffness evaluated by Equation (1); and (GA) eff , is the effective shear stiffness evaluated by Equation (2).

Validation of the Proposed Factor
The proposed modification co-efficient was used to estimate the deflection of panels under four-point bending from literature. The validation was based on results from literature. All the individual mechanical properties required to determine deflection were often not present in those works and so values for those mechanical properties were obtained from different sources. This is summarized in Table 6. It should be noted that, the experimental results for uniformly distributed loading is scarce in literature and therefore, results from four-point bending has been used in this study. The test schematic is shown in Figure 16.  [29] 1 The mean ratio of G LT /E L for hardwood is 0.067 according Green et al. [29]; 2 The mean ratio of G RT /E L for hardwood is 0.017 according Green et al. [29]; 3 The mean ratio of G RT /E L for softwood is 0.0625 according to Green et al. [29] and Crovell et al. [38]. The deflection at 20 kN load for the panels deduced from the mean apparent bending stiffness recorded in the literature [38][39][40] under four point bending panels made from these species are estimated and compared against those estimated from Timoshenko method and modified Timoshenko method. The results are tabulated in Table 7. The corresponding geometry of the panels from the literature and the reference of the work is also included in the table. From the table, it is observed that modified Timoshenko method and Timoshenko method yields similar outcomes for panels with high span to depth ratio, but modified Timoshenko method seems to be more effective for panels with lower span to depth ratio. This is seen for Irish Sitka with dimension 1440 × 584 × 120, where the modified Timoshenko method outperforms Timoshenko method by a significant margin. Modified Timoshenko method is also found to perform better for the CLT made from Eucalyptus urophylla (hardwood) as it has significantly higher rolling shear modulus. For CLT panels made from Red Maple and White ash, which are also hardwoods, both Timoshenko and modified Timoshenko method yield similar results. This can be explained by the higher length of these specimens which was found to have negative correlation with the modification co-efficient as explained in Section 5.1.

Conclusions
The presented numerical study compared the accuracy of Timoshenko method for CLT panels under out-of-plane bending over a range of rolling shear modulus considering softwood-hardwood hybrid compositions. It is worth noting that the scope of the work was limited to three layered CLT panels subjected to uniformly distributed load with simply supported boundary conditions. It was observed that when the rolling shear modulus of the cross-layer of a CLT is considerably higher than that of the longitudinal layers, i.e., hardwood middle layer combined with softwood outer layers, the deflection predicted using Timoshenko method was significantly different from that obtained from FE models. This trend was significant in panels of lower span-to-depth ratio, with deviation in deflection varying from 12 percent to 22 percent for panels of span-to-depth ratio of 11.4 over the range of rolling shear modulus in between 50 to 250 MPa. A comprehensive parametric study was conducted using FE models, and it was observed that the error in deflection prediction was primarily due to inaccurate estimation of shear stiffness by Timoshenko method. Span length, panel depth, and rolling shear modulus of the cross-layer were identified to be the parameters affecting the shear stiffness of the models. A modification co-efficient was introduced as a function of identified parameters, i.e., length, depth, and rolling shear modulus to calculate the effective shear stiffness more accurately while using Timoshenko method. The proposed technique was shown to produce good agreement with relevant tests results.

Conflicts of Interest:
The authors declare no conflict of interest.