A new method of determination of the angle of attack on rotating wind turbine blades

: The angle of attack (AoA) is the key parameter when extracting the aerodynamic polar from the rotating blade sections of a wind turbine. However, the determination of AoA is not straightforward using computational ﬂuid dynamics (CFD) or measurement. Since the incoming streamlines are bent because of the complex inductions of the rotor, discrepancies exist between various existing determination methods, especially in the tip region. In the present study, ﬂow characteristics in the region near wind turbine blades are analyzed in detail using CFD results of ﬂows past the NREL UAE Phase VI rotor. It is found that the local ﬂow determining AOA changes rapidly in the vicinity of the blade. Based on this ﬁnding, the concepts of e ﬀ ective AoA as well as nominal AoA are introduced, leading to a new method of AOA determination. The new method has 5 steps: (1) Find the distributed vortices on the blade surface; (2) select two monitoring points per cross-section close to the aerodynamic center on both pressure and suction sides with an equal distance from the rotor plane; (3) subtract the blade self-induction from the velocity at each monitoring point; (4) average the velocity of the two monitoring points obtained in Step 3; (5) determine the AoA using the velocity obtained in Step 4. Since the monitoring points for the ﬁrst time can be set very close to the aerodynamic center, leading to an excellent estimation of AoA. The aerodynamic polar extracted through determination of the e ﬀ ective AoA exhibits a consistent regularity for both the mid-board and tip sections, which has never been obtained by the existing determination methods.


