Evaluating Global Blockage engineering parametrizations with LES

Under the term global blockage, the cumulative induction of wind turbines in a wind farm has been recently suspected to be responsible for observed overestimations of the energy yield in large-size wind farms. In this paper, the practice of modeling this effect after linear superposition of single turbine inductions, calculated with three of the most recent analytic models, is compared to Large-Eddy-Simulations of wind farms. We compare the models across two different farms, composed of 9 and 49 turbines, with two different heights of the atmospheric boundary layer, 300 and 500 m. The results show hat the differences between the analytical models are negligible while they substantially differ from the LES results. The linear superposition of induction consistently underestimates the velocity deficit in front of the farm with an error that increases as the wind farm size grows and the ABL height decreases. Also, when calculating the power output at the turbines of the farm, all the analytical models considered do not agree with the LES. These comparisons reveal that the farm interactions with the atmospheric boundary layer may highly outclass the turbine induction in determining the extent of the global blockage effect. Therefore, we present a first dimensional approach to the problem based on LES, aimed at simplifying its characterization.


Introduction
Despite being already conceived in the very early theory of Betz, it is only in recent years that the induction of a wind turbine was finally given a representation in engineering models that aim at modelling the flow inside of wind farms. There is, in fact, a growing awareness that upstream interactions between wind and wind turbines should no longer be neglected [1]. Supporting this theory, there is consistent evidence of wind speed deceleration upstream of wind farms suspected to cause efficiency losses, especially of the first turbines row [2]. This effect, discussed under the name of wind farm blockage or global blockage, is at our day and age still debated in science and industrial applications, and different postulates on its origin exist. An early study on wind farm power losses caused by inflow wind speed reduction can be found in [3], in there the wind farm is suspected to generate atmospheric gravity waves, similarly to mountains, capable of reducing the inflow speed through pressure perturbation. The same idea has been more recently extensively investigated by Allaerts and Meyers in [4,5,6], where by means of Large-Eddy Simulations (LES) the authors were able to better describe the extent of the phenomena and to connect it to the stratification of the atmospheric boundary layer (ABL).
A more general discussion on the existence of global blockage and its effect on power production can be found in [1]. In this work, the authors compared the wind speed measurements collected by different met-masts before and after a nearby wind farm was installed. The wind speed was observed mostly to decrease in the second period even when the met-masts measured 2÷ 10 rotor diameters upstream the wind farm. Other, more specific, in situ evidence of the phenomena can be found in the measurements with remote sensing technique of Asimakopoulos et al. [7] and especially Schneemann et al. [8]. In this novel work, the wind speed at transition piece height (≈25 m above mean sea level) is found to decrease up to 4.5 % between 40 and 10 rotor diameters upstream of the Global Tech 1 wind farm when stable stratification characterizes the ABL. As the pool of evidence for global blockage is growing larger, the discussion on whether and how to consider its effect in wind farm energy yield assessment is still quite open. An important contribution in this direction is given by the work of Segalini and Dahlberg [9]. They derived an empirical correlation describing velocity deficit for the first row at varying farm layout, from a series of wind tunnel measurements. A different approach can be found in [10,11]. In both publications, the authors propose to model the global blockage effect by linear superposition of the induction effect of the single turbines, in a similar fashion to what is currently done for wake modeling.
Nygaard et al. [11] also attempted a validation of the method. In there, the power output of the first row at the UK-based Gunfleet Sands wind farm is simulated for two particular wind direction bins at a fixed velocity. By comparing the results with averaged Supervisory control and data acquisition (SCADA) data, the authors observe that "in its present form the model is likely missing some element of the relevant physics". This important conclusion cannot, unfortunately, be considered definitive. The wind farm chosen is surrounded either by the coast or by nearby farms, in this condition the blockage effect is not the only factor determining wind velocity distributions, and phenomena such as external wind farm wakes or coastal wind speed gradients may even play a greater role. Within this work, we want to contribute in assessing how much the linear superposition of single turbine induction should be expected to suffice in the modeling of wind speed deceleration and hence wind farm power reduction connected to the effect of global blockage. As reference for the comparison we use LES capable of describing the development of the atmospheric boundary layer under physical sound conditions. With the numerical environment we remove parasite perturbations that may compromise the fairness of the comparison in measurement-based validation approaches. After a brief description of the induction models, their implementation in the Fraunhofer IWES in-house software flappy (Farm Layout Program in Python), and the LES set-up in section 2, we present results for the comparison models/LES on wind velocity and power production in section 3. As a closing argument, in section 4, we propose a simple dimensional study to characterize the wind velocity deficit the LES predict in front of the farm, seeing the necessity of improving the current models.

