A pressure-driven atmospheric boundary layer model satisfying Rossby and Reynolds number similarity

Idealized models of the atmospheric boundary layer (ABL) can be used to leverage understanding of the interaction between the ABL and wind farms towards the improvement of wind farm flow modeling. We propose a pressure-driven one-dimensional ABL model without wind veer, which can be used as an inflow model for three-dimensional wind farm simulations to separately demonstrate the impact of wind veer and ABL depth. The model is derived from the horizontal momentum equations and follows both Rossby and Reynolds number similarity; use of such similarity reduces computation time and allows rational comparison between different conditions. The proposed ABL model compares well with solutions of the mean momentum equations that include wind veer if the forcing variable is employed as a free parameter.

1 Introduction

The interaction between the atmospheric boundary layer (ABL) and wind farms is important for wind energy, because it influences the energy yield and wind turbine lifetime. Many models of the ABL exist; these range from mesoscale models like the Weather Research and Forecasting (WRF) Model (Skamarock et al., 2019)¹ to microscale models such as large-eddy simulation (LES) (Stoll et al., 2020) and Reynolds-averaged Navier–Stokes (RANS) (Apsley and Castro, 1997; Blackadar, 1962; van der Laan et al., 2020b). LES is a transient method that resolves large-scale turbulence, while RANS is a steady-state method that models all turbulence scales. A RANS turbulence model that can handle all relevant turbulence scales does currently not exist. Despite this fact, RANS is our method of choice because it is roughly 3 orders of magnitude faster compared to LES, and it can be used to study trends of atmospheric wind farm flows. For example, it is possible to get good results of wind turbine wake losses in a wind farm subjected to a neutral atmospheric surface layer when using a RANS solver with modified two-equation turbulence models (van der Laan et al., 2015). In addition, RANS can simulate wind turbine interaction (van der Laan et al., 2015; Bleeg et al., 2018) – meaning both wake and blockage effects – which is not trivial for engineering wind farm flow models that rely on a predefined wake (and induction) shape and wake superposition. Finally, RANS can leverage the understanding of the interaction between the ABL and wind farms, because one can add or remove components of ABL physics (representing atmospheric stability, Coriolis forces, etc.) by including or deleting the corresponding physical terms in the RANS equations. For example, a RANS model can be used to model the effect of non-neutral atmospheric stability (van der Laan et al., 2020a). When the model operates in neutral mode, all model components that represent non-neutral conditions are switched off.

The use of higher-fidelity ABL inflow models in RANS for wind farm flows is a research area of both practical and academic interest. One can include the effects of surface layer atmospheric stability on a wind turbine wake using analytical profiles following Monin–Obukhov similarity theory (MOST, Monin and Obukhov, 1954), as shown by, e.g., Doubrawa et al. (2020), although it is not expected to be accurate for large wind turbines that operate outside the atmospheric surface layer. An idealized model of the ABL can also be employed in RANS; for example the two-equation turbulence model of Apsley and Castro (1997) includes Coriolis forces, a constant pressure gradient, and a turbulence length-scale limiter that determines the ABL depth. Neither MOST nor ABL inflow models require a temperature equation or a buoyancy contribution in the vector momentum equation; the effect of atmospheric stability can be modeled, e.g., by source terms in the turbulence model equations that only depend on velocity gradients (MOST) or via limitation of the turbulence length scale (ABL). An advantage of MOST is that for a prescribed (fixed) Obukhov length, the shape of the inflow profile is independent of wind speed. As a consequence, the simulated normalized wake losses in a wind farm using a MOST inflow follow Reynolds number similarity, as shown in van der Laan et al. (2020a). This is because the viscous forces can be neglected due to the high Reynolds number of atmospheric wind farm flows, and all external forces (wind turbine forces) scale as ${\mathcal{U}}^{\mathrm{2}}/\mathcal{L}$, with 𝒰 and ℒ as characteristic velocity and length scales, respectively. The wind speed independence can be exploited when calculating the wake effects of a wind farm for different wind speeds in a single wind farm simulation: different wind speed flow cases are run consecutively by scaling the wind turbine controller without changing the inflow profile, as shown in van der Laan et al. (2019). This method reduces the total number of iterations required to simulate multiple wind speed cases by a factor 2–3, because only local changes in the flow field need to be solved for since the global inflow is kept constant. The ABL inflow model does not follow Reynolds number similarity because the Coriolis force in the momentum equations scales linearly with 𝒰, instead of ${\mathcal{U}}^{\mathrm{2}}/\mathcal{L}$. However, the ABL inflow model does follow Rossby similarity, where the ABL profiles are only dependent on two Rossby numbers if the height z is normalized as $(z+z_{\mathrm{0}})f_{\mathrm{c}}/G$ or $(z+z_{\mathrm{0}})/z_{\mathrm{0}}$, with z₀ as the roughness length, f_c as the Coriolis parameter, and G as the geostrophic wind speed. The downside of Rossby similarity is that it cannot be used to speed up wind farm simulations as was done for MOST inflow profiles obeying Reynolds similarity, because the wind turbine size (hub height and rotor diameter) does not scale by $f_{\mathrm{c}}/G$ or z₀.

In this article, we present a new pressure-driven ABL model, which can be employed to follow both Reynolds and Rossby similarity. The ABL profiles of the proposed model are very similar to the ABL model of Apsley and Castro (1997) including ABL depth, by applying a momentum source term that represents the balance between Coriolis force and a fixed pressure gradient but without turning of the wind with height (veer). Several authors (Wilson et al., 1998; Parente et al., 2011; Cindori et al., 2020) have developed unidirectional atmospheric inflow models for two- or three-dimensional simulations of complex terrain, urban areas, and forests, mainly for the purpose of wind tunnel validation where the turbulent kinetic energy varies with height. These models do not include an ABL depth and could be interpreted as atmospheric surface layer (ASL) models. Our proposed ABL model can be used as an inflow model for three-dimensional wind farm simulations to isolate the effects of wind veer or ABL depth, when the results are compared with wind farm simulations using an inflow based on the ABL model of Apsley and Castro (1997) (with wind veer) or a neutral surface layer, respectively. Isolating the effect of wind veer can be of interest for wake steering control studies, where wind veer can have a significant impact as discussed by Brugger et al. (2020). In addition, it is possible to use the model to obtain Reynolds similarity and employ it as an inflow model for wind farm simulations where wind speed flow cases are simulated consecutively to reduce the number of required iterations. The pressure-driven ABL model is not equivalent to a common pressure-driven half-channel flow, because it includes a source term representing a balance between a constant pressure gradient and a type of Coriolis force, and the wind veer is removed by considering the scalar momentum equation in the direction of the mean wind (in effect swapping the U and V momentum source terms and applying the same sign). While the resulting equations may not initially appear to make sense physically, the simulated ABL profiles are very similar to the ABL profiles from the ABL model of Apsley and Castro (1997) where the correct equations are employed. Furthermore, the proposed ABL model has a physical and mathematical basis derived from a scalar momentum equation (in the mean wind direction), as shown in Sect. 3 (following similar ideas as seen in, e.g., Sogachev et al., 2005). Section 2 first gives an introduction of idealized ABL models in RANS, which is needed to understand the derivation of the pressure-driven ABL model. Rossby and Reynolds number similarity of the proposed ABL model is discussed and numerically demonstrated in Sects. 5 and 6, respectively. A comparison and application of the ABL models with and without veer is presented in Sect. 7.

