On the high fracture toughness of wood and polymer-filled wood composites – Crack deflection analysis for materials design

Cracks oriented in the toughest direction across the grain of wood (0 ◦ ) tend to deflect at 90 ◦ to the precrack rather than extending in 0 ◦ direction. Fracture toughness data across the grain are therefore difficult to interpret. Crack growth mechanisms and effects from replacing wood pore space with a polymer are investigated. Crack growth is analyzed in four-point bending fracture mechanics specimens of birch and two different polymer-filled birch composites using strain-field measurements and finite element analysis (FEA). Calibrated cohesive zone models in both pre-crack and 90 ◦ -directions describe fracture process zone properties in orthotropic FEA-models. Conditions for 0 ◦ crack penetration versus 90 ◦ crack deflection are analyzed based on cohesive zone properties. Stable, subcritical crack deflection takes place at low load, reduces crack tip stress concentration, and contributes to high structural toughness, provided the 90 ◦ toughness is not too low. Polymer-filled neat birch composites have the best structural toughness properties in the present investigation, since 90 ◦ toughness is not compromised by any chemical treatment.


Introduction
Trees have evolved by competing with other plants for sunlight, and height is an important advantage.During storms, the tree stem needs to provide structural support for the plant.The stem is subjected to bending loads, for which the structure is well adapted (axial grain direction, high microfibril angle in fibers).Fracture surfaces of broken stems offer interesting insight with respect to wood failure mechanisms.Although ultimate fracture is by crack growth across the grain, the cracking pattern is complex with local transverse cracking, jagged fracture surfaces, high roughness at multiple scales and very high energy absorption for the tree stem structure.If only partial failure takes place the tree may still survive to reproduce for many years through its seeds.
In wood, opening mode I cracks initially oriented perpendicular to the fiber direction (LT (longitudinal-tangential) and LR (longitudinal-radial)), see Fig. 1, often deviate in the perpendicular direction (90 • ).This crack deflection phenomenon may inspire materials design of strong and tough man-made composites, including polymer-filled wood composites."Weak interfaces" have been suggested to improve toughness also in other biological composites (nacre [1]) by related mechanisms.The proposed mechanisms for wood are seldom quantified and backed up by continuum mechanics models, verifying strain field measurements etc., as in the present investigation.
Here we consider neat wood, but also polymer-filled wood composites, since they offer many possibilities in terms of microstructural, nanostructural, and even molecular scale tailoring [2,3].In these composites, native wood in original form, or chemically modified is used as the reinforcement phase.The porous wood reinforcement is impregnated by monomer so that after polymerization, the pore space is replaced by the polymer.
In previous wood fracture toughness investigations, most studies were carried out on neat wood ("clear wood" without knots).Wood is a porous, orthotropic composite, where mechanical properties scale with relative density (volume fraction of solid cell wall tissue) [4].Ashby (1985) [5] showed that fracture toughness, K c , in LT direction scales quadratically with relative density (volume fraction solid tissue).In an earlier study, Jeronimidis [6] investigated LT fracture toughness and crack growth in sitka spruce and teak.A major problem is that the induced LT precrack tends to deviate 90 • during growth.Ashby et al. [5] pointed out that K C (determined at maximum specimen load) along and across the fiber direction differs by a factor of ten between the tough LT (across the fiber direction) and weaker TL (parallel to the fiber direction) fracture planes.For LT oriented precracks, strain-field measurements have suggested that there are transverse tensile strains at the crack tip [7].The so-called Cook-Gordon mechanism for general crack deflection at interfaces depend on transverse tensile stresses in front of the crack tip, even in isotropic materials [8].
There are quite a few experimental studies reporting LT fracture toughness of wood for crack growth initiation, and the crack path nearly always deviates 90 • during crack growth.Typical values of critical strain energy release rate, G c , range 1-20 kJ/m 2 [4,[9][10][11], and fracture toughness, K c , ranges from 0.1 to 10 MN/m 3/2 .Linear elastic fracture mechanics (LEFM) is the dominating approach, although several common LEFM assumptions are invalid for wood and polymer-filled wood, e.g., small scale yielding and isotropic properties.Tattersall and Tappin did not use fracture mechanics, but developed a specific specimen geometry were the total energy to break a specimen with a pre-crack was measured [8].Jeronimidis [7] successfully used this specimen for wood and reported a work of fracture for sitka spruce of 16 kJ/m 2 .In a previous study, we used edge-grooved fracture mechanics specimens to force actual mode-I LT crack growth penetration.G c for neat birch was 8 kJ/m 2 [12] For polymer-filled wood composites, there are very few reports on LT fracture toughness [12,13], with G c for birch/PMMA reported to be 16 kJ/m 2 .Toughness from impact bending tests have been reported for polymer-filled wood composites [14][15][16], which is not a well-defined material property but strongly dependent on specimen geometry and test conditions.
Bao and Zhou discussed limitations with the LEFM-approach; they specifically pointed out that the powerful small-scale yielding assumption in LEFM means that mechanisms of deformation in the "yield" zone tend to be ignored [17].For composites this can mean missed opportunities for materials design, such as microstructural organization and tailoring of constituent properties.Furthermore, a common simplification in wood fracture toughness analyses is that orthotropic wood materials are assumed to be isotropic.Orthotropic analysis indeed contributes a better fundamental understanding of wood material fracture toughness.An important orthotropic analysis result for wood is that the 90 • deflection of an LT oriented crack correlates with large transverse (or tangential) tensile strains, ε T , ahead of the crack tip and perpendicular to precrack orientation [7].Crack initiation occurs already at ≈40-50 % of the specimen peak load at which the specimen becomes more compliant due to several cracks growing at 90 • , parallel to the grain direction [18].There are several examples of orthotropic models developed for crack growth in wood [19,20].Landis reviews finite element-based methods, including crack bridging approaches [19] where properties of a fracture process zone ahead of the crack tip is included.Matsumoto and Nairn are discussing limitations of LEFM-approaches [20] and presenting experimental data.Both refs [19,20] focus on crack growth with tensile load in tangential or radial direction (weakest directions); more recent work on low toughness crack growth by Nairn and coworkers combines crack bridging models with strain-field measurements and FEM [21].Recently, general and very sophisticated methods for crack growth analysis were presented [22,23].Füssl and coworkers used variational analysis and minimization of energy to predict the crack path [22].Again, weak direction crack growth was the example selected; realistic zig-zag crack patterns influenced by the growth rings were predicted.The reason for the focus on weak direction crack growth is its importance in engineering design.In a previous investigation, we addressed this problem for the materials in the present study [24].Growth of LT oriented precracks, however, is not well understood.Better understanding may contribute to materials design guidelines for wood Fig. 1.Illustration of a cross section from a crown cut wood veneer with an arrow indicating the fiber direction, longitudinal direction (L).Sections a) and b), indicated with dashed lines on the veneer, illustrates a crack oriented perpendicular to the fiber direction for longitudinal-radial (LR) and longitudinal-tangential (LT), respectively.
E. Jungstedt et al. composites, and possibly improved toughness.
The main objective is to investigate crack growth conditions for 90 • crack deflection of LT-precracks in wood and measure effects from the polymer phase in polymer-filled wood composites.The phenomenon is favorable for the work of fracture of tree stems but may under certain conditions reduce work of fracture in other wood and wood composite structures.A better understanding of the phenomenon may explain scale effects and could contribute to tailoring of new wood composites.Three materials were included in the study: native birch, native birch impregnated with poly(methyl methacrylate) (Birch/PMMA), and delignified birch impregnated with PMMA (D-Birch/PMMA).The delignified wood reinforcement has chromophores removed to make transparent wood biocomposites [25].
We will use cohesive zone models (CZM) combined with orthotropic continuum models for wood, finite element modeling (FEM) and strain field measurements around the crack tip by digital image correlation (DIC) techniques.These CZMs represent fracture properties of the materials expressed as softening function, critical cohesive stress, limiting opening displacement, and fracture energy.Such models may be used to analyze development of the fracture process zone (FZP), crack growth onset and the crack path.The present model aims to analyze 90 • crack deflection.It therefore includes an LT-precrack and two CZMs positioned in two possible crack growth directions (continued crack propagation along LT and crack deflection parallel to grain along TL).For accuracy, FEM strain fields are fitted to DIC strain fields by a FEM updating (FEMU) routine, including orthotropic linear elastic constants for the materials and CZM fracture properties.For wood, FEMU has been used to identify orthotropic linear elastic properties [26], and we previously used it to investigate fracture properties in a weak transverse crack growth direction [24].

