LES Modelling of Highly Transient Wind Speed Ramps in Wind Farms

Large eddy simulations of wind farms are often performed with canonical atmospheric conditions, where the background flow is based on precursor simulations with idealized model setups yielding statistically stationary turbulent flows. However, precursor simulations can only handle gradually changing flow conditions, and are not capable of modelling highly transient and statistically non-stationary flows, e.g. frontal passages or large gusts. Such flows frequently occur in nature, and can influence the operation—and potentially design—of wind turbines. It is generally not possible to impose non-stationary features through inlet boundary conditions, if the imposed flow violates the most fundamental assumption of micro-scale flow simulations, namely conservation of mass. This work presents a method for modelling highly transient wind speed ramps by extending and adapting the method of applying body forces to achieve specific flow scenarios, where the wind speed ramps are embedded as constrained turbulent boxes. Several scenarios with significant increases in the streamwise wind speed are simulated. Analyses of the transient wake dynamics, as the wind speed ramps propagate through large wind farms are performed to show how well the momentum is maintained throughout the numerical domain and the influence and operation of turbines during the ramp passages.


Introduction
Large Eddy Simulations (LES) have been used extensively to model large wind farms e.g. [1,2,3,4], and such simulations are generally performed with statistically stationary background flow conditions. Precursor simulations can generate such statistically stationary turbulent flow scenarios based on simplifying assumptions, e.g. wall-modelling based on Monin-Obukhov similarity theory with the assumptions of homogeneous terrain and constant geostrophic wind. The resulting atmospheric flow exhibits a fully developed balance between the driving geostrophic wind and the shear stress on the ground. The flow will be turbulent and contain variations and localized gusts, but temporal and spatial changes in the flow will be limited.
However, atmospheric flows are typically not in balance as they continuously develop, and scientific focus is increasing on describing the experimental statistics involved with non-canonical scenarios, including extreme wind shear events [5], storms [6], and/or frontal passages, which can give rise to large increases in wind speeds and large directional changes [7]. Yet changes, e.g. [8]. Highly transient events can significantly influence the performance of large wind farms, as recently shown using LES to model large directional changes [9]. Additionally, the extreme events can also be design driving for the turbines [10]. However, previous studies have predominantly focused on the influence of extreme events on single turbines with constrained turbulence simulations [11,12], and not on the aerodynamic developments as highly transient events propagate through wind farms and potentially are altered through the highly dynamic wake interaction.
This work aims to remedy these shortcoming by presenting a method for modeling highly transient events. Here, the focus is particularly on wind speed ramps with rapid wind speed changes, although the proposed methodology can be applied for other extreme events as well.

Flow Solver
The flow is simulated using EllipSys3D, which a 3D incompressible Navier-Stokes solver [13,14]. The governing equations are solved using a finite volume formulation in general curvilinear coordinates on collocated grids in a multigrid and multiblock framework. Turbulence closure is achieved by employing Large Eddy Simulations (LES), which solves the governing equations for the largest scales in time and space, but model the smallest scales through a subgrid-scale model. The mixed length scale model of Phuoc et al. [15] is used in the present work.

Turbine Modelling
The turbines are modelled using the actuator line method, which imposed body forces directly into the flow solver along rotating lines [16]. The actuator lines are here directly coupled to Flex5, which is an aeroelastic tool and computes deflections and forces on the blades based on extracted velocities along the rotating lines, while forces and deflections are transferred from Flex5 to the flow. The aeroelastic coupling also includes the Flex5 controller, which is essential as the operational regime will change during the ramp and hence realistic control modelling is paramount. See [17,18] for further details on Flex5 and [19] for details on the coupling between EllipSys3D and Flex5.

Ramp Modelling
Stochastic turbulence can be generated using the Mann model [20,21], where the spectral tensor is modelled using rapid distortion theory. The Mann model generates turbulent velocity fields, which are Gaussian, anisotropic, homogeneous, yet captures second order statistics for a neutral atmosphere. However, ramps are not homogeneous, so the generated turbulent velocity field has to be corrected by imposing an underlying ramp, where the streamwise velocity changes significantly as done by Hannesdóttir et al. [7]. Alternatively, turbulent flowfields could be generated using precursor simulations of pre-ramp and post-ramp conditions, and subsequently combined.