2 Idealized ABL modeling in RANS

A steady-state idealized atmospheric boundary layer can be modeled by the incompressible RANS equations of momentum when considering homogeneous terrain and neglecting mesoscale effects:

$$\begin{array}{}\text{(1)}& \begin{array}{rl}& {\displaystyle \frac{\mathrm{D}U}{\mathrm{D}t}}=f_{\mathrm{c}}(V-V_{\mathrm{G}})+{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}U}{\mathrm{d}z}}\right)=\mathrm{0},\\ & {\displaystyle \frac{\mathrm{D}V}{\mathrm{D}t}}=-f_{\mathrm{c}}(U-U_{\mathrm{G}})+{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}V}{\mathrm{d}z}}\right)=\mathrm{0},\end{array}\end{array}$$

where U and V are the streamwise and lateral horizontal velocity components, U_G and V_G are the geostrophic velocities which represent constant pressure gradients, f_c is the Coriolis parameter dependent on latitude, z is the height, and t is the time. In addition, ν_T is the turbulent eddy viscosity, which is a result from employing the linear relationship of the Reynolds stresses and strain-rate tensor following Boussinesq (1897). The boundary conditions are $U=V=\mathrm{0}$ at z=z₀, with z₀ as the roughness length, and $\mathit{\{}U,V\mathit{\}}=\mathit{\{}{U}_{\mathrm{G}},V_{\mathrm{G}}\mathit{\}}$ for z→∞. Analytic solutions of Eq. (1) exist if the turbulent eddy viscosity ν_T is set as a constant (Ekman, 1905) or defined by a linearly increasing function with height (Ellison, 1956). Such solutions tend to not compare well with observations (e.g., Jensen et al., 1984; Hess and Garratt, 2002); this motivates the use of higher-fidelity turbulence models for the eddy viscosity. Blackadar (1962) applied the mixing-length model of Prandtl employing a prescribed turbulence length scale ℓ including a maximum ℓ_max:

$$\begin{array}{}\text{(2)}& {\mathit{\nu }}_T={\mathrm{\ell }}^{\mathrm{2}}\mathcal{S},\mathrm{\ell }={\displaystyle \frac{\mathit{\kappa }z}{\mathrm{1}+\frac{\mathit{\kappa }z}{{\mathrm{\ell }}_{max}}}},\end{array}$$

where κ is the von Kármán constant (we use κ=0.4). In addition, 𝒮 is the magnitude of the strain-rate tensor. For $\mathrm{\ell }\ll {\mathrm{\ell }}_{max}$, the neutral surface layer solution is obtained. The parameter ℓ_max is a proxy for ABL depth and can be used to model either neutral or stable atmospheric conditions, as discussed by Apsley and Castro (1997). For ${\mathrm{\ell }}_{max}\to \mathrm{\infty }$, neutral conditions are obtained that correspond to the analytic solution of Ellison (1956), as discussed in van der Laan et al. (2020b). For small values of ℓ_max, an ABL profile that has the characteristics of stable conditions is obtained – a shallow ABL, a strong shear and wind veer, and a small eddy viscosity, all with respect to neutral conditions. Hence, a potential temperature equation is not necessary, and the effects of a stable ABL are solely modeled by a limitation of the turbulence length scale. One can also model the eddy viscosity by a two-equation turbulence model including a turbulence length-scale limiter, for example the k–ε model of Apsley and Castro (1997):

$$\begin{array}{}\text{(3)}& \begin{array}{rl}& {\mathit{\nu }}_T=C_{\mathit{\mu }}{\displaystyle \frac{k^{\mathrm{2}}}{\mathit{\epsilon }}},{\displaystyle \frac{Dk}{\mathrm{D}t}}={\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\displaystyle \frac{{\mathit{\nu }}_T}{{\mathit{\sigma }}_k}}{\displaystyle \frac{\mathrm{d}k}{\mathrm{d}z}}\right)+\mathcal{P}-\mathit{\epsilon },\\ & {\displaystyle \frac{D\mathit{\epsilon }}{\mathrm{D}t}}={\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\displaystyle \frac{{\mathit{\nu }}_T}{{\mathit{\sigma }}_{\mathit{\epsilon }}}}{\displaystyle \frac{\mathrm{d}\mathit{\epsilon }}{\mathrm{d}z}}\right)\\ & +(\left[C_{\mathit{\epsilon },\mathrm{1}}+\left(C_{\mathit{\epsilon },\mathrm{2}}-C_{\mathit{\epsilon },\mathrm{1}}\right){\displaystyle \frac{\mathrm{\ell }}{{\mathrm{\ell }}_{max}}}\right]\mathcal{P}-C_{\mathit{\epsilon },\mathrm{2}}\mathit{\epsilon }){\displaystyle \frac{\mathit{\epsilon }}{k}},\end{array}\end{array}$$

where C_μ is a constant, k is the turbulent kinetic energy, and ε is its dissipation. Both k and ε are modeled by a transport equation, where 𝒫 is the mechanical production of turbulence; $\mathrm{\ell }=C_{\mathit{\mu }}^{\mathrm{3}/\mathrm{4}}{k}^{\mathrm{3}/\mathrm{2}}/\mathit{\epsilon }$ is the model-based local turbulence length scale; and C_ε,1, C_ε,2, σ_k, and σ_ε are model constants. The two-equation turbulence model also provides results of the turbulence intensity, which is not the case for the mixing-length model. In previous work (van der Laan et al., 2020b), we have shown that the analytic solutions of Ekman (1905) and Ellison (1956) are bounds of the mixing-length model of Blackadar (1962) and the two-equation model of Apsley and Castro (1997) for ${\mathrm{\ell }}_{max}\to \mathrm{0}$ and ${\mathrm{\ell }}_{max}\to \mathrm{\infty }$. In addition, the results of the numerical models follow a Rossby similarity, and all possible solutions of the ABL can be defined by two Rossby numbers based on different length scales (van der Laan et al., 2020b):

$$\begin{array}{}\text{(4)}& {\mathit{Ro}}_{\mathrm{0}}\equiv {\displaystyle \frac{G}{\left|f_{\mathrm{c}}\right|z_{\mathrm{0}}}},{\mathit{Ro}}_{\mathrm{\ell }}\equiv {\displaystyle \frac{G}{\left|f_{\mathrm{c}}\right|{\mathrm{\ell }}_{max}}}.\end{array}$$

