Building blocks towards progressive fatigue damage modeling of a whole wind turbine blade

This paper presents a concise summary of the recent fatigue damage modeling research in the Mechanics of Materials and Structures (MMS) research group at Ghent University, which utilizes physics-based models and microscopic experiments to develop components for a progressive fatigue damage model of laminated composite structures. The paper discusses the evolution of ply cracking and delamination under cyclic loading and highlights distinct modeling elements for stress analysis and fatigue damage evolution characterization. It also explores the experiments required for model calibration and validation, including the use of digital image correlation to quantify damage mechanisms. Additionally, the paper discusses the use of neural networks for efficient and robust modeling of damage effects and outlines the remaining steps required to develop a comprehensive fatigue design tool. Finally, the recent developments with the BladeMesher software of the research group provide a platform to integrate this fatigue design tool for complete wind turbine blades.


Introduction
The use of continuous fiber composite laminates in wind turbine design has grown in popularity due to their superior mechanical properties.These laminates, made of high-strength fibers embedded and oriented in multiple directions in a polymer matrix, offer a unique combination of strength, stiffness, and lightweight characteristics, ideal for wind turbine blades.Continuous fiber laminates are also highly resistant to fatigue and environmental degradation [1], allowing them to withstand the extreme loads and harsh weather conditions that wind turbines face.However, their use presents challenges such as the expensive design iterations, anisotropic material properties, and extensive experimental testing.While modeling can provide insights into the material's behavior and reduce the cost and time associated with testing, it should be validated with experimental data.Computer modeling can simulate the behavior of composite materials under different loading conditions and provide insight into their mechanical properties without the need for extensive experimental testing [2].Additionally, modeling can provide a better understanding of the anisotropic behavior of composite materials [3], allowing designers to optimize the layup to achieve the desired mechanical properties.However, while modeling has several advantages, it should be noted that it is only as accurate as the underlying material models and input parameters used.Therefore, it is crucial to ensure that the modeling results are validated with experimental data to ensure their reliability and accuracy.
Estimating the life time of a wind turbine blade or any composite structural element when subjected to fatigue or creep deformations is a difficult undertaking.While metal failure is defined as the 1293 (2023) 012001 IOP Publishing doi:10.1088/1757-899X/1293/1/012001 2 separation of the metal into multiple pieces, composite failure follows a series of localized damage effects such as ply cracking, delamination, fiber/matrix debonding, and fiber fracture.These damage mechanisms are interactive and involve statistical uncertainty, leading to a slow deterioration of composite performance that ultimately results in catastrophic failure.To create a useful design methodology for predicting the lifespan of composite structures, it is essential to be able to predict the initiation and progression of damage leading to failure with confidence.Over the years, the Mechanics of Materials and Structures (MMS) research group at Ghent University has undertaken both theoretical and experimental activities to develop the necessary components for a progressive fatigue damage model of laminated composite structures.In parallel, the group has developed a BladeMesher software, that can generate shell, solid and hybrid shell/solid meshes of a full wind turbine blade in a few minutes, without any manual operations.Such high-quality meshes, providing converged local stress fields, are essential to be combined with fatigue damage models.
This paper primarily focuses on some of the physics-based developed models, along with the experiments performed that were supported by microscopic observation of damage mechanisms.Additionally, the paper outlines the remaining steps required to increase confidence in implementing the models for larger structural components while maintaining computational efficiency.

Modeling of ply cracking under multiaxial fatigue loading
The occurrence of stress concentration caused by geometrical discontinuities (such as holes or free edges) or manufacturing defects can initiate and interact with various damage mechanisms.However, in areas where these discontinuities are absent, ply level damage progression in laminated composites typically begins with the formation of ply cracks.Ply cracking degrades all the in-plane, out-of-plane and bending laminate stiffnesses and changes its coefficients of thermal expansion, moisture absorption, and the natural vibration frequency.Although ply cracking may not be a critical damage mechanism in most applications, it can significantly degrade the stiffness of laminates, such as Glass/Epoxy laminates used in wind turbine blades, and increase deformation.In extreme cases, this can lead to blade failure and impact with the tower.Ply cracking can also trigger the development of other damage mechanisms that are critical from the perspective of final fracture.
Three distinct modeling ingredients are typically required for prediction of ply cracking initiation, propagation and its effect on laminate properties under multiaxial fatigue loading conditions.It is highly challenging to predict the number of cycles required for the initiation of a measurable-sized ply crack, due to the involvement of intricate microscopic mechanisms that contribute to the formation of small ply cracks.These mechanisms operate at the interfaces of fiber and matrix, causing debonded areas that can grow in various directions, including through-thickness, width, and between plies, leading to ply cracks and delaminations.Additionally, the variable quality of the fiber/matrix interface, uneven spacing of fibers together with non-uniform distribution of residual stresses contribute to significant scatter in ply crack initiation times, requiring a statistical approach.To address this problem, a practical solution is to consider the ply crack initiation phase of fatigue damage as a known parameter (e.g. a negligibly small crack density if the laminate is not damaged), allowing predictive modeling to focus on ply crack growth and catastrophic failure onset.

