Finite element analysis of fibre reinforced timber beams under flexural loading

ABSTRACT

Laminated Veneer Lumber (LVL) is a popular product used in the building and construction industry because of its high strength, serviceability, aesthetic characteristics, and cost-effectiveness. However, LVL has traditionally been limited to residential and small light commercial buildings. To expand the use of LVL and reduce the need for steel and reinforced concrete structural members, a comprehensive computational study is conducted to understand the material’s structural performance and find an optimal strengthening ratio. In this paper, a new 3-dimensional Finite Element model was introduced and evaluated to predict the structural behaviour, strength properties, and failure damage modes of LVL timber beams reinforced with a Carbon Fibre Reinforcement Polymer (CFRP) sheet. The model is validated with three experimental tests from the literature, with the reinforced timber model showing a 1.8% error compared to the plain timber model’s 7% error. Additionally, a parametric study is conducted to understand the effects of reinforcement length on the structural behaviour of the reinforced beam, with mathematical equations established to express the relationships. The study found that reinforcing 50% of the beam span generally resulted in the most significant strength improvement. This investigation is made for various strength and serviceability parameters.

KEYWORDS:

• Engineered wood product
• Laminated veneer lumber
• Carbon fibre reinforcement polymer
• Finite element modelling
• Conventional shell composite layup
• Hashin damage criteria

1. Introduction

Timber has been a significant source of construction material in the building industry for many decades (Subhani et al. 2017). Timber was primarily used to construct small, lightweight buildings and short-span bridges due to its strength characteristics. Timber is widely used in structural systems because of its high strength-to-weight ratio and serviceability properties, including its aesthetic appeal (Rafn Thorhallsson, Ingi Hinriksson, and Thor Snæbjörnsson 2017; Ramage et al. 2017; Subhani et al. 2017; Yeoh et al. 2009; Yusof and Rahman 2017). These properties make timber a preferred choice over concrete or steel when designing structural systems. In terms of construction, timber-framing members are well-regarded for their ease of handling during construction. Unlike other materials, timber-framing members do not require heavy machinery to construct, and they can be constructed in less time, contributing to more economically efficient structural developments.

In recent times, the development of engineered wood products, such as Laminated Veneer Lumber (LVL), glued laminated wood (Glulam), Oriented Strand Board (OSB), plywood, and Cross-Laminated Timber (CLT), has been significant in producing large and sturdy framing members. This development has encouraged architects and designers to use these engineered wood products in framing mid-rise to high-rise buildings (Chun, Balen, and Pan 2016; Fossetti, Minafò, and Papia 2015; Janssens 2017; Khorsandnia and Crews 2015; Ramage et al. 2017; Subhani et al. 2017). These products offer economical and sustainable solutions as an alternative to hardwood products. However, engineers face challenges in designing long-span timber beams that meet acceptable serviceable limit states. For instance, designing a 300 mm deep timber beam spanning 6 m is typically not feasible, as it would exceed limit state requirements for deflection limits. Consequently, the timber framing option is replaced with a steel member, thus forfeiting its economic and sustainable advantages.

Several studies have explored the use of Fibre Reinforced Polymer (FRP) plates, wraps, or rebar composites for strengthening structural members, and they have shown significant opportunities for strength improvement. Various approaches have been introduced and applied to retrofit existing timber members, aiming to improve their mechanical properties and strength and stiffness behaviours. FRP composites are considered a practical strengthening component due to their proven high tensile and shear resistance capacities. They can be used to reinforce both existing structural timber members, such as those found in heritage buildings, and new timber members. Besides their high-strength properties, the ease of on-site or factory application further validates the effectiveness of using FRP in strengthening timber structures. These benefits have been highlighted in multiple studies (Fossetti, Minafò, and Papia 2015; Hamid Hashemi, Akbar Maghsoudi, and Rahgozar 2008; Janssens 2017; Jimmy Kim and Harries 2010; Raftery and Harte 2011; Ramage et al. 2017; Razavian et al. 2020; Subhani et al. 2017; Tareq Khshain Al-Saadi et al. 2019; Raftery and Harte 2013).

Researchers have investigated different types of reinforcing fibres in FRP composites, such as carbon, glass, and basalt fibres. Each of these materials has distinct mechanical properties. For instance, carbon fibres have the highest strength and stiffness capacities, as reported by Šedivka et al. (2015). However, in the case of reinforcing low-grade timber members, glass fibre is a more cost-effective option than carbon fibre (Raftery and Harte 2011). If sustainability and fireproof properties are of utmost importance for reinforcing timber structures, basalt fibre is the preferable choice to achieve serviceability requirements (Soares et al. 2016).

Researchers have explored various installation methods of FRP composites on timber members to identify the most effective reinforcement technique for achieving economically and architecturally acceptable sizes of framing members. Apart from the FRP installation method, the reinforcement performance is heavily influenced by several other factors, such as the timber type and grade, the mechanical properties of the FRP component, and the adhesive bonding agent.

Table 1 represents the findings of nine previous studies of reinforcement cases aimed at investigating the influence of applying FRP composite to timber beams (as listed in column 1). The table classifies external and internal strengthening approaches used to reinforce the timber beams, indicating the timber’s physical and mechanical properties and reinforcement components (as shown in columns 2 to 12). The findings of these studies regarding strength properties reported for each reinforcement approach are listed in column 13.

Table 1. The nine reviewed studies on reinforcement timber beams using FRP materials.

