Characterization of Wind Turbine Wakes with Nacelle-Mounted Doppler LiDARs and Model Validation in the Presence of Wind Veer

Accurate prediction of wind turbine wakes is important for more efficient design and operation of wind parks. Volumetric wake measurements of nacelle-mounted Doppler lidars are used to characterize the wake of a full-scale wind turbine and to validate an analytical wake model that incorporates the effect of wind veer. Both, measurements and model prediction, show an elliptical and tilted spanwise cross-section of the wake in the presence of wind veer. The error between model and measurements is reduced compared to a model without the effect of wind veer. The characterization of the downwind velocity deficit development and wake growth is robust. The wake tilt angle can only be determined for elliptical wakes.


Introduction
Renewable energy sources are important for satisfying future energy demands and mitigating climate change. Limitations in the production capacity of a single wind turbine have led to wind parks with multiple turbines clustered together. However, the interaction between the wind field and the rotor blades of the wind turbine leads to a flow region of reduced kinetic energy and increased turbulence levels downstream of the wind turbine [1]. This flow region, called the wake, can affect other turbines negatively by decreasing their energy production and increasing their mechanical wear [2,3]. Therefore, accurate prediction of the wake is of interest to design and operate wind parks.
Wind turbines are located within the atmospheric boundary layer (ABL), which is the lowest layer of the atmosphere and is directly influenced by the underlying surface leading to larger temporal and spatial variation of flow variables compared to the overlying free atmosphere [4]. In the context of wind turbines, the vertical gradient of the wind speed (wind shear) and wind direction (wind veer) within the ABL and their effect on the wind turbine wake is of special interest [5][6][7][8][9][10][11]. The effect of wind veer resulting from the interaction between the synoptic scale pressure gradients, the Coriolis forces, and surface friction can affect the wake recovery of large wind turbines [12,13]. Stable conditions can enhance this wind veer [14] and strong wind veer can also result from katabatic flows in nocturnal or stable ABLs (e.g., [15,16]). Field measurements with active remote sensing instruments showed qualitatively [17,18] and quantitatively [19] that the wake shape is skewed in the direction of wind veer.
Low-order models describing the effect of wind veer on the wake are useful for optimizing wind park layouts and developing turbine steering algorithms due to their short computation time [20].
Validating those models at full-scale turbines with in-situ measurements is challenging due to the variable wake position, its height above the ground and the considerable spatial volume to be covered. Doppler lidars are actively scanning remote sensing instruments, which can sample the line-of-sight with sufficient spatial and temporal resolution to resolve turbulence scales within the inertial subrange of the ABL [21]. Several studies showed them to be viable instruments for investigating the wakes of full scale wind turbines under a wide range of atmospheric conditions (e.g., [18,[22][23][24][25][26]). Especially nacelle-mounted Doppler lidars have the advantage of alignment with rotor axis independent of the wind direction and favorable errors compared to a ground-based setup [27].
The objectives of this study are (i) to characterize the three-dimensional structure of the longitudinal velocity field in the wind turbine wake from nacelle-mounted Doppler lidar measurements, and (ii) to validate an analytical wake model that includes the effect of veer using those measurements.

Methods
The measurement site and instrument setup described in this section is the same as the one presented by Fuertes et al. in [26], but it is repeated here for completeness, and follows the post-processing and wake-analysis methods employed. The measurements were conducted from 20 August 2017 to 16 October 2017.

Measurement Site
The wind turbine is located at the Kirkwood Community College campus in Cedar Rapids, Iowa, United States. It is a 2.5 MW Liberty C96 model from Clipper Windpower with a rotor diameter of D = 96 m and a hub height of z hub = 79 m. The optimal operating range of the turbine is at velocities between 5 m/s and 10 m/s for which a power coefficient C p = 0.37 was computed from the SCADA data and a thrust coefficient C t = 0.82 is assumed in absence of manufacturer data. The surrounding area of the wind turbine is characterized by gentle rolling hills that do not vary more than 30 m from the height of the turbine base. The two main wind directions (NW and SSE) cover the urban area of Cedar Rapids to the northwest and agricultural fields to the south and east.

Meteorological Tower
A meteorological tower was located at a distance of 900 m south of the wind turbine. The tower was equipped with cup anemometers and wind vanes at heights of 20 m, 32 m, 80 m, and 106 m above the ground [28]. The cup anemometers used were A100LK from Campbell Scientific (Logan, UT, USA) and the wind vanes used were NRG 200P from NRG Systems (Hinesburg, VT, USA). All instruments were newly installed on the tower before the campaign. The wind speed profile (u mt (z)) and the wind direction profile (dir mt (z)) were used from their measurements. The data of the instruments at 32 m boom is missing from the 23 September 2017-1010 UTC until the end of the campaign.

Doppler Lidar Setup
Two Stream Line Doppler lidars from Halo Photonics (Worcestershire, UK) were deployed on the nacelle of the Kirkwood wind turbine. The steerable scanner head of this instrument emits a laser pulse, which is scattered and Doppler shifted by aerosol in atmosphere. The backscattered light is received by the instrument and the radial (line-of-sight) velocity is estimated from the Doppler shift. The travel time of the signal provides the distance information. The laser wavelength (1.5 μm) is in the near infrared with a pulse repetition frequency of 10 kHz. A range gate length of 18 m was used with 3000 (5000) samples per estimate for the scans (stares). One Doppler lidar was measuring the inflow upwind of the turbine and the second the wake downwind of the turbine.

Inflow Scanning
The inflow scanning Doppler lidar performed a repeating scan schedule of 30-min period that consisted of four scan types (a detailed description is provided in [26]): 1. A Plan Position Indicator (PPI) scan at hub height covering azimuth angles ±60 • from the rotor axis with 28 sweeps in 6 min. The yaw angle (γ ppi ) and the incoming wind speed (u hub ) at hub height were estimated from these scans. 2. A Range Height Indicator (RHI) scan in the rotor axis covering elevation angles ±15 • from the horizontal with 28 sweeps in 6 min. The vertical profiles of the inflow wind speed (u in (z)) and longitudinal turbulence intensity (TI RHI (z)) were derived from these scans. 3. A staring beam in longitudinal and transversal direction at hub height for 9 min each. They yielded the yaw angle (γ st ), the incoming wind speed at hub height (u st ), and the longitudinal and transverse turbulence intensity (TI st,x and TI st,y ).
The redundancy in the determination of the yaw angle (γ ppi and γ st ) and the incoming wind speed at hub height (u hub , u in (z hub ) and u st ) is used for stationarity tests in the quality assurance (Section 2.6).

Wake Scanning
The wake scanning Doppler lidar scanned successive PPI scans with an azimuth range of 165 • to 195 • with an azimuth step of 2 • and an elevation range of −15 • to 15 • with an elevation step of 3 • (Figure 1). This scan pattern was repeated 24 times within a 30-min period.  The measured line-of-sight velocities were mapped on a regular spherical coordinate system with a grid resolution corresponding to the azimuth and elevation steps of the scans and a temporal average was applied to gain the mean radial velocity (rv(el, az, r)). Then a provisional mean longitudinal velocity was computed with u p (el, az, r) = rv(el, az, r) cos(az − γ ppi ) cos(el) (1) and interpolated on a regular Cartesian grid with a resolution of 10 m yielding u p (x, y, z). The subscripted "p" indicates that it is a provisional mean longitudinal velocity, because the yaw angle of the turbine and the wind veer will introduce a bias that is corrected in the following. The resolution of the spherical grid also mitigates in part the effects of turbine vibrations and rocking on the beam positions of the Doppler lidar. Assuming horizontal homogeneity of the flow outside of the wake, the effect of wind shear and wind veer on the velocity deficit is corrected in approximation with a linear profile of the inflow at each height. The longitudinal velocity deficit is then given by with and y right (z) and y le f t (z) the average y coordinate of the data points. The range 6 ≤ x/D ≤ 10 was chosen, because measurements with |y/D| > 1.5 only exist for x/D ≥ 6 due to the maximum azimuth angle of the scan pattern. Using the minimum velocity of a spanwise cross-section as done by Iungo and Porté-Agel in [29] is not sufficient for analyzing the whole wake region. Correcting the yaw angle and subtracting u p (x, y, z) from the inflow wind speed at hub height to compute the velocity deficit as done for planar scans at hub height by Fuertes et al. in [26] would introduce a bias at z = z hub due to wind shear and wind veer ( Figure 2).  (1)) and neglecting the effect of wind veer (black line-i.e., applying Equation (2) on the radial velocities) will introduce a bias compared to Equation (2) (blue line). The effect of wind shear would manifest as a constant offset.