Modelling Challenges and Methodology
The generated turbulence fields of highly transient events can not be imposed as an inflow boundary condition. It conflicts with the most fundamental concept of fluid dynamics, namely conservation of mass. Any change in mass flux on the inflow boundary has to propagate instantaneously throughout the numerical domain, to fulfill conservation of mass. This will result in an nonphysical acceleration throughout the domain. This effect will also occur in compressible flow solvers as the increased mass flow can not, and should not, be fully absorbed in density changes. An alternative could be to make adaptive boundary conditions in the flow direction, but it would require very large domains in order not to enforce the local flow 3 development in the area of interest. Mesoscale models (e.g. WRF) can cover large areas, which enables the modelling of transient weather events, such as thunderstorms [22], but it is not capable of resolving the smaller turbulent scales governing the wake aerodynamics within wind farms due to computational constraints.
Hence, a different methodology is required. The turbulent velocity fields can be imposed in the numerical domain using body forces, see [23,24] and recently, gusts have been modelled using a similar approach [23]. Typically, the forces are applied in a plane upstream the object of interest, in this case wind turbines operating in a wind farms.
Here, the method of imposing the constrained turbulent velocity fields within the domain is applied, but extended in order to address two significant challenges involved when modelling highly transient events using this approach: Additional steps beyond the standard approach are required to reduced the implications of numerical blockage and the flow degradation, and numerous setups were tested to derive the following remedies.

Flow degradation
The imposed ramp is advected downstream from the plane by the flow solver, which takes over and governs the flow development. Eventually, the imposed highly transient wind speed ramp will dissipate until it is restored to be in balance with the surrounding background flow. However, this flow degradation needs to be minimized in order to model the aerodynamic wake interaction within farms as good as possible over as long distance as necessary.
The required forcing to achieve large rapid changes in velocity is very large and a significant disturbance to the surrounding flow. The surrounding background flow is maintained throughout the simulation and exists outside the ramp region. Therefore, it is beneficial to choose a background flow, which minimize the required forcing relative to the background flow. Hence, the background flow is chosen to be the average of the wind speed before the ramp and after the ramp. Thereby, the required forcing is divided approximately evenly to perturb the background flow in the ramp region to achieve both pre-ramp and post-ramp conditions, so the forcing is initially reducing the flow compared to the background wind speed (U 0 ) and subsequently accelerating the flow from the background wind speed (U 0 ) in the ramp zone, respectively. Figure 1 shows an example of the wind speed (U ) with an imposed ramp as well as the mean background wind speed (U 0 ). Additionally, the degradation process will be enhanced as the large forcing required to impose the wind speed ramp will form a sharp and artificial shear layer relative to the surrounding background flow, where large shear stresses will enhance mixing. The artificial shear layer can be reduced significantly by gradually tapering the forcing to zero outside of the imposed planes in a buffering region. Here, the forcing is tapered off following a half-normal distribution in both lateral and vertical direction. The streamwise force distribution is shown in Figure 2 for both the full domain (left panel) and in the vicinity around the turbine (right panel), where the rotor area is given for reference.

Numerical blockage
Numerical blockage occurs if the numerical domain is not sufficiently large, so the presence of an obstruction, for instance a cylinder or a wind turbine, will cause an acceleration of the flow around the obstacle. Here, the imposed ramp plane upstream acts as a virtual obstacle, essentially behaving like an entire actuator plane with large body forces. The ramp plane needs to be significantly larger than the turbines, to properly mimic the entrainment of turbulent kinetic energy from above, see [25], but the blockage ratio, i.e. size of the turbulence plane relative to the size of the numerical domain, should also be kept low. Here, the extend of the entire imposed ramp plane is 13.5R × 7.25R and the central region with full forcing is 8R × 5.25R as shown in Figure 2. This corresponds to a blockage ratio of approximately 6.1% and 2.6%, respectively, which roughly match the general recommendations of Baetke et al. [26]) of a blockage ratio of 3%.

Methodology Summary
The methodology with recommendations is summarizes as follows:   Table 1: Summary of simulated scenarios in terms of pre-ramp wind speed (U * pre−ramp ), ramp velocity increase (∆U * ramp ), mean background velocity (U 0 ), pre-ramp wind shear (dU/dz pre−ramp ), rise time (∆t), and turbulent length scale (L) for Mann box generation.
The following will assess the simulated flows, but the suggested criteria here can obviously be changed to simulate even more challenging flow scenarios or larger wind farms.

