Large-eddy simulation of multiple wakes in offshore wind farms

In this study, high-resolution large-eddy simulations of two German offshore wind farms are performed with the LES model PALM in which the turbines are parameterized with an actuator disk approach. The simulation of a single wind farm is realized by a stationary model domain with a turbulent inflow. Different atmospheric stratifications from slightly stable to unstable have been simulated and their effect on the wake characteristics is investigated. The results show a clear development of the flow within the wind farms with large differences between single, double and triple wakes. The wake deficit is increasing up to the second wake for the neutral and stable boundary layer and up to the third wake in the convective case before reaching a nearly constant value. This is explained with the turbulence intensity being much higher behind the second and subsequent turbines compared to the wake of the first turbine. With increasing atmospheric stability the wake deficits are shown to be stronger due to reduced turbulence.


Introduction
Many offshore wind projects have been integrated into the transmission grid during the last few years. According to the European Wind Energy Association (EWEA), Europe's wind energy capacity is expected to grow up to a total of 400 GW in 2030, out of which 150 GW is expected to be located offshore [1]. Most of the projects are concentrated in clusters in rather small regions, mainly in the Southern North Sea and the Western Baltic Sea. For finding the optimal layout of wind farms, a deep knowledge of inter and intra wind farm effects is required. Thus, the modeling of wind turbine wakes is a key task for the energy yield prediction of wind farms as well as for the design of new wind farms. So far, most studies used Reynolds-Averaged-Navier-Stokes (RANS) models to investigate wake flow within wind farms. In this type of model the turbulence cannot be resolved but has to be parameterized completely. Only mean quantities are computed. Current high performance computing (HPC) clusters allow the use of large-eddy simulation (LES) models for wake simulations. With turbulence resolving LES it is possible to study the interaction between the turbulent atmospheric boundary layer and the turbulent wakes as well as the interaction between the individual wakes in a wind farm. Whereas the large scales of turbulence are resolved, the turbulent flow around the rotor blades requires much finer model resolutions in the order of 10 −3 -10 −2 m. Such a high resolution is not feasible for wind farm simulations, even with the currently most powerful supercomputers. Thus, the effect of the wind turbine on the flow has to be parameterized. During the last few years, first LES studies of wind farms have been published. Jimenez et al. [2] presented first results of an LES study with an actuator disk model which was refined by Calaf et al. [3] who simulated an infinite wind farm by using periodic boundary conditions. Similar studies were performed by Meyers and Meneveau [4]. Porté-Agel et al. [5] compared different wind turbine models with wind-tunnel measurements and simulated the flow through a subset of five turbines of an operational onshore wind farm. Recently, Churchfield et al. [6] simulated the offshore wind farm Lillgrund with 48 turbines using non-periodic boundary conditions with a turbulent inflow. In Churchfield et al. [7] the flow through two turbines was investigated for different stratifications and roughness lengths. Most of these studies considered only neutral stratification and simplified wind turbine arrays. In this study we use the LES model PALM to investigate the flow through existing offshore wind farms for different atmospheric stabilities which has, to our knowledge, not been done before.

Model and wind farm description
The simulations in this study were performed with the parallelized LES model PALM [8]. The wind turbines are parameterized with an actuator disk model (ADM) [3]. The ADM represents the impact of the rotor on the flow as a porous disk which acts as a homogeneous momentum sink with a thrust force F T acting against the mean flow: with the thrust coefficient C T , the rotor area A, the axial induction factor a and the temporally and rotor disk averaged velocity in direction of the mean flow u r . This model with its uniformly loaded actuator disk represents a strong simplification and neglects rotational effects. It is, however, widely used for RANS and LES due to its simplicity and capability to deliver reasonable results with rather coarse grids. A much more complex actuator line model (ALM) is implemented in PALM as well. Though the ALM delivers more realistic results than the ADM, concerning the structure of the near wake, it is not feasible for simulations of large wind farms, as it requires several times more computing time compared to the ADM. For the results shown in this paper, the effect of towers and nacelles are neglected. We use non-periodic inflow boundary conditions with a turbulent inflow. Initial boundary layer turbulence is created by a cyclic precursor run without any turbines. The results of the precursor run are used for the initialization of the main run during which the turbulence is constantly recycled.
In this study, we investigate the flow within the two German offshore wind farms alpha ventus and EnBW Baltic 1. Germany's first offshore wind farm alpha ventus (in the following AV ) is The results of the AV simulations will only be briefly discussed, the focus is on the EB1 results.

