Implementation of the BEM skewed-wake model within the multibody aero-elastic solver Cp-Lambda

To account for the problem of an azimuthally constant induction in the BEM method, which influences on incorrectly predicted aerodynamic loads in the yawed flow, a skewed-wake model implementation to the BEM method has been performed. The numerical aerodynamic loads have been compared with the wind tunnel data of the NREL Phase VI and against another numerical campaign. At first, the model predictions have been validated against experimental data performed with aligned flow conditions, showing a reasonable match. Successively, the model predictions are validated against experimental results obtained with the wind turbine yawed. Results show, a possible better prediction of loads at yawed flow with Skewed-Wake correction, however the method does not overall correlate better, compared to the BEM method with implemented local variability of the induction factor.


Introduction
The turbines in a wind farm are often experiencing yawed operating conditions or highly sheared inflows, which could lead to a consistent reduction of the power production and to an increase of the fatigue loads [1]. Moreover, these operating conditions will become more common in the next future, since it has been found that yawing a wind turbine away from the wind is an effective way of altering the path of the wake and thus boosting the overall wind farm power [2,3,4,5,6].
To evaluate the performance/loads of a yawed wind turbine, the most widely used approach is coupling the Blade Element Momentum (BEM) method to an aeroelastic solver. Such a solution is indeed preferred to more advanced approach, like computational fluid-structure interaction and flow simulation, thanks to its low computational cost and its relatively good accuracy. However, one of the main disadvantages of the method is that its underlying hypothesis are violated when the flow is not aligned with the rotor disk. Hence, there is a need to improve the classical BEM method in order to enable its use in yawed conditions. In the present work, the BEM Skewed-Wake model developed in [7] has been implemented within the multibody aero-elastic solver Cp-Lambda [8], and the numerical results are compared with the data from a well-known experimental campaign.

Objectives
One of the principal limitations of the BEM is the assumption of an azimuthally constant induction. Although most formulations assume that the induction can vary radially, the induction is assumed to be the same along each annulus of the rotor disk. In reality, it is well known that in general the induction can be far from azimuthally constant, although for the flow aligned with the rotor axis experimental examination has shown that the assumption of radial independence is acceptable [1].
In yawed conditions, i.e. when the rotor is not pointing straight into the wind, the angle of attack of the rotating blade is constantly changing and the induced velocity varies azimuthally and radially, which affects the efficiency of the energy conversion process [1]. In order to account for these phenomena, a modification to BEM method implemented within the multibody aero-elastic solver Cp-Lambda has been performed. The code already considered an azimuthally-radial dependent induction factor, but was not accounting for the effects produced by a skewed wake on the axial and tangential component of the inflow velocity at the rotor disk, which was found to improve the accuracy of some simulation results against experimental data [7]. Hence, the BEM Skewed-Wake model developed by [7] has been implemented within Cp-Lambda software, and the numerical predictions have been compared with the wind tunnel data of the NREL Phase VI [9], which provides vast accessible datasets of yawed inflow conditions.

