4.1. Flow physics

We start our analysis with figure 5, which shows a side view of turbulent kinetic energy (TKE) perturbation with respect to the precursor simulation averaged in the $y$ direction along the width of the farm, together with the base and top of the inversion layer computed by fitting the LES data with the Rampanelli & Zardi (Reference Rampanelli and Zardi2004) model. In the H150 and H300 cases, the vicinity of the capping inversion to the turbine-tip height limits the flow development in the vertical direction. Moreover, the capping inversion behaves as a pliant surface, limiting flow entrainment from the free atmosphere. Further, the high stability attained within the inversion layer dampens turbulence, therefore severely limiting the vertical development of the internal boundary layer. Consequently, figure 5(a,b) display relatively high values of TKE only within the farm, which then decay downwind of the last row of turbines. Moreover, the very low vertical mixing of momentum due to the vicinity of the inversion layer significantly reduces the wind shear at hub height in the wake region (see also figure 6 a), which further causes the TKE to drop below the values attained in the precursor simulations. The combination of these effects leads to very strong velocity deficit in the farm wake and a slow recovery rate, as shown in figure 6(a,b). The increase in the $H/z_{{tip}}$ ratio for cases H500 and H1000 allows for higher TKE values and accounts for a significant vertical development of the internal boundary layer, as shown in figure 5(c,d). The latter enhance vertical turbulent transport of momentum, therefore giving rise to lower wake deficit and a faster recovery rate, as visible in figure 6(c,d). Moreover, case H1000 shows the formation of a high speed channel between the tip height and the capping-inversion base, which further increases the TKE and flow mixing. Stieren & Stevens (Reference Stieren and Stevens2022) mentioned that a wind farm operating in the wake of an upstream one is negatively affected in terms of power output, but has a higher wake efficiency due to the increased turbulence level in the incoming inflow. Figure 5 shows that this is the case only in deep boundary layers and only when the farm is less distant than $10$ – $20$ km, as the TKE decays much faster than the velocity deficit. The latter, shown in figure 6, illustrates that the wake is not fully replenished in all cases, meaning that it would impact on downwind farms located more than $70$ km downstream. Finally, above the inversion layer, the internal gravity waves generate a slanted pattern of alternating positive and negative perturbations in all fields (Allaerts & Meyers Reference Allaerts and Meyers2017). The phase of these perturbations in each field is governed by the polarisation equation, which is derived from linear gravity-wave theory (Nappo Reference Nappo2002). This equation shows that the perturbations in the streamwise velocity component are out of phase with respect to the pressure perturbation (see Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ) Appendix C for a simplified derivation). Since internal waves generate a negative pressure perturbation in the free atmosphere and downstream of the maximum capping-inversion displacement (i.e. above the farm), this results in a flow speed-up in the same region, particularly accentuated in shallow-boundary-layer cases.

Figure 5.

Contours of the time-averaged TKE perturbation with respect to the precursor simulation taken in an $x$ – $z$ plane further averaged along the farm width in the spanwise direction for cases (a) H150, (b) H300, (c) H500 and (d) H1000. The black lines represent the bottom and top of the inversion layer computed with the Rampanelli & Zardi (Reference Rampanelli and Zardi2004) model. The locations of the turbine-rotor disks are indicated with vertical white lines.

Figure 6.

Contours of the time-averaged horizontal velocity magnitude in an $x$ – $z$ plane further averaged along the farm width in the spanwise direction for cases (a) H150, (b) H300, (c) H500 and (d) H1000. The black lines represent the bottom and top of the inversion layer computed with the Rampanelli & Zardi (Reference Rampanelli and Zardi2004) model. The locations of the turbine-rotor disks are indicated with vertical white lines.

Figure 7.

Contours of the time-averaged (a–d) horizontal velocity magnitude, (e–h) flow angle, (i–l) vertical velocity, (m–p) pressure perturbation with respect to the inflow, and (q–t) TKE perturbation with respect to the precursor simulation in an $x$ – $y$ plane, taken at hub height for cases (a,e,i,m,q) H150, (b,f,j,n,r) H300, (c,g,k,o,s) H500 and (d,h,l,p,t) H1000. The location of the wind farm is indicated with the black rectangle. The last 10 km of domain in the streamwise direction are not displayed as the solution is influenced by the fringe forcing in that region.

Next, figure 7(a–d) show a top view taken at hub height of horizontal velocity magnitude. First, we observe that in all cases, the velocity deficit does not seem to spread horizontally. A similar behaviour can be observed in numerical simulations and SAR images, such as figure 1 (Baas & Verzijlbergh Reference Baas and Verzijlbergh2022; Maas & Raasch Reference Maas and Raasch2022; Finserås et al. Reference Finserås, Herrera Anchustegui, Cheynet, Gebhardt and Reuder2024). Therefore, the assumption that the velocity deficit spreads linearly, which is typically made in analytical wake models, seems to not hold for wind-farm wakes. As shown previously, the vicinity of the capping inversion to turbine-tip height in cases H150 and H300 limits the flow development in the vertical direction. Therefore, to conserve mass, high-speed regions form at the sides of the wake. The same phenomenon was also observed by Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ). The strong flow deceleration attained in these two cases decreases the Coriolis force magnitude, which scales linearly with the wind speed. Hence the wake turns towards the direction of the background pressure gradient, therefore undergoing an anticlockwise flow rotation.

