849830 WIE0010.1177/0309524X19849830Wind EngineeringBoudounit et al. research-article2019 Research Article Structural analysis of offshore wind turbine blades using finite element method Wind Engineering 1–13 © The Author(s) 2019 Article reuse guidelines: sagepub.com/journals-permissions https://doi.org/10.1177/0309524X19849830 DOI: 10.1177/0309524X19849830 journals.sagepub.com/home/wie Hicham Boudounit1,2 , Mostapha Tarfaoui2,3, Dennoun Saifaoui1 and Mourad Nachtane1,2 Abstract Wind energy is one among the most promising renewable energy sources, and hence there is fast growth of wind energy farm implantation over the last decade, which is expected to be even faster in the coming years. Wind turbine blades are complex structures considering the different scientific fields involved in their study. Indeed, the study of blade performance involves fluid mechanics (aerodynamic study), solids mechanics (the nature of materials, the type of solicitations …), and the fluid coupling structure (IFS). The scope of the present work is to investigate the mechanical performances and structural integrity of a large offshore wind turbine blade under critical loads using blade element momentum. The resulting pressure was applied to the blade by the use of a user subroutine “DLOAD” implemented in ABAQUS finite element analysis software. The main objective is to identify and predict the zones which are sensitive to damage and failure as well as to evaluate the potential of composite materials (carbon fiber and glass fiber) and their effect on reduction of rotor’s weight, as well as the increase of resistance to wear, and stiffness. Keywords Wind turbine blade, composite materials, FEA, BEM, damage, Hashin’s criterion Introduction Electricity production sector is one of the most emitting sectors of carbon dioxide and greenhouse gas, due to the use of fossil fuels and coal. In order to reduce those emissions, industrialized countries have agreed on common objectives; renewable energy sources are an inevitable part of the solution, particularly wind energy that presents itself as a real reservoir of energy, and has been identified as one of the most promising green energy. It does not require fuel and its operating cost is low, it is even one of the cheapest among the renewable energies when installed in a well-ventilated site. Wind turbine blades’ (WTBs) structures are designed in a way to recuperate the kinetic energy from wind and transform it to mechanical energy, which is converted to electrical energy with the help of a generator. WTB are subjected to aerodynamic, centrifugal, and gravity loads, which impose large fatigue stresses on the moving rotor, notably at the transition region of the blade (Tarfaoui et al., 2017). Considering the maintenance constraints for renewable marine energies, durability and stiffness are essential for the offshore WTBs, and within this framework, composite materials will play the main role in the development of wind energy conversion systems (Schubel and Crossley, 2012; Schubel et al., 2013). Due to their advantageous mechanical properties, composites materials are used in all field such as energy, aeronautics, and automotive, but the lack of knowledge about the impact and fatigue characteristics of those materials leads to over-dimension the structure and subsequently the manufacturing will be expensive (Arbaoui et al., 2016a, 2016b, 2016c; Laurens et al., 2015; Shah and Tarfaoui, 2014). WTBs are made of composite materials with glass fibers and/or carbon fibers. Due to their large size, they are assembled using adhesive bonding. The thickness of the adhesive can affect the structural behavior. This is due to their very 1Laboratory for Renewable Energy and Dynamic Systems, FSAC–UH2C, Casablanca, Morocco CNRS 6027, IRDL, ENSTA Bretagne, Brest, France 3University of Dayton, Dayton, OH, USA 2UMR Corresponding author: Hicham Boudounit, UMR CNRS 6027, IRDL, ENSTA Bretagne, F-29200 Brest, France. Email: mostapha.tarfaoui@ensta-bretagne.fr 2 Wind Engineering 00(0) Table 1. Mechanical properties of the materials used (Nachtane et al., 2017; Rappolt, 2015). Material ρ ( kg / m3 ) E1 (GPa) E2=E3 (GPa) ν12 ν13=ν23 G12=G13 (GPa) G23 (GPa) Comp-G Comp-C 1900 1600 48.16 123.52 11.21 6.52 0.274 0.321 0.096 0.35 4.42 2.5 9 2.3 Material X T (MPa) X C (MPa) YT (MPa) YC (MPa) SLT = STT’ (MPa) Comp-G Comp-C 1021.3 1429 978 530 29.5 41 171.8 145 35.3 83.4 Material Gft (N/mm) Gfc (N/mm) Gmt (N/mm) Gmc (N/mm) Comp-G Comp-C 0.3 12.5 0.3 12.5 0.6 1 0.6 1 low stiffness, which deforms before showing a visible crack. Thus, the higher deformation can “hide” the crack where it is not visible, but the adhesive is already plastically deformed enough to not be able to resist the crack propagation (Shah and Tarfaoui, 2016a). Another study carried out by Shah and Tarfaoui (2014) shows that there is a relationship between damage on long fiber polymer composite material and the heat generated. They concluded that any sharp decrease in rigidity or increase in heat generation is signs of the material fiber-matrix interfacial crack growth, in addition to hysteresis in the polymer matrix and which was experimentally validated. The main objective of this project is to develop an effective process for the design of structurally optimized turbine blades, reduce the risk of blade failures, and maximize aerodynamic efficiency and structural robustness, and thereby reducing the blade mass and cost of total production. This article deals with determining the mechanical performances and structural integrity of offshore WTBs using ABAQUS finite element analysis software, evaluating the potential of composite materials and the effect of layup parameters, as well as defining the size and shape of the damaged area and the maximum displacement at the tip of the blade. Structure and materials properties Materials Due to their mechanical performance, composite materials offer new prospects for the renewable energy industry and specifically wind energy. However, the variability of their behavior, due in particular to the presence of initial or induced microscopic defects in service, is an important constraint to their development (Tarfaoui and Akesbi, 2001a, 2001b; Tarfaoui et al., 2001). The majority of wind turbine components are made of GFRP (glass fiber-reinforced polymer), due to its ease of use, nonrequirement of high technology, reasonable cost, and elastic qualities. Nowadays, carbon fibers have taken over, as they are more robust and lighter. Moreover, its resistance to fatigue is much more important (Kong et al., 2005; Tarfaoui et al., 2007), but the lack of knowledge about the behavior of composite materials induces too many additional costs in the fabrication of structures. The structural design requires high material strength, high material stiffness, high fatigue resistance, and low density. For the rest of this study, the following names will be adopted for the two composites: Comp-G: glass/epoxy composite, Comp-C: carbon/epoxy composite. The materials used in the simulations have the properties classified in Table 1. Failure mode: Hashin’s criterion The damage in composites materials occurs in two phases: damage initiation and damage evolution. Damage response. The model suggested by Matzenmiller et al. (1995) and Warren et al. (2016) has been adjusted the way to compute with the degradation of the stiffness matrix coefficients σ = M σ (1) 3 Boudounit et al. where M is the damage operator that has the following form (Lapczyk and Hurtado, 2007)  1 1 − d f   M = 0    0   0    0   1   1 − d s  0 1 1 − dm 0 (2) d f , d m , d s are internal variables that characterize the fiber damage, matrix damage, and shear damage, respectively. The damage compliance matrix is given by equation (3; Lapczyk and Hurtado, 2007) 1    1 − d f E1  v H =  − 12 E1    0   ( − )      0   1  (1 − d s ) G  v21 E2 0 1 (1 − d m ) E2 0 (3) The damage stiffness matrix, C, has the following form (Lapczyk and Hurtado, 2007) ( )  1 − d f E1  C =  1 − d f (1 − d m ) E2 v12   0  ( (1 − d ) (1 − d ) E v f ) m 1 21 (1 − d m ) E2 0    0  D (1 − d s ) G   0 (4) where ( D = 1− 1− d f ( )( ) (1 − d ) v m ) 21v12 >0 d s = 1 − 1 − d ft 1 − d fc (1 − d mt ) (1 − d mc ) (5) (6) E , v are properties of materials classified in Table 1. Damage initiation. Damage initiation refers to the occurrence of damage in a composite material (Camanho, 2002; Nachtane et al., 2017). Hashin’s criterion uses more than one constraint component to evaluate the different modes of failure. The failure modes in Hashin’s (1980) criterion are linked to fiber and matrix failure modes, and implicate four different modes of damage initiation for composite material: fiber tension (HSNFTCRT), fiber compression (HSNFCCRT), matrix tension (HSNMTCRT), and matrix compression (HSNMCCRT). According to Nachtane et al. (2017), the damage initiation criteria “Hashin’s criteria” for fiber-reinforced composites that are integrated on ABAQUS are based on Hashin (1980) theory, and Hashin and Rotem (1973) theory. f1, f2, f3, and f4 represent the four modes of failure included in Hashin’s criteria (Table 1), and the strength values, are referred to the Hashin’s failure criteria in Table 2. Damage evolution. Damage evolution describes the rate of degradation of the material stiffness once the corresponding initial criterion has been reached. The evolution of damage variable use the critical energy release rates, which depend on the materials (Matzenmiller et al., 1995; Warren et al., 2016), see Figure 1. The modeling approach is a generalization of that used to model cohesive elements, which is based on the work of Davila and Camanho (2002), and the evolution law is based on the energy dissipated during the process, and we assume that materials are linear softening. 4 Wind Engineering 00(0) Table 2. Hashin’s failure criteria (Laulusa et al., 2006). Failure mode and condition Failure criteria Tensile fiber failure for σ11⩾0 f1 ≥ 1 2   2    σ 11   + α ×  σ 12  = f1 where 0 ≤ α ≤ 1  SLT   XT        Compressive fiber failure for σ11<0 f2 ≥ 1 2     σ 11  = f 2  XC    Tensile matrix failure for σ22+σ33>0 f 3 ≥ 1    σ 22  YT  Compressive matrix failure for σ22+σ33<0 f 4 ≥1    σ 22  2STT  2     +  σ 12   SLT   2 2   = f3   2        +  YC  −1 × σ 22 +  σ 12   2STT   YC  SLT    2   = f4   α determines the contribution of shear stress on fiber tensile, and σ 11 , σ 12 , σ 22 determines the contribution of shear stress on the stress tensor. Figure 1. Damage evolution. Table 3. Blade properties. Length (m) 48 Power (MW) Number of blades Weight (tons) 5 3 12.18 (glass fibers)–9.54 (carbon fibers) The damage evolution is done by equation (7) d= ( (δ δ eqf δ eq − δ eq0 δ eq f eq − δ eq0 ) ) (7) Structure and design The model studied in this article is a 48-m-long industrial WTB modeled with shell elements, deformable, with box-side spars. Table 3 gives the properties of the composite blade. By using carbon fiber, we reduce the blade mass by 2.64 tons, so for a three-blade wind turbine, the rotor must be lighter by 7.92 tons, Table 3. Then the use of carbon fiber allows to increase the size of the blade, as well as to economize on materials when dimensioning the main shaft of the wind turbine also on the tower, which serves as a support for the entire structure. Boudounit et al. 5 Figure 2. Profile form used (NACA 4424). Figure 3. Evolution of blade chord. Figure 4. Views of the blade. We used NACA 4424 profile for its well-known aerodynamic properties (Figure 2). It has a maximum camber of 4%–40% from the leading edge, with a maximum thickness of 24%. Figure 3 shows the evolution of the chord with the blade length. Figure 4 shows different views of the blade with various chord lengths and different twist angles. Figure 5 shows the stratification model followed for the blade, and Figure 6 shows the distribution of materials for the blade surface shell (a) and for the spar (b). Simulation of the behavior of an in-service WTB Boundary conditions The boundary conditions of “Encastre” type is applied to the composite blade at the root region, where the blade is fixed (Figure 7). 6 Wind Engineering 00(0) Figure 5. Stratification model. Figure 6. Distribution of materials and layup thickness: (a) blade surface shell, (b) spar. Figure 7. Encastre boundary condition. Figure 8. Distribution of aerodynamic load. Loads Aerodynamic loading. When the wind is moving and stopped by a surface, it transforms into a pressure, which acts on a surface of intrados and extrados. Figure 8 shows the aerodynamic loads, which are calculated by BEM (blade element momentum) for wind speeds varying from 10 to 75m/s (Table 4). Current wind turbines are designed to harness the wind to produce the maximum amount of energy at a lower cost throughout the year. The wind turbines start to turn when the wind reaches a starting speed, generally 14km/h (3m/s), and stop when the winds exceed the speed of stop, namely, 90km/h (25m/s). It would be possible to build small wind turbines capable of exploiting the slightest breeze, 7 Boudounit et al. Table 4. Wind speed cases. In service V (m/s) 10 15 20 25 Maximum pressure (Pa) 1,93,028 2,90,909 4,00,423 5,23,291 Off service V (m/s) 35 45 55 65 75 Maximum pressure (Pa) 8,15,665 11,78,139 16,16,680 21,34,665 27,34,019 Figure 9. Centrifugal load on the wind turbine blade. but they would produce little energy. Indeed, the amount of available wind energy is proportional to the cube of the wind speed, while the power produced is proportional to the square of the length of the blades. If the exploitation of low wind speeds is, therefore, not interesting (little energy, little power), this is not the case for winds that blow more than 90km/h. Prototypes are trying to convert the energy of these winds, but for the moment, the conventional wind turbines just put their blades in the flag (in the alignment of the wind) beyond 90km/h. There is indeed a risk of premature wear and accident in this case. The operating range of wind turbines, therefore, depends on a compromise between wind strength and technical operating constraints. The loads that vary with position were applied as pressure using a DLoad subroutine in ABAQUS software. Gravity loading. Once the materials have been defined, ABAQUS finite element analysis software will calculate the load according to layup parameters. Centrifugal load. The centrifugal force depends on the mass and the rotational velocity of the blade (Figure 9) and given by equation (8) Fc = 1 mω 2r 2 (8) “ ω ” is the rotational velocity, “r” is the distance of the blade element from the center of rotation. Equation (9) allows for calculating the centrifugal force at each of the elements along the blade’s length Fc = ω2 ms rdr 2 ∫ (9) Figure 10 shows the wind turbine in service under all loads, where the green arrows represent the centrifugal force and the violet, the aerodynamic force. In this article, one studied the mechanical behavior of an offshore WTB under different wind speed, and we are interested in following two loading cases: in service with a wind speed of 25m/s, in a storm, when the blades are in flag mode (do not turn), with a wind speed of 75m/s. 8 Wind Engineering 00(0) Figure 10. Wind turbine blade in service. Figure 11. Mesh size study. Mesh study In order to determine the size of the mesh or the optimal mesh that will be used during the simulations, we carried out a mesh study to determine the size interval for the convergent results. To consider all possible cases, we performed simulations with the same loads, but with different mesh sizes. This will allow us to establish an interval of optimal mesh size, which gives more precise results with less calculation time (Figure 11). Referring to this figure, one can notice that the mesh convergence is obtained for an element size of 300 mm, see Figure 12. According to Laulusa et al. (2006), shell structures are widely encountered in many engineering applications and are commonly used in the aerospace, aeronautics, and automobile industries. Shell finite elements are complex compared with other structural elements because they involve a dimensional reduction procedure and treat the initial curvatures. The shell element (the mesh is of S4R type for the entire blade except at the level of curvature at the tip of the blade, which is S3R type of small size to eliminate the mesh defects) gives the same level of precision as the solid elements with less substantial computational costs. It is for these reasons that it is better to use the shell element than the solid elements (Laulusa et al., 2006; Shah and Tarfaoui, 2016b). Results and discussion In this study, the Hashin failure criterion has been used to model the damage development (occurrence) in the blade. The blade laminated with glass fibers shows very important values of displacement at its tip, unlike the carbon fibers blade that shows a high resistance to damage and displacement, which is obvious when comparing the results of the two composite blades with different loads (2.613 ;m for Comp-G and only 0.813 ;m for Comp-C; Figure 13 and Table 5). Figures 14 and 15 show the damage occurred on Comp-G and Comp-C blades in the storm conditions (V=75m/s), respectively. The damage modes are localized and appear on the area near the root, and their values depend on the composite material used as well as on wind speed (Table 6). Boudounit et al. 9 Figure 12. Optimal mesh for the blade. Figure 13. Displacement at the tip of the blade (mm). The damaged areas are similar for both materials. However, a comparison of the results of Hashin’s criterion for fiber and matrix failure modes shows that at low velocity the damage values for both materials are almost equal, and are decreased due to the absence of centrifugal force, when the pitch system turns the blades out of the wind to keep the rotor from turning when wind speed reaches 25m/s (Table 6). It is important to note here that the damage condition is not satisfied for carbon fiber blade. The model shows very high stiffness and good resistance to damage in the fiber and matrix modes (Figure 13, Tables 5 and 6). Despite the high wind velocity, the Hashin failure criterion is satisfied for Comp-G model only after reaching 55m/s, which means cracking of matrix occurs, followed by breaking of fiber–matrix interface, and ends with fiber failure. This leads to glass fiber model failure (Figure 14 and Table 6). The zones with high-stress concentration, where the damage develops, are located in the internal plies near the root of the blade and at the level of spars. We conclude that the following: •• The damages that appear near the area of the assembly of the blades with the hub are the result of a transition in thickness between those zones. •• The damages that appear on the spar are due to the adhesive, which is not very resistant. Using carbon fiber with the suggested layup parameters, one reduced the weight of a three blades rotor by 7.92tons, therefore, reduced transport and installation cost. Furthermore, the blade gains 5 times in stiffness, 10 Wind Engineering 00(0) Table 5. Displacement at the tip of the blade. Wind velocity (m/s) Comp-G U (mm) Comp-C U (mm) 10 15 20 25 35 45 55 65 75 630.3 819.4 1019 1309 1571 2390 2430 2583 2613 125.1 178.4 217.2 324.8 327.5 497.6 549.5 670.6 813.1 Figure 14. Damage on wind turbine blade Comp-G, V=75m/s. reduced displacement at the tip of the blades and the damages, then provide to the structure a good resistance and durability (Table 7). Conclusion This study shows that numerical simulation can be considered as a powerful tool to reproduce the experimental approach; it offers the possibility of predicting failure modes and identifies the sensitive zones on a WTB by using finite element method. The design methodology developed to take into account the aerodynamic loads on the blade, which are estimated using blade element momentum theory that gives a good initial design capable of being refined in later iterations. The results of the simulations allow to identify and predict the sensitive zones to damage and failure, as well as to determine the effect of the layup parameters. Generally, the damage begins at microscopic scale first with matrix cracking, then breaking of fibers–matrix interface, delamination, and ends with fibers failure. 11 Boudounit et al. Figure 15. Damage on wind turbine blade Comp-C, V=75m/s. Table 6. Value of damage given by Hashin’s criterion. V (m/s) 10 15 20 25 35 45 55 65 75 HSNFCCRT HSNFTCRT HSNMCCRT HSNMTCRT Comp-G Comp-C Comp-G Comp-C Comp-G Comp-C Comp-G Comp-C 0.000562 0.001161 0.001746 0.004818 0.002572 0.005265 0.007932 0.01225 0.01870 0.0001918 0.0002920 0.0004648 0.001139 0.001275 0.002598 0.003478 0.003751 0.004208 0.0002431 0.0006893 0.003119 0.004585 0.001387 0.002772 0.003557 0.003945 0.004446 0.0001017 0.0002105 0.0007513 0.001585 0.0004280 0.0008580 0.001258 0.001917 0.002914 0.005735 0.01848 0.02831 0.07195 0.07375 0.1698 0.2561 0.2878 0.3387 0.00253 0.004171 0.006116 0.0167 0.01428 0.02678 0.03834 0.05754 0.09145 0.1158 0.1753 0.7153 0.9693 0.3909 0.7755 1 1.030 1.128 0.0135 0.02784 0.04181 0.1140 0.06398 0.1287 0.189 0.2901 0.4427 HSNFCCRT: Hashin’s criterion for fiber in compression; HSNFTCRT: Hashin’s criterion for fiber in tension; HSNMCCRT: Hashin’s criterion for matrix in compression; HSNMTCRT: Hashin’s criterion for matrix in tension. Table 7. Comparative of Comp-G and Comp-C blades, main parameters. Comp-G Comp-C Weight of the blade (tons) Maximum displacement (mm) at 25m/s Maximum displacement (mm) at 75m/s Stiffness (kN/m) 12.18 9.54 1309 324.8 2613 813.1 7.065 33.14 12 Wind Engineering 00(0) Each of the composite materials used has its qualities and limits, glass fibers are economically very attractive and are an interesting choice for medium-sized WTBs, but composites with carbon fibers offer excellent properties and good resistance to fracture under quasi-static and dynamic loadings. The finite element modeling of the proposed design shows that from a strength point of view, the suggested lamination model allows to both blades to perform adequately in normal service conditions “less than 25m/s,” furthermore, the blade laminated in carbon fibers can withstand a storm “75m/s” if the wind turbine is off. ORCID iDs Hicham Boudounit Mourad Nachtane https://orcid.org/0000-0002-7191-0242 https://orcid.org/0000-0001-7381-2331 Declaration of conflicting interests The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article. Funding The author(s) received no financial support for the research, authorship, and/or publication of this article. References Arbaoui J, Tarfaoui M, Bouery C, et al. (2016a) Comparison study of mechanical properties and damage kinetics of 2d and 3d woven composites under high-strain rate dynamic compressive loading. International Journal of Damage Mechanics 25(6): 878–899. Arbaoui J, Tarfaoui M and El Malki Alaoui A (2016b) Dynamical characterisation and damage mechanisms of E-glass/vinylester woven composites at high strain rates compression. Journal of Composite Materials 50(24): 3313–3323. Arbaoui J, Tarfaoui M, El Malki and Alaoui A (2016c) Mechanical behavior and damage kinetics of woven E-glass/vinyl ester laminate composites under high strain rate dynamic compressive loading: Experimental and numerical investigation. International Journal of Impact Engineering 87: 44–54. Camanho PP (2002) Failure Criteria for Fiber-Reinforced Polymer Composites. Porto: Department of Mechanical Engineering and Industrial Management Faculty of Engineering, University of Porto. Hashin Z (1980) Failure criteria for unidirectional fiber composites. Journal of Applied Mechanics 47(2): 329–334 Hashin Z and Rotem A (1973) A fatigue failure criterion for fibre reinforced materials. Journal of Composite Materials 7: 448–464. Kong C, Bang J and Sugiyama Y (2005) Structural investigation of composite wind turbine blade considering various load cases and fatigue life. Energy 30(11): 2101–2114. Lapczyk I and Hurtado JA (2007) Progressive damage modeling in fiber-reinforced materials. Composites Part A: Applied Science and Manufacturing 38(11): 2333–2341. Laulusa A, Bauchau OA, Choi JY, et al. (2006) Evaluation of some shear deformable shell elements. International Journal of Solids and Structures 43(17): 5033–5054. Laurens J-M, Ait Mohammed M and Tarfaoui M (2015) Design of bare and ducted axial marine current turbines. Renewable Energy 89: 181–187. Matzenmiller A, Lubliner J and Taylor RL (1995) A constitutive model for anisotropic damage in fiber-composites. Mechanics of Materials 20(2): 125–152. Nachtane M, Tarfaoui M, El Moumen A, et al. (2017) Damage prediction of horizontal axis marine current turbines under hydrodynamic, hydrostatic and impacts loads. Composite Structures 170: 146–157. Rappolt JT (2015) Analysis of a carbon fiber reinforced polymer impact attenuator for a formula SAE vehicle using finite element analysis. Available at: https://digitalcommons.calpoly.edu/theses/1426/ Schubel PJ and Crossley RJ (2012) Wind turbine blade design. Energies 5(9): 3425–3449. Schubel PJ, Crossley RJ, Boateng EKG, et al. (2013) Review of structural health and cure monitoring techniques for large wind turbine blades. Renewable Energy 51: 113–123. Shah OR and Tarfaoui M (2014) Effect of damage progression on the heat generation and final failure of a polyester–glass fiber composite under tension-tension cyclic loading. Composites Part B: Engineering 62: 121–125. Shah OR and Tarfaoui M (2016a) Effect of adhesive thickness on the mode I and II strain energy release rates—Comparative study between different approaches for the calculation of mode I & II Serr’s. Composites Part B: Engineering 96: 354–363. Shah OR and Tarfaoui M (2016b) The identification of structurally sensitive zones subject to failure in a wind turbine blade using nodal displacement based finite element sub-modeling. Renewable Energy 87: 168–181. Tarfaoui M and Akesbi S (2001a) A finite element model of mechanical properties of plain weave. Colloids and Surfaces A: Physicochemical and Engineering Aspects 187: 439–448. Tarfaoui M and Akesbi S (2001b) Numerical study of the mechanical behavior of the textile structures. International Journal of Clothing Science and Technology 13(3–4): 166–175. Boudounit et al. 13 Tarfaoui M, Drean JY and Akesbi S (2001) Predicting the stress-strain behavior of woven fabrics using the finite element method. Textile Research Journal 71(9): 790–795. Tarfaoui M, Gning PB and Collombet F (2007) Residual strength of damaged glass/epoxy tubular structures. Journal of Composite Materials 41(18): 2165–2182. Tarfaoui M, Nachtane M, Khadimallah H, et al. (2017) Simulation of mechanical behavior and damage of a large composite wind turbine blade under critical loads. Applied Composite Materials 25: 237–254. Warren KC, Lopez-Anido RA, Vel SS, et al. (2016) Progressive failure analysis of three-dimensional woven carbon composites in single-bolt, double-shear bearing. Composites Part B: Engineering 84: 266–276. Appendix 1 Notation E1: longitudinal Young modulus E2: transversal Young modulus E3: Young modulus along the thickness G12: shear modulus in 1-2 plane G13: shear modulus in 1-3 plane G23: shear modulus in 2-3 plane Gft, Gfc, Gmt, Gmc: damage evolution coefficients HSNFCCRT: Hashin’s criterion for fiber in compression HSNFTCRT: Hashin’s criterion for fiber in tension HSNMCCRT: Hashin’s criterion for matrix in compression HSNMTCRT: Hashin’s criterion for matrix in tension S LT : longitudinal shear strength STT ’ : transverse shear strength X T : longitudinal tensile strength X C : longitudinal compressive strength YT : transverse tensile strength YC : transverse compressive strength ν12, ν13, ν23: Poisson’s ratio