Study on composite bend-twist coupled wind turbine blade for passive load mitigation

Composite Structures 213 (2019) 173–189

Contents lists available at ScienceDirect

Composite Structures journal homepage: www.elsevier.com/locate/compstruct

Study on composite bend-twist coupled wind turbine blade for passive load mitigation

T

⁎

Jinge Chen, Xin Shen, Xiaocheng Zhu , Zhaohui Du School of Mechanical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China

A R T I C LE I N FO

A B S T R A C T

Keywords: Load alleviation Bend-twist coupling Composite blade VABS Wind turbine Shear wind

The wind turbine size is nowadays becoming increasingly larger as an effective approach to reduce cost of energy. Passive load control technique is introduced to alleviate the increasing loads on an upscaling wind turbine blade to improve its fatigue life. In this paper, the composite bend-twist coupled blade is investigated and utilized to mitigate the cyclic fluctuating loads in shear wind. The NREL 5-MW blade is firstly inversely redesigned of its composite layup configuration. The fiber on the spar cap is rotated away from the blade axis to implement the bend-twist coupling (BTC). An advanced aeroelastic model of the wind turbine blade is developed to conduct the study. The geometrically exact beam theory is used to account for the large deformation effects, while the tool VABS is utilized to generate the cross-section stiffness matrices. Meanwhile, the free wake lifting surface model is adopted to calculate the aerodynamic loads, which is physically more realistic than the traditionally used BEM method. Based the established model, the influences of BTC on the rotor aerodynamic performance are investigated. Load mitigation effects of the blades with various spar cap fiber angles are evaluated. The influences of coupling region and different wind speeds are also discussed.

1. Introduction The wind energy industry nowadays has growing interests in designing increasingly larger wind turbines due to economic considerations, because larger wind turbines can capture more wind energy and help to reduce the cost of energy. The diameter of modern utility-scale wind turbine rotors nowadays has been over 160 m, and even larger rotors with the size of 10–20 MW are under development to serve as the next-generation power units. However, as the rotor size increases, loads on the blades increase dramatically. Load control techniques, which enable to achieve better optimized material usage and blade fatigue life improvement, are therefore essential for large and super-large wind turbine design. There exist various load control methods and they can be generally divided into two categories: active and passive load control techniques. The active control methods are proved effective, but come at a price of actuator requirements which add cost and induce further uncertainty to the system reliability [1]. Hence, it is attractive to investigate passive load alleviation technologies that are free of actuators due to their simplicity, relatively quick response and cost-effectiveness [2]. The passive load control, which is also described as aeroelastic tailoring, is typically implemented by using bend-twist coupling (BTC) blades. This kind of self-adaptive blades are designed to passively twist under ⁎

bending loads in a way that changes the angle of attack and alleviates the load variations. The desired bend-twist coupling can be achieved either through the use of anisotropic composite material properties (material BTC) or by creating an offset between the local elastic axis and blade pitch axis (geometric BTC) [3–5]. The present work only focuses on the employment of material-based adaptive blades. Numerous studies have been conducted on the subject of material BTC wind turbine blades. The concept of using material-based adaptive blades to passively mitigate loads on wind turbines is firstly investigated by Lobitz and Veers [6–9] in Sandia National Laboratory (SNL). Significant reductions in fatigue damage under turbulent wind are observed. The BTC effect is considered in these works by introducing a coupling term to the stress-strain relation in the beam model, without relating the BTC magnitude to a specified composite layup. Experimental tests are afterward conducted [10–12] and prove that the BTC can be achieved by rotating the fiber angle of the unidirectional (UD) laminates in spar cap and/or skin away from the beam axis. Ong et al. [10] show that the off-axis ply orientation, the laminate material and the proportional volume of the anisotropy layers in the laminate are the three key parameters that influence the coupling magnitude. By numerical method, Hayat and Ha [2] further find out that the coupling magnitude and fatigue load reduction can be increased if the plythickness unbalance for a bi-axial (BX) laminate is added in the skin

Corresponding author.at: 800 Dongchuan Road, Minhang district, Shanghai 200240, China. E-mail address: [email protected] (X. Zhu).

https://doi.org/10.1016/j.compstruct.2019.01.086 Received 5 July 2018; Received in revised form 24 November 2018; Accepted 25 January 2019 Available online 30 January 2019 0263-8223/ © 2019 Elsevier Ltd. All rights reserved.

Composite Structures 213 (2019) 173–189

J. Chen et al.

structure problem to a combination of 1-D beam dynamics and 2-D cross-section analysis through equivalence of strain energy. Exact strain-displacement relations are used in the 1-D beam equations to account for large deflections. Correspondingly, the related 2-D finite element based cross-section analysis tool VABS is developed by Yu [32], which provides a 6 × 6 stiffness matrix for a given composite section of the equivalent beam. Blasques and his coauthors [33–36] conduct similar works and develop the sectional analysis software called BECAS. Similar nonlinear beam models are also developed by Bauchau et al. [37,38]. By using nonlinear beam theory and/or composite cross-section analysis method, several works are carried out to structurally model a composite wind turbine blade. Kim et al. [39] at DTU implement a newly developed anisotropic beam element into their nonlinear multibody code HAWC2 to enhance its capability to calculate composite blade dynamics and bend-twist coupling effects. Otero and Ponta [40] develop their own code based on a modified implementation of the VABS and nonlinear beam theory. A 40-m long composite test blade is analyzed of the vibrational modes and the stress loads in operations. A similar study is performed by Fleming and Luscher [41], in which they develop the nonlinear beam code NLBeam and use VABS tool to provide section stiffness. The CX-100 wind turbine blade is analyzed and both static and dynamic behaviors are obtained. The material BTC blade technique is also further investigated in wind turbine industry in recent years, by using nonlinear beam theories and/or cross-section analysis methods. Bottasso et al. [42] conduct a multilevel constrained structural optimization design for BTC blades considering various design requirements. The combination of passive coupled blade and active individual pitch control is also investigated. The code Cp-Max [43] is used to conduct the work, in which the code ANBA similar to BECAS/VABS is utilized to calculate sectional stiffness and the multibody simulator Cp-Lambda similar to HAWC2 is adopted to compute the dynamic responses. Another multibody aeroelastic code hGAST, combined with BECAS, is used by Bagherpour et al. [44] to investigate the material BTC blades. The bend-twist coupling coefficient distributions along the span are assessed and associated with the ply offset angles. Blade root flapwise bending fatigue and ultimate load reduction are finally achieved. Scott et al. [45] further compare the material BTC blade and the material/geometric combined BTC blade, in which work the software PreComp and GH-Bladed are used. It shows that the blade with both couplings displays scope for potential increases in energy yield. Sener et al. [46] recently conduct a study to evaluate the effect of spar cap fiber orientation of composite BTC blades on the wind turbine aero-structure performance, where the E-glass/epoxy blade and the hybrid E-glass/carbon/epoxy blade are compared. The tool VABS and nonlinear-beam-based aeroelastic software PHATAS are utilized to conduct the simulations. Results reveal that certain hybrid BTC blade designs are more effective in fatigue load mitigation. In the present work, a nonlinear structural model for analyzing composite BTC blades is also developed, by implementing the displacement-based geometrically exact beam theory proposed by Bauchau [29]. The software VABS is used to obtained the sectional stiffness matrices. However, it is worth noting that all the above works are based on BEM method to compute aerodynamic loads on wind turbines. BEM method is a simplified engineering model whose fidelity highly depends on various correction models, leading to limitations in several situations [47,23]. For example, it cannot account for the 3-D shape effects of the largely deformed blade [48]. The more physically accurate vortex method is an alternative, and is gradually taking the place of BEM to predict the wind turbine’s aerodynamic performance [49–55]. Therefore, in this paper, a free wake lifting surface model is adopted by combination with the nonlinear beam model to set up the aeroelastic model for analyzing the composite BTC blades. Differing from the previous works, the BTC blade is further explored of its capability of alleviating cyclic fluctuating loads under shear wind inflow instead of typically studied turbulent condition. To get a

layup, combined with ply-angle and ply-material unbalances. Vesel and McNamara [13] go further to investigate the blade optimization design for minimum COE by a simultaneous consideration of aerodynamics and BTC effects. In their work, the degree of BTC in the blade is taken as one of the design variables along with airfoil shapes, chord and twist distributions. The reductions in flapwise bending loads are assumed to decrease rotor cost through reduced material requirements. The material BTC adaptive blade technique is also applied in marine propellers and tidal turbines. Studies on marine propellers [14–16] have demonstrated that the composite adaptive blades can help to increase the energy efficiency and delay cavitation inception, particularly in spatially and/or temporally varying inflows and in off-design conditions. Similar passive control strategies are then introduced to marine turbines. Researches have demonstrated the potential of composite BTC blades for load reductions and power regulation of marine turbines by both numerical [17–19] and experimental methods [20–22]. By using a vortex method-finite element method solver, Nicholls-Lee et al. [17] find that a decrease of up to 12% in thrust and an increase of up to 5% in power capture can be achieved through using properly designed BTC blades in a tidal turbine. By adopting a boundary element method coupled with a finite element method, Motley and Barber [18] further study the capabilities of BTC adaptive blades under both instantaneous and long-term variable amplitude loading while considering practical design and operational restrictions. It is shown that passively controlled blades that increase energy capture would result in higher blade loads and increased required active control range and vice versa. Higher unsteady stress profiles due to the increased flexibility of the blades and the off-axis orientation of the fibers are also observed. Barber and Motley [19] afterwards use the same solver to show that the passively adaptive pitch-to-feather marine turbine blades can be used both to delay cavitation and to reduce cavitation volume over the blade surface. Compared with the marine propeller or tidal turbine blades, the modern wind turbine blade usually has a much larger aspect ratio and therefore is more often treated with a beam theory, which can provide enough accurate results if constructed properly [23]. The above mentioned numerical works in wind turbine applications calculate the blade responses based on traditional linear beam theory and account for the BTC effects by simply using a coupling coefficient as done in Lobitz’s method [6]. However, as the composite wind turbine blade typically has a complex internal structure, the section stiffness property changes due to fiber offset are difficult to be well captured in Lobitz’s method. Furthermore, the BTC adaptive blades are likely be characterized by relatively high flexibility and undergo large deflections in operations, which linear beam theory cannot account. To accurately represent the material BTC blade behaviors, advanced techniques must be introduced to calculate composite beam’s sectional stiffness properties, as well as the highly flexible blade’s dynamic responses. The 3-D finite element method has a high level of fidelity and is usually coupled with an advanced fluid solver to fully take effects in fluid–structure interaction (FSI) simulations of wind turbines [24–27]. However, many advanced fluid solvers, e.g. the Navier-Stokes equations based computational fluid dynamics (CFD) solver, are still currently impractical to be used in FSI simulations especially in primary design stage due to high computational resource requirements. Nonlinear 1-D beam models combined with 2-D sectional analysis methods are therefore introduced to fill the gap between model accuracy and computational efficiency, which is well suited to be coupled with a moderately improved and computationally practical fluid solver, for instance a vortex method solver. Nonlinear beam theories have gone through a long-term developing process in helicopter industry, which is well reviewed in [28,29]. A latest nonlinear beam technique is developed by Hodges [30,31] and is called geometrically exact beam theory or generalized Timoshenko beam theory, which is able to compute any initially twisted and curved, inhomogeneous, anisotropic beam with arbitrary cross-sectional geometries. This beam theory reduces the 3-D 174