Methodology
The multibody aero-servo-elastic solver Cp-Lambda is a complex algorithm where aerodynamic loads on the numerical bodies (eg. Blades, tower, …) are accounted for in the FEM analysis of the elastic structure. At every time step the set convergence norm must be achieved in order to complete the analysis. Aerodynamic loads are evaluated through the routine of extracting the data of the aerodynamical model, updating and interpolating the local velocities at every computational point, accounting for the rotation and the deformation of the part and performing the computation with updated data. All the bodies that contribute to the aerodynamic loads, are represented by a lifting line, which is composed from a prescribed number of airstations (AST). These are the computational nodes in the aerodynamical model, where the algorithm for loads based on a 2D strip theory is performed. The BEM method is one of the possibilities used here, to acquire the local flow velocities, which affect the aerodynamic angles and therefore forces. It is a combined blade element theory and momentum theory and is as such only applicable for cases with radially and azimuthally constant induction factor [1].
The implementation of the method is shown in Fig 1. Firstly, the AST properties from the previous step are extracted (pitch and twist in Fig 2, induction factor ,…) in order to evaluate the effective angle of attack ( shown in Fig 2) from the local velocities and deformations in the AST. If the Prandtl tip and hub loss model is included, the Prandtl coefficient is set. With effective angle of attack computed, lift drag and moment coefficients are interpolated from the lookup tables. If axial induction factor is above the threshold, Glauert's empirical thrust coefficient factor is evaluated. With all of the data extracted, induction factors are computed and residual is compared with the set convergence norm, then AST properties are updated. As the BEM method is fast and robust, modifications to account for skewed wake correction, have been proposed as presented in [10]. The correction on the blade element method can be made, by including for local velocities at airstation. This is already accounted for in the Cp-Lambda software. Many corrections on the momentum side of equations are based on the Glauert's correction [11], where it was proposed, that the induction factor would be ( ( ) ) , where is the average induction factor, K is a function of yaw angle , specific spanwise location and azimuthal angle (Fig 2). In the current work, the modified Pitt and Peters correction is used, as described in [7] ( ) (1) Where a is the local induction factor and the wake skew angle ( ) (2) In the case where the turbine is not tilted, the angle between the normal of the rotor disk and wind is the same as yaw angle . The correction for the axial induction factor is implemented at the last step of the BEM method.
For the verification of the method, the simulation outputs are compared with experimental data of the NREL Phase VI [9], where a 10m diameter stallregulated wind turbine with a power rating of 20kW has been used. The experimental campaign chosen for the validation (Sequence H), regards data obtained by testing an upwind wind turbine, equipped with two twisted and tapered blades. Aerodynamic forces were derived from the pressure sensors, which were positioned at 5 spanwise positions. The numerical model of the experimental wind turbine is built based on information from [7], [9] and [12], where geometry, controller, airfoil lookup tables and material properties are listed. Simulations within the Cp-Lambda software were performed with the standard BEM model and the new implemented Skewed-Wake model. At first, the model predictions have been validated against experimental data performed with aligned flow conditions at wind velocity of 7m/s and small shear exponent of 0.05, showing a reasonable match. Successively, the models predictions are validated against experimental results obtained with the wind turbine yawed of 20° and 40° at the same wind velocity. Aerodynamic loads in the in-plane and out-of-plane directions ( and as shown in Fig 2) are derived from lift and drag loads ( and ) for comparison with numerical [7] and experimental [9] values.

Results
The comparison of the model with the aligned flow condition (Fig3), shows a reasonable agreement with experimental values. The values of both methods coincide as is expected, since the wind is aligned with rotor axis and hence the wake skew angle and the correction is 0. The difference between azimuthal positions is also negligible and therefore not plotted.
As the largest difference between the numerical methods is at azimuthal values of 90°and 270°, only these plots are selected to be plotted in Fig 4 -6. At higher yaw Skewed-Wake model predicts the same trend as it is seen from the measurements, that the loads are higher at 270°, although the values are over predicted. Similar results were also obtained in [7], which is why the denominator of the Skewed-Wake model is larger than originally proposed by Pitt and Peters [13], to better correlate with the measurements. Over prediction at the inboard sections could be the influence of the inadequate correction/extension of the lifting line lookup tables for the effects of deep stall and rotational stall delay.      In Fig. 8. and 9., in-plane and out-of-plane load distributions are plotted against the azimuthal angle of the blade, at specific spanwise locations * + . Inboard, the difference between the models and measurements is even better seen due to the different axis values, however this only contributes smaller portion of the overall aerodynamic load. The asymmetrical trend of experimental values is here better correlated with the Skewed-Wake model at both yaw angles as the increase of loads at is better captured, however the loads are over predicted. Part of the overestimate is due to the discrepancies in the lookup tables from [7] and [9], with different corrections applied, which directly affect the computed aerodynamic loads.

Conclusions
An implementation and validation of the skewed wake model within the Cp-Lambda solver with respect to available experimental results was performed. Results show a reasonable agreement with experimental values for the aligned flow and are well in line by what was observed in the literature [7]. At some spanwise locations, a better trend of load dependence to azimuthal blade position of the skewed-wake method is seen. However, the method does not overall correlate better with the measurements, compared to the already implemented BEM method, which accounts for the local variability of the induction factor and is therefore already modified. Part of the discrepancies could be due to the fact of the different lookup tables in [7] and [9], with different corrections applied, which directly influence on the computed loads. The results imply to the fact, that a better model might be necessary to implement.
One of the drawbacks of a vast experimental campaign with many sensors, is the influence of the supported structure on the flow. It was shown that in yawed flow in some cases, the instrumentation boom and enclosures in front of the wind turbine, affect the pressure measurements due to the wake [9]. Moreover, the local flow angle at yawed wind turbine at the pressure probes, might exceed the design range for the probes. Therefore, the quantities based on the pressure measurements at azimuthal angle from 0° to 120° might be affected. Due to the measurement uncertainties, different dataset should also be investigated, to find some other trends.
All in all, the results of the implemented skewedwake model show some improvement at the aerodynamic load assessment, but for a better model prediction, further research is necessary.