J integral strain-curvature approach for multilayer fracture specimens with unknown thickness and stiffness of layers

In the coming years, many wind turbines will reach their planned design life, and it is of great interest to investigate if their service life can be extended so that the turbines can produce a larger amount of electricity before decommissioning. Such an assessment can include the measurement of the fracture mechanics properties of materials interfaces, e.g., for trailing edge bondlines, web foot bondlines and interfaces between layers in the load-carrying main spars. It would then be preferable to conduct the fracture mechanics testing and data analysis by approaches that give accurate results with the smallest amount of testing. In the present work we propose an approach for measurement of mixed mode fracture resistance (and mixed cohesive laws) that does not require knowledge of the stiffness properties and the elastic centre of the laminates of the test specimens. This is advantageous since the layer stiffness, layer thickness and the lay-up of a laminate cut from an old blade may not be known. Furthermore, the elastic properties of polymer resins of bondlines and matrix materials in fibre composites can have changed (aging) due to the blades long-time environmental exposure.


Introduction
When the use-time of a wind turbine in a wind farm is close to the projected service life (typically, 20 to 25 years), it is relevant to consider if the service life of the turbines can be extended in order to get more energy produced before the turbine is taken out of service.The process of assessing the extension of the use-life of old wind turbine rotor blades in a wind turbine park can involve the characterization of the mechanical properties of the materials and interfaces by testing specimens cut from a single representative blade.The fracture mechanics properties of the blade can then be used to assess the growth of damages in the blades of the remaining turbines.The primary load-carrying parts of a wind turbine blade is made of layers of composites materials (typically, made of glass and/or carbon fibres in a polymer matrix) and therefore it would be relevant to make fracture mechanics testing of interfaces between layers in laminates (e.g., within the main spar) and laminate/adhesive interfaces (e.g., web-foot and trailing edge bondlines).It would be convenient to cut fracture mechanics specimens of regular shape from a blade, so that established fracture mechanics testing procedures can be used.There are well-established fracture mechanics testing methods for linear elastic fracture mechanics characterisation, both for pure Mode I and mixed mode [1,2] (valid for small scale fracture process zone) and J integral test specimens (valid for both small-scale fracture process zones and large scale bridging problems) [3,4].
However, the elastic properties, thickness and orientation of each layer may not be known for fracture mechanics test specimens cut from an old wind turbine blades.Therefore, it would be convenient to use fracture mechanics testing methods for which the calculation of the fracture mechanics properties (mixed mode fracture resistance or mixed mode cohesive laws) can be done without the need for knowing the elastic properties and thickness of the individual layers.For Mode I Double Cantilever Beam (DCB) specimens loaded with wedge forces or bending moments, there are known J integral equations that are independent of the elastic properties [5,6,7].These approaches have been generalized 1293 (2023) 012024 IOP Publishing doi:10.1088/1757-899X/1293/1/012024 2 to DCB specimen mixed mode specimens loaded by transverse forces [8] and DCB specimen loaded with axial forces and bending moments [9].For beam-like specimens loaded by axial forces and bending moments the J integral solution can be expressed in terms of the applied axial forces, bending moments as well as the curvature and the axial strain of the elastic centres of the specimen beams [9].Neither the elastic properties of the individual layers nor the elastic properties of the beams are required.The approach is fairly general; it is applicable for multilayered and graded specimens and, as a J integral solution, remains valid for both small scale fracture (linear elastic fracture mechanics conditions) and for large-scale bridging problems.The approach requires, however, knowledge of the elastic centres for each of the three beams of the DCB-specimen.Although this can be determined experimentally [9], it would be preferable if these additional measurements were not required.Only in the absence of axial forces (i.e., pure bending), the approach does not require knowledge of the elastic centre.
The purpose of the present paper is to establish a J integral equation where forces, moments, curvatures and axial strains are defined from the geometric mid-plane instead of the elastic centre.It is found that while the bending moments and axial strains are different from those defined from the elastic centres, the J integral equation takes the exact same mathematical form.This enables a novel fracture mechanics testing approach, which only requires the overall dimensions (thickness and width) of each beam as well as continuous measurements of the applied axial forces, moments, curvatures and axial strains of the geometric mid-plane of each beam.The strains and curvatures can be measured experimentally by the use of the digital image correlation technique, e.g. by measurement of axial strains at the top and bottom of the three beams, i.e. six strain measurements.For the measurement of the fracture resistance to be accurate, it is necessary that these strains are measured with sufficient accuracy.It is therefore also important to consider how these strains can be measured.Strain gauges are known to be very accurate but could be difficult to place on beams running parallel to each other.Therefore, the use of digital image correlation (DIC) [10] should be considered.However, since strain measurements by DIC are not as accurate as strain measurements by strain gauges, this paper considers how to check, if the accuracy is sufficient to enable J integral measurement from curvature measurements obtained from DIC.
The paper is organised as follows.First, in Section 2 we define the problem.In Section 3, we prove mathematically that the J equation of Toftegaard and Sørensen [9] holds true irrespective of what coordinate system is used to define forces and moments, and set up equations for estimation of the uncertainty of the J integral value accounting for uncertainties associated with DIC strain measurements.Section 4 presents results in the form of tables for the estimated uncertainties.Section 5 discusses the advantages and disadvantages of the proposed approach, and Section 6 contains the main conclusion of the study.