Here, Ro₀ is the well-known surface Rossby number, and Ro_ℓ is a Rossby number based on the maximum turbulence length scale. The Rossby similarity applies to the normalized ABL profiles, where the height z is normalized as $(z+z_{\mathrm{0}})f_{\mathrm{c}}/G$ and the flow variables are normalized by G and ℓ_max. The two-equation turbulence model of Apsley and Castro (1997) is only applicable to flat terrain, but it can be used as an inflow model for atmospheric wind farm flows in homogeneous terrain and roughness, as performed in previous work (van der Laan and Sørensen, 2017b).

3 A pressure-driven model of the ABL without wind veer

Our goal is to develop a pressure-driven one-dimensional model of the idealized ABL in terms of wind speed but without wind veer. One can derive such a model by combining the momentum equations of U and V, i.e., Eq. (1), and rewriting them as a single equation in terms of the magnitude of geostrophic deficit (Wyngaard, 2010): $\widehat{S}\equiv \sqrt{(U-U_{\mathrm{G}}{)}^{\mathrm{2}}+(V-V_{\mathrm{G}}{)}^{\mathrm{2}}}$. Sogachev et al. (2005) also derived a momentum equation of wind speed using $S\equiv \sqrt{U^{\mathrm{2}}+V^{\mathrm{2}}}$ for 2-D flows, but we will use normalized velocity variables for derivation and transform the final result back to the common velocity variables. An equation for $\widehat{S}$ can be derived in a number of ways. While a textbook method is to write Eq. (1) in the complex form (Wyngaard, 2010)

$$\begin{array}{}\text{(5)}& {\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}W}{\mathrm{d}z}}\right)=i{f}_{\mathrm{c}}\left(W-W_{\mathrm{G}}\right),\end{array}$$

where $W\equiv U+iV$, $W_{\mathrm{G}}\equiv U_{\mathrm{G}}+i{V}_{\mathrm{G}}$, and $i^{\mathrm{2}}=-\mathrm{1}$, we instead use the components in order to keep our result clear. Taking the sum of $(\widehat{U}+\widehat{V})\mathrm{D}U/\mathrm{D}t$ and $(-\widehat{U}+\widehat{V})\mathrm{D}V/\mathrm{D}t$ from Eq. (1), using the normalized variables $\widehat{U}\equiv U-U_{\mathrm{G}}$ and $\widehat{V}\equiv V-V_{\mathrm{G}}$,² and defining the wind direction as $\widehat{\mathit{\phi }}=\arctan (\widehat{V}/\widehat{U})$, we get

$$\begin{array}{}\text{(6)}& \begin{array}{rl}\mathrm{0}& =f_{\mathrm{c}}\left({\widehat{U}}^{\mathrm{2}}+{\widehat{V}}^{\mathrm{2}}\right)+\widehat{U}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}\right)+\widehat{V}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}\right)\\ & -\widehat{U}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}\right)+\widehat{V}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}\right),\\ & =f_{\mathrm{c}}\left({\widehat{U}}^{\mathrm{2}}+{\widehat{V}}^{\mathrm{2}}\right)+{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T\left[\widehat{U}{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}+\widehat{V}{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}\right]\right)\\ & -{\mathit{\nu }}_T\left({\left[{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}\right]}^{\mathrm{2}}+{\left[{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}\right]}^{\mathrm{2}}\right)-{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T\left[\widehat{U}{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}-\widehat{V}{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}\right]\right),\\ & =f_{\mathrm{c}}{\widehat{S}}^{\mathrm{2}}+{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\widehat{S}}^{\mathrm{2}}\left[{\displaystyle \frac{\mathrm{1}}{\widehat{S}}}{\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}-{\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right]\right)\\ & -{\mathit{\nu }}_T{\widehat{S}}^{\mathrm{2}}\left({\left[{\displaystyle \frac{\mathrm{1}}{\widehat{S}}}{\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}\right]}^{\mathrm{2}}+{\left[{\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right]}^{\mathrm{2}}\right),\\ & =f_{\mathrm{c}}\widehat{S}+{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}\right)-{\mathit{\nu }}_T\widehat{S}{\left({\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right)}^{\mathrm{2}}\\ & -{\displaystyle \frac{\mathrm{1}}{\widehat{S}}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\widehat{S}}^{\mathrm{2}}{\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right).\end{array}\end{array}$$

Here, we have applied the chain and product rules of differentiation, assumed a zero geostrophic shear ($\mathrm{d}G/\mathrm{d}z=\mathrm{0}$), and the following relations (Kelly and van der Laan, 2021) for the wind veer $\mathrm{d}\widehat{\mathit{\phi }}/\mathrm{d}z$ and wind shear $\mathrm{d}\widehat{S}/\mathrm{d}z$ are employed:

$$\begin{array}{}\text{(7)}& \begin{array}{rl}{\widehat{S}}^{\mathrm{2}}{\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}& =\widehat{U}{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}-\widehat{V}{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}},\widehat{S}{\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}=\widehat{U}{\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}+\widehat{V}{\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}},\\ & {\left({\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}\right)}^{\mathrm{2}}+{\left(\widehat{S}{\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right)}^{\mathrm{2}}={\left({\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}\right)}^{\mathrm{2}}+{\left({\displaystyle \frac{\mathrm{d}\widehat{V}}{\mathrm{d}z}}\right)}^{\mathrm{2}}.\end{array}\end{array}$$

If we take the wind veer to be much less than $(\mathrm{1}/\widehat{S})\mathrm{d}\widehat{S}/\mathrm{d}z$ in Eq. (6) then we recover an equation for geostrophic deficit $\widehat{S}$, which looks identical to Eq. (5) for the magnitude of the ageostrophic wind vector (again assuming that the geostrophic shear $\mathrm{d}G/\mathrm{d}z$ is zero):

$$\begin{array}{}\text{(8)}& {\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}\right)=-f_{\mathrm{c}}\widehat{S}.\end{array}$$

The assumption $\mathrm{d}\widehat{\mathit{\phi }}/\mathrm{d}z\ll (\mathrm{1}/\widehat{S})\mathrm{d}\widehat{S}/\mathrm{d}z$ could be considered a weak assumption since it in principle allows for some veer. However, for some cases this is violated; e.g., in the Ekman (1905) solution (constant eddy viscosity) $\mathrm{d}\widehat{\mathit{\phi }}/\mathrm{d}z$ equals $(\mathrm{1}/\widehat{S})\mathrm{d}\widehat{S}/\mathrm{d}z$. A stronger and simpler assumption, which also leads to Eq. (8), is to simply take $\mathrm{d}\widehat{\mathit{\phi }}/\mathrm{d}z=\mathrm{0}$ in Eq. (6). Note that neglecting veer gives $\widehat{S}=|S-G|$, so then we are not really dealing with a Coriolis force, per se. In addition, solving Eq. (8) will result in a solution for S where its magnitude cannot be larger than G for all z (which will be further elucidated at the end of this section), and we can use $\widehat{S}=-(S-G)$. Thus, we rewrite Eq. (8) for the wind speed as