Introduction
Wind energy has become an important and growing source of electric power.CO 2 emission reduction targets and rapidly falling costs per kilowatt-hour have made wind energy the priority option for new power generating capacity in many places of the world.The country with the largest new and total wind power capacities is China (23 GW and 211 GW in 2018, respectively [1]).The region with the highest proportion of wind power in electricity consumption is Europe.In 2018, wind energy covered an estimated 14% of the EU's annual electricity consumption and a much higher share of 40.8% was achieved in its Member State Denmark [2].The utilization of wind energy is expected to grow steadily in the future [1], though as a natural energy it still has some unfavorable factors such as unstable output and its related demand for control [3,4].
The great achievements of wind power were driven by the continuous improvement of wind turbine technologies, among which the aerodynamics [5][6][7]  essential for the design of advanced wind turbines with high power efficiency and low structural loads.There are mainly three general methods in wind turbine aerodynamic computation: the blade-element momentum (BEM) theory [8][9][10], the vortex wake theory (VWT) [11,12], and the computational fluid dynamics (CFD) [13,14].CFD is the most general one with a minimum of hypotheses, giving the most complete flow information in the whole simulation field.On the other hand, it is the most computationally expensive method.It is meaningful to relate the results of CFD to BEM theory or VWT and validate some key issues such as interference factors [15], tip loss [16,17], rotational effect [18], etc.Many of these investigations highly depend on the determination of the angle of attack (AoA).For example, researchers have to know the AoA in order to compare the lift and drag of a blade section in CFD with those in the BEM theory.However, the determination of AoA in CFD is not straightforward.
Several methods have been introduced for determining AoA using CFD or experimental data.These can be divided into three categories: inverse BEM or VWT, averaging technique (AT), and self-induction subtraction (SIS).
Inverse BEM was developed and used by Snel et al. [19], Bruining et al. [20], Laino et al. [21], Lindenburg [22], and Bak et al. [23], applying the equations of BEM theory to solve the AoA with given forces on blade sections.Inverse VWT was developed by Sant et al. [24,25] and Bretton et al. [26], in which the bound circulation of a free-wake vortex model is computed with section forces by the Kutta-Joukowski law, and AoA is consequently determined from the converged wake structure.An advantage of inverse BEM/VWT is its simplicity: it is easy to convert the normal BEM/VWT code into an inverse one.Nevertheless, the accuracy is restricted to the reliability of BEM/VWT method itself.
The AT method determines AoA by averaging the velocity at multiple locations.Hansen et al. [27] and Johansen et al. [28] performed an annular average along certain rings upstream and downstream of the rotor plane.This is usually called the azimuthal average technical (AAT).The average of Jost et al. [29] is along a symmetric, closed line (usually a circle) around the rotor blade.Rahimi et al. [30] used three monitoring points along the chord length on each side of the rotor.The above two ways of averaging are called the line AT and 3-point AT, respectively.The main function of various types of AT is the same: to cancel out the blade self-induction and estimate the required velocity in rotor plane.Inspired by the idea of cancelling out, Herráez et al. [31] employed only one monitoring point in the rotor plane for each cross-section.The point is located on the bisectrix of the angle between two arbitrary blades where the influence of the self-induction of each blade is cancelled out by the other blades.
Instead of cancelling out blade self-induction, SIS developed by Shen et al. [32,33] directly subtracts self-induction from the velocity of a monitoring point in front of the blade.The Biot-Savart law is used to compute self-induction from the bound circulation.The bound circulation is simplified to a concentrated vortex in [32], which only requires section forces as an input, but leads to a shortcoming that the monitoring point must be kept at a certain distance from the leading edge of the blade.As an improvement, the bound circulation is represented by distributed vortices in [33], which requires the pressure distribution over blade sections as an input.In this situation, the monitoring point can be chosen closer to the blade and a distance between 0.5c to 2.0c (c is the chord length) is recommended.For the sake of simplicity, the technique using a concentrated vortex [32] is hereinafter referred to as SIS1 and that using distributed vortices [33] is referred to as SIS2.
The above methods have played significant roles in the research of wind turbine aerodynamics.The accuracy of BEM method was improved by using the so-called three-dimensional (3-D) lift and drag coefficients extracted through the determined AoA.Such kinds of work were made by Yang et al. [34,35] using SIS1, Schneider et al. [36] using inverse BEM as well as AAT, and Syed Ahmed Kabir et al. [37] using inverse BEM.There exists several other applications with AoA determination, such as improving the actuator line/Navier-Stokes (AL/NS) simulation (Wimshurst et al. [38]), investigating stall delay as well as dynamic stall (Zhu et al. [18]), discussing unsteady phenomena under yawed conditions (Elgammi et al. [39] and Wen et al. [40]), and analyzing the measurements of a wind turbine in the field (Wu et al. [41]).As a key parameter, accurate determination of AoA is of superior importance.However, discrepancies between the methods were observed by many researchers.In the IEA Task 29 project [42], several methods were attempted to determine the AoA and the corresponding aerodynamic polar of blade sections, including inverse BEM, inverse VWT, and SIS1, but they resulted in remarkable difference (greater than 25%) in lift coefficient for some cases.Guntur et al. [43] observed a discrepancy between inverse BEM and AAT even in the mid-board sections, and attributed it to the effect of trailing vortices.Schneider et al. [36] found that the AoA determined by AAT produced a lower tangential force than the CFD values.Rahimi et al. [30] implemented four different methods of inverse BEM, AAT, 3-point AT, and SIS1.In their results, AAT and SIS1 lead to about 15% lower axial interference factor.Jost et al. [29] analyzed several AT and SIS methods, pointing out that the different considerations of trailing vortices lead to the discrepancies in the tip region.A recent comprehensive evaluation of different methods, including inverse BEM, AAT, line AT, 3-point AT, SIS1, and SIS2, was made by Rahimi et al. [44] on two different 10 MW reference wind turbines.The results are consistent with each other in the mid-board area but discrepancies arise and increase as approaching the tip.In addition, the chosen locations of monitoring points of AT methods were found to be a significant influence on the resulting induced velocity.The above investigations indicate that a deeper understanding of the determination of AoA is still necessary, especially for the tip region.
The present study presents a new method of AoA determination using the CFD velocity field.It combines the advantages of SIS2 and AT methods and determines the effective AoA by five steps: (1) Find the distributed vortices on the blade surface; (2) select two monitoring points per cross-section close to the aerodynamic center on both pressure and suction sides with an equal distance; (3) subtract the blade self-induction from the velocity at the monitoring points; (4) average the velocity obtained in Step 3; (5) determine the AoA using the velocity obtained in Step 4. In addition to the new method itself, the innovation of this study lies in the following items: (a) Based on an analysis of the induction variation in the vicinity of a rotor blade, the definitions of effective AoA and nominal AoA are proposed; (b) effective AoA is found to be different from nominal AoA in the tip region, which to some extent explains the observed discrepancies between the existing methods; (c) the lift and drag coefficients extracted through the effective AoA are largely consistent with those of the 2-D airfoil for both the mid-board and tip sections, which has never been obtained by the existing AT or SIS methods.
The CFD results used in the present study are from the simulations of the NREL UAE Phase VI rotor [45] solving the Reynolds averaged Navier-Stokes (RANS) equations using the k-ω SST turbulence model.The detailed computational settings and validations can be found in [13].
The paper is organized as follows.In Section 2, an analysis of AoA is made and the concepts of effective AoA and nominal AoA are introduced.In Section 3, the method of computing self-induction is derived.In Section 4, the velocity field after self-induction subtraction is analyzed.In Section 5, the effective AoA is determined by averaging the velocity on two monitoring points adjacent to the aerodynamic center.In Section 6, a comparison is performed between the determined effective AoA and nominal AoA, and the aerodynamic coefficients extracted through them are compared.In Section 7, conclusions of the present study are given.

