{'fraction_non_alphanumeric': 0.035266191929467616, 'fraction_numerical': 0.020881297852974247, 'mean_word_length': 4.360887868350555, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Two-scale interaction of wake and blockage effects in large wind farms

Turbine wake and farm blockage effects may significantly impact the power produced by large wind farms. In this study, we perform Large-Eddy Simulations (LES) of 50 infinitely large offshore wind farms with different turbine layouts and wind directions. The LES results are combined with the two-scale momentum theory (Nishino&Dunstan 2020, J. Fluid Mech. 894, A2) to investigate the aerodynamic performance of large but finite-sized farms as well. The power of infinitely large farms is found to be a strong function of the array density, whereas the power of large finite-sized farms depends on both the array density and turbine layout. An analytical model derived from the two-scale momentum theory predicts the impact of array density very well for all 50 farms investigated and can therefore be used as an upper limit to farm performance. We also propose a new method to quantify turbine-scale losses (due to turbine-wake interactions) and farm-scale losses (due to the reduction of farm-average wind speed). They both depend on the strength of atmospheric response to the farm, and our results suggest that, for large offshore wind farms, the farm-scale losses are typically more than twice as large as the turbine-scale losses. This is found to be due to a two-scale interaction between turbine wake and farm induction effects, explaining why the impact of turbine layout on farm power varies with the strength of atmospheric response.


Introduction
The global wind energy production is set to rise in the next few decades. To achieve this, wind farm clusters are expected to be built which are an order of magnitude larger than existing ones (Maas & Raasch 2022). To optimise their design, it is important to accurately predict the total farm power as well as aerodynamic loads on each turbine. However, it is very difficult to model the aerodynamics of wind farms because of the multi-scale nature of the wind farm flows (Porté-Agel et al. 2020). In 2019 Ørsted, one of the largest offshore wind farm developers, announced that its wind farms were producing less power than expected (Ørsted 2019). The underprediction of farm power was attributed to two effects: wake and farm blockage effects.
Behind every turbine there is a turbulent wake which has a reduced wind speed. When turbine wakes interact with other turbines within a farm this can cause substantial power losses. This is known as the 'wake effect' and has been measured to reduce the power of turbines in an existing wind farm by up to 40% in the worst case (Barthelmie et al. 2010). This effect has been extensively investigated in the literature using Large-Eddy Simulations (LES) (Porté-Agel et al. 2013;Wu & Porté-Agel 2015;Stevens et al. 2016a). The 'farm blockage' is a recently observed effect in large wind farms (Bleeg et al. 2018). The flow resistance caused by a wind farm reduces the wind speed upstream as well inside the farm. Hence, the total power produced by the wind farm is reduced compared to the ideal situation where the upstream wind speed is not affected by the farm (Nishino & Dunstan 2020).
Wind farm aerodynamics have been traditionally modelled using 'wake' models, which predict the velocity behind a turbine (e.g., Jensen 1983; Bastankhah & Porté-Agel 2014). To account for interactions between multiple turbines, the wake velocity deficits are superposed (e.g., Katic et al. 1986;Zong & Porté-Agel 2020). However, these models do not account for the response of the atmospheric boundary layer (ABL) to the farm. As such they tend to perform poorly for large wind farms (Stevens et al. 2016a). A different approach to modelling wind farms is to use 'top-down' models (e.g., Frandsen 1992;Frandsen et al. 2006;Calaf et al. 2010). They model the response of an idealised ABL (which follows a logarithmic law) to an infinitely large wind farm but cannot take into account the details of turbine-wake interactions explicitly. Hence, these wind farm models cannot correctly capture the two-way interactions of turbine-scale and farm-scale flow effects which determine the performance of large wind farms.
To capture this two-way interaction, it has been proposed to couple wake and top-down models; for example, by adjusting parameters in both models to match the hub-heightaveraged velocities (Stevens et al. 2016b;Starke et al. 2021). Examples of two-way coupling can also be found in existing engineering software, e.g., the 'Deep Array Wake Model' (Brower & Robinson 2012). However, these models involve the coupling of low-order flow models and are therefore limited by the assumptions made by the constituent models, e.g., wake superposition or a log-law wind profile. To account for the effects of more realistic flow physics, it would be beneficial to use an approach based on more fundamental laws of fluid mechanics.
The optimal design of a large wind farm under realistic atmospheric conditions remains a challenge as it requires consideration of the complex atmospheric response to the farm. It is often too expensive to run a large number of simulations which resolve both individual turbines and the atmospheric response to the farm. In addition, to find an optimal design for the long-term performance of a wind farm, such as the annual energy production (AEP), the range of timescales we would need to consider is too wide. As such, Nishino & Dunstan (2020) proposed the 'two-scale momentum theory' to split the multi-scale flow problem into 'internal' turbine/array-scale and 'external' farm/atmospheric-scale problems. West & Lele (2020) performed LES of infinitely-large wind farms and the results showed a good agreement with the two-scale momentum theory. However, their LES study was limited to 'fully aligned' turbine layouts, and their discussion was also limited to the special case where the momentum supplied by the atmosphere to the wind farm site was fixed. In reality, the strength of atmospheric response to the wind farm resistance depends on mesoscale weather patterns (Patel et al. 2021) as well as atmospheric stability and gravity waves (Allaerts & Meyers 2017, 2018, 2019.
The aim of the present study is to better understand the fluid mechanics processes which determine the power production of large wind farms, using a combined theoretical and computational approach. First, we will perform a large suite of LES of infinitely large wind farms with different turbine layouts and wind directions. We then combine the results of LES with the two-scale momentum theory to investigate and explain expected performance of large finite-sized wind farms with a realistic range of atmospheric response strengths. Using this approach allows the combined effects of turbine-scale and farm-scale flow characteristics on wind farm power to be determined.
In section 2 we summarise the definitions of key wind farm parameters in the two-scale momentum theory (Nishino & Dunstan 2020). Section 3 details the methodology of the LES and wind turbine implementation. In section 4 we present the results including validation of the LES code. These results are discussed in section 5 and concluding remarks are given in section 6. Figure 1 shows a pair of control volumes for a given farm site. We first consider the momentum balance of the control volume without the turbines (figure 1a). For this scenario the following equation can be derived for a 'short-time-averaged' flow:

Two-scale momentum theory
is the velocity in the hub-height wind direction (i.e. streamwise direction) , is the net streamwise momentum injection through the top and side boundaries of the control volume (due to advection and Reynolds stress), is the streamwise component of the Coriolis force averaged over the control volume, / is the pressure gradient in the direction , is the bottom shear stress, is the wind farm area and is the volume of the control volume. The subscript 0 refers to values without the turbines present, [ ] refers to control-volume-averaged values and refers to farm-area-averaged values. By considering the control volume with the turbines present (figure 1b) the following equation can be derived: where is the thrust of turbine in the farm and is the total number of turbines in the farm. Equations 2.1 and 2.2 can be combined to obtain the non-dimensional farm momentum where is the farm wind-speed reduction factor defined as ≡ / 0 (with defined as the average wind speed in the nominal farm-layer of height , and 0 is the farmlayer-averaged speed without the presence of the turbines); is the array density defined as ≡ / (where is the rotor swept area); * is the (farm-averaged) 'local' or 'internal' turbine thrust coefficient defined as * ≡ =1 / 1 2 2 ; 0 is the natural friction coefficient of the surface defined as 0 ≡ 0 / 1 2 2 0 ; is the bottom friction exponent defined as ≡ log ( / 0 ); and is the momentum availability factor given by: . (2.4) Note that the first term in the left-hand-side of equation 2.3 can be extended to include the impact of support structure drag (see Ma et al. 2019) but their impact is usually small and is therefore neglected in this study for simplicity. The height of the farm-layer, , is used to define the reference velocities and 0 . is typically between 2 ℎ and 3 ℎ (where ℎ is the turbine hub-height); the NDFM equation (2.3) is valid as long as the same value is used for both 'internal' and 'external' problems. In this study we use a fixed definition of = 2.5 ℎ . This is discussed further in appendix A.
The equation 2.3 helps the analysis of large wind farm aerodynamics. This is because * and on the left hand side are expected to depend primarily on the turbine/array-scale flow physics or 'internal' conditions, for example the turbine layout, operating conditions and local wind conditions, whereas on the right hand side is expected to depend largely on 'external' conditions. Following Nishino & Dunstan (2020), in this study we assume that the 'internal' problem (to be modelled using LES in section 3 to calculate * and ) can be modelled without explicitly considering the effects of 'external' conditions such as wind farm size and location, and the response strength of the atmosphere.
The 'external' problem is to determine the parameter , which represents how much the amount of momentum available to the farm site differs from its 'natural' value. This problem is largely independent of the small-scale flow features and can be modelled using a numerical weather prediction (NWP) model with a wind farm parameterisation, i.e., without resolving individual turbines. Patel et al. (2021) used such an NWP model to demonstrate that, for most cases, varied almost linearly with (for a realistic range of between 1 and 0.8). Therefore, can be approximated by (2.5) where is called the 'momentum response' factor or 'wind extractability' factor. Patel et al. (2021) found to vary between 5 and 25 for a typical offshore wind farm site. Note that = 0 corresponds to the case where momentum available to the farm site is assumed to be fixed, i.e., = 1. An infinitely-large value of corresponds to the case where there is Focus on Fluids articles must not exceed this page length no wind speed reduction in the farm-layer, i.e., = 1. Preliminary results of an extended study from Patel et al. (2021) show that changes according to atmospheric conditions and decreases exponentially with increasing farm size (see appendix B). Although the details of how changes with weather conditions are still unclear and need to be clarified in future studies, in the present study we take 0 < < 25 as a typical range of the wind extractability for large offshore wind farms. Using obtained from equation (2.3) for a set of * , and , the following equation can be used to calculate the power coefficient of the turbines within the farm, where is the (farm-averaged) turbine power coefficient defined as ≡ =1 / 1 2 3 0 ( is power of turbine in the farm) and * is the (farm-averaged) 'local' or 'internal' turbine power coefficient defined as * ≡ =1 / 1 2 3 .

Analytical model of ideal wind farm performance
In this study we consider arrays of actuator discs (or aerodynamically ideal turbines operating below the rated wind speed). For an actuator disc * = * where is the turbinescale wind speed reduction factor defined as ≡ / ( is the streamwise velocity averaged over the rotor swept area). We can estimate using the expression = √︁ * / where ≡ / 1 2 2 is a turbine resistance coefficient describing the turbine operating conditions (noting that this is strictly valid only for infinitely large regular arrays of turbines where the farm-averaged turbine thrust is identical to the thrust of each individual turbine). The theoretical of an actuator disc is therefore given by where * may be predicted using a simple analytical model (Nishino 2016) given by * = 4 (1 − ) = 16 (4 + ) 2 (2.8) using the expression = * / 2 to express * as a function of . The model predicts * as a function of turbine-scale wind-speed reduction by using an analogy to the classical actuator disc theory. This simple analytical model will be compared with LES results later in section 4. Using the analytical model of * (equation 2.8), and the linear approximation of (equation 2.5), equations 2.3 and 2.7 can be solved to give a theoretical prediction of , which we will call , ℎ . Note that West & Lele (2020) also introduced this but only for the special case with = 0. As shown by Nishino & Dunstan (2020), , ℎ is sensitive to but much less sensitive to . If we assume that = 2.0, then we can obtain an analytical expression for , ℎ , i.e., (2.9) It is worth noting that the power coefficient of an isolated turbine, , , is given by (2.10) which takes the well-known maximum value of 16/27 at = 2. This equation can be obtained by substituting 2.8 into 2.7 with = 1 (i.e., assuming flow mechanisms as described by the classical actuator disc theory and no farm-scale wind speed reduction). Note that this is the same as solving equation 2.9 for two special cases: (1) with / 0 = 0; and (2) with an infinitely large value of .

Governing equations of the flow
We performed LES of flow over periodic turbine arrays using the MetOffice/NERC Cloud (MONC) Model (Brown et al. 2018). The flow is driven by an imposed pressure gradient and is neutrally stratified. The flow is governed by the incompressible Navier-Stokes equations, i.e., = 0 (3.1) where is the resolved velocity in the direction, is the pressure perturbation from the reference state, is the reference density, is the subgrid stress term, is the force added to model the wind turbines and ∞ / is the imposed pressure gradient. The subgrid stress model is a standard Smagorinsky model (Smagorinsky 1963) given by = where is the subgrid-scale eddy viscosity and is the rate of strain tensor. The eddy viscosity is given by a mixing length model = 2 where is the mixing length scale and is the modulus of the rate of strain tensor = / √ 2. Near the bottom boundary, the mixing length scale is damped using the function 1/ 2 = 1/( 0 ) 2 + 1/[ ( + 0 )] 2 described in Brown et al. (1994) where 0 is the basic mixing length scale. 0 is given by 0 = Δ with a coefficient of = 0.23 and a grid spacing given by Δ = max(Δ , Δ ) (Brown et al. 1994).
All velocity components are set to zero at the bottom boundary. The shear stress at the surface is parameterised by specifying following the classical Monin-Obukhov similarity theory. The horizontal boundary conditions are periodic for all prognostic quantities. The top boundary has a zero vertical velocity boundary condition and a damping layer for the top 200 metres of the domain (which was not necessary in the present study for neutrally stratified flows but still included for future studies to explore the effect of atmospheric stability).

Actuator disc implementation
We model individual turbines as actuator discs following the methodology used by the KULeuven code described in Calaf et al. (2010). The approach uses a Gaussian convolution filter to apply the turbine force from the rotor plane onto the LES grid. This allows the position and orientation of turbines to be changed easily. The thrust force exerted by a single turbine is given by where is the time-filtered disc-averaged velocity and is the turbine diameter. This turbine thrust force is spatially distributed using a normalised indicator function R ( ), defined as where ( ) is a filtering kernel. This integral is calculated over the surface of the disc. We divide the disc area into 10 segments in the radial and angular directions, respectively. This was sufficient for R ( ) to be independent of the number of segments. MONC uses a staggered grid where the and velocities are evaluated at different points. As such, two different indicator functions, R ( ) and R ( ), are calculated for the and directions. We use the same filtering kernel as described in Shapiro et al. (2019), where is the filter width, which following the approach of Shapiro et al. (2019) is given by = 1.5 √︁ Δ 2 + Δ 2 + Δ 2 . The force per unit density at a given grid point in the direction is given by The disc-averaged turbine velocity is calculated using the indicator function R ( ) as a weighting function, where is the wind direction relative to the direction. We use a constant value for which is the direction of the pressure gradient forcing. Note that refers to the velocity in the direction whereas describes velocities in the wind direction. The spatially-averaged velocity is then temporally averaged using a one-sided exponential time filter with a time window of 10 minutes to calculate . To calculate * from the LES we use the following relationship, * = 2 ∑︁ =1 2 (3.8) noting that turbine velocity is time filtered before being squared and then averaged over all discs.
is calculated by integrating the streamwise velocity ( cos + sin ) between the surface and 2.5 ℎ across the entire domain. Unlike , no time filter is used to calculate . * varies with time during the LES so is time-averaged over a long period to give a single value of * .
To calculate the (farm-averaged) turbine power coefficient from the LES we use the expression, where 0 is the farm-layer-averaged velocity in an LES without turbines. varies with time so is calculated by time averaging over a long period.

LES code validation
We firstly validate our LES framework with the new actuator disc implementation by comparing with the benchmark cases reported in Calaf et al. (2010). We then investigate the sensitivity of our results to horizontal resolution, domain size and pressure solver.
For the validation cases summarised in table 1 we use a surface roughness length of 0 = 0.1m and a pressure gradient of (1/ ) ∞ / = 1 × 10 −3 m/s 2 . The turbines all have a hub height of 100m and a diameter of 100m. We use a turbine resistance of = 1.33 and the same turbine spacing as for Case A1 in Calaf et al. (2010) ( = 7.85 and = 5.23 ). The surface roughness length, pressure gradient, turbine design and resistance are chosen to match the values used by Calaf et al. (2010). Validation cases V-1, V-2 and V-3 use a FFT pressure solver whereas V-4 uses an iterative pressure solver. Cases V-1, V-3, V-4 have a domain size × × of 3.14 × 3.14 × 1 km (with 24 turbines) and case V-2 has a domain size of 6.28 × 6.28 × 1 km (with 96 turbines). All validation cases were run for 100,000 seconds and flow data averaged between = 30,000 and 100,000 seconds. The convergence of two flow statistics for case V-1 are shown in figure 2. Figure 3 shows the instantaneous streamwise velocity plotted on a cross-streamwise plane 2.5 behind a row of turbines in the validation case with a double horizontal resolution. Figure 4a shows the time and horizontally averaged streamwise velocity for the four validation cases. The velocity profiles agree well with the results in Calaf et al. (2010). Figure 4b shows the profiles of total shear stress and the results reported in Calaf et al. (2010). The shear stress profiles match well except for the region near the bottom surface. Note that the shear stress in our LES does not approach zero at the bottom surface as it includes both modelled and resolved components; the latter of which approaches zero as in the results of Calaf et al. (2010). The velocity profiles are insensitive to the horizontal domain size and the pressure solver.
In wind farm LES using the actuator disc method, the disc-averaged velocity is usually overpredicted (Shapiro et al. 2019). This is because at coarse resolutions, the vorticity shed from the disc edge is not fully captured. This can be seen in our results by considering cases V-1 and V-3 in figure 4a. In the coarse grid case V-1, the disc-averaged velocity is overpredicted so the turbine thrust applied is greater. This results in a slightly lower throughout the entire domain. This effect is also seen in the higher uncorrected * value in case V-1 compared to V-3 (table 1), which is due to the overprediction of the disc-averaged velocity.
To correct the overprediction of disc velocity in coarse LES the following correction factor was proposed in Shapiro et al. (2019), (4.1) We apply this correction factor to the disc-averaged velocity by multiplying our uncorrected * values by 2 . The correction is applied after the simulation and not during. After correction, the horizontal resolution used in cases V-1 and V-3 only had a small impact on * , suggesting that this correction factor can be successfully applied to a periodic array of actuator discs. For the value of used here the analytical model of * in (2.8) gives a  Case Δ / Δ / Δ / × × (km) Pressure solver Uncorrected * * V-1 0.245 0.245 0.0787 3.14 × 3.14 × 1 FFT 0.8517 0.6845 V-2 0.245 0.245 0.0787 6.28 × 6.28 × 1 FFT 0.8587 0.6902 V-3 0.1225 0.1225 0.0787 3.14 × 3.14 × 1 FFT 0.7844 0.6957 V-4 0.245 0.245 0.0787 3.14 × 3.14 × 1 Iterative 0.8603 0.6914 value of 0.75. All the validation cases in table 1 have a lower * than this because of wake interactions between turbines. We also consider the effect of resolution on the wake velocity deficit. Figure 5 shows the average wake profiles for each of our validation cases in table 1. The wake velocity profiles are normalised by the farm-averaged velocity for each validation case. This shows the far wake velocity deficit does not vary with the domain size, horizontal resolution and pressure solver. Comparing the wake profiles for V-1 and V-3 at 2 downstream (figure 5a) shows a small difference in the velocity deficit at the centre. This is because V-3 uses a different filter size for the projection of the turbine area (see section 3.2). When comparing the wake profiles 4 and 6 downstream this difference is negligible. This shows that the far wake velocity profile is insensitive to the filter size used for the turbine area projection. 0.75 1.00

LES results
The LES cases for this study have the same setup as the validation case V-4 (see section 4.1) except for the domain size, the surface roughness length and the streamwise pressure gradient. To model offshore wind farms, we use a surface roughness length of 1 × 10 −4 m and a pressure gradient of (1/ ) ∞ / = 8.38 × 10 −5 m/s 2 in the wind direction , which results in 0 = 10.103m/s and 0 = 1.607 × 10 −3 for a fixed nominal farm-layer height of = 250m (both obtained from LES with no turbines). All cases were run for 100,000 seconds and flow data averaged between = 60,000 and 100,000 seconds. Note that we adopted a different spin-up and averaging period (compared to section 4.1) because of the different pressure gradient forcing.
We performed a suite of 50 simulations with different turbine layouts which are described by the parameters , the turbine spacing in the direction, , the turbine spacing in the direction and , the wind direction relative to the direction (see figure 7a). The turbine operating conditions are the same for all simulations and is given by = 1.33. We consider a realistic range of turbine layouts and wind directions: ∈ [5 , 10 ], ∈ [5 , 10 ] and ∈ [0 , 45 ]. We only consider regular arrays and so by symmetry we only need to consider wind directions up to 45 . We adopt the minimum possible horizontal domain size ( and ) within the range between 3.14km and 6.28km (depending on and ) as the validation results presented in section 4.1 suggest that the results would be insensitive to the domain size within this range.
We use a space filling maximin design (Johnson et al. 1990;Santner et al. 2018) to select different turbine layouts in the parameter space ( , , ). The maximin algorithm iteratively selects a point which maximises the minimum distance to other points and to the boundaries of the parameter space. Figure 7b shows the 50 different turbine layouts selected in the parameter space. Figure 8 shows the time-averaged flow fields from 4 of 50 cases. Figure 8a is for a case where the wind direction is almost perfectly aligned with a relatively small streamwise spacing between turbines of 5.76 . This case gives a low * value of 0.585 due to strong wake effects. High speed regions between rows of turbines are formed because of the large cross- streamwise turbine spacings and aligned wind direction. Figure 8b is for a case with a high turbine density and the wind direction almost aligned along the diagonal. This arrangement is similar to a staggered layout. The streamwise spacings between turbines is larger than for figure 8a so the * has a higher value of 0.669 because of the increased wake recovery between turbines. The flow field for a case with an intermediate wind angle is shown in figure 8c. The turbine wakes are mostly misaligned with downstream turbines which minimises wake effects and gives a high * value of 0.752. This result agrees qualitatively with Stevens et al. (2014) in which it was found that the maximum farm power was produced by an intermediate wind direction. The results also give further evidence that the analytical model of * proposed by Nishino (2016) can be used to predict an upper bound to wind farm performance as it gives * = 0.75 in this case. Figure 8d shows the streamwise velocity for a partially waked turbine layout. The partial wake effects cause the * to be reduced slightly to 0.713. Figure 8 also shows the effect of the turbine layout on the farm-averaged wind speed . The farm shown in figure 8a has a low array density and so has a high farm-averaged wind speed of = 0.347 (shown by the brighter colour). Figure 8b shows a farm with a high array density which resulted in a low farm-averaged speed of = 0.248. Figures 8c and 8d Figure 9 shows that * was not a strong function of or . * was found to be a much stronger function of (see figure 10). The lowest * values were for small values of because of the high degree of turbine-wake interactions. When was very small * was also sensitive to the value of (see figure 10a). As increases * increases rapidly until the maximum value around 15 . As increases further * slowly decreases. This is because turbines start to become aligned along the diagonal (similar to the layout shown in figure 8b). When is greater than 15 , the minimum * value tends to be observed when tan( )/ is close to 1 (figure 10b). This corresponds to layouts where turbines are aligned along the diagonal (similar to the layout shown in figure 8b). Figure 11 shows that the range of from the wind farm LES is small (varying between 1.7 and 1.8). Nishino (2016) suggested that would be slightly less than 2 because the presence of turbines would increase the turbulence intensity in the ABL. These results, along with the     Figure 11 shows that there is a slight variation of with wind direction and effective array density.

Prediction of wind farm performance
Now we compare the (farm-averaged) turbine power coefficient from the wind farm LES with the analytical model derived from the two-scale momentum theory, , ℎ (equation 2.9).
To make a fair comparison of between the theory and the LES, we need to consider the fact that the coarse LES resolution caused the turbine thrust to be overpredicted. For * this was corrected for in post-processing using the correction factor (equation 4.1).
However, for there are two simultaneous factors that need to be corrected. The first is that the disc velocity (relative to the farm-layer velocity) has been overpredicted, which increases the turbine power. The second is that the farm-layer velocity has been underpredicted (as shown earlier in figure 4a), which reduces the turbine power. The first effect can be corrected by using the correction factor (i.e., multiplying the raw from the LES by 3 ) but the second effect cannot be corrected in this manner.
To adjust for the second effect we estimate the farm wind-speed reduction that would be obtained if a sufficiently fine resolution was used for the LES,  is calculated by solving equation 2.3 (with = 1) using * from the LES and the corresponding array density / 0 . To calculate , ℎ we do the same but with the uncorrected turbine thrust (i.e., * / 2 ). Since equation 2.3 is derived from the law of momentum conservation, the only assumption we are making in equation 4.2 is that the flow is independent of Reynolds number (which is a reasonable assumption as the change in Reynolds number between the fine and the coarse resolution cases is only of the order of 10%). We can therefore estimate the turbine power from a fine resolution LES, using the expression The green symbols marked by =0 in figure 12 show the corrected for each of the 50 wind farm LES runs. These are plotted against the effective array density / 0 . The blue line shows the prediction using the theoretical model (equation 2.9). The theory matches remarkably well with the LES results, despite that the theory does not account for the turbine layout parameters, namely , and . The reason for this excellent agreement for = 0 will be discussed later in section 5.
Next, we use the results from LES of infinitely large wind farms to estimate the average power coefficients for large but finite-sized wind farms, following the concept of the twoscale momentum theory. We combine the infinite wind farm LES results with the simple linear model of the momentum availability factor (equation 2.5). This does not capture the finite-size effects observed near the edge of the farm but reveals general trends of turbine layout effects with different atmospheric responses. We assume that the farms are still sufficiently large that the flow over the farm is mostly fully developed (or more specifically, the dependency of * on , and is still approximately the same as that in the corresponding infinitely-large farm). 'Sufficiently large' depends on atmospheric conditions (Wu & Porté-Agel 2017) but it is likely to be on the order of 10km.
Firstly, we consider the balance of the pressure gradient forcing (PGF) with surface stress and turbine thrust in our LES, where * is the 'local' or 'internal' friction coefficient * ≡ / 1 2 2 and Δ /Δ is the PGF applied to the LES of an infinitely large wind farm. If we assume * and * are independent of the PGF (as the Reynolds number of the flow is very high), then equation 4.4 shows that 2 should be proportional to Δ /Δ ; hence where is the farm-layer-averaged wind speed in the LES and Δ /Δ is the corresponding PGF, , is the new farm-layer-averaged wind speed that would result from an unknown PGF for a finite-sized wind farm, Δ /Δ . As the flow is quasi-steady and horizontally periodic with no Coriolis force in this study (meaning that we ignore the contributions of all terms except the PGF in equation 2.4), the right hand side of equation 4.5 is equivalent to the momentum availability factor in equation 2.4. This will not be true for real wind farms under real atmospheric conditions. However, this is valid for our simplified analysis where the flow across the farm is mostly fully developed and is driven purely by a PGF. As noted earlier, it was found in Patel et al. (4.6) Since and 0 are both known from LES results, , can be calculated analytically for a given . Finally we assume (again based on the Reynolds number independency) that the disc-averaged velocity scales with and as such , = , 3 . (4.7) where , is the average turbine power coefficient of the large finite-sized wind farm. Figure 12 also shows the average power coefficients estimated for three different values. This shows that the atmospheric response can significantly impact the farm power. for = 25 are roughly an order of magnitude higher than for the same layout with = 0. Note that these results are with a constant turbine operating condition of = 1.33 which may not be optimal for a given . Figure 4 in Nishino & Dunstan (2020) shows how the theoretically optimal turbine power coefficient, , varies with the strength of atmospheric response.
Of particular interest in figure 12 is that, for a given / 0 , the variation of (due to different turbine layouts) increases with . To better understand this trend, the values for the 50 turbine layouts are presented together with their * values in figure 13 for 6 different values separately. Figure 13a shows the results for infinitely large farms with =0. This shows that the average turbine power coefficient is a strong function of the effective array density and insensitive to * meaning that turbine-wake interactions have little effect on the farm power. However, for finite-sized farms (figures 13b-f) the layouts with high wake interactions and a low * value tend to produce less power. These are shown by the darker plots which fall well below the theoretical prediction (blue line). Interestingly, some of the layouts seem to produce slightly higher power than predicted by the two-scale momentum theory. These layouts have a * slightly higher than 0.75 (the value for an isolated turbine) and this seems to be due to locally accelerated flow caused by the local blockage effect (Nishino & Draper 2015;Ouro & Nishino 2021). These results suggest that both the array density and turbine-wake interactions are important for the performance for large finite-sized wind farms. Figure 13 also suggests that finite-sized wind farms are less sensitive to the effective array density than infinitely-large farms. Figure 13a shows a roughly inverse relationship between and / 0 whereas figures 13b-f show a more linear decrease of . Overall, the analytical model (equation 2.9) predicts the variation of with the effective array density well for all of the atmospheric responses. The theoretical model does not predict the effect of turbine-wake interactions because the model of * (equation 2.8) is a function of turbine operating conditions only. These results support the argument that the two-scale momentum theory can be used to provide an approximate upper limit on wind farm performance.

Discussion
The analytical model of * (equation 2.8) has been shown in this study to provide an approximate upper bound to wind farm performance. The model is based on the classical actuator disc theory with an upstream velocity of . Equation 2.8 would provide accurate predictions of * if the following two conditions are met: (1) the wind speed upstream of each turbine in the farm is ; and (2) the mechanism of the flow around each turbine is the same as that around an isolated actuator disc. In reality, the wind speed upstream of each turbine is often lower than due to wake effects, reducing * (and thus ). Conversely, local blockage effects (Nishino & Draper 2015) may increase * because (1) it creates locally accelerated flows which may allow the upstream velocity of most turbines to be higher than (see e.g., Ouro & Nishino 2021); and (2) it changes the mechanism of the flow around each turbine, allowing for a higher (for a given ) than that predicted by the classical actuator disc theory. However, such a positive effect of local blockage can be exploited only when the layout is carefully optimised for a specific wind direction. As such, equation 2.8 can be used to predict the upper bound of farm performance, or the performance of an ideal wind farm without the negative effect of turbine-wake interactions.
The results from this study also show that the of infinitely large wind farms depends mainly on the effective array density (see figure 13a), i.e., the farm power is insensitive to the turbine-scale flow interactions. A closer look at the results suggest that the limited impact of turbine-scale interactions is due to the fact that the turbine drag is typically much greater than the surface drag. The ratio of total turbine drag to surface drag is given by: We measured the time-averaged value of * / * for all 50 turbine layouts simulated. The mean value was 5.22, the minimum value was 2.93 and the maximum was 8.59. Therefore, the turbine drag was typically 5 times greater than the surface drag. Consider an offshore farm where * / * = 5, / 0 = 10 (both typical values) and = 0. This means the turbine drag is 5 times greater than the surface drag and the momentum supplied by the atmosphere does not change in response to the farm. The composition of the total drag (normalised by 0 ) for this scenario is shown by the bar 1 in figure 14a. This state 1 is an equilibrium state, i.e., the normalised drag is balanced by as in equation 2.3. If there is a sudden small change in wind direction which increases the degree of wake interactions then this will decrease * (or the ratio of to 2 ) for the farm. In this example, the turbine drag is now only 4 times greater than the surface drag (bar 2 in figure 14a). This corresponds to * decreasing from 0.75 to 0.6 (which is close to the largest difference observed in the 50 wind farm LES). This state 2 is a non-equilibrium state, i.e., the normalised drag is not balanced by (due to the sudden small change of wind direction). However, if the momentum supplied by the atmosphere is unchanged the wind speed in the farm will eventually increase to compensate for the reduced turbine drag. The surface and turbine drag both scale with 2 so they are both expected to increase at the same rate. Bar 3 in figure 14a shows the new equilibrium state with a new composition of turbine and surface drag after the increase in . Comparing states 1 and 3, the total turbine drag has been reduced only slightly because of the constant amount of momentum supplied by the atmosphere to the farm site. Therefore, for aerodynamically ideal turbines (or turbines operating below the rated wind speed), the average turbine power coefficient is also insensitive to turbine-scale flow interactions when the momentum supplied by the atmosphere is constant ( = 0). Figure 14b explains why turbine-scale flow interactions become more important as increases. Now the momentum supplied by the atmosphere changes with the farm windspeed reduction factor according to = 1 + (1 − ). The ratios between turbine and surface drag for states 1 and 2 are exactly the same as for the = 0 case in figure 14a. However, now when the wind speed increases in response to the reduced turbine drag the momentum supplied by the atmosphere to the farm site changes. As the wind speed increases, increases so the momentum supplied by the atmosphere decreases (see figure  14c). Therefore the sum of the turbine and surface drag for state 3 decreases compared to state 1. This explains why there is a much greater reduction in turbine drag (and thus power) for the non-zero case than for =0. As increases there will be a larger decrease in the total drag in response to an increased wind speed. Therefore the power losses due to turbine-scale flow interactions also increase as increases.
To better understand the factors which determine the power output of wind farms we propose three power loss factors. Firstly, the turbine-scale loss factor Π which is defined by where , ℎ is given by equation 2.9. Π represents the power losses due to turbinescale (or internal) flow interactions only, separate from the losses due to the farm-scale atmospheric response. For real turbines, Π would also include turbine design losses (i.e., power losses due to a non-ideal rotor design for a given ). Figure 15a shows that the losses caused by turbine interactions are small for infinitely large farms, typically less than 5%. For some turbine layouts Π is negative, meaning that exceeds , ℎ . These layouts have * values greater than 0.75 which is the value given by equation 2.8. As discussed earlier, this is likely to be due to local blockage effects increasing the turbine incident velocity above the farm-layer velocity . Figures 15b-f show that the turbine-scale losses are greater in finite-sized farms with the same turbine layout. Under different atmospheric conditions, the same turbine layout can give different turbine-scale losses. Across a realistic range of wind extractability factors, the turbine-scale losses from the same layout can vary significantly. As an example, the losses from one layout varies from 13% to 22% as changes from 5 to 25. Figure 15 suggests the maximum losses due to turbine-scale flow interactions is likely to be about 20% for large offshore wind farms.
The farm-scale loss factor, Π , represents the power loss due to the atmospheric response to the whole farm, and is defined by

c)
Turbine drag Surface drag Figure 14: Examples of offshore wind farm drag composition at 3 different states: 1) equilibrium state before wind direction change, 2) non-equilibrium state immediately after a small change of wind direction, and 3) new equilibrium state after atmospheric response, for a) = 0 and b) = 15. c) schematic of the three states on the vs (1 − ) plot.
where , is given by equation 2.10. Note that in this study = 1.33 and hence , = 0.563 ( = 2 would give the optimal performance for an isolated turbine of , = 16/27). Π represents the power loss accompanied by the reduction of the farmaverage wind speed. Figure 16 shows that the farm-scale losses are typically more than twice the turbine-scale losses (i.e., Π /Π is generally less than 0.5). This suggests that for large offshore wind farms the atmospheric response to the array density is more important than the turbine-scale interactions (i.e., wake effects). Similarly the total power loss factor, Π, can also be defined as It is important to note that the turbine-scale loss factor discussed above is smaller than what is typically referred to as 'wake losses' (see figure 17). 'Wake losses' are traditionally evaluated by comparing the farm power with the power produced by the first row of turbines ( ,1 in figure 17). However, this includes not only the losses due to wake interactions between turbines but also the atmospheric response to the array density. In contrast, Π represents the power losses solely due to the interactions between turbines. For the same reason, the farm-scale loss factor, Π , is larger than what is often referred to as 'farm blockage losses'. These two different classifications of power losses (one using , ℎ and the other using ,1 as a point of reference) are both useful in different ways. The latter classification is straightforward when the value of ,1 is known; however, sometimes ,1 cannot be defined unambiguously (e.g., when there is no regular 'front row' which is perpendicular to the wind direction). The merit of using , ℎ is that this can be predicted analytically using equation 2.9.
It should be noted that figure 17 is for actuator discs (i.e., ideal turbines). Real turbines will experience additional power losses due to practical (non-ideal) turbine design. As demonstrated earlier in figures 13, 15 and 16, the strength of atmospheric response alters wind farm performance. The importance of array density and turbine layout vary with the large-scale atmospheric conditions. The results shown in appendix B suggests that the value of tends to decrease as the wind farm size increases (at least within the range between 10km to 30km; note that the values of reported in Patel et al. (2021) were for a fixed wind farm size of 20km). These results seem to suggest that larger wind farms will be less sensitive to turbine-scale wake effects and more sensitive to farm-scale losses. This could have significant implications for the design of future wind farms.
It is also worth noting that similar trends of large wind farm aerodynamics have already been reported in the literature, e.g., in the wind farm LES performed by Wu & Porté-Agel (2017) in which they also used a relatively low surface roughness length of 0.05m. They performed LES of finite and infinite wind farms under a weak and strong free-atmosphere stratification, using a rotational actuator disc model (representing Vestas V-80 2MW turbine) with an aligned and a staggered turbine layout with the same array density. The varying free-atmosphere stratifications changes the strength of atmospheric response, which affects the wind farm blockage. Figure 11 in Wu & Porté-Agel (2017) shows the power production of different turbine rows for the different wind farms. For both atmospheric conditions, the layout of the infinite wind farm did not change the farm power. This agrees with our finding that the power output of infinite wind farms are insensitive to turbine interactions. The stronger free-atmosphere stratification induced a larger pressure gradient across the farm which implies a larger value. The power output of the finite farm under the strong freeatmosphere stratification was more sensitive to the turbine layout. Under a weak stratification the layout had a smaller impact on farm power. This is further confirmation that the largescale atmospheric response changes the importance of turbine-scale interactions. Other LES studies have also generally found that the power output of finite wind farms (e.g., Porté-Agel et al. 2013;Stevens et al. 2014;Archer et al. 2013) is more sensitive to turbine layout than infinite farms (e.g., Yang et al. 2012;Yang & Sotiropoulos 2014;Abkar & Porté-Agel 2013). However, these previous studies used surface roughness lengths typical of rough onshore locations rather than offshore, meaning that the ratio of total turbine drag to surface drag (equation 5.1) was smaller than the present study.
The impact of turbine layout on wind farm performance is typically investigated by reporting the normalised turbine power, where the turbine power is normalised by either the power of a turbine in the first row or that of a standalone turbine. The normalised power is often used as a measure of turbine-wake interactions within the wind farm. However, this could be misleading for a large wind farm where the power is also reduced due to the ABL response to the farm resistance. The turbine-scale loss factor, Π , gives the power losses due to wake interactions separated from farm-scale effects. Our results suggest that the power losses in large farms with the same layout can change with the strength of large-scale atmospheric response. Hence, the relative performance of different turbine layouts can also vary with the strength of atmospheric response. Layouts which produce more power with a specific large-scale atmospheric response may not perform as well under different responses.
A limitation of this study is that the wind farm LES consider only neutrally stratified atmospheric boundary layers. One uncertainty is how different stratifications would affect wake recovery within large wind farms, which could affect the * value for a given turbine layout. However, since the two-scale momentum theory has been derived from the principle of momentum conservation which is valid for all atmospheric conditions, we expect that the trends found in this study (on how the importance of turbine-scale interactions changes with the strength of the atmospheric response) would be observed generally. The study also used a fixed turbine operating condition given by = 1.33. This is typical for the current wind turbines operating below their rated wind speed (at which the farm blockage effects are most significant). Therefore the results of this study give the trends for large farms if the operating conditions are similar to those currently used. As discussed earlier, the trends found in this study have also been observed in LES of finite wind farms with different atmospheric conditions and a different turbine model, supporting the generality of our findings.