Composite Structures 213 (2019) 173–189

J. Chen et al.

The 3-D strain energy of the cross-section or the strain energy density of the beam is now written as

comprehensive understanding of the BTC blades, influences of BTC effects on rotor aerodynamic performance across the full operating wind speeds are also investigated. The NREL 5-MW wind turbine blade is selected as the baseline. As this blade originally has no specified composite layup design, a composite layup configuration for the NREL 5-MW blade is also developed. This paper is organized as follows. The theoretical model, including the basic principle of VABS, the geometrically exact beam theory and the free wake lifting surface model, is briefly described in Section 2. The detailed composite layup design of the NREL 5-MW blade and the implementation of bend-twist coupling via spar cap fiber angles are introduced in Section 3, followed by the assessment of influences of BTC on rotor aerodynamic performance. The load mitigation effects by using the BTC adaptive blades under shear wind are discussed in detail in Section 5 and finally some conclusions are drawn.

U=

The aeroelastic model based on free wake lifting surface model and geometrically exact beam theory has been established and validated in our previous works [56–58]. A brief outline of the theoretical model is described in order to make this paper self-contained.

2U = εcT Aεc + 2εcT Bε′c + εc′ T Eε′c + 2εcT Dεc ″

(4)

These newly introduced matrices carry information on the material properties as well as the geometry of a given cross-section, and are obtained by 2-D finite element discretization of the section. The classical strain measures εc now should be transformed into Timoshenko-like strains. By substituting the transformation relations and dropping higher-order terms, the strain energy is finally expressed in terms of Timoshenko-like strain measures:

2.1. Variational asymptotic beam sectional analysis The sectional analysis starts from the blade kinematics representation. The kinematics of the rotating blade, including both the undeformed reference configuration and the deformed one, are described based on the 1-D curve line and the 2-D cross-sections in the fixed inertia coordinate system as shown in Fig. 1. The position of an arbitrary point within a particular cross-section in the deformed beam can be specified as

2U1 = ∊T A ∊ + 2 ∊T AQs′ + 2 ∊T APs + 2 ∊T B ∊′ + ∊′ T E ∊′ + 2 ∊T D ∊″ (5)

]T