In the H500 and H1000 cases, the wake recovery rate is higher, enhanced by the vertical turbulent entrainment of momentum. The latter adds clockwise-turning flow from aloft into the wake region. Therefore, a clockwise wake deflection is observed, which is dictated by the wind veer present in the background flow. This behaviour is also illustrated in figure 7(e–h), which display the flow angle at hub height. These plots clearly show that effects of Coriolis forces and wind veer are opposed. Whether one dominates the other depends upon the capping-inversion height. Moreover, we can also conclude that in the absence of wind veer in the background flow, the wind-farm wake would be always deflected anticlockwise in the northern hemisphere. We note that similar results have been observed for both single-turbine and wind-farm wake simulations (van der Laan & Sørensen Reference van der Laan and Sørensen2017; Howland, Ghate & Lele Reference Howland, Ghate and Lele2018; Heck & Howland Reference Heck and Howland2025). The typical signature of trapped waves is visible in the vertical velocity field shown in figure 7(i–l), particularly in the shallow-boundary-layer cases where the capping inversion is closer to hub height. The V-shaped pattern attained in sub-critical flows, also visible in the pressure perturbation fields shown in figure 7(m–p), has been already observed and analysed in detail by Allaerts & Meyers (Reference Allaerts and Meyers2019) and Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ). However, we notice the presence of slanted waves with a high wavenumber that cross the full domain in the streamwise direction, particularly visible in cases H150 and H300. The same behaviour has been observed by Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ), who showed that these slanted waves originate from the perturbation introduced by the turbine spacing, and their angle with respect to the $x$ -axis can be predicted with linear theory for deep-water waves. We refer the reader to Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ) for more details. Finally, we note in figure 7(i,j) the presence of a region of positive vertical velocity on the left edge of the wake. This denotes a region of flow convergence, which will be discussed later in this section. To conclude, figure 7(q–t) illustrate the TKE perturbation with respect to the precursor simulation at hub height. As expected, a higher TKE level is obtained in deeper boundary layers. For shallow boundary layers, the TKE at hub height in the wake of the farm is smaller than the precursor value, as also shown in figure 5(a,b). This is a consequence of the reduced velocity shear in the farm wake in the region below $H$ , which therefore diminishes the production of turbulence. Similar results have been observed in wind tunnel experiments and mesoscale models (Chamorro & Porté-Agel Reference Chamorro and Porté-Agel2009; Fitch et al. Reference Fitch, Olson, Lundquist, Dudhia, Gupta, Michalakes and Barstad2012). Similar to Maas (Reference Maas2023b ), we notice a higher level of TKE on the left side of the wake, i.e. at $y=55$ km.

Figure 8.

Contours of the time-averaged streamwise velocity in an $y$ – $z$ plane focused on the wake region and its sides taken at (a) $x=50$ km, (b) $x=70$ km and (c) $x=100$ km for case H300. The arrows represent the in-plane velocity vector given by the spanwise and vertical velocity components. The vector is normalised with its magnitude $(\bar {v}^2+\bar {w}^2)^{1/2}$ , so that all arrows have the same length. The blue lines represent the bottom and top of the inversion layer computed with the Rampanelli & Zardi (Reference Rampanelli and Zardi2004) model. The vertical black dashed lines denote the centreline of the mixing layers located at the sides of the wake – see the fitting model introduced in § 4.2.

Next, we investigate the flow behaviour in the wake region and at its sides. To do so, we illustrate $y$ – $z$ planes of streamwise velocity taken at $x=50$ , $70$ and $100$ km for case H300 in figure 8. Moreover, we add vectors that denote the flow direction based on the spanwise and vertical velocity components, while their colours indicate the in-plane velocity magnitude (i.e. $(\bar {v}^2+\bar {w}^2)^{1/2}$ ). Further, we remark that due to the convention used, the right side of the wake corresponds to the left side of the figure, and vice versa – see figure 2.

At $x=50$ km (i.e. $20$ km downstream of the farm), figure 8(a) shows a strong streamwise velocity deficit, which spreads over the whole ABL height. Within the wake and on its right side, the flow direction follows the wind veer profile of the background flow (i.e. a positive and negative spanwise velocity below and above turbine-hub height). However, the zero flow-angle height is shifted upwards, close to turbine-tip height. This is a consequence of the reduced Coriolis force and absence of vertical mixing, which causes the flow to turn towards the background pressure gradient. Moreover, we observe the presence of anticlockwise flow recirculation within the wake region, induced by the trapped waves that also corrugate the base and top of the capping inversion. The flow near the left edge of the wake behaves differently. Here, we observe an anticlockwise motion of the wake, with a positive spanwise component below tip height, an upward motion at the wake edge, and a negative spanwise component near the base of the capping inversion. As a result, the spanwise velocity component attains values close to zero at this edge (i.e. at $y\approx 56$ km), significantly limiting the advection of mean flow – see § 4.3. On the left side, we observe the presence of trapped waves which leads to an anticlockwise flow motion.

Figure 8(b,c) show the flow behaviour at $x=70$ km and $x=100$ km, respectively. As we move downstream, the streamwise velocity component recovers at a slow rate along the downstream direction, while the wake gradually narrows. The zero flow-angle height moves downwards towards turbine-hub height, with the wind veer profile converging towards that of the background flow within the wake and on its right side. However, the anticlockwise motion at the left edge of the wake does not die out and is still clearly visible $50$ km downstream of the farm. Contrarily, the trapped waves at the left side of the wake are gradually dampened along the streamwise direction, so that the anticlockwise flow rotation at the left side of the wake vanishes at $x=100$ km, and the flow assumes a negative spanwise component at each height. Finally, figure 8(b,c) also illustrate how the strong wind veer attained on the right edge of the wake causes the streamwise velocity deficit to spread horizontally.

Figure 9.

Contours of the time-averaged streamwise velocity in an $y$ – $z$ plane focused on the wake region and its sides taken at (a) $x=50$ km, (b) $x=70$ km and (c) $x=100$ km for case H1000. The arrows represent the in-plane velocity vector given by the spanwise and vertical velocity components. The vector is normalised with its magnitude $(\bar {v}^2+\bar {w}^2)^{1/2}$ , so that all arrows have the same length. The vertical black dashed lines denote the centreline of the mixing layers located at the sides of the wake – see the fitting model introduced in § 4.2.

The same analysis for case H1000 is shown in figure 9. Here, the flow behaviour is significantly different from that obtained in shallow boundary layers. In fact, the strong flow mixing along the vertical direction adds clockwise-turning flow from aloft into the wake region. Consequently, figure 9(a) illustrates that the streamwise velocity deficit recovers at a much faster rate, with the spanwise velocity component being negative at each height within the wake region. Moreover, the wind veer is significantly reduced if compared to that of the background flow. Further, we also observe a much stronger in-plane velocity magnitude compared to case H300. At the sides of the wake, the wind veer follows the same orientation as that attained within the wake region. Moreover, flow recirculation induced by trapped waves remains very limited. In fact, the influence of trapped waves on the flow behaviour within the rotor height is very limited in deep boundary layers – see figure 7(l). Further downstream, figure 9(b,c) show that the streamwise velocity deficit is mostly recovered and the wake deflects clockwise. Moreover, the streamwise velocity deficit does not spread along the horizontal direction as a consequence of the reduced wind veer.

4.2. Wake properties