Simulation Setup
The numerical domain for the simulations has a total size of 136R×14R×100R in the streamwise (X), lateral (Y ), and vertical direction (Z). Periodic boundary conditions are applied on the lateral sides, no-slip on bottom, symmetry on top boundary, and inlet and outlet on the streamwise boundaries. The grid is stretched towards the boundaries, and the grid points are equidistantly spaced in the vicinity of the turbines from Y = 4R − 10R and Z = 0R − 5R corresponding to approximately 10 cells per turbine blade. The same grid resolution is used in the streamwise direction from X = 5R −130R. This is quite coarse for actuator line simulations, but it is expected to only give a difference of approximately 1% in C T [27], which governs the wake dynamics of interest here for assessing the methodolgy. The mesh consists of a total of 6 · 10 6 mesh points.
Eight different ramp scenarios have been modelled as an ensemble of various parameters governing ramp events, namely pre-ramp wind speed (U * pre−ramp ), ramp velocity increase (∆U ramp ), mean background velocity at hub height (U 0 ), background wind shear (dU/dz pre−ramp ), rise time (∆t), and turbulent length scale (L) for Mann box generation. The scenarios are summarized in Table 1. Note, the values reported here corresponds to the resulting values of U * pre−ramp and ∆U * ramp calculated as the mean velocity at hub height 2R upstream the first turbine for t = [100s; 700s] for pre-ramp values and for t = [900s; 1600s] for postramp. These values are somewhat different than the design values described in [28,29], which performed the statistical analyses of these ramps. The required forcing to attain a certain level can be achieved through additional calibration.
A total of 9 turbines have been modelled in a row for the different scenarios, where the turbines are spaced 14R apart. The first turbine is located in (X, Y, Z) = (10R, 7R, 2R). The modelled turbine is an uprated version of NM80, which has a rated power of 2750M W and rated velocity U = 14m/s as well as a radius of R = 40m [30].

Flow Description
Simulating the ramps with the presence of the wind turbines affect the flow as the turbines extract power. Figure 3 shows an instantaneous snapshot of the streamwise velocity contour through a row of nine wind turbines, which indicate how the flow behaves as it propagate through the wind farm. The bottom plot shows the streamwise wind speed at hub height, which is clearly different for the first four turbines compared to the last four, while turbine five is experiencing the ramp passage at this instance. The wake regime is also clearly different before and after the passage of the ramp as the turbine operations transitions from below rated wind speed to above rated wind speed, which show significantly lower wake deficits, as expected.

Assessment of Flow Degradation
It is important to assess how well the ramps are simulated throughout the domain, so the flow degradation is assessed by computing the streamwise momentum flow rate through planes at different locations in the numerical domain in the absence of wind turbines for scenario 1. The local streamwise momentum flow rate is computed aṡ where X 1 , X 2 , Y 1 , and Y 2 are integration limits, ρ is density, and U is the streamwise velocity. The momentum flow rate has been computed for square areas of ±1R and ±2R in both directions around the wind turbines to show the sensitivity of the area size, see green boxes shown in Figure  2. Figure 4 show the normalized momentum flow rate as function of time and downstream distance for the two different areas shown in green in Figure 2, where it is normalized by the initial momentum flow rate throughout the domain before the forcing is imposed to show the reduction (blue) and increase (red) in momentum flow rate pre-and post-ramp, respectively. The ramp is initiated at t = 800s and the propagation front of the ramp is clearly seen in white. The momentum flow rate is computed instantaneously, so there are some local changes in time and space as larger turbulent structures propagate downstream. The level is generally quite constant and within the same contour levels of ±5%, which verifies the proposed methodology as the flow degradation is small. The difference between the two areas show that the sensitivity of the momentum flow rate is higher for the smaller areas, whereas it is more constant for the large area. The constant level of momentum flow rate is essential to ensure the correct entrainment once the turbines are introduced, so the imposed momentum is not drained by the turbines, but that there is sufficient momentum available in the vicinity of the turbines to enable proper wake recovery.

