66 American Academic Scientific Research Journal for Engineering, Technology, and Sciences ISSN (Print) 2313-4410, ISSN (Online) 2313-4402 http://asrjetsjournal.org/ Rotor Blade Aerodynamics Forces Modelling for Dynamic Analysis Mohamed Abdalla Almheriegh * Associate Professor, Department of civil Engineering, Faculty of Engineering - Tripoli University, P O Box 82677, Tripoli, Libya Email: malmherigh@gmail.com Abstract In order to establish a rational method for structural analysis of wind turbine rotor with respect to failure in ultimate loading, either in strength or fatigue fashion. Two famous methods in predicting the complicated nature of the two main lift and drag forces acting ‘normal and tangential’ forces on the rotor blade are discussed and numerical solution ‘iteration’ using the coefficients involved in establishing these forces is presented and, therefore are readily introduced in the analysis stage typically done by dynamics analysis codes. Keywords: ultimate loads on rotor blades; lift and drag forces; normal and tangential forces on blades; structural analysis of wind turbine blades. 1. Introduction Various methods are used to calculate the aerodynamic forces acting on the blades of a wind turbine. The most advanced are numerical methods solving Navier-Stokes equations for the compressible flow as well as the flow near the blades. The two major approaches to calculating the forces are the Actuator Disk Model and the Blade Element Model. In the paper to follow a brief introduction to these methods will be presented which extracted from Freris [1,2,3,4,5,6]. ------------------------------------------------------------------------ * Corresponding author. http://asrjetsjournal.org/ American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 67 2. Actuator disc model Figure 1: Actuator disk Based on Bernoulli’s equation and energy balances [5,6,4]. It is assuming that, the rotor is replaced by an actuator disc, through which the static pressure decreases discontinuously. By examining the flow through a control volume the extractable power from the turbine can be calculated Figure (1). The stream tube has a cross-sectional area larger than the cross-sectional area for the upstream disc and a smaller area than the downstream disc. Within the stream tube, continuity is required and the rate of the mass flow must be constant.   m  U  =   = www  1 By introducing an axial interference factor, a, as the fractional decrease in wind velocity between the free stream and the rotor plane represented by    v a 2 It is found that )1( a  3 The air, which passes through the disk, undergoes an overall change in velocity. The velocity multiplied by the flow rate gives the rate of change of momentum, more known as a force )( wmT L     4 Combining the equations above with the fact that the change of momentum comes entirely from the pressure American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 68 difference across the actuator disc, it is obtained that     )1()( aAwpp   5 To obtain the pressure difference, the Bernoulli’s equation is applied separately to the upstream and downstream sections of the stream tube. For the upstream section it becomes     pp 22 2 1 2 1  6 Similarly, downstream pUpU w      22 2 1 2 1 7 Subtracting equation (6) from equation (7) yields )( 2 1 22 UUpp w      8 Equations (8) and (5) give UU a w   )21( 9 Substituting (9), (3), and (1) into (4) obtains the force, T, which gives )1(2 2 aaT UA    10 Combining (3), and (9) and the rate of work done by the force; UTP   the power extraction from the air is obtained as: 23 )1(2 aaP UA    11 Or, by introducing the dimensionless power coefficient, 2)1(4 aaCp  pCP UA 3 2 1     12 The power coefficient represents the efficiency of the turbine, which depends on variables like the wind speed, American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 69 the rotor speed and the pitch angle. The coefficient shows how much of the kinetic energy in the air stream that is transformed into mechanical energy. The maximum pC as a function of a is 59.0 27 16 pC , at 3 1a (obtained by taking the first derivative of the power coefficient 2)1(4 aaCp  with respect to ‘a’ and equating it to zero). The maximum value of pC = 0.59 is called the Betz limit and applies to all types of wind turbines. Intuitively there must be a wind speed change at which the conversion efficiency is maximum. If there were no change in wind speed, no energy would be extracted, and the power of the wind turbine would be zero. If the air were brought completely to rest, all its energy would dissipate. However, a rotating wind turbine will not completely prevent the flow of air, so it can only extract a proportion of the kinetic energy in the wind. Hence in terms of force exerted by wind on rotor area based on this assumption pressure is underestimated on blades. Therefore, its use for predicting these forces is in decline. Modern wind turbines operate at performance coefficient of about (0.4) 3. Blade Element Theory For the use of aeroelastic codes in design calculations, the aerodynamic method has to be very time efficient. The Blade Element Momentum (BEM) theory has been shown to give good accuracy with respect to time cost. In this method, the turbine blades are divided into a number of independent Figure 2: Load forces on blade section elements along the length of the blade. At each section, a force balance is applied involving 2-D section lift and drag with the thrust and torque produced by the section. At the same time, a balance of axial and angular American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 70 momentum is applied. This produces a set of non-linear equations, which can be solved numerically for each blade section, the presented discussion follows Andres [5] and Det Norske [4], only the force in the flow direction was regarded. The BEM theory, also takes notice of the tangential force due to the torque in the shaft. The left force L per unit length is perpendicular to the relative speed V rel of the wind and equals: CV Lrel c L 2 2   13 Where c is the blade chord length. The drag force D per unit length, which is parallel to relV is given by CV Drel c D 2 2   14 Since the interest only in the forces, normal to and tangential to the rotor-plane, the lift and drag are projected on these directions, Figure (2)  sincos DLFN  15 And  cossin DLFT  16 The application of this theory requires information about the lift and drag air foil coefficients LC and DC . Those coefficients are generally given as functions of the angle of incidence, Figure (3)   17 Further, it is seen that ra a U   )'1( )1( tan     18 In practice, the coefficients are obtained from a 2D wind tunnel tests. If  exceeds about 015 , the blade will stall. This means that the boundary layer on the upper surface becomes turbulent, which will result in a radical increase of drag and a decrease of lift. The lift and drag coefficients need to be projected onto the normal and tangential directions.  sincos DLN CCC  19 And American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 71  cossin DLT CCC  20 Figure 3: Velocities at rotor plane Rotor Blade Chosen RISΦ-1 airfoil Further, a solidify  is defined as the fraction of the annular area in the control volume, which is covered by the blades r Bc r r   2 )(  21 Where B denotes the number of blades. The normal force and the torque on the control volume of thickness dr since TN andFF are forces per length drcC a NdrNFdT NN U   2 22 sin )1( 2 1    22 And rdrcC ara NdrrNFdQ LT U    cossin )'1()1( 2 1    23 Finally, the two influence factors are declared by American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 72 1 sin4 1 2   NC a   24 And 1 cossin4 1 '   TC a   25 These two factors are the key to establish a value for the forces normal and tangential to rotor plane using blade element method ‘BEM’. These factors are partially empirical and solution could be attained via iteration for each value of r/R. To increase accuracy tip loss correction factor need be applied, this is to allow for the velocities and forces not being circumferentially uniform due to the rotor having a finite number of blades. This factor is expressed as: )) sin2 (arccos(exp 2  r rRB F   26 This reduction factor is called Prandtl’s tip loss factor Det Norske [4]. This is yields: 1 sin4 1 2   NC F a   27         1 cossin4 1 ' Cr F a   28 All terms as defined before. More practical implementation of this method is further detailed in the literature, while results of the iteration done using the MathCAD software to arrive at linear pressure profile for typical one case blade parameters is shown next to this paragraph: American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 73 4. MathCAD sheet for blade wind load typical iteration Input parameters Iteration No 1  2.31rad sec 1  c 0.604 m r 10 m  10.6 deg Re 1.6 10 6  CN 1.125 CT 0.464 Tangential Force Coefficient Pitch Angle Using alpha, calculate Lift and Drag coefficients: Normal Force Coefficient  atan 1 a( ) V0 1 a1   r         28.419deg    lpha1 3 deg lpha2 5 deg lpha3 10 deg lpha4 15 deg lpha5 20 deg lpha6 25 deg CL1 0.6 CL2 0.94 CL3 1.21 CL4 1.25 CL5 1.18 CL6 0.98 CD1 0 CD2 0 CD3 0.03 CD4 0.085 CD5 0.16 CD6 0.32 CL linterp lpha CL   CL 1.211 CD linterp lpha CD   CD 0.127 CN CL cos   CD sin   CT CL sin   CD cos   r c B 2  r  r 0.029 F 2  acos e B R r 2 r sin          F 0.999 a 1 4 F sin  2 r CN 1  a1 1 4 F sin   cos   r CT 1  a 0.035 a1 8.069 10 3   17.819deg American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 74 Iteration No 2 Pitch Angle Using alpha, calculate Lift and Drag coefficients Normal Force Coeffiecient Tangential Force Coefficient  atan 1 a( ) V0 1 a1   r         27.394deg     16.794deg lpha1 3 deg lpha2 5 deg lpha3 10 deg lpha4 15 deg lpha5 20 deg lpha6 25 deg CL1 0.55 CL2 0.9 CL3 1.2 CL4 1.25 CL5 1.2 CL6 1.0 CD1 0 CD2 0 CD3 0.03 CD4 0.085 CD5 0.18 CD6 0.31 CL linterp lpha CL   CL 1.232 CD linterp lpha CD   CD 0.119 CN CL cos   CD sin   CN 1.149 CT CL sin   CD cos   CT 0.461 r c B 2  r  r 0.029 F 2  acos e B R r 2 r sin          F 0.999 a 1 4 F sin  2 r CN 1  a1 1 4 F sin   cos   r CT 1  a 0.038 a1 8.213 10 3  American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 75 Iteration No 3 Pitch Angle Using alpha, calculate Lift and Drag coefficients: Normal Force Coeffiecient Tangential Force Coefficient  atan 1 a( ) V0 1 a1   r         27.316deg     16.716deg lpha1 3 deg lpha2 5 deg lpha3 10 deg lpha4 15 deg lpha5 20 deg lpha6 25 deg CL1 0.55 CL2 0.9 CL3 1.2 CL4 1.25 CL5 1.2 CL6 1.0 CD1 0 CD2 0 CD3 0.03 CD4 0.085 CD5 0.18 CD6 0.31 CL linterp lpha CL   CL 1.233 CD linterp lpha CD   CD 0.118 CN CL cos   CD sin   CN 1.149 CT CL sin   CD cos   CT 0.461 r c B 2  r  r 0.029 F 2  acos e B R r 2 r sin          F 0.999 a 1 4 F sin  2 r CN 1  a1 1 4 F sin   cos   r CT 1  a 0.038 a1 8.231 10 3  American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 76 FINAL FORCES Resolve Forces K 4 F sin  2 r CN  ac 0.2 a if a 0.2 a 0.5 2 K 1 2 ac  K 1 2 ac  2  2 4 K ac 2  1             1.025 kg m 3  r 10 m FN 0.5  V0 2 1 a( ) 2  sin  2  c CN  FN 244.45N m 1  FT 0.5  V0 1 a( )  r 1 a1  sin   cos    c CT FT 98.102N m 1   10 deg Fres FN 2 FT 2  Fres 263.4N m 1  FORCEN FN cos   FT sin   FORCEN 257.771N m 1  FORCET FN sin    FT cos   FORCET 54.163N m 1  Fres FORCEN 2 FORCET 2  Fres 263.4N m 1  Forces0 FORCEN Forces1 FORCET Forces 0 0 1 257.771 54.163 N m 1  American Academic Scientific Research Journal for Engineering, Technology, and Sciences (ASRJETS) (2021) Volume 82, No 1, pp 66-77 77 5. Conclusions Theoretical investigation based on blade element method theory (BEM) to derive left and drag forces thus calculating wind pressure acting on rotor blade of wind turbine, the mathematical formula developed is solved via iteration using developed Mathcad sheet and forces need be applied to rotor during finite element analysis are readily calculated. References [1] Freris L. L “Wind Energy Conversion Systems” Printice Hall 1996. [2] Walker F. J and Jenkins N “Wind Energy Technology” John Wiley & Sons UK 1997. [3] DET NORSKE VERITAS (DNV), Offshore Standard, DNV-OS-E301, Position Mooring, 2004, http://www.dnv.com, October 2004. [4] Det Norske Veritas “Guidelines for Design of Wind Turbines” second edition 2002, printed by Judsk Centraltrykkeri, Denmark ISBN 87-550-2870-5 [5] Anders Ahlström “Simulating Dynamic Behaviour of Wind Power Structures” Licentiate Thesis Royal Institute of Technology, Department of Mechanics Stockholm2002. http://www2.mech.kth.se/~andersa/thesis/Licthesis.pdf, January 2003. [6] David M. Eggleston and Forrest S. Stoddard “Wind Turbine Engineering Design” 1987 Van Nostrand Inc. New York, ISBN 0-442-22195-9. http://www.dnv.com/ http://www2.mech.kth.se/~andersa/thesis/Licthesis.pdf