Materials and methods
Crown cut silver birch veneers (Betula Pendula), free from knots, were bought from Holm Trävaror AB (Sweden), with an oven-dried weight of around 590 kg/m 3 and a moisture content of around 8 %, both were measured after one week of conditioning in an oven at 105 • C. The thickness of the veneers was around 2.5 mm, measured with a digital dial indicator (Mitutoyo S112SB, Japan), and the samples herein are made from one and the same veneer sheet.A wood veneer, cut far away from the pith of the wood stem, have three material directions due to the microstructure, see Fig. 1a; longitudinal material direction (L) follows the fibers, and tangential (T) and radial (R) directions are perpendicular to the fiber direction.Here, the samples were cut so that the loading axis follows the L-direction, with the T-direction transverse to the load.This was to ensure that the wood microstructural directions coincided with the orthotropic material directions L, T and R. For studying fracture mechanics of native wood, six fracture planes are commonly postulated to describe possible crack paths [27], denoted by two indices of the three material directions L, T, and R; the first index denotes direction normal to the fracture plane and the second index denotes the crack propagation direction, (see Fig. 1a-b).The initial crack orientation of the specimens was cut to ensure that a machined crack was within the LT fracture plane, see Fig. 1b.

Transparent wood biocomposite preparation
Briefly, the biocomposites were prepared by filling the porous native wood substrate with the monomer methyl methacrylate and cured by in-situ heat-treatment.The transparent wood-polymer biocomposites were made similarly, but the wood substrates were first delignified (removal of lignin and chromophores).Details of the biocomposite manufacturing is explained in previous work [25].The transparency of the biocomposites helped to identify defects and unsatisfactory specimens.A thickness of 2.8 mm of a 10 cm x 10 cm sheet gave satisfactory results, without any visually observable defects.Thicker samples would be very challenging to prepare.