Momentum Analysis
The momentum flow rate is also computed for the eight wind farm scenarios, and shown in Figure 5. Overall, the same trends are seen with decreases in forcing prior to the ramp and increases afterwards. However, streaks now appear in both the blue pre-ramp wind speeds as well as in the red post-ramp regions, and clearly shows how the presence of the turbines extract momentum from the flow. The turbines operate below rated before the passage of the ramp and hence impose a relatively large thrust force, i.e. thrust coefficient. Conversely, the relative influence of the turbines is reduced after the ramp passage as the resulting wind speed is above rated for several of the scenarios, while the post-ramp influence of the turbines is significant in scenarios 2, 6, 7, and 8 as the resulting wind speed inside the wind farm is again below rated wind speed, i.e. high thrust coefficient. Figure 6 shows the difference between the momentum of scenario 1 without and with the turbines, i.e. difference between Figure 4b) and Figure 5a), and normalized by the initial momentum of the scenario without turbines. The normalized difference shows how the momentum through the domain is always larger without the turbines, and mostly in the order of 0 − 30% larger. The difference is clearly largest during the ramp passage due to the presence of the turbines, which clearly corresponds to significant changes in the flow.
The total thrust forcing of the nine turbines is opposing the momentum, but the question is whether the spatial change is dominated by the general flow degradation or the turbine forcing, i.e. the physically correct wake aerodynamics as turbine extract power. Figure 7 shows the average momentum flow rate for scenario 1 with and without turbines at the nine turbine locations as well as the average thrust force exerted by the turbines. The average has been performed for t = [350s; 750] (pre-ramp) and t = [1150s; 1550] (post-ramp) as indicated by the black boxes in Figure 4. The average thrust force of the individual turbines are continuously cumulated to provide direct comparison of the magnitude of turbine forcing relative to the momentum throughout the domain. The momentum decrease for both pre-ramp and post-ramp   Figure 6: Normalized difference in momentum flow rate for scenario 1 with and without wind turbines, see Figure 4b) and Figure 5a).
due to the presence of the turbines, but the turbine forcing corresponds to Tn ∆Ṁ = 104% of the momentum decrease during the pre-ramp, while it is only Tn ∆Ṁ = 58% for the post-ramp. This indicates that the turbines dominate the flow degradation for below rated operation, which is physically correct, although additional momentum will also be entrained from the surroundings. Part of the imposed momentum is essentially dissipated to the surrounding background flow for above rated conditions, i.e. a negative entrainment from the turbine region to the surrounding as the wake deficits are small, but here the effect is assessed on average after the ramp passage, and deemed minor during the ramp passage itself. Note, that the degradation could be reduced further if required for simulating for instance larger wind farms or different transient events, e.g. large changes in wind direction or to model specific flow scenarios observed in nature [31]. The degradation can be decreased directly by increasing both the size of the numerical domain and the size of the imposed ramp region.  Figure 8 shows the power production of the nine turbines during the ramp passage of scenario 1. The wake effects clearly reduce the power production of the downstream turbines before the ramp passage. The ramp only have minor influence on the power production of the first turbine, but yields significant increases for the downstream turbines as they reach rated production. The ramp propagation is clearly seen in the power signal, although additional dynamics appear to be introduced, particular for turbines 5, 6, 7, and 9, which show a intermediate peak in power during the ramp passage. For more detailed analyses of the operation and aeroelastic loads on wind turbines operating in wind farms, please refer to [28,29].

Conclusion and Outlook
Highly transient flow scenarios are inherently difficult to model in LES, as it violates mass conservation. A methodology is presented here, which enables the simulation of such highly 10 transient events to study the wake aerodynamic flow and wind farm operation. Here, the simulated cases focus on large sudden increases in velocity, i.e. ramps, but the methodology can be applied for other types of extreme scenarios, such as large wind directional changes or storms. These scenarios will generally be influenced by the atmospheric stability, and the ability of the methodology to also mimic rapid transitions in atmospheric stability requires further investigation. Here, the methodology is assessed in terms of the ability to maintain the momentum flow throughout the domain in the vicinity of the operating wind turbines.
Generally, the verification show that the momentum is maintained well throughout the domain to mimic the ramp propagation, although the presence of the turbines counteract the imposed ramps. For below rated, this is physically correct, while part of the flow degradation for above rated ramp velocities appear to be caused by increased mixing with the surrounding. However, the ramp passage itself is expected to be less influenced by this. It should be noted that the general flow degradation can be reduced to achieve a desired level for a given scenario by increasing the size of both the numerical domain and the forcing region. The results also indicate significant changes to the wake aerodynamics during ramp passages, which requires further investigations as they can have significant impact on the turbine operation and be design driving.