Wind Shear and Wind Veer
Similar to the computation of the longitudinal velocity deficit, the wind shear and wind veer can be estimated in approximation from the measurements of the wake scanning Doppler lidar with and For the 30-min periods that passed quality assurance (see Section 2.6), a comparisons with u in (z) and u mt (z) showed that u w (z) had an absolute mean error of less than 0.5 m s −1 across all heights where comparison was possible. Comparison of dir w (z hub ) with γ ppi showed a mean absolute error of 2.6 • , which translates to an error of 0.027 • m −1 for the wind veer across the rotor region and comparison of dir w (z) with dir mt (z) showed a mean absolute error of 3.6 • or 0.037 • m −1 . However, those errors include errors arising from the spatial separation of the measurements and temporal averaging difference (the wake measurements are averaged over 30 min and the inflow measurements only 6 min). In the following, the wind veer (ϕ w ) is defined by the slope of a linear fit to dir w (z) in the region 0 ≤ z/D ≤ 1.

Volumetric Wake Analysis
The aim of the wake analysis is to characterize the wake structure with a small set of parameters. To this end, the method of fitting a Gaussian to the velocity deficit introduced by Fuertes et al. in [26] for planar scans was extended for volumetric measurements by fitting a 2D-Gaussian function to spanwise (y-z) cross-sections of the wake at downwind distances of 2 ≤ x/D ≤ 10. The 2D-Gaussian function is given by and with the amplitude C (m s −1 ) as the maximum velocity deficit, the standard deviations σ y and σ z (m) giving the wake width along the principle axes, the wake deviation from the rotor centerline y 0 and z 0 (m), and the wake tilt angle α (deg) between the vertical and the major axis of the wake ( Figure 3). The trigonometric functions have the argument α − 90 to agree with the geometry shown in Figure 3. The fit uses a weighted non-linear least square regression with a weighting function that corresponds to 50% increased standard deviations. Initial values at x/D = 2 were C = 3, σ y = σ z = 50, y 0 = z 0 = 0, and α = 90 and for x/D > 2 the results of the previous cross-section were used as initial values. The coordinate system is a right-hand Cartesian system with the x-axes pointing from the wind turbine in the downwind direction. The cross-section view has positive y values to the right and corresponds to an observer standing downstream of the turbine looking towards the turbine (i.e., looking against the direction of the x-axis). The schematic illustrates a tilted wake for clock-wise rotation of the wind direction with height, which is common in the ABL.