Models Description and set-up
Within the last years a good number of engineering models describing the velocity deficit in the induction zone of a wind turbine have been proposed. Branlard et al. [12] review most of them, offering a comparison of their performance against an actuator disk (AD) Reynolds-Averaged Navier-Stokes (RANS) simulation. The Vortex Cylinder (VC) model of [10] appears as the one offering the best accuracy. In this manuscript, we consider the VC model and two additional models with lower computational effort: the Self-Similar model of Troldborg et al. [13] and the potential flow model of Gribben and Hawkes [14], here as Rankine Half-Body.

Vortex Model (VC)
Inspired by the cylinder model for impellers of Joukowski, the Vortex Cylinder model [10] computes the perturbation velocity components by linear superposition of the effects from the three main vortexes generated in the interaction between wind and turbine blades: the straight root vortex, the straight bound lifting lines at the blades, and helical vortex trailed at the blades tip. As an infinite number of blades is assumed, the lifting lines at the rotor disk merge into a bound vortex sheet and the helical vortex can be decomposed in line vortex filaments and a circular component responsible for a tangential circulation, γ t [15]. The computation of velocity perturbation happens thanks to the Biot-Savart law integration applied to any of the vortex systems defined. Multiple concentric cylinders supporting γ t may be considered to describe the wake expansion, but, as done in [10,12], only the single cylinder matching the rotor dimension is considered here. Under the previous assumptions, only the tangential circulation is the capable of providing the axial induction. Therefore it is possible to connect γ t to the induction factor,ā, thanks to the relation γ t = −2 U refā , where U ref is the velocity of the free stream approaching the turbine, as demonstrated in [15]. In the present work,ā for the model is computed according to the polynomial expansion on the global thrust coefficient C t , , found in [10] to provide a better agreement with an AD RANS throughout the inflow region x ∈ [−2.5 D,−0.5 D], with D the rotor diameter, with respect to the one-dimension momentum theory. The VC model has been introduced in our wind farm modeling tool implementing the approach followed for FLORIS [16] available at https://github.com/ebranlard/wiz.

Self-Similar model (SS)
The axial velocity deficit at hub height, u x (x, r = 0), in the VC model takes a very simplified formulation that reads as, with x and r being the coordinates of a cylindrical frame or reference centered at the rotor hub height (z HH ). Troldborg et al., [13], have observed the induction profiles for different turbines simulated with RANS to exhibit self-similar behaviour when x > R. They hence proposed to compute the velocity deficit using (2) and describing the radial decay as Differently to what chosen for the VC model, the induction factorā is based here on the 1 D momentum theory,ā = 0.5 (1 − √ 1 − νC t ) with the constant ν = 1.1 providing better agreement with RANS simulations.

Rankine Half-Body (RHB)
Gribben and Hawks in [14] suggest that the flow around a wind turbine resembles the response of a uniform flow field to a constant momentum source. In this typical problem of potential flow theory, the main flow is inviscidly deflected around a Rankine Half-Body, the imaginary surface described by the main stagnation streamline, described in a spherical frame of coordinates (r, θ) centered at turbine hub height as with σ = 2ā A rotor U ref being the momentum source, andā again derived from the onedimensional momentum theory. Any of the cartesian velocity perturbation components is given by 2.4. Empirical wind speed reduction by Segalini and Dahlberg (Seg) As additional model for comparison, we consider the empirical wind speed reduction derived by Segalini and Dahlberg [9] (Seg) Here, u x 0 is the wind speed deficit at the front-most row middle turbine, ∆ x and ∆ y the streamwise and spanwise distances between the turbines and N rows the number of streamwise rows of the wind farm. While not explicitly stated in [9], we use here u x 0 as initialization velocity at any turbine of the farm.