1	2	3	4	5	6	7	8	9	10	11	12	13
Studies	Reinforcement approach	Reinforcement application	FRP composites	EFRP (MPa)	FRP Dimension (mm)	Timber Type	Specimen number	Timber beam’s cross section (HxW) (mm)	Total length (mm)	Total Span (mm)	ETimber (MPa)	Improvement %
1. (Subhani et al. 2017)	Externally	Flat plates	CFRP	216000	L 1500 x W 45	LVL-F17	3	240 x 45	2400	2100	14000	K = 14%
		U-Wrap sheets			L 1500 x H 120 (2 sides) x W 45							K = 30%
2. (Corradi et al. 2017)	Externally	Sheets	CFRP	417600	W 100 x A 0.082 x2 (layers)	Fir-Wood	17	200 x 200	4000	N/A	N/A	M = 58%
			GFRP	78650	W 100 x A 0.052 x2 (layers)							M = 51%
3. (Thorhallsson et al., 2017)	Externally	Sheets	BFRP	84000	T = 0.8	Glulam	24	65 x167	3200	N/A	N/A	M = 43%
			GFRP	73000	T = 0.8							M = 57%
4. (de la Rosa García, Cobo Escamilla, and Nieves González García 2013)	Externally	U-Wrap sheets	BFRP	84000	T= 0.103	Pine timber	26	78 x 155	1090	1000	N/A	M = 24.5%
			CFRP	230000	T= 0.167							M = 3.46%
5. (Yusof and Rahman 2017)	Externally	Plates	CFRP	165000	W 30 x T 1.4	Pine timber	5	100 x 200	3000	N/A	N/A	M MULT. = 35%
					W 60 x T 1.4							M = 44%
6. (Fossetti, Minafò, and Papia 2015)	Internally	rods	BFRP	N/A	Ø8	Glulam	3	80 x 160	3400	3150	11,000	M = 26%
7. (Chun, Balen, and Pan 2016)	Internally -NSM	Rods	CFRP	90000	Ø6	Pine Wood	20	100 x 150	1700	1500	12204.5	MOR = 13% ~37.4%
					Ø8	Fir-wood					10263	MOR = 13.6%
		Plates		180000	L 30 x T 1.4	Pine Wood					12204.5	MOR = 22%
					L 30 x T 2.8	Fir-wood					10263	MOR = 29%
8. (Nowak, Jasieńko, and Czepiżak 2013)	Internally	plates	CFRP	165000	L 50 x T 1.2	Hardwood	*	120 x 200	4000	N\A	8000	MOR = 8%
9. (Bal 2014)	Internally	Sheets	WGF	N\A	N\A	LVL	72	600 x 600	3000	N\A	N\A	MOR = 10.25%

K: Bending stiffness.

M: Ultimate Bending Moment Capacity.

MOR: Modulus of Rupture.

L: Length.

W: Width

H: Height

A: Area

T: Thickness

* Existing members

The studies listed from 1 to 5, which considered external reinforcement forms such as flat plates, U-warp sheets, and plates, demonstrated significant increases in various strength properties. Promising improvements of 30% and 58% in bending stiffness (K) and ultimate bending moment capacity (M), respectively, were observed using this reinforcement approach. Among these studies, the study by Corradi et al. (2017) using CFRP sheet resulted in a remarkable improvement of 58% in M. Following closely, the experiment conducted by Thorhallsson et al, 2017) achieved a significant enhancement of 57% in M using GFRP sheets on Glulam timber beams. These findings indicate that the external reinforcement application offers substantial improvement in the load-bearing capacity of the structural timber beams.

Conversely, the findings of the studies from 6 to 9 using internal reinforcement approaches, although effective, generally reveal lower improvements in strength properties. The highest improvement using this approach was reported by Chun et al. (2016), who obtained an improvements ranging from 13% − 37.45% in MOR property.

Based on the findings observed in the Table, the external reinforcement approach appears to be an optimal option for enhancing the structural performance of the timber structural beams. Despite this, limited attempts are made to analyse and predict the behaviour of reinforced timber beams using FRP composite material, considering their grain orientations and the simultaneous occurrence of shear, tension, and compression failure modes.

This paper presents the development and application for an efficient and accurate FE simulation of the load-deflection behaviours and damage prediction in timber beams reinforced with an external reinforcement sheet of CFRP composite materials utilising a well exercise reinforcement application method in the industry, which is the U-wrap form. The developed models are validated using empirical data reported by (Subhani et al. 2017; Ardalany et al. 2013) Subhani et al. (2017) and Ardalany et al. (2013) (discussed in the next section). In addition, the proposed model is used to conduct a novel parametric study to evaluate the influence of various lengths of U-wrap CFRP sheets applied to the tension side of LVL beams.

2. Observations from experimental investigations

This section presents a discussion of three literature-surveyed experiment programmes that were conducted for plain and reinforced LVL timber beams. The details of the tests, including their setup, preparation, and observed results, are critically identified for careful consideration in developing a reliable FE model demonstrated in Section 3. Test 1 for one of the unreinforced timber samples examined in this study is provided in Table 2. The experiment was conducted by Ardalany (2012) to evaluate the damage behaviour in the vicinity of holes within LVL timber beams. The specimens presented in Tests 2 and 3 were carried out by Subhani et al. (2017) for plain and reinforced LVL timber beams to evaluate the bending strength improvement gained from reinforcing timber beams with U-shaped CFRP sheet arrangements. Table 2 shows the dimensions and span tests for the reviewed tests, while Table 3 depicts the material properties of the specimens, such as the density of the timber ($\mathit{\rho }$), the modulus of elasticity in bending (E), the bending strength ($f_b^{\prime }$), the tensile strength parallel to the grain ($f_t^{\prime }$), the shear strength of the beam ($f_s^{\prime }$), and the compressive strength parallel to the grain ($f_c^{\prime }$).