We have seen in the previous subsection that the wake region is heavily dependent upon the capping-inversion height. In this subsection, we will investigate the vertical and horizontal wake development. Moreover, we will rigorously define the wake edges and centre, which we use to evaluate the wake width, strength and deflection. To do so, we propose a simple function that we use to fit the velocity magnitude profiles along the spanwise direction:

(4.1) \begin{equation} F(y) = \begin{cases} M_{{side}}& \text{if } y_{b,r} \leqslant y \lt y_r,\\ M_{{side}} - \Delta M (y-y_r)/\delta _r & \text{if } y_r \leqslant y \lt y_r + \delta _r, \\ M_{{wake}} & \text{if } y_r + \delta _r \leqslant y \lt y_l - \delta _l, \\ M_{{wake}} + \Delta M (y-y_l+\delta _l)/\delta _l& \text{if } y_l - \delta _l \leqslant y \lt y_l, \\ M_{{side}}& \text{if } y_l \leqslant y \lt y_{b,l}, \end{cases} \end{equation}

where $M_{{side}}$ and $M_{{wake}}$ denote the velocity magnitudes attained at the sides and in the wake of the farm, while $y_r$ and $y_l$ denote the right- and left wake-edge locations, located in the region around $y=45$ km and $y=55$ km, respectively. We remark that due to the convention used, the right side of the wake corresponds to the left side of the figure, and vice versa. Moreover, we note that the velocity $M_{{side}}$ is measured in a region that extends $10$ km sideways of the farm location, so that its value does not depend upon the spanwise domain dimension. Hence $y_{b,r}$ and $y_{b,l}$ assume values $35$ km and $65$ km, respectively. Further, $\delta _r$ and $\delta _l$ represent the distances over which the velocity decays from $M_{{side}}$ to $M_{{wake}}$ . Further, $\Delta M = M_{{side}} - M_{{wake}}$ indicates the velocity deficit with respect to the wind speed at the wake sides. Note that this is also the definition of velocity deficit typically used in real case scenarios (Djath et al. Reference Djath, Schulz-Stellenfleth and Cañadillas2018). Following Pope (Reference Pope2000), we define the farm wake width as the distance between the centreline of the regions $\delta _l$ and $\delta _r$ , i.e. $\delta _{{wake}} = y_l - \delta _l/2 - (y_r + \delta _r/2)$ . Further, we define the wake centre as the middle point of the region where $M_{{wake}}$ is attained, which is formally defined as $y_c=(y_l-\delta _l+y_r+\delta _r)/2$ . The fitting function $F(y)$ reported in (4.1) together with a visual representation of its parameters is shown in figure 10.

Figure 10.

Sketch of the fitting function $F(y)$ and its parameters – see (4.1).

Figure 11.

Time-averaged horizontal velocity magnitude as a function of the spanwise direction obtained at $x=x_0$ (brightest yellow) and $x=x_0+70$ (darkest blue), with increments of $5$ km in between, for cases (a) H150, (b) H300, (c) H500 and (d) H1000. Each symbol denotes a spanwise average over one turbine spacing (i.e. $5D$ ). The red dashed lines illustrate the results of our fitting model for each streamwise location. The vertical dashed black lines denote the location of the first- and last-column turbines. We define the regions around $y=45$ km and $y=55$ km as the right- and left edges, respectively, and $x_0=30$ km denotes the beginning of the wake region.

Figure 12.

(a,b) Wake width, wake centre and (c,d) right- and left mixing layer widths normalised with the farm width. (e,f) Velocity deficit with respect to $M_{{side}}$ , and velocity magnitude at the sides of the wake, normalised with the velocity magnitude obtained in the precursor domain at hub height. All quantities are predicted by the fitting model shown in (4.1). The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

The only purpose of the model reported in (4.1) is to identify the wake location and shape in a horizontal plane given a velocity profile along a direction normal to the farm wake. The profiles obtained for this study are shown in figure 11 at various locations downwind of the farm (i.e. from $x=x_0$ to $x=x_0+70$ with increments of $5$ km). figure 11(a,b) show the profiles for cases H150 and H300. Here, we can observe strong velocity deficits within the first $10$ km, where the profiles assume a top-hat shape with peaks in the location of the first and last column of turbines. This behaviour may be caused by the strong flow redirection at the farm sides, which yaws the first and last columns of turbines up to 8 $^\circ$ in case H150, allowing for a higher energy extraction and therefore causing higher velocity deficit at $y=45$ km and $y=55$ km (Lanzilao & Meyers Reference Lanzilao and Meyers2024b ). In these cases, the wake turns anticlockwise. Therefore, the velocity gradient within the right edge of the wake (i.e. at $y=45$ km) diminishes along the streamwise direction. Moreover, the vertical velocity profiles at the lateral edges of the wake are only slightly dependent on the vertical direction for $z\lt H$ and the ratio $\delta _{{wake}}/H \gg 1$ – see figure 12(a). Consequently, the regions $\delta _l$ and $\delta _r$ can be interpreted as two distinct mixing layers with width varying along the streamwise direction. In figure 11(a,b), we observe that the mixing layer located at the right edge of the wake expands along the streamwise direction, with its centreline turning towards the low-velocity side. On the contrary, the mixing layer on the left edge of the wake (i.e. at $y=55$ km) has a limited growth, possibly due to the fact that the wake turns in the direction opposite to the wind veer of the background flow. Next, figure 11(c,d) show the results for the deep-boundary-layer cases, i.e. H500 and H1000. The velocity profiles fully resemble a top-hat shape with smooth edges. The wake turns clockwise in these cases, and the velocity speed-up at the wake sides is less pronounced. Finally, for case H1000, with $\delta _{{wake}}/H \approx 10$ , we find that the mixing layer analogy no longer works, with velocity profiles that depend significantly on height (see figure 15), and wake recovery that is dominated by vertical entrainment – see § 4.3. The velocity magnitude profiles discussed above are then fitted using the model reported in (4.1). The fitting is performed using the curve_fit method provided by the Scipy library, adopting the trust region reflective algorithm (Virtanen et al. Reference Virtanen2020). The results are shown in figure 11, where the red dashed lines illustrate the fitted velocity profiles. Here, we can observe that the model captures the wake location and strength very well.