$$\begin{array}{}\text{(9)}& {\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}S}{\mathrm{d}z}}\right)=f_{\mathrm{pg}}\left(S-G\right),\end{array}$$

where we have replaced f_c by f_pg; in lieu of f_c(S−G), the replacement f_pg(S−G) represents the magnitude of the pressure gradient and Coriolis effects. One needs to employ a different f_pg than f_c, in order to get similar ABL profiles of wind speed, turbulence intensity, and turbulence length scale, when disregarding the Coriolis-induced wind veer; this will be shown numerically in Sect. 7 and further explained below at the end of this section. Upon neglecting the wind veer, $\mathrm{D}U/\mathrm{D}t$ and $\mathrm{D}V/\mathrm{D}t$ become decoupled; then we can write the momentum equations by taking the product of Eq. (9) and cos (φ) or sin (φ), since φ is constant ($\mathrm{d}\mathit{\phi }/\mathrm{d}z=\mathrm{0}$). Preserving the relationship between magnitudes as evoked by Eqs. (8) and (9), we then have a 1-D pressure-driven ABL model:

$$\begin{array}{}\text{(10)}& \begin{array}{rl}& {\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}U}{\mathrm{d}z}}\right)=f_{\mathrm{pg}}\left(U-U_{\mathrm{G}}\right),\\ & {\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\displaystyle \frac{\mathrm{d}V}{\mathrm{d}z}}\right)=f_{\mathrm{pg}}\left(V-V_{\mathrm{G}}\right).\end{array}\end{array}$$

Equation (10) is the basis of the proposed ABL model without wind veer. It is the same as the original set of momentum equations that describe an idealized ABL, Eq. (1); however, through solving the scalar (decoupled in x and y) equation, the source terms are swapped and have the same sign. Furthermore, one could interpret the forcing of the ABL model as a pressure gradient, hence the subscript pg in f_pg. However, the ABL model does include a balance between a fixed pressure gradient and a type of Coriolis force, while a standard pressure-driven ABL model does not.

With the neglect of $\mathrm{d}\widehat{\mathit{\phi }}/\mathrm{d}z$, Eq. (7) implies

$$\begin{array}{}\text{(11)}& {\displaystyle \frac{\mathrm{d}{\widehat{S}}^{\mathrm{2}}}{\mathrm{d}z}}={\displaystyle \frac{\mathrm{d}{\widehat{U}}^{\mathrm{2}}}{\mathrm{d}z}}\left(\mathrm{1}+{\displaystyle \frac{{\widehat{V}}^{\mathrm{2}}}{{\widehat{U}}^{\mathrm{2}}}}\right)\text{and}{\displaystyle \frac{\mathrm{d}\widehat{S}}{\mathrm{d}z}}={\displaystyle \frac{\mathrm{d}\widehat{U}}{\mathrm{d}z}}{\left(\mathrm{1}+{\displaystyle \frac{{\widehat{V}}^{\mathrm{2}}}{{\widehat{U}}^{\mathrm{2}}}}\right)}^{\mathrm{1}/\mathrm{2}}.\end{array}$$

This can be seen as an approximation whereby the minor effect of lateral winds provides a perturbation to the streamwise gradients, when considering the full wind shear $\mathrm{d}S/\mathrm{d}z$ and gradient of mean kinetic energy $\mathrm{d}{S}^{\mathrm{2}}/\mathrm{d}z$. Ghannam and Bou-Zeid (2020) considered the effect of veer on ABL profiles and derived an approximate model for such; the neglected terms in Eq. (6), with the above equation, can be compared to magnitudes implied by their model.

Figure 1 Analytic solutions with and without wind veer. (a, b) Constant ν_T. (c, d) Linear ν_T and Ro₀=10⁵. (a, c) Wind speed. (b, d) Ageostrophic wind speed.

Analytic solutions of the wind speed profile can be derived from Eq. (10) using a constant or a linearly increasing eddy viscosity similarly to Ekman (1905) and Ellison (1956), respectively, using the original equation including wind veer, Eq. (5). The constant and linear eddy viscosity solutions of the decoupled “veerless” ABL model become

$$\begin{array}{}\text{(12)}& \begin{array}{ll}\mathrm{Constant}{\mathit{\nu }}_T& S\left(\mathit{\xi }\right)=G\left[\mathrm{1}-\exp \left(-\mathit{\xi }\right)\right],\\ \mathrm{Linear}{\mathit{\nu }}_T=\mathit{\kappa }{u}_{*\mathrm{0}}z& S\left(\mathit{\eta }\right)=G\left[\mathrm{1}-c{K}_{\mathrm{0}}\left(\mathit{\eta }\right)\right],\end{array}\end{array}$$

with $\mathit{\xi }=z\sqrt{f_{\mathrm{pg}}/{\mathit{\nu }}_T}$ and $\mathit{\eta }=\mathrm{2}\sqrt{z{f}_{\mathrm{pg}}/\left(\mathit{\kappa }{u}_{*\mathrm{0}}\right)}$ as normalized heights, K₀ as the zero-order modified Bessel function of the second kind, $c=\mathrm{2}{u}_{*\mathrm{0}}/\left(\mathit{\kappa }G\right)=-{\left[{\mathit{\gamma }}_e+\frac{\mathrm{1}}{\mathrm{2}}\ln \left(z_{\mathrm{0}}{f}_{\mathrm{pg}}/\left(\mathit{\kappa }{u}_{*\mathrm{0}}\right)\right)\right]}^{-\mathrm{1}}$ as a constant, u_*0 as the friction velocity at the surface, and γ_e as the Euler–Mascheroni constant. The derivation of the constant ν_T solution in Eq. (12) is identical to the classical textbook Ekman (1905) solution (e.g., Wyngaard, 2010), taking $i\left|f_{\mathrm{c}}\right|\to f_{\mathrm{pg}}$, which also indicates that $\sqrt{f_{\mathrm{pg}}}=\mathrm{\Re }\mathit{\left\{}\sqrt{i\left|f_{\mathrm{c}}\right|}\mathit{\right\}}=\sqrt{\left|f_{\mathrm{c}}\right|/\mathrm{2}}$. Note that the friction velocity in the linear ν_T solutions with and without wind veer is solved from an implicit relation of $u_{*\mathrm{0}}/G$, through the constant c derived both via $\mathrm{d}S/\mathrm{d}z=u_{\mathrm{0}*}/\left(\mathit{\kappa }z\right)$ and S=0 for z→z₀. The latter is effectively a form of the geostrophic drag law, which generally arises when G is used as a boundary condition in wind profile forms which include the surface stress and z₀ (Kelly and Troen, 2016). The solution including veer is further discussed in van der Laan et al. (2020b), also based on Ellison (1956) and Krishna (1980).