Problem description
The problem we consider is a beam-like fracture specimen consisting of three beams, as shown in figure 1.The specimen has a constant width,  and the beam height is denoted  1 ,  2 and  3 , for Beam 1, Beam 2 and Beam 3, respectively.The specimen can consist of a single, homogenous, linear-elastic material or of multiple layers made of orthotropic materials oriented parallel to the plane of the fracture process zone (the  1 - 3 -plane) or the elastic properties can be graded in the  2 -direction, but must be independent of the  1 -position.For a layered structure, the elastic properties and the layer thicknesses must be independent of the  1 -coordinate.The specimen is loaded with a combination of axial forces and bending moments.The axial forces, acting at the elastic centres, are denoted by  1 ,  2 and  3 (subscript indicates beam number), and the bending moments are denoted  1 ,  2 and  3 .The forces and moments are taken to be positive as indicated in figure 1.The distance from the ends (where the forces and moments apply) to the fracture process zone is several times the beam thickness.The fracture process zone may not necessary be small, it can include large-scale fibre bridging.Therefore, we will analyse the specimen by the path-independent J integral [11], not by linear elastic fracture mechanics.
As mentioned above, this problem has been analysed by Toftegaard and Sørensen [9].The J integral evaluation around the external boundaries can be written as where  1 0 ,  2 0 and  3 0 are the axial strain component,  11 , at the elastic centres of Beam 1, Beam 2 and Beam 3, and  1 ,  2 and  3 are the (elastic centre axis) curvatures of Beam 1, Beam 2 and Beam 3. Thus, by measuring the strain at the elastic centres and the curvatures along with the applied axial forces and moments, the J integral value can be determined without needing to know anything about the elastic properties of the beams, except that they are linearly elastic.
But this approach requires the knowledge of the position of the elastic centre for each of the three beams.In the following, we will express (1) in terms of axial forces and moments associated with the geometric mid-height of each beam instead of the elastic centre of each beam.