where ∊ = [γ11 κ1 κ2 κ3 is the 1-D generalized Timoshenko-like strains associated with the Bi system, and s = [2γ12 2γ13 ]T is the column matrix of Timoshenko transverse shear strain measures also associated with the Bi system. In detail, γ11 is the sectional axial strain, γ12 and γ13 are transverse shearing strains, κ1 is the sectional twist rate, and κ2 and κ3 are bending curvatures. On the other hand, by definition, the Timoshenko-like strain energy should be written as

R (x1, x2 , x3) = r (x1) + u (x1) + x2 B2 + x3 B3 + wi Bi (1)

where r is initial position vector of the reference point, u is sectional translation vector, CBb is the sectional rotation tensor defining the elastic deformation along with u , C bI is the known rotation tensor given by initial twist and/or curvature, x2 , x3 are material coordinates in local frame Bi , and wi is warping displacement. The first step is to derive the 3-D strain energy expression for the deformed beam. Starting from the above kinematics description and based on the concept of decomposition of the rotation tensor [59], the 3-D strain field of the beam structure can be described in terms of classical 1-D generalized strain measure εc = [γ11 κ c1 κ c 2 κ c3]T and warping function (including the transverse shear deformation) and its derivative [32]:

Γ = Γh w + Γε εc + ΓR w + Γl w′

(3)

where the notation 〈·〉 = ∫s · g dx2 dx3 means integration over the cross section, and D is the 6 × 6 symmetric material matrix in the bi basis. The next step is to reproduce the above 3-D energy in terms of 1-D quantities, by using the variational asymptotic method (VAM). The VAM is a powerful mathematical method that can be used to find the 1D energy which approximates the 3-D energy as closely as possible. Firstly, the characteristic size of the section h is identified as the small parameter because h l < 1, h R < 1, where l is the characteristic wavelength of deformation along the beam axial coordinate, and R is the characteristic radius of initial curvature/twist of the beam. By expanding the unknown warping field in asymptotic series of h and discarding unnecessary high-order terms, the first-order approximation of warping is solved and the second-order asymptotically correct energy can then be written as

2. Theory

= r (x1) + u (x1) + CBb·(x2 b2 + x3 b3 + w (x1, x2 , x3))

1 T 〈Γ D Γ〉 2

X Z⎤ ε 2U = [ εT sT ] ⎡ T ⎡ ⎤ = εT Xε + 2εT Zs + sYsT ⎣Z Y ⎦⎣s⎦

(6)

where the matrices X , Y and Z compose the 6 × 6 cross-sectional stiffness matrix consistent with the equivalent 1-D generalized Timoshenko beam. The final step is to find out these matrices in such a way that the strain energy in Eq. can asymptotically approximates the energy in Eq. up to second order. The nonlinear 1-D equilibrium equations by Hodges [30] are used to achieve this goal. Detailed derivations can be found in the reference [32].

(2)

2.2. Geometrically exact beam theory As have discussed in last subsection, the 3-D strain energy of the deformed beam is rewritten in terms of six 1-D generalized strains, which can be re-expressed as force strains γ = [γ11 γ12 γ13 ]T and moment strains κ = [κ1 κ2 κ3 ]T . The 3-D structure like a wind turbine blade is in this way reduced to a 1-D beam problem. The remaining issue is to derive a set of motion equations to describe the beam dynamics. The geometrically exact beam theory takes the full geometrical nonlinearity effects into account and therefore is suitable for a beam undergoes arbitrarily large displacements and rotations. In this approach, the beam deformations are described by six variables including 3 rotation displacements and 3 translation displacements at each reference point. The rotation of the cross-section in three directions is further represented by a rotation tensor, as shown in Fig. 1 and Eq. (1).

Fig. 1. Curved beam in the reference and deformed configurations. 175

Composite Structures 213 (2019) 173–189

J. Chen et al.

The kinematics of geometrically beams was first proposed by Reissner [60]. Simo and Vu-Quoc [61] extended Reissner’s work to deal with 3D dynamic problems. Similar developments were proposed by Borri and Merlini [62], Danielson and Hodges [59], or Bauchau and Hong [38]. Based on different representations of finite rotation and numerical solving methods, many formulations have been proposed. The key issue is the fidelity of the strain-displacement relations. Danielson and Hodges [59] derived the geometrically exact kinematics and generalized strain expressions of an initially curved and twisted beam in 1987, which forms a basis for consequent nonlinear beam models. For a beam with shallow initial curvature, these generalized strains can be calculated by [29,63]

γ = r ′ + u′ − CBbC bIe1

(7)

∼ κ = (CBbC bI)′ (CBbC bI)T

(8)

Y

Far wake Lagrangian markers Centroid point Near wake Control Point

where u is the column matrix containing the measure numbers of u in the inertial basis, so that ui = u·ii ; similarly, ri = r ·ii . Note that e1 = [1 0 0]T , and ( )' denotes the derivative with respect to x1, while ∼ the tilde symbol ( ) for a column matrix represents a second-order, skew-symmetric tensor. Based on the above displacement and strain field, the kinetic energy K and strain energy U of the deformed beam can then be derived and substituted into the Hamilton’s principle to get the final governing equations: t2

∫t ∫0

l

[δ (K − U ) + δW ] dxdt = 0

1

Bound vortex X Z Fig. 2. Schematic of rotor coordinates and the lifting surface model.

2.3. Free-wake vortex lattice model

(9)

The lifting surface model assumes the flow field around the rotor as incompressible potential flow. The blade is divided into chord-wise and span-wise panels, with bound vortices set at the 1/4 length and control points at the 3/4 length at the middle of each panel. The free vortices trailed from the ends of the bound vortices at each panel extend downstream, forming sets of horseshoe filaments. The trailing vortices comprise the near wake sheet and far wake represented by a single tip vortex. The near wake extends from the bound vortices, and then rolls up into the far wake. The bound vortices on the blade and wake structures are illustrated in Fig. 2. The circulation strength of each horseshoe filaments is determined by the boundary condition on the control point at each panel, which implies that the incident velocity normal to the span-wise segment at the control point should be zero to prevent flow penetrate the blade:

where δW is the virtual work of the externally applied loads. In the present work, the displacement-based formulation developed by Bauchau [29,63] is applied. The detailed derivation can be found in the reference. Here we directly give the final governing equations of motion:

h ̇ − F ′ = fext ∼ g ̇ + u ̇ h − T ′ − (∼ r′ + ∼ u ′) F = mext

(10)

where fext and mext are externally applied force and moment. The no.

tation ( ) denotes a derivative with respect to time. The linear and angular momenta (h and g , respectively), and the sectional force and moment (F and T , respectively) are defined by the constitutive equations:

⎛ h ⎞ = M u̇ ω ⎝g⎠

()

(11)

(TF ) = C ( γκ )

(12)

(V∞ + Vind + Vrlv )·n = 0

(14)

where V∞ is the upstream far-field flow speed, n the lifting surface’s unit normal vector at the control points, Vind the induced velocity, and Vrlv the relative speed caused by blade rotation, etc.. The flow on each control point is influenced by the bound vorticity of each segment and all the vortices trailed behind the blade. In other words, the induced velocity equals to the sum of the induced velocities from the bound vortices, the near wake and the far wake:

where M and C are the beam’s 6 × 6 sectional mass matrix and stiffness matrix, respectively. ω is the section angular velocity determined by CBb . In the present work, the finite rotations (CBb and C bI ) are represented by means of vectorial parameterization, where a specified trielement vector p = [p1 p2 p3 ]T is used to describe rotation, e.g. CBb = CBb (p) . Hence we have six unknown quantities at each point on the beam, represented as q = [u p]T . The nonlinear governing equations of motion are solved by Newton-Raphson method, resulting in linearized equations in incremental formulations. These incremental equations are then discretized by finite element method, finally leading to the following discretized incremental equations:

 Δv + K  Δa + G  Δq = F  M

,

Vind = VB + VNW + VFW

(15)

where the Biot-Savart law is employed to determine the induced velocities. A time-marching free wake model, which is suitable for transient conditions, is used in this paper for the determination of the far-wake shape and thereby the velocity VFW . The vortex line is assumed as material lines convected through the flow field at the local velocity. A Lagrangian fluid particle is then used to describe the motion of a point on the far-wake filament:

(13)

drw = V (rw ) dt

, K  and F  are generalized mass matrix, gyroscopic matrix, , G where M stiffness matrix and generalized force matrix, respectively. And Δa , Δv and Δq are incremental acceleration, incremental velocity and incremental displacement, respectively. Generalized-α scheme is used to solve the derived ODEs.

(16)

where rw is the position of the control point on the far-wake, and V (rw ) is the local flow-field velocity at that point. The local velocity consists of the free stream velocity, the induced velocity and other sources: 176

Composite Structures 213 (2019) 173–189

J. Chen et al.

Vy Vx

N

L

Vrel D

T

Rotor Plane

rȍ

ș V’

Įe

Įe Fig. 4. Illustration of the cross-section internal geometry.

Vrel Vrel

Wind

aV’

only distributed structural properties of the blades are provided, lacking detailed composite layup information. A structural layup concept for the blade is needed to investigate the composite structure dynamics and bend-twist coupling effects via unbalanced fibers. Hayat [2] developed layups for a 5-MW blade by down-scaling the geometry and layup of the all-glass 13.2 MW SNL 100-00 wind turbine blade designed by Sandia National Laboratories [65]. Clearly, the equivalent sectional bending and torsion stiffness of the scaled blade are not necessarily consistent with NREL 5-MW blade. Specifically, a basic layup design for the NREL 5-MW blade is developed by Resor [66], based on the composite layup concept from the UpWind program [67] and layup and materials information from the SNL 100-00 blade. Although the design result is not a fully vetted blade design which is ready for manufacture, it provides a good starting point for more detailed and targeted investigations. Based on the design by Resor, the composite layup for the NREL 5MW blade is established in the present work. The typical internal geometry, as shown in Fig. 4, is used to build up the blade structure. The spar box is centered at the chordwise location generally with the maximum airfoil thickness. The spar cap is strengthened mainly by carbon fiber laminate, with a cap width of 600 mm. Foam laminate constitutes the core of shear webs, leading-edge (LE) panel, and trailing-edge (TE) panel. No core is applied at the two narrow regions with a width of 100 mm at the LE and TE connection positions. The thin trailing edge region is further reinforced with glass fiber laminate to avoid structure failure. Sandwich constructions are then formed at spar cap, LE and TE panels by bonding internal and external triax skins, while saertex skins are utilized in shear webs. Triax (TX) laminates with increased thickness are used to reinforce the blade root. The whole blade is coated with gelcoat at the outside surface. The material properties are listed in Table 2, where longitudinal, transverse and shear stiffness properties are denoted by E1, Eα and Gij , respectively, and Poisson ratio is denoted as vij . The layer thicknesses or ply numbers are then determined to build a fully composite blade. Our first objective is to design the layup such that the section stiffness distribution of the composite NREL 5-MW blade approximates the original design [64] as closely as possible, which can serve as a baseline blade for further study. For simplicity, the unidirectional carbon fibers in spar cap and glass fibers in TE reinforcement segment are oriented at 0 ° with respect to the span axis. Compared with the layer design in [66], the ply numbers are readjusted to better match the original stiffness distribution. Specifically, the ply numbers are reduced in carbon laminate of spar cap and foam cores of shear web and LE/TE panels. Fig. 5 shows the readjusted layer thicknesses along the span. This layup concept leads to a structurally redesigned NREL 5-MW composite blade which is referred as the baseline blade. Based on the above layup concept, the blade’s 1-D equivalent crosssectional stiffness distributions along the span are then obtained by using VABS. The 2-D meshes for airfoil sections at different spanwise positions are displayed in Fig. 6. These meshes are produced by preVABS code developed by Su et a. [68] as a pre-processor for VABS. The NREL 5-MW blade uses different airfoils shapes over the whole span. The TU-Delft family of airfoils is used at approximately two-thirds span, while NACA 64-series airfoils are used in the final one-third blade span.

rȍa’

Fig. 3. Schematic of velocity and force decompositions of a blade element.

(17)

V = V∞ + Vindv + Vex

where the induced velocity Vindv is comprised of itself and mutually induced velocities in the wake as well as the bound vortex on the blade. And Vex is other source of the velocity such as the unsteady term of the velocity in the wind. A rotor fixed coordinates is used for the definition of the positionrw , which uses the blade azimuthal location ψ as the temporal coordinate, and vortex age ς as the spatial coordinate. Then the left hand of Eq. can be written as follow:

drw ∂r dψ ∂r dς = w + w dt ∂ψ dt ∂ς dt

(18)

Once the circulation strength on the blade segments are solved, the aerodynamic forces and moments are calculated by table look-up according to effective angle of attack (AOA). The forces on the airfoil are illustrated in Fig. 3, where φ is the angle between the plane of rotation and the total relative velocity Vrel , θ is the sum of pitch angle and blade pre-twist angle, and the effective AOA is represented as α e . 3. Adaptive blade design In this section, the structural design for adaptive blade is detailed. The composite layup for the NREL 5-MW blade is firstly developed, based on which the bend-torsion coupling is implemented by offsetting the spar cap fiber angle to build the adaptive blade for further study in the next sections. 3.1. Composite layup for the NREL 5-MW blade The NREL 5-MW conceptual reference wind turbine is developed by US. National Renewable Energy Laboratory (NREL), and is used extensively by researchers in various topics as a system that represents the current state of the art in an offshore horizontal wind turbine system [64]. Some key parameters and settings in this paper are listed in Table 1. Detailed specifications of this utility-scale 5 MW wind turbine, including geometry (chord, twist, airfoils and pitch-axis locations), mass and stiffness properties, are well documented in [64]. However, Table 1 Main parameters and settings of the NREL 5-MW reference wind turbine. Parameters

Values

Rating Rotor configuration Control Rotor, Hub diameter Cut-in, Rated, Cut-out wind speed Cut-in, Rated rotor speed Rated tip speed Shaft tilt, Precone

5 MW Upwind, 3 Blades Variable speed & pitch control 126 m, 3 m 3 m/s, 11.4 m/s, 25 m/s 6.9 rpm, 12.1 rpm 80 m/s 0° , 0°

177

Composite Structures 213 (2019) 173–189

J. Chen et al.

Table 2 Material properties, reproduced from [66]. Material

Type (–)

Density (kg/m3)

E1 (MPa)

E2 = E3 (MPa)

G12 = G13 = G23 (MPa)

v12 = v13 = v23 (–)

Carbon (UD) Foam Eglass (UD) SNL (TX) Saertex (DB) Gelcoat

orthotropic orthotropic orthotropic orthotropic orthotropic isotropic

1220 200 1920 1850 1780 1235

114,500 256 41,800 27,700 13,600 3440

8390 256 14,000 13,650 13,300 –

5990 22 2630 7200 11,800 1380

0.27 0.3 0.28 0.39 0.49 0.3

3.2. Implementation of bend-twist coupling via spar cap fiber angle

Similar to Resor [66], transition section shapes are developed by interpolation between the root circle and airfoil at maximum chord. The shear webs begin at the span of 1.35 m and ends at the blade tip. An example of the calculated sectional 6 × 6 stiffness matrix at the span of 53 m is shown in Eq. (19). It consists of the diagonal “classical” stiffness and the non-diagonal “coupling” stiffness. The blade extension, edgewise bending, flapwise bending and torsion stiffnesses are denoted as c33, c44 , c55 and c66 , respectively.

− 0.08⎤ 0 0 0 ⎡3.74 0.20 18.38 0 0 0 1.48 ⎥ ⎢ 168.39 24.01 − 10.03 0 ⎥ ⎢ 7 [cij] = 10 × ⎢ 42.42 − 1.64 0 ⎥ ⎢ 5.66 0 ⎥ sym . ⎥ ⎢ 1.62 ⎦ ⎣

In the present work, bend-twist coupling is introduced into the blade by locating off-axis plies with constant fiber angle in the spar cap. Mirrored layup with respect to the mid-plane of the section are located on both suction side and pressure side of the blade. This concept of material-based coupled adaptive blade is illustrated in Fig. 8. By changing the fiber angle, both twist-to-feather and twist-to-stall BTC blades can be achieved. Many literatures [2,6–9,13,42–46] have proved that the twist-to-feather design can help reduce the fatigue loads of the blade in turbulence wind. Thus, positive fiber angle as shown in the figure is used in the present analysis. Moreover, there can be full-span (100%) or partial-span (less than 100%) coupled BTC blade designs, which means off-axis orientation of the fibers starts at the root or at a specified span location. The fully coupled BTC blade is discussed in the present study except where noted. A total of 7 BTC blades are developed with a series of fiber angles of 5° to 30° with a step of 5 degree, and a extreme case of 45° fiber angle. It should also be noted that the off-axis fibers in spar cap will cause decreased stiffness in the blade. Though the load carrying capability of the BTC blade might thus be reduced, but in the present work we do not restore the stiffness of the blade by adding new materials. Because, as pointed out by Bagherpour [44], restoring the stiffness would be a rather conservative approach considering that through BTC not only the stiffness but also the loading of the blade will decrease. A full assessment of the load carrying capability variations of the BTC blades requires detailed structural strength analysis, which is out of the scope of this paper. However, an assessment of the stiffness variations of the BTC blades is conducted. In fact, the stiffness variation due to fiber off-axis rotation will strongly influence the BTC blade design and evaluation, which therefore must be taken into account in reality. Fig. 9 shows the variations of the sectional stiffness with respect to fiber angle, including different span stations. The stiffness reduction of the BTC blade is calculated by comparing to the baseline blade which is

(19)

The computed sectional stiffness results, for both blades with ply numbers of Resor’s design and reduced ply numbers used in the present work, are shown in Fig. 7. The results are compared with the original stiffness of the NREL 5-MW blade released in Jonkman’s report. Overall, the inversely designed stiffness distributions, except the tension stiffness, are in relatively good agreement with the original design. Significant differences can be observed in the tension stiffness between the computed value and original design value. This is because the original tension stiffness value is not exact but estimated to be 107 times the average mass moment of inertia at each blade station, due to lack of reference blade information [64]. By reducing the ply numbers in composite laminates of spar cap, LE/TE panels and shear web, the bending stiffness of the inversely designed blade is in much better agreement with the original design value, while the torsion stiffness is less influenced by the adjusted laminate thicknesses. It should also be stated that the stiffness differences in the root region (less than 10 m span), though relatively obvious, do not significantly influence the blade structural dynamic behaviors because of high stiffness and small deflections therein.

TE reinf., Eglass TE panel, Foam LE panel, Foam

40 Thickness (mm)

40 Thickness (mm)

50

Root reinf., Triax Spar cap, Carbon Shear web, Foam LE/TE/Spar cap, Gelcoat LE/TE/Spar cap, Triax skin Shear web, Saertex skin

50

30 20 10

30

20

10

0

0 0

10

20

30

40

50

60

0

10

20

30 Radius (m)

Radius (m) Fig. 5. Layup of the NREL 5-MW composite blade. 178

40

50

60

Composite Structures 213 (2019) 173–189

J. Chen et al.

Fig. 6. Cross-sectional 2-D meshes used in VABS code to calculate the blade’s stiffness properties, generated by preVABS.

stemming from the off-axis fibers in spar cap. The edgewise bend-twist coupling factor can also be defined in a similar way. Fig. 10 shows the variation of the coupling factors at different span locations as a function of off-axis fiber angle. It is seen that a strong flapwise bend-twist coupling relation is established as the fiber in spar cap is orientated at nonzero angles. The maximum BTC is obtained at approximately 20° fiber angle. It is also observed that the coupling magnitude increases as the spanwise location moves to the blade tip. That is because the ratio of spar cap in the whole section area increases along with the span, as the spar box width remains constant while the blade chord decreases in the region from the maximum chord location (approximately at 12.5 m span) to the tip. The drop of coupling magnitude of the section at 60 m span is caused by the small thickness of carbon laminates in spar cap therein (see Fig. 5). The edgewise bendtwist coupling factors are also displayed as a comparison, which shows negligible coupling magnitude, indicating that the off-axis fiber in spar cap does not introduce edgewise bend-twist coupling to the blade. The BTC blades are numerically tested by applying a static uniform distributed load of 3000 N/m along the span length in flapwise direction. The steady deflections are calculated by the developed finite element beam code and presented in Fig. 11. As bending stiffness decreases with fiber angle, the tip flapwise deflection increases along with the fiber angle. Meanwhile, due to the bend-twist coupling effect, the blade twists towards feather under the thrust loading. The maximum twist deformation occurs at 20° − 25° fiber angle, which is in accordance with the coupling factor distribution.

developed with zero fiber angle in the previous section. It is observed that the tension and flapwise bending stiffness are significantly reduced due to fiber orientation rotation. The edgewise bending stiffness, though also slightly reduced, is little influenced by the fiber orientation. The torsional stiffness, on the other hand, is increased along with fiber angle in the region of 5° to approximately 30°, and then is slightly decreased when fiber angle increasingly gets larger. Overall, it can be expected that the BTC blade will undergo larger aeroelastic deflections in operation due to decreased stiffness compared to the uncoupled baseline blade. An example of the sectional stiffness matrix at the span of 53 m of the BTC blade with 10° fiber angle is given in Eq. . Compared with the stiffness matrix in Eq. for the baseline blade, it is clearly seen that the bend-twist coupling terms herein are nonzero.

0.01 − 0.07 ⎤ ⎡3.52 0.23 − 0.22 0.01 19.63 − 9.23 − 0.73 0.61 1.58 ⎥ ⎢ 153.48 23.21 − 9.04 − 0.74 ⎥ ⎢ 7 [cij] = 10 × ⎢ 41.59 − 1.57 − 0.06 ⎥ ⎢ sym . 5.06 0.59 ⎥ ⎢ ⎥ 1.74 ⎦ ⎣

(20)

The magnitude of flapwise bend-twist coupling of the BTC blade is assessed by introducing a factor expressed as

α=

c56 c55·c66

(21)

where c55 is the sectional flapwise bending stiffness, c66 is the sectional torsional stiffness, and c56 is the bend-twist coupling stiffness term 179

Composite Structures 213 (2019) 173–189

J. Chen et al.

27

2

9

9

12 9 6 3 0 0

5

4

2

Jonkman Present with Resor's ply number Present with reduced ply number

15 Tension stiffness (x10 Nm )

Flapwise bending stiffness (x10 Nm )

18

24 21

2 18 15

0

12

15

9

0 5

6

Flapwise bending stiffness 0.5

15

2

Edgewise bending stiffness (x10 Nm )

35

0.4

4

25

2

20

0

15

15

30

45

60

Jonkman Present with Resor's ply number Present with reduced ply number

10

12

0.3 0.2

9

2

Torsional stiffness (x10Nm )

9

30

10 15 20 25 30 35 40 45 50 55 60 Span (m)

Tension stiffness 8

60

3

Span (m)

40

45

Jonkman Present with Resor's ply number Present with reduced ply number

6

0

10 15 20 25 30 35 40 45 50 55 60

30

5

9 0.1 0.0

6

15 30 45 60 Jonkman Present with Resor's ply number Present with reduced ply number

3 0

0 0

5

0

10 15 20 25 30 35 40 45 50 55 60

5

10 15 20 25 30 35 40 45 50 55 60 Span (m)

Span (m)

Edgewise bending stiffness

Torsional stiffness

Fig. 7. Structural properties of the redesigned NREL 5-MW composite blade.

present work, the pre-twist processing is conducted at rated condition. Due to the aero-elastic coupling effect, the aerodynamic pre-twist for a BTC blade should be determined in an iterative way. In each iterative step, aeroelastic simulation of the wind turbine in uniformed inflow is performed to get the steady-state blade torsional deformation and rotor output power. The flexible blade is then modified of its aerodynamic twist based on the obtained rotor power. This iterative loop is continued until the steady-state output power matches with that of the reference rotor by a given tolerance. The reference rotor is equipped with rigid blades that have original aerodynamic twist distributions. Fig. 12 shows an example of redesigned aerodynamic twist for the BTC blade with 20° fiber angle, compared with the original twist distribution of the uncoupled rigid blade. The twist is set unchanged in the root sections (less than 20% span) where the cylinder and transition sections mainly provide supporting forces and generate no lift forces for power production. It is clearly shown that pre-twist toward stall is demanded for the BTC blade to compensate the torsional deformation. Fig. 13 compares the steady-state aeroelastic tip deflections of the BTC blades before and after pretwist processing. In both situations, the maximum blade tip twist due to BTC is achieved at 20° − 25° fiber angle, which is in consistency with the previous numerical static test result although the blades now are applied with more physically realistic distributed aerodynamic loads. It therefore also verifies that one

Fig. 8. Spar cap fiber orientation of the BTC blade.

3.3. Pretwist of the bend-twist coupled blade The BTC blade is designed to alleviate the unsteady load fluctuations. Due to the BTC effect, however, the blade would significantly twist toward feather under the mean aerodynamic loads in operations, leading to decreased effective angle of attacks and therefore decreased rotor power production. Clearly, the aerodynamic twist of the BTC blade must be readjusted in order to restore the power loss due to nonoptimum twist distribution of the deformed blade. More specifically, aerodynamic pre-twist toward stall is required for a BTC blade to achieve rotor power comparable with that of uncoupled blades. In the 180

Composite Structures 213 (2019) 173–189

Flapwise bending stiff ness reduction (-

J. Chen et al.

Tension stiffness reduction (-)

0.1 0.0 -0.1 -0.2 -0.3 -0.4

span = 20 m span = 30 m span = 40 m span = 50 m span = 60 m

-0.5 -0.6 -0.7

0

5

10

15

20

25

30

35

40

45

50

0.1 0.0 -0.1 -0.2 -0.3 -0.4 -0.5 -0.6 -0.7 -0.8 0

5

span span span span span

= = = = =

10

15

Fiber angle (deg)

Torsional stiffness reduction (-)

Edegewise bending stiffness reduction (-)

0.30

0.00 -0.02 -0.04 -0.06

-0.12

span = 20 m span = 30 m span = 40 m span = 50 m span = 60 m 0

5

10

15

20

25

30

20

25

30

35

40

45

50

(b) Flapwise bending stiffness

0.02

-0.10

m m m m m

Fiber angle (deg)

(a) Tension stiffness

-0.08

20 30 40 50 60

35

40

45

span span span span span

0.25 0.20

20 30 40 50 60

m m m m m

0.15 0.10 0.05 0.00 -0.05 0

50

= = = = =

5

10

15

20

25

30

35

40

45

= = = = =

m m m m m

50

Fiber angle (deg)

Fiber angle (deg)

(c) Edgewise bending stiffness

(d) Torsion stiffness

Fig. 9. Variations of sectional stiffness of BTC blades with different fiber angles.

0.003

Edge-Twist coupling factor

Flap-Twist coupling factor

0.30 0.25 0.20 0.15 0.10 span = 20 m span = 30 m span = 40 m span = 50 m span = 60 m

0.05 0.00 -0.05 -5

0

5

10

15

20

25

30

0.000 -0.003 -0.006 span span span span span

-0.009 -0.012

35

40

45

50

-5

0

5

10

15

20

25

30

35

40

20 30 40 50 60

45

50

Fiber anle (deg)

Fiber anle (deg)

(a) flapwise bend-twist coupling factor

(b) edgewise bend-twist coupling factor

Fig. 10. Variation of bend-twist coupling factors at different span locations with respect to spar cap fiber angle.

4. Influence of bend-twist coupling on rotor aerodynamic performance

can analyze the blade BTC characteristics by performing static bending tests in practice. As the pre-twist angle is introduced into the blade, the aerodynamic load becomes heavier and thus the blade undergoes larger deflections.

The NREL 5-MW wind turbine is variable speed and pitch controlled 181

Composite Structures 213 (2019) 173–189

J. Chen et al.

5

8

flapwise

torsion

Power (MW)

6 Tip deflection (m, deg)

0 deg 10 deg 25 deg Rigid

4

4 2

3 2 1 0

Tip flapwise deflection Tip torsion deflection

0

2

3

4

5

6 7 8 Wind speed (m/s)

9

10

11

12

(a) Power output

-2 0.6

0 deg

10 deg

25 deg

Rigid

0.5

0

5

10

15

20

25

30

35

40

45

50

0.4 Cp (-)

-4 -5

CAP fiber off-axis angle (deg) Fig. 11. Tip deflections of BTC blades with different fiber angles under distributed 3000 N/m in flapwise direction.

0.3 0.2 0.1 0.0

Aerodynamic twist (deg)

15

2

4

5

6 7 8 Wind speed (m/s)

9

10

11

12

(b) Power coefficient (Cp)

10

Fig. 14. Comparison of rotor power output between rigid reference blade and BTC blades with different spar cap fiber angles, at below-rated wind speeds.

5

over its full operating range. As the BTC blade is pretwisted at rated condition, it still remains a problem how it behaves at other wind speeds. In this section, the influences of adopting BTC blades on rotor aerodynamic performance across the operating range are assessed. The main objective is to prove that the BTC blade (with pretwist) would not significantly affect the rotor power production while it can be used to alleviate load fluctuations. The BTC blades with 0°, 10° and 25° fiber angles are selected for the study in this section. Firstly, aeroelastic simulations are performed for these BTC blades at steady below-rated wind speeds, with the same rotor speeds as the rigid reference wind turbine. The computed steadystate power outputs are presented in Fig. 14, along with the corresponding power coefficients. It shows that the rotor with flexible uncoupled blades of 0° fiber angle produces almost the same power as the rigid reference rotor, indicating that uncoupled blade produces nearly no twist deformation under below-rated wind speeds. As for the BTC blades with nonzero fiber angles, the rotor power is also not reduced and even a little bit higher than the reference one under these wind speeds, thanks to the pretwist processing. As we have discussed, the BTC blades are pretwisted based on the torsions at rated wind condition. Considering that

0 original 20 deg BTC blade with pre-twist

-5

0

10

20

30

40

50

60

Radius (m) Fig. 12. Aerodynamic pre-twist for the BTC blade with 20° fiber angle.

10 8 6 4 without pretwist, flapwise without pretwist, torsion with pretwist, flapwise with pretwist, torsion

2 0

25

-2

pitch setting (deg)

Tip deflection (m, deg)

3

-4 -6 -5

0

5

10

15

20

25

30

35

40

45

50

CAP fiber off-axis angle (deg) Fig. 13. Tip bending and torsion deflections of the BTC blades with and without pretwist, under uniformed wind at rated condition.

20

rigid 0 deg 10 deg 25 deg

15 10 5 0

15

20 wind speed (m/s)

25

Fig. 15. Comparison of pitch settings between rigid reference blade and BTC blades with different spar cap fiber angles, at above-rated wind speeds. 182

Composite Structures 213 (2019) 173–189

J. Chen et al.

1.5

0 deg 10 deg 25 deg

6

0 deg 10 deg 25 deg

1.2 Tip Ip deflection (m)

Tip Oop deflection (m)

8

4

2

0.9 0.6 0.3 0.0

0 2

4

6

8

2

10 12 14 16 18 20 22 24 26 Wind speed (m/s)

4

6

8

10 12 14 16 18 20 22 24 26 Wind speed (m/s)

(a) out-of-plane deflection

(b) in-plane deflection

Tip twist deflection (deg)

0 -1 -2 -3 0 deg 10 deg 25 deg

-4 -5

2

4

6

8

10 12 14 16 18 20 22 24 26 Wind speed (m/s)

(c) twist deflection Fig. 16. Tip deflections of BTC blades across the operating wind speeds.

settings are required for BTC blades to limit loads and keep the rotor power at the rated value. Fig. 16 presents the steady-state aeroelastic tip deflections of the three BTC blades across the full operating wind speed range varying from 3 m/s to 25 m/s. It is seen that the bending deflections in out-ofplane and in-plane directions increase with wind speed in below-rated region and achieve the largest at the rated wind speed. In the aboverated region, though the wind speed becomes increasingly larger, the bending deflections decrease with wind speed due to the active pitch control system. Obviously, the out-of-plane deflections are relatively larger than in-plane deflections across the operating range, because of larger bending stiffness in the later direction. As the fiber angle increases from 0 to 25 degree, the bending stiffness of BTC blade decreases and so the blade bending deflections become larger. In this situation, as mentioned previously, further study would be required to balance the reduction of stiffness and the bend-twist coupling magnitude, such as by adding new materials and/or by adjusting the ratio of off-axis layer thickness to non-off-axis one. The tip twist deflections of the blades over the varying wind speeds are also presented in Fig. 16(c). All the blades are found twisted toward feather. It is clearly observed that the uncoupled blade with 0° fiber angle differs from the other two BTC blades. For the uncoupled blade, very slight and almost constant torsion deformations are produced at below-rated wind speeds, while in above-rated region relatively significant torsion deformations are achieved and increase along with the increasing wind speed. This is caused by the resultant torsional

the torsions are smaller under below-rated wind speeds (later in Fig. 16), the blades are actually operating under smaller aerodynamic twists and so higher angle of attacks at these wind speeds compared with the rated situation. In this situation, the rotor speeds should be readjusted to capture the optimal power coefficients under below-rated wind speeds. However, it requires an aerodynamic-structural coupled optimization procedure which is out of the scope of the present study. In above-rated wind speed situations, the rotor speed is kept constant and the blade is pitched toward feather to limit the loads and keep the rotor power at the rated value. Due to the different torsion deformations under these wind speeds, the pitch angles for BTC blades with fiber angles are expected to differ from the original settings. A binary search program is developed to determine the corresponding pitch angles for those BTC blades with different fiber angles. The resulting pitch angles in cases of 15 m/s, 20 m/s, and 25 m/s wind speed are shown in Fig. 15, by comparing the blades of 0°, 10° and 25° fiber angles to the rigid reference one. It is seen that the uncoupled blades of 0° fiber angle needs a little bit lower pitch setting than the reference blade to achieve the rated power, which is because of the larger nose-down torsion deformations at above-rated speeds (later in Fig. 16). Similar to the below-rated speed situations, the BTC blades undergo smaller nose-down torsion deformations at above-rated wind speeds than the rated one. Considering that the BTC blades are pretwisted at rated condition, it would result in increased effective angle of attacks and above-rated rotor power if the original reference pitch settings are applied. Therefore, larger pitch

183

Composite Structures 213 (2019) 173–189

J. Chen et al.

9

6

Root out-of-plane moment (x10 Nm)

moments applied on the blade, which consist of the aerodynamic torsional moment and the aerodynamic lift and drag multiply by the offset distance of the aerodynamic center relative to the elastic axis of the NREL 5-MW blade, which is detailly discussed in our previous work [57]. As for the BTC blades with nonzero fiber angles, on the other hand, the torsion deformations are mainly determined by the blade bending deflections due to coupling effects. Significant torsion deformations are achieved and the largest twist occurs at rated wind, which is consistent with the bending deflections. 5. Load alleviation of wind turbine in shear flow The existence of the atmospheric shear layer leads to vertical wind shear, in which the wind speed at the surface of the ground is zero and increases with the height. As the rotor diameter increases, the variation of wind speed over the rotor height becomes increasingly considerable. In the present work, the wind shear profile is defined by the power law as follows:

y ⎞ V (y ) = Vhub ⎜⎛ ⎟ ⎝ yhub ⎠

5

0 deg 20 deg

4 60

120

180

240

300

360

Fig. 18. Comparison of root bending moments between the uncoupled blade and the BTC blade with 20° fiber angle, under shear flow at rated wind speed.

α

(22)

1.2

180

Tip twist variation (deg)

-0.8

0

60

120

180

240

300

360

Fig. 19. Comparison of tip twist cyclic variations between the uncoupled blade and the BTC blade with 20° fiber angle in shear wind.

reduced compared to the uncoupled blade, which is beneficial to increase the blade fatigue life. This result proves that the BTC blade with nonzero fiber angle can potentially be used to mitigate the fatigue loads of wind turbines in shear wind. To explain the load mitigation procedure, the cyclic tip twist variations for both the uncoupled blade and the BTC blade are presented in Fig. 19, with the average twist value subtracted for comparison convenience. It is seen that the twist deformation behavior of the BTC blade significantly differs from the uncoupled one. Under the aerodynamic torsion moments, the uncoupled blade is twisted toward stall at 0° azimuth and toward feather at 180° azimuth, which can even further increase the load imbalance due to wind shear. As for the BTC blade, however, the twist deformation variation is totally opposite as illustrated in Fig. 19. The BTC blade is twisted toward feather at 0° azimuth because the aerodynamic load therein is higher and so the bend-twist coupled effect leads to increased nose-down torsion deformation, and vice versa when the blade is at 180° azimuthal position. The twist deformation toward stall at 180° azimuth helps alleviate the reduction of effective AOA due to decreased inflow speed, while the twist toward feather at 0° azimuth attributes to mitigate the increase of AOA due to increased inflow speed. In this way, the AOA variation amplitude is reduced effectively, as shown in Fig. 20, and therefore alleviates the load variations.

hub height Vhub=8m/s Vhub=11.4m/s Vhub=20m/s

0 20 V(y), wind speed (m/s)

-0.4

Azimuth (deg)

120

10

0.0

toward feather

rotor plane

30

0 deg 20 deg

0.4

-1.2

150

60

toward stall

0.8

Firstly, the BTC blade with 20° fiber angle is studied as an example to prove the capability of BTC blade to mitigate the load imbalance due to wind shear. Aeroelastic simulations are performed for both the uncoupled and coupled blades, under rated wind speed at hub height. For the convenience of discussion, the blade azimuthal position is defined 0° when the blade at the highest position pointing upward and 180° when it rotates to the lowest positon pointing downward. The results of root bending moments are compared between the uncoupled blade and the BTC blade with 20° fiber angle, as shown in Fig. 18. Clearly, the blade undergoes cyclic varying loads in shear wind inflow. The bending moment reaches its largest value at around 0° azimuthal position and gradually decreases to the lowest at around 180° position. It is obviously observed that, by using BTC blade, the amplitude of root out-of-plane bending moment fluctuation is significantly

y, height (m)

6

Azimuth (deg)

5.1. Load variation using bend-twist coupling blade with 20° fiber angle

0

7

0

where Vhub is the wind speed at the hub height yhub , y is the vertical distance from the ground, and α is the wind shear exponent and is set to 0.4 in the present work. Fig. 17 shows the vertical sheared wind speed profiles at different hub wind speeds of 8 m/s, 11.4 m/s (rated) and 20 m/s. When the rotor is operating under the shear flow, each blade experiences varying inflow speed at different azimuthal positions within a revolution, resulting in cyclic fluctuations of aerodynamic loads. These cyclic varying loads during the long-term operations of a wind turbine would produce significant fatigue loads on the blades, as well as the drive system and tower. In this section, the designed BTC blades with nonzero fiber angles are utilized to mitigate the azimuthal load fluctuations due to wind shear.

90

8

30

Fig. 17. The shear wind velocity profiles at different hub wind speeds. 184

Composite Structures 213 (2019) 173–189

J. Chen et al.

there must be another influence factor that affects the load reduction. The blade vibrating velocity is found to account for the increasing reduction of root moment amplitude in the fiber angle region of 20°–45°. As discussed previously, the stiffness of BTC blade decreases with the fiber angle, leading to increasingly larger bending deflections and higher vibrating velocity of the blade in shear wind (see Fig. 22d). The blade vibrating velocity in downwind direction at 0 degrees azimuth helps to further decrease the effective angle of attack, while the vibrating velocity in upwind direction at 180 degrees azimuth helps to further increase the effective angle of attack. Therefore, though the variation amplitude of blade twist deformation due to BTC effects is indeed decreased in the fiber angel region of 20° to 45° (Fig. 22c), the variation amplitudes of AOA (Fig. 22b) and root moment (Fig. 22a) are reduced under the combined influences of blade vibrating velocity and twist deformation.

8

AOA (deg)

7 6 5 0 deg 20 deg

4 3

0

60

120

180

240

300

360 5.3. Influence of different bend-twist coupling regions

Azimuth (deg)

The BTC blade can be designed as either full-span or partial-span coupled. For partial-span coupled blade, the off-axis fiber placement starts from somewhere along the span and ends at the blade tip. Three blades with different coupling regions are discussed in this subsection, including 100%-, 75%- and 50%-span coupled blades, with the same spar cap off-axis fiber angle of 25 degree. Fig. 23 presents the steady-state twist deformations of the three blades in rated steady wind. It shows that the twist deformation mainly occurs at the bend-twist coupled region. As the coupling region becomes smaller, the torsion deformation also decreases. Due to the high stiffness in the root section, the bending and torsion deformations are relatively small in the root section. Therefore, the 75%-span coupled blade produces torsion deformation similar to the 100%-span coupled one. Aeroelastic simulations under the same wind shear are then performed for the three BTC blades. The amplitude reductions of the root bending moments are also obtained by comparison to the uncoupled blade. Fig. 24 presents the comparison of the reduction ratios between the three blades. Overall, the load mitigation effects for both out-ofplane and in-plane bending moments are weakened as the coupling region becomes smaller. However, the root out-of-plane bending moment amplitude reduction percentage is decreased by 19.6% if switching from the 75%-span coupled blade to the 50%-span coupled one, while the percentage is only decreased by 4.6% if switching from the 100%-span coupled blade to the 75%-span coupled one. This result is consistent with the torsion distributions in Fig. 23. If one expects to keep the root section stiffness unchanged to help resist loads, the 75%span coupled blade can be optional because it achieves a balance between the load mitigation capability and the root section strength.

Fig. 20. Comparison of effective angle of attack (AOA) variations between the uncoupled blade and the BTC blade with 20° fiber angle in shear wind.

5.2. Influence of different spar cap fiber angles

Reduction of root moment amplitude (%)

In this subsection, the influences of spar cap fiber angle on the load mitigation effects of BTC blades are investigated. Aeroelastic simulations are conducted for all the BTC blades with seven different fiber angles, under the shear flow with rated wind speed at hub height. The root bending moments of these blades are obtained from the simulations. Then the peak-to-peak amplitude reduction of root moment on each BTC blade is calculated by comparison to that of the uncoupled blade with 0° fiber angle. Fig. 21 presents the variation of root bending moment amplitude reduction with respect to the spar cap fiber angle. Both the out-of-plane and in-plane bending moments are included. The figure shows that as the fiber angle increases from 5° to 45°, the load mitigation magnitude increases though the increasing rate gradually decreases. The result is as expected in the fiber angle region of 5° to 20°, considering that the blade achieves the maximum BTC magnitude at approximately 20° fiber angle. The gradually enhanced load mitigation capability of the BTC blade with fiber angle in this region is mainly due to the increased BTC magnitude. However, though the BTC magnitude decreases in the region of 20° to 45°, the reduction percentage of root moment amplitude does not decrease but instead keeps increasing slightly. In this situation,

50

40

5.4. Cases in other wind speeds All the above works are conducted under the rated wind speed. It remains a problem how the BTC blade behaves under other wind speeds because of the varying aerodynamic loads. In this subsection, the different load mitigation effects of the BTC blade under below-rated and above-rated wind speeds are discussed. To illustrate the load distributions, aeroelastic simulations are firstly conducted for the uncoupled blade under the wind shear with different wind speeds at hub height. Two below-rated conditions (5 m/s, 8 m/s) and three above-rated conditions (15 m/s, 20 m/s and 25 m/s) are considered. Fig. 25 presents the cyclic fluctuations of the root bending moments applied on the blade in these shear wind conditions, by plotting the mean values and peak-to-peak amplitudes, respectively. The mean root out-of-plane bending moment reaches its largest value at approximately the rated wind speed and decreases at both below- and above-rated speeds. However, the fluctuation amplitude of the root moment gets

30

20

10

0

Out-of-plane In-plane 0

5

10

15

20

25

30

35

40

45

50

Spar cap fiber angle (deg) Fig. 21. Amplitude reduction of blade root bending moments in shear wind by using different BTC blades with varying spar cap fiber angles. 185

Composite Structures 213 (2019) 173–189

J. Chen et al.

8.5 8.0

8

7.5 7.0

7

AOA (deg)

6

Root out-of-plane moment (x10 Nm)

9

6 0 deg 10 deg 20 deg 30 deg 45 deg

5 4 0

60

6.5 6.0 5.5 5.0

0 deg 10 deg 20 deg 30 deg 45 deg

4.5 4.0 3.5

120

180

240

300

3.0

360

0

60

120

Azimuth (deg)

(a) Root out-of-plane bending moment 1.0

0.4 0.2 0.0 -0.2 -0.4 -0.6 60

120

180

300

360

3 Tip out-of-plane velocity (m/s)

Tip twist variation (deg)

0.6

0

240

(b) Angle of attack at r = 50 m

0 deg 10 deg 20 deg 30 deg 45 deg

0.8

-0.8

180

Azimuth (deg)

240

300

2 1 0

-2 -3

360

0 deg 10 deg 20 deg 30 deg 45 deg

-1

0

60

120

Azimuth (deg)

180

240

300

360

Azimuth (deg)

(c) Tip twist variation

(d) Tip out-of-plane vibrating velocity

Amplitude reduction of root moments (%)

Fig. 22. Variation of blade aerodynamic and structural responses along with azimuth angle under shear flow.

1

Torsion deformation (deg)

0 -1 -2 -3 100% span coupled 75% span coupled 50% span coupled

-4 -5

0

10

20

30

40

50

60

50

Root out-of-plane bending moment Root in-plane bending moment

40

4.6% 19.6%

30

20

10

0

100% span

75% span

50% span

BTC region

Radius (m) Fig. 24. Comparison of amplitude reduction ratios between the three BTC blades with different coupling regions under shear wind with rated speed at hub.

Fig. 23. Steady-state twist deformations of the three BTC blades with different coupling regions under rated steady wind.

186

Composite Structures 213 (2019) 173–189

peak-to-peak root moment (x10 Nm)

J. Chen et al.

6

out-of-plane in-plane

7 6 5 4 3 2 1 0

5 m/s

8 m/s 11.4 m/s 15 m/s

12 out-of-plane in-plane

6

Mean root moment (x10 Nm)

8

20 m/s

25 m/s

10 8 6 4 2 0

5 m/s

8 m/s 11.4 m/s 15 m/s

20 m/s

25 m/s

Wind speed (m/s)

Wind speed (m/s)

(a) Mean value

(b) Peak-to-peak amplitude

Amplitude reduction of root moments (%)

Fig. 25. Root moments on the uncoupled blade under shear wind with different hub wind speeds: (a) mean value; (b) peak-to-peak amplitude.

50

model is established to conduct the study. The unbalanced laminate induced BTC is represented by the cross-section stiffness matrices, which are computed by utilizing the tool VABS. The main conclusions are drawn as follows.

out-of-plane in-plane 40

1) Both the cross-sectional stiffness and the BTC magnitude of the composite adaptive blade are related with the spar cap off-axis fiber angle. The tension and flapwise bending stiffnesses are remarkably reduced as the fiber angle increases, while the edgewise bending stiffness is little influenced and the torsional stiffness are increased. Meanwhile, as the fiber angle increases, the coupling magnitude firstly increases and then decreases. The highest coupling magnitude is achieved at approximately 20° fiber angle for the designed adaptive blade. 2) With the pretwist processing applied at the rated condition, the BTC adaptive blade is shown capable to produce adequate aerodynamic rotor power. The output power of the BTC blades at below-rated wind speeds is even a little bit higher than the original rigid design, which indicates a further optimization potential of the BTC blade. It is also found that higher pitch angle settings at above-rated wind speeds are required for the BTC blade to maintain the constant rotor power production when compared with the original rigid design. 3) The composite material BTC blade is proved capable to effectively alleviate the cyclic fluctuating loads under shear wind. The reduction of the cyclic load variation amplitude is highly influenced by the spar cap fiber angle which induces BTC effects. As the fiber angle increases from 0° to 20°, the amplitude reduction of the out-of-plane and in-plane root bending moments are gradually increased up to 38% and 16%, respectively. As the fiber angle continues to increase, the load reduction increases slightly though the BTC magnitude is decreased, which is because of the increased vibration velocity due to reduction of blade bending flexibility. 4) The load mitigation effects are weakened for partial-span coupled adaptive blades. The load mitigation capability is less decreased for the 75%-span coupled blade compared with the 50%-span coupled one, due to the relatively high stiffness distribution in blade root region. 5) The blade root bending moments vary across over the full operating hub wind speed range. The mean values of both out-of-plane and inplane root moments reach the highest at approximately the rated wind speed, while the peak-to-peak amplitudes continue to increase along with the hub wind speed. The BTC blade is further proved capable to alleviate the loads over the full operating range and