The analytical solutions corresponding to a constant and a linear ν_T are depicted in Fig. 1, where results of the wind speed, S, and the ageostrophic wind speed, $\widehat{S}$, are plotted. One solution with wind veer and two solutions without wind veer are shown using $f_{\mathrm{pg}}=\left|f_{\mathrm{c}}\right|$ and $f_{\mathrm{pg}}=\left|f_{\mathrm{c}}\right|/\mathrm{2}$. We have chosen to plot the results against $z\sqrt{\left|f_{\mathrm{c}}\right|/{\mathit{\nu }}_T}$ and $\mathrm{2}\sqrt{z\left|f_{\mathrm{c}}\right|/\left(\mathit{\kappa }{u}_{*\mathrm{0}}\right)}$ to depict the differences between different values of f_pg. The latter would not be clear if the results are plotted against ξ, since the results for different f_pg would collapse for the Ekman solutions without wind veer. When $f_{\mathrm{pg}}=\left|f_{\mathrm{c}}\right|/\mathrm{2}$, the ageostrophic wind speed for both Ekman solutions (with and without wind veer) is the exactly the same (Fig. 1b): $\widehat{S}=-G\exp (-\mathit{\xi })$, while the wind speed compares better if $f_{\mathrm{pg}}=\left|f_{\mathrm{c}}\right|$ (Fig. 1a). This result may seem counter intuitive, but it simply follows from the analytic solution of Eq. (12). A similar conclusion can be made for the linear ν_T solutions with and without wind veer (Fig. 1c and d), although $\widehat{S}$ is not exactly the same when $f_{\mathrm{pg}}=\left|f_{\mathrm{c}}\right|/\mathrm{2}$. Figure 1 clearly shows that the wind speed of the analytic solutions without wind veer (Eq. 12) cannot exceed the geostrophic wind speed. We also find this for the higher-fidelity turbulence model closures since their solutions are bounded by the two analytic solutions (van der Laan et al., 2020b). The ABL model with wind veer includes the supergeostrophic wind speed (jet), which typically occurs below the ABL top, as predicted by the Ekman equations (Blackadar, 1957). The jet is a consequence of the Coriolis-induced interaction of alternate horizontal momentum and stress components; it does not exist in an idealized ABL model when the wind veer – and more importantly the oscillating part of the solution (which results from the coupling of the u and v equations) – is removed. This is not a new insight because one could deduce it from a text book (e.g., chap. 10 of Wyngaard, 2010) employing Eq. (5); however, we have shown the relation between the Coriolis-induced wind veer and the jet more explicitly by deriving two analytic solutions without wind veer (Eq. 12), as depicted in Fig. 1.

4 Methodology of numerical simulations

The methodology of the numerical one-dimensional simulations of the present article is very similar to that performed in van der Laan et al. (2020b), and a brief summary is presented here. The RANS simulations are carried out with a one-dimensional version of EllipSys (van der Laan and Sørensen, 2017a), which is an in-house incompressible finite-volume flow solver initially developed by Michelsen (1992) and Sørensen (1994). The numerical grid represents a 10⁵ m line with 384 cells that increase with height using an expansion ratio of 1.2 and a first cell height of 10⁻² m. The number of cells is a conservative choice based on a grid refinement study, as performed in previous work (van der Laan et al., 2020b). The bottom and top boundary conditions are set as a rough wall (Sørensen et al., 2007) and symmetry boundaries, respectively. Ambient source terms in transport equations of the k–ε model are employed in order to prevent zero values (van der Laan et al., 2020b). The RANS simulations are solved transiently with a fixed large time step set to $\mathrm{1}/\left|f_{\mathrm{c}}\right|$ or $\mathrm{1}/f_{\mathrm{pg}}$ s and converge to a steady-state solution. The following turbulence model constants are employed: $(C_{\mathit{\mu }},C_{\mathit{\epsilon },\mathrm{1}}$, $C_{\mathit{\epsilon },\mathrm{2}},{\mathit{\sigma }}_k,{\mathit{\sigma }}_{\mathit{\epsilon }},\mathit{\kappa })=(\mathrm{0.03},\mathrm{1.21},\mathrm{1.92},\mathrm{1.0},\mathrm{1.3},\mathrm{0.4})$ following Sørensen (1994).

Figure 2 Rossby number similarity of the proposed ABL model without wind veer for two turbulence closures. (a–c) Mixing-length model. (d–g) k–ε model.

5 Rossby number similarity

The ABL model including wind veer follows a Rossby similarity as shown in previous work (van der Laan et al., 2020b). As a consequence, all possible normalized solutions of the ABL are only dependent on two Rossby numbers, each with a different length scale (see Eq. 4). The proposed ABL model without wind veer, as derived in Sect. 3, also follows a Rossby similarity. One can show this by writing the equation of wind speed, Eq. (9), in non-dimensional form:

$$\begin{array}{}\text{(13)}& {\displaystyle }{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}{z}^{\prime }}}\left({\left[{\displaystyle \frac{\mathit{\kappa }{z}^{\prime }}{\mathrm{1}+\mathit{\kappa }{z}^{\prime }{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}/{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}}}\right]}^{\mathrm{2}}{\left({\displaystyle \frac{\mathrm{d}{S}^{\prime }}{\mathrm{d}{z}^{\prime }}}\right)}^{\mathrm{2}}\right)=S^{\prime }-\mathrm{1},\end{array}$$

where $z^{\prime }=z/z_{\mathrm{0}}$, $S^{\prime }=S/G$, and ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}/{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}=z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$. Here, we have used the mixing-length turbulence model and prescribed the turbulence length scale of Blackadar (1962) (Eq. 2), but the same Rossby similarity applies to the two-equation turbulence model of Apsley and Castro (1997) (Eq. 3). In addition, the two Rossby numbers are defined as

$$\begin{array}{}\text{(14)}& {\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}\equiv {\displaystyle \frac{G}{f_{\mathrm{pg}}{z}_{\mathrm{0}}}},{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}\equiv {\displaystyle \frac{G}{f_{\mathrm{pg}}{\mathrm{\ell }}_{max}}}.\end{array}$$

A numerical proof of the Rossby similarity of the ABL model without wind veer for both the mixing-length and two-equation turbulence models is depicted in Fig. 2. Four sets of Rossby numbers are used, and each set is simulated by four cases using two different values of f_pg and G. Normalized results of wind speed, wind direction, eddy viscosity, and turbulence intensity, $I\equiv \sqrt{\mathrm{2}/\mathrm{3}k}/S$, collapse and are only dependent on the two Rossby numbers.