Experimental setup
A four-point bending (4 PB) fracture tests were used to measure the crack growth of small wood veneer-based biocomposites, Fig. 2.This specimen was developed by Sorensen and coworkers, and provides stable crack growth [28].We used a similar experimental set up in fracture tests of D-Birch/PMMA for TL fracture properties [25].Three to five samples were tested for each present material.
A deep edge notch was introduced to the specimens with LT orientation, see Fig. 1b.The cut was done with an IR laser (Redsail CM1490, China) at 20 % maximal power, a laser head speed of 20 mm/s, and a cutting diameter of around 0.2 mm.To avoid generating pre-cracks in the weaker TL fracture plane, potentially caused by heating from the laser, the last millimeter of the notch was made with a 0.2 mm thin blade jigsaw at low speed.A deep and straight notch made with the jigsaw was challenging, therefore, the laser cut was useful.Finally, a sharp pre-crack in the LT fracture plane was made with a 100 μm thick razor blade.The crack length of the specimens was around a 0 = 7.5 mm.The notch geometry was inspected by using an optical microscope (Leica GmbH, M205 FA, Germany).
For the specimens in Fig. 2a-b, dimensions and notch length were used following previous work [24].The specimens were cut with a laser to a length of 60 mm, and a height b = 15 mm.To support stability for the thin specimen geometry and to avoid damage from the rollers of the 4 PB fracture test, two steel supports at each end of the specimens were glued and firmly fixed with screws.The remaining specimen length was w = 20 mm, see Fig. 2b.The steel plates had a length of 50 mm, a height of 15 mm, and a thickness of 3 mm.The 4 PB fracture tests were done in an electromechanical material testing machine (ElectroPulse E1000, Instron, U.S.A), with a 2 kN load cell attached (Dynacell 2 kN, Instron, U.S.A).The load points were controlled with a monotonic prescribed displacement, with a rate of 0.015 mm/min.Images of the specimens were continuously recorded during the 4 PB fracture tests and the experiments took about 30 min to reach the peak load.The crack opening displacement, δ COD , is measured at the lower end of the notch, using digital image correlation (DIC), see Fig. 2b.

Digital image correlation
Stereo-vision camera system (isi-sys GmbH, Germany) was calibrated for a field of view (FOV) of 40 mm with two CMOS 8.9 Mpixel cameras (IMX267, Sony, Japan) and a pixel size of 3.45 μm.The average pixel to mm ratio was 71 pixel/mm.Cameras had Rodagon 5.6/135 mm lenses (Qioptic, UK) with 1.1 x magnification polarized filters.The camera stage was mounted with one camera perpendicular to the specimen and the other tilted above the former camera.The focus was adjusted at f/5.6, thereafter, set to f/16 for camera system calibration, following the manufacturer's instructions.The cameras were connected to a data acquisition unit (DAQ, T8D-16, isi-sys GmbH, Germany) and strobe-light setup (Blue-x-focus, isi-sys GmbH).Load and displacement data from the testing machine were synchronized with the DAQ unit.Displacement and strain fields were analyzed using the DIC software MatchID (MatchID stereo, Belgium).
A white and black speckle pattern was applied on the specimen surface with an airbrush.Local DIC (subset-based) algorithm was used with a zero-mean normalized sum-of-squared differences (ZNSSD) matching criterion.To quantify the static noise of the DIC setup, deformation and strain field resolutions were determined by the spatial standard deviation from ten static images of each sample.The subset size was set to 25 pixels, step size 7 pixels, and strain window size was 11 data points, which satisfied strain resolution ≤ 100μ strain.The in-plane and out-of-plane deformation resolutions were approximately 0.5 μm and 2 μm, respectively.