The wake width and centre, the velocity deficit and the velocity at the wake sides, together with the mixing layers width, all evaluated at hub height and estimated using the fitting model reported in (4.1), are shown in figure 12. The wake width, which is defined as the distance between the left- and right mixing layers, is shown in figure 12(a). Here, we can observe a width $10$ – $20$ % higher than $L_y^f$ at the beginning of the wake region. Further downstream, the mixing layer centrelines turn towards the low-velocity side, resulting in a reduction of $\delta _{{wake}}$ . At $x-x_0=65$ km, cases H150, H300 and H500 show a wake-to-farm width ratio of roughly 1, while case H1000 displays a narrower wake. The wake narrowing along the streamwise direction can also be observed in SAR images and in the results shown by Baas & Verzijlbergh (Reference Baas and Verzijlbergh2022), Maas & Raasch (Reference Maas and Raasch2022) and Baas et al. (Reference Baas, Verzijlbergh, van Dorp and Jonker2023). In particular, Maas (Reference Maas2023b ) performed a wind-farm simulation with inflow conditions similar to our H300 case, and observed very similar trends in both wake width and deflection. The wake centre is displayed in figure 12(b). Here, we observe how the capping-inversion height influences the wake deflection, generating a clear distinction between shallow and deep boundary layers. The strongest wake deflection is observed in case H1000, where the wake centre deviates by approximately $45$ % of the wake width $65$ km downstream of the farm. Next, figure 12(c,d) show the widths of the regions at the side of the farm over which the velocity decays from $M_{{side}}$ to $M_{{wake}}$ . The region $\delta _r$ becomes three times wider while going from $x_0$ km to $x_0+70$ km in the H150 case, while it shows a more limited growth in the other cases. Except for case H1000, this growth can be interpreted as an expanding mixing layer along the streamwise direction. On the left side, the region $\delta _l$ preserves approximately its width, as illustrated in figure 12(d). A growth is observed only in the H1000 cases, although it happens mostly within the first $10$ km downwind of the farm. These behaviours can also be observed in figure 11.

Next, figure 12(e) illustrates the velocity deficit with respect to the velocity measured at the sides of the wake and further normalised with the precursor velocity at hub height. In all cases, the wake shows a very fast recovery rate within the first $5$ km. Afterwards, the velocity deficit attains a very slow recovery, except for case H1000 where $\Delta M$ keeps decreasing up to $40$ km downwind of the farm. We note that the ratio $\Delta M/M_{{prec}}$ measures approximately $0.25$ at the end of the wake region in cases H150 and H300, confirming that the wake is not replenished. Moreover, we remark the presence of oscillations in the $\Delta M$ profiles with decreasing amplitude for increasing values of $H$ . These oscillations are related to the flow divergence–convergence generated by interfacial waves (see figure 7 i–l), and therefore are more accentuated in shallow boundary layers. Figure 12(f) shows the velocity at the sides of the wake normalised with $M_{{prec}}$ . Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ) have shown that shallow boundary layers cause the flow to speed up at the farm sides. The same phenomenon is visible here. For instance, in case H150, $M_{{side}}/M_{{prec}}$ goes from $0.99$ to $1.07$ in the first $20$ km downwind of the farm. A high capping inversion allows the internal boundary layer to grow vertically, resulting in lower $M_{{side}}/M_{{prec}}$ ratios. We also note that the maximum $M_{{side}}/M_{{prec}}$ value shifts downstream as $H$ decreases.

Figure 13.

Mixing layers centreline and wake centre measured at $z_h-D/2$ , $z_h$ and $z_h+D/2$ for cases (a) H150, (b) H300, (c) H500 and (d) H1000 predicted by the fitting model shown in (4.1). The black horizontal dashed lines denote the locations of the first and last turbine columns. The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

Figure 14.

Contours of the time-averaged horizontal velocity magnitude downwind of the farm for cases (a) H150, (b) H300, (c) H500 and (d) H1000 in an $x$ – $y$ plane taken at hub height. The black lines denote the centreline of the mixing layers, while the red dashed lines denote their width. The pink line indicates the wake centre. All quantities are predicted by the fitting model shown in (4.1). The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

Up to this point, we have investigated only the wake properties at hub height. Looking at the wake evolution at different heights, figure 13 shows the centreline of the mixing layers $y_r+\delta _r/2$ and $y_l-\delta _l/2$ together with the wake centre at the heights $z_h-D/2$ , $z_h$ and $z_h+D/2$ , for all cases. While the wake turns anticlockwise at hub height for case H150, we can observe a strong clockwise rotation of the velocity deficit at tip height. This effect is caused by the strong wind veer attained within the capping inversion, which is very near to tip height in this case. Moreover, the noise observed in the blue lines around $x=90$ km is due to the vicinity of the capping inversion, where the velocity deficit assumes a different profile than the one that our model can fit. In all other cases, the wake edges and centre do not differ substantially across the rotor height. Moreover, we also note that the wake recovery rate does not vary across the rotor height (not shown).

Next, we apply a quadratic fit along the streamwise direction to the profiles obtained with the fitting model in order to average out small oscillations, and we superimpose them to the hub-height horizontal velocity magnitude field. The results for all cases are shown in figure 14. Here, we can observe that we correctly track the wake evolution along the streamwise direction. Moreover, we can observe the evolution of the mixing layers width at the sides of the wake. This makes the fitting model reported in (4.1) a useful tool for inferring wake properties whenever a horizontal plane of velocity in the wake of a farm is available, which could be provided by numerical simulations, lidar measurements or SAR images.

Figure 15.

Vertical profiles of velocity deficit $\Delta M = M_{{side}} - M_{{wake}}$ obtained at $x=x_0$ (brightest yellow) and $x=x_0+70$ (darkest blue), with increments of $5$ km in between, for cases (a) H150, (b) H300, (c) H500 and (d) H1000. Here, $M_{{side}}$ and $M_{{wake}}$ are obtained by averaging the vertical profiles over the wake and side regions defined by the fitting model reported in (4.1). The grey dashed line denotes the turbine-hub height, while the black dashed lines are representative of the rotor dimensions.