Case study setup and LES
The framework of our comparison is to perform simulations for the very same wind farm layout either in AD LES and in the Farm Layout Program in Python 1 (flappy) developed at Fraunhofer IWES [17], the successor of flapFOAM (Farm Layout Program coupled to OpenFOAM) [18]. The wind turbines considered in the simulations have a constant C t of 0.86 and a constant power coefficient. Thus, they are characterized by a rather high, but not unusual, induction factor. The turbine's geometry is typical for offshore application with a rotor diameter (D) of 160 m, resulting in a rated power of approximately 8 MW, and a hub height of 110 m. Two basic wind farm setups were simulated, one with three turbines in streamwise and spanwise direction, respectively (3 × 3), and the other with seven turbines in each direction (7 × 7). The distance between the turbines in the streamwise and in spanwise directions is 7 D and 5 D, respectively. Additional simulations with a single turbine and just the first turbine row were conducted to serve as references.
The LES data was created with the Parallelized Large Eddy Simulations Model (PALM) [19]. Further details about the used model setup can be found in [20].  In flappy the final flow field is derived in an iterative approach. Any velocity perturbation computed by induction or wake models are linearly superimposed. No induction model is allowed to provide a contribution in the projected cylinder at the back of the rotor at which the induction is computed. The presence of the ground is mimicked by also modelling a second turbine mirrored at the ground surface. Wake effects are modelled with a modified version of the Bastankhah wake model as described in [21]. However, as the wake modelling only has a marginal influence on the modification of the flow field by the induction models, we will not explore the representation of the wakes further.

Results
In this section, the predictions of the engineering models are compared against the LES results. First we focus mostly on the axial velocity deficit that we define as , the reference velocity is taken as the spanwise averaged value from the LES at the inlet of the domain. Later, we look at the models effect when applied to wind farm power prediction. In this section (3.2), the empirical correlation of Segalini and Dahlberg [9] will be discussed as well.
It is important to mention that the induction models have no dependence on ABL height. Thus, the comparison with the LES of different ABL height serves as an illustration of an effect that is not reproducible with the implemented models. Figure 1 presents the distribution of the velocity deficit in the vertical plane upstream of a wind turbine. We observe a good agreement between models and LES when the single turbine  is simulated in isolation, Fig. 1(top), confirming the results of Branlard et al. [12] in their comparison with AD-RANS simulations. Qualitatively, the trend of the velocity deficit match up quite well in the region x > −2 D, especially for the VC model. For further upstream distances the models tend to overestimate the deficit slightly.

Velocity field prediction
However, repeating the same comparison for the middle turbine of the front row of the 7×7 500 case, Fig. 1(bottom), the previous difference are exacerbated, with significant underestimations of the wind speed reduction for all three models. While the wind farm in the LES induces wind speed reductions of up to 2 % at distances of 5 D upstream, the reductions in the models do not exceed 1 %. Interestingly, a good agreement remains for the near-rotor region in the VC model simulation. Hence we deduce that the linear superposition method fails to describe the cumulative effect of the wind turbines on the flow. Figure 2 provides a different visualisation of the upstream velocity reduction with a display of the mean axial wind speed at hub height. By comparison to the LES, the models clearly underestimate the velocity deficit throughout the considered domain and compared to the LES results the model difference become marginal. From the 3×3 to the 7×7 wind farm, the velocity deficit at x = −10 D grows of 3 times larger in the LES simulations, while the models predict only a 70 % increase. Similarly, varying the height of the boundary layer strongly affects the inflow velocity deficit at both the farms in the LES. The complete lack of sensitivity on this parameter proves to be critical for the model.
Given the general good agreement observed for the single turbine case, the results obtained until here point out that there are significant overhead effects contributing to wind speed decelerations in front of the farm that turbine induction cannot justify.
It is by looking at the vertical distribution of the velocity deficit, Figure 3, that it becomes even more evident that the inductions effects have too local traits and they cannot fully motivate the global blockage effect. The LES shows the velocity perturbation upstream of the first row to not only be more intense with respect to the models, but also to interest the whole vertical dimension of the ABL. Figure 3 also displays properties of the flow within the wind farm. The vertical profiles of the velocity deficit at the middle and last row in the LES reveal that the initial differences between the two ABL simulated are recovered within the wind farm, and above the rotors the wind speed may even exceed the velocity of the free flow far upstream. This effect is completely outside the analytical models' capability, even the one describing the speed-ups around the rotors.
To our understanding, the flow patterns observed up to now can only be motivated by pressure perturbations in the lower boundary layer originating, as explained by Allaerts and Meyers in [4], as a consequence of a displacement of air masses with different density aloft. Clearly, such an effect cannot be replicated with the engineering model approaches where, to derive the perturbation of the velocity field, pressure is considered to be in equilibrium.
As a conclusion to our analysis of the velocity field representation, Fig. 4 visualises the horizontal velocity contours at hub height for the 7×7 500 case. We observe that the use of an induction model enhances the similarity with the results of the LES, especially if the model accounts for lateral speed-ups like the VC model used in the comparison. This underlines once more that the models are capable of reproducing the local turbine blockage, but their linear superposition cannot adequately represent the global blockage effect, as the wind farm interactions with the atmospheric boundary layer makes the phenomenon strongly non-linear.