Wake Model with the Effect of Wind Veer
Abkar et al. developed an analytical model for the mean longitudinal velocity deficit in the wake region in the presence of vertical wind veer [30]. The basis of this model is a mass-and momentum-conserving model [31], which was extended to incorporate a skewing of the wake due to a lateral velocity component caused by the wind veer. The model is given by with and where ϕ(z) = ϕ w z is the incoming wind angle across the rotor region, C t is the thrust coefficient introduced in Section 2.2, and k * y and k * z are wake growth rates for the lateral and vertical direction, x 0 is the length of the near wake, and α * = 2.32 and β * = 0.308 in analogy with jet flows [32]. Similar wake growth rates in lateral and vertical directions that follow a linear relationship with TI x,st are assumed following common practice in analytical modeling [33]. The growth rates are then given by with the constant parameter a that will be determined later as a = 0.30 in Section 3. In order to isolate the effect of wind veer for the comparison with the measurements, we use y 0 and z 0 from the 2D-Gaussian fit to remove effects of the wake deviation from the centerline.

Quality Assurance
Quality assurance followed a two-step procedure. The first part is a programmatic filtering following the criteria outlined in Fuertes et al. (2018) that rejects data if the inflow has a non-stationary wind speed (by comparing u hub , u st , u in (z hub ), and u mt (80 m)), a non-stationary wind direction (by comparing γ ppi and γ st or a change in the orientation of the turbine), a non-stationary turbulence intensity (by comparing TI st,x and TI RHI (z hub )), or a Doppler lidar signal-to-noise ratio (SNR) below −16 dB . Further rejected were 30-min periods with a wind speed outside of the optimal operating range of the turbine or if the rotor axis was not aligned with the wind direction (γ ppi > 5 • or γ st > 5 • ). The successful 30-min periods were then subjected to a visual inspection of the raw data to reject intervals with (i) an horizontally inhomogeneous flow outside of the wake, (ii) an influence of the wake on the regions |y/D| > 1.5, (iii) hard targets within the scan area, or (vi) a too short range of the wake scanning Doppler lidar. Lastly, 30-min periods were rejected if the coefficient of determination of 2D-Gaussian fit was below 0.96 throughout the fitted volume.