Figure 3 Reynolds number similarity of the proposed ABL model without wind veer for two turbulence closures. (a–c) Mixing-length model. (d–g) k–ε model.

6 Reynolds number similarity

If all external forces in the momentum equations scale by ${\mathcal{U}}^{\mathrm{2}}/\mathcal{L}$, with 𝒰 and ℒ as characteristic velocity and length scales, respectively, then one can obtain Reynolds number similarity. The Reynolds number can be defined as

$$\begin{array}{}\text{(15)}& \mathit{Re}\equiv {\displaystyle \frac{\mathcal{U}\mathcal{L}}{\mathit{\nu }}}={\displaystyle \frac{f_{\mathrm{c}}{\mathcal{L}}^{\mathrm{2}}}{\mathit{\nu }}}\mathit{Ro},\end{array}$$

with ν as the molecular viscosity. The Reynolds number can be related to the Rossby number $\mathit{Ro}\equiv \mathcal{U}/\left(f_{\mathrm{c}}\mathcal{L}\right)$. Wind turbine wake simulations where only pressure rotor forces are considered (and viscous rotor forces are neglected) follow Reynolds number similarity if the inflow model does as well (van der Laan et al., 2020a). For high Reynolds numbers, the wind turbine wake simulations then become independent of inflow wind speed and wind turbine size (rotor diameter D and hub height z_H, as long as $D/z_{\mathrm{H}}$ is kept constant), which can be employed to reduce the total number of required iterations of parametric studies or annual energy calculations of wind farms using RANS (van der Laan et al., 2019). The Coriolis force in the ABL model with wind veer scales by 𝒰, which means that Reynolds number similarity cannot be obtained. However, the f_pg parameter in the pressure-driven ABL model can be redefined to obtain Reynolds number similarity:

$$\begin{array}{}\text{(16)}& f_{\mathrm{pg}}\equiv C{\displaystyle \frac{\mathcal{U}}{\mathcal{L}}}={\displaystyle \frac{\mathrm{1}}{{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}}}{\displaystyle \frac{G}{z_{\mathrm{0}}}},\end{array}$$

with C as a constant, which is equivalent to setting ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}=\mathrm{1}/C$ and ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}=z_{\mathrm{0}}/{\mathrm{\ell }}_{max}/C$. We can allow ourselves to redefine f_pg because it does not directly represent a physical Coriolis parameter as f_c does in the ABL model with wind veer. Substitution of Eq. (16) and employing the mixing-length model of Eq. (2) in Eq. (9) lead to same results as given in Eq. (13). For a constant ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$, the only parameter that changes the ABL profile shape is $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$. Hence, the ABL model becomes wind speed independent, and one can use it as an inflow model for Reynolds-number-independent wind turbine wake simulations for a fixed $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$, which represents a fixed ABL depth or a fixed atmospheric stability. Note that a fixed atmospheric stability refers to a set value of a stability parameter, as opposed to a calculated stability condition that one could obtain from an ABL solution, which is typically the case for a transient ABL model including a potential temperature equation, as shown by Sogachev et al. (2012). The Reynolds number independence is shown in Fig. 3 for two turbulence closures. Two values of ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ are used, representing onshore and offshore roughness for latitudes around $\pm {\mathrm{45}}^{\circ }$. Each Rossby number is simulated with two different values of $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$. Note that we have chosen to use a different set of $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$ for each Rossby number because the range of meaningful $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$ is dependent on the choice of ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$. The four resulting ABL cases are then simulated with two geostrophic wind speeds and roughness lengths. Figure 3 shows that the normalized ABL profiles are independent of G, as long as ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ and $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$ are kept constant. Figure 3 can also be interpreted as Rossby number similarity as depicted in Fig. 2; however, the difference is that f_pg in Fig. 3 is now used to keep ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ constant for different values of G and z₀. In addition, it is clear that $z_{\mathrm{0}}/{\mathrm{\ell }}_{max}$ is a proxy for a normalized ABL depth, when comparing the profiles in pairs for a constant ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ (${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}={\mathrm{10}}^{\mathrm{6}},z_{\mathrm{0}}/{\mathrm{\ell }}_{max}={\mathrm{10}}^{-\mathrm{1}}$ compared to ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}={\mathrm{10}}^{\mathrm{6}},z_{\mathrm{0}}/{\mathrm{\ell }}_{max}={\mathrm{10}}^{-\mathrm{3}}$ and ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}={\mathrm{10}}^{\mathrm{9}},z_{\mathrm{0}}/{\mathrm{\ell }}_{max}={\mathrm{10}}^{-\mathrm{4}}$ compared to ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}={\mathrm{10}}^{\mathrm{9}},z_{\mathrm{0}}/{\mathrm{\ell }}_{max}={\mathrm{10}}^{-\mathrm{6}}$).

Figure 4 Comparison of k–ε ABL models (with and without wind veer) for stable (${\mathrm{\ell }}_{max}=\mathrm{5}$ m) and neutral (${\mathrm{\ell }}_{max}=\mathrm{30}$ m) conditions offshore using G=10 ms⁻¹, $f_{\mathrm{c}}=\mathrm{2}{f}_{\mathrm{pg}}={\mathrm{10}}^{-\mathrm{4}}$ s⁻¹ and $z_{\mathrm{0}}={\mathrm{10}}^{-\mathrm{4}}$ m. Bottom plots (e–h) are a magnified view of the top plots (a–d), focused around a 5 MW wind turbine rotor area, which is depicted as black dashed lines. (a, e) Wind speed. (b, f) Wind direction. (c, g) Turbulence intensity. (d, h) Turbulence length scale.

7 Comparison of ABL models and application to inflow profiles

In this section, one-dimensional RANS simulations are performed to compare the proposed ABL model without wind veer to the ABL model including wind veer, and we investigate the application to use the model as an inflow model. Three-dimensional RANS simulations are not performed in this article and will be carried out in future work. In addition, the k–ε ABL model of Apsley and Castro (1997) is used because it also provides an estimate of the turbulence intensity, while the mixing-length model of Blackadar (1962) does not.