30

20

10

0

5 m/s

8 m/s 11.4 m/s 15 m/s

20 m/s

25 m/s

wind speed Fig. 26. Amplitude reduction of root moments under shear wind with different wind speeds at hub, using the BTC blade with 25 deg fiber angle.

increasingly larger as the wind speed increases, which may result in higher fatigue loads on the blades. The load mitigation effects are then investigated by simulating the BTC blade under these hub wind speeds. The 100%-span coupled BTC blade with 25 degree fiber angle is selected to conduct the study. Fig. 26 shows the amplitude reduction results. It is seen that the BTC blade plays positive roles in all the wind speed cases. Both the out-of-plane and in-plane root moments are reduced of the fluctuation amplitudes. It is worth noting that the load mitigation effects are enhanced at high wind speed region, where the load variations are also higher. 6. Conclusion In the present work, the composite material bend-twist coupled blade is studied of its capability to alleviate the cyclic fluctuating loads in shear wind. The NREL 5-MW blade is inversely redesigned of its composite layup configuration, which is basically in consistence with the original bending and torsion stiffness distributions. The BTC adaptive blade is implemented by rotating spar cap fibers away from the blade axis. An aeroelastic model of the wind turbine blade based on nonlinear geometrically exact beam theory and free wake lifting surface 187