Table 2. Characteristics of tests 1, 2, and 3 specimens.

Test	Specimen’s index	Details	Specimens number	Characteristic (Length × Width × Depth (mm))	Span (mm)	References
1	LVLC	LVL Timber Beam	1	1500 × 200 × 45	1300	(Ardalany 2012; Ardalany et al. 2013)
2	LVLC	LVL Timber Beam	3	2400 × 240 × 45	2100	(Subhani et al. 2017)
3	LVLU	LVL Reinforced Timber Beam	3	2400 × 240 × 45	2100

Table 3. Material properties of tests 1, 2, and 3 specimens.

Test	Specimen’s index	$\mathit{\rho }$ (kg/m^3)	$\mathit{E}$ (GPa)	$f_b^{\prime }$ (MPa)	$f_t^{\prime }$ (MPa)	$f_s^{\prime }$ (MPa)	$f_c^{\prime }$ (MPa)
1	LVLC	573	12	48	33	6	12
2	LVLC	650	14	62	34	5.3	47
3	LVLU	650	14	62	34	5.3	47

The test specimens were bent using universal loading frames and acquisition systems to capture load-displacement measurements. The three-point bending test setup was used. Figure 1 shows the test setups and cross-section preparation for the tests reviewed. Test 1 followed the AS/NZS 4063.2 standard (Australia and New Zealand Standards 2010), while Tests 2 and 3 were performed in linewith the Australia Standard (AS 1720.1 2010).

The beams were loaded until they failed at room temperature. Sections 2.1 and 2.2 present the beams’ strength and serviceability performances and their failure modes, including crack damage initiation and development. The loading system used a load cell with an accuracy of 0.1 kN, and midspan deflection was measured with an accuracy of ±0.1 mm. Test 1 used lateral supports at the beam ends to prevent lateral buckling behaviour. In Tests 2 and 3, two L-shaped steel braces were installed at the midspan to prevent movement or local buckling during load application, as shown in Figs. 1(a)-(b).

Figure 1. Experimental setup of the three test specimens: (a) tests 1 (Ardalany 2012) and 2 (Subhani et al. 2017) plain timber beams set up, (b) test 3 reinforced timber beam set up (Subhani et al. 2017).

2.1. Bending performance of plain LVL timber beam

The plain timber specimens tested in bending Tests 1 and 2 showed identical crack development and failure modes. Three common failure mechanisms were observed, which are illustrated in a schematic diagram in Figure 2. These mechanisms consist of bottom tensile damage, diagonal shear cracks, and top compressive crushing damage, as reported by Ardalany (2012) and Subhani et al. (2017). The load-deflection responses of the beams tested in Tests 1 and 2 are presented in Figures 3 and 5, respectively. Both tests show predominantly linear behaviour in the plotted results up to the point where the beams fracture. At this load, the beams typically show the initiation of crack damage in the tension region at mid-span. These failure modes and behaviours are important for evaluating the performance of the FE model presented in Section 3.

Figure 2. Failure modes observed in the tested beams: (a) three failure modes found in the structural beam of test 1 (Ardalany 2012), (b) splintering tensile failure found in the structural beam of test 2 (Subhani et al. 2017), and (c) Cross-grain tensile failure (b) splintering tensile failure found in the structural beam of test 3 (Subhani et al. 2017).

Figure 3. Load-deflection curve of tests 1 (Ardalany 2012).

2.2. Bending performance for reinforced timber beams

Test 3, in which timber specimens were reinforced with a U-shaped CFRP sheet, demonstrated a significant increase in ultimate load-carrying capacity, stiffness, and displacement ductility when compared to the non-reinforced beams in Test 2. Table 4 shows the improvement percentages of these characteristics when applying a U-shaped CFRP sheet.

Table 4. Improvement percentages in ultimate load carrying capacity, stiffness, and displacement ductility of tests 2 and 3 compared to non-reinforced beam.

Parameters	Unit	Test 2	Test 3	Improvement %
Ultimate load-carrying capacity	kN	42.03^1	52.6^1	25%
Stiffness	kN	1.48^1	1.78^1	19.8%
Displacement Ductility	mm	1.18^1	1.58^1	39%

¹The estimated values were found based on the data provided by Subhani et al. (2017).

Table 5 displays the material properties of the CFRP sheet used for reinforcement, including the Modulus of Elasticity (E) and the mean tensile strength (T_s). The rupture strain (R_S) is also provided, along with the length and thickness of the CFRP sheet, represented by L and T, respectively. Sikadur 330 adhesive was used to fix the CFRP sheet onto the timber beam. Table 6 illustrates its key features, where F_T represents tensile strength, E_F denotes flexural modulus of elasticity, and E_T stands for tensile modulus of elasticity. The reinforced beam displayed a failure in cross-grain tension at its mid-span, as illustrated in Figure 4. The failure took place perpendicular to both the direction of the LVL timber beam’s grains and the direction of the CFRP fibres, occurring along the weak plane of the LVL timber beam (Subhani et al. 2017).

Table 5. CFRP characteristic.