Figure 4 shows a comparison between the original k–ε ABL model with wind veer and the proposed k–ε ABL without wind veer. Results of a stable (${\mathrm{\ell }}_{max}=\mathrm{5}$ m) and a neutral ABL (${\mathrm{\ell }}_{max}=\mathrm{30}$ m) are depicted for an offshore roughness length of 10⁻⁴ m, a geostrophic wind speed of 10 ms^-1, and a Coriolis parameter of 10⁻⁴ s⁻¹. In addition, the parameter f_pg of the ABL model without veer is set as $f_{\mathrm{pg}}=\left|f_{\mathrm{c}}\right|/\mathrm{2}$, as suggested by the classic Ekman (1905) solution discussed in Sect. 3. The ABL profiles of Fig 4 could be used as inflow profiles for offshore wind farm simulations, and we have chosen to normalize the results by a reference height, z_ref=90 m, which corresponds to the hub height of the NREL 5 MW reference wind turbine (Jonkman et al., 2009). In addition, the swept rotor area is marked as dashed black lines in Fig. 4, and the bottom plots are magnified views of the top plots, with a focus on the ABL around the fictitious wind turbine. Figure 4a shows the wind speed, and it is clear that the supergeostrophic jet (Blackadar, 1957) is a feature that cannot be predicted by the ABL model without wind veer, as also found for the analytic solutions from Eq. (12). As a consequence, the wind shear in the pressure-driven ABL model is smaller than the original ABL model including wind veer, as seen in Fig. 4e. Figure 4b and f show that the proposed ABL model predicts a zero wind veer as intended. Finally, the turbulence intensity and length scale are similar between the ABL models around the wind turbine rotor area but are different for higher altitudes, as shown in Fig. 4c and d.

It is possible to use f_pg as a free parameter in the pressure-driven ABL model, for a given height, to give results approximating those from the original ABL model considering veer around a reference height. One can justify this choice because f_pg does not simply represent the Coriolis parameter f_c; while it is a proxy for f_c with use of the scalar wind speed equation, it also contains the effects of the neglected wind veer, since we can write

$$\begin{array}{}\text{(17)}& {\displaystyle }{f}_{\mathrm{pg}}=f_{\mathrm{c}}-{\mathit{\nu }}_T{\left({\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right)}^{\mathrm{2}}-{\displaystyle \frac{\mathrm{1}}{{\widehat{S}}^{\mathrm{2}}}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}z}}\left({\mathit{\nu }}_T{\widehat{S}}^{\mathrm{2}}{\displaystyle \frac{\mathrm{d}\widehat{\mathit{\phi }}}{\mathrm{d}z}}\right).\end{array}$$

Here, we have set Eq. (9) equal to the final result of Eq. (6). It is not trivial to solve for Eq. (17), and we choose to obtain f_pg from a library of pre-calculated ABL profiles, which are only dependent on the two Rossby numbers ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ and ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}$. For wind farm simulations, one would like to obtain an inlet ABL profile for a desired reference wind speed S_ref and turbulence intensity I_ref, specified at a reference height z_ref for a given site where z₀ and f_c are known. In Appendix A, a procedure is presented for how to obtain a desired ABL profile using pre-calculated libraries of normalized ABL profiles based on Rossby number similarity. Examples of neutral and stable ABL profiles are made by using the chosen and derived values listed in Table 1, and the results of both ABL models are depicted in Fig. 5. Here, we have used f_pg as a free parameter, and the geostrophic wind speed is derived differently for both ABL models. In general, we find that f_pg<f_c (as discussed previously) and G_pg>G. A higher geostrophic wind speed is required in the ABL model without wind veer to compensate for a reduced wind shear at z=z_ref caused by the lack of the supergeostrophic jet. Figure 5 shows that a close match between the ABL models can be achieved in terms of wind speed, turbulence intensity, and turbulence length scale, especially around the rotor area. Hence, the pressure-driven ABL model can be used as an inflow model in wind farm simulations to isolate the effect of wind veer when it is compared to an inflow model based on the original ABL model including wind veer. The main differences between the ABL model can be found in terms of wind speed near the supergeostrophic jet and above. For very stable conditions, the location of the jet could approach the wind turbine rotor, and the effect of wind veer cannot be isolated from the effect of the jet unless one considers the jet to be part of the effect of wind veer (as previously described).

Figure 5 Comparison of k–ε ABL models (with and without wind veer) for stable and neutral conditions offshore using f_c≠2f_pg and parameters from Table 1. Bottom plots (e-h) are a magnified view of the top plots (a–d), focused around a 5 MW wind turbine rotor area, which is depicted as black dashed lines. (a, e) Wind speed. (b, f) Wind direction. (c, g) Turbulence intensity. (d, h) Turbulence length scale.

Table 1 Summary of input and derived parameters for ABL models.

In previous work, the k–ε ABL model with wind veer has been applied as an inflow to wind farm simulations using RANS to investigate the effect of Coriolis forces (van der Laan and Sørensen, 2017b), and the k–ε ABL was coupled with a k–ε developed for wake simulations under neutral ASL conditions. A similar coupling could be made with the proposed ABL model without wind veer, which will be investigated in future work.

It should be noted that both k–ε ABL models (with and without wind veer) cannot be used as an inflow model for complex terrain simulations in RANS because the employed global length-scale limiter of Apsley and Castro (1997) does not perform well when the turbulence length scale associated with the terrain is larger than ℓ_max. We plan to modify the length-scale limiter in future work to overcome this issue. However, the proposed momentum source from Eq. (10) can still be used in combination with an alternative turbulence model suited for complex terrain.

8 Conclusions

We have proposed a pressure-driven model of the mean ABL without wind veer, based on the streamwise (scalar) momentum equation. One-dimensional RANS simulations of the pressure-driven ABL model are performed to show that the model follows both Rossby and Reynolds number similarity. The similarities can be employed to quickly find a desired ABL profile based on a pre-calculated library of ABL profiles, which can be used as an inflow profile for three-dimensional RANS simulations of wind farms. The pressure-driven ABL model compares well with an ABL model including wind veer if the forcing variable f_pg is used a free parameter. The largest differences between the models are found near the location of the supergeostrophic jet and above (i.e., around the top of the ABL), because the pressure-driven ABL model cannot represent wind speeds exceeding the geostrophic wind. The absence of the geostrophic jet in the pressure-driven ABL model is related to the lack of Coriolis-induced wind veer (lack of coupling between the equations for $\mathrm{D}U/\mathrm{D}t$ and $\mathrm{D}V/\mathrm{D}t$), as explicitly shown by analytic solutions without wind veer for constant and linearly increasing eddy viscosity. The difference between the ABL models can become important for shallow boundary layers representing very stable atmospheric conditions if the pressure-driven ABL model is used as an inflow profile for three-dimensional RANS simulations of large wind turbines. Despite this challenge, one can employ the pressure-driven ABL model to isolate the effect of wind veer or ABL depth, when it is compared to an ABL model including wind veer or an ASL model, respectively. In addition, the Reynolds number similarity of the pressure-driven ABL model can be used to perform parametric studies of the effect of wind speed on wind farm flow simulations more quickly, similar to using an ASL inflow model (van der Laan et al., 2019). The proposed ABL model in combination with the turbulence model of Apsley and Castro (1997) cannot yet be used for complex terrain simulations, since the length-scale limiter does not behave appropriately over such terrain. However, it is possible to use the momentum source term of the ABL model in combination with a turbulence model suited for complex terrain; such work is still under development. Furthermore, softening of the effect of the length-scale limiter, to capture the influence of the strength of ABL-capping inversion (Kelly et al., 2019b), is also needed (and underway) to better capture the top-down effects entraining momentum into the wind farm (e.g., Kelly et al., 2019a).