Definitions of Effective AoA and Nominal AoA
The classic definition of AoA for a two-dimensional (2-D) airfoil is the angle between the chord line of the airfoil and the undisturbed uniform inflow.The direction of the undisturbed inflow can be obtained by setting a monitoring point far away from the airfoil where the self-induction of the airfoil is negligible.If we have to define the AoA using a monitoring point close to the airfoil, the self-induction at the monitoring point must be subtracted from the measured velocity.In the present paper, self-induction is defined as the velocity induced by the bound circulation of an airfoil/wing/blade.
For a fixed blade with a finite span, an additional induction generated by the tip and root vortices (as well as the trailing vortices before rolling up into the tip and root vortices) arises, which is called tip-root-induction in the present paper.Accordingly, there are two kinds of AoA definition: the geometric AoA and the effective AoA.The geometric AoA is the angle between the local chord and the undisturbed inflow, while the effective AoA is defined by the inflow with tip-root-induction included.The difference between the geometric AoA and the effective AoA of a blade section is an angle determined by the downwash because of the tip-root-induction.It is clear that the effective AoA varies along the wingspan, since the strength of downwash changes at different wing sections.
The inflow experienced by a blade section of a rotor is more complex than described above since the tip and root vortices become helical and there are three kinds of induction.The first induction, which is the azimuthally uniform part of the tip-root induction, corresponds to the induced velocity in which the rotor is simplified to a theoretical disc with an infinite number of blades.It is hereinafter referred to as disc-induction.The second is the blade self-induction.The third is the azimuthally non-uniform part of the tip-root induction, which leads to tip/root loss and is hereinafter referred to as 3D-induction because of its 3-D nature.The disc-induction leads to a continuous increase of the axial interference factor a as well as the tangential inference factor a in the axial direction, causing bended streamlines shown in Figure 1a.In lifting line theory, the most reasonable definition of AoA in this situation should employ the velocity at the aerodynamic center (point C in Figure 1) of the blade section.In practice, however, it is impossible to detect any velocity at the aerodynamic center since it is inside the blade (that is, outside the flow field).As an alternative choice, a monitoring point (point P in Figure 1) in the rotor plane can be set in front of the blade (along the azimuthal direction).The velocity of point P could be regarded to be equal to that of point C (supposing there is a theoretical permeable flow), if the self-induction illustrated in Figure 1b and the downwash (due to 3D-induction) illustrated in Figure 1c were excluded, according to the azimuthally uniform assumption of BEM theory.Table 1.List of the angle of attack (AoA) determination methods and their corresponding kinds of AoA.

Representing Airfoil by Distributed Vortices
According to the Kutta-Joukowski law, the lift of an airfoil in a potential flow can be related to a bound circulation of the airfoil by where L is lift, ρ is air density, V ∞ is inflow speed, and Γ is the circulation of a concentrated vortex which typically is located at the quarter chord point of the airfoil.The concentrated vortex of an airfoil at a positive AoA is illustrated in Figure 2. The circulation is defined as positive when the induced velocity of the vortex is clockwise in the present paper, meaning that the illustrated Let us think about the effects of the three inductions to the aerodynamic center.The disc-induction and 3D-induction make realistic contributions to the (theoretical) velocity of the aerodynamic center, while the self-induction makes no such a contribution because the aerodynamic center is exactly the center of the bound vortex.Therefore, excluding the self-induction is always necessary for any monitoring point away from the aerodynamic center, but excluding the 3D-induction leads to a different physical meaning of the determined AoA.In the present paper, after the self-induction is subtracted, the AoA defined by the velocity in a certain front of a blade section (where the 3D-induction is negligible) is called the nominal AoA, while the AoA defined at the aerodynamic center is called the effective AoA.
The existing methods of AoA determination lead to AoA with different physical meanings.The inverse BEM (with a tip / root loss correction) or inverse VWT method determines the velocity at the blade axis (which usually is very close to, even coincident with, the aerodynamic center), leading to an estimation of effective AoA.The AAT method estimates the annularly averaged velocity in the rotor plane, which significantly weakens the effect of 3D-induction and thus results in an AoA close to the nominal AoA.The line AT and 3-point AT methods attempt to estimate the velocity of the aerodynamic center but their monitoring points have to be kept within a certain distance (usually about 1 chord length) from the aerodynamic center, leading to a blend of effective AoA and nominal AoA.The SIS1 and SIS2 methods employ a monitoring point in front of the blade where the 3D-induction is small and thus results in a nominal AoA.
The methods as well as their corresponding kinds of determined AoA are listed in Table 1.It is obvious that different methods should be employed for different purposes, depending on which kind of AoA is required.For the blade sections of mid-board, the effect of 3D-induction is usually small and thus the nominal AoA and the effective AoA tend to be consistent.
Table 1.List of the angle of attack (AoA) determination methods and their corresponding kinds of AoA.

Representing Airfoil by Distributed Vortices
According to the Kutta-Joukowski law, the lift of an airfoil in a potential flow can be related to a bound circulation of the airfoil by where L is lift, ρ is air density, V ∞ is inflow speed, and Γ is the circulation of a concentrated vortex which typically is located at the quarter chord point of the airfoil.The concentrated vortex of an airfoil at a positive AoA is illustrated in Figure 2. The circulation is defined as positive when the induced velocity of the vortex is clockwise in the present paper, meaning that the illustrated circulation is positive and leads to a positive lift according to Equation (1).A shortcoming of representing the airfoil by a concentrated vortex is that the induced velocity of the vortex, e.g.