Composite Structures 213 (2019) 173–189

J. Chen et al.

shows larger amplitude reduction at higher wind speed. [27]

Acknowledgements

[28]

This work was supported by the National Natural Science Foundation of China [No. 11872249], China.

[29] [30]

References [31] [1] Kragh KA, Hansen MH. Load alleviation of wind turbines by yaw misalignment. Wind Energy 2014;17:971–82. https://doi.org/10.1002/we.1612. [2] Hayat K, Ha SK. Load mitigation of wind turbine blade by aeroelastic tailoring via unbalanced laminates composites. Compos Struct 2015;128:122–33. https://doi. org/10.1016/j.compstruct.2015.03.042. [3] Larwood S, Van Dam CP, Schow D. Design studies of swept wind turbine blades. Renew Energy 2014;71:563–71. https://doi.org/10.1016/j.renene.2014.05.050. [4] Liebst BS. Wind turbine gust load alleviation utilizing curved blades. J Propuls Power 1986;2:371–7. https://doi.org/10.2514/3.22897. [5] Ashwill T, Kanaby G, Jackson K, Zuteck M. Development of the Sweep-Twist Adaptive Rotor (STAR) Blade. 48th AIAA Aerosp Sci Meet Incl New Horizons Forum Aerosp Expo 2010. p. 1–14. doi:10.2514/6.2010-1582. [6] Lobitz DW, Veers PS. Aeroelastic behavior of twist-coupled HAWT blades. ASME/ AIAA Wind Energy Symp Reno, NV 1998:75–83. https://doi.org/10.2514/6. 1998-29. [7] Veers P, Laboratories SN, Bir G, Renewable N, Wind N, Lobitz D. Aeroelastic tailoring in wind-turbine blade applications. Wind Energy 1998. [8] Lobitz DW, Veers PS, Laino DJ. Performance of twist-coupled blades on variable speed rotors 1999. p. 11. doi:10.2514/6.2000-62. [9] Lobitz DW, Veers PS. Load mitigation with bending/twist-coupled blades on rotors using modern control strategies. Wind Energy 2003;6:105–17. https://doi.org/10. 1002/we.74. [10] Ong CH, Wang J, Tsai SW. Design, manufacture and testing of a bend-twist D-spar. AIAA Pap. No. 99-0025; 1999, p. 43–52. doi:10.2514/6.1999-25. [11] Berring P, Branner K, Berggreen C, Knudsen HW. Torsional performance of wind turbine blades – part I: Experimental investigation. 16th Int Conf Compos Mater July 8-13 2007; Kyto Japan 2007. p. 1–10. [12] Luczak M, Manzato S, Peeters B, Branner K, Berring P, Kahsin M. Dynamic investigation of twist-bend coupling in a wind turbine blade. J Theor Appl Mech 2011;49:765–89. [13] Vesel RW, McNamara JJ. Performance enhancement and load reduction of a 5 MW wind turbine blade. Renew Energy 2014;66:391–401. https://doi.org/10.1016/j. renene.2013.12.019. [14] Motley MR, Liu Z, Young YL. Utilizing fluid-structure interactions to improve energy efficiency of composite marine propellers in spatially varying wake. Compos Struct 2009;90:304–13. https://doi.org/10.1016/j.compstruct.2009.03.011. [15] Liu Z, Young YL. Utilization of bend-twist coupling for performance enhancement of composite marine propellers. J Fluids Struct 2009;25:1102–16. https://doi.org/10. 1016/j.jfluidstructs.2009.04.005. [16] Young YL, Motley MR, Barber R, Chae EJ, Garg N. Adaptive composite marine propulsors and turbines: progress and challenges. Appl Mech Rev 2016;68:62001. https://doi.org/10.1115/1.4034659. [17] Nicholls-Lee RF, Turnock SR, Boyd SW. Application of bend-twist coupled blades for horizontal axis tidal turbines. Renew Energy 2013;50:541–50. https://doi.org/10. 1016/j.renene.2012.06.043. [18] Motley MR, Barber RB. Passive control of marine hydrokinetic turbine blades. Compos Struct 2014;110:133–9. https://doi.org/10.1016/j.compstruct.2013.11. 026. [19] Barber RB, Motley MR. Cavitating response of passively controlled tidal turbines. J Fluids Struct 2016;66:462–75. https://doi.org/10.1016/j.jfluidstructs.2016.08. 006. [20] Young YL, Garg N, Brandner PA, Pearce BW, Butler D, Clarke D, et al. Load-dependent bend-twist coupling effects on the steady-state hydroelastic response of composite hydrofoils. Compos Struct J 2018;189:398–418. https://doi.org/10. 1016/j.compstruct.2017.09.112. [21] Barber RB, Hill CS, Babuska PF, Wiebe R, Aliseda A, Motley MR. Flume-scale testing of an adaptive pitch marine hydrokinetic turbine. Compos Struct 2017;168:465–73. https://doi.org/10.1016/j.compstruct.2017.02.051. [22] Murray RE, Ordonez-Sanchez S, Porter KE, Doman DA, Pegg MJ, Johnstone CM. Towing tank testing of passively adaptive composite tidal turbine blades and comparison to design tool. Renew Energy 2018;116:202–14. https://doi.org/10. 1016/j.renene.2017.09.062. [23] Hansen MOL, Sørensen JN, Voutsinas S, Sørensen N, Madsen HA. State of the art in wind turbine aerodynamics and aeroelasticity. Prog Aerosp Sci 2006;42:285–330. https://doi.org/10.1016/j.paerosci.2006.10.002. [24] Bazilevs Y, Hsu M-C, Akkerman I, Wright S, Takizawa K, Henicke B, et al. 3D simulation of wind turbine rotors at full scale. Part I: Geometry modeling and aerodynamics. Int J Numer Meth Fluids 2010. https://doi.org/10.1002/fld. [25] Bazilevs Y, Hsu MC, Kiendl J, Wuchner R, Bletzinger KU. 3D simulation of wind turbine rotors at full scale. Part II: Fluid – structure interaction modeling with composite blades. Int J Numer Meth Fluids 2011;65:236–53. https://doi.org/10. 1002/fld. [26] Carrion M, Steijl R, Woodgate M, Barakos GN, Munduate X, Gomez-Iradi S. Aeroelastic analysis of wind turbines using a tightly coupled CFD-CSD method. J