Simulations of alpha ventus
For the two AV simulations, a flow from about 270 degrees has been simulated so that a full wake situation arises. Two cases with different atmospheric stabilities have been simulated: a near neutral situation with a surface heat flux of w ′ θ ′ s = 0.005 Kms −1 and a slightly unstable situation with w ′ θ ′ s = 0.03 Kms −1 .
The wake deficits of the 12 turbines can clearly be seen in Fig. 2. They seem to be larger in the less unstable case, where the recovery of the wake deficit needs a longer distance. A roll-like structure exists in the less unstable simulation which makes the comparison of the two cases difficult. These organized structures are known as roll convection and commonly occur in neutral to moderately unstable atmospheric boundary layers [10]. The stability parameter z i /L has a value of -1.7 in the simulation with w ′ θ ′ s = 0.005 Kms −1 which is in a range where roll convection has been observed. In the simulation with w ′ θ ′ s = 0.03 Kms −1 , z i /L is -12.2, which is probably too large for roll convection. The roll-like structure comes along with a pronounced updraft region that moves very slow (due to the very small v-component) and results in a convergence of the ambient flow in the central part of the model domain. This leads to a deflection of the wakes in positive y-direction in the lower part of Fig. 2 and in negative y-direction in the upper part.
To allow for a better comparison of the wake effects between the two simulations, the velocities Figure 3. Cross-section of the minimum wake velocities at hub height for all four turbine rows and both AV simulations. Continuous lines: of the wake center lines have been compared for each row of turbines (Fig. 3). The wake center line was defined as position of the minimum velocity along the x-direction. It can clearly be seen that the wake deficits are similar behind the first wind turbine in both simulations. The recovery of the wake is however faster in the case of the larger near-surface heat flux, so that the wake deficits behind the second and third turbine are larger in the less unstable case. Another interesting effect is that the wake recovery after the second turbine is much faster than the recovery after the first turbine. Thus, the wind speed in front of the third turbine is about the same as that in front of the second turbine, although the wake deficit behind the first turbine was much smaller than that behind the second turbine. This observation is in agreement with reports from operating wind farms where often a huge decline in the power output is observed between the first and the second row of turbines, while there is no huge difference between the power output of the second and third row (see e.g. [11]). Furthermore, it was found that the turbulence intensity is higher behind the second turbine than behind the first turbine which can explain the differences in wake recovery. This will be analyzed in more detail in the following section for the EB1 results.
The results of the AV simulation point out that there is a development of the flow within the wind farm. The situation between the first row of wind turbines is different from the situation behind the second row of turbines. However, the recovered velocities reaching the second and third row are similar. Therefore, it is interesting to study the wind farm flow in larger wind farms with multiple wakes of four or more turbines in a row as in the EB1 simulations presented in the following section. With a flow from the North, a variety of wake situations from a single wake up to a sixfold wake can be studied in one simulation (see: EB1 map in Fig. 1).