Understanding the behaviour of wind-farm wakes along the horizontal directions is of fundamental importance for assessing the impact on the downwind farm. However, the vertical spread of the wake is also of interest. Therefore, we show in figure 15 vertical profiles taken at different locations downwind of the farm for all cases. To be consistent with previous results, we define $\Delta M = M_{{side}} - M_{{wake}}$ , where $M_{{side}}$ and $M_{{wake}}$ are obtained by averaging the vertical profiles over the wake and side regions defined by the fitting model reported in (4.1). The results are shown in figure 15 for all cases. At $x=x_0$ (i.e. in close proximity to the last row of turbines), the velocity deficit resembles a Gaussian profile in all cases. However, the profiles quickly evolve towards a top-hat distribution, with the velocity deficit spreading up to the capping-inversion height for the whole streamwise extent of the wake. Only case H1000 shows a much more gradual wake growth downwind of the farm before also converging towards a velocity deficit profile with a top-hat shape.

To conclude, we show in figure 16 the turbulence intensity measured at hub height and as a function of the streamwise direction. We note that these profiles are further averaged along the $y$ direction over the wake width evaluated with the model reported in (4.1). For case H150, the turbulence intensity attains a constant value $20$ km downwind of the farm, while $50$ km is necessary in case H1000 for the farm-added turbulence intensity to decay. While the turbulence intensity tends to the precursor value in cases H500 and H1000, a lower value is attained in shallow-boundary-layer cases due to the reduced wind shear within the turbine rotor disk with respect to the background flow. Similarly to Maas (Reference Maas2023a ), we observe that the decay of the turbulence intensity can be well approximated by the expression $a_{t} \exp(-b_{t} (x-x_0)^{c_t} ) +d_{t}$ , where $\Phi _{{TI}} = [a_{t} ,b_{t} ,c_{t} ,d_{t}]$ represents a set of fitting parameters, and $x_0=30$ km denotes the beginning of the wake region. The dashed lines in figure 16 represent the fitting function. The fitting parameters for each capping-inversion height are reported in the caption of figure 16. These parameters can support the development and calibration of turbulence intensity parametrisations for engineering wake models (Crespo & Hernandez Reference Crespo and Hernandez1996).

Figure 16.

Time-averaged turbulence intensity further averaged over the wake width as a function of the streamwise direction further normalised with the turbulent intensity value attained in precursor simulation – see table 1. The dashed lines are obtained using the expression $a_t \exp(-b_t(x-x_0)^{c_t}) +d_t$ , with $x_0=30$ km, and $\Phi _{\textit{TI}} = [a_t,b_t,c_t,d_t]$ is the set of fitting parameters, with $\Phi _{\textit{TI}}^{H150}=[8.79,0.34,0.88,2.94]$ , $\Phi _{\textit{TI}}^{H300}=[10.44,0.29,0.92,3.30]$ , $\Phi _{\textit{TI}}^{H500}=[12.72,0.25,0.81,3.61]$ and $\Phi _{\textit{TI}}^{H1000}=[16.21,0.39,0.52,3.59]$ . The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

4.3. Streamwise momentum budget analysis

To highlight the dominant recovery mechanisms for the different cases, we perform a streamwise momentum budget analysis. A mass balance analysis is also performed, with results reported in Appendix C. We remark that in both analyses, we focus solely on the region downwind of the farm. We refer to Lanzilao (Reference Lanzilao2023) and Lanzilao & Meyers (Reference Lanzilao and Meyers2024b ) for an in-depth analysis of the momentum and energy budget upstream and within the farm.

The momentum analysis adopts a control volume $\Omega$ that measures $s_x D/2$ along the streamwise direction. The vertical dimension of the control volume coincides with the turbine rotor height, i.e. from $z_1 = z_h - D/2$ to $z_2 = z_h + D/2$ . Since we focus on the momentum balance within the farm wake, the lateral extension of $\Omega$ is dictated by the wake width provided by the fitting model discussed in § 4.2. Hence the lateral faces of the control volume change location along the streamwise direction and are defined as $y_1 = y_r + \delta _r/2$ and $y_2=y_l-\delta _l/2$ (i.e. the centreline of the mixing layers). We note that we assume the lateral edges of the control volume to be constant with height, meaning that $y_1$ and $y_2$ are not functions of $z$ . Figure 13 shows that this is a valid assumption, apart from case H150 where the farm wake abruptly deviates clockwise at the turbine-tip height. In summary, a generic control volume $\Omega _i$ has dimensions $s_x D/2 \times \delta _{{wake}}(x_i) \times D$ . We remark that we denote with $x_1$ , $x_2$ , $y_1$ , $y_2$ , $z_1$ and $z_2$ the boundaries of the control volume, while the $y$ – $z$ , $x$ – $z$ and $x$ – $y$ faces are denoted with $\Gamma _x$ , $\Gamma _y$ and $\Gamma _z$ .

Figure 17.

Streamwise variation of momentum sources and sinks for case H300 normalised by the time-averaged total wind-farm thrust force and the control volume width: (a) all terms of (4.2), including $\mathcal{S}=\Delta \mathcal{A}_{uu,x} + \Delta \mathcal{A}_{uw,z} + \mathcal{P}_x$ ; (b) the contribution of each added in the $\mathcal{S}$ term. Only the region downwind of the farm is considered, so the term $\mathcal{F}_{t,x}$ has no contribution. The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

The budget equation is obtained by taking a time average of the streamwise momentum equation, further integrating it over the control volume $\Omega$ to average out local oscillations. Additionally, we apply the divergence theorem, enabling us to eliminate the divergence operator and transition from a volume to a surface integral. As a result, the streamwise momentum equation for any control volume $\Omega$ reads as