CFRP - Type	$\mathit{E}$(GPa)	Ts (MPa)	RS (MPa)	L (mm)	T (mm)	Reference
Unidirectional sheet - Sika-Wrap 230 C	216	3176	1.8%	1500	0.131	(Subhani et al. 2017)

Table 6. Adhesive characteristics.

Adhesive - Type	FT (MPa)	EF (MPa)	ET (MPa)	Reference
Sikadur 330	30	3800	4500	(Subhani et al. 2017)

Figure 4. Failure mode observed in the reinforced beam (Subhani et al. 2017).

The plot in Figure 5 illustrates the load-deflection curve of the beam, which exhibits a linear pattern up to the yield point of 38.5 kN. Beyond this point, the beam demonstrates non-linear behaviour until it reaches the average fracture point at a load of 51 kN, where cross-grain tensile cracks appear.

Figure 5. Load-deflection curves of tests 2 and 3 (Subhani et al. 2017).

In Section 3, three FE models have been created for LVL timber beams using the testing arrangements and material properties listed in Tables 2–3. These models were developed using a sophisticated FE simulation approach that can accurately determine various properties of the beams, such as ultimate load capacity, fracture load, stiffness, ultimate bending moment capacity, and displacement ductility. Additionally, the models are capable of predicting the damage failure modes that may occur in timber structural elements.

3. A refined finite element model for reinforced timber

The literature presents several FE models that can replicate the strength and serviceability behaviour of structural timber members under different loading conditions. However, advanced FE simulations of structural timber elements require complex modelling arrangements (Fragiacomo 2005), such as the creation of user-subroutines or numerous computational processes, particularly when producing three-dimensional solid models (Ardalany 2012; Janssens 2017; Khelifa et al. 2015; Mirianon, Fortino, and Toratti 2008; Sandhaas and Van de Kuilen 2013). To illustrate, Figure 6 depicts a 3D solid model developed using a cohesive element approach to simulate material behaviour under damage. This modelling process is time-consuming and demands a large number of elements and partitions to ensure an effective simulation. Moreover, selecting the appropriate mesh size for the analysis is among the difficulties encountered when utilising this modelling approach (Camanho and Hallett 2015).

Figure 6. The developed 3D solid model using a cohesive element approach.

This paper presents a computationally efficient three-dimensional FE model, which employs shell-type elements available in the ABAQUS standard software 2018 to simulate complex materials such as timber in less time. The shell model used Hashin’s damage model, which is a widely used simulator in the research area for composite structures (Al Abadi et al. 2018; Chybiński and Polus 2021; Khorsandnia, Valipour, and Crews 2013; Ramage et al. 2017). This theory has recently been employed to model timber elements due to their anisotropic behaviour (Chybiński and Polus 2021; Khorsandnia, Valipour, and Crews 2013; Stenrud Nilsen 2015).

The refined FE model is utilised to simulate the three testing arrangements discussed in Section 2 (Table 2) by simulating the applied bending actions for the three tests. The predicted load-deflection behaviours for each test are then critically evaluated to determine the validity of the proposed refined model. Figure 7(a) illustrates the schematic representation of the developed FE models for the non-reinforced timber Tests 1 and 2, while Figure 7(b) depicts the modelling arrangements for applying the U-shaped CFRP component. The geometrical properties of the three tested models are presented in Table 7.

Figure 7. (a) Schematic presentation of FE model for pure timber tests 1 and 2, (b) geometric figures of the reinforced timber test model.

The timber and CFRP components are modelled using a 4-node general-purpose shell-type S4R. This element type is a deformable planar element that utilises reduced integration and hourglass control for minor membrane strains. According to the Abaqus manual ABAQUS Manual (2012), using reduced integration and hourglass capabilities is more effective when simulating three-dimensional FE models instead of utilising solid models. It is also an appropriate option for extensive strain analysis, explicit finite membrane strains, and arbitrarily large rotations, all of which are expected in this investigation (ABAQUS Manual 2012).

Both parts, namely timber and CFRP, are represented using the geometrical element shape of Quad, which includes thick shells utilising the hourglass control formulation with reduced integration. Hence, this allows using three integration points element direction due to the type of the geometrical element shape for the node, Quad shape (ABAQUS Manual 2012). This formulation enhances the accuracy of the material’s mesh and non-linear response at significant deformation levels (ABAQUS Manual 2012).

Figure 8 shows a simple presentation of the geometrical shape nodes used in this model for both parts.

Figure 8. 4-nodes geometrical shape for both timber and CFRP parts.

Table 7. Geometric properties of test models.

Parts	Geometric characteristics	Model 1 (mm)	Model 2 (mm)	Model 3 (mm)
Timber	Span	1300	2100	2100
	Total Length	1500	2400	2400
	Width	200	240	240
	Depth	45	45	45
U-Shape CFRP Sheet	Total Length	N/A	N/A	1500
	Width			45
	Edge Length			120
	Thickness			0.13

This method was used to ensure the accuracy of the analysis results. In addition, a fine mesh consisting of 1140 elements with an element size of 30 mm was selected for both parts. Sensitivity analyses for the element size were conducted for the adopted FE model across a range of element sizes from 20 mm to 40 mm (consisting of 2514 to 626 elements, respectively). Over the evaluated range, the difference in predicted values of ultimate strength was below 0.0037%. Consequently, a 30 mm element size was adopted for this modelling work.