Ply cracking fatigue evolution law
Expanding on the idea proposed in reference [4,5], an energy balance approach can be applied to predict unstable fatigue crack growth [6].To begin with, an energy balance formulated based on effective laminate properties needs to be considered for the region around the crack tip during cyclic deformation.Although most of the cyclic deformation generates localized heat, some of it also contributes to the degradation of the material in the fracture process zone.As fatigue crack growth can be attributed to the dissipation of energy in the failure process zone at the crack tip, energy concepts can be used to define the range of effective (laminate level) energy release rate, ΔG, as follows: , where where GI, GII and GIII are the effective strain energy release rates at the maximum load of the cycle for different deformation modes and Gmin is the effective energy release rate at the minimum load of the stress cycle.The parameter ρ is the effective ply crack density.It is essential to define an effective ply crack density for fatigue loading because, unlike quasi-static loading, ply cracks may not propagate unstably through the width or ply thickness.This definition suggests that the gradual growth of numerous small ply defects is equal to a consistent density ρ of fully developed ply cracks through the thickness and width, from a macroscopic perspective.The definition provided above for effective energy release rate range has two valuable attributes.Firstly, it considers the impact of R-ratio or the maximum and minimum applied stresses.Another desirable feature of this definition is that G0 →as min  →  .Finally, a fatigue crack growth law can be assumed of the form where CI, CII, CIII, nI, nII and nIII are material constants that should be calibrated using fatigue experiments on laminates with microscopic ply crack quantifications.Furthermore, the fatigue energy release rate threshold, denoted by th G  , represents the point at which fatigue crack growth ceases.This value may also vary depending on the deformation mode (shown as a single value for simplicity) and is reached when the effective energy release rate falls below it, due to laminate stiffness reduction.Experimental evidence, such as [7], suggests that cracks in off-axis plies do not grow continuously and reach a limit in number at a certain cyclic load level.This number of cracks remains relatively constant for the rest of the specimen's lifespan.The equation for fatigue ply crack growth in Eq. (2) bears resemblance to the widely-used empirical Paris Law for fatigue crack growth in metals.However, it is important to note that there is no clear physical justification for choosing this particular form of the fatigue ply crack growth law in Eq. (2), other than the fact that it is expressed in terms of energy parameters that take into account the complex triaxial stress fields present at the tips of ply cracks and laminate anisotropy.

Stress analysis in the presence of ply cracking
A crucial modeling tool is yet required to predict the stress distribution in the existence of ply cracks, which aids in computing the effective energy release rate range in various deformation modes and estimating the influence of fatigue-induced ply cracks on the properties of the laminate.Ply cracks tend to manifest in large numbers, often interacting with other cracks in the same or different plies (Figure 1).Considering the complex patterns of ply cracks that appear during the fatigue loading, the use of Finite Element Method (FEM) is impractical.Additionally, the stress distribution must be continuously updated throughout the fatigue loading process as a result of crack propagation.Furthermore, numerous calculations are required to calibrate the parameters in Eq. ( 2) and during the design process, different laminate lay-ups and ply thicknesses must be assessed.As a result, a robust stress analysis modeling tool is required that not only matches the accuracy of FEM but also provides a feasible run-time for the stress field, contingent upon the crack patterns and laminate specifications.To meet this demand, the variational approach [5] has been extensively developed over the years to handle intricate crack patterns in both symmetric and non-symmetric laminates under multiaxial in-plane[6], out-of-plane[7], and bending [8] loading conditions.It considers the impact of residual stresses and non-uniformity in crack distributions, which could be caused by the presence of manufacturing defects.The variational model involves (Figure 2) incorporating the geometry of the damaged state as a perturbation to the undamaged geometry.A limited set of assumptions are then made about the 3D stress field, and equilibrium equations are solved under periodic boundary conditions that correctly represent the stress states.The reader should refer to Ref. [4,5] to find detailed information about the calculation of ply cracking energy release rate in terms of laminate properties and the assumptions made.The run time of these models is a split second and extensive validation with highly refined FEM models has been done by UGent-MMS.

