# Universal Functions with Ekman Spiral and Monin-Obukhov Surface Layers

Martin Ferrand, Romain Pennel, Eric Dupont

► To cite this version:

Martin Ferrand, Romain Pennel, Eric Dupont. Universal Functions with Ekman Spiral and Monin-Obukhov Surface Layers. Boundary-Layer Meteorology, 2025, 191 (9), pp.43. 10.1007/s10546-025-00938-5. hal-05320518

HAL Id: hal-05320518

https://hal.science/hal-05320518v1

Submitted on 18 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Universal Functions with Ekman Spiral and Monin–Obukhov Surface Layers

Martin Ferrand<fnref n="1" />,<fnref n="2" /> · Romain Pennel<fnref n="2" /> · Eric Dupont<fnref n="1" />,<fnref n="2" />

Received: 19 June 2025 / Accepted: 13 August 2025 / Published online: 9 September 2025
© The Author(s) 2025

## Abstract

Today, computational fluid dynamics (CFD) is widely used for atmospheric dispersion at building scale. This type of simulation requires boundary conditions for wind, turbulence, and temperature. A classical way of imposing representative flow at open boundaries is to derive (universal) functions corresponding to idealised situations (e.g. constant shear stress, constant heat flux, such as in the Monin–Obukhov theory). In this paper, we first propose an analysis of universal functions for the surface layer and compare them to a 5-year data base of measurement on the SIRTA observatory (France). The conclusion is that these functions are in good agreement with the experimental data for moderate ratios of ζ (distance to the ground divided by Obukhov length), but most of them do not fulfil asymptotic behaviours in the very stable or convective limit. We also compare the measurements with the extensions proposed by Gryning et al. (Bound Layer Meteorol 124(2):251–268, 2007), which take into account the height of the atmospheric boundary layer. Then, following the work of Nieuwstadt (Noct Bound Layer J Atmos Sci 41(14):2202–2216, 1984. https://doi.org/10.1175/1520-0469(1984)041$lt$2202:TTSOTS$gt$2.0.CO;2), we propose universal functions able to reproduce the Ekman spiral and consistent with moment-turbulence closures for stably stratified atmospheres.

**Keywords** Universal functions · Stable condition · Ekman spiral · Boundary layer profiles · SIRTA

## 1 Introduction

Idealized wind and temperature profiles are useful for boundary conditions of computational fluid dynamics (CFD) simulations of the atmospheric boundary layer (ABL) e.g. to estimate wind energy production of a farm or to study local-scale pollutant dispersion (see

<footnote>Martin Ferrand
martin.ferrand@edf.fr</footnote>

<footnote n="1">¹ EDF R&D, MFEE, Chatou, France</footnote>

<footnote n="2">² CEREA, EDF, École nationale des ponts et chaussées, Chatou, France</footnote>
e.g., Kelly et al. 2019; Kelly and Van Der Laan 2023). For years, the Monin–Obukhov similarity theory (MOST) has given the theoretical wind speed and potential temperature profiles within the atmospheric surface boundary layer (ASL) using the so-called “universal functions” (Monin and Obukhov 1954) that are the subject of a very vast well-established literature (see e.g. Dyer 1974, for a review). However, in spite of the various expressions given in the literature (see, e.g., Carl et al. 1973; Högström 1988; Cheng and Brutsaert 2005; Hartogensis and Bruin 2005, which will be used in the comparative work in this article), none is reaching consensus, as pointed out by Foken (2006). Furthermore, this theory gives the meteorological profiles only within the ASL, also called constant flux layer, which is about 10% of the whole ABL, as pointed out in neutral cases by Richards and Norris (2019). However, CFD is expected to simulate the flow over a depth of several hundreds of metres, even a few kilometres, which is far beyond the surface layer depth. An extension has been proposed in Gryning et al. (2007), which gives a correct reproduction of the velocity magnitude for the entire atmospheric boundary layer, but it does not take the Coriolis forces into account, and thus the Ekman spiral, first introduced in Ekman (1905), is not reproduced. However, modelling the wind veer is very important for many applications such as atmospheric dispersion of pollutants or wind farms production assessment. Nieuwstadt (1984) derived profiles for stable ABL starting from second-order moment turbulent closures and accounting for Coriolis force, and thus recovered the Ekman spiral. However, this theory is only valid above the surface layer for stable situations where the Richardson flux number reaches a constant critical value. Recent work of Narasimhan et al. (2024) combines MOST and Nieuwstadt (1984), taking up the hypothesis of a 3/2 power law for the total turbulent horizontal stress and compare large eddy simulations with the proposed velocity profiles with good match. The value of the exponent 3/2 in Nieuwstadt theory is in fact intrinsically linked to the value of the thermal turbulent flux decay.

The main two objectives of this paper are, first, to combine Monin–Obukhov theory with a generalization of the Nieuwstadt proposition in the case of stable stratification with an arbitrary power law exponent in the thermal turbulent flux and, second, to validate these profiles against long-term experimental data.

In Sect. 2, we recall some theoretical aspects of the ASL and the ABL, such as MOST and asymptotic behaviours including second order statistics, we present the representative universal functions proposed in the literature and considered in this article, then we propose an extension of (Nieuwstadt 1984) that combines MOST and Ekman layer. In Sect. 4, the proposed profiles are compared with experimental data of *Site Instrumental de Recherche par Télédétection Atmosphérique* (SIRTA – https://sirta.ipsl.fr/), and with various proposit ions of Monin–Obukhov universal functions and the proposition of Gryning et al. (2007).

## 2 Atmospheric Boundary Layer Theories

### 2.1 The Monin–Obukhov Similarity Theory Within the Surface Layer

#### 2.1.1 Wind and Temperature Vertical Gradients

The surface universal profile derivatives for momentum and energy are derived by assuming an atmospheric flow in a steady regime, neglecting the coriolis forces, the large-scale pres-
sure gradients and heat source terms. We assume here that the mean flow is aligned with the x-axis. It follows from momentum and energy balances that the shear stress $\tau_{xz,0} = -\rho u'w'$ and the heat flux $q_0 = -\rho C_p w'\theta'$ are constants (and so are the turbulent fluxes), where u and w are the horizontal and vertical components of the wind speed respectively, with $\epsilon_z$ orientated upward, $\theta$ is the potential temperature, $\langle \cdot \rangle$ is the reynolds averaging operator, and $\langle \cdot \rangle'$ is its fluctuation. The friction velocity is written using the turbulent momentum flux: $u_{*0} := \sqrt{-u'w'}$. In the same way, a friction temperature is derived from the turbulent heat flux: $\theta_{*0} := -\frac{w'\theta'}{u_{*0}}$. The turbulent kinetic energy (TKE) budget under the same assumptions gives:

$$P + B - \epsilon \approx 0, \quad (1)$$

where $P(= -\overline{u'w'} \frac{\partial \overline{u}}{\partial z}$ here), $B(= \frac{g}{\theta_0} \overline{w'\theta'}$ here), and $\epsilon$ are, respectively, the TKE production by the shear, the TKE production by buoyancy and the TKE dissipation rate, $g$ is the gravity acceleration, and $\theta_0$ is the average potential temperature within the surface layer. This is formally rewritten as $\frac{\partial \overline{u}}{\partial z} = F(u_{*0}, \frac{g}{\theta_0} \overline{w'\theta'}, z)$, where $z$ is the distance to the ground, as $\epsilon$ is itself function of $P$, $B$ and the distance to the ground.

Applying Vaschy-Buckingham theorem gives the well-known universal function for the momentum derivative:

$$\frac{\partial \overline{u}}{\partial z} = \frac{u_{*0}}{\kappa z} \varphi_m(\zeta), \quad (2)$$

where $\kappa$ is the Von Kàrmàn constant and $\zeta := z/L$, with $L$ being the Obukhov length defined by $L := -\frac{u_{*0}^3}{\kappa \frac{g}{\theta_0} \overline{w'\theta'}} = \frac{u_{*0}^2}{\kappa \frac{g}{\theta_0} \overline{w'\theta'}}, z$. Note that $\kappa$ is used in (2) by convention, so that

$\varphi_m(0) = 1$ (see Yaglom (1977)).

Similarly, using temperature variance one gets:

$$\frac{\partial \overline{\theta}}{\partial z} = \frac{\theta_{*0}}{\kappa z} \varphi_h(\zeta). \quad (3)$$

$\varphi_m$ and $\varphi_h$ are the so-called momentum and heat universal functions. Note that $\varphi_h(0)$ is formally the turbulent Prandtl number in neutral condition, denoted by $Pr_{T,0}$ (see e.g. Jacobson (2005); Foken and Mauder (2024)). As for $\kappa$, the value of $Pr_{T,0}$ is subject to debate and their range vary among the authors. Businger et al. (1971) proposed $\kappa \approx 0.35$ and $Pr_{T,0} \approx 0.74$, but had to be corrected because the shell cup anemometers had too much over-speeding (Wieringa (1980)). We retain in this article the value of $\kappa \approx 0.4$ which reads in the classical range of validity (see Högström (1996); Foken and Mauder (2024)). For a correct comparison of $\varphi_h$ functions, also the same value of $Pr_{T,0}$ must be considered, as explained by Yaglom (1977).

Classically, we define the following dimensionless numbers which are of importance in atmospheric boundary layer theories. First, the gradient Richardson number $Ri_g$ is written as: $Ri_g := \frac{g \frac{d\overline{\theta}}{dz}}{(\frac{d\overline{u}}{dz})^2}$, which can be rewritten using the universal functions:
$$Ri_g = \zeta \frac{\varphi_h(\zeta)}{\varphi_m^2(\zeta)}. \qquad (4)$$

Then, the flux Richardson number is defined by: $Ri_f := -\frac{B}{P}$. In the ASL it leads to:

$$Ri_f = \frac{\zeta}{\varphi_m(\zeta)}. \qquad (5)$$

## 2.1.2 Asymptotic Behaviours

The objective of this section is to recall some important behaviours presented in Monin and Yaglom (2007) that must verify universal functions. Although these limits are not directly evaluated in the present work, they can provide objective criteria for choosing one empirical function over another when implementing the generalized profiles proposed in 2.3 and 2.4.

### *Unstable Stratification:*

First, we look at the limiting case $\zeta \to -\infty$, corresponding for a given height $z$ to strongly unstable conditions ($L \to 0^-$), which may happen in the case of purely thermal turbulence with very low wind. The magnitudes, different from the wind speed, no longer explicitly depend on $u_{*0}$. The potential temperature depends on $z$, $q_0/(\rho_0 C_p)$ and $g/\theta_0$. Using a dimensional analysis, the potential temperature can be written as:

$$\bar{\theta}(z) = \theta_{\infty} + C_1 \left( \frac{q_0}{\rho_0 C_p} \right)^{2/3} \left( \frac{gz}{\theta_0} \right)^{-1/3},$$

where $\theta_{\infty}$ is a constant temperature and $C_1$ is a dimensionless constant. Differentiating the potential temperature profile against $z$, gives:

$$\frac{\partial \bar{\theta}}{\partial z} \underset{\zeta \to -\infty}{\sim} -\frac{C_1}{3} \left( \frac{q_0}{\rho_0 C_p} \right)^{2/3} \left( \frac{g}{\theta_0} \right)^{-1/3} z^{-4/3}.$$

Using Eq. (3), we obtain:

$$\varphi_h(\zeta) \underset{\zeta \to -\infty}{\sim} \frac{C_1 \kappa^{4/3}}{3} |\zeta|^{-1/3}.$$

Since the turbulent Prandtl number tends to a constant value $Pr_t^{-\infty} = \lim_{\zeta \to -\infty} Pr_t$ for very unstable stratification (see Monin and Yaglom (2007), Section 7.3, equation (7.38)), we get:

$$\varphi_m(\zeta) \underset{\zeta \to -\infty}{\sim} \frac{C_2}{3} |\zeta|^{-1/3},$$

where $C_2 = C_1 \kappa^{4/3} / Pr_t^{-\infty}$.

From this last expression and Eq. (5), we deduce the behaviour of the flux Richardson number:
$$Ri_f \underset{\zeta \to -\infty}{\sim} -\frac{3}{C_2} |\zeta|^{4/3}. \qquad (6)$$

This behaviour is coherent with O'KEYPS equation (see Panofsky (1963)):

$$\varphi_m^4 - \gamma \zeta \varphi_m^3 = 1, \qquad (7)$$

$\gamma$ a constant.

In the universal functions presented in Appendix 7, only Carl et al. (1973) verifies asymptotic behaviour Eq. (5), which was an objective of this proposition compared to Businger et al. (1971), where attention was paid to small negative values of $\zeta$ (also in agreement with Eq. (7) for $\zeta \to 0^-$).

*Stable Stratification:*

Now, we will look at the limit case $\zeta \to +\infty$. First, the flux Richardson number increases as the stability increases. In fact, when stability increases, the turbulence production term $P$ within the TKE balance equation decreases; hence the flux Richardson number increases. Furthermore, according to Eq. (1), the flux Richardson number is lower than 1. We deduce that the flux Richardson number tends toward a finite critical limit $Ri_f^c$ when $\zeta$ tends toward plus infinity. From Eq. (5), it follows that:

$$\varphi_m(\zeta) \underset{\zeta \to +\infty}{\sim} \frac{\zeta}{Ri_f^c}. \qquad (8)$$

The critical flux Richardson number value is about 0.2 (see, e. g., Nieuwstadt 1984).

This is a well-known issue, which Holtslag and Nieuwstadt (1986) referred to as the so-called z-less scaling. Some authors (e.g., Foken and Mauder (2024)) argue that the universal function can only be applied for $z/L$ values not exceeding approximately one. Nevertheless, with a view toward applying these profiles in CFD, it is important to ensure that $\varphi_m \le \frac{\zeta}{Ri_f^c}$.

In the universal functions presented in Appendix 7, only the one proposed by Hartogensis and Bruin (2005) has an asymptotic behaviour proportional to $\zeta$ under stable conditions as shown in Fig. 1.

### 2.1.3 Universal Functions of Turbulent Kinetic Energies

In this section, we recall the asymptotic behaviours of turbulent kinetic energy for all components presented in Monin and Yaglom (2007), so that they can be compared with the experimental data in Sect. 4.2. For direction $i \in \{1, 2, 3\}$ (with by convention $u_1 = u$, $u_2 = v$, $u_3 = w$), we denote by $\varphi_{u_i}$ the universal functions of velocity fluctuation in direction $i$, which varies with $\zeta$:

$$\varphi_{u_i}(\zeta) = \frac{\sigma_{u_i}}{u_{*0}} = \frac{\sqrt{u_i'^2}}{u_{*0}}, \qquad (9)$$
<figure><img src="image_3.png" /><figcaption>Fig. 1 Universal functions $\varphi_m$ considered in the present article (the acronyms used in the legend are given in Sect. 2.2)</figcaption></figure>

where $\sigma_\alpha$ is the standard deviation of the quantity $\alpha$. Note that in this section, we neglect the dependence on $z$, which is only valid in the surface layer.

Therefore, the turbulent kinetic energy writes $k/u_{*0}^2 = \frac{1}{2}(\varphi_u^2 + \varphi_v^2 + \varphi_w^2)$.

First, for fairly small $|\zeta|$ in nearly neutral conditions, the fluctuating velocity universal functions could be developed in power series:

$$ \varphi_{u_i}(\zeta) \sim_{\zeta \to 0} A_i + B_i\zeta, \qquad (10) $$

where $A_i$ and $B_i$ are dimensionless constants.

Second, when $\zeta \to -\infty$, that is, in the limiting case of purely thermal turbulence, according to Monin and Yaglom (2007) there is no horizontal mean velocity, therefore the horizontal components of the velocity fluctuations become statistically equal, i.e., $\varphi_{u_1} \approx \varphi_{u_2}$. Moreover they vary as $|\zeta|^{1/3}$, so:

$$ \varphi_{u_i} \sim_{\zeta \to -\infty} C_i |\zeta|^{1/3}, \qquad (11) $$

where $C_1 = C_2$ and $C_3$ are dimensionless constants.

Finally, when $\zeta \to +\infty$, all universal functions must tend toward constants, so:

$$ \lim_{\zeta \to +\infty} \varphi_{u_i} = D_i, \qquad (12) $$

with $D_i$ dimensionless constants. As presented in Sect. 4.2, experimental data are in good agreement with this theory.

## 2.2 Some Universal Functions Proposed in the Literature

In this paper, we limit the study to a selection of universal functions, presented in more detail in the Appendix 7. First, Businger et al. (1971) universal functions (denoted B71 hereafter) are quite popular. But, we chose Högström (1988) universal functions (denoted
HO88), who proposed a corrected version of the universal functions of the same type as B71, with same value of $κ ≃ 0.4$ taken here. Carl et al. (1973) universal functions (denoted C73) under unstable conditions are also considered because $φ_m$ satisfies the theoretical behaviour of the limiting case $ζ → −∞$. Finally, Hartogensis and Bruin (2005) (denoted HB05) and Cheng and Brutsaert (2005) (denoted CB05) are finally considered under stable conditions.

## 2.3 Gryning Proposition for the Atmospheric Boundary Layer Profiles

Strictly speaking, MOST is only valid within the ASL, which represents about 10% of the ABL, and not above. At maximum, the ASL is a few hundred metres high. However, in CFD simulations, it is often important to model the atmosphere from several hundreds metres to a few kilometres high. Gryning et al. (2007), proposed an extension of the wind speed profile up to the ABL height, and was complemented in Peña et al. (2010). However, this proposition makes no distinction between the horizontal components of the velocity, and reads:

$$ \frac{\partial \bar{u}}{\partial z} = \frac{u_*}{\kappa l}, \qquad (13) $$

where the turbulence length scale $l$ is expressed as Gryning et al. (2007):

$$ \frac{1}{l} = \frac{\varphi_m(z/L)}{L_{SL}} + \frac{1}{L_{MBL}} + \frac{1}{z_i - z}, \qquad (14) $$

where $L_{SL} = z$ is the length scale of the surface layer, $L_{MBL}$ is the length scale of the middle boundary layer, and $L_{UBL} = z_i - z$ is the length scale of the upper boundary layer, $z_i$ being the depth of the boundary layer.

For the middle boundary layer, Blackadar (1962) showed that $L_{MBL}$ is independent of height. Furthermore, Gryning et al. (2007) give a parametrization for $L_{MBL}$:

$$ \frac{u_{*0}}{f L_{MBL}} = \left( -2 \ln \left( \frac{u_{*0}}{f z_0} \right) + 55 \right) \exp \left( - \frac{(u_{*0}/(f L))^2}{400} \right), \qquad (15) $$

where $f$ is the Coriolis parameter. Above the surface layer, considering a large-scale pressure gradient leads to momentum fluxes dependent on $z$, which leads to a friction velocity profile, taken as a power law of $(1 - z/z_i)$. Furthermore, as the turbulence vanishes at the top of the atmospheric boundary layer, the friction velocity writes:

$$ u_*(z) = u_{*0} \left( 1 - \frac{z}{z_i} \right)^p, \qquad (16) $$

with $p = 1$ in Gryning et al. (2007).

### 2.3.1 Neutral Conditions

Integrating Eq. (12) from $z_0$ to $z$ leads to:
<figure><img src="image_4.png" /><figcaption>Fig. 2 Turbulent length scale profiles (left) and wind profiles (right) in neutral conditions for $z_i = 1000$ m and $L_{MBL,N} = 150$ m, the influence of each term of $l^{-1}$ is displayed</figcaption></figure>

$$ \frac{\bar{u}(z)}{u_{*0}} = \frac{1}{\kappa} \left( \ln \left( \frac{z}{z_0} \right) + \frac{z - z_0}{L_{MBL,N}} \left( 1 - \frac{z + z_0}{2z_i} \right) \right). \quad (17) $$

Figure 2 reproduces the profile of the turbulence length scale and the profile of the wind speed module and shows the influence of each term on the expression of the length scale. We can note that the theory of Gryning et al. (2007) provides a higher wind speed at the top of the ABL compared to assuming that the Monin–Obukhov theory is valid for the entire atmospheric boundary layer, that is corroborated by profiles presented in Sect. 4.

### 2.3.2 Stable Conditions

For stable conditions, Gryning et al. (2007) took<fnref n="1" /> $\varphi_m(\zeta) = 1 + b\zeta$. In the following, we take $b = 4.8$.

Integrating Eq. (12) from the roughness length to the height $z$ gives:

$$ \frac{\bar{u}(z)}{u_{*0}} = \frac{1}{\kappa} \left( \ln \left( \frac{z}{z_0} \right) + (z - z_0) \left( \frac{b}{L} + \frac{1}{L_{MBL}} \right) \left( 1 - \frac{z + z_0}{2z_i} \right) \right). \quad (18) $$

### 2.3.3 Unstable Conditions

For unstable conditions, Gryning et al. (2007) took the proposition of Carl et al. (1973)<fnref n="2" />

$$ \varphi_m(\zeta) = (1 - 12\zeta)^{-1/3}. $$

Integrating Eq. (12) from the roughness length to the height $z$ gives:

$$ \frac{\bar{u}(z)}{u_{*0}} = \frac{1}{\kappa} \left( \ln \left( \frac{z}{z_0} \right) - \psi \left( \frac{z}{L} \right) + \frac{z}{L_{MBL}} - \frac{z}{z_i} \left( \frac{z}{2L_{MBL}} \right) \right), \quad (19) $$

<footnote n="1">¹As Businger et al. (1971) with $b = 4.7$ or Högström (1988) with $b = 4.8$ for $\zeta < 0.5$.</footnote>

<footnote n="2">²Note that this choice is compatible with the asymptotic behaviour presented in Sect. 2.1.2</footnote>
where:

$$ \psi\left(\frac{z}{L}\right) = \frac{3}{2} \ln\left(\frac{1+X+X^2}{3}\right) - \sqrt{3}\arctan\left(\frac{1+2X}{\sqrt{3}}\right) + \frac{\pi}{\sqrt{3}}, \quad (20) $$

and $X = (1 - 12z/L)^{1/3}$.

In the following, the profiles (15), (16) and (17) presented in this section are denoted by G07. Note that, under stable conditions, G07 is equivalent to HO88 linear in the surface layer.

## 2.4 Generalized Nieuwstadt Theory

In this section, we propose an extension of the ABL profiles derived in Nieuwstadt (1984) for stable conditions, which account for the flux variation with height. This section is independent of and alternative to Gryning et al. (2007). The atmospheric flow is still assumed to have reached a steady regime, but the Coriolis force is taken into account, and a large-scale pressure gradient is considered (or equivalently, a geostrophic wind is imposed). Under stable conditions, above the surface layer, the gradient and flux Richardson numbers tend to constant values (see Monin and Yaglom (2007), Section 7.3). From equilibrium assumptions on second-order moment closures, Nieuwstadt (1984) derived profiles for the velocity components and temperature.

These profiles are generalized in two ways. First, following Nieuwstadt (1984), a radiative source term $S$ is added to the Reynolds averaged energy equation, and is proportional to a power law of $1 - \frac{z}{z_i}$ (see also Yokoyama et al. (1979)), but the exponent is kept unspecified. Under these assumptions, Reynolds averaged Navier-Stokes equations and energy conservation write:

$$ 0 = f(\bar{v} - v_g) - \frac{\partial \overline{u'w'}}{\partial z}, \quad (21) $$

$$ 0 = -f(\bar{u} - u_g) - \frac{\partial \overline{v'w'}}{\partial z}, \quad (22) $$

$$ 0 = -\frac{\partial \overline{w'\theta'}}{\partial z} + S, \quad (23) $$

where $f$ is the Coriolis parameter and $u_g$ and $v_g$ are the geostrophic wind speed components. The gradient and flux Richardson numbers rewrite:

$$ Ri_g = \frac{\frac{g}{\theta_0} \frac{\partial \bar{\theta}}{\partial z}}{\left(\frac{\partial \bar{u}}{\partial z}\right)^2 + \left(\frac{\partial \bar{v}}{\partial z}\right)^2}, \quad (24) $$
$$Ri_f = - \frac{\frac{g}{\theta_0} \overline{w'\theta'}}{-\overline{u'w'} \frac{\partial \overline{u}}{\partial z} - \overline{v'w'} \frac{\partial \overline{v}}{\partial z}}. \quad (25)$$

Above the surface layer, under stable conditions, they tend towards constant critical values, respectively $Ri_g^c$ and $Ri_f^c$.

As the sink term $S$ follows a power law of the height, the turbulent heat flux writes:

$$\overline{w'\theta'} = \overline{w'\theta'}_0 \left(1 - \frac{z}{z_i}\right)^\alpha, \quad (26)$$

where $\overline{w'\theta'}_0$ is the turbulent heat flux value at the surface, $z_i$ is the atmospheric boundary-layer height and $\alpha$ is the dimensionless exponent. In Nieuwstadt (1984), $\alpha$ is taken to be 1, which means that the source term is constant with height. In Deardorff (1972), the turbulent heat flux follows a quadratic law with height, that is, $\alpha = 2$, which means that the source term is linear. In order to be as general as possible and to solve the problem of singular profiles at the junction with the free atmosphere, we let $\alpha$ without specifying a value yet.

Then, Eqs. (18) and (19) have a solution if the height of the atmospheric boundary layer satisfies the following criterion (see Appendix 6):

$$z_i = \left( \frac{(\alpha + 1)\sqrt{\alpha(\alpha + 2)}}{2} \kappa R_i^c \right)^{1/2} \sqrt{\frac{u_{*0}L}{f}}. \quad (27)$$

Equation (24) gives a parameterization of $z_i$ as a function of the flow characteristics. We note that this expression of $z_i$ has the same form as the one given by Zilitinkevich (1972).

For $\alpha = 1$, the constant is $\sqrt{3\kappa R_i^c} \approx 0.37$ while it is 0.4 in Zilitinkevich (1972) parameterisation. Furthermore, the friction velocity writes:

$$u_* = \sqrt[4]{\overline{u'w'}^2 + \overline{v'w'}^2} = u_{*0} \left(1 - \frac{z}{z_i}\right)^{(\alpha+2)/4}. \quad (28)$$

Finally, the equations for the vertical profiles of each mean horizontal wind component (aligning the x component with wind direction near the ground) are:

$$\frac{\partial \overline{u}}{\partial z} = \frac{u_{*0}}{\kappa L R_i^c} \left(1 - \frac{z}{z_i}\right)^{\alpha/2-1} \cos \left( \frac{\sqrt{\alpha(\alpha+2)}}{2} \ln \left(1 - \frac{z}{z_i}\right) \right), \quad (29)$$

$$\frac{\partial \overline{v}}{\partial z} = \frac{u_{*0}}{\kappa L R_i^c} \left(1 - \frac{z}{z_i}\right)^{\alpha/2-1} \sin \left( \frac{\sqrt{\alpha(\alpha+2)}}{2} \ln \left(1 - \frac{z}{z_i}\right) \right), \quad (30)$$

and the mean potential temperature vertical profile equation is:

$$\frac{\partial \overline{\theta}}{\partial z} = \frac{\theta_{*0} R_i^c}{\kappa L (R_i^c)^2} \left(1 - \frac{z}{z_i}\right)^{\alpha-2}. \quad (31)$$
One can notice in Eq. (28), that the temperature gradient profile is continuous between the free atmosphere and ABL if and only if $\alpha > 2$. In Nieuwstadt (1984), the Richardson numbers are assumed to be constant. The second ingredient to extend this theory is to blend surface layer profiles with (26) (27) using (4) and (5):

$$ \frac{\partial \bar{u}}{\partial z} = \frac{u_{*0}}{\kappa z} \varphi_m \left(\frac{z}{L}\right) \left(1 - \frac{z}{z_i}\right)^{\alpha/2-1} \cos \left( \frac{\sqrt{\alpha(\alpha+2)}}{2} \ln \left(1 - \frac{z}{z_i}\right) \right), \quad (32) $$

$$ \frac{\partial \bar{v}}{\partial z} = \frac{u_{*0}}{\kappa z} \varphi_m \left(\frac{z}{L}\right) \left(1 - \frac{z}{z_i}\right)^{\alpha/2-1} \sin \left( \frac{\sqrt{\alpha(\alpha+2)}}{2} \ln \left(1 - \frac{z}{z_i}\right) \right), \quad (33) $$

and the proposed potential temperature vertical profile used is:

$$ \frac{\partial \bar{\theta}}{\partial z} = \frac{\theta_{*0}}{\kappa z} \varphi_h \left(\frac{z}{L}\right) \left(1 - \frac{z}{z_i}\right)^{\alpha-2}. \quad (34) $$

We can note that, in the lower part of the atmospheric boundary-layer, for $z \ll z_i$, $1 - z/z_i \approx 1$ such that:

$$ \begin{aligned} \frac{\partial \bar{u}}{\partial z} &\approx \frac{u_{*0}}{\kappa z} \varphi_m \left(\frac{z}{L}\right), \\ \frac{\partial \bar{v}}{\partial z} &\approx 0, \\ \frac{\partial \bar{\theta}}{\partial z} &\approx \frac{\theta_{*0}}{\kappa z} \varphi_h \left(\frac{z}{L}\right), \end{aligned} \quad (35) $$

and the vertical profiles of MOST are recovered; whereas, above the surface layer, the generalized Nieuwstadt profiles are recovered.

# 3 Data Description and Analysis

## 3.1 Data Description

For this study, 5 years of measurements taken at the French national atmospheric observatory SIRTA were used, from 2012 to 2017. The SIRTA observatory is located on the campus of École Polytechnique, about 20 km south of Paris, France (Haeffelin et al. 2005). These data were 10-min runs obtained with two ultrasonic anemometers (Metek USA-1) placed on the top of a 10 m mast and on the top of a 30 m mast, respectively. These anemometers measured the 3 wind components and "sonic" temperature (which corresponds approximately to virtual temperature), at 10 Hz frequency. Measurements of the 3 wind components at higher heights were obtained with two pulsed Doppler wind Lidars (light detection and ranging). The first one is a Vaisala Windcube WLS7 v2 measuring from 40 m to a maximum of 290 m, whereas the second one is a Vaisala Windcube WLS70 measuring from 100 m to a maximum of 2000 m, both at about 1 Hz frequency. All these instruments were installed on
the main measurement area of SIRTA. As for other similar studies (e.g., Peña et al. (2010)), 10-min averages are used, which is more than ten times the turbulence integral time-scale for $z \le 10m$.

To account for missing data in the temporal measurements, only 10-min periods with at least 80% completeness were retained across all measuring devices. Finally, there were several obstacles around the instrumented area that can disrupt the data, such as forests, so only the measurements for the wind coming from sector 265°–295° at 30 m height were kept, because in this sector the terrain is flat and there is no upstream obstacle. After these filterings, around 12000 10-min samples were left.

## 3.2 Data Analysis

First, all the parameters governing the ABL physics had to be computed. The friction velocity at the surface was computed by:

$$u_{*0} = \left( \overline{u'w'}^2 + \overline{v'w'}^2 \right)^{1/4} . \qquad (36)$$

Figure 3 shows that the trend between the turbulent heat flux measurements at 10 m and at 30 m is much closer (the slope of the linear regression is 0.98). For the following, since the profiles are normalized by ground friction velocity and friction temperature, we take the ones obtained using the measurement at 10 m.

Then, the Obukhov length was computed by:

$$\frac{1}{L} = - \frac{\kappa \frac{g}{\theta_0} \overline{w'\theta'_0}}{u_{*0}^3}, \qquad (37)$$

taking measurements at 10 m for $\theta_0$ and $\overline{w'\theta'_0}$ for the same reason as above.

In order to avoid low winds, which are not well suited for wind and temperature profiles reconstruction testing, we chose two thresholds on the friction velocity, $u_{*s}$ and $u_{*u}$, and we

<figure><img src="image_5.png" /><figcaption>(a) Friction velocity. The linear regression equation is $y = 1.14x + 0.03$. The coefficient of determination R2 is 0.8687</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>(b) Turbulent heat flux. The linear regression equation is $y = 0.98x + 0.00$. The coefficient of determination R2 is 0.9030</figcaption></figure>

<caption>**Fig. 3** Friction velocity and turbulent heat flux derived from ultrasonic anemometers measurements at 30 m against 10 m. +: measurements, red line: $y = x$, blue line: linear regression</caption>
kept only data for which $u_{*0} > u_{*s}$ for stable conditions and $u_{*0} > u_{*u}$ for unstable conditions. As the wind speed normalized by the friction velocity is often lower under unstable conditions than under stable conditions, the friction velocity threshold must be more restrictive under unstable conditions than under stable conditions. That is why the thresholds such as $u_{*u} > u_{*s}$ were chosen. For the following, we chose $u_{*s} = 0.2\text{m/s}$ and $u_{*u} = 0.3\text{m/s}$. After filtering on the friction velocity, about 6500 samples were left.

Figure 4 shows the occurrences of the wind speed measured at 10 m before and after filtering on the friction velocity. We see that the filtering mainly removes wind speed less than 3 m/s.

Then, the wind speed has been normalized by $u_{*0}$.

In order to compare theoretical values to measurements for similar stability conditions, the data was classified into 6 categories. First, a threshold on $L$ was chosen to classify the measurements under neutral conditions: $|L| > L_N = 1000\text{m}$. Then, with measurements left, inspired by Pasquill's classification proposed in Pasquill (1974), three categories under unstable conditions and two categories under stable conditions were established. For each of these conditions, the borders between the categories were chosen so that there was the same number of measurements within each category.

Then an estimation of the roughness length $z_0$ has been derived from the wind measurements. Using the dimensionless wind speed, in neutral conditions ($1/L \approx 0\text{m}^{-1}$) the roughness length is given by:

$$z_0 = z \exp \left( -\kappa \frac{\| \bar{u}(z) \|}{u_{*0}} \right). \qquad (38)$$

To compute the roughness length with this method, measurements at 10 m height within the neutral category ($|L| > L_N = 1000\text{m}$) were used and the median value of the roughness lengths was kept. Figure 5 shows the probability density function (PDF) of the computed values of the roughness length. We see that the probability is maximum for roughness length of about 5 cm and the probability becomes very weak after 20 cm, which is perfectly consistent with the dominant land use (grass) of the considered area. For the following, we will use the median value which is $z_0 = 6.6$ cm.

<figure><img src="image_7.png" /><figcaption>Fig. 4 Occurrences of wind speed at 10 m before and after filtering on $u_{*0}$</figcaption></figure>
Finally, we made the height dimensionless with the ABL height to avoid dependence on $z_i$. Once the ABL parameters calculated, we computed the mean, median, standard deviation and quantile 5% and 95% of the horizontal wind components and temperature, within each stability category for the same dimensionless height $z/z_i$.

### 3.3 Estimation of the ABL Height

The ABL height is a very important parameter as its estimation is required for the wind and temperature vertical profiles reconstruction. However, $z_i$ is often not directly accessible from the measurements, therefore a method is needed to derive it from other parameters.

The first way to estimate $z_i$ is the use of parameterizations given in the literature:

$$z_i = c_N \frac{u_{*0}}{f}, \qquad (39)$$

under neutral stratification with $c_N \approx 0.2$ (Duynkerke 1988) and:

$$z_i = c_S \sqrt{\frac{u_{*0}L}{f}}, \qquad (40)$$

under stable stratification with $c_S \approx 0.4$ (Zilitinkevich 1972).

The second way is to use the parameterization (24) derived from Nieuwstadt (1984) theory.

Finally, a last method consists in using the vertical profile of the wind components. Denote $z_{max}$ the height at which the maximum value of the wind component $u$ is reached. From Eq. (26), it follows that:

$$\cos \left( \frac{\sqrt{\alpha(\alpha + 2)}}{2} \ln \left( 1 - \frac{z_{max}}{z_i} \right) \right) = 0. \qquad (41)$$

<figure><img src="image_8.png" /><figcaption>Fig. 5 Probability density function of the local roughness height $z_0$</figcaption></figure>
As $\ln(1 - z/z_i) < 0$, the value within the cos will be negative and so, the last equation will be fulfilled for the first time when the magnitude within the cos will be equal to $-\pi/2$:

$$ \frac{\sqrt{\alpha(\alpha + 2)}}{2} \ln \left( 1 - \frac{z_{\max}}{z_i} \right) = -\frac{\pi}{2}, \quad (42) $$

and then:

$$ z_i = \frac{z_{\max}}{1 - \exp\left(-\frac{\pi}{\sqrt{\alpha(\alpha + 2)}}\right)}. \quad (43) $$

In practice, Stull (1988) criterion is used on wind speed module to keep only data which showed a low-level jet:

$$ \|\bar{u}\|_{\max} - \|\bar{u}_g\| \ge 2 \text{ m/s}. \quad (44) $$

The relations (34), (24), and (35) are valid only under stable conditions. For unstable conditions, from statistical analysis of $z_i$ measured at SIRTA (see Haeffelin et al. (2005)), we take 2000 m, 1600 m and 1200 m, respectively from the most to the least unstable category.

## 3.4 Estimation of the Turbulent Heat Flux Profile Parameter

Finally, the last parameter introduced with the generalized Nieuwstadt theory is the power of the turbulent heat flux profile $\alpha$. This parameter does have a high impact on the wind speed distribution profiles. It is then important to know a good estimation of it for each data sample. This parameter depends on the meteorological characteristics, particularly on the source/sink term of the radiative heat flux. It can be estimated by equating the expression of $z_i$ Eq. (24) and the one obtained with estimation (35), this yields:

$$ \left( \frac{(\alpha + 1)\sqrt{\alpha(\alpha + 2)}}{2} \kappa R i_f^c \right)^{1/2} \left( 1 - \exp\left(-\frac{\pi}{\sqrt{\alpha(\alpha + 2)}}\right) \right) = z_{\max} \sqrt{\frac{f}{u_{*0}L}}. \quad (45) $$

# 4 Results

## 4.1 Monin-Obukhov Similarity Theory Results

In this section, a comparison between the MOST mean wind speed profiles and the measurements is presented for some universal functions proposed in the literature. For that, the height of the atmospheric boundary layer $z_i$ is estimated using parameterisation (33) under neutral conditions and parameterisation (34) under stable conditions.
Figure 6 shows the theoretical (solid lines) and measured (dots) wind speed profiles. For each graph, the values of $L$ and $z_i$ represent the median values derived from the measurements within the respective stability category and are the values used to compute the theoretical profiles. The number of samples within each category is indicated. As mentioned previously, we kept only timestamps for which there were measurements with both anemometers and low-range LIDAR at all heights; therefore, the number of samples corresponds to the number of simultaneous validated measurements with anemometers and low-range LIDAR. Black bars represent the mean values plus or minus the standard deviation, i.e. about 68% of the measurements are within the bars. Red bars represent the range between quantiles 5% and 95%, i.e. 90% of the measurements are within the bars.

This figure shows that for some stability classes the theoretical profiles follow the averaged measurements relatively well up to heights larger than the ASL height (estimated as around 10% of the ABL) although MOST is valid only within this layer. However, large differences between MOST profiles and averaged measurements are often observed, especially under weakly unstable (c) and neutral (d) stratification, for which a change of universal function has no or a weak effect.

One can also notice in Fig. 6 and the following ones, that the standard deviation is quite large compared to the average wind value, which is similar to what is obtained in other work (e.g. see Peña et al. (2010), Figure 6). This is an illustration that steadiness is almost never reached for the whole ABL, but ensemble average of observations do verify good agreement compared with steady analytical profiles.

## 4.2 Wind Fluctuations Profiles

In this section, we verify the asymptotic behaviours discussed in Sect. 2.1.3. We used wind speed fluctuations measured with both anemometers and the high-range LIDAR from 100 m. The median of the standard deviation of each wind component, normalized by the friction velocity, has been computed as a function of $ζ$, for each stability condition (unstable, neutral, and stable). Furthermore, to validate the asymptotic behaviour, we computed the coefficients $A_i, B_i, C_i$ and $D_i$ defined in Eqs. (9), (10) and (11) using the least squares method.

Figures 7, 8 and 9 show that all the fluctuations follow the 1/3 power law remarkably accurately for the unstable stratification. Under stable stratifications, all dimensionless fluctuations tend to a constant value. Under neutral stratifications, both horizontal components fluctuations verify the negative slope we expected, and tend to the constant value observed in stable conditions, whereas the vertical component standard deviation is already approximately constant. Note that the value of $\phi_w(0) \simeq 1.24$ is with the 1.25 value proposed by Panofsky and McCormick (1960), while values from 0.8 to 1.3 can be found in the literature (Wilson (2008)). For $\phi_u(0) \simeq 2.51$, this is slightly higher than the values commonly cited from Panofsky et al. (1977), which range from $\phi_u(0) \simeq 2$ to 2.3 according to their equations (5) and (6).

### 4.2.1 Gryning Proposition for the Atmospheric Boundary Layer

For the computation of the vertical profile of wind speed according to Gryning et al. (2007) theory, the height of the atmospheric boundary layer has been computed using the parameterisation Zilitinkevich (1972).
<figure><img src="image_9.png" /><figcaption>**Fig. 6** Plot of the wind speed module using Monin and Obukhov (1954) theory and Duynkerke (1988) and Zilitinkevich (1972) parameterizations (33) and (34) for $z_i$. The acronyms used in the legend are explicited in Sect. 2.2</figcaption></figure>

Figure 10 shows that Gryning et al. (2007) theory provides theoretical profiles closer to the measurements than MOST, especially for neutral and weakly unstable conditions. This is due to the two additional length scales, as illustrated on Fig. 2. In Fig. 10, other $\phi_m$ functions than the one used by Gryning et al. (2007) in (13) are also plotted, revealing only small differences.
<figure><img src="image_10.png" /><figcaption>Fig. 7 Wind speed normalized fluctuations $\varphi_u$ against $\zeta$</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 8 Wind speed normalized fluctuations $\varphi_v$ against $\zeta$</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 9 Wind speed fluctuations $\varphi_w$ against $\zeta$</figcaption></figure>

## 4.3 Generalised Nieuwstadt Theory Results

This section presents a comparison between the theoretical wind components profiles from a generalisation of the Nieuwstadt theory and the measurements. The blending with the MOST profiles has been performed with the different universal functions listed in Sect. 2.2. First, as the generalised Nieuwstadt theory models a low-level jet (because of the trigonometric functions in the wind speed components vertical profiles equations), only data that show a low-level jet are kept. For that, the Stull (1988) criterion (36) is used. Furthermore, this theory is valid only under stable stratification. Finally, around 260 10-min samples are left.
<figure><img src="image_13.png" /><figcaption>(a) $L = -43\text{m}; L_{\text{MBL}} \to +\infty$; $z_i = 2000\text{m}$; (b) $L = -109\text{m}; L_{\text{MBL}} = 2688\text{m}$; $z_i = 911$ values</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>(c) $L = -327\text{m}; L_{\text{MBL}} = 185\text{m}$; $z_i = 1200\text{m}$; (d) $L = 8841\text{m}; L_{\text{MBL}} = 134\text{m}$; $z_i = 867\text{m}$; 910 values</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(e) $L = 331\text{m}; L_{\text{MBL}} = 132\text{m}$; $z_i = 420\text{m}$; (f) $L = 91\text{m}; L_{\text{MBL}} = 474\text{m}$; $z_i = 185\text{m}$; 1586 values</figcaption></figure>

<caption>**Fig. 10** Comparison with the SIRTA measurements of the wind speed vertical profile computed using Gryning et al. (2007) theory and Duynderke (1988) and Zilitinkevich (1972) parameterizations (33) and (34) for $z_i$</caption>

Profiles (29) and (30) contain two major parameters: $z_i$ and the power $\alpha$.

For the first case, the ABL height $z_i$ is estimated using the parameterisation of Zilitinkevich (1972) Eq. (34). The turbulent heat-flux profile power $\alpha$ has been set to 2 (as given by Deardorff (1972)). The first comparison (Fig. 11) considers the vertical profiles of the wind components under weakly and strongly stable stratification. In the weakly stable case (panel a), the profiles are generally well reproduced, although the gradient of $\bar{v}$ is underestimated
near the surface. According to (30), this suggests that the estimates of $\alpha$ and/or $z_i$ may need refinement. Under strongly stable stratification (panel b), larger discrepancies appear, particularly close to the ABL top.

To address this, the second configuration (Fig. 12) uses the theoretical parameterisation of $z_i$ given by (24), with $\alpha = 2$ as in the previous case. This approach improves agreement with measurements in the strongly stable regime but slightly degrades performance for weak stability.

In the third configuration (Fig. 13), $z_i$ is estimated from the height at which the u-component reaches its maximum value (35), again with $\alpha = 2$. Because this method uses the same criterion for both measured and theoretical profiles, the maximum of $u$ occurs at the same height in both. The resulting $u$ and $v$ profiles agree well with measurements in both stability regimes-except when using the HO88 linear function. However, as with the previous method, specifying $\alpha$ remains necessary and introduces uncertainty.

We now study the use the estimation of $\alpha$ provided by Eq. (37) (using the least squares method) and the corresponding $z_i$ from (35).

Figure 14 shows the histogram of the computed values of $\alpha$ for each stability category. We see that under weakly stable stratification (a), most of the values are between 0 and 4 (97 out of 112) and the number of occurrences decreases quickly when $\alpha$ increases. In contrast, under strongly stable stratification (b), the number of values occurrences decreases much slowly when $\alpha$ increases and only half of the values are between 0 and 15 (77 out of 150)

<figure><img src="image_16.png" /><figcaption>**Fig. 11** Comparison with the SIRTA measurements of the horizontal wind components vertical profile (left: u; right: v) computed using blended vertical wind speed distributions (29) and (30) and Zilitinkevich (1972) parameterization (34) for $z_i$</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Fig. 12** Plot of the wind speed (left: *u*; right: *v*) using blended vertical wind speed distributions (29) and (30) and *z*ᵢ parameterization derived from Nieuwstadt (1984) theory (24)</figcaption></figure>

and the distribution is very spread out. This is coherent with the literature<fnref n="3" /> For subsequent profile comparisons, we retain only cases with 0 ≤ α ≤ 4.

With this refined α estimation, Fig. 15 shows that the height of the wind-component maxima aligns well with observations. This agreement arises from matching the theoretical *z*ᵢ parameterisation with the *z*ᵢ inferred from the maximum of *u*. The amplitudes of these maxima are also close to measured values–particularly in the strongly stable category–making this approach more effective than the previous three.

# 5 Conclusions and Perspectives

In this article, we have used a five-year data set obtained at the SIRTA site to evaluate different wind theoretical profiles between the surface and the ABL top. First, a comparison between the wind speed vertical profiles calculated using the framework of the Monin–Obukhov theory with the ensemble averaged profiles measured at SIRTA has shown a relatively good agreement in the surface layer, and even up to higher heights for some stability classes. No clear differences between the studied universal functions except for more stable

<footnote n="3">³Holtslag and Nieuwstadt (1986) claimed that there is no theoretical justification for profile such as $w'/\theta'/w'/\theta'_0 = (1-z/z_i)^\alpha$ and, therefore, they should only be interpreted as convenient approximations. Yokoyama et al. (1979) also reported values from 1 to 3, though they acknowledged the number of usable data were very limited.</footnote>
<figure><img src="image_18.png" /><figcaption>Fig. 13 Plot of the wind speed (left: $u$; right: $v$) using blended vertical wind speed distributions (29) and (30) and $z_i$ estimation from the maximum value of $u$-component profile (35)</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. 14 Histograms of the parameter $\alpha$ used in the turbulent heat flux vertical profile</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>**Fig. 15** Plot of the wind speed (left: *u*; right: *v*) using blended vertical wind speed distributions (29), *α* estimation from the maximum value of *u*-component profile (37) and the corresponding atmospheric boundary layer height *z*ᵢ.</figcaption></figure>

or unstable cases, where asymptotic behaviours may give some argument to select Carl et al. (1973) for the unstable case, and Hartogensis and Bruin (2005) for the stable case. The experimental data also gives some insights on the variability around the ensemble averaged profiles, which is bigger more the flow is unstable.

Second-order statistics also agree globally well with the theory of Monin–Obukhov. However, large discrepancies are also noted when this theory is applied up to the ABL top, especially in neutral and weakly unstable situations. Then two proposals of the literature to extend the vertical profiles in the whole ABL have been studied. The proposition of Gryning et al. (2007), simple to implement, improves significantly the wind speed profiles, especially for the neutral conditions. However it does not allow to reproduce the Ekman spiral. The proposition of Nieuwstadt (1984) is appealing because it is derived from second-order moment turbulence closures and takes into account Coriolis forces. However it can be applied only above the surface layer, and only in stable conditions. This article presents an extension of this formulation in which the exponent governing the power law followed by the turbulent heat flux profile can be adjusted. We also propose a continuous blending of the surface layer and the rest of the ABL in order to get a formulation valid from the surface to the ABL top. With this theory it is possible to derive in stable conditions the vertical profiles of the two wind horizontal components separately. The comparison with the measurements has been presented for different methods of ABL height *z*ᵢ estimations and different values
of the heat flux profile power coefficient $\alpha$. Even though they require an estimation of $z_i$ (as in Gryning et al. (2007); Peña et al. (2010)) and the value of $\alpha$ is subject of debate (with either a value around 2 or a parametrization proposed in this study), the proposed normalized wind-component profiles have been shown to agree reasonably well with the ensemble-averaged normalized wind components from the SIRTA measurements.

As derived from second-order turbulence closures, these profiles offer a good candidate for inlet boundary conditions of CFD models with this type of turbulence models. In this context, it would be interesting to complete this study with a numerical validation of the derived profiles using them as inlet conditions for CFD simulations, and to test their conservation in time and space. In addition, it is also important to use these parameterizations to simulate atmospheric flows at local scale (less than 30 km), in order to assess their ability to reproduce the transport and diffusion processes of atmospheric pollutants. Finally, further work is needed to extend generalized Nieuwstadt profiles in unstable conditions.

## Derivation of the Generalized Nieuwstadt Theory

We recall the derivation of the profiles following Nieuwstadt (1984), and we generalize it to other profiles of the thermal source term $S$ than linear ones. As in Nieuwstadt (1984), we add a closure hypothesis in this section that is the Richardson numbers have reached their critical values $Ri_1^q$ and $Ri_2^q$. Starting with first-moment (18), (19), (20), supplemented with the equations on the second-order moments proposed, e.g., by Brost and Wyngaard (1978),<fnref n="4" /> which read in steady and horizontally homogeneous flows for respectively $q^2 = \overline{u'^2} + \overline{v'^2} + \overline{w'^2}$, $\overline{w'^2}$, $\overline{u'^w}$, $\overline{v'^w}$, $\overline{w'^\theta}$, $\overline{\theta'^2}$, $\overline{u'^\theta}$ and $\overline{v'^\theta}$:

$$0 = \left( -\overline{u'w'} \frac{\partial \overline{u}}{\partial z} - \overline{v'w'} \frac{\partial \overline{v}}{\partial z} \right) + \frac{g}{\theta_0} \overline{w'\theta'} - C_\epsilon \frac{q^3}{l}, \quad (46)$$

$$0 = \frac{2}{3} C_2 \left( -\overline{u'w'} \frac{\partial \overline{u}}{\partial z} - \overline{v'w'} \frac{\partial \overline{v}}{\partial z} \right) + \left( 2 + \frac{4}{3} C_3 \right) \frac{g}{\theta_0} \overline{w'\theta'} - \frac{2}{3} C_\epsilon \frac{q^3}{l} - C C_\epsilon \frac{q}{l} \left( \overline{w'^2} - \frac{q^2}{3} \right), \quad (47)$$

$$0 = \left( (1 - C_2) \overline{w'^2} - \frac{C_1}{2} q^2 \right) \frac{\partial \overline{u}}{\partial z} - (1 + C_3) \frac{g}{\theta_0} \overline{u'\theta'} + C C_\epsilon \frac{q}{l} \overline{u'w'}, \quad (48)$$

$$0 = \left( (1 - C_2) \overline{w'^2} - \frac{C_1}{2} q^2 \right) \frac{\partial \overline{v}}{\partial z} - (1 + C_3) \frac{g}{\theta_0} \overline{v'\theta'} + C C_\epsilon \frac{q}{l} \overline{v'w'}, \quad (49)$$

$$0 = - \overline{w'^2} \frac{\partial \overline{\theta}}{\partial z} + (1 - a_1) \frac{g}{\theta_0} \overline{\theta'^2} - d C_\epsilon \frac{q}{l} \overline{w'\theta'}, \quad (50)$$

$$0 = - \overline{w'\theta'} \frac{\partial \overline{\theta}}{\partial z} - C_{\theta} C_{\epsilon} \frac{q}{l} \overline{\theta'^2}, \quad (51)$$

<footnote n="4">⁴There is a typo error in Nieuwstadt (1984), a minus sign was missing in front of $a_2$.</footnote>
$$0 = \left( -\overline{u'w'} \frac{\partial \bar{\theta}}{\partial z} - (1 - a_2) \overline{w'\theta'} \frac{\partial \overline{u}}{\partial z} \right) - dC_{\epsilon} \frac{q}{l} \overline{u'\theta'}, \quad (52)$$

$$0 = \left( -\overline{v'w'} \frac{\partial \bar{\theta}}{\partial z} - (1 - a_2) \overline{w'\theta'} \frac{\partial \overline{v}}{\partial z} \right) - dC_{\epsilon} \frac{q}{l} \overline{v'\theta'}, \quad (53)$$

where $l$ denotes the turbulence length scale and with the following constants:

$$\left\{ \begin{array}{lll} C_1 = 0.207, & C_2 = 0.1, & C_3 = -0.8, \\ C = 6, & C_{\epsilon} = 0.139, & d = 9.7, \\ a_2 = 0.5, & C_{\theta} = 1.4, & a_1 = 0.5 + 1.5Ri_g^2 - Ri_g^3. \end{array} \right.$$

One can check from Eqs. (40), (41), (44) and (45) that the following relationship holds:

$$\overline{u'w'} \frac{\partial \overline{v}}{\partial z} - \overline{v'w'} \frac{\partial \overline{u}}{\partial z} = 0. \quad (54)$$

Therefore, Richardson numbers (21) and (22) can be rewritten as:

$$Ri_g = \frac{g \frac{\partial \bar{\theta}}{\partial z}}{\frac{\bar{\theta}_0}{\Phi^2}}, \quad (55)$$

$$Ri_f = \frac{g \overline{w'\theta'}}{\frac{\theta_0}{u_*^2 \Phi}}, \quad (56)$$

where $u_*^2$ is the kinematic shear stress and $\Phi$ is the norm of the derivative of the wind speed with respect to the height defined by:

$$\begin{aligned} u_*^2 &= \left( \overline{u'w'}^2 + \overline{v'w'}^2 \right)^{1/2}, \\ \Phi &= \left( \left( \frac{\partial \overline{u}}{\partial z} \right)^2 + \left( \frac{\partial \overline{v}}{\partial z} \right)^2 \right)^{1/2}. \end{aligned}$$

Note that both $u_*$ and $\Phi$ vary with $z$.

To solve system (18), (19) and (20), it is easier to use complex numbers. We denote the following dimension-less variables:

$$\Psi = - \frac{\overline{u'w'}}{u_{*0}^2} - i \frac{\overline{v'w'}}{u_{*0}^2}, \quad (57)$$

$$\Upsilon = \frac{\overline{w'\theta'}}{\overline{w'\theta'}_0}, \quad (58)$$
$$ \eta = \frac{z}{z_i}, \qquad (59) $$

$$ s = \frac{L}{u_{*0}} \left( \frac{\partial \bar{u}}{\partial z} + i \frac{\partial \bar{v}}{\partial z} \right), \qquad (60) $$

$$ \Theta = \frac{L}{\theta_{*0}} \frac{\partial \bar{\theta}}{\partial z}, \qquad (61) $$

where $u_{*0}$ is still the friction velocity within the surface layer (constant), $\theta_{*0}$ is the friction temperature within the surface layer (constant) and $w'\theta'_{0}$ is the turbulent heat flux within the surface layer (constant). Using (49), Eq. (20) is rewritten:

$$ \frac{d\Upsilon}{d\eta} = \frac{S z_i}{w'\theta'_{0}}. \qquad (62) $$

Differentiating (18)+i(19) with respect to z leads to:

$$ \frac{d^2\Psi}{d\eta^2} - ic^2s = 0, \qquad (63) $$

where c is a dimensionless number defined by:

$$ c = z_i \sqrt{\frac{f}{u_{*0}L}}. \qquad (64) $$

Finally, using (49) and (46), the Richardson numbers can be rewritten:

$$ Ri_g = \frac{\Theta}{\kappa |s|^2}, \qquad (65) $$

$$ Ri_f = \frac{\Upsilon}{\kappa \Psi^* s}, \qquad (66) $$

where $|\cdot|$ is the norm of the complex number and $(\cdot)^*$ is the complex conjugate.

The goal is now to solve the four equation system (54)–(58). The boundary conditions are as follows:

$$ \Psi = 1 \quad \Upsilon = 1 \quad \text{for } \eta = 0, $$

$$ \Psi = 0 \quad \Upsilon = 0 \quad \text{for } \eta = 1. $$

The first condition allows us to find the surface layer properties at the ground and the second one allows to let the turbulence vanishes at the top of the atmospheric boundary layer. First, since we made the assumption that S follows a power law with height, Eq. (54) easily leads to:
$$ \Upsilon = (1 - \eta)^{\alpha}, \qquad (67) $$

which is equivalent to $\frac{w'\theta'}{w'\theta'_0} = (1 - \frac{z}{z_i})^\alpha$.

Then, after eliminating $s$ from (55) with the expression of the flux Richardson number (58) and using the expression on $\Upsilon$ (59), we get the new differential equation on $\Psi$:

$$ \Psi^* \frac{d^2\Psi}{d\eta^2} - i \frac{c^2}{\kappa R i_f} (1 - \eta)^\alpha = 0. \qquad (68) $$

We search for the solution of (60) under the form:

$$ \Psi = (1 - \eta)^{a+i b}. \qquad (69) $$

where $a$ and $b$ are real constants. The differential equation implies the following three relations between $a$, $b$, $c$ and $\alpha$:

$$ \begin{aligned} a &= \frac{\alpha + 2}{2}, \\ b &= \frac{\sqrt{\alpha(\alpha + 2)}}{2}, \\ c &= \left( \frac{(\alpha + 1)\sqrt{\alpha(\alpha + 2)}}{2} \kappa R i_f^c \right)^{1/2}. \end{aligned} \qquad (70) $$

Then, we deduce that $\Psi$ is the solution of the problem if $c$ fulfils the third relation, which gives a relation on the atmospheric boundary layer height:

$$ z_i = \left( \frac{(\alpha + 1)\sqrt{\alpha(\alpha + 2)}}{2} \kappa R i_f^c \right)^{1/2} \sqrt{\frac{u_{*0}L}{f}}. \qquad (71) $$

To continue, taking the norm of $\Psi$, we find the expression of the friction velocity:

$$ u_*^2(z) = u_{*0}^2 \left(1 - \frac{z}{z_i}\right)^{(\alpha+2)/2}. \qquad (72) $$

Using the expression of the flux Richardson number (58), we deduce the expression of $s$:

$$ s = \frac{1}{\kappa R i_f^c} (1 - \eta)^{\frac{\alpha}{2} - 1 + i \frac{\sqrt{\alpha(\alpha + 2)}}{2}}. \qquad (73) $$

With gradient Richardson number expression (57), it follows the expression of $\Theta$:

$$ \Theta = \frac{1}{\kappa^2 (R i_f^c)^2} (1 - \eta)^{\alpha - 2}. \qquad (74) $$
Taking the real part and the imaginary part of s, we find the expressions of the vertical distributions of each wind speed components and the expression of Θ gives the vertical distribution of the temperature:

$$ 
\begin{aligned}
\frac{\partial \bar{u}}{\partial z} &= \frac{u_{*0}}{\kappa L R i_f^c} \left(1 - \frac{z}{z_i}\right)^{\alpha/2-1} \cos\left(\frac{\sqrt{\alpha(\alpha+2)}}{2} \ln\left(1 - \frac{z}{z_i}\right)\right), \\
\frac{\partial \bar{v}}{\partial z} &= \frac{u_{*0}}{\kappa L R i_f^c} \left(1 - \frac{z}{z_i}\right)^{\alpha/2-1} \sin\left(\frac{\sqrt{\alpha(\alpha+2)}}{2} \ln\left(1 - \frac{z}{z_i}\right)\right), \\
\frac{\partial \bar{\theta}}{\partial z} &= \frac{\theta_{*0} R i_f^c}{\kappa L (R i_f^c)^2} \left(1 - \frac{z}{z_i}\right)^{\alpha-2}.
\end{aligned}
\quad (75) $$

Presentation of Some Universal Functions Used in this Study

For the sake of completeness, here we recall the formulae of the universal functions used in this article. $\varphi_h$ functions are given only for information, but are not used in the comparison to experiment part of this article. If one wants to compare them in a correct manner, same $Pr_{T,0}$ should be used (see Yaglom (1977); Foken and Mauder (2024)). Carl et al. (1973) (C73):

$$ 
\begin{aligned}
\varphi_m(\zeta) &= (1 - 16\zeta)^{-1/3}, \quad \zeta < 0 \\
\varphi_h(\zeta) &= 0.74(1 - 16\zeta)^{-1/2}, \quad \zeta < 0
\end{aligned}
\quad (76) $$

Högström (1988) (HO88):

$$ 
\begin{aligned}
\varphi_m(\zeta) &= \begin{cases} (1 - 19.3\zeta)^{-1/4}, & \zeta < 0 \\ 1 + 4.8\zeta, & 0 \le \zeta < 0.5 \\ 7.9 - \frac{4.25}{\zeta} + \frac{1}{\zeta^2}, & 0.5 \le \zeta < 10 \\ 0.7485\zeta, & 10 \le \zeta \end{cases} \\
\varphi_h(\zeta) &= \begin{cases} 0.95(1 - 11.6\zeta)^{-1/2}, & \zeta < 0 \\ 0.95 + 7.8\zeta, & 0 \le \zeta \end{cases}
\end{aligned}
\quad (77) $$

Hartogensis and Bruin (2005) (HB05):

$$ 
\begin{aligned}
\varphi_m(\zeta) &= 1 + \zeta \left( a + b \exp(-d\zeta) - bd \left(\zeta - \frac{c}{d}\right) \exp(-d\zeta) \right), & 0 \le \zeta \\
\varphi_h(\zeta) &= 1 + \zeta \left( a \left(1 + \frac{2}{3}a\zeta\right)^{1/2} + b \exp(-d\zeta) - bd \left(\zeta - \frac{c}{d}\right) \exp(-d\zeta) \right), & 0 \le \zeta
\end{aligned}
\quad (78) $$

where $a = 1, b = 2/3, c = 5, d = 0.35$.

Cheng and Brutsaert (2005) (CB05):
$$ 
\begin{aligned}
\varphi_m(\zeta) &= 1 + a \left( \frac{\zeta + \zeta^b (1 + \zeta^b)^{(1-b)/b}}{\zeta + (1 + \zeta^b)^{1/b}} \right), \quad 0 \le \zeta \\
\varphi_h(\zeta) &= 1 + c \left( \frac{\zeta + \zeta^d (1 + \zeta^d)^{(1-d)/d}}{\zeta + (1 + \zeta^d)^{1/d}} \right), \quad 0 \le \zeta
\end{aligned}
\qquad (79) $$

where $a = 6.1, b = 2.5, c = 5.3, d = 1.1$.

With these universal functions, integrating Eqs. (2) and (3) from $z_0$ to $z$ lead to the expressions of wind speed and temperature:

$$ 
\begin{aligned}
\bar{u}(z) &= \frac{u_*}{\kappa} \Psi_m(\zeta_0, \zeta) \\
\bar{\theta}(z) &= \bar{\theta}(z=0) + \frac{\theta_*}{\kappa} \Psi_h(\zeta_0, \zeta)
\end{aligned}
\qquad (80) $$

where

$$ \Psi_{\alpha}(\zeta_0, \zeta) = \int_{z_0}^{z} \varphi_{\alpha}\left(\frac{x}{L}\right) \frac{d^{\frac{x}{L}}}{\frac{x}{L}} \qquad (81) $$

is the integrated form of the universal function $\varphi_{\alpha}$, $\alpha = m, h$, where $\zeta_0 = z_0/L$. HO88 gives:

$$ \Psi_m(\zeta, \zeta_0) = \begin{cases} \ln\left(\frac{\zeta}{\zeta_0}\right) - 2\ln\left(\frac{1+X}{1+X_0}\right) - \ln\left(\frac{1+X^2}{1+X_0^2}\right) + 2\arctan X - 2\arctan X_0, & \zeta < 0 \\ \ln\left(\frac{\zeta}{\zeta_0}\right) + 4.8(\zeta - \zeta_0), & 0 \le \zeta < 0.5 \\ 7.9\ln(2\zeta) + \frac{4.25}{\zeta} - \frac{0.5}{\zeta^2} - \ln(2\zeta_0) - 4.8\zeta_0 - 4.1, & 0.5 \le \zeta < 10 \\ 0.7485\zeta + 7.9\ln 20 - 11.165 - \ln(2\zeta_0) - 4.8\zeta_0, & 10 \le \zeta \end{cases} \qquad (82) $$

where $X = (1 - 19.3\zeta)^{1/4}$ and $X_0 = (1 - 19.3\zeta_0)^{1/4}$

$$ \Psi_h(\zeta, \zeta_0) = \begin{cases} 0.95 \left( \ln \left( \frac{\zeta}{\zeta_0} \right) - 2 \ln \left( \frac{1+X}{1+X_0} \right) \right), & \zeta < 0 \\ 0.95 \ln \left( \frac{\zeta}{\zeta_0} \right) + 7.8(\zeta - \zeta_0), & 0 \le \zeta \end{cases} \qquad (83) $$

where $X = (1 - 11.6\zeta)^{1/2}$ and $X_0 = (1 - 11.6\zeta_0)^{1/2}$ C73 gives:

$$ \Psi_m(\zeta, \zeta_0) = \ln\left(\frac{\zeta}{\zeta_0}\right) - 2\ln\left(\frac{1+X}{1+X_0}\right) - \ln\left(\frac{1+X^2}{1+X_0^2}\right) + 2\arctan X - 2\arctan X_0, \quad \zeta < 0 \qquad (84) $$

where $X = (1 - 16\zeta)^{1/3}$ and $X_0 = (1 - 16\zeta_0)^{1/3}$

$$ \Psi_h(\zeta, \zeta_0) = 0.74 \left( \ln \left( \frac{\zeta}{\zeta_0} \right) - 2 \ln \left( \frac{1+X}{1+X_0} \right) \right), \quad \zeta < 0 \qquad (85) $$

where $X = (1 - 16\zeta)^{1/2}$ and $X_0 = (1 - 16\zeta_0)^{1/2}$ HB05 gives:

$$ \Psi_m(\zeta, \zeta_0) = \ln\left(\frac{\zeta}{\zeta_0}\right) + a(\zeta - \zeta_0) + b\left(\zeta - \frac{c}{d}\right) \exp(-d\zeta) - b\left(\zeta_0 - \frac{c}{d}\right) \exp(-d\zeta_0), \quad 0 \le \zeta \qquad (86) $$
$$ \Psi_h(\zeta, \zeta_0) = \ln\left(\frac{\zeta}{\zeta_0}\right) + \left(1 + \frac{2}{3}a\zeta\right)^{3/2} + b\left(\zeta - \frac{c}{d}\right)\exp(-d\zeta) - \left(1 + \frac{2}{3}a\zeta_0\right)^{3/2} - b\left(\zeta_0 - \frac{c}{d}\right)\exp(-d\zeta_0), \quad 0 \le \zeta \le \zeta_0 \quad (87) $$

where $a = 1, b = 2/3, c = 5, d = 0.35$. CB05 gives:

$$ \Psi_m(\zeta, \zeta_0) = \ln\left(\frac{\zeta}{\zeta_0}\right) + a \ln\left(\zeta + (1 + \zeta^b)^{1/b}\right) - a \ln\left(\zeta_0 + (1 + \zeta_0^b)^{1/b}\right) \quad (88) $$

$$ \Psi_h(\zeta, \zeta_0) = \ln\left(\frac{\zeta}{\zeta_0}\right) + c \ln\left(\zeta + (1 + \zeta^d)^{1/d}\right) - c \ln\left(\zeta_0 + (1 + \zeta_0^d)^{1/d}\right) \quad (89) $$

where $a = 6.1, b = 2.5, c = 5.3, d = 1.1$.

**Acknowledgements** The authors would like to acknowledge SIRTA for providing the measurements used in this study.

**Author Contributions** M.F. and R.P. wrote the main manuscript text and R.P. prepared figures. M.F. and R.P. developed the theory. All authors reviewed the manuscript.

**Funding** Open access funding provided by Ecole Nationale des Ponts et Chaussees.

**Data Availability** No datasets were generated or analysed during the current study.

## Declarations

**Conflict of Interest** The authors declare no conflict of interest.

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

## References

Blackadar AK (1962) The vertical distribution of wind and turbulent exchange in a neutral atmosphere. J Geophys Res 67(8):3095-3102. https://doi.org/10.1029/JZ067i008p03095

Brost RA, Wyngaard JC (1978) A model study of the stably stratified planetary boundary layer. J Atmos Sci 35(8):1427-1440. https://doi.org/10.1175/1520-0469(1978)035(1427:AMSOTS)2.0.CO;2

Businger JA, Wyngaard JC, Izumi Y, Bradley EF (1971) Flux-profile relationships in the atmospheric surface layer. J Atmos Sci 28(2):181-189. https://doi.org/10.1175/1520-0469(1971)028$lt$0181:FPRITA$gt$2.0.CO;2

Carl DM, Tarbell TC, Panofsky HA (1973) Profiles of wind and temperature from towers over homogeneous terrain. J Atmos Sci 30(5):788-794. https://doi.org/10.1175/1520-0469(1973)030$lt$0788:POWATF$gt$2.0.CO;2

Cheng Y, Brutsaert W (2005) Flux-profile relationships for wind speed and temperature in the stable atmospheric boundary layer. Bound-Layer Meteorol 114(3):519-538. https://doi.org/10.1007/s10546-004-1425-4
Deardorff JW (1972) Numerical investigation of neutral and unstable planetary boundary layers. J Atmos Sci 29(1):91–115

Duynkerke PG (1988) Application of the E - $\epsilon$ turbulence closure model to the neutral and stable atmospheric boundary layer. J Atmos Sci 45(5):865–880. https://doi.org/10.1175/1520-0469(1988)045$lt$0865:AO TTCM$gt$2.0.CO;2

Dyer A (1974) A review of flux-profile relationships. Bound-Layer Meteorol 7:363–372

Ekman VW (1905) On the influence of the earth's rotation on ocean-currents. Arkiv för matematik, astronomi och fysik

Foken T (2006) 50 years of the Monin-Obukhov similarity theory. Bound-Layer Meteorol 119:431–447

Foken T, Mauder M (2024) Microclimatology. In: Micrometeorology. Springer, pp 331–351

Gryning SE, Batchvarova E, Brümmer B, Jørgensen H, Larsen S (2007) On the extension of the wind profile over homogeneous terrain beyond the surface boundary layer. Bound-Layer Meteorol 124(2):251–268

Haeffelin M, Barthès L, Bock O, Boitel C, Bony S, Bouniol D, Chepfer H, Chiriaco M, Cuesta J, Delanoë J et al (2005) Sirta, a ground-based atmospheric observatory for cloud and aerosol research. Annales Geophysicae Copernicus GmbH 23:253–275. https://doi.org/10.5194/angeo-23-253-2005

Hartogensis OK, Bruin HARD (2005) Monin-Obukhov similarity functions of the structure parameter of temperature and turbulent kinetic energy dissipation rate in the stable boundary layer. Bound-Layer Meteorol 116(2):253–276. https://doi.org/10.1007/s10546-004-2817-1

Högström ULF (1988) Non-dimensional wind and temperature profiles in the atmospheric surface layer: a re-evaluation. In: hicks BB (ed) topics in micrometeorology. A festschrift for arch dyer. Springer, Dor- drecht, pp 55–78. https://doi.org/10.1007/978-94-009-2935-7-6

Högström U (1996) Review of some basic characteristics of the atmospheric surface layer. Bound-Layer Meteorol 78(3):215–246

Holtslag AA, Nieuwstadt FT (1986) Scaling the atmospheric boundary layer. Bound-Layer Meteorol 36(1):201–209

Jacobson MZ (2005) Fundam atmos modeling, 2nd edn. Cambridge University Press

Kelly M, Van Der Laan MP (2023) From shear to veer: theory, statistics, and practical application. Wind Energy Sci Discuss 2023:1–37

Kelly M, Cersosimo RA, Berg J (2019) A universal wind profile for the inversion-capped neutral atmospheric boundary layer. Q J R Meteorol Soc 145(720):982–992

Monin AS, Yaglom AM (2007) Statistical fluid mechanics: mechanics of turbulence. Courier Corporation, google-Books-ID: EtTyy14 DvIC

Monin AS, Obukhov AM (1954) Basic laws of turbulent mixing in the surface layer of the atmosphere. Con- trib Geophys Inst Acad Sci USSR 151(163):e187

Narasimhan G, Gayme DF, Meneveau C (2024) Analytical model coupling ekman and surface layer structure in atmospheric boundary layer flows. Bound-Layer Meteorol 190(4):16

Nieuwstadt FTM (1984) The turbulent structure of the stable. Noct Bound Layer J Atmos Sci 41(14):2202–2216. https://doi.org/10.1175/1520-0469(1984)041$lt$2202:TTSOTS$gt$2.0.CO;2

Panofsky HA (1963) Determination of stress from wind and temperature measurements. Q J R Meteorol Soc 89(379):85–94

Panofsky HA, McCormick RA (1960) The spectrum of vertical velocity near the surface. Q J R Meteorol Soc 86(370):495–503

Panofsky HA, Tennekes H, Lenschow DH, Wyngaard J (1977) The characteristics of turbulent velocity com- ponents in the surface layer under convective conditions. Bound-Layer Meteorol 11(3):355–361

Pasquill F (1974) Atmospheric diffusion: the dispersion of windborne material from industrial and other sources. Wiley, New York

Peña A, Gryning SE, Hasager CB (2010) Comparing mixing-length models of the diabatic wind profile over homogeneous terrain. Theoret Appl Climatol 100(3):325–335

Richards PJ, Norris SE (2019) Appropriate boundary conditions for computational wind engineering: still an issue after 25 years. J Wind Eng Ind Aerodyn 190:245–255

Stull RB (1988) An introduction to boundary layer meteorology. Springer, Dordrich. https://doi.org/10.1007/978-94-009-3027-8

Wieringa J (1980) A revaluation of the kansas mast influence on measurements of stress and cup anemometer overspeeding. Bound-Layer Meteorol 18(4):411–430

Wilson J (2008) Monin-Obukhov functions for standard deviations of velocity. Bound-Layer Meteorol 129(3):353–369

Yaglom A (1977) Comments on wind and temperature flux-profile relationships. Bound-Layer Meteorol 11(1):89–102

Yokoyama O, Gamo M, Yamamoto S (1979) The vertical profiles of the turbulence quantities in the atmo- spheric boundary layer. J Meteorol Soc Jpn Ser II 57(3):264–272
Zilitinkevich SS (1972) On the determination of the height of the Ekman boundary layer. Bound-Layer Meteorol 3(2):141–145. https://doi.org/10.1007/BF02033914

**Publisher's Note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.