The orientational characteristics of the timber grain and CFRP fibres necessitated the use of a composite layup tool in ABAQUS for modelling each part. The timber component was modelled using a single 45 mm thick ply (to represent timber beam width), while the CFRP component was modelled using a 0.13 mm thick ply, which is the thickness of the reinforcement sheet (Subhani et al. 2017).

The beam model with tie constraints, which was used to replicate the testing environment, is illustrated in Figure 9. Restraints were enforced at the beam’s supporting ends to limit rotation and transition degrees of freedom. Additionally, a steel bracket was set up at the loading point (mid-span) to provide extra tie constraints and prevent rotation during loading, as shown in Figure 9.

Figure 9. Boundary conditions of the reinforced model and the application of “tie” constraints.

The elastic properties of the model were determined using the elastic lamina material type found in (ABAQUS Manual 2012), which is a commonly used method in research to achieve suitable anisotropic elastic behaviour for composite materials and timber (Al Abadi et al. 2018; Chybiński and Polus 2021). Table 8 provides the elastic properties used for both the timber and CFRP parts in the model.

Table 8. Elastic properties of the test models.

Parts	ELASTIC PROPERTIES	Units	Symbol	MODEL 1	MODEL 2	MODEL 3
Timber	Elastic modulus parallel to grain direction	MPa	E1	12000^1	14000^4	14000^4
	Elastic modulus perpendicular to grain direction	MPa	E2	485^1	485^1	485^1
	Longitudinal Poisson’s Ratio	-	V12	0.38^1	0.4^4	0.4^4
	Longitudinal shear modulus	MPa	G12	518^1	160^2	160^2
	Shear modulus	MPa	G13	518^1	160^2	160^2
	Rolling shear modulus	MPa	G23	96^3	25^1	25^1
CFRP	Elastic modulus parallel to grain direction	MPa	E1	N/A	N/A	216000^5
	Elastic modulus perpendicular to grain direction	MPa	E2			26000^3
	Longitudinal Poisson’s Ratio	-	V12			0.3^1
	Longitudinal shear modulus	MPa	G12			160^2
	Shear modulus	MPa	G13			160^2
	Rolling shear modulus	MPa	G23			24^1

Adopted resources.

¹Based on (Ardalany 2012)

²Based on assumptions and estimations made from multiple resources to match the simulation and the experiment results (Christoforo et al. 2014; Sandhaas and Van de Kuilen 2013; Ardalany et al. 2013).

³Based on (Al Abadi et al. 2018)

⁴Based on (Subhani et al. 2017)

⁵Based on the manufacturer datasheet (Stenrud Nilsen 2015)

To simulate the interlaminar damage mechanisms leading to the ultimate collapse, Hashin’s damage initiation theory (ABAQUS Manual 2012; Chybiński and Polus 2021) is employed. In Abaqus, the HSNFTCRT, HSNFCCRT, HSNMTCRT, and HSNMCCRT parameters are utilised to depict the Hashin damage initiation criterion for the fibre in tension, fibre in compression, matrix in tension, and matrix in compression, respectively (ABAQUS Manual 2012). The expressions for these four criteria, adopted from (ABAQUS Manual 2012), are given as follows:

Criterion for fibre tensile initiation (HSNFTCRT) $\left({\sigma }_{11}\ge 0\right)$:

(1) $F_f^t={\left(\frac{{\mathit{\sigma }}_{11}}{X^T}\right)}^2+\mathit{\alpha }{\left(\frac{{\mathit{\tau }}_{12}}{S^L}\right)}^2$(1)

Fibre compressive initiation criterion (HSNFCCRT) $\left({\sigma }_{11}<0\right)$:

(2) $F_f^c={\left(\frac{{\mathit{\sigma }}_{11}}{X^C}\right)}^2$(2)

Matrix tensile initiation criterion (HSNMTCRT) $\left({\sigma }_{22}\ge 0\right)$:

(3) $F_m^t={\left(\frac{{\mathit{\sigma }}_{22}}{Y^T}\right)}^2+{\left(\frac{{\mathit{\tau }}_{12}}{S^L}\right)}^2$(3)

Matrix compressive initiation criterion (HSNMCCRT) $\left({\sigma }_{22}<0\right)$:

(4) $F_m^c={\left(\frac{{\mathit{\sigma }}_{22}}{2S^T}\right)}^2+\left[{\left(\frac{Y^c}{2S^T}\right)}^2-1\right]\frac{{\mathit{\sigma }}_{22}}{Y^C}+{\left(\frac{{\mathit{\tau }}_{12}}{S^L}\right)}^2$(4)

Where X^T, X^C, Y^T, and Y^C refer to longitudinal tensile strength along grain direction, longitudinal compressive strength perpendicular to grain direction, transverse tensile strength along grain direction and transverse compressive strength perpendicular to grain direction, respectively. S^L and S^T represent the longitudinal shear strength along the grain direction and transverse shear strength perpendicular to the grain direction, respectively (ABAQUS Manual 2012; Stenrud Nilsen 2015). The coefficient α, which is 1 in this study, is a factor that determines the contribution of shear stress to the fibre tensile initiation criterion. It can be 0 or 1 (ABAQUS Manual 2012). The effective stresses for the timber and CFRP parts in the primary directions, σ₁₁, σ₂₂, and σ₁₂, are utilised to evaluate the damage initiation criteria. It is noteworthy that these stresses act on the undamaged material (ABAQUS Manual 2012)