Results and Discussion
The data set consists of 847 30-min periods with volumetric wake measurements between the 20 August 2017 and the 16 October 2017 of which 88 passed the first step of the quality assurance. From those 88 30-min periods, we accepted 25 in the second step of the quality assurance. This results in a validity rate of 3%, which is lower than the validity rate of 6% found by Fuertes et al. (2018) for 2D scans during the same time. This is primarily due to the more complex volumetric scan pattern, which required more strict data filtering criteria. The reduced SNR due to the usage of 3000 pulses per ray instead of 5000 might have an effect on the validity rate as well. The remaining intervals cover weak to medium turbulence cases (TI st,x < 0.15). Figure 4 shows two examples of the longitudinal velocity deficit retrieved from the measurements. The first example features a case with weak wind veer (Figure 4a-c) and the second example a case with strong wind veer (Figure 4d-f). An increased tilt angle and a more elliptical shape of the wake can be visually observed from the measurement data and the 2D-Gaussian fit in the strong wind veer case compared to the weak wind veer case (Figure 4c,f). These qualitative observations will be investigated systematically and extended in the following paragraphs. Further, the model given by Equation (1), whose predictions are also shown in Figure 4c,f, will be evaluated, too. For this purpose, we will concentrate on the region of 4 ≤ x/D ≤ 8 downwind of the turbine. This region was chosen because for x/D < 4 the field of view of the scan is too narrow to cover whole wake leading to an uncertain fit of the Gaussian (which we defined as an R-squared value below 0.96) and for x/D > 8 the rejection of measurement data with a SNR < −16 dB led to gaps in the data in some cases.
First, the wake growth and the velocity-deficit recovery with downwind distance from the turbine are analyzed. The equivalent wake width σ eq (x) = σ y (x)σ z (x) is used to estimate the wake growth for an elliptical wake from the 2D-Gaussian fit [12]. The equivalent wake growth rate k * eq was estimated from σ eq (x) with a linear fit for 4 ≤ x/D ≤ 8. The found relationship k * eq = 0.30TI x,st agrees with the findings for planar horizontal scans during the same period (Figure 5a). The faster wake growth at high turbulence intensities due to the stronger entrainment of momentum into the wake is well established in literature (e.g., [22]). The slight underestimation compared to the planar horizontal scans is within the 95% confidence bounds of the linear fit. Previous numerical [33] and experimental [26] studies have shown that k * eq increases approximately linearly with TI x,st for moderate to high TI x,st . However, the linear relationship between wake growth and ambient turbulence levels might not be valid at very low turbulence levels (TI x < 0.06 according to [33]), because it only accounts for wake growth from ambient turbulence, but not for wake growth from turbulence created by the wind turbine itself. This interpretation is supported by the positive residual of the linear fit at small TI x,st (Figure 5b) and is in line with findings from large-eddy simulations (Figure 5a). Notwithstanding, it is possible that k * eq = 0.30TI x,st holds for sufficiently large TI x and at small TI x the wake growth rate obeys to a different regime (e.g., levelling off to a constant value or depending on wind turbine parameters).   (7)) and the black line a contour line of the prediction by the wake model (Equation (11)).
Based on the relationship of k * eq = 0.30TI x,st , the velocity deficit is predicted with the far-wake model and compared with the measurements (Figure 6). The model agrees well with the measurement data for the weak wind veer case, but shows major deviations for the strong wind veer case (Figure 6a,b). The case with strong wind veer had a turbulence intensity of only TI x,st = 0.03, for which the model presumably underestimated the wake growth due to turbulence created by the wind turbine as mentioned above. Using the wake growth parameters and near-wake length estimated for this period by the 2D-Gaussian fit shows better agreement (Figure 6b). The measured velocity deficit agrees on average for all 30-min periods with the model, but for individual cases the model can exhibit deviations of up to 15% of the incoming velocity at hub height from the measurements (Figure 6c).
In Figure 4, it was qualitatively observed that the wake appears more elliptical for the case with strong wind veer compared to the case with weak wind veer. This is quantified by investigating the relationship between the wind veer and the ratio of the minor wake width (σ s ) to the major wake width (σ l ) as a measure of the wakes circularity. The ratio of the wake widths was determined from the 2D-Gaussian fit with σ s σ l = σ y /σ z , for σ y < σ z σ z /σ y , for σ y > σ z and averaged over 4 ≤ x/D ≤ 8. The distinction of the two cases is necessary, because it was observed that the 2D-Gaussian fit could rotate by 90 • and swap its principle axes if σ x ≈ σ y (an example is shown in Appendix A). Figure 7 shows that the ratio of two principle axes of the Gaussian fit approaches unity for weak wind veer indicating a more circular shape of the wake and decreases for strong wind veer indicating a more elliptical shape of the wake. The analytical model prediction for the principle axes ratio shows agreement with the measurements (Figure 7). This effect can be explained by the fact that the transversal wind component associated with wind veer advects the wake in opposite directions below and above the hub height [30]. Next, it is investigated whether the tilt angle of the wake increases with increasing wind veer and how it develops with the downwind distance from the turbine. The 95%-confidence interval of the wake tilt angle increases for decreasing wind veer (Figure 8a), because the wake tilt angle is ambiguous for a circular wake. For the cases with a confidence interval smaller than 25 • , the mean wake tilt angle increases with the wind veer and the data points scatter around the model prediction (Figure 8b). Large eddy simulations of a stable boundary layer with a wind veer of 0.165 • m −1 across the ABL had an average α of 59 • for 4 ≤ x/D ≤ 8 [12], which is lower than the closest measurement of this data set (α = 62 • for a wind veer of 0.14 • m −1 ), but within the confidence interval. Further, an increase of α with increasing downwind distance from the turbine is expected [12]. The downwind development of the wake tilt angle was estimated from the slope a linear fit to the wake tilt angle for 4 ≤ x/D ≤ 8. The expected increase is not recovered from the measurement data (Figure 8c), because the aforementioned uncertainty of the fitted α leads to errors that can be more than twice as high as the expected trend. This serves to illustrate the limits of the wake characterization with this method for the case of the wake tilt angle. shows the (equivalent) wake growth rate k * eq as a function of the streamwise turbulence intensity TI x,st . The blue crosses show results from the volumetric measurements and the red dashed line a linear fit with k * eq = 0.30TI x,st . The black dashed line shows the result from planar wakes scans during the same experiment presented by Fuertes et al. [26] and the black squares show the results of large-eddy simulations by Wu and Porté-Agel [34]. The lower panel (b) shows the residuals of the linear fit (blue crosses) as a function of the turbulence intensity together with a moving mean spanning 7 data points (red line). Model with k eq * from meas.
Model with k eq * and x 0 from meas. Lastly, it is investigated how including the effect of wind veer in the analytical model affects the error of the model. Therefore, the root-mean-square error (RMSE) between model and measurement was computed for the volume 4 ≤ x/D ≤ 8, −0.75 ≤ y/D ≤ 0.75 and −0.5 ≤ z/D ≤ 0.75. The model without the effect of wind veer (Equation (11) with α = 0) has an average RMSE of 0.95 m s −1 or 9% based on the incoming velocity at hub height compared to the measurement data. For comparison, the 2D-Gaussian fit has an average RMSE of 0.28 m s −1 or 2.5% of the incoming velocity at hub height. The model including the effect of wind veer has a slightly lower average RMSE of 0.85 m s −1 or 8%. However, Figure 9 shows that the error reduction by including the wind veer in the model increases if the wake is more elliptical. This behaviour is expected from Equation (11) given that circular wakes are connected to weak wind veer ( Figure 7) and it can be concluded, that only strong wind veer has a pronounced effect on the wake. Moreover, the analytical model that includes the effect of wind veer could reduce the errors for strongly elliptical wakes to the same level as for circular wakes. The shown results are subject to errors that are quantified and discussed in the following. The velocity deficit field derived from the Doppler lidars can be subject to averaging and interpolation errors and errors arising from a violation of the stationarity assumption. The averaging and interpolation errors were estimated to be below 1.9% of the inflow wind speed at hub height on average, with a maximum of 7.1% [27]. Non-stationarity was minimized with rigorous quality assurance as described in Section 2.6. The degree of the remaining non-stationarity is estimated from the range of mean values during sub-periods of the investigated 30-min periods. The mean range between γ ppi and γ st is 1.9 • (maximum 2.9 • ), 0.8% (maximum 3.2%) for TI st,x and TI RHI (z hub ), and 0.36 m s −1 (maximum 0.85 m s −1 ) for u hub , u st , u in (z hub ), and u mt (80 m). Non-stationarity-especially of the wind direction-can bias the estimated wake growth rates and recovery of the velocity deficit with downwind distance. We assume that the effect of the remaining non-stationarity on the results is small, because in the case of u hub the variability is considerably smaller than the smallest absolute values (see Section 2.6) and in the case of TI st,x the variability is small compared to the range of values used to estimate the relationship with k * eq (see Figure 5). In the case of γ ppi the variation is mitigated by the grid resolution of the wake characterization. On the model side, errors of u hub will directly lead to an error of the maximum velocity deficit (Equation (12)) affecting the whole wake region, errors of TI st,x will result in an over-or underestimation of the velocity deficit with increasing downwind distance, and errors of the wind veer will lead to a wrong tilt angle and subsequently a wrong spatial position of parts of the wake (this does not imply that input parameters without errors would provide a perfect prediction). All of those effects will lead to an increased RMSE in Figure 9.