Modeling of delaminations at the tip of ply cracks under multiaxial fatigue loading
Ply cracks typically trigger delaminations due to stress concentrations at their tips.Delamination is a significant form of damage that worsens the effects of ply cracks, leading to largely reduced stiffness and stress redistribution.This damage is particularly problematic in multiaxial loading conditions where out-of-plane deformations occur, as it significantly impacts out-of-plane laminate properties.Delamination can initiate ply cracks in adjacent plies and cause fiber breakages, leading to catastrophic failure.Its detection during inspection is difficult, making it a challenging damage mechanism.
As previously mentioned, cracks in a ply reach saturation at a certain cyclic load level, but this does not imply that delamination only occurs after ply crack saturation.Recent experiments on glass/epoxy cross-ply laminates have revealed that micro-scale delamination (100 μm) occurs immediately after an off-axis crack fully penetrates the ply thickness [9].This is attributed to high stress fields caused by the presence of ply cracks at the interface, which include a significant peeling component, near the crack tip.Therefore, considering a very small preexisting delamination length at the tips of ply cracks can be helpful to focus on fatigue delamination growth.

Delamination fatigue evolution law
After initiating at the tip of a ply crack, a delamination propagates during cyclic loading.It is reasonable to assume that the delamination spreads perpendicular to the off-axis crack faces in a self-similar manner, at least when far from free edges.In comparison to fatigue ply crack growth, the development of delaminations at the ply crack tips is more akin to the growth of a fracture mechanics crack, where the delamination increment during each loading cycle is infinitesimal.Thus, it is reasonable to consider a new delamination growth law that is analogous to the empirical Paris Law.This would involve characterizing the rate of delamination growth in terms of delamination energy release rate range where CD, and nD are material constants that should be calibrated using fatigue experiments on laminates with microscopic quantification of delamination and where GI, GII, GIII concern different deformation modes.Equation (3) suggests that the total energy release rate is utilized to characterize delamination growth primarily to minimize the number of calibration parameters.However, it is important to consider a similar discussion as that for the case of ply cracking fatigue growth.As long as it is feasible to calculate the energy release rate for delamination under various deformation modes, it is possible to adopt any form of delamination fatigue evolution law, provided there are sufficient well-monitored experiments to calibrate the necessary inputs.

Stress analysis in the presence of delamination at the tips of ply cracks
It is exceedingly challenging to calculate the stress field, energy release rate range, and subsequent degradation of laminate properties in the interactive presence of both ply cracks and delaminations.The use of finite element method (FEM) for this purpose would be computationally inefficient due to the necessity of mesh refinements at the tips of both ply cracking and delamination.This problem would have been less severe if only a few ply cracks and delaminations existed in a single ply and/or interface.However, this complex damage pattern becomes increasingly intricate during fatigue loading as more cracks and delaminations appear even in adjacent plies.The variational models previously introduced for predicting the stress field in the presence of ply cracking only have recently been significantly expanded to accurately account for the impact of delaminations.This model [10] is applicable to general symmetric laminates subjected to general triaxial loading (see Figure 3).It is important to recognize that while Eqs.( 2) and (3) display the energy release rates as a function of effective crack density and delamination length, they are interdependent.In high cycles, both damage mechanisms require stress transfer to occur simultaneously.The reader should refer to Ref. [10] for detailed information about the model and its assumptions.