Simulations of EnBW Baltic 1
The size of the model domain, especially of the precursor run, is crucial for simulations with turbulent inflow. If the domain size is too small, the largest turbulent structures cannot be resolved and the turbulence spectrum is truncated. To estimate the necessary domain sizes for the wind farm simulations, a sensitivity study concerning the domain size of the precursor run has been performed for both neutral and convective boundary layers. A model domain of 2560 m x 2560 m was sufficient to cover the largest turbulent structures in the neutral boundary layer simulation whereas even a model domain of 10240 m x 5120 m was too small for the convective case. The simulated convective boundary layer was typical for a situation over land with strong heating from the ground. A marine boundary layer, as in the offshore wind farm simulations, is usually characterized by weaker and smaller convection due to a smaller surface heat flux and a shallower boundary layer. It can be expected that the convective structures will be smaller, so that a domain size of about 12000 m x 6000 m, as it is planned for the EB1 precursor run, is sufficient. Based on these findings, the precursor runs for the three EB1  Figure 4 shows the flow structure at the end of the three pre-runs (after 24 hours of simulated time). In the convective case, elongated large structures are visible with a size of up to 10 km in flow direction and 1-2 km normal to the flow direction. This is again roll convection as it has been observed in the AV simulations. Here, z i /L is -2.13. The structures are quite persistent, so that even the one hour time-averaged flow field features distinctive bands aligned with the mean flow. No roll convection is visible in the vertical velocity field of the neutral case. However, the flow field in mean flow direction shows coherent regions of high velocity, much smaller in size than the convective rolls. These structures are known as streaks which originate from the surface layer [12]. They are not to be mistaken for convective rolls which span throughout the whole boundary layer, whereas the much smaller streaks only occur in the lower boundary layer and are not necessarily accompanied by coherent up-and downdrafts. Similar but smaller and weaker streaks can also be found in the stable case. The model domain of the pre-runs is in any case large enough to cover all scales of turbulence. Again, the need for a much larger domain in the convective boundary layer simulation is demonstrated.  behind the first turbine which is in agreement with the results of the AV simulations. Also shown in Fig. 5 is the mean turbulence intensity T I u,v,w = σ u,v,w /u inflow which is greatly enhanced in the wake of the turbines. Behind the second and subsequent turbines the turbulence level is significantly higher than in the wake of the first turbine. This can be seen more clearly in Fig. 6 (right). In Fig. 6, some differences can be observed between the three simulated atmospheric stabilities. The convective case features much weaker wake deficits compared to the neutral and stable case. The deficits are strongest for the stable case. The turbulence intensity within the wakes is similar for all cases. The undisturbed flow has a slightly enhanced turbulence intensity in the convective case. The slight deflection of the wakes in the stable case (about 1.7 • ) occurs due to an inflow not exactly aligned with the x-axis. This happens as the same geostrophic wind was prescribed for the neutral and stable cases but the turning of the wind direction within the Ekman layer is stronger in the stable case. The minimum wake velocities along x for the five-turbines row are compared in Fig. 7. For all cases, a maximum wake deficit is reached behind the second turbine. Behind the subsequent turbines, the minimum velocity is nearly constant or even slightly increasing. The largest wake deficit is found for the stable case with a minimum velocity of about 40-45 % of the inflow value. The convective case shows the weakest wake deficit with a minimum of 50-55 % of the inflow velocity. Interestingly, there is no wake recovery behind the first turbine for the neutral and stable cases and very little recovery for the convective case. Figure 6 (left) reveals that there is indeed a wake recovery in the outer parts of the wake due to mixing with the surrounding flow but the distance to the second turbine is too small that the recovery can penetrate to the wake center. Behind the subsequent turbines the wake recovery is present for all cases and very pronounced. Hence, the inflow velocities for the third and subsequent turbines are similar to that of the second turbine. Consequently, the wake deficit is not increasing further. Downstream of the wind farm, the 90 % recovery level is reached about 20 to 25 D downstream of the last turbine. The 95 % level is reached not until about 40 D downstream of the last turbine for all cases, and a full wake recovery requires at least 80 D.