[32]

[33] [34]

[35]

[36]

[37] [38] [39]

[40]

[41]

[42]

[43]

[44]

[45]

[46]

[47]

[48] [49]

[50]

[51]

[52] [53] [54]

[55]

[56]

[57]

188

Fluids Struct 2014;50:392–415. https://doi.org/10.1016/j.jfluidstructs.2014.06. 029. Yu DO, Kwon OJ. Predicting wind turbine blade loads and aeroelastic response using a coupled CFD-CSD method. Renew Energy 2014;70:184–96. https://doi.org/ 10.1016/j.renene.2014.03.033. Hodges DH. Nonlinear composite beam theory. Reston: American Institute of Aeronautics and Astronautics, Inc.; 2006. Bauchau OA. Flexible multibody dynamics. Dordrecht: Springer; 2013. Hodges DH. A mixed variational formulation based on exact intrinsic equations for dynamics of moving beams. Int J Solids Struct 1990;26:1253–73. https://doi.org/ 10.1016/0020-7683(90)90060-9. Hodges DH. Geometrically exact, intrinsic theory for dynamics of curved and twisted anisotropic beams. AIAA J 2003;41:1131–7. https://doi.org/10.2514/2. 2054. Yu W, Hodges DH, Volovoi V, Cesnik CES. On Timoshenko-like modeling of initially curved and twisted composite beams. Int J Solids Struct 2002;39:5101–21. https:// doi.org/10.1016/S0020-7683(02)00399-2. Blasques JP, Bitsche RD. BECAS – an Open-Source Cross Section Analysis Tool Motivation & Overview; 2012. Blasques JP, Stolpe M. Multi-material topology optimization of laminated composite beam cross sections. Compos Struct 2012;94:3278–89. https://doi.org/10. 1016/j.compstruct.2012.05.002. Blasques JP, Bitsche RD, Fedorov V, Lazarov BS. Accuracy of an efficient framework for structural analysis of wind turbine blades. Wind Energy 2016;19:1603–21. https://doi.org/10.1002/we. Blasques JP. Multi-material topology optimization of laminated composite beams with eigenfrequency constraints. Compos Struct 2014;111:45–55. https://doi.org/ 10.1016/j.compstruct.2013.12.021. Bauchau OA, Hong CH. large displacement analysis of naturally curved and twisted composite beams. AIAA J 1987;25:1469–75. Bauchau OA, Hong CH. Nonlinear composite beam theory. J Appl Mech 1988;55:156–63. Kim T, Hansen AM, Branner K. Development of an anisotropic beam finite element for composite wind turbine blades in multibody system. Renew Energy 2013;59:172–83. https://doi.org/10.1016/j.renene.2013.03.033. Otero AD, Ponta FL. Structural analysis of wind-turbine blades by a generalized timoshenko beam model. J Sol Energy Eng 2010;132:11015. https://doi.org/10. 1115/1.4000596. Fleming I, Luscher DJ. A model for the structural dynamic response of the CX-100 wind turbine blade. Wind Energy 2014;17:877–900. https://doi.org/10.1002/we. 1603. Bottasso CL, Campagnolo F, Croce A, Tibaldi C. Optimization-based study of bendtwist coupled rotor blades for passive and integrated passive/active load alleviation. Wind Energy 2013;16:1149–66. https://doi.org/10.1002/we. Bottasso CL, Campagnolo F, Croce A. Multi-disciplinary constrained optimization of wind turbines. Multibody Syst Dyn 2012;27:21–53. https://doi.org/10.1007/ s11044-011-9271-x. Bagherpour T, Li XM, Manolas DI, Riziotis VA. Modeling of material bend-twist coupling on wind turbine blades. Compos Struct 2018;193:237–46. https://doi.org/ 10.1016/j.compstruct.2018.03.071. Scott S, Capuzzi M, Langston D, Bossanyi E, McCann G, Weaver PM, et al. Effects of aeroelastic tailoring on performance characteristics of wind turbine systems. Renew Energy 2017;114:887–903. https://doi.org/10.1016/j.renene.2017.06.048. Şener Ö, Farsadi T, Gözcü MO, Kayran A. Evaluation of the effect of spar cap fiber angle of bending-torsion coupled blades on the aero-structural performance of wind turbines. J Sol Energy Eng 2018. https://doi.org/10.1115/1.4032391. Hauptmann S, Bülk M, Schön L, Erbslöh S, Boorsma K, Grasso F, et al. Comparison of the lifting-line free vortex wake method and the blade-element-momentum theory regarding the simulated loads of multi-MW wind turbines. J Phys Conf Ser 2014;555:12050. https://doi.org/10.1088/1742-6596/555/1/012050. Ahlström A. Influence of wind turbine flexibility on loads and power production. Wind Energy 2006;9:237–49. https://doi.org/10.1002/we.167. Shen X, Chen J, Zhu X, Liu P, Du Z. Multi-objective optimization of wind turbine blades using lifting surface method. Energy 2015;90:1111–21. https://doi.org/10. 1016/j.energy.2015.06.062. Jeong MS, Kim SW, Lee I, Yoo SJ, Park KC. The impact of yaw error on aeroelastic characteristics of a horizontal axis wind turbine blade. Renew Energy 2013;60:256–68. https://doi.org/10.1016/j.renene.2013.05.014. Saverin J, Peukert J, Marten D, Pechlivanoglou G, Paschereit CO, Greenblatt D. Aeroelastic simulation of multi-MW wind turbines using a free vortex model coupled to a geometrically exact beam model. J Phys Conf Ser 2016;753:1–10. https:// doi.org/10.1088/1742-6596/753/8/082015. Chattot J-J. Helicoidal vortex model for steady and unsteady flows. Comput Fluids 2006;35:733–41. Chattot JJ. Tower shadow modelization with helicoidal vortex method. Comput Fluids 2008;37:499–504. https://doi.org/10.1016/j.compfluid.2007.07.006. Sebastian T, Lackner MA. Development of a free vortex wake method code for offshore floating wind turbines. Renew Energy 2012;46:269–75. https://doi.org/ 10.1016/j.renene.2012.03.033. Gupta S, Leishman J. Comparison of Momentum and Vortex Methods for the Aerodynamic Analysis of Wind Turbines. 43rd AIAA Aerosp. Sci. Meet. Exhib., Reno, Nevada: AIAA; 2005. p. 1–24. doi:10.2514/6.2005-594. Chen J, Shen X, Liu P, Zhu X, Du Z. Design Tool for Aeroelastic Analysis of Wind Turbine Blades Based on Geometrically Exact Beam Theory and Lifting Surface Method. 35th Wind Energy Symp., 2017, p. 1–18. doi:10.2514/6.2017-0450. Chen J, Shen X, Zhu X, Du Z. Impact of Blade Flexibility on Wind Turbine Loads and