Appendix A: Obtaining a desired ABL profile using Rossby number similarity

The results of the k–ε ABL model of Apsley and Castro (1997) (with wind veer) and the proposed ABL model without wind veer from Sect. 3 can both be used as inflow profiles for three-dimensional wind farm simulations for homogeneous terrain. Wind farm simulations are often run as flow cases, where a desired turbulence intensity I_ref and wind speed S_ref is set at a reference height z_ref. In this section, a methodology is presented on how to use the Rossby number similarity (Sect. 5) to quickly find the geostrophic wind speed and maximum turbulence length scale that corresponds to the desired inflow profile from a pre-calculated library of all possible normalized ABL profiles based on two Rossby numbers (Eqs. 4 or 14). In order to find a unidirectional ABL profile that matches an ABL profile with veer we perform the following steps.

1. Simulate non-dimensional libraries of ABL profiles for both models (with and without wind veer) based on the two Rossby numbers (i.e., both ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ and ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}$ as well as Ro₀ and Ro_ℓ). These Rossby numbers can be made up of any combination of G, f_c, z₀, and ℓ_max as shown in Fig. 2. For example, choose $G=G_{\mathrm{lib}}=\mathrm{10}\mathrm{m}/\mathrm{s}$, $f_{\mathrm{pg}}=f_{\mathrm{c}}=f_{\mathrm{c},\mathrm{lib}}={\mathrm{10}}^{-\mathrm{4}}{\mathrm{s}}^{-\mathrm{1}}$, ${\mathrm{\ell }}_{\max ,\mathrm{lib}}=G_{\mathrm{lib}}/\left(f_{\mathrm{c},\mathrm{lib}}{\mathit{Ro}}_{\mathrm{\ell }}\right)$, and $z_{\mathrm{0}}=G_{\mathrm{lib}}/\left(f_{\mathrm{c},\mathrm{lib}}{\mathit{Ro}}_{\mathrm{0}}\right)$. We have chosen to use a parametric study of two Rossby numbers in logarithmic space, e.g., $({\mathit{Ro}}_{\mathrm{0}},{\mathit{Ro}}_{\mathrm{\ell }})=({\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}},{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }})=({\mathrm{10}}^a,{\mathrm{10}}^b)$, with $a=[\mathrm{5},\mathrm{10}]$ using a spacing of 0.2 and $b=[\mathrm{2},\mathrm{4.5}]$ with a spacing of 0.1 for b<3.5 and a finer spacing of 0.05 for b>3.5. The results are stored as a function of a normalized height $z_{\mathrm{norm}}=(z+z_{\mathrm{0}})f_{\mathrm{c},\mathrm{lib}}/G_{\mathrm{lib}}$.

2. Set z_0,ref, f_c,ref, S_ref, I_ref, and z_ref.

3. Find Ro₀ and Ro_ℓ from the ABL library with veer that satisfy the reference values of Step 2 and calculate G and ℓ_max:

   • a.

     For each (Ro₀, Ro_ℓ) pair interpolate z_norm, where I=I_ref is obtained and calculate the corresponding geostrophic wind speed $G({\mathit{Ro}}_{\mathrm{0}},{\mathit{Ro}}_{\mathrm{\ell }})=(z_{\mathrm{ref}}+z_{\mathrm{0},\mathrm{ref}})f_{\mathrm{c},\mathrm{ref}}/z_{\mathrm{norm}}$ and wind speed S(Ro₀,Ro_ℓ).

   • b.

     Curve A in Fig. A1 is a set of points satisfying $\left(\mathit{\right\{}{\mathit{Ro}}_{\mathrm{\ell },A}\mathit{\}},\mathit{\{}{\mathit{Ro}}_{\mathrm{0},A}\mathit{\left\}}\right)=\mathit{\left\{}\right({\mathit{Ro}}_{\mathrm{\ell }},{\mathit{Ro}}_{\mathrm{0}}\left)\right|S({\mathit{Ro}}_{\mathrm{\ell }},{\mathit{Ro}}_{\mathrm{0}})=S_{\mathrm{ref}}\mathit{\}}$.

   • c.

     Curve B in Fig. A1 represents $\left(\mathit{\right\{}{\mathit{Ro}}_{l,B}\mathit{\}},\mathit{\{}{\mathit{Ro}}_{\mathrm{0},B}\mathit{\left\}}\right)=\left(\mathit{\right\{}{\mathit{Ro}}_{l,A}\mathit{\}},\mathit{\{}{G}_A\mathit{\}}/(f_{\mathrm{c},\mathrm{ref}}{z}_{\mathrm{0},\mathrm{ref}}\left)\right)$, where G_A corresponds to extracted G values from curve A.

   • d.

     Ro₀ and Ro_ℓ can be obtained by the intersection of curves A and B.

   • e.

     Calculate ${\mathrm{\ell }}_{max}=z_{\mathrm{0},\mathrm{ref}}{\mathit{Ro}}_{\mathrm{0}}/{\mathit{Ro}}_{\mathrm{\ell }}$ and $G=z_{\mathrm{0},\mathrm{ref}}{f}_{\mathrm{c},\mathrm{ref}}{\mathit{Ro}}_{\mathrm{0}}$.

4. Find ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$ and ${\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{\ell }}$ from the ABL library without wind veer similar to steps 3a–d, calculate $f_{\mathrm{pg}}=f_{\mathrm{c},\mathrm{ref}}{\mathit{Ro}}_{\mathrm{0}}/{\stackrel{\mathrm{̃}}{\mathit{Ro}}}_{\mathrm{0}}$, and then correct the geostrophic wind speed as $G_{\mathrm{pg}}=G{S}_{\mathrm{ref}}/S_{\mathrm{pg}}$ using Reynolds number similarity.

Here, the subscript pg is used for the parameters corresponding to the ABL model without wind veer, where S_pg is the obtained wind speed at z_ref before correcting the geostrophic wind speed G to G_pg. In addition, we use Ro₀ and Ro_ℓ in logarithmic space in steps 3a–d. If the ABL libraries only contain a few profiles for different Rossby numbers, then obtaining the desired Rossby numbers may lead to errors. In case this one could use the ABL libraries as an initial guess for a numerical optimization that results in the required G and ℓ_max and f_pg and G_pg. An example of obtaining a set of Rossby number from an ABL library with wind veer (Step 3c) is depicted in Fig. A1. The example represents the neutral case as listed in Table 1.

Figure A1 Example of obtaining a set of Rossby numbers from an ABL library with wind veer, Step 3d.