(4.2) \begin{equation} \begin{aligned} &\underbrace {-\Big [ \int _{\Gamma _x} \bar {u}\bar {u}\, {\rm d}\Gamma _x\Big ]_{x_1}^{x_2}}_{\Delta \mathcal{A}_{uu,x}} \underbrace {-\Big [ \int _{\Gamma _y} \bar {u}\bar {v}\, {\rm d}\Gamma _y\Big ]_{y_1}^{y_2}}_{\Delta \mathcal{A}_{uv,y}} \underbrace {-\Big [ \int _{\Gamma _z} \bar {u}\bar {w}\, {\rm d}\Gamma _z\Big ]_{z_1}^{z_2}}_{\Delta \mathcal{A}_{uw,z}} \underbrace{{-}\Big [ \int _{\Gamma _x} \Big (\overline {u'u'}^r+\overline {u'u'}^{{sgs}} \Big )\, {\rm d} \Gamma _x \Big ]_{x_1}^{x_2}}_{\Delta \mathcal{F}_{uu,x}}\\[5pt] & \underbrace {-\Big [ \int _{\Gamma _y} \Big (\overline {u'v'}^r+\overline {u'v'}^{{sgs}} \Big )\, {\rm d} \Gamma _y \Big ]_{y_1}^{y_2}}_{\Delta \mathcal{F}_{uv,y}} \underbrace {-\Big [ \int _{\Gamma _z} \Big (\overline {u'w'}^r+\overline {u'w'}^{{sgs}} \Big )\, {\rm d} \Gamma _z \Big ]_{z_1}^{z_2}}_{\Delta \mathcal{F}_{uw,z}} \\[5pt] &{}+\underbrace {\int _\Omega f_c \overline {v}\, {\rm d}\Omega }_{\mathcal{C}_x} \underbrace {{}-\int _\Omega \frac {1}{\rho _o} \frac {\partial \overline {p}^\ast }{\partial x}\, {\rm d}\Omega }_{\mathcal{P}_x} \underbrace {{}-\int _\Omega \frac {1}{\rho _o} \frac {\partial \overline {p}_{\infty }}{\partial x}\, {\rm d}\Omega }_{\mathcal{P}_{\infty ,x}} + \underbrace {\int _\Omega \overline {f}_x\, {\rm d}\Omega }_{\mathcal{F}_{t,x}}= 0, \end{aligned} \end{equation}

where $u'=u-\bar {u}$ . The terms $\Delta \mathcal{A}_{uu,x}$ , $\Delta \mathcal{A}_{uv,y}$ and $\Delta \mathcal{A}_{uw,z}$ represent the advection of streamwise momentum by the streamwise, spanwise and vertical velocity. The terms $\Delta \mathcal{F}_{uu,x}$ , $\Delta \mathcal{F}_{uv,y}$ and $\Delta \mathcal{F}_{uw,x}$ denote the total (i.e. resolved + modelled) turbulent transport of streamwise momentum along the streamwise, spanwise and vertical directions. Further, $\mathcal{C}_x$ denotes the contribution of the Coriolis force, $\mathcal{P}_x$ and $\mathcal{P}_{\infty ,x}$ represent the contributions of the streamwise pressure gradients induced by the wind-farm forcing and the background pressure gradient, and $\mathcal{F}_{t,x}$ indicates the turbine forcing term. The latter term is zero out of the wind-farm region and therefore is not included in the analysis below. We note that due to the sign convention chosen, the first six terms in (4.2) are positive when the mean flow or turbulence transports more momentum in than out of the control volume, and negative when the opposite occurs. For example, the term $\Delta \mathcal{A}_{uu,x}$ is positive when $\bar {u} \bar {u}$ integrated over $\Gamma _{x_1}$ is higher than when integrated over $\Gamma _{x_2}$ , meaning that an increase in streamwise velocity along the $x$ direction causes this term to be negative.

Figure 18.

Streamwise variation of momentum sources and sinks for case H1000 normalised by the time-averaged total wind-farm thrust force and the control volume width: (a) all terms of (4.2), including $\mathcal{S}=\Delta \mathcal{A}_{uu,x} + \Delta \mathcal{A}_{uw,z} + \mathcal{P}_x$ ; (b) the contribution of each added in the $\mathcal{S}$ term. We note that only the region downwind of the farm is considered, so the term $\mathcal{F}_{t,x}$ has no contribution. The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

Figure 17(a) shows the streamwise evolution of all terms of (4.2) for case H300. The close vicinity of the capping inversion to turbine-tip height dampens the vertical turbulent entrainment of momentum. In fact, $\Delta \mathcal{F}_{uw,z}$ drops abruptly in the first $5$ km, attaining a similar value to that observed in the precursor simulation. This term contributes to the replenishment of the wake only in the near-wake region, where a fast recovery rate is observed – see e.g. Figure 12(c). Further, $\Delta \mathcal{F}_{uu,x}$ also attains a positive value in the near-wake region. However, this term is about three times smaller than $\Delta \mathcal{F}_{uw,z}$ , and it becomes negligible $5$ km downstream of the farm. The spanwise advection of mean flow $\Delta \mathcal{A}_{uv,y}$ is negligible in the first $10$ km of the wake, but becomes the dominant term contributing to the wake recovery in the far-wake region. However, flow advection along the spanwise direction is far less efficient than the vertical transfer of momentum in the wake recovery process. As a consequence, the velocity deficit shows a very slow recovery rate. Similarly to Bastankhah et al. (Reference Bastankhah, Mohammadi, Lees, Diaz, Buxton and Ivanell2024), we note a negligible contribution of the Coriolis term $\mathcal{C}_x$ , together with the spanwise turbulent transport of momentum. The background pressure gradient term $\mathcal{P}_\infty$ brings a positive contribution, although it speeds up the flow throughout the whole domain. Finally, we group terms that show an oscillatory behaviour into $\mathcal{S}=\Delta \mathcal{A}_{uu,x} + \Delta \mathcal{A}_{uw,z} + \mathcal{P}_x$ , which is the sum of mean advection along the streamwise and vertical directions and pressure perturbation. The contribution of each term is shown in figure 17(b). The strong oscillatory behaviour is due to the presence of trapped waves, which considerably alter the flow behaviour within the ABL when $H$ is low. In fact, the gravity-wave horizontal wavelength predicted using linear theory measures approximately $4.5$ km (see Vosper Reference Vosper2004; Sachsperger, Serafin & Grubišić Reference Sachsperger, Serafin and Grubišić2015), which is in line with the period of the oscillations shown in figure 17(b). We note that this mechanism in not visible in the results of Bastankhah et al. (Reference Bastankhah, Mohammadi, Lees, Diaz, Buxton and Ivanell2024) as they fixed the capping-inversion height. Here, we can also observe that $\Delta \mathcal{A}_{uu,x}$ is out of phase with respect to $\Delta \mathcal{A}_{uw,z}$ , while the latter is in phase with $\mathcal{P}_x$ . This behaviour can be explained through continuity. In fact, as the flow accelerates ( $\Delta \mathcal{A}_{uu,x}$ negative), it also moves downwards ( $\Delta \mathcal{A}_{uw,z}$ positive), and vice versa. Moreover, a positive $\mathcal{P}_x$ value corresponds to a favourable pressure gradient, which causes the flow to locally accelerate, therefore generating this out-of-phase behaviour with $\Delta \mathcal{A}_{uu,x}$ . The fact that both $\Delta \mathcal{A}_{uu,x}$ and $\Delta \mathcal{A}_{uw,z}$ do not attain a zero value in the far-wake region underlines that the wake does not fully recover.