Composite Structures 213 (2019) 173–189

J. Chen et al.

[58] [59]

[60] [61] [62] [63]

analysis. J Am Helicopter Soc 1993;38:3–14. https://doi.org/10.4050/JAHS.38.3. [64] Jonkman J, Butterfield S, Musial W, Scott G. Definition of a 5-MW reference wind turbine for offshore system development. National Renewable Energy Laboratory; 2009. doi:10.1002/ajmg.10175. [65] Griffith DT, Ashwill TD. The Sandia 100-meter All-glass Baseline Wind Turbine Blade : SNL100-00. Albuquerque, New Mexico: 2011. doi: Report No. SAND20113779. [66] Resor BR. Definition of a 5MW/61.5m Wind Turbine Blade Reference Model. Sandia National Laboratories; 2013. [67] Nijssen R. WMC5MW laminate lay-out of reference blade for WP 3. Project UpWind Report; 2007. [68] Tian S, Liu X, Yu W. PreVABS Manual 2018. Available from: https://cdmhub.org/ resources/1597.

Pitch Settings. Proc. ASME 2018 Turbomach. Tech. Conf. Expo., Oslo, Norway: ASME; 2018, p. 1–15. Shen X, Yang H, Chen J, Zhu X, Du Z. Aerodynamic shape optimization of nonstraight small wind turbine blades. Energy Convers Manage 2016;119:266–78. Danielson DA, Hodges DH. Nonlinear beam kinematics by decomposition of the rotation tensor. J Appl Mech 1987;54:258–62. https://doi.org/10.1115/1. 3173004. Reissner E. On one-dimensional large-displacement finite-strain beam theory. Stud Appl Math 1973(L):87–95. Simo JC, Vu-Quoc L. A three-dimensional finite-strain rod model. Part II Computational aspects. Comput Methods Appl Mech Eng 1986;58:79–116. Borri M, Merlini T. A large displacement formulation for anisotropic beam analysis. Meccanica 1986;21:30–7. https://doi.org/10.1007/BF01556314. Bauchau OA, Kang NK. A multibody formulation for helicopter structural dynamic

189