Summary and Conclusions
It was demonstrated that volumetric measurements of the longitudinal velocity deficit with nacelle-mounted Doppler lidars are feasible in practice and methods for wake characterization were presented. Fitting a 2D-Gaussian function to spanwise cross-sections of the wake characterized it in terms of amplitude, widths, and tilt angle as well as the development of those parameters with the downwind distance from the turbine. The wake characterization required temporal stationary and horizontally homogeneous conditions, which was satisfied for 3% of the data in this study. The characterization of the amplitude and widths was reliable, but the wake tilt angle could only be determined for an elliptical wake shape and its downwind development is prone to errors. Further, it is emphasized that the vertical structure of the inflow, especially the wind shear and wind veer, is essential for retrieving the velocity deficit from volumetric measurements and should be accounted for in the scan design or by using additional instruments.
In agreement with simulations and wind tunnel experiments, we observed a skewed and tilted wake in the presence of wind veer. It was shown that the wake becomes more elliptical and tilted for increasing wind veer as expected from theoretical considerations. Improved agreement between measured velocity deficit and an analytical wake model was found if the effect of wind veer is included in the model. If strong wind veer is present, the RMSE between model and measurements is reduced by up to 34% compared to a model without the effect of wind veer. Including the effect of wind veer in the analytical model reduced the errors for elliptical wakes to the same level as for circular wakes.
An example of a 2D-Gaussian fit to the longitudinal velocity deficit with approximately equal major and minor axis is shown in Figure A1. For such cases, the 2D-Gaussian fit can arbitrarily rotate by 90 • and α has large confidence bounds indicating the weak certainty. These cases are connected to weak wind veer and highlight a problem with fitting Equation (7) in cases without a pronounced elliptical shape. However, fit parameters other than α were unaffected by this problem (e.g., wake depth or wake width).