Elastic behavior
The bulk material was assumed to be homogeneous and in-plane orthotropic linear elastic with four independent elastic parameters: E L , E T , ν LT , G LT, and a plane stress condition.The orthotropic linear elastic properties were known from previous tensile tests [25] and here used as initial input parameters.The elastic properties were validated by comparing DIC, and finite element method (FEM) results, as explained in the following sections.

Cohesive zone model
For the non-linear softening behavior, a cohesive zone model (CZM) was used, allowing for fracture surface separation (here interpreted as a fracture process zone).The crack propagation path was known from the experimental observations, and postulated crack paths were introduced in the model, see Fig. 2a.The crack propagation path was configured to two fracture planes: LT and TL.Within the LT fracture plane with a CZM (CZM LT ), only mode I opening fracture was assumed since the 4 PB load case was symmetric, as shown in Fig. 2a-b.For the TL fracture plane with a CZM (CZM TL ), two cohesive zone formulations were used.The first CZM TL only allowed normal opening, mode I fracture (CZM TL -mode I).In the second CZM TL (CZM TL -PPR), a mixed-mode formulation, developed by Park et al. (PPR) [29], was used, assumed that the fracture process zone exhibit both normal opening, mode I, and in-plane fracture surface sliding by shear, mode II.
To incorporate mixed-mode fracture in the cohesive zone, a potential-based CZM developed by Park et al. was available for implementation into Abaqus through a user-defined element subroutine [29].Details of the PPR cohesive model are found in [30].However, a brief overview of the cohesive stress-displacement relation is included here for clarity.The cohesive zone deformation had two degrees of freedom: opening displacement in the normal direction of the fracture surface, δ n , and in-plane tangential (parallel to fracture surface), δ t , displacements.The opening stress (normal to the fracture surface), σ n , and in-plane shear stress, τ, are defined as where α and β are shape parameters for the opening and shear stresses, respectively; the corresponding fracture energies φ n and φ t are . Final fracture occurs when any of the critical cohesive displacements (opening or shear deformation), δ cn and δ ct are reached [30].The softening damage behavior initiates when the critical cohesive stress is reached σ n = σ ccs or τ = τ ccs , and complete separation of the fracture surfaces occurs when σ n (δ cn ,δ t ) = 0, or τ(δ n , δ ct ) = 0.The limiting opening displacements are defined as, ( The advantage of the PPR cohesive model is that the softening relation is easily changed by the two functional shape parameters α and β.The value α = β = 2 results in almost linear softening behavior.Larger values than 2 result in a convex functional shape, similar to more brittle fracture behavior, and lower values than 2 means a concave shape.The initial linear elastic behavior of the cohesive zone to the specimen deformation was assumed negligible; hence, the cohesive zone mainly governs softening damage behavior.

FEM model
The FEM analysis was processed with an implicit solver in batch mode on a desktop with a six-core processor (Intel Xeon E5-1650 v3, 3.5 GHz), 16 GB RAM, and NVIDIA Quadro M2000 (4 GB RAM).A quad mesh was generated in Abaqus with a mesh concentrated along the crack paths with an element side length of about 0.05 mm (4 pixels).Four-noded elements with reduced integration (CPS4R) and hourglass control were used for the orthotropic linear elastic bulk material.The PPR cohesive zone user-defined elements were attached between the CPS4R elements at the experimentally observed and assumed crack paths.Boundary conditions mimicked the 4 PB experimental setup with four rigid roller supports (see Fig. 2a).The specimen steel supports at the specimen's ends, illustrated in Fig. 2b, were modeled as isotropic linear elastic with Young's modulus E = 205 GPa, and Poisson's ratio ν = 0.3.

Inverse identification of material model parameters
The material model parameters were identified by minimizing the difference in strain fields determined from DIC and FEM, using MatchID and Abaqus 2020.To compare the DIC and FEM displacement field data, the FEM nodal displacements were converted to pixel positions by using the MatchID software.This procedure is described in previous work [24], and details that deviate are briefly introduced.The Nelder-Mead (N-M) optimization algorithm [31] was used to minimize the objective function, J, defined here as where ε i is the strain fields with strain components ε L , ε T , and shear γ LT from DIC and FEM from the ith evaluated load point, λ is a vector of model parameters, n is the number of parameters to be determined, and N is the total number of load points evaluated.A penalty was added to the objective function if trial points violated non-allowable values, e.g., infeasible negative material model parameters.The material model parameters were identified in two sequences: first, the elastic properties were determined, then the CZM parameters.Three load points (N = 3) within the elastic regime were used for elastic property identification and four load points (N = 4) within the damage development and crack growth regime for the CZM parameters.The minimization was terminated when the standard deviation of J i was equal to 10 − 4 .

Table 1
Orthotropic linear elastic properties with previous experimental data [25] and present data from material identification using FEMU.

Results and discussion
This section is divided in three subsections.First section deals with identification of orthotropic linear elastic properties using FEMU and strain field analysis.In the second, the mechanisms of crack growth are analyzed.The fracture process zone is predicted from previously estimated CZM parameters for mode I fracture in the LT and TL fracture planes.Strain fields from DIC (measured) and FEM (predicted) are compared, and the CZM parameters are adjusted by minimizing the DIC-FEM strain field discrepancy.In the last section, the conditions for crack growth by 90 • deflection or penetration are analyzed.

Identification of elastic properties from FEMU
Orthotropic linear elastic properties are presented in Table 1, with experimental results from previous tensile tests and from the present FEMU analysis.The precracked specimens (Fig. 2a-b) are subjected to four-point bending (4 PB) in the elastic regime, and the strain fields are measured with DIC.Initial values for the FEMU procedure was based on previous tensile data for orthotropic linear elastic properties [25].After about 30 iterations, the minimization of the objective function in Eq. ( 3) converged, and the elastic constants were lower than the starting values, except for the Poisson's ratio ν LT for D-Birch/PMMA.The anisotropy ratio, E L /E T , is lower for Birch/PMMA and D-Birch/PMMA biocomposites, compared with neat Birch due to increased E T from the PMMA phase.Note that the transverse T-direction corresponds to the tangential direction in wood.Comparing the starting values with the FEMU results, the Birch/PMMA biocomposites have relatively little change in E L , a decrease of ≈3 %, compared to E L of native birch, D-Birch/PMMA, and pure PMMA (isotropic linear elastic), which decreased by more than 30 %.The time scale of the present experiments is longer than the tensile tests in [25], and creep deformation is likely, as previously observed for beech wood [32].The PMMA matrix reduces creep deformation in native wood but delignification (D-Birch/PMMA) increases the biocomposites compliance.The elastic parameters should be considered as apparent elastic properties fitted for the present load conditions.

Strain field analysis and FEMU data
In Fig. 3, principal strain fields, ε 1 , from DIC are presented at low δ COD ≈ 0.1 mm, measured from the lower end of the pre-crack tip, see Fig. 2. The colored-coded contour plots show strain fields from DIC and discrepancy-plots between DIC and FEM (based on FEMU identified elastic parameters from Table 1).The ε 1 is the maximum strain in a data point, with the principal directions indicated with white arrows for improved visualization.Strain fields for the PMMA specimen is included in the Appendix for comparison, see Fig. A1.
Native birch has an "elliptical" appearance in the DIC strain field around the notch region with considerable tensile (positive) ε 1 values in front and beside the crack tip, and a narrow neck near the crack tip.These strain fields result from the large tensile strains in the T-Fig.3. First principal strain field (ε 1 ) plots from DIC and difference between DIC and FEM for a single specimen of native Birch, Birch/PMMA, and D-Birch/PMMA, at the same δ COD ≈ 0.1 mm.δ COD is defined in Fig. 2. The principal strain ε 1 direction within the field is indicated with white arrows to visualize the strain field gradient.
E. Jungstedt et al. direction, of critical importance in the present study, which are related to the orthotropic nature of the material (see T-direction strain in Fig. A2 in Appendix) since T-direction stiffness is low.The peculiar experimental shape of the large strain region for native Birch is also captured by the FE-model.In Fig. 4, strain fields from each strain component are shown; an important contribution to the strainfield shape at the neat birch crack tip is from large strains in T-direction influenced by the precrack and the orthotropy.Similar shapes of the DIC strain fields are observed for Birch/PMMA and D-Birch/PMMA, although the strain levels are lower compared to native Birch.The reason is that data are from the same opening displacement of the precrack in the beam specimen, δ COD , see Fig. 2b.The presence of PMMA in native wood significantly reduces local strain field heterogeneity and anisotropy.From Fig. 3, we note that the strain difference between measured (DIC) and predicted from FEM (based on FEMU predictions) is the largest for native birch, whereas agreement is significantly better for both Birch/PMMA and D-Birch-PMMA.Since the continuum-analysis from FEM can predict measured DIC strain field shapes, the shapes are not an effect from microstructural heterogeneity but caused by the orthotropic linear elastic properties.
The discrepancy between strain fields from DIC and FEM can be further analyzed to learn about model limitations, e.g., due to microstructural heterogeneities.The standard deviation (SD) of a strain field describes the strain variation around the average strain field value.SD is used here to indicate microstructural heterogeneities by comparing the difference in strain field SD between DIC and FEM (denoted by SD DIFF ); a larger SD DIFF would indicate larger heterogeneity.Only the bulk material is considered, to avoid geometrical complexity effects in the near precrack region.For native birch, the data in Fig. 3 results in SD DIFF ≈ 2•10 − 4 , which is more than an order of magnitude larger than for Birch/PMMA with SD DIFF ≈ 8•10 − 6 .The upper compression side of the native birch sample has the largest discrepancy, see DIC-FEM strain fields in Fig. 3.For D-Birch/PMMA, SD DIFF ≈ 3•10 − 5 and for PMMA SD DIFF ≈ 5•10 − 5 , which are slightly larger values compared to Birch/PMMA.Strain inhomogeneities are much larger for native birch compared to the biocomposites and PMMA, thus indicating the limitation of a continuum-based linear elastic material model for small native birch specimens.

Mechanisms of crack propagation and 90 • crack deflection
In Fig. 5a, early 90 • crack deflection along the grain occurs for D-Birch/PMMA before peak load is reached.Similar observations were made for native birch and Birch/PMMA.Previous investigations on wood and wood biocomposites report the same phenomenon, although peak load is still used for LT (or LR) fracture toughness data [5,13].For accurate mode I LT crack growth fracture toughness, the crack must grow in the LT fracture plane.This can be promoted by introducing side-grooves to force crack growth in same direction as the LT-oriented pre-crack, as we did previously [12].Here, all specimens without side-groove showed crack growth initiation before Fig. 4. Strain field data for all in-plane strain components ε L , ε T , and shear, γ LT , from the linear elastic FE-model of Birch, Birch/PMMA, and D-Birch/PMMA, using the FEMU material properties from Table 1.Strains are estimated at the same opening displacement of the precrack in the beam specimen, δ COD ≈ 0.1 mm.The crack tip location is indicated by the black star in the middle.peak load by immediate 90 • deflection into the TL fracture plane, see Fig. 5a.A qualitative explanation is transverse failure of wood tissue [6] due to the stress-state in the crack tip region.Fractography suggests subcritical damage of the wood fiber cell wall [9].

Cohesive zone model analysis and fracture properties from FEMU
CZMs are estimated here to predict fracture process zone development and analyze the load-δ COD of a D-Birch/PMMA specimen, see Fig. 5b-d (note that δ COD is measured at global specimen level, see Fig. 1, in contrast to displacements in the CZM).Previously determined mode I CZM parameters for crack growth in TL [24] and LT are used initially [12] (see D-Birch/PMMA CZM TL and CZM LT parameters in Fig. 5b).For the CZM TL in Fig. 5c, the shear properties were arbitrarily set to be identical as for opening mode, e.g., φ TL n =φ TL t and σ TL ccs = τ TL ccs .The poor agreement between the model and experimental data in Fig. 5c is addressed by FEM updating using DIC strain fields.Note that original CZM TL data were obtained on specimens with pure opening mode I crack growth [24], in contrast to the present case.
The present CZM TL parameters are then improved by minimizing the difference in strain fields between DIC and FEM by FEMU.The agreement between DIC data and the FEMU model is then quite good, as is verified in Fig. 6.Note that data in show the presence of subcritical 90 • cracks for all three materials already at small δ COD for the specimen.This reduces stress concentrations at the LT  [24] and LT [12] fracture tests, presented in the image, and the cohesive shear properties are initially and arbitrarily chosen to be equal to the cohesive tensile properties in TL. c) Load-plot from the model c) compared to experimental results of a D-Birch/PMMA sample.Load-plots with calibrated CZM parameters, using FEMU, from Tables 2-3 with the model shown in Fig. 5b) for single specimens of Birch, Birch/PMMA, and D-Birch/PMMA.Red dashed line indicates when unstable crack growth is initiated.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)precrack tip.
The new CZM TL parameters are presented in Table 2. CZM LT parameters are found Table 3.We can see from Fig. 5d that the combination of CZM LT and CZM TL leads to good agreement between FEM and experimental load-δ COD data for all materials.The values of σ TL ccs and φ TL n identified by FEMU are larger compared to previous work on TL precracks [24].The reason is that present data are from an LT precrack showing mixed-mode TL crack growth.
Comparing the CZM TL properties for the three materials, Birch/PMMA shows the largest values for φ TL t and for limiting opening displacements (δ TL cn and δ TL ct ), see Table 2. High TL fracture energy leads to short 90 • cracks at a given load, and experimental load-δ COD curves in Fig. 5d show the highest load values for Birch/PMMA.The lower CZM TL properties for D-Birch/PMMA (from FEMU) is the major reason for lower performance compared to native Birch/PMMA.Delignification removes the lignin-rich middle lamella bondline between fibers, which reduces CZM TL properties.Interestingly, for Birch/PMMA, φ TL t and δ TL ct are significantly increased compared to φ TL t and δ TL ct for native birch.In general, the PMMA matrix improves CZM properties.It is interesting to consider the conditions for LT crack growth in the direction of the original crack.In FEM, we can extract nodal data at the crack tip from the CZMs and analyze each degree of freedom in the cohesive zone.In Fig. 7, the stress state at the crack tip of the FE-model is illustrated, using CZM properties from Table 2 and Table 3.The stress components σ L , σ T , and τ LT are normalized to the corresponding cohesive strengths σ TL ccs , τ TL ccs , and σ LT ccs respectively.When the local stress in σ T at the crack tip in the TL cohesive zone reaches its critical value σ TL ccs , the other ratios are lower, so TL crack growth is favored.The exception is native birch where σ L /σ LT ccs ≈ 1 when σ T /σ TL ccs = 1.For all materials, τ LT is always lower than τ TL ccs , meaning that 90 • crack deflection is initiated by crack opening from σ T , when σ T /σ TL ccs = 1.Fig. 7 shows the specimen displacement δ COD at which stress-induced damage is initiated in the 90 • CZM TL for the three materials.
In all cases, the mechanism is indeed 90 • crack deflection in the TL fracture plane (highest values are for σ T /σ TL ccs ).A conservative estimate for damage initiation in the composites (peak value for σ T /σ TL ccs ) is quite low, δ COD ≈ 0.05 mm, whereas δ COD is larger for  2 and Table 3. maximum σ T /σ TL ccs ratio in neat birch (lower transverse modulus).In Fig. 6, experimental DIC data show that substantial TL cracks did not form until much later, at specimen δ COD ≈ 0.3 mm.The corresponding FEMU predictions confirm the 90 • crack deflection mechanism, see lower part of Fig. 6.Exact experimental size of the crack lengths and strain fields (DIC, Fig. 6) differs somewhat from FEMU predictions in Fig. 6, because of microstructural heterogeneities.The displacement δ COD ≈ 0.3 mm is marked in load-displacement curves of Fig. 5d and this figure illustrates that 90 • cracks form well ahead of peak load for all materials.
By further analysis of DIC and FEM results, TL crack surface separation can be measured by the deformation of a "virtual" extensometer across the crack path of the 90 • crack, see Fig. 8.Although 90 • TL cracks are formed early, they are growing slowly.For Birch/PMMA and D-Birch/PMMA at δ COD ≈ 0.145 mm, there is a rapid increase in deformation of the virtual extensometer, due to crack growth, and it levels out at δ COD ≈ 0.15 mm.For native birch, see Fig. 8, the increased deformation of the "virtual" extensometer from subcritical crack growth is detected by DIC already before δ COD = 0.05 mm.The slope in Fig. 8 increases even more for native birch at δ COD ≈ 0.25 mm, when further crack growth is known to take place from changes in DIC strain field shape (see Fig. 6).From FEM of D-Birch/PMMA, the slope of virtual extensometer deformation is lower than for DIC data even after damage is initiated in the CZM (δ COD ≈ 0.05 mm, see Fig. 7).It is possible that there is a dynamic crack growth effect not captured by the FE-analysis.* Rounded values are from experimental data [12], and an assumed linear cohesive softening behavior with α = 2.