The shear stresses incorporated in the HSNMTCRT and HSNMCCRT failure criteria are utilised to anticipate matrix tension and compression cracking that could result in matrix shear failure (Al Abadi et al. 2018). Table 9 presents a collection of parameters applied in Abaqus to calculate the damage variables for the modelled components mentioned above.

Table 9. Damage characteristics of the reinforced timber beam model.

	Damage characteristics		Unit	MODEL 1	MODEL 2	MODEL 3
Timber Part	Hashin Damage	longitudinal Tensile strength	MPa	41.9^1	27^5	27^5
		longitudinal compressive strength	MPa	50.3^1	50.3^1	50.3^1
		Transverse tensile strength	MPa	10^1	10^1	10^1
		Transverse compressive strength	MPa	15^1	15^1	15^1
		longitudinal shear strength	MPa	10^1	10^1	10^1
		Transverse shear strength	MPa	5^1	5^1	5^1
	Damage Evolution	longitudinal Tensile Fracture Energy	N/mm	45^1	120^2	120^2
		longitudinal Compressive Fracture Energy	N/mm	45^1	130^2	130^2
		Transverse Tensile Fracture Energy	N/mm	0.1^1	0.1^1	0.1^1
		Transverse Compressive Fracture Energy	N/mm	0.1^1	0.1^1	0.1^1
	Damage Stabilization	Viscosity - longitudinal Tensile	MPa	0.01^2	0.01^2	0.01^2
		Viscosity - longitudinal Compressive	MPa	0^2	0^2	0^2
		Viscosity - Transverse Tensile	MPa	0^2	0^2	0^2
		Viscosity - Transverse Compressive	MPa	0^2	0^2	0^2
CFRP part	Hashin Damage	longitudinal tensile strength	MPa	N/A	N/A	400^2
		longitudinal compressive strength	MPa			570^3
		Transverse tensile strength	MPa			400^2
		Transverse compressive strength	MPa			570^3
		longitudinal shear strength	MPa			29^4
		Transverse shear strength	MPa			29^4
	Damage Evolution	longitudinal Tensile Fracture Energy	N/mm			1200^2
		longitudinal Compressive Fracture Energy	N/mm			1300^2
		Transverse Tensile Fracture Energy	N/mm			0.1^1
		Transverse Compressive Fracture Energy	N/.mm			0.1^1
	Damage Stabilization	Viscosity - longitudinal Tensile	MPa			0.01^2
		Viscosity - longitudinal Compressive	MPa			0^2
		Viscosity - Transverse Tensile	MPa			0^2
		Viscosity - Transverse Compressive	MPa			0^2

Adopted resources.

¹Based on (Zhang et al. 2022)

²Based on assumptions and estimations made from multiple resources to match the simulation and the experiment results (Christoforo et al. 2014; Janssens 2017; Sandhaas and Van de Kuilen 2013;and Ardalany2012).

³Based on the manufacturer’s datasheet adopted from (Subhani et al. 2017)

⁴Based on (Al Abadi et al. 2018)

⁵Based on (Futurebuild 2020).

The DAMAGEFT, DAMAGEFC, DAMAGEMT, and DAMAGEMC parameters in ABAQUS depict the damage modes resulting from tensile fibre damage, compressive fibre damage, matrix tensile damage, and matrix compressive damage, respectively. These parameters vary in value from zero, indicating no damage, to 1.0, representing complete damage (Al Abadi et al. 2018).

3.1. Discussion of FE simulation results and comparison with experiments

The implemented FE analyses and predicted results are used to plot the load-deflection behaviour and damage development modes. In the subsequent sections, the predicted results are examined and assessed in comparison to the experimental observations listed in Table 2 for the tested plain and reinforced timber cases.

3.1.1. Plain timber beam models

The developed FE model predicted various failure damages for the plain timber beam case (Tests 1 and 2). The first type of damage involved a tensile crack that emerged at the midspan of the beam’s soffit when the load reached around 41.57 kN and 21.85 kN for Tests 1 and 2, respectively (Depicted in Figure 10). Figs 11(a)-(b) illustrate the predicted failure modes of the FE models for Tests 1 and 2 at tension areas, respectively. The models effectively simulated the initial tension cracks observed in the plain timber beams. Furthermore, the FE model can simulate the compressive failure damage that occurs at the beam’s top, as shown in Figure 12. The simulation predicted that the compressive failure damage for Test 1 started at a load of 48.31 kN, which corresponded well with the experimental result. However, the compressive failure damage for Test 2 occurred at a load of 35.41 kN, slightly lower than the experimental result. Nevertheless, these two predictions are in good agreement with the experimental data plotted in Figure 10, which showed the crack initiation at the end of the load-displacement curves’ elastic range.

Figure 10. Load vs. deflection for test 1 and test 2 models.

Figure 11. (a) Tensile failure damage observed in FE beam model for test 1 (b) tensile failure damage observed in FE beam model for test 2.

Figure 12. (a) Compressive failure observed in FE beam model for test 1 (b) compressive failure observed in FE beam model for test 2.

It is worth mentioning that the created models have effectively captured the bending characteristics of the beam during the phase of elastic deformation, displaying a small indication of non-linear behaviour before the load of fracture is reached (Figure 10). In general, the minor variations in the outcomes arising from the use of mechanical properties in the experiment that were marginally distinct from those presumed in the FE model.

3.1.2. FE reinforced timber beam model