J integral result
In the following, we will analyse the J integral contribution from Beam 2 with reference to two coordinate systems, one located at the elastic centre and another coordinate system located elsewhere (e.g., at the geometric mid-height), as shown in figure 2. The J integral is evaluated away from the fracture process zone, where the stress and strain fields are uniform, i.e., independent of the  1 -position [9].Define the local coordinate system with origin in the elastic centre by (̅ 2 ,  � 2 , ̅ 2 ), and the other coordinate system by ( 2 * ,  2 * ,  2 * ) (for both coordinate systems, subscript 2 indicates Beam 2).In the global   -coordinate system, the  2 -coordinate of the elastic centre is denoted by  2 and the  2coordinate of the origin of the other coordinate system is denoted  2 (again, for both symbols, subscript 2 is used to denote Beam 2).The axial force and moment defined with reference to the ( With respect to the original  � 2 -coordinate, the normal strain in Beam 2 can be written as where  2 0 =  11 ( 2 =  2 ).In the ( 2 * ,  2 * ,  2 * )-coordinate system the strain distribution would be written as where  2 0 * =  11 ( 2 =  2 ).Note that the curvature  2 is the same for both descriptions.Since We proceed by expressing  2 0 in terms of  2 0 * .By setting  2 * = 0 (and thus by (5),  � 2 =  2 −  2 ), we get from (3) and (4) using (5): Next, from force equilibrium we find: Finally, we insert  2 from (7),  2 0 from (6) and  2 from ( 7) into (2): and tidying up we notice that the terms with  2 and  2 cancel out, so that eq. ( 8) reduces to The J integral equation expressed in the ( 2 * ,  2 * ,  2 * )-coordinate system, eq. ( 9), takes exactly the same mathematical form as the J integral expression with the normal force and moment defined from the elastic centre, eq. ( 2), with  2 replaced by  2 * ,  2 0 replaced by  2 0 * and  2 replaced by  2 * .It follows, that we do not need to define the axial force and moment from the elastic centre of the beam.Irrespective of what coordinate system is used when defining the axial force and moment, the same J integral equation applies.
Similar results can be found for Beam 1 and Beam 3 (since the integration path of Beam 3 is opposite of the integration path of Beam 1 and Beam 2, the J integral contributions become negative), so that the J integral evaluated around the external boundaries of the specimen can be written as 3.2.Measurements of  2 0 * and  2 from surface strains In practice, the variables  1 0 * , κ 1 ,  2 0 * , κ 2 as well as  3 0 * and κ 3 must be measured during fracture experiments.Focussing again on Beam 2, we will now describe an approach from which we can determine  2 0 * and κ 2 from measurements of the axial strain at the top and bottom of the beam.Denote the axial strain at the top of the beam by  11 2 and the axial strain at the bottom of the beam by  11 2 (see figure 2).From (4) we get the strains as: and Solving the two equations ( 11) and ( 12) for  2 0 * and κ 2 gives and Note that ( 14) is independent of the choice of  2 .
For the special case that the coordinate system is chosen at the geometric mid-plane of the beam,  2 = − 2 2 ⁄ , eq. ( 13) reduces to the expected result Similar equations apply for Beam 1 and Beam 3, i.e., for  1 0 * , κ 1 ,  3 0 * and κ 3 .To sum up, the width  should be measured before the experiment, and the following variables should be measured during the experiment:

Equations for sensitivity analysis
Before using the proposed approach experimentally, it is useful to contemplate the accuracy of the approach of measuring axial strains and curvatures instead of measuring the extension and bending stiffness properties of the three beams, Beam 1, Beam 2 and Beam 3. In the following we will present some simple estimates of the uncertainties on the J integral contributions.Our focus will be on the sensitivity of the measurement of curvature from strain measurements in accordance with (14).For simplicity we will consider a symmetric (i.e., identical beam heights,  1 =  2 ), homogenous DCBspecimen loaded with symmetric moments,  1 = − 2 (with  2 > 0), with  1 =  2 =  3 = 0 and  3 = 0.In an experiment, the relevant moments and strains ( 2 ,  11 2 and  11 2 ) would be measured.In this sensitivity analysis we instead assume the material to have a certain stiffness (Young's modulus) IOP Publishing doi:10.1088/1757-899X/1293/1/0120246 and set up some equations to establish consistent sets of moments and strains as a function of   and  2 .
For a symmetric, homogenous DCB-specimen loaded with symmetric moments (Mode I), the J integral result (10) (or equivalently (1)) reduces to Assuming that Beam 1 and Beam 2 are homogeneous beams, it follows that With the moment of inertia of the rectangular cross section of Beam 1 and Beam 2 given by [12] and the Young's modulus in the  1 -direction denoted by  11 , the moment is given by [12] we find by inserting  2 (expressed in terms of  2 ) from ( 18) into ( 16) From ( 19), we can find the moment required for a given J integral value: Eq. ( 20) gives  2 as a function of  2 (since  depends on  2 , see (17)).With  2 determined, we can derive equations for  11 2 and  11 2 .With  2 = 0,  2 0 * must be zero for a homogenous specimen.An equation for Beam 2, similar to (14), gives so that an equation for Beam 2 similar to (15) gives or Now we can set up calculations to assess the uncertainly calculations.Since, in some test configurations, the bending moment is generated by a force couple [3], we will consider the uncertainty of the moment in terms of both the applied force, , and the moment arm,  2 , where the moment is given by Then, by ( 22) and ( 24), the J integral equation ( 16) becomes (dropping the subscript of   ) The uncertainty of , denoted by ∆, can now be estimated using the general formula for error propagation [13], as follows: IOP Publishing doi:10.1088/1757-899X/1293/1/012024 where ∆, ∆ 2 , ∆, ∆ 2 and ∆ 11 2 are the uncertain of measured values of the applied force, , the moment arm,  2 , the specimen width, , the specimen height,  2 and the strain at the top of the beam,  11 2 , respectively.Partial differentiation of (25) gives ( Inserting ( 27) into (26) gives: Taking the square root of (28) gives the relative uncertainty of , Δ  ⁄ .

Input parameters for sensitivity study
Although our analysis, Eq. ( 28), includes the sensitivity of all parameters, we will in the following focus on the effect of strain accuracy for curvature measurements.This is motivated by the fact that different strain measurement techniques have different accuracy, and it is important to choose a strain measurement techniques that lead to a sufficient accuracy result for the J integral value.The uncertainty of strain measurements by foil strain gauges bonded to the surface is estimated to be ± 5 µstrain by Carlsson et al. [14].Hannah and Reed [15] discuss many of the parameters that can influence the output from strain gauge measurement but wrote that, with reasonably care, one can realize strain readings accurate to ±2 % with a threshold strain sensitivity of 5 µstrain.Similar values can be found in other text books on strain measurements [16].
The uncertainty of strain measurements by DIC has been estimated to be about 10 -4 (100 µstrain) [17].In cases where the strain is uniform, the strain values for 100 images can be averaged and accuracy of the averaged value can be about 3 x 10 -6 (3 µstrain) [18].However, Reu et al. [19] found, comparing results from 11 different DIC codes, that for a varying strain field (Sample 14 in that study) that by averaging 50 rows of data, the standard deviation in the strain was between 172 x 10 -6 to 875 x 10 -6 (172 -875 µstrain), with the majority in the range of 300 -600 x 10 -6 (300 -600 µstrain).
In the following, we use these values of typical DCB test specimens for glass fibre composites:  = 30.0mm,  2 in the range of 2 to 16 mm and  11 = 40 GPa [4].We use a moment arm with  2 = 86.5 mm.We investigate the situation of characterising a rather weak interface (fracture energy of 10 J/m 2 ), as well as a value corresponding to onset of crack growth of typical glass fibre composites (200 J/m 2 ) and a relative high fully-developed (steady-state) fracture resistance value (1000 J/m 2 ) [4].
Table 1 lists the squared uncertainty terms that are kept fixed in the present study.They constitute some of the parameter to go into (28).The remaining terms are explored next.

Table 1. Fixed uncertainty terms.
Specifying the fracture resistance, and thus the applied J integral value,   , the stiffness,  11 , as well as the width and height of the beams,  and  2 , eq. ( 20) then gives the required moment,  2 .Next,  2 is calculated from (18), and from (23) we get  11 2 .Inserting these values together with the estimated uncertainty terms into (28) gives the relative uncertainty of .

Estimated uncertainties
Results for uncertainty terms are shown in Table 2 for  2 = 2 mm (rather thin DCB specimens for wind turbine blade materials [4]), Table 3 for  2 = 8 mm (a representative DCB test height [4]) and Table 4 for  2 = 16 mm (thicker specimens).The squared terms shown under the light green background colour enter eq. ( 28).Insight can be gained by comparing the magnitude of the squared uncertainty terms (showed under the light green colors in Tables 1 to 4) associated with the various parameters.The contributions from ,  2 and  are in the order of 10 -5 to 10 -4 (Table 1).The squared uncertainty terms associated with  2 and  11 2 are (in nearly all cases) larger.The squared uncertainty term from  2 is in the range of approximately 10 -4 to 0.01, decreasing with increasing  2 .The squared uncertainty contribution from  11 2 is in the order of magnitude of 10 -4 to 1, increasing with increasing  2 .Comparing the tables, it is seen that in most cases, the uncertainty contribution from  11 2 is the largest contributor to the combined uncertainly of .Optimization of specimen geometry could include reducing the term (∆ 11 2  11 2 ⁄ ) 2 to the same magnitude as the other terms.

Limitation of the approach
Tables 2 to 4 show that for  = 10 J/m 2 , a strain resolution of 300 µstrain gives an uncertainty of  that is larger than 20%.This suggests that DIC would not be suitable for measurement of strains for curvature for materials and interfaces having low fracture energy values.Measurement using strain gauges would most likely be needed.
In practice, it would be difficult to bond foil strain gauges to the bottom of Beam 1 and top of Beam 2 (the cracking plane).It is therefore highly desirable to be able to obtain the strain measurements from the frontal face of the beams using DIC.A broader sensitivity study could be made to determine the best specimen dimensions following the lines of the small example presented here.
The approach is not limited to the type of results given in Table 2 to 4. The formulas in Section 3.3 can for example be used to determine the strain resolution necessary to obtain a desired relative uncertainty (e.g., ∆/ = 20% or maybe ∆/ = 10%) for different levels of  2 and .

Practical measurement of 𝜀𝜀 11
2 and  11

2𝑏𝑏
As described earlier, for a beam segment loaded by an axial force and a bending moment, the strain component  11 varies linearly across the beam height in accordance with (3) and (4), being independent of the  1 -position.Therefore, the averaging approach proposed by Wang and Tong [18] can be used for the present specimen remote from the fracture process zone.More precisely, the strain component  11 could be averaged along the length for the same  2 -position (or as an example,  2 * -position for Beam 2).Furthermore, the fact that the strain varies linearly across the beam in accordance with (4), enables precise measurement of  11 2 and  11 2 (or  2 and  2 0 * ) as follows.First, as shown in figure 3a), the strain component  11 is determined for many points along a line at a given  2 -position (the approach of Wang and Tong).These strain values are then averaged (see Fig. 3b).Likewise, the average value of  11 is determined for points along another line at a smaller  2 -position.This can be done for several lines along the  1 -direction, creating a data set of  11 =  11 ( 2 ).A straight line can be fitted to the  11 ( 2 )data and the parameters  2 and  2 0 * can be obtained by linear regression in accordance with (4).Such an approach might give a strain measurement accuracy better than the range of 300 -600 µstrain obtained by Reu et.al. [19].

Conclusions
The J integral contributions for beams of beam-like specimens subjected to axial forces and bending moments can be calculated from measurement of applied loads and resulting deformation (strains and curvature) of three beams.No stiffness data are required.The approach eliminates the need to make additional tests to determine the tensile and bending stiffnesses.
Such a testing approach is particularly useful for fracture mechanic characterisation of materials with unknown elastic properties, and specifically for analysis of substructure cut from old blades where the elastic properties and thickness of layers are not readily known.The proposed approach of using digital image correlation for strain and curvature measurements for fracture mechanical characterisation by  may be associated with large uncertainties for materials and interfaces that have a low fracture energy (here exemplified by a fracture energy of 10 J/m 2 ).But for materials and interfaces with higher fracture resistance (fracture energy in the order of 200 J/m 2 or higher), the proposed approach seems feasible.

Figure 1 .
Figure 1.The problem under investigation: A beam-like specimen loaded with axial forces and bending moments remote from the fracture process zone.The elastic properties of the beams are independent of the  1 coordinate but can vary layer wise or graded in the  2 -direction.

Figure 2 .
Figure 2. The definition of the axial force, moment, strain and curvature with respect to the two coordinate axes,  � 2 , and  2 * .

Figure 3 .
Figure 3.An approach for measurement of strain variation across a beam segment subjected to uniform extension and bending by the DIC technique.(a) measurement of strains of many points along several lines parallel to the beam (the  1 -directions), (b) for each line average points are anticipated to lie on a straight line, which allows extrapolation to identify  11 2 and  11 2 .

Table 2 .
Estimated uncertainty terms for beam height  2 = 2 mm calculated as a function of applied  value and various uncertainties of strain measurements.The strain uncertainty values in the fourth column from the left are absolute values (Δ 11 2 ) except for relative values (Δ 11 2 / 11 2 ) marked with #.