Effects of cohesive fracture energy and critical stress on crack penetration and 90 • crack deflection
In general, LT crack growth by 90 • deflection or penetration depends on material properties, both the fracture energies and cohesive strengths [12], as well as local stress state.In Fig. 9, the alternative mechanisms of 90 • crack deflection (TL) or crack penetration (LT) are analyzed.This is based on FEM data of the stress state and parametric studies of variations in the CZM TL and CZM LT parameters.The dimensionless ratios σ LT ccs /σ TL ccs and φ LT n /φ TL n are used to illustrate conditions for crack penetration or 90 • deflection.The cohesive zones in the FE-model can be considered as interphases between elastic regions, allowing for crack growth along the TL and/or LT fracture planes with corresponding CZM TL and CZM LT parameters.For D-Birch/PMMA analyzed in Fig. 9, it is apparent that the high σ LT ccs /σ TL ccs ratio for D-Birch/PMMA results in 90 • deflection for the present material and load case.Crack penetration will only occur provided σ LT ccs /σ TL ccs ≤ 5.0, where φ LT n /φ TL n = 80 for the present composite.At even larger values of φ LT n /φ TL n , a lower limit is reached for the strength ratio σ LT ccs /σ TL ccs ≈ 5.0 and further increase of φ LT n /φ TL n will not affect the mechanism of crack growth.This is in agreement with a general study on the competition of crack deflection and penetration at interfaces [33], although the elastic orthotropy and the load case are different.Note that the condition for crack deflection and penetration is influenced by a characteristic fracture length scale, L = Eφ/(σ 2 ) [34].Here, a set value on L was used based on the results obtained by FEMU in Tables 1-3.From a materials design perspective, fracture mechanisms such as crack deflection at interphases are crucial for composite toughness [35].