Wind Farm power
In the previous section, the comparison with the LES has shown the linear superposition of induction models to be far from sufficient to describe the velocity deficit caused by global blockage. The enhancement of the flow visualization may not be a sufficient reason justifying the method deployment in context such as Annual Energy Production (AEP) assessment or wind farm optimization algorithms. The many degrees of freedom involved in the solution of these problems lead to iterative procedures that require fast calculation of the farm output at a given state. Any computational cost added to this task must be compensated by clear improvements in the accuracy of the prediction. In the following section, we debate the benefits of considering an induction model in the assessment of wind farm production.
In Fig. 5(a,b) we show the percentage of power reduction at the first row and to the whole farm, by the application of the induction models and the global wind speed reduction at any of the farm's turbine, according to Eq. 5. Because induction in the LES cannot be switched off, Figure 5: Reduction of the average power production of the first row and the whole wind farm caused by modeling the induction effect, for farm 3×3 a) and 7×7 b). The references are the same simulations without modelling any induction effects. For the LES, the reference of the first row power is a simulation of the isolated first row. For the whole farm power no reference is possible for the LES case. c) The reduction of the average power production of the first row and the whole wind farm for the two different wind farms in the simulation setups with different boundary layer height. The reference case is the setup ABL 500.
we define its reference as the results for power output of a simulation where a single row of the farm is considered. This practice is identical to the reference taken for the analytical models, as removing induction in this models is equivalent to removing the effects of the downwind turbines. Such a normalization for the LES is also beneficial to remove row effects that lead to higher power production as shown in [22]. Unfortunately, it is not possible to define a similar reference for the whole farm in the LES.
The trend of velocity reduction observed in Fig. 2, can also be observed in the first row power predicted by the LES, that decreases evidently as the farm size grows and the ABL height decreases. From the 3×3 to the 7×7 the power reduction grows of 14 times in the LES at ABL 500, while in the models it merely double. A further evidence of a strongly non linear phenomena. The models also consistently underestimate the power reduction at first row for all the studied cases but the 3×3 500, in which the effect of flow blockage are minimal. This result shows that the inlet velocity deficit may be partially recovered, another sign of the importance of the pressure effects.
The power production of the whole farm is less affected by the global blockage effect than the first row in all the simulations with the induction models. This result is also confirmed by the LES. The relative power reduction across the two ABL simulated is much greater for the first row than for the whole farm, Fig. 5c, indicating as in Fig. 2 that the inflow deficit is recovered at the downwind turbines. We therefore urge caution with use of approaches based on reducing the velocity at all the turbines in the farm according to a correction measured upstream, e.g. Eq. 5. Despite we could not prove it here, we expect the methods to be biased towards underestimating the whole farm power.
To conclude the results' section, we finally discuss the impact induction models have on the whole farm power production. In our simulations, the use of neither the VC model nor the RHB model significantly affect this prediction. To our understanding, this is connected to the speed up around the rotor described in these models, that despite being subtle it is sufficient to balance the negative effect of induction on power. To this extent, only the Self-Similar model provides a prediction of the order of magnitude typically assumed for the global blockage effect. More research in this direction is still necessary but for the other two models, the significant computational cost we experienced in our simulations, especially with the VC, and the marginal influence on power production are results that do not encourage their use in AEP calculation.