Figure 19.

(a) Streamwise variation of $\Delta \mathcal{F}_{uw,z}$ (solid line) and $\Delta \mathcal{A}_{uv,y}$ (dashed line) normalised by the time-averaged total wind-farm thrust force of the respective cases and the control volume width. (b) Sum of the advection terms $\mathcal{A} = \Delta \mathcal{A}_{uu,x} + \Delta \mathcal{A}_{uv,y} + \Delta \mathcal{A}_{uw,z}$ normalised by the time-averaged total wind-farm thrust force of the respective cases and the control volume width. The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

The same analysis for case H1000 is shown in figure 18. The main difference from case H300 is the strong vertical turbulent entrainment of momentum in the first $20$ km of the wake, which causes the velocity deficit to rapidly recover. Since the wake is partially recovered after $20$ km, the advection of spanwise momentum by the streamwise velocity remains negligible. Similarly to case H300, $\Delta \mathcal{F}_{uu,x}$ , $\Delta \mathcal{F}_{uv,y}$ and $\mathcal{C}_x$ have negligible contributions. Figure 18(b) illustrates the components of the $\mathcal{S}$ term. Here, we can observe that the flow is accelerating and moving downwards (i.e. negative $\Delta \mathcal{A}_{uu,x}$ and positive $\Delta \mathcal{A}_{uw,z}$ ). Similarly to Bastankhah et al. (Reference Bastankhah, Mohammadi, Lees, Diaz, Buxton and Ivanell2024), we note that the strong flow acceleration is attained approximately where the vertical turbulent transport of momentum is also high, resulting in a strong correlation between the terms $\Delta \mathcal{A}_{uu,x}$ and $\Delta \mathcal{F}_{uw,z}$ . Moreover, although $\Delta \mathcal{A}_{uw,z}$ oscillates along the streamwise direction, its magnitude is on average slightly higher than $\Delta \mathcal{F}_{uw,z}$ , meaning that the vertical advection of mean flow plays an important role in the wake recovery process. We note that the H1000 case has the same oscillatory frequency as case H300. This is due to the fact that the horizontal trapped-wave wavelength depends on the capping-inversion strength and free-atmosphere lapse rate, which are equal for the two cases. We also observe the same in-phase and out-of-phase behaviour, although the magnitude of the oscillations caused by interfacial gravity waves is smaller than case H300 due to the higher distance between the turbine region and the capping-inversion height. Moreover, both $\Delta \mathcal{A}_{uu,x}$ and $\Delta \mathcal{A}_{uw,z}$ tend to zero at the end of the wake region, denoting that the wake has mostly recovered.

Excluding the terms heavily influenced by oscillation due to the presence of gravity waves, the two dominant terms that contribute to wake recovery are $\Delta \mathcal{F}_{uw,z}$ and $\Delta \mathcal{A}_{uv,y}$ , which are shown for all cases in figure 19(a). We can see that the vertical turbulent entrainment of momentum is directly proportional to $H$ . In fact, in deep-boundary-layer cases, the wake is mostly recovered within the first $20$ – $30$ km downwind of the farm, leaving a negligible role to the flow advection in the spanwise direction. Contrarily, $\Delta \mathcal{F}_{uw,z}$ dies out abruptly in shallow boundary layers, therefore causing very long wind-farm wakes, which are slowly replenished by the advection of spanwise momentum by the streamwise velocity. Interestingly, we note that the $\Delta \mathcal{F}_{uw,z}$ and $\Delta \mathcal{A}_{uv,y}$ profiles intersect in proximity of the zero line, except for case H150. Next, figure 19(b) shows $\mathcal{A} = \Delta \mathcal{A}_{uu,x} + \Delta \mathcal{A}_{uv,y} + \Delta \mathcal{A}_{uw,z}$ , i.e. the sum of the advection terms. We can observe oscillations that gradually decrease in amplitude along the streamwise direction, which denote a net acceleration and deceleration of the flow induced by trapped waves. In the far wake, $\mathcal{A}$ tends to a negative value for shallow-boundary-layer cases, denoting that the wake is still accelerating. This is in accordance with the results shown in figure 12(c). However, in deep-boundary-layer cases, $\mathcal{A}$ tends to small negative values in the far wake. In fact, the wake is mostly recovered by the vertical turbulent entrainment within the first 20–30 km downwind of the farm, leaving a negligible role to the flow advection. We remark that all cases share the same potential temperature profile above the ABL, resulting in in-phase profiles of $\mathcal{A}$ up to approximately 40 km downstream of the farm.

Figure 20.

Time-averaged advection of (a,b) spanwise momentum by the streamwise velocity and (c,d) total spanwise turbulent transport of momentum, as functions of the spanwise direction obtained at $x=x_0$ (brightest yellow) and $x=x_0+70$ (darkest blue), with increments of $5$ km in between, for cases (a,c) H300 and (b,d) H1000. Each symbol denotes a spanwise average over one turbine spacing (i.e. $5D$ ), while the vertical yellow and blue lines denote the location of the wake edges at $x=x_0$ and $x=x_0+70$ , respectively, where $x_0=30$ km indicates the beginning of the wake region.

To further investigate the spanwise advection of mean flow and spanwise turbulent transport of momentum, we analyse the quantities $\bar {u} \bar {v}$ and $\overline {u' v'}$ along the $y$ direction at hub height, taken at locations downwind of the farm (i.e. from $x=x_0$ to $x=x_0+70$ with increments of $5$ km). figure 20(a,b) display the spanwise advection of mean flow for cases H300 and H1000, respectively. In the shallow-boundary-layer case, we observe a high value of $\bar {u} \bar {v}$ at the location of the wake right edge (i.e. $y\approx 45$ km), which further increases along the downstream direction. Within the wake region, we observe a constant and positive spanwise advection that decreases significantly in proximity to the left edge of the wake (i.e. $y\approx 55$ km). In the H1000 case, the spanwise advection of mean flow attains a very similar magnitude at the wake edges, resulting in low $\Delta \mathcal{A}_{uv,y}$ values. In the far-wake region, we observe a higher spanwise advection of mean flow on the left edge. However, the wake deficit is almost fully recovered at these streamwise locations. Next, figure 20(c,d) illustrate the spanwise turbulent transport of momentum for cases H300 and H1000. In the shallow-boundary-layer case, this term is close to zero within the wake, attaining higher values in proximity of the left edge of the wake, where a steep spanwise gradient of streamwise velocity takes place – see figure 11(b). Since the latter increases in magnitude along the streamwise direction, we observe the peak in $\overline {u' v'}$ increasing as well. A similar behaviour is observed in the H1000 case, although higher fluctuations of $\overline {u' v'}$ are attained within the wake region. However, the role of the spanwise turbulent transport of momentum in the wake recovery process remains negligible since this term is approximately two orders of magnitude smaller than $\bar {u} \bar {v}$ .