Conclusions
In this study we performed a combined theoretical and computational analysis of large wind farms, using the two-scale momentum theory (Nishino & Dunstan 2020) and new LES of infinitely large wind farms with different turbine layouts and wind directions. To consider a range of wind directions, we used a new implementation of the actuator disc model which is not aligned with the structured grid. Firstly, we validated the LES against the results published by Calaf et al. (2010) and found an excellent agreement. We also confirmed that, when the correction factor proposed by Shapiro et al. (2019) is applied, the average internal turbine thrust coefficient * (representing the turbine thrust relative to the square of the average wind speed across the wind farm layer) is insensitive to the grid resolution. We then performed 50 wind farm LES with different turbine spacings and wind directions. * was found to be a strong function of wind direction with a weaker dependence on turbine spacings.
This study adopted an approach to adjust the average power coefficient of an infinitely large wind farm to estimate the power of large but finite-sized farms with the same layout, following the concept of the two-scale momentum theory. A full validation of this approach is currently infeasible as it would require a large set of finite-size wind farm LES, which will be the subject of future studies. However, we have provided detailed theoretical arguments on why we can estimate the power of large finite-size farms from the results of infinitely large farm LES. For infinitely large farms, was found to be a strong function of the array density and agree remarkably well with the analytical model (equation 2.9). The power produced by large finite-sized farms was found to depend on both the array density and turbine layout (including the effect of wind direction), although the analytical model still provides an approximate upper limit of the power (for a given array density). These results confirm that, to model large but finite-sized wind farms it is important to consider both the effect of array density and the turbine layout.
The impact of array density and turbine layout on farm power varies with the strength of atmospheric response or the 'wind extractability'. The analytical model seems to predict very well the impact of array density for a given farm-scale atmospheric response. As such, we propose a new classification of power losses in large wind farms by introducing new metrics, namely a turbine-scale loss factor Π , a farm-scale loss factor Π and a total loss factor Π. The turbine-scale loss factor describes the power losses due to turbine-scale wake interactions only (which is smaller than what is commonly called 'wake losses'). Importantly, the turbine-scale loss factor varies with the strength of large-scale atmospheric response. As an example, Π varied from 0.13 to 0.22 for a given layout across a realistic range of the extractability factor (from = 5 to 25 in this study). Although further studies are required to better quantify the extractability factor, the results obtained to date suggest that as wind farms get larger the proportion of power losses due to turbine-wake interactions will decrease. This is because (1) having a better turbine layout (with less wake interactions and higher * ) does not substantially increase the momentum available to the farm site when is small ( figure 14) and (2) seems to decrease as the wind farm size increases. The farm-scale loss factor describes the power loss due to the atmospheric response to the array density (which is larger than what is often called 'farm blockage losses'). Π was typically more than twice as large as Π for the 50 offshore wind farms considered in this study with a realistic range of atmospheric responses. This suggests that farm-scale flow effects have a greater impact than turbine-scale flow effects on the performance of large offshore wind farms. = 3 ℎ . Figure 18 compares the values obtained for = 2.5 ℎ and 3 ℎ . Changing the value of has no effect on the theoretical predictions of for a given . Increasing decreases the LES-based results of for each farm because the reference velocity 0 slightly increases. The increase in 0 is the same for all wind farm cases. Therefore, increasing to 3 ℎ reduces the of all wind farm cases by the same factor of 1.043 in this study. The absolute values of are affected by the value of but the trends reported in this study remain unchanged. The * values obtained from the LES are also reduced when the farm-layer height is increased. This is because the reference velocity slightly increases. Nevertheless, the two-scale momentum theory still captures the trend well and can provide an upper bound estimate to wind farm performance.
Nishino & Dunstan (2020) proposed a definition for based on the undisturbed velocity profile 0 ( ).
is defined as the farm-layer height with which the farm-layer average of the undisturbed velocity is equal to that of the rotor average, i.e., Using this definition, the exact value of will vary with turbine design and (undisturbed) ABL profile. Figure 19 shows the value of according to A 1 for a wide range of velocity profiles and turbine designs. We assume that the velocity profile is given by the logarithmic law for a wide range of surface roughness lengths that includes onshore and offshore locations. For the turbine design, we assume that the distance between the surface and bottom of the rotor is fixed at 50m. The range of / ℎ used corresponds to turbine diameters from 100m up to 300m. Across a wide range of turbine design and velocity profiles the value of obtained from equation A 1 only varies between 2.5 ℎ and 3 ℎ . The average turbine power coefficient is changed only slightly when is increased from 2.5 ℎ to 3 ℎ . The trends reported in this study are unaffected by a change in . is relatively insensitive to velocity profiles and turbine design and is expected to be in the range 2.5 ℎ to 3 ℎ . Therefore, the results of the present study are insensitive to the exact value of . Figure 19: Variation of normalised wind farm-layer height / ℎ with surface roughness length 0 and turbine design parameter / ℎ for a log law wind profile.

Appendix B. Effect of wind farm size on the wind extractability factor
To explore the effect of wind farm size on the wind extractability factor , we performed additional twin NWP simulations using the same methodology as Patel et al. (2021). Within the simulations, the diameter of a hypothetical circular wind farm (located off the east coast of Scotland) was varied from 10km to 30km. Example results for a 24-hour period with a relatively constant wind speed across the farm site (corresponding to Case B of figure  9 of Patel et al. (2021)) are shown in figure 20. It can be seen that the value of tends to decrease with increasing the wind farm size for all the time. Further simulations are required in future studies to derive a correction function for the effect of wind farm size on , but our preliminary results suggest that an exponential relationship may exist between the wind farm size and . Such an exponential correction would satisfy the expected asymptotic behaviour of (i.e., approaches to +∞ and 0 as the farm size approaches to 0 and +∞, respectively). Figure 20 also shows that varies significantly with time due to changing atmospheric conditions. These preliminary results show that depends both on wind farm size and atmospheric conditions.