A simple dimensional analysis of global blockage effect
According to our analysis, we find that wind farm size, as well as atmospheric boundary layer height, have a significant impact on the extent of velocity reduction in front of the wind farm, see Fig. 2. This motivates us to perform a simple dimensional analysis to find out whether a universal behaviour of the velocity deficit can be identified.
Our LES simulations suggest that the wind speed reduction upstream the farm, here defined as u x,loss = U x,ref − U x,LES , can be represented as a firt approach by the functional relation In Eq. 7, x denotes the axial coordinate with x = 0 at the farm first row, U x,ref the undisturbed velocity, L f arm the length of the farm, and H abl the height of the ABL. We decide to neglect the width of the farm, as having fixed the aspect ratio in the two farms simulated, L f arm is sufficient to describe their geometry change.
Casting the problem in the context of the Buckingham-Vaschy π-theorem [23,24], there are five total variables, and two fundamental dimensions, length and time, therefore three π parameters must be defined.
The choice of the first two groups appears as the most convenient, the one for the last group deserves a further comment. From our understanding, the wind can adapt to the presence of the farm by deflecting upward, downward, and to the sides of the area occupied by the rotors. As any deflection to the sides comes virtually with no consequence, the deflection up and down can cause extra shear as the flow is bounded by the ground and the capping inversion. The parameter π h is therefore representing the height of the channels upwards and downwards the wind farm rotors, whose cumulative height is H abl − D. In terms of the π-theorem, we have incorporated the variable D (the rotor diameter) only here to avoid an additional π-term. This choice, taken for sake of simplicity, is also justified by the fact that D is constant across the LES, therefore we could have not investigated the dependence on this parameter.
According to the new non-dimensional π-terms, Eq. (7) can be rewritten as, Mathematically, we now propose that F factorises to F (π x , π h ) = G (π x ) g(π h ). As the LES results provide us the trend of π u (π x ), we search for the functional relation g(π h ), so that with G(π x ) being the universal curve describing the velocity deficit. As demonstrated in the right figure of (6), choosing gains the four curves predicted by the LES, Fig. 6 (left), to collapse to one in the inflow region. It must be underlined that this is only a first and early attempt at characterising the interactions of wind farm and ABL we believe to originate the GBE. It is interesting that this simple application of the π-theorem provided a rather promising result. Lacking a more substantial set of observations to test the model on, we avoid trying to interpret to much out of it. But we believe the whole modelling of global blockage would benefit from more work in this direction.

Conclusion
By means of comparison with LES, the analytic induction models are shown to correctly describe the induction generated by a single wind turbine in isolation and to improve the flow representation within the wind farm, especially if the model accounts for lateral speed-ups. Under the conditions simulated, these results do not carry over when these models are linearly superimposed to describe the wind farm inflow deceleration due to global blockage. We observe a substantial underestimation in the model predictions that increases as the wind farm dimension grows and the ABL height decreases. These results show the poor sensitivity, or complete lack thereof, to these parameters to generate large inaccuracies in the models predictions.
We deduce that the deceleration of the inflow is mainly determined by the interactions and the physics of the ABL rather than by the turbines' induction. Compelling evidence at support of this thesis are also found in the production data. We observe the inflow deficit to be partially recovered at first row and especially at the row downstream in the farm. Suggesting the global blockage effect to be not correlated to the turbines induction but rather to exogenous pressure perturbations, likely due to the displacement of air masses with different density aloft [4], capable of redistributing the available energy. We conclude that methods to predict the GBE at wind farm based on reducing the whole farm speed velocity according to the measured inflow deficit, e.g. as in Eq. (5), are not advisable. We also believe the linear superposition of induction models to be inadequate to the modeling of GBE in its present form where no dynamics of the ABL are considered.
Eventually, we perform a dimensional analysis of the inlet velocity deficit predicted by the LES aimed at simplifying the description of the interactions between wind farm and atmospheric boundary layer. Despite being very early, the results are promising and they could give birth to a new approach for the modeling of the GBE.
Future work is necessary to more broadly characterize the extent of the GBE at varying conditions of the ABL, and the wind farm size and shape. This will gain a better testing ground for the presented dimensional model, and nourish further development of the approach. In parallel, we foresee the modeling of the lower atmosphere's vertical structure, such the one presented in [6], as a possible step forward of the current status.