The load-deflection behaviour obtained from the reinforced timber FE model accurately represents the results of Test 3, as illustrated in Figure 13. The curves closely resemble the linear elastic behaviour, which records the yielding load at 40.2 kN with approximately 23 mm deflection and accurately depicts the non-linear behaviour up to the fracture point at 50 kN and 41.9 mm deflection. The failure modes predicted by the FE model are shown in Figs. 14(a)-(b). As indicated by the experimental data, the FE model anticipated the initial tensile failure damage at the midspan of the beam (shown in Figure 14(a)) at a load of 38.83 kN, which is believed to contribute to initiating the non-linear behaviour. Furthermore, the model presents the development of damage at the top, compressed side of the beam, where the damage is predicted to initiate at 43.12 kN, the point where the beam starts to deflect non-linearly, and ultimately predicts the fracture failure damage at the 50 kN load.

Figure 13. Predicted load-deflection curve of FE test 3 model.

Figure 14. (a) Tensile failure damage (b) compressive failure damage.

3.2. Evaluation of FE model predictions

In order to assess the precision of the FE model in replicating the load-deflection behaviour of the three beams that were tested, five strength properties and displacement ductility values were calculated from the experimental results and compared to those estimated from the FE model results. The strength properties evaluated are listed in Table 10, including bending moment capacity (M), stiffness (K), fracture load (F), maximum deflection at fracture load (D_Max), and displacement ductility (µ). The results of the predicted FE Test 1 beam were marginally less than the experimental results for M, K, F, D_Max, and µ, with an error percentage of less than 10%. The reason for this discrepancy was due to certain assumptions made in the FE analysis. These assumptions included the longitudinal shear modulus (G₁₂), shear modulus (G₁₃), and damage stabilisation parameters, such as the viscosity coefficient in longitudinal tensile, longitudinal compressive, transverse tensile, and transverse compressive directions. These factors likely contributed to the higher error percentage for F.

The M, K, F, D_Max, and µ, properties of the predicted FE Test 2 model were slightly lower than the experimental results, with an error percentage of about 7%, indicating precision in simulating the plain timber beam. In the case of the reinforcement beam (Test 3), the experimental and FE model predictions were in good agreement, with an error percentage reported as low as 1.8%, except for µ, which had an error percentage of around 13% due to certain assumed properties used in the FE analysis, such as timber and CFRP shear modulus, timber longitudinal tensile and compressive fracture energy, damage stabilisation parameters (for both timber and CFRP parts), and CFRP longitudinal tensile and compressive strength. With the presented FE model’s high level of accuracy in simulating the performance of fibre-reinforced timber beams, it was utilised for a parametric study on timber beams reinforced with varying lengths, as presented in Section 4.

Table 10. Strength properties of the tested models and test results.

1	2	3	4	5	6	7	8	9	10
Parameters (Unit)	Test 1	FE Test 1 model	Error %	Test 2	FE Test 2 model	Error %	Test 3	FE Test 3 model	Error %
M (kN.m)	18.14	16.58	−8.6%	21.02	19.90	−5.37%	26.72	26.25	−1.8%
K (kN/mm)	3.68	4.0	−8%	1.48	1.5	1.53%	1.78	1.77	0.4%
F (kN)	55.8	51.0	−9%	40.05	37.9	−5.37%	50.9	50	−1.8%
DMax (mm)	0.007	0.006	−8.6%	0.010	0.010	0%	0.013	0.013	0%
µ (mm)	1.17	1.22	4.2%	1.18	1.1	7%	1.58	1.82	−13%

4. Parametric study

The study involved a parametric FE analysis of timber beams reinforced with different lengths of U-shaped CFRP to determine the effects on the beam’s strength and serviceability. Six lengths were considered for a 2400 mm span LVL timber beam demonstrated in Figure 15 and listed in Table 11. The study used seven FE beam models developed from a refined finite element model, including an unreinforced timber beam specimen referred to as TB in Table 11. The other models were reinforced with varying lengths of U-shaped CFRP attached to the tension (lower) side of the beams, denoted as TCB_L_CFRP in Table 11, where L_CFRP indicates the adopted reinforcement length. The geometric properties of the timber and CFRP parts are presented in Table 11, and the FE test arrangements used in the study are depicted in Figs. 15(a)-(b).

Figure 15. (a) Geometry of TB model, (b) geometry of TCB_L_CFRP model, and cross-section of U-Shape CFRP.

Figure 16 depicts the load-deflection curves obtained from the simulations. These curves were used to determine five parameters for evaluating the performance of composite beams with progressively longer reinforcement lengths. These parameters include the ultimate load capacity (Fy) at the yield point, fracture load (F), stiffness (K), ultimate bending moment capacity (M), and ductility parameter (µ). To compare the results with the unreinforced specimen case (TB-FE model), these parameters were also predicted for that scenario.

Figure 16. Load-midspan deflection curve of the simulated beams used in the parametric study.

Table 11. Geometric properties for the FE beam Model.

FE Beam Model	Total length (mm)	Span (mm)	Width (mm)	Deep (mm)	CFRP Lengths (mm)	CFRP Side’s Length (mm)
TB	2400	2100	240	45	N/A
TCB_300CFRP					300	120
TCB_600CFRP					600
TCB_900CFRP					900
TCB_1200CFRP					1200
TCB_1500CFRP					1500
TCB_1800CFRP					1800

4.1. CFRP reinforcement length impact on strength and serviceability

In Figure 17, the Fy values are graphed to display the maximum load capacity of the beams’ linear load-deflection behaviour as a function of the reinforcement length and the improvement ratio (IR%) function.

