Fracture predictions in impact three-point bending test of European beech

Kubíketal. Journal of Wood Science (2024) 70:42 https://doi.org/10.1186/s10086-024-02157-x ORIGINAL ARTICLE Open Access © The Author(s) 2024. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/. Journal of Wood Science Fracture predictions inimpact three-point bending test ofEuropean beech Petr Kubík1, František Šebek1* , Mojtaba Hassan Vand2, Martin Brabec2 and Jan Tippner2 Abstract Hardwood has become widespread in European forests. The strongest factor is climate change and damage to conifers by the bark beetle. The effort to study hardwoods grows with increasing volume of applications. Therefore, European beech wood was investigated under two impact loads in two material directions, resulting in four unique combinations supplemented by the measurement of the friction coefficient. Then, it was computationally simulated to reproduce the cracking, while the material model reflected the orthotropic behaviour in elasticity, plasticity and failure. The model was coded using the user subroutine in Abaqus to initiate and propagate the crack using the element deletion. The resulting reaction forces were in good agreement with those from the experiments. Cracking was numerically simulated in three of four cases as experimentally observed, however, upon larger deflections. Therefore, the model is applicable for further investigations. Keywords Anisotropy, Barlat, Dynamic, Explicit, Flexure, Rupture Introduction Wood is a widespread natural composite, which has been used as a construction material for centuries because it is renewable, biologically degradable, environmentally friendly, strong, lightweight, easy to manufacture, electrically resistant, absorbing noise and aesthetic. Therefore, models capable of predicting deformation and failure under operational or random loading are needed. Bending is one of the common tests of wood. Jansson [1] investigated failure modes and stresses using static and impact bending tests showing decreasing impact bending strength with decreasing time to failure for the Spruce Pine Fir (SPF). Yoshihara etal. [2] used a static three-point bending test to determine the shear modulus with the help of a correction function and a modified Timoshenko beam for six wood species: Sitka spruce, western hemlock, akamatsu, yellow poplar, shioji and balsa. Yoshihara etal. [3] conducted a similar study for asymmetric four-point bending tests. Yoshihara and Oka [4] carried out the compression bending test to determine the elastic modulus, proportional limit and bending strength on specimens of Japanese fir of various lengthto-thickness ratios. Then they compared the results with the conventional bending test to confirm that the correct flexural properties can be obtained from the compression bending test for a large length-to-thickness ratio. Kubojima etal. [5] used the impact bending test to estimate the elastic modulus for Japanese cedar, hondo spruce, hiba arborvitae, Japanese red pine, paulownia, Manchurian ash and Japanese evergreen oak. Polocoşer etal. [6] investigated the effect of low-velocity impact on failure stresses and stiffness using the bending test consisting of a pendulum and specimens having a length of 650mm, width of 50mm and thicknesses of 20, 30 and 40mm to find that failure was significantly different compared to quasi-static tests for beech, larch and pine. Polocoşer et al. [7] found higher energy absorption in pine and spruce reinforced with E-glass on the surface loaded in *Correspondence: František Šebek [email protected].cz 1 Institute of Solid Mechanics, Mechatronics and Biomechanics, Faculty of Mechanical Engineering, Brno University of Technology, Technická 2896/2, 616 69 Brno, Czech Republic 2 Department of Wood Science and Technology, Faculty of Forestry and Wood Technology, Mendel University in Brno, Zemědělská 810/3, 613 00 Brno, Czech Republic Page 2 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 tension on a pendulum compared to specimens without reinforcement. The factor was 1.4 for pine and 2.5 for spruce. Jacques etal. [8] carried out full-scale testing on individual light frame lumbers made from SPF. Then, the four-point bending tests with strain rates of 6 × 10–6 to 4 × 10–1 s–1 served to propose the stress–strain relationship applicable to the modelling of blast loadings. Olmedo etal. [9] followed the dynamic response of fresh stems under impact loading using the Mouton–Charpy pendulum. Then they analysed the load-bearing capacity of wooden constructions made of felled trees serving as protection against falling rocks. Brancheriau etal. [10] proposed an analytical relationship, which converts the elastic modulus estimated from the three-point bending test to the elastic modulus appearing in the fourpoint bending test. The results were experimentally validated on specimens from six wood species with different densities. Different elastic moduli of spruce and oak in threeand four-point bending tests were also analysed by Babiak etal. [11]. Gaff etal. [12] examined the effect of thermal modification of European oak and Norway spruce using the Charpy pendulum to find that changes in basic chemical components and impact bending strength were less affected for oak than for spruce. The effect of thermal modification on impact bending strength was also investigated by Hassan Vand and Tippner [13] to find that thermal modification caused a decrease in deflection and maximum longitudinal strain up to approximately 50% according to digital image correlation (DIC) for three-point bending tests of five wood species: ash, beech, larch, oak and spruce. Hassan Vand etal. [14] also tested the effect of moisture content on the behaviour of three wood species of beech, oak and spruce. They evaluated the work required to initiate the crack and break the specimens and found, using DIC, that the maximum deflection and longitudinal tensile strain increased with increasing moisture content. The minimum of the literature deals with the computational modelling of wood under impact. Therefore, this work focuses on the dynamic behaviour of European beech (Fagus sylvatica L.) under the three-point bending test. All tests were carried out using a drop-weight impact testing machine with two hammers to achieve various impact velocities. Then, numerical simulations followed to develop a material model with good fracture predictability compared to experiments. It should be noted that many works [15, 16] use the orthotropic yield criterion according to Hill [17] for the description of plastic behaviour. However, some material parameters (three out of six) for the equivalent stress can be negative due to significantly distinct yield stresses in respective wood directions [18–20]. The negative material parameters can then cause the equivalent stress to be undefined. Therefore, the orthotropic yield criterion according to Barlat etal. [21] will be used to overcome this issue. Experiments Material The material studied was European beech belonging to hardwoods. The specimens were made from tree grown near Brno (Czech Republic). The wood was considered orthotropic with a direction parallel with three primary (height) growth called longitudinal (L), a direction perpendicular to the tree primary growth passing through the pith of the tree called radial (R) and a direction perpendicular to tree growth and simultaneously as a tangent to annual rings called tangential (T) as depicted in Fig.1. All specimens were stored in a climatic chamber at a temperature of 20°C and a relative humidity of 65% to produce a uniform equilibrium moisture content of 12%. Three‑point impact bending test All specimens were manufactured with fibres oriented along the length of 300mm having a square cross section with an edge of 20mm. Before testing, specimens were weighed in analytical balance with 1mg readability, which resulted in an arithmetic mean of density of 727kg/m3 with a standard deviation of 33kg/m3. The symmetric three-point bending was carried out in two material directions R and T (Fig. 1), respectively. The span-to-depth ratio was 12, resulting in a support span of 240mm. Both the supports and the hammer had a radius of 15mm and a surface roughness Ra of 0.2μm. The tests were performed on the drop-weight impact testing machine DPFest 400 from Labortech (right Fig.2) at room temperature. Two different hammers were dropped from two initial heights to obtain two impact velocities as given in Table1 to provide a diverse material for numerical simulations. The testing machine automatically measured the true impact velocities (Table1). For each of these two conditions, 6 specimens were tested in the R and T directions, resulting in 4 datasets. The hammer displacement was measured with 0.01 mm precision, while the force was measured using a piezoelectric force transducer CFT + 50 kN from HBM attached to the hammer. The tests were recorded with Fastcam SA-X2 from Photron. The high-speed camera had a cell size of 20μm and was equipped with a Nikon Micro-Nikkor G lens with a focal length of 105mm and a Nikon Z TC-2.0 × . It was placed approximately 0.9m from the lateral specimen surface, which was parallel to the camera sensor. Additional lighting was ensured by two MultiLed QT standalone lamps (left Fig.2). The field of view was fitted to a central part of lateral specimen’s surface with a hammer. Images with a resolution of 1024 × 672px were Page 3 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 captured with a frame rate of 20,000fps. The image postprocessing was performed in Vic-2D 2010 software from Correlated Solutions. In order to determine the conversion factor, a simple calibration was done using the one-dimensional scale defined by the known distance. The displacement field was calculated from the field of 3 × 3pt. The Lagrange strain field was determined with the lowest possible strain filter size of 5 × 5pt so that the maximum possible spatial resolution was achieved. In 20 T LR L R T 20 240 300 R15 R15 R15 Fig. 1 Loading in the R (top) and T directions (bottom) with all material directions during impact three-point bending tests (all dimensions in mm) Impact test i ng mac hi n e H igh -spee d camera L ED ligh ts D ata ac q u i s i t i on s y ste m Fig. 2 High-speed camera with accessories (left) and drop-weight impact testing machine (right) Page 4 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 addition, the optical measurement served for the analysis of the fracture (if any). Measurement offriction coefficient Tribological experiments were performed to collect input for numerical modelling. The tests were carried out using the Universal Mechanical Tester (UMT) TriboLab from Bruker (Fig.3). The pin of 5mm in diameter was made of AISI 52100 alloy steel having a surface roughness Ra of 0.2μm, the same as the supports and the hammer in the three-point impact bending test. It was forced against wooden specimens by 196 N, which was measured by a dual force sensor DFH-100 with a sampling frequency of 100 Hz and corresponded to a contact pressure of 10MPa. This contact pressure was chosen as it reached tens of MPa in the following numerical simulations. The pin oscillated with a frequency of 2Hz and amplitude of 5mm for 120s. Six specimens with dimensions of 40 × 20 × 10mm were tested. These specimens were manufactured from the same tree as those for the three-point impact bending test and stored in the climatic chamber to reach the same uniform equilibrium moisture content of 12%. The pin was forced against three LR and three LT surfaces of 40 × 20mm moving in the R and L directions, respectively. Therefore, the measurement was within one annual ring for the L direction and approximately 5 annual rings for the R direction with an average annual ring width of 2.5mm. The pin trajectory was 10mm, but the results were processed in MATLAB R2024a for a 6-mm-long portion of the specimen in order to omit the dead ends where the pin decelerated and accelerated. Then, the considered friction was kinetic and area for evaluation approximately 50 mm2. Also, the evaluation was not conducted for the beginnings of the tests where a ploughing was present. Numerical simulations The calculations were performed in Abaqus/Explicit commercial code based on the explicit formulation of the finite element method. The material model was implemented using the VUMAT user subroutine, as it is not a standard one. The failure was modelled by deleting elements that reached critical damage, which is a simple technique that does not require remeshing, contrary to the node separation method. Model ofmaterial The hammer and supports were made of steel, which is much stiffer than the tested wood specimens. Therefore, the steel parts were modelled as rigid, saving some computational time. On the contrary, wood was modelled as a homogeneous orthotropic continuum. The heterogeneity was neglected to save some computational time again through material model simplification, including its calibration. The orthotropic elasticity was described by the generalised Hooke’s law as where εi , Ei and σi are the normal strain, elastic modulus and normal stress, respectively, εij , νij , Gij and σij are the tensorial shear strain, Poisson’s ratio, shear modulus and shear stress for i,j=L,R,T , respectively. All elastic constants (Table2) were taken from [20], where the disintegration of the same wood was carried out under various strain rates in various material directions. The orthotropic plasticity was described by the model of Barlat etal. [21], who proposed the yield condition as (1)        εL εR εT εLR εRT εTL        =           1 EL − νLR EL − νLT EL 000 − νRL ER 1 ER − νRT ER 000 − νTL ET − νTR ET 1 ET 000 0001 2GLR 00 0 0 001 2GRT 0 0 0 0001 2G TL                  σL σR σT σLR σRT σTL        (2) σB−σy=0, Table 1 Configuration of the impact three-point bending tests Hammer weight (kg) Initial height (mm) Average true impact velocity (m/s) Number of specimens (–) Material direction (–) 9.05 458.87 2.70 6 R 9.05 458.87 2.70 6 T 4.55 815.77 3.25 6 R 4.55 815.77 3.25 6 T Drivin g mechanism Sp ec i men h o ld e r Sp ecimen Pin F o r ce se n sor L oa di n g mec h an i sm Fig. 3 UMT TriboLab used for measurement of friction coefficient Page 5 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 where σy is the yield stress, which is dependent on the equivalent plastic strain (Fig.4), while σB is the equivalent stress according to Barlat etal. [21] as (3) σ B =1 m √2 m | K1 − K2 | m +| K2 − K3 | m +| K3 − K2 | m , where m is the exponent (influencing the shape of the yield surface), K1,K2 and K3 are the principal values of the linearly transformed deviatoric stress tensor where a , b , c , f , g , and h are the material parameters, which were fitted (Table3; Fig.5) so that the yield surface corresponded to that in [20], where the model of Hill [17] was used with a negative material parameter (H<0) . A change in the yield criterion was sought due to the elimination of a negative material parameter resulting from a significant difference in the yield stresses in the L direction compared to others (R and T in Fig.5), which can be problematic in some cases. The yield criterion was taken from [20], where it was not necessary to model different yield stresses in tension and compression in finite elements. Then, the yield criterion was recalibrated for European beech, still neglecting the strength differential effect. However, more sophisticated yield criteria should be considered when significant differences in tensile and compressive yield stresses were observed. (4) K =    c (σL − σR) −b (σT − σL) 3hσLR gσTL hσLR a(σR−σT)−c(σL−σR) 3fσRT gσTL fσRT b(σT−σL)−a(σR−σT) 3   , Table 2 Elastic constants for European beech wood [20] Material direction L R T LR RT TL Elastic modulus (MPa) 13,000 3500 3000 – – – Poisson’s ratio (–) – – – 0.302 0.362 0.318 Shear modulus (MPa) – – – 1608 460 1059 00.2 0.4 0.6 0.8 1.0 Equivalent plastic strain [–] 0 10 20 30 40 50 Y i e ld stress [MPa] Fig. 4 Flow curve of the European beech wood [20] Table 3 Plasticity-related material parameters for European beech m (–) a (–) b (–) c (–) f (–) g (–) h (–) 2 2.357 0.862 0.009 0.456 0.642 0.525 15 30 37.5 22.5 7.5 –7.5 –22.5 –37.5 –75 –45 –15 15 45 75 –75 –45 –15 15 45 75 R[MPa] 2 3 ( ( 0,5 T[MPa] 2 3 ( ( 0,5 R[MPa] L[MPa] 2 3 ( ( 0,5 L [MPa ] LR [MPa] Fig. 5 Yield locus in the Haigh–Westergaard space (left) and the space of two normal and one shear stresses (right) for the European beech wood Page 6 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 The orthotropic cumulative damage was described by the damage parameter as: where εf i , εD i and εp i are the fracture strain, plastic strain for a given loading path and plastic strain, respectively. As mentioned above, the element is removed when the critical value is reached, specifically when max (Di)=1 . Furthermore, fracture strains are dependent on stress triaxiality (apart from the rate dependence introduced later) to incorporate the tension/compression failure asymmetry (which is simpler than the fracture models in [15, 22–24]) where σm is the mean stress and σ is the equivalent stress according to von Mises The dependence of fracture strains on stress triaxiality was based on [20] and was further refined using the trial and error method using tens of numerical simulations so that the results matched experimental observation in impact three-point bending tests. Fracture strains dependent on the stress triaxiality, εf i (η ) , are shown in Fig.6 for the reference equivalent plastic strain rate of 1 s–1. Finally, fracture strains were additionally dependent on the equivalent plastic strain rate as experiments were conducted at various impact velocities. Contrary to the linear rate dependence in [20] for significantly greater strain rates, an exponential dependence on the equivalent plastic strain rate based on the equation proposed by Johnson and Cook [25] was recalibrated to fit the experiments as follows (using εf i (η ) that has already been calibrated previously for the reference equivalent plastic strain rate of 1 s–1) where C is the material parameter, calibrated as 0.1 after several additional numerical simulations for reference (5) D i= ε D i  0 1 εf i  dεp i , (6) η= σm σ , (7) σ m = σ L +σ R +σ T 3, (8) σ= 1 √2 (σL − σR)2 + (σR − σT)2 + (σT − σL)2 + 6σ2 LR + 6σ2 RT + 6σ2 TL . (9) ε f i=εf i(η)  1+Cln  ˙ ε p ˙ ε0 , equivalent plastic strain rate ˙ ε 0=1s −1 . Finally, ˙ εp is the equivalent plastic strain rate where ˙εp is the plastic strain rate tensor. Model ofgeometry andboundary conditions As mentioned above, the advantage of the element deletion technique is its simple implementation. However, it has its drawbacks as it depends on the size of the element, which is usually kept as small as possible to realistically propagate the crack using one or two elements. Therefore, the elements had a size of 0.1mm in areas of wood failure as well as in contact regions (Fig.7). Other areas were meshed with elements of a size of 1.82mm. Only a 0.1-mm-width specimen was modelled as a plane strain to represent the inner layer of wood instead of using 20mm to save considerable computational time. Geometry was meshed with 8-node linear brick finite elements with reduced integration and hourglass control (labelled C3D8R in Abaqus). The mesh consisted of 128,774 nodes and 63,850 elements of which 4 were additionally 6-node linear wedge finite elements with reduced integration and hourglass control (labelled C3D6R in Abaqus) near the contact regions where the fine mesh changed into coarse mesh (highlighted in red in Fig.7). As mentioned earlier, the hammer and supports were modelled as rigid, therefore, as surfaces which come into contact with the specimen. The width of the rigid (10) ˙ ε p = 2 3˙ εp :˙ εp , Stress triaxiality [–] –0.5 00.5 1.0 0 0.5 1.0 1.5 2.0 Fracture stra i n [–] 0 L R T Fig. 6 Fracture strain for all material directions and the reference equivalent plastic strain rate of 1 s–1 Page 7 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 bodies, 0.3mm, was greater than the width of the specimen, 0.1mm, and were centred with each other. The hammer and supports were meshed with 4-node bilinear quadrilateral rigid finite elements (labelled R3D4 in Abaqus) having a size of 0.1mm, therefore, 3 elements per width. One support, modelled the same as a hammer, had 1892 nodes and 1416 elements. Moreover, an element was added to the reference point of the hammer (highlighted by a blue pentagram in Fig.7) in order to incorporate the volume load (labelled MASS in Abaqus), equivalent to the hammer weight. Displacements and rotations of reference points of the supports were restricted. The same applies to the hammer except for the vertical displacement. The initial velocity was prescribed in this direction, corresponding to the average true impact velocity of the experiment (Table1). The specimen had constrained displacements only in the width direction to simulate the plane strain condition (the geometry of the specimen was discretised by 1 element along the width). The acceleration due to gravity of 9.807 m/s2 was applied on the whole geometry in the direction of impact. In addition to the definition of the contact between the specimen, hammer and supports, self-contact was applied to all element edges in the area of expected failure (where the elements with 0.1mm edges were). It ensured the contact of new free surfaces, which were not initially present in the geometry, emerging after the element deletion that simulated the cracking. The normal behaviour of the contacts was set as ‘hard’ to allow any pressure when the surfaces are in contact, while the friction coefficient of 0.15 obtained from the experiments in the previous section was set in the tangential direction. Results anddiscussion Figures8 and 9 summarise all the force responses against the hammer displacements obtained from the impact testing machine for hammers that weighed 9.05 and 4.55kg, respectively. The DIC did not serve to obtain the displacements. The dotted lines highlight the tests in which the specimens remained intact. The corresponding specimens bended forward and back (the hammer returned approximately to the reference position corresponding to the displacement of 0mm—it is zero deflection). The dashed-dotted lines highlight the tests in which the specimens partially broke (did not break through the whole height). The corresponding specimens exhibited significant springback, therefore, the hammer returned approximately to the reference position. The solid lines highlight the tests in which the specimens failed (broke through the entire height). The corresponding specimens remained bent after the test (all in Fig.8). Moreover, the moment of crack initiation is highlighted by a circle. The average value of the hammer displacement corresponding to the moment of crack initiation is highlighted by a thick vertical dashed grey line. This line is not plotted in left Fig.9, where only 2 out of 6 specimens partially broke and none failed for the hammer weighing 4.55kg loading in the R direction. Only 1 specimen failed for the hammer that weighed 4.55kg in the T direction (right Fig.9). Other specimens partially broke, while the specimen corresponding to Test 6 in right Fig.9 was cracked along almost the whole height under these conditions, therefore, the hammer started returning later. The arithmetic mean friction coefficient of 0.15 with a standard deviation of 0.005 was measured irrespective of the material direction. Then, the experimentally and computationally obtained results are compared in Figs.10 and 11. The predicted forces were in good Fig. 7 Finite element mesh with details Page 8 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 agreement with the experiments before the maximum deflection. Forces oscillated more in numerical simulations than in impact three-point bending tests. However, this is an intrinsic feature of explicit dynamics. Cracking was predicted later in numerical simulation than in experiments for a 9.05kg hammer in Fig.10. No cracking 0246810 12 14 Displacement [mm] 0 1000 2000 3000 4000 Force [N] Test 1 Test 2 Test 3 Test 4 Test 5 Test 6 0246810 12 14 Displacement [mm] 0 1000 2000 3000 4000 Force [N] Test 1 Test 2 Test 3 Test 4 Test 5 Test 6 R, 9.05 kg T, 9.05 kg Fig. 8 Force–displacement responses for a 9.05 kg hammer in the R (left) and T directions (right), respectively (for interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article) 0246810 12 14 Dis p lacement [mm] 0 1000 2000 3000 4000 Force [N] T, 4.55 kg Test 1 Test 2 Test 3 Test 4 Test 5 Test 6 0246810 12 14 Dis p lacement [mm] 0 1000 2000 3000 4000 Force [N] R, 4.55 kg Test 1 Test 2 Test 3 Test 4 Test 5 Test 6 Fig. 9 Force–displacement responses for a 4.55 kg hammer in the R (left) and T directions (right), respectively (for interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article) 0246810 12 14 Dis p lacement [mm] 0 1000 2000 3000 4000 Force [N] R, 9.05 kg Experiments Simulation 0246810 12 14 Dis p lacement [mm] 0 1000 2000 3000 4000 Force [N] T, 9.05 kg Experiments Simulation Fig. 10 Comparison of force–displacement responses from numerical simulations and experiments for a 9.05 kg hammer in the R (left) and T (right) directions, respectively Page 9 of 12 Kubíketal. Journal of Wood Science (2024) 70:42 was achieved for the R direction with a 4.55kg hammer (left Fig. 9) computationally to reproduce the experiments (left Fig.11). However, partial cracking was not achieved computationally for the T direction with a 4.55kg hammer as in the experiments (right Fig.9). The elastic strain energy was so high in numerical simulations for that case that led to the prediction of crack propagation along the whole height, contrary to the impact three-point bending tests where specimens broke just partially (right Fig. 9). The homogeneity of the model also contributed to the lack of prediction of crack arrest, which occurred in experiments after some delamination because, among others, the wood is strongly inhomogeneous. This should be taken into account by modelling the geometry based on data from computed tomography specifically for each specimen, which might not be practical for industrial applications. The predicted contours of the damage parameter are displayed in Figs.12 and 13 after the impact three-point bending tests compared to experiments approximately in scale. However, it should be noted that the predicted deflections are greater (approximately 1.5mm), as is evident from Figs.10 and 11. The damage parameters DR and DT reached their critical value in respective R and T 0246810 12 14 Dis p lacement [mm] 0 1000 2000 3000 4000 Force [N] R, 4.55 kg Experiments Simulation 0246810 12 14 Dis p lacement [mm] 0 1000 2000 3000 4000 Force [N] T, 4.55 kg Experiments Simulation Fig. 11 Comparison of force–displacement responses from numerical simulations and experiments for a 4.55 kg hammer in the R (left) and T (right) directions, respectively R, 9.05 kg T, 9.05 kg 38.029.000.185.076.057.033.024.005.080.071.052.00.00 Fig. 12 Experiments (top) compared to predicted contours of the damage parameter (bottom) for the 9.05 kg hammer and the R (left) and T (right) directions (approximately in scale)