Modeling of ply cracking and delamination in multiple off-axis plies and interfaces
As previously discussed, calculating the energy release through stress transfer analysis becomes increasingly challenging when dealing with complex patterns of ply cracking and delamination, particularly when they occur in multiple plies and interfaces.Construction of a unit cell with several cracks in different orientations is almost impossible due to the requirement of applying periodic boundary conditions at the same time in different orientations.However, the evolution laws presented in Eq. ( 2) and (3) offer a valuable advantage in that they utilize laminate effective properties to accurately capture all the details of these defects without relying on the stress field of each individual crack.
In order to handle complex damage patterns involving ply cracking and delaminations in various orientations, the variational models described earlier are utilized in tandem with a homogenization approach.Through homogenization, a cracked ply with or without its delaminated interfaces can be replaced with an uncracked ply with perfect interfaces possessing adjusted elastic properties.As a result, the actual properties of the damaged laminate and the homogenized laminate remain the same.This technique of modified effective properties of a cracked ply can expand the scope of stiffness reduction models to tackle ply cracks and delaminations occurring in multiple orientations.In reference [4], it is demonstrated that to effectively homogenize ply cracking effects in symmetric laminates, changes in four parameters are necessary: the in-plane ply transverse stiffness E22, in-plane ply shear stiffness G12, out-of-plane transverse shear stiffness G23, and a coupling term φ26, specifically for ply cracking.Analytical solutions [11] can be derived to express these parameters in terms of laminate geometry, average crack spacing considering uniform and non-uniform distributions [12].When dealing with more intricate damage patterns that involve delamination, and when subjected to complex bending loads, a greater number of parameters are required for homogenization.In such cases, deriving simple analytical solutions is not feasible.To address these scenarios, neural networks are utilized and trained using data generated by running variational models [3, 10,13] for multiple cases.The outcome of the trained neural 1293 (2023) 012001 IOP Publishing doi:10.1088/1757-899X/1293/1/0120017 networks is then used to calculate both the laminate properties and energy release rate range as a function of effective crack density and delamination length.

Fatigue experiments with in-situ microscopic observation of damage mechanisms
In-situ fatigue damage quantification in composite laminates is cumbersome specially in non-transparent thermoplastic or carbon fiber laminates, however, it provides vital information for developing reliable fatigue crack and delamination growth laws and for their validation.Recently, the advances in optical measurements, namely Digital Image Correlation (DIC), are utilized to detect automatically crack density in glass/epoxy [14] and glass/propylene [8] multidirectional laminates with different lay-ups under quasi-static [15] and fatigue loading conditions.
In order to measure full-field strains and assess the extent of damage mechanisms in specimens, two 3D-DIC (Stereo Digital Image Correlation) systems are utilized, as shown in Figure 4.The experimental setup involves capturing in-situ images synchronized with data acquired from the testing machine during quasi-static and fatigue experiments.DIC analysis using strain and displacement fields is performed to detect cracks and obtain optimal results.Subsequently, an automated program is developed to evaluate crack density and weighted crack density from the edge and top of the specimen, respectively.The details can be found in Refs.[8,14,15].

Results and discussions
The first step for validation of the approach is to assess its quality in predicting the laminate effective properties as a function of fatigue driven ply cracks.For comparison with experiments, a [02/902]s laminate made of Glass/Propylene, which was tested in UGent-MMS [8], is first considered.This laminate was tested under uniaxial fatigue loading with R=0.1, f = 2 Hz, at two different amplitude loads σA=115 MPa, 150 MPa.For each case and at specific loading cycles, the in-plane axial stiffness EA and in-plane axial Poisson's ratio υA were measured together with the average value of effective crack density in the off-axis ply while delamination was reported to be negligible.We have used exactly the reported average crack densities [8] together with the provided input material properties [8] to predict elastic constants as a function of loading cycles.Figs.5a and b, show, respectively, the variation of the in-plane axial stiffness and Poisson's ratio.It is now useful to focus on the prediction of effective ply cracking density evolution as a function of applied cycle using the fatigue crack growth law in Eq. (2).To accomplish this, experimental data from the [02/902]s laminate tested under 115MPa amplitude loading are utilized for model calibration.It is important to note that in a cross-ply laminate under uniaxial loading, ply cracking growth at the homogenized level is exclusively in mode I. Therefore, a single experiment is adequate for determining CI and nI in Eq. ( 2).The calibrated evolution law is subsequently employed to predict the fatigue evolution of ply cracks under higher amplitude loading conditions, as depicted in Figure 6.To further compare with experimental data, particularly when both ply cracking and delamination occur during fatigue loading, two glass/epoxy laminates -a [02/904]s and a [02/454]s -previously tested in Ref. [9,16], were examined.These laminates were subjected to uniaxial fatigue loading with various amplitude loads.Normalized in-plane axial stiffness (EA/EA0) and normalized in-plane axial Poisson's ratio (υA/υA0) were measured for each sample at specific loading cycles, as well as the average crack density in the off-axis ply and the length of delamination at the interface between the off-axis ply and the 0 o ply.The reported average crack densities and delamination lengths [9,16] , along with the provided input material properties, were used to forecast elastic constants as a function of loading cycles. .As seen in Figure 7a, the variational approach predictions for the in-plane axial stiffness of the [02/904]s samples are in good agreement with the experimental results, despite some scattering between the experimental results for different samples of the same type.Additionally, the predictions in Figure 7b correspond well with the experimental outcomes for the [02/454]s laminate.
In Ref. [17], experimental characterization of ply cracking propagation was conducted in the off-axis plies of different laminates such as [0/902]s, [0/752]s, [0/602]s, and [0/753]s, subjected to uniaxial cyclic loading with varying amplitudes.Consequently, the crack density in the off-axis plies was quantified with respect to the applied cycles.For unbalanced laminates, ply cracks grow under mixed mode conditions, while for cross-ply laminates, cracks form under mode I deformation during uniaxial loading.Hence, in order to fully calibrate the model, experiments on both the cross-ply laminate and an arbitrary laminate with an off-axis ply other than 90-degree are required.The experimental crack density measurements on the [0/902]s and [0/752]s samples were utilized to calibrate CI, nI and CII, nII in Eq. (2), respectively.Once the model is calibrated, it can predict the crack density evolution in other laminates.Ply crack modeling in these laminates is illustrated in Figure 8, showing that the model can accurately predict the crack density evolution under cyclic loading when cracks propagate under mixed mode loading conditions.
Further validation of the model is required, particularly in predicting delamination fatigue growth and ply cracking fatigue threshold together with capability of neural networks in prediction of complex damage patterns.