4.4. Spanwise momentum budget analysis

In this analysis we adopt the same control volumes $\Omega$ used in § 4.3, but we investigate the budget terms for the spanwise momentum equation. The latter reads as

(4.3) \begin{equation} \begin{aligned} &\underbrace {-\bigg [ \int _{\Gamma _x} \bar {u}\bar {v}\, {\rm d}\Gamma _x\bigg ]_{x_1}^{x_2}}_{\Delta \mathcal{A}_{uv,x}} \underbrace {{}-\bigg [ \int _{\Gamma _y} \bar {v}\bar {v}\, {\rm d}\Gamma _y\bigg ]_{y_1}^{y_2}}_{\Delta \mathcal{A}_{vv,y}} \underbrace {{}-\bigg [ \int _{\Gamma _z} \bar {v}\bar {w}\, {\rm d}\Gamma _z\bigg ]_{z_1}^{z_2}}_{\Delta \mathcal{A}_{vw,z}} \\[5pt] &\underbrace {{}-\bigg [ \int _{\Gamma _x} \bigl (\overline {u'v'}+\overline {\tau }_{12}^{{sgs}} \bigl )\, {\rm d} \Gamma _x \bigg ]_{x_1}^{x_2}}_{\Delta \mathcal{F}_{uv,x}} \underbrace {-\bigg [ \int _{\Gamma _y} \bigl (\overline {v'v'}+\overline {\tau }_{22}^{{sgs}} \bigl )\, {\rm d} \Gamma _y \bigg ]_{y_1}^{y_2}}_{\Delta \mathcal{F}_{vv,y}} \underbrace {-\bigg [ \int _{\Gamma _z} \bigl (\overline {v'w'}+\overline {\tau }_{23}^{{sgs}} \bigl )\, {\rm d} \Gamma _z \bigg ]_{z_1}^{z_2}}_{\Delta \mathcal{F}_{vw,z}} \\[5pt] &\underbrace {{}-\int _\Omega f_c \overline {u}\, {\rm d}\Omega }_{\mathcal{C}_y} \underbrace {{}-\int _\Omega \frac {1}{\rho _o} \frac {\partial \overline {p}^\ast }{\partial y}\, {\rm d}\Omega }_{\mathcal{P}_y^\ast } \underbrace {{}-\int _\Omega \frac {1}{\rho _o} \frac {\partial \overline {p}_\infty }{\partial y} \,{\rm d}\Omega }_{\mathcal{P}_{\infty ,y}} \underbrace {{}+\int _\Omega \overline {f}_y\, {\rm d}\Omega }_{\mathcal{F}_{t,y}}= 0. \end{aligned} \end{equation}

The terms $\Delta \mathcal{A}_{uv,x}$ , $\Delta \mathcal{A}_{uv,y}$ and $\Delta \mathcal{A}_{vw,z}$ represent the advection of spanwise momentum by the streamwise, spanwise and vertical velocity. The terms $\Delta \mathcal{F}_{uv,x}$ , $\Delta \mathcal{F}_{vv,y}$ and $\Delta \mathcal{F}_{vw,z}$ denote the total (i.e. resolved + modelled) turbulent transport of spanwise momentum along the streamwise, spanwise and vertical directions. Further, $\mathcal{C}_y$ denotes the contribution of the Coriolis force, $\mathcal{P}_y$ and $\mathcal{P}_{\infty ,y}$ represent the contributions of the spanwise pressure gradients induced by the wind-farm forcing and the background pressure gradient, and $\mathcal{F}_{t,y}$ indicates the turbine forcing term. The latter term is zero out of the wind-farm region and therefore is not included in the analysis below. Moreover, we adopt the same sign convention as used in § 4.3.

Figure 21.

Spanwise variation of momentum sources and sinks for cases (a) H150, (b) H300, (c) H500 and (d) H1000 normalised by the time-averaged total wind-farm thrust force and the control volume width. Only the region downwind of the farm is considered, so the term $\mathcal{F}_{t,y}$ has no contribution. The $x$ -axis is rescaled with $x_0=30$ km, which denotes the beginning of the wake region.

Figure 21 illustrates the spanwise momentum sources and sinks for all cases. The spanwise background pressure gradient is related to $f_c G \cos\alpha$ and is therefore constant along the streamwise direction. The contribution of the Coriolis force is proportional to the streamwise velocity, and it increases in magnitude along the streamwise direction as a results of the wake recovery. In the free atmosphere, $\mathcal{P}_{\infty ,y} + \mathcal{C}_y = 0$ . However, the momentum sink generated by the farm decreases the magnitude of the Coriolis force, so that $\mathcal{P}_{\infty ,y} + \mathcal{C}_y \gt 0$ in the turbine region. This imbalance is compensated by $\mathcal{P}_y$ in cases H150, H300 and H500, which assumes a negative value denoting the presence of an opposing pressure gradient along the spanwise direction, as also observed in figure 7(m–o). The magnitude of this opposing spanwise pressure gradient is inversely related to $H$ and may explain the lack of growth in the mixing layer located on the left side of the farm (i.e. $\delta _l$ ), as shown in figure 12(d). The oscillatory behaviour in $\mathcal{P}_y$ is related to the presence of atmospheric gravity waves. In case H1000, $\mathcal{P}_y \approx 0$ and the imbalance between background pressure gradient and Coriolis force is compensated by the turbulent transport of spanwise momentum along the vertical direction. We note that the term $\Delta \mathcal{F}_{vw,z}$ is non-negligible also in case H500, but contributes only in the first 10 km downwind of the farm, while it becomes close to zero in shallow boundary layers due to the vicinity of the inversion layer to turbine-tip height. All other terms reported in figure 21 are negligible, assuming values of the order of the residual.