→
v at a point P(x,y) in Figure 2, is different from that of the airfoil entity in the near field, because of the ignorance of the airfoil's geometry.Table 1.List of the angle of attack (AoA) determination methods and their corresponding kinds of AoA.

Representing Airfoil by Distributed Vortices
According to the Kutta-Joukowski law, the lift of an airfoil in a potential flow can be related to a bound circulation of the airfoil by where L is lift, ρ is air density, V ∞ is inflow speed, and Γ is the circulation of a concentrated vortex which typically is located at the quarter chord point of the airfoil.The concentrated vortex of an airfoil at a positive AoA is illustrated in Figure 2. The circulation is defined as positive when the induced velocity of the vortex is clockwise in the present paper, meaning that the illustrated circulation is positive and leads to a positive lift according to Equation (1).A shortcoming of representing the airfoil by a concentrated vortex is that the induced velocity of the vortex, e.g.v′  at a point P(x,y) in Figure 2, is different from that of the airfoil entity in the near field, because of the ignorance of the airfoil's geometry.Another choice is to properly distribute vortices over the airfoil, as illustrated in Figure 3a.The vortex distribution needs to meet two requirements.First, the sum of the circulations of the distributed vortices should be equal to the circulation of the concentrated vortex, where γ k is the circulation of the kth vortex, and n is the total number of the distributed vortices.Second, at any point P(x,y), the velocity induced by the distributed vortices should be consistent with that induced by the airfoil entity, Energies 2019, 12, x; doi: FOR PEER REVIEW www.mdpi.com/journal/energies where ds  is a microelement along the clockwise and closed circuit A-B-C-D-A, and v  is the corresponding local velocity.The segment BC is located on the edge of the boundary-layer.
According to the boundary-layer theory [46], the velocity along DA (on the wall) is zero, and the velocity on AB or CD is almost perpendicular to the segment itself.As a result, the integrations along the segments AB, CD, and DA are zero, and the equation becomes where k v  is the averaged velocity on BC, and k l  is a vector pointing from B to C with Although the velocity at any location in the flow field can be obtained from CFD simulations, The above two conditions will be met only when γ k is determined by the actual velocity distribution over the airfoil.
In order to determine γ k , the airfoil surface is divided into n small segments, each corresponding to one of the distributed vortices.Figure 3b illustrates the viscous boundary-layer flow on the kth small segment.The definition of circulation leads to → s is a microelement along the clockwise and closed circuit A-B-C-D-A, and → v is the corresponding local velocity.The segment BC is located on the edge of the boundary-layer.According to the boundary-layer theory [46], the velocity along DA (on the wall) is zero, and the velocity on AB or CD is almost perpendicular to the segment itself.As a result, the integrations along the segments AB, CD, and DA are zero, and the equation becomes where → v k is the averaged velocity on BC, and → l k is a vector pointing from B to C with Although the velocity at any location in the flow field can be obtained from CFD simulations, | → v k | is difficult to be directly determined because the thickness of boundary-layer varies along the airfoil surface.Instead of calculating the thickness of boundary layer, Shen et al. [33] suggested determining v k from the local pressure by applying the Bernoulli equation.For an airfoil with fully attached flow, the viscous effect is substantially limited within the boundary-layer that is usually very thin.Therefore, the Bernoulli equation can be regarded as valid outside the boundary-layer, giving the relation between the velocity and pressure, where p k is the pressure at the location of the kth vortex, which can be taken from the airfoil surface since pressure is almost unchanged across the boundary-layer.
The circulation of the kth vortex, γ k , can now be determined by Equations ( 5) and (7).The filaments of the distributed vortices can be regarded infinitely long in the direction perpendicular to the plane of the airfoil, since the concept of airfoil is 2-D.The induced velocity at any point P(x,y) can then be computed by the following formula according to the Biot-Savart law for infinitely long vortex filaments, where → r k is the vector pointing from the kth vortex to point P(x, y), and → e is the unit vector in the direction of the vortex (perpendicular to the airfoil plane and points inside in Figure 2.
The velocity field around the airfoil can now be reproduced by the distributed vortices, through a simple summation of In order to evaluate the method introduced above, 4, where → v air f oil is the induced velocity of the airfoil entity, in which → v (x, y) is the velocity at point P(x,y), obtained from CFD in the present paper.As shown in Figure 4, the two contour maps look largely similar to each other.There is no difference observed except in a slender region of the viscous wake originated from the trailing edge.The vortices give a result with no viscous wake because of the potential assumption.

Representing Blade by Distributed Vortices
As in the airfoil case, a blade can also be represented by bound vortices distributed over its surface.Part of a rotor blade is illustrated in Figure 6 where the surface is divided into small faces.A relative error can further be defined as which is contoured in Figure 5.It is observed that: (1) The relative error is generally small in the flow domain, and can be considered negligible (δ < 0.5%) beyond a certain distance from the airfoil; (2) the error is negligible in the region close to the front half of the upper surface of the airfoil; (3) the error increases as approaching the trailing edge, but is almost less than 5% except in the viscous wake.

Representing Blade by Distributed Vortices
As in the airfoil case, a blade can also be represented by bound vortices distributed over its surface.Part of a rotor blade is illustrated in Figure 6 where the surface is divided into small faces.Supposing that the pressure distributions of the two neighboring blade sections (section 1 and Section 2 in Figure 6) are known, the bound vortices between the two sections can be determined.
Let us consider a small face ABCD in Figure 6, which is indexed as the kth face in the next derivation.Points A and B are located at the equal dimensionless locations on their respective sections, as well as points C and D. The vortex of this face can be simplified to a line vortex starting from point P1(x1,y1,z1) and ending at point P2(x2,y2,z2) if the lengths of AD and BC are small enough.Points P1 and P2 are the midpoints of AD and BC, respectively.
In a rotational reference frame fixed with the rotor, supposing the pressure values of points P1 and P2 are known, the Bernoulli equation leads to ( ) ( ) resulting in In general, the results shown in Figures 4 and 5 indicate a good representation of the airfoil entity by the distributed vortices.

Representing Blade by Distributed Vortices
As in the airfoil case, a blade can also be represented by bound vortices distributed over its surface.Part of a rotor blade is illustrated in Figure 6 where the surface is divided into small faces.Supposing that the pressure distributions of the two neighboring blade sections (Section 1 and Section 2 in Figure 6) are known, the bound vortices between the two sections can be determined.

Subtraction of Blade Self-Induction
One basic assumption of BEM theory is that the flow remains uniform in the azimuthal direction.However, this assumption cannot be directly validated by the velocity field simulated by CFD, because of the existence of the blade self-induction.Taking a mid-board section of r/R = 47% as an example, Figure 7 shows the contours of axial interference factor defined by using a CFD velocity field before the self-induction is subtracted.It is clearly shown that a strong variation exists around the blade, which is obviously far from the azimuthally uniform assumption of BEM theory and can only be caused by the blade self-induction.In order to make the CFD and BEM results comparable, the blade self-induction must be subtracted.
The method of representing the blade entity by distributed vortices provides a way for subtracting the blade self-induction.The velocity after the subtraction, denoted by v   , can be computed by: The axial interference factor defined by the velocity v   is shown in Figure 8.The thick dashed Let us consider a small face ABCD in Figure 6, which is indexed as the kth face in the next derivation.Points A and B are located at the equal dimensionless locations on their respective sections, as well as points C and D. The vortex of this face can be simplified to a line vortex starting from point P1(x 1 ,y 1 ,z 1 ) and ending at point P2(x 2 ,y 2 ,z 2 ) if the lengths of AD and BC are small enough.Points P1 and P2 are the midpoints of AD and BC, respectively.
In a rotational reference frame fixed with the rotor, supposing the pressure values of points P1 and P2 are known, the Bernoulli equation leads to Energies 2019, 12, 4012 where → v 1 and → v 2 are the velocity of points P1 and P2, respectively; p 1 and p 2 are the pressure of the points P1 and P2, respectively; r S1 and r S2 are the radial location of the Sections 1 and 2, respectively.The circulation of the line vortex can then be estimated by the following average: where → AD and → BC are the path vectors ( → AD = AD, → BC = BC, the vector directions are consistent with the clockwise direction in Figure 6).At point P(x,y,z) in the flow domain, the velocity induced by the kth line vortex can then be determined by the Biot-Savart law, k are unit vectors in the x, y, z directions.The velocity induced by all the distributed vortices can be computed by a summation, where B is the number of rotor blades, and n is the number of small faces on each blade.
→ ν vortices (x, y, z) is exactly the blade self-induction.It is clear that the result will be more accurate if the faces are smaller.

Subtraction of Blade Self-Induction
One basic assumption of BEM theory is that the flow remains uniform in the azimuthal direction.However, this assumption cannot be directly validated by the velocity field simulated by CFD, because of the existence of the blade self-induction.Taking a mid-board section of r/R = 47% as an example, Figure 7 shows the contours of axial interference factor defined by using a CFD velocity field before the self-induction is subtracted.It is clearly shown that a strong variation exists around the blade, which is obviously far from the azimuthally uniform assumption of BEM theory and can only be caused by the blade self-induction.In order to make the CFD and BEM results comparable, the blade self-induction must be subtracted.

Determination of Effective AoA
Figure 10 shows the contours of axial interference factor in the flattened section of r/R = 95%, defined by the velocity after the self-induction subtraction.The 3D-induction dominates the flow field of this tip section, which is very different from the result of r/R = 47% (Figure 8). Figure 11 shows the local perspectives of the two sections.In the r/R = 47% section, the influence of viscous trailing wake is limited to an area after the trailing edge.In the r/R = 95% section, the 3D-induction The method of representing the blade entity by distributed vortices provides a way for subtracting the blade self-induction.The velocity after the subtraction, denoted by → v , can be computed by: The axial interference factor defined by the velocity → v is shown in Figure 8.The thick dashed line represents the rotor plane with many fewer contours than that in Figure 7, implying much less variation in the azimuthal direction.

Determination of Effective AoA
Figure 10 shows the contours of axial interference factor in the flattened section of r/R = 95%, The azimuthal (ψ) variation of axial interference factor in the rotor plane is plotted in Figure 9, with the results before and after the self-induction subtraction.In front of the blade (left of the blade axis), the curve of "After subtraction" is approximately horizontal in a wide range of azimuthal angle, while the curve of "Before subtraction" presents a significant decline.Obviously, the result of "After subtraction," rather than that of "Before subtraction," is comparable with the BEM theory.As for the sharp fluctuation behind the blade (right of the blade axis), it is a result of the viscous effect of trailing wake.

Determination of Effective AoA
Figure 10 shows the contours of axial interference factor in the flattened section of r/R = 95%, defined by the velocity after the self-induction subtraction.The 3D-induction dominates the flow field of this tip section, which is very different from the result of r/R = 47% (Figure 8). Figure 11 shows the local perspectives of the two sections.In the r/R = 47% section, the influence of viscous trailing wake is limited to an area after the trailing edge.In the r/R = 95% section, the 3D-induction

Determination of Effective AoA
Figure 10 shows the contours of axial interference factor in the flattened section of r/R = 95%, defined by the velocity after the self-induction subtraction.The 3D-induction dominates the flow field of this tip section, which is very different from the result of r/R = 47% (Figure 8). Figure 11 shows the local perspectives of the two sections.In the r/R = 47% section, the influence of viscous trailing wake is limited to an area after the trailing edge.In the r/R = 95% section, the 3D-induction covers the whole blade section including the aerodynamic center, implying a significant influence of the effective AoA.
Energies 2019, 12, x; doi: FOR PEER REVIEW www.mdpi.com/journal/energies the effective AoA needs to be determined.The contours of AoA determined by the local velocity (after the self-induction subtraction) are shown in Figure 13.The monitoring points of the 3-point AT and the monitoring circle of the Line AT are also illustrated.It is clear that the determination would be more accurate if the monitoring points (circle) were closer to the airfoil.However, the two techniques have to be set at a certain distance between the monitoring points (circle) and the blade in order to avoid the error introduced by the blade self-induction, because no self-induction subtraction was performed.Once the blade self-induction is subtracted, much closer monitoring points can be set.
In the present determination, two monitoring points (A and A' in Figure 13) are set adjacent to the airfoil surface, being above and below the aerodynamic center and with a same distance to the rotor plane, respectively.In the present study, the aerodynamic center is assumed to be at the 25% chord and the specific locations of the points A and A' are defined as: their azimuthal coordinates are identical and are equal to that of the 25% chord, and their axial distances from the rotor plane are also identical and are determined by a minimum range covering the inside area and the boundary-layer (which can be easily identified by the velocity contours) of the blade.The two locations are very close to the aerodynamic center and thus can be expected as the best choice for determining the effective AoA.The final determined effective AoA is where a v  and t v  are the axial and tangential velocity after the self-induction subtraction.
To summarize, the new method of determination of AOA has five steps: (1) Find the distributed vortices on the blade surface (Section 3); (2) select two monitoring points per cross-section close to the aerodynamic center on both pressure and suction sides with an equal distance from the rotor plane; (3) subtract the blade induction from the velocity at each monitoring point (Section 4); (4) average the velocity of the two monitoring points obtained in Step 3; (5) determine the AoA using the velocity obtained in Step 4 (Equation ( 20)).The 3D-induction causes an azimuthally non-uniform flow in the rotor plane, leading to a dependence of the determined AoA on the locations of monitoring points.Figure 12 illustrates the definition of nominal AoA and effective AoA.If the AoA is determined at a monitoring point that is a certain distance in front of the blade, it should be regarded as nominal AoA.The closer the monitoring point is to the aerodynamic center the closer the determined AoA is to the effective AoA.However, it is never possible to get an exact effective AoA since no velocity can be detected if the monitoring point is placed inside the blade.A technique has to be introduced to estimate the virtual velocity at the aerodynamic center, if the effective AoA needs to be determined.The contours of AoA determined by the local velocity (after the self-induction subtraction) are shown in Figure 13.The monitoring points of the 3-point AT and the monitoring circle of the Line AT are also illustrated.It is clear that the determination would be more accurate if the monitoring points (circle) were closer to the airfoil.However, the two techniques have to be set at a certain distance between the monitoring points (circle) and the blade in order to avoid the error introduced by the blade self-induction, because no self-induction subtraction was performed.Once the blade self-induction is subtracted, much closer monitoring points can be set.In the present determination, two monitoring points (A and A' in Figure 13) are set adjacent to the airfoil surface, being above and below the aerodynamic center and with a same distance to the rotor plane, respectively.In the present study, the aerodynamic center is assumed to be at the 25% chord and the specific locations of the points A and A' are defined as: their azimuthal coordinates are identical and are equal to that of the 25% chord, and their axial distances from the rotor plane are also identical and are determined by a minimum range covering the inside area and the boundary-layer (which can be easily identified by the velocity contours) of the blade.The two locations are very close to the aerodynamic center and thus can be expected as the best choice for determining the effective AoA.The final determined effective AoA is where v a and v t are the axial and tangential velocity after the self-induction subtraction.
To summarize, the new method of determination of AOA has five steps: (1) Find the distributed vortices on the blade surface (Section 3); (2) select two monitoring points per cross-section close to the aerodynamic center on both pressure and suction sides with an equal distance from the rotor plane; (3) subtract the blade induction from the velocity at each monitoring point (Section 4); (4) average the velocity of the two monitoring points obtained in Step 3; (5) determine the AoA using the velocity obtained in Step 4 (Equation ( 20)).

Comparison Between Nominal AoA and Effective AoA
A comparison between the nominal AoA (α n ) and the effective AoA (α e ) is shown in Figure 14.The nominal AoA is determined by setting a monitoring point in the 0.5c front of the blade in the velocity field after the self-induction subtraction.From the mid-board to the tip, the nominal AoA goes down first and then goes up, while the effective AoA presents a downward trend and drops faster in the tip region.This result answers the question raised from the existing evaluations [43,44]; why different AoA determination methods lead to different trends in the tip region.The reason is that different determinations may lead to different definitions of AoA (nominal, effective, or blended).

Comparison Between Nominal AoA and Effective AoA
A comparison between the nominal AoA ( n α ) and the effective AoA ( e α ) is shown in Figure 14.The nominal AoA is determined by setting a monitoring point in the 0.5c front of the blade in the velocity field after the self-induction subtraction.From the mid-board to the tip, the nominal AoA goes down first and then goes up, while the effective AoA presents a downward trend and drops faster in the tip region.This result answers the question raised from the existing evaluations [43,44]; why different AoA determination methods lead to different trends in the tip region.The reason is that different determinations may lead to different definitions of AoA (nominal, effective, or blended).
The difference between the nominal AoA and the effective AoA is the induced angle because of downwash which is called downwash angle in the present paper,   The difference between the nominal AoA and the effective AoA is the induced angle because of downwash which is called downwash angle in the present paper, The trend of α i along r/R is shown in Figure 15a.The curves rise faster as r/R approaches the tip.In addition, a higher wind speed leads to a larger value of α i at a certain location of r/R.The ratio of α i to α n , which is shown in Figure 15b, looks interesting: all the curves generally coincide with each other.That implies a rule exists in the relationship between α i and α n .The value of α i might be estimated from α n that can be more easily determined, if a model of their relationship could be established in the future work.

Validation of Extracted Lift and Drag Coefficients
As AoA is determined, the lift and drag coefficients of blade elements can be extracted through the following relation, in which φ = α + θ is the flow angle, C n and C t are the normal and tangential force coefficients determined by where F n and F t are the normal and tangential forces of the blade section, and V rel is the relative velocity.Figure 16a,b shows the lift coefficients extracted through α n and α e , respectively.The forces used to determine the coefficients are obtained from pressure integration over the airfoil/blade sections.The 2-D experimental data [47] and 2-D CFD data are for the airfoil S809 at Reynolds number of 1 × 10 6 .Remarkable differences between C l (α n ) and C l (α e ) are seen.For the sections of r/R = 63% and r/R = 85%, the curves of C l (α n ) is slightly lower than the 2-D experimental/CFD data, while the curves of C l (α e ) are in good agreement with each other as well as the 2-D experimental/CFD data.For the sections closer to the tip (r/R = 95% and r/R = 99%), the curves of C l (α n ) are observed to decline significantly, being far away from the 2-D data.In contrast, the curves of C l (α e ) do not show a sustained decline.The curves of r/R = 95% and r/R = 99% are consistent with each other.They are only slightly lower than the 2-D curves, which could be caused by some phenomena, e.g., the de-cambering effect [48], which are not taken into account in the present study.Figure 17a,b shows the drag coefficients corresponding to α n and α e , respectively.The curves of C d (α n ) are observed much higher than the 2-D experimental/CFD data.In contrast, the curves of C d (α e ) match much better with each other and the 2-D experimental/CFD data, except the curve of r/R = 99%.The exception implies some physical difference in the mechanism of drag between the section of extreme tip and the 2-D airfoil.

Conclusions
A new method has been developed and presented for determining the AoA of rotating wind turbine blades.For AT methods, there must be a sufficient distance between the monitoring points and the blade in order to avoid the errors caused by the blade self-induction.The present method completely overcomes this shortcoming by representing the blade entity with distributed vortices and subtracting their induction from the velocity at monitoring points.The monitoring points are for the first time set very close to the aerodynamic center of the blade, leading to an excellent estimation of the effective AoA.The success of this method is based on a combination of the advantages of the SIS2 and AT methods.The method may be regarded as an extension of the SIS2, and thus can be called SIS3.The following conclusions are drawn from the present study:

•
The determination of AoA depends on the consideration of different inductions (disc-induction, blade self-induction, and 3D-induction) experienced by the rotor blade.The blade self-induction should always be excluded, disc-induction leads to the nominal AoA, and the sum of disc-induction and 3D-induction, i.e., the tip-root-induction, leads to the effective AoA.The discrepancies between existing methods observed by other researchers are to some extent caused by the unclassified comparisons.The above results have definitely shown that the lift and drag coefficients extracted through the effective AoA are largely consistent with those of the 2-D airfoil for both the mid-board and tip sections.This regularity of consistence has never been obtained by existing AT or SIS methods.

Conclusions
A new method has been developed and presented for determining the AoA of rotating wind turbine blades.For AT methods, there must be a sufficient distance between the monitoring points and the blade in order to avoid the errors caused by the blade self-induction.The present method completely overcomes this shortcoming by representing the blade entity with distributed vortices and subtracting their induction from the velocity at monitoring points.The monitoring points are for the first time set very close to the aerodynamic center of the blade, leading to an excellent estimation of the effective AoA.The success of this method is based on a combination of the advantages of the SIS2 and

Figure 1 .
Figure 1.Different kinds of inductions in the flattened section of a certain radial location: (a) Bent streamlines because of disc-induction; (b) self-induction; (c) downwash because of 3D-induction.

Figure 1 .
Figure 1.Different kinds of inductions in the flattened section of a certain radial location: (a) Bent streamlines because of disc-induction; (b) self-induction; (c) downwash because of 3D-induction.

Energies 2019 ,Figure 1 .
Figure 1.Different kinds of inductions in the flattened section of a certain radial location: (a) Bent streamlines because of disc-induction; (b) self-induction; (c) downwash because of 3D-induction.

Figure 2 .
Figure 2. Illustration of the concentrated vortex of an airfoil.

Figure 3 .
Figure 3. Illustration of the distributed vortices on an airfoil: (a) vortices; (b) viscous flow on a small segment.

Figure 3 .
Figure 3. Illustration of the distributed vortices on an airfoil: (a) vortices; (b) viscous flow on a small segment.

21 Figure 6 .
Figure 6.Illustration of the vortex at a small face of blade surface: the faces are intentionally enlarged in order to be displayed clearly.

Figure 6 .
Figure 6.Illustration of the vortex at a small face of blade surface: the faces are intentionally enlarged in order to be displayed clearly.

Energies 2019 ,Figure 11 .
Figure 11.Contours of axial interference factor after the self-induction subtraction around the blade the dashed black line is a typical contour line ( 7m/s V ∞ = , 72rpm ω =).

Figure 12 .
Figure 12.Values of AoA determined by the local velocity after the self-induction subtraction at various azimuthal locations on the r/R = 95% ring in the rotor plane; the dashed black curves near the blade axis as well as the effective AoA are estimated ( 7m/s V ∞ = , 72rpm ω =).

Figure 13 .Figure 12 .
Figure 13.Contours of AoA determined by the local velocity after the self-induction subtraction in the flattened section of r/R = 95%; the solid points P1, P2, P3, P1', P2', P3' are the monitoring points of 3-point AT; the dashed circle are used for Line AT; the star points A and A' are the monitoring points of the present determination ( 7m/s V ∞ = , 72rpm ω = ).

Figure 12 .
Figure 12.Values of AoA determined by the local velocity after the self-induction subtraction at various azimuthal locations on the r/R = 95% ring in the rotor plane; the dashed black curves near the blade axis as well as the effective AoA are estimated ( 7m/s V ∞ = , 72rpm ω = ).

Figure 13 .
Figure 13.Contours of AoA determined by the local velocity after the self-induction subtraction in the flattened section of r/R = 95%; the solid points P1, P2, P3, P1', P2', P3' are the monitoring points of 3-point AT; the dashed circle are used for Line AT; the star points A and A' are the monitoring points of the present determination ( 7m/s V ∞ = , 72rpm ω = ).

Figure 13 .
Figure 13.Contours of AoA determined by the local velocity after the self-induction subtraction in the flattened section of r/R = 95%; the solid points P1, P2, P3, P1', P2', P3' are the monitoring points of 3-point AT; the dashed circle are used for Line AT; the star points A and A' are the monitoring points of the present determination (V ∞ = 7m/s, ω = 72rpm).

Figure 15 .
Figure 15.Variation of the downwash angle along the blade span: (a) downwash angle; (b) the ratio of downwash angle to nominal AoA.
is primary.Accurate aerodynamic computation is