Figure 17. Fy values and IR% of the simulated beams.

The study found that the ultimate load to failure (Fy) of the reinforced beam increased progressively as the reinforcement length increased from 300 mm to 1500 mm. The Fy of the beam reinforced with a 300 mm CFRP increased only by 0.3% compared to the un-reinforced beam (TB) but increased by 6%, 12%, 15%, and 17% for reinforcement lengths of 600 mm, 900 mm, 1200 mm, and 1500 mm, respectively. Although the Fy capacity only showed slight improvement when the reinforcement length was increased to 1800 mm compared to 1500 mm, the fracture load (F) capacity of the composite beams improved in line with the findings for Fy. The F value is the point at which damage initiates in the composite section materials, as determined from the load-deflection curves shown in Figure 16. The F influences by the level of reinforcement lengths are plotted in Figure 18, which again indicates that the reinforcement length ranging from 300 mm to 1200 mm mainly improves the F capacity in line with the findings for Fy.

Figure 18. F values and IR% of the simulated beams.

The values of K predicted through analytical methods by analysing the elastic region of the load-deflection curves. Figure 19 demonstrates how the reinforcement lengths affect these parameters. The graph provide evidence that reinforcement lengths between 300 mm to 1500 mm play an important role in enhancing the K of the beam.

Figure 19. Flexural stiffness and stiffness improvement ratio of the simulated models.

Furthermore, the values of M were predicated through additional analysis of the load-deflection curves using the ultimate strength value, which represents the maximum load that the simulated beams can reach before fracturing (Depicted in Figure 16). These M values are displayed in Figure 20, which also indicates that reinforcement lengths between 300 mm to 1500 mm effectively enhance the M value of the beam.

Figure 20. M values and IR% of the simulated models.

Figure 21 depicts the relationship between the influence ratio (IR) and displacement ductility parameter (µ) based on the reinforcement lengths. The graph indicates that the improvement in the influence ratio is mainly dominated by reinforcement lengths between 300 mm and 900 mm.

Figure 21. Μ values and IR% of the simulated beam models.

After analysing the performance of different parameters, it has been concluded that the Fy , F, K, and M parameters are not significantly affected by reinforcement lengths below 300 mm or above 1200 mm. To present the Reinforcement Ratio (RR) of the composite beam under investigation, the reinforcement lengths have been normalised to the overall length of the beam, which is 2400 mm. The non-effective RR falls within the ranges of 0 to 12.5% and above 50%. Figure 21 shows that despite limited improvements in other parameters at a reinforcement length of 300 mm, there is a significant 95% improvement in the ductility parameter. This indicates the importance of reinforcement in increasing the ductility behaviour of the beams by more than two folds.

Table 12 provides mathematical expressions for the five evaluated strength parameters as a function of the Reinforcement Ratio (RR) using the refined FE simulation method presented in this paper. The expressed relationships listed in Table 12 are limited to the timber beam size, range of reinforcement lengths (300 mm and 1200 mm) and the bending boundary conditions adopted in the conducted FE simulations (Section 3). The expressions can be refined to predict the strength properties comprehensively.

Table 12. Mathematical expression of strength properties to RR relationships.

Strength property	Mathematical expression over 300–1200 reinforcement range	RR range	Equation number
Fy	Fy = 13 RR + 31	where 0.125 ≤ RR ≤ 0.5	5
F	F = 21.5 RR + 39	where 0.125 ≤ RR ≤ 0.5	6
K	K = 2 RR +1.3	where 0.125 ≤ RR ≤ 0.5	7
M	M = 6.8 RR + 16	where 0.125 ≤ RR ≤ 0.5	8
µ	µ = 0.6 RR + 1.9	where 0.125 ≤ RR ˂ 0.5	9

5. Conclusions

This paper has introduced a novel computational FE modelling approach to develop a design tool for predicting the strength and serviceability behaviour of plain and reinforced timber structural members with high accuracy and reduced computational time. The model employed shell-type elements with a composite layup tool and Hashin’s damage model available in ABAQUS to simulate the damage mechanisms that cause a fracture in structural members. To validate the performance of the model, it was critically evaluated with three experiments reported in the literature, and a parametric study was conducted to optimise the reinforcement ratio of LVL timber beams reinforced with various lengths of U-shaped CFRP.

The following conclusions can be drawn from this computational modelling work:

• The FE model predicts the load-deflection curves and strength properties of reinforced timber beams with strong agreement with experimental results, with differences attributed to some assumed properties used in the simulation process.

• The developed model successfully predicts various failure modes experienced under flexure loading, such as tensile failure, compressive failure, splintering tensile fracture, and shear damage. The model is also found to accurately estimate each failure mode’s initiation with its specific time occurrence.

• The parametric study shows that the reinforcement lengths between no-reinforcement and 300 mm and above 1200 mm have minimal impact on the strength parameters presented by Fy, F, K, and M. Therefore, the non-effective RR is defined within 0 to 12.5% and more than 50%. However, the ductility parameter shows a significant improvement with a reinforcement ratio of only 12.5%, which differs from the improvement percentages found in the other parameters using the same RR.

• The paper also provides mathematical expressions to approximately estimate the Fy, F, K, M, and µ based on the RR ranging between 12.5% − 50%, which can significantly assist in designing timber structural members strengthened with a U-shaped CFRP sheet.

Acknowledgments

The research described in this paper was financially supported by La Trobe University, Australia.

Disclosure statement

No potential conflict of interest was reported by the author(s).