Outlook
The prediction of damage mechanisms in complex geometries presents a significant challenge.One way to address this challenge is by using finite element techniques to analyze complex structures, such as wind turbine blades.This paper's developments can serve as a useful test data at the coupon level, providing inputs for anisotropic continuum damage models that can be integrated into FEM analyses.The benefits of using the methods described in this paper over actual coupon level experiments are not limited to cost and time efficiencies.Semi-analytical approaches can provide information on the damage state of laminates under various loading conditions, including those that are difficult or impossible to simulate experimentally.To develop appropriate constitutive laws that can be used as an anisotropic continuum damage model in FEM analyses at the structural level, it is also important to generate highquality meshes, that can accurately resolve the local stress geometries.
Over the past years, UGent-MMS has developed a BladeMesher software [18] (Figure 9), that can generate in a fully automated way high fidelity FEM models of full wind turbine blades.The models can be fully meshed with (linear or quadratic) solid elements (Figure 10), shell elements or a hybrid combination.Cohesive elements can be embedded in the adhesive bondlines to simulate delamination/debonding in the bondlines.
The BladeMesher software [18] has been successfully used for (i) simulation of static certification tests, [19], (ii) fluid-structure-interaction (FSI) simulations with a CFD solver [20,21], (iii) accurate torsional stiffness prediction of large blades [22], and (iv) global-local analysis of adhesive and bolted joints.This software platform opens up the opportunity to couple the advancements of physics-based, but extremely computationally efficient fatigue models with the structural analysis of large wind turbine blades.

Figure 4 .
Figure 4. Fatigue setup and images taken by the cameras on the front and edge surfaces [8].
Figure7aillustrates the reduction of the normalized in-plane axial stiffness (EA/EA0) for two [02/904]s laminate samples tested under uniaxial cycling load with amplitude ,max 120 results are also demonstrated for each sample separately due to slight variations in crack densities and delamination lengths reported for each sample at different cycles.Figure7billustrates the variation of the normalized in-plane axial stiffness (EA/EA0) and normalized in-plane axial Poisson's ratio (υA/υA0) for the [02/454]s laminate tested under uniaxial cycling load with

Figure 7 .
Variation of normalized (a) in-plane axial stiffness EA/EA0 in two samples with [02/904]s lay-up (b) In-plane axial stiffness EA/EA0 and in-plane Poisson's ratio υA/υA0 in one sample with [02/454]s lay-up, as a function of cycles which all were tested under a cyclic load

Figure 8 .
Figure 8. Variation of crack density as a function of cycles in a carbon/epoxy laminates with different lay-ups under a cyclic loading with different amplitude at R=0 .