Figure 7.
Cross-section of the minimum wake velocities at hub height for the five-turbines row. The turbine positions are indicated by arrows.
In Fig. 8, the turbulence intensity profile along x is shown for the three velocity components. The turbulence intensities are averaged over a 500 m section along y which covers most of the wake. Again, it can be clearly seen that the turbulence intensity is much higher in the wake of the second and subsequent turbines than in the wake of the first turbine. This can be partly explained by the weak velocity deficit behind the first turbine leading to much a weaker wind shear compared to the subsequent wakes. Besides, the wake is rather homogeneous behind the first turbine before instabilities start to deform the wake and meandering sets in. This happens several diameters downstream of the first turbine, only shortly before the second turbine is reached which then generates additional turbulence. This is also the explanation for the missing wake recovery behind the first turbine and the much stronger recovery behind the second turbine. The turbulence intensities further increase behind each subsequent turbine with the maximum occurring in the wake of the last turbine of the row. This leads to a faster recovery of the velocity deficit and the slight increase of the velocity minimum which can be seen in Fig. 7. The turbulence intensities of the u-component are significantly higher within the wake for unstable stratification compared to the neutral and stable cases. Contrary to what one would expect, the v-and w-components show little differences between the cases, except for a higher background turbulence level in the convective case. An explanation may be that the vertical wind shear in the range of the turbine is significantly larger in the neutral and stable cases than in the convective case. Thus, additional turbulence can be generated by shear in the neutral and stable cases. Figure 9 demonstrates the fluctuation of the inflow wind speed, the turbines are exposed to, for the convective case. For a single location, the wind speed varies significantly on very short time scales of seconds. However, longer-scale fluctuations of several minutes are visible as well. The amplitude of fluctuation increases for turbines in the wake which is due to enhanced turbulence intensity and wake meandering. The thrust generated by the wind turbines in the model is based on the wind speed averaged over the rotor area which fluctuates, of course, much less. However, the thrust is proportional to u 3 which causes, in combination with the variable C T , very large fluctuations of the generated thrust. Figure 9 shows that sometimes a turbine in the wake can generate even a larger thrust than a turbine in the undisturbed flow whereas most of the time the thrust is much lower for turbines in the wake.

Summary and outlook
We investigated the effect of atmospheric stability on wake characteristics by simulating the flow around and within the two German offshore wind farms alpha ventus and EnBW Baltic 1 with the LES model PALM. The results show a clear development of the flow within wind farms with an increase of the wake deficit from the wake of the first turbine towards the wake of the second turbine of each row of turbines aligned with the mean flow. For the subsequent turbines, the minimum wake velocity remains constant or is even slightly increasing. These results are in agreement with measurements (see e.g. [11]) and other modeling studies (see e.g. [6]) and can be explained with the development of the turbulence intensities within the wind farm. Behind the first turbine of each row the turbulence intensity is barely enhanced compared to the turbulence intensity of the incoming flow, so there occurs very little wake-induced turbulence. Thus, the recovery of the wake deficit is rather small, and the inflow velocity for the second turbine is significantly lower than for the first turbine which leads to a stronger wake deficit behind the second turbine. The turbulence intensity is still increasing considerably behind the second turbine, resulting in a much more pronounced wake recovery and inflow velocities for the third turbine similar to the second turbine. Hence, the wake deficits behind the subsequent turbines do not further increase. Atmospheric stability affects the wind farm flow, as a comparison of three simulations with different stabilities shows. For a convective (unstable) boundary layer, the wake turbulence is found to be higher compared to a neutral or stable boundary layer. Hence, the wake deficit is strongest for stable stratification and weakest in the convective case. The wakes are persisting over a long distance downstream of the wind farm. For the rows with two or more turbines, a value of 90 percent of the inflow velocity is not reached before 20 to 25 rotor diameters downstream of the last turbine. For a full recovery of the wake at least 80 rotor diameters are required. It has also been shown that the inflow a turbine is exposed to, and hence the generated thrust, is subject to strong fluctuations on very different scales.
The results demonstrate that atmospheric stability significantly affects the wake flow and power output of offshore wind farms. It is also shown that atmospheric turbulence is crucial for the wake recovery and the development of the flow within wind farms. Thus, wake models should include adequate parameterizations of both atmospheric stability and turbulence. We are currently implementing an enhanced non-uniformly loaded actuator disk model in PALM which includes rotational effects and the effect of the turbine towers. The EB1 simulations will then be repeated and validated with measurements within the wind farm and from the offshore met mast FINO 2. Different flow directions with varying distances between the turbines will be investigated and validated with LiDAR measurements within the wind farm.