Conclusions
Previous literature for wood and polymer-filled wood composite fracture toughness across the grain (LT) is frequently based on a peak load criterion and LEFM analysis.The present cohesive zone description can better capture mechanisms of crack growth.Although the LT precrack is perpendicular to the grain, there are deflecting 90 • TL cracks along the grain formed at very low load.The reason is high opening tensile stress σ T at the crack tip perpendicular to grain and low fracture energy for 90 • TL cracks.These cracks, however, grow very slowly so that maximum specimen load becomes much higher.The TL cracks lower the strain concentration at the tip of the LT precrack.
The fracture energy values φ LT n estimated from peak load for toughness across the grain are misleadingly low, since the growth of deflecting 90 • cracks is controlling peak load.For transparent polymer-filled composites from delignified wood, the chemical pretreatment is reducing fracture energy for 90 • cracks so that peak load is low as for neat birch wood.The peak load is about twice as high for wood/polymer composites without chemical treatment, due to much higher fracture energy for 90 • cracks and contributions from the polymer phase.
The results show that "weak interfaces" is not enough as a guideline for tough composites.Wood-inspired composites, however, with delicately balanced fracture energies (φ TL n , φ TL t , and φ LT n ) and critical cohesive stresses (σ TL n , τ TL , and σ LT n ) could be tailored for structures combining strength with high fracture energy.Based on these parameters, it would be possible to design composites for failure with controlled extent of cracking across and along the grain.strain field ε y gradients can be found for the PMMA specimen (which can be considered isotropic), but with significantly lower strain values ε y compared to native birch.

Fig. 2 .
Fig. 2. A) Schematic of the 4 PB fracture test with hypothetical crack propagation paths along LT (crack propagation) and TL fracture plane (crack deflection).Material axes and fiber direction are shown in left part of the specimen.b) 4 PB specimen geometry and cross-section, with steel plates at the ends.The crack opening displacement, δ COD is measured at the lower end of the precrack using DIC.

Fig. 5 .
Fig. 5. A) Speckled d-birch/pmma sample showing immediate 90 • crack deflection along the T-direction.b) FE-model with two CZMs (CZM LT and CZM TL ) based on a mixed-mode CZM (PPR).Properties are based on previous TL[24] and LT[12] fracture tests, presented in the image, and the cohesive shear properties are initially and arbitrarily chosen to be equal to the cohesive tensile properties in TL. c) Load-plot from the model c) compared to experimental results of a D-Birch/PMMA sample.Load-plots with calibrated CZM parameters, using FEMU, from Tables 2-3 with the model shown in Fig.5b) for single specimens of Birch, Birch/PMMA, and D-Birch/PMMA.Red dashed line indicates when unstable crack growth is initiated.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 6 .
Fig. 6.Principal strain field ε 1 from DIC and adjusted FEMU model at δ COD ≈ 0.3 mm, during stable 90 • crack deflection.The white arrows are used to indicate the principal ε 1 direction.The CZM properties in FEMU results are from Table2 and Table 3.

Fig. 7 .
Fig. 7. Stress ratios at the crack tip of the cohesive zone (local stress divided by critical stress) as a function of δ COD in the specimen.A value of unity for a specific ratio means that the fracture criterion is fulfilled.Data show development of stress ratios before, during, and after crack growth in Birch, Birch/PMMA, and D-Birch/PMMA specimens.

Fig. 8 .
Fig. 8. Experimental DIC data and predicted FEM crack separation (deformation of virtual strain gauge, see text) in the TL fracture plane, generated by growth of the 90 • crack.a) Illustration of how crack separation (deformation of virtual strain gauge) was measured in DIC, similarly done in FEM.b) Crack separation from DIC for birch, Birch/PMMA, and D-Birch/PMMA, also from FEM for D-Birch/PMMA, measured as illustrated in a).FEmodel parameters are from Tables 1 to 3.

E
.Jungstedt et al.

Fig. A2 .
Fig. A2.Strain field data ε y from DIC in the y-direction (vertical direction) of native birch, Birch/PMMA, D-Birch/PMMA, and PMMA specimen in the elastic regime at δ COD ≈ 0.1 mm.Large tensile strains in the vertical direction of the specimen are affected by the material anisotropy.Similar

Table 2
CZM TL PPR (mixed-mode) cohesive model properties determined with FEMU for opening properties (φ TL

Table 3
[12]LT parameters: (φ LT n , σ LT ccs , and δ LT n ) from our previous study on LT fracture toughness of native Birch, Birch/PMMA, and D-Birch/PMMA[12], determined by minimized load-displacement differences between experiment and FEM.δ LT cn is incorporated into the PPR model and determined by Equation (2).