Characterization of Turbulence in the Neutral and Stable Surface Layer at Jang Bogo Station, Antarctica

Abstract

The availability of 5-year time series of velocity and temperature data from two sonic anemometers installed at Jang Bogo Station, Antarctica, allowed a systematic investigation of the turbulence features in a stable layer affected by submeso motions and characterized by the vertical divergence of some second-order moments for a large fraction of time (quite a non-ideal surface layer). The investigation of the effect of the averaging time interval on the statistics of the second-order moments showed that this is greater for the variances of the velocity components with respect to that for the vertical fluxes. This corresponds to a greater contribution from low-frequency motions. The turbulence statistics were investigated and compared with current literature results in terms of vertical structure, share of energy between horizontal and vertical components, skewness of the vertical velocity and turbulent velocities. As a general result, all the normalized second-order moments show a clear change passing from neutral to stable conditions, passing through the range of bulk Richardson number equal to 0.1–1.

1. Introduction

The stable boundary layers call for a number of exceptions in the application of the Monin–Obukhov similarity theory (MOST). This is based on the hypothesis of the driving effect of vertical fluxes of momentum and heat, and of the empirical validity of the flux–gradient relationship. MOST can be understood as a local theory [1], i.e., eddies enable the local exchange of properties of the flow. Grachev et al. [2] noted that beyond a critical value of the gradient Richardson number, the inertial subrange associated with the Richardson–Kolmogorov cascade dies out, and that this is a limit in the application of MOST.

Recently, the hockey-stick transition (HOST) paradigm has been developed to interpret the dependence of turbulence velocity scales on local wind velocity [3,4,5], i.e., on the bulk shear.

The increasing amount of turbulence observations in the atmospheric boundary layer challenges our skill in the understanding of a variety of phenomena, and stimulates the development of new approaches and new parameterizations.

The observations give evidence of the absence of a critical value for the gradient of the bulk Richardson number, i.e., turbulence does not disappear beyond such a threshold. The share of the turbulent kinetic energy among vertical and horizontal components (or the ratio between the vertical and horizontal velocity variances) at the inertial scales also changes with respect to the typical values in shear-dominated turbulence, the horizontal components being dominant with respect to the vertical one. This topic is widely discussed for instance by Zilitinkevich et al. [6] (Figures 4–6) and by Cheng et al. [7] (Sect. 6) who examined many different datasets and numerical simulations in the frame of the formulations of closure models without a critical Richardson number.

Beyond the classical assumption of constant fluxes in the surface layer, the vertical divergence of second-order moments has been discussed by many authors. Högström et al. [8] and Mahrt and Vickers [9] studied the variance of the vertical component of velocity, while Mahrt et al. [10] studied that of the friction velocity, i.e., the vertical momentum flux.

Last but not least, the effect of the averaging time interval, that allows to include different ranges of frequencies in the sample used for evaluating (co)variances of the turbulent quantities and thus to account for different scales of motion, is one of the open topics.

In this work, we will deal with neutral to stable conditions, in line with a number of contributions in the literature, and give an overall characterization of the surface layer turbulence in a coastal site at the Korean Jang Bogo Station, Antarctica. In particular, in this paper the vertical divergence of the second-order moments and the dependence of some second-order moments properly normalized on stability (measured by the Richardson number) are considered. An observation of the conditions under which the local similarity (in basic terms, MOST) works, and those under which large eddies are driving (with HOST as a paradigm) is given. Moreover, comparing the results obtained using different averaging times, the influence of low-frequency (submeso) on the statistics was investigated.

2. Background

2.1. Definitions

The capital letters U and $\Theta$ indicate the averaged values of the wind velocity and the potential temperature, while the small letters u, v, w or $u_i$ and $\theta$ indicate the fluctuations around these averages (by a proper choice of axes, $V=W=0$). Brackets indicate the time averaging.

The bulk Richardson number (e.g., [11]) reads:

$$R_B\equiv \frac{(z_2-z_1)\beta \Delta \Theta (z_1,z_2)}{\Delta U^2(z_1,z_2)}$$

(1)

where $z_1=2m$, $z_2=9m$ are the measurement heights, and $\beta =g/{\Theta }_0$, where g is the gravity acceleration and ${\Theta }_0$ is a reference value for the potential temperature. In the following analysis, we shall refer when appropriate to two ranges of the bulk Richardson number: $0.001<R_B<0.01$, near-neutral conditions, and $1<R_B<10$, stable conditions.

The Obukhov length reads:

$$L=-\frac{{\tau }^{3/2}}{\kappa \beta \langle w\theta \rangle }$$

(2)

where $\tau ={[{\langle uw\rangle }^2+{\langle vw\rangle }^2]}^{1/2}$ is the absolute value of the vertical momentum flux and $\langle w\theta \rangle$ is the vertical heat flux; $\kappa$ is the von Kármán constant.

2.2. Remarks about MOST

According to MOST (Monin–Obukhov similarity theory, [12]), the non-dimensional gradients of mean velocity and potential temperature are expressed as a function of L, but the actual form of the function is not specified. The common log-linear relationships read:

$$\frac{1}{{\tau }^{1/2}}\frac{\mathrm{d}U}{\mathrm{d}z}=\frac{1}{\kappa z}+\frac{\alpha }{L}$$

(3)

$$\frac{{\tau }^{1/2}}{-\langle w\theta \rangle }\frac{\mathrm{d}\Theta }{\mathrm{d}z}=\frac{{\beta }_0}{\kappa z}+\frac{{\beta }_1}{L}$$

(4)

According to Högström [13] (see also [11], pages 94–96), we used $\alpha =6$, ${\beta }_0=0.95$, ${\beta }_1=8$ and $\kappa =0.37$. Note that, in the literature, other values are proposed: for instance $\alpha =5$ [14] and ${\beta }_1=5$ [15]; Yagüe et al. [16] also found better agreement for the non-dimensional temperature gradient using ${\beta }_1\simeq 5$. However, no relevant differences arise in the present context.

These gradients are local [1], and can be integrated to give the velocity and temperature differences in a layer, prescribing some form for the flux profiles. It can be shown that in a shallow layer, the effect of changing the fluxes is quite small—well below the experimental variability. Thus, we can use the integrated profiles assuming fluxes’ constant with height and thus derive a relationship between $R_B$ and L.

It is well known that the expressions in Equations (3) and (4) show a larger increase than the observations for large stability, and predict an upper limit to the bulk Richardson number (a critical value) which depends on the values of the coefficients and on the heights used in the computations (note that if the gradient Richardson number is concerned, the critical value is of course independent from height). Note also that formulations that allow an unlimited Richardson number have been proposed (see for instance, [17], Equations 32 and 33).

In order to anticipate some results from HOST (see Section 2.3 below), it is useful to design a paradigm relating the turbulent kinetic energy and the wind velocity. The formula in Equation (3) states that the square root of the vertical turbulent momentum flux ${\tau }^{1/2}\left(z\right)$ is proportional to the velocity gradient $\mathrm{d}U/\mathrm{d}z$, the proportionality being the constant function of the distance from the surface z and the Obukhov length L. Under near-neutral conditions $z\ll L$, and for a constant momentum flux, the equation can be integrated on z and results $u_{*}\left(z\right)\propto U\left(z\right)$:

$${\tau }^{1/2}\left(z\right)\equiv u_{*}\left(z\right)=\frac{U\left(z\right)\kappa }{ln(z/z_0)}$$

(5)

where $z_0$ is the roughness length.

In simplified form (steadiness, horizontal homogeneity, negligible divergence of third-order moments), the budget equation for the mean turbulent kinetic energy (TKE) $E_K$ [11] reads:

$$\begin{array}{c}\hfill \tau \frac{\mathrm{d}U}{\mathrm{d}z}\left(1+\frac{\beta \langle w\theta \rangle }{\tau \mathrm{d}U/\mathrm{d}z}\right)-ϵ=0\end{array}$$

(6)

In general, this budget does not state the amount of TKE, but only the equilibrium between production/destruction (the first term) and dissipation (the second term). An approximate step can be made considering that common parameterizations of the dissipation involve TKE itself. Here, we use the results of Basu et al. [18]: for near-neutral or moderately stable conditions, they suggest:

$$ϵ={\alpha }_B{E}_K\frac{\mathrm{d}U}{\mathrm{d}z}$$

(7)

with ${\alpha }_B\approx 0.23$ ([18], Equation 18a).

Using the log-linear velocity gradient Equation (3) and the dissipation Equation (7) the relation between $E_K$ and $\tau$ reads:

$$E_K=\frac{\tau }{{\alpha }_B}\left[1-\frac{z}{L}{\left(1+\alpha \frac{z}{L}\right)}^{-1}\right]$$

(8)

which means that $E_K^{1/2}\propto U\left(z\right)$, using Equation (5), under near neutral conditions.

By summarizing, MOST and the simplified TKE budget predict that $E_K^{1/2}$ or ${\tau }^{1/2}$ (called turbulence scale velocities) are proportional to the velocity gradient (a local relationship).

2.3. Some HOST Suggestions

Based on experimental observations, Sun et al. [3,4,5] formulated the so-called HOST paradigm, suggesting that Equations (5) and (8) approximately hold for wind velocity less than a threshold $V_S$, above which the turbulence scale velocities are a function of the local wind velocity, i.e., of the velocity bulk gradient $U\left(z\right)/z$: see for instance Sun et al. [5], Equations 4 and 5). Thus, for the square root of the vertical turbulent momentum flux:

$$\begin{array}{c}\hfill u_{*}\left(z\right)={\alpha }_{u1}\left(z\right)U\left(z\right);U\left(z\right)<V_S\left(z\right)\end{array}$$

(9)

$$\begin{array}{c}\hfill u_{*}\left(z\right)={\alpha }_{u2}\left(z\right)U\left(z\right)+{\beta }_u\left(z\right);U\left(z\right)>V_S\left(z\right)\end{array}$$

(10)

where ${\alpha }_{u1}$ is determined from Equation (3) and ${\alpha }_{u2}\approx 0.1$ (see [5], Figure 6): above the threshold velocity (which depends on the height, and also on the local features of the surface), the relation among scale velocities and wind velocity is determined by the large vortices of scale z and by the bulk shear $U\left(z\right)/z$—instead of the local one.

HOST states that similar relationships hold for TKE (let substitute $E_K^{1/2}$ to $u_{*}$, and the subscript $E_K$ to u in Equations (9) and (10)): below the threshold velocity, an ‘equilibrium’ relation (a relation consistent with MOST, see Equations (5) and (8)) holds; however, for a wind velocity larger than the threshold velocity, the TKE increases in a similar way to Equation (10): ${\alpha }_{E2}$ is approximately $0.28$ (from [3], Figure 1). Moreover, Sun et al. [3] noted that TKE does not tend towards zero as the velocity approaches zero, and identified as ‘regime 3’ a condition under which TKE is larger than the ‘equilibrium’ relation (see their Figure 2). Actually, this regime can be interpreted as the overwhelming effect of low-frequency motions. From the same paper, the threshold velocity can be estimated to be approximately 2 ms ${}^{-1}$ at $z=2$ m and 4 ms ${}^{-1}$ at $z=9$ m.

2.4. Remarks about Second-Order Moments

According to MOST, the normalized second-order moments can be expressed as functions of L (see, e.g., [11,19]). For instance, under near-neutral conditions (and under stable conditions, within the so-called z-less limit (see for instance [20])), the ratio between the velocity component variances and the momentum flux read:

$$\frac{\langle u_i^2\rangle }{\tau }={\gamma }_i$$

(11)

with ${\gamma }_i\approx 5.3,4,1.6$ for $i=1,2,3$, respectively, [11,19].

Under stable conditions, the normalized moments are expected to be dependent on stability, but because of self-correlation, the relations fitted on the data may give ambiguous results. To overcome the problem, Mauritsen and Svensson [21] used the bulk Richardson number as a stability parameter, and observed that $\tau /E$ and $\langle w\theta \rangle /{\left[E\langle {\theta }^2\rangle \right]}^{1/2}$ are both critically dependent on the Richardson number, showing a well-defined change in the Richardson number from 0.1 and 1. This behavior has been treated by Zilitinkevich et al. [6] and reviewed by Li et al. [22] in the frame of models without a critical Richardson number. If instead of taking under exam $E/\tau$ the normalized horizontal and vertical variances of the velocity are considered, the behavior may be different, as suggested for instance by Juang et al. [23].

3. Site and Dataset

The Jang Bogo station, managed by the Korea Polar Research Institute (KOPRI), was opened in February 2014 and is located in the coastal region of Terra Nova Bay, Antarctica (74${}^{\circ }$37${}^{\prime }$26${}^{″}$ S, 164${}^{\circ }$13${}^{\prime }$43${}^{″}$ E). The station is surrounded by hills up to ∼670 m high from the west to the northeast within a 2 km distance and most of land gently slopes close to the coast. This area is mainly composed of exposed bedrock and glacial moraines. The station is operative all year-round. During summer, with slightly positive temperatures, even for short periods, some of the area of the station grounds is exposed (i.e., free of snow), while all the area is covered with snow and ice during non-summer seasons.

An automated weather station (AWS) was installed when the station was opened, at approximately 150 m east from the main building, acting as a reference for the meteorological conditions at the station. Based on seasonal statistics obtained from single 10 min averages for the years 2016–2019, we found that temperatures are always well below freezing all year-round, except during the austral summer, when the seasonal averages are negative by a few degrees. On the contrary, the winter is very harsh, with average values lower than −20 ${}^{\circ }$C. Prevailing winds are from the north/west to the south/west sector during the whole year, while during summer this sector expand towards the south and even the east. Seasonal average speed ranges between 3 and 4 m/s, but the 90th percentile can be higher by as much as 10 m/s. Many other atmospheric observations are routinely carried out at the station, including for radiation (shortwave, longwave, UV), green-house gas and aerosols. A meteorological radiosonde is routinely launched once per day at noon (UTC), while an ozonesonde is routinely launched once per week.

A 9 m tower for meteorological and turbulence measurements was also installed in 2014 upwind from the station, in the west to northeast direction. Turbulence statistics and mean values for three wind components and temperature (corrected for crosswind and humidity effects) were obtained from two sonic anemometers (Campbell scientific model CSAT3), which were mounted facing west at a height of 2 and 9 m of the tower. As declared by the manufacturer, the resolution for the wind components u and v is 1 mm/s (RMS), 0.5 mm/s for w and 0.002 ${}^{\circ }$C for temperature (the last at 25 ${}^{\circ }$C, might increase the value at very low temperatures). The accuracy is only known for wind speed: <±8.0 cm/s (u, v), <±4.0 cm/s (w). Since air temperature does not fall below −40 ${}^{\circ }$C, the sonic anemometers worked year-round. Flow distortion by the tower is not likely significant for a 9 m-high sonic anemometer in any direction but could be significant for a 2 m-high one from the northeast to southeast. Therefore, we analyzed data for wind direction from 180${}^{\circ }$ to 360${}^{\circ }$.

Raw data at a sampling rate of 20 Hz from both the two sonic anemometers were recorded on a data logger, and then regularly retrieved to a computer in a laboratory. For this work, the data for 2016–2019 were analyzed over time intervals of 30, 10 and 2 min, in order to investigate the different behavior of fluctuations at different frequencies. As a general statement, shorter time intervals tend to filter out low-frequencies phenomena, such as submeso motions. These generally include a rich variety of motions on horizontal scales larger than the turbulence but smaller than a few kilometers. The 30 min results will include a large fraction of the low-frequency fluctuations, although not all of them. For the main part of this work, the two extreme cases (30 and 2 min intervals) were used in order to distinguish which are the quantities affected by submeso motions. A total of 20,329 and 407,003 cases were obtained for the two cases, respectively.

4. Results

4.1. Effect of the Averaging Time Interval

As a preliminary analysis, we discuss in short the effect of the length of the averaging time interval on the statistics we are concerned with.

As introduced above, using different averaging time intervals on a time series is equivalent to applying a box filter, which selects (even approximately) a frequency range in the spectrum ([24], Equation 13.1). Thus, using different intervals corresponds to select different frequency ranges, and we shall use this simple tool to investigate the presence and the effects of the low-frequency (submeso) motions on the turbulence statistics.

In

Table 1. Median values for horizontal and vertical velocity variances (in m${}^2$ s${}^{-2}$) measured at the two heights, for different stability conditions and different time intervals.

	30 min	10 min	2 min
$0.001<R_B<0.01$
$\langle u^2\rangle +\langle v^2\rangle$ (2 m)	4.2	3.8	3.0
$\langle u^2\rangle +\langle v^2\rangle$ (9 m)	4.5	3.9	2.6
$\langle w^2\rangle$ (2 m)	0.28	0.29	0.29
$\langle w^2\rangle$ (9 m)	0.33	0.33	0.33
$1<R_B<10$
$\langle u^2\rangle +\langle v^2\rangle$ (2 m)	0.80	0.40	0.13
$\langle u^2\rangle +\langle v^2\rangle$ (9 m)	1.04	0.48	0.14
$\langle w^2\rangle$ (2 m)	0.015	0.014	0.012
$\langle w^2\rangle$ (9 m)	0.035	0.029	0.019

, the median values for horizontal ($\langle u^2\rangle +\langle v^2\rangle$) and vertical ($\langle w^2\rangle$) velocity variances are reported for different averaging time intervals, for different stability conditions and for the two measuring levels. In passing, note that reading the table by row gives evidence of the effects of averaging time, while reading by column allows to guess the vertical divergence of the variances (see Section 4.2). Similarly, in

Table 2. As in Table 1, but for the vertical fluxes of momentum $\tau$ (in m${}^2$ s${}^{-2}$) and of heat $\langle w\theta \rangle$ (in K m s${}^{-1}$).

	30 min	10 min	2 min
$0.001<R_B<0.01$
$\tau$ (2 m)	0.22	0.22	0.22
$\tau$ (9 m)	0.19	0.19	0.18
$\langle w\theta \rangle$ (2 m)	−0.016	−0.016	−0.018
$\langle w\theta \rangle$ (9 m)	−0.028	−0.028	−0.029
$1<R_B<10$
$\tau$ (2 m)	0.009	0.008	0.006
$\tau$ (9 m)	0.024	0.018	0.011
$\langle w\theta \rangle$ (2 m)	−0.008	−0.008	−0.007
$\langle w\theta \rangle$ (9 m)	−0.010	−0.009	−0.006

, the median values for the vertical fluxes of momentum $\tau$ and heat $\langle w\theta \rangle$ are reported.

We already noticed that the increase in the length of the averaging time interval allows the statistics to include low-frequency motions. Under near-neutral conditions, $\langle w^2\rangle$ is practically unaffected by the length of the interval, i.e., it is mainly affected by the high-frequency (small-scale) eddies, while $\langle u^2\rangle +\langle v^2\rangle$ is reduced by approximately $30÷40$% (at least for the measurements at 2 m) from 30 min to 2 min averaging times, indicating a contribution by such eddies. Under stable conditions, the reduction is quite evident both for the horizontal and the vertical velocity variances, suggesting that low-frequency motions affect both, to a larger extent for $\langle u^2\rangle +\langle v^2\rangle$ with respect to $\langle w^2\rangle$, and for 9 m measurements with respect to 2 m.

As regards the fluxes, under near neutral conditions both are practically unaffected by the averaging time interval length at both heights. Under more stable conditions, again at both heights, the interval length affects the momentum flux, and to a smaller extent, the heat flux; as for overall observations, fluxes are less affected than the horizontal variance by the averaging time interval length. As a consequence, the estimate of the Obukhov length is less affected by the choice of the averaging time than the velocity variances, at least for this specific dataset (see Section 4.3).

4.2. Vertical Structure (Divergence) of the Layer

As noted above, an inspection by columns of the Table 1 and Table 2 suggests that the vertical divergence of the second-order moments is a common feature of this surface layer.

A quantitative evaluation of the divergence of a generic quantity x is given by the ratio between the difference of the values of x at the upper and lower levels of the observation and its mean value in the layer between the two levels:

$$I_x=\Delta x/\overline{x}.$$

(12)

The indexes for the variances of the horizontal and vertical components of the velocity and for the vertical fluxes as a function of the bulk Richardson number and for 30 min and 2 min averages are reported in Figure 1.

From Figure 1a (see also Table 1), the results show that, for 30 min averages, under almost neutral conditions, the vertical divergence is negligible on the median, although a 10% increment can be observed for $\langle w^2\rangle$, while under stable conditions, there is an increase in the horizontal velocity variance of approximately 20% in height and by more than 100% for the vertical velocity variance.

The different behavior of the index I for horizontal and vertical velocity variances under stable conditions is related to the buoyancy damping of the vertical motions, while horizontal motions are less affected, as can be seen [25].

From Figure 1b (see also Table 2), it is interesting to note that the vertical momentum flux $\tau$ slightly decreases in height under near neutral conditions, and in median (but the distribution is quite large) increases of more that 100% in stable cases. A different behavior is found for the heat flux, which shows an increase in height under near neutral conditions (though it is quite small in these cases) while it is almost constant on average under stable conditions.

With regard to the 2 min averages, as can be seen in Figure 1c,d, the results show that the general behavior is similar to that of the 30 min averages, but with values of divergence that are smaller on average and at the same time more dispersed, so much so that in the case of horizontal velocity, they pass from purely positive to a distribution with almost half the negative values. In other words, the vertical divergence is mainly associated with the low-frequency motions.

It may be worth noticing that Högström et al. [8] discussed the increase in vertical velocity variance in the surface layer as related to the blocking due to the surface.

4.3. The Obukhov Length and the Richardson Number

The inverse of the Obukhov length, computed using the mean values of the fluxes in the $2÷9$ m layer, is reported as a function of the bulk Richardson number in Figure 2. In order to present a non-dimensional variable, in the figure, the quantity $z/L$ with $z=2.5$ m is plotted; the height is quite arbitrarily chosen in the layer of interest and it must be noted that any choice does not affect the physical meaning of the observations. It results that $R_B$ steadily increases while $z/L$ diverges for $R_B\approx 0.2$. The empirical dependence of $z/L$ from $R_B$ is well described by the relation derived from the mean profiles for $R_B\le 0.1$; the log-linear profiles (by integration of Equations (3) and (4)) predict a critical value of approximately $0.21$, while the profiles suggested by Gryanik et al. [17] exhibit an increase for $R_B$ at $1<z/L<10$—in any case, their validity is limited to $z/L<100$.

Moreover, the different time averaging intervals almost do not affect the distribution of data, consistently with the observations that the fluxes show a minor influence of the length of the averaging interval.

This result supports and justifies the choice of Richardson number as a measure of stability: in fact, many observations would collapse in a small interval of $z/L$ around $0.25$, while being well distributed in the range $R_B>0.2$.

4.4. Turbulence Characterization

In this subsection, we analyze the behavior of some characteristic quantities as a function of the bulk Richardson number.

The anisotropy ratio A:

$$A=\frac{\langle u^2\rangle +\langle v^2\rangle }{\langle w^2\rangle }$$

(13)

is a measure of the share of the TKE between the horizontal and vertical components of velocity and is characterized by an almost universal value $A\approx 5÷6$ for the shear produced by turbulence (i.e., under near-neutral or weakly stable conditions), while it presents a large variability under strong stability conditions $A\approx 15÷30$ (as can be seen, for instance, in Zilitinkevich et al. ([6], Figure 6), Tampieri ([11], Figure 3.29)).

In the frame of this analysis, we limited ourselves to the examination of second-order moments, and did not use the spectra. Thus, the anisotropy ratio has a two-fold meaning: it increases its value from neutral to stable conditions, and also increases its value in the presence of submeso motions (which typically affect the horizontal more than the vertical velocity variances). In the following discussion, we shall outline how to separate (at least statistically) the two aspects.

The anisotropy ratio A is shown in Figure 3 as a function of $R_B$.

Under near neutral conditions (say, $R_B<0.02$), A is quite similar at the two levels, and its values for the 2 min averages are only slightly smaller than the ones for the 30 min averages. Broadly, the data are quite in agreement with the theoretical value from [6], with modified coefficients: Zilitinkevich and coworkers (private communication) did re-evaluate the coefficients of the 2013 model using further atmospheric observations and LES simulations, so the coefficients we use are slightly different from the published ones.

Under stable conditions (e.g., $R_B>0.5$), the ratio A for 30 min averages increases, and its variability also increases; the increase and the spread around the median are much more evident at the lowest level. This difference between the two levels occurs because $\langle w^2\rangle$ increases with height while $\langle u^2\rangle +\langle v^2\rangle$ much less, as discussed before. The A ratios for 2 min averaged data are quite similar to the ones of the near neutral case, thus practically independent on stability.

The theoretical values predicted by Zilitinkevich et al. [6] show fairly good agreement for near neutral conditions, and at larger $R_B$ only for the observations at $z=9$ m averaged over 30 min. A major change related to the increase in stability occurs at $z=2$ m for the 30 min averaged data, while the 2 min averaged data appear to be almost independent from stability.

In Figure 3, the wind velocity binned for $R_B$ intervals is also reported (black symbols). The wind velocity reaches the highest values under near neutral conditions and drops to very small values with increasing stability, thus supporting the empirical evidence of minor effects of submeso motions (which are a reason for the increase in the ratio A) under strong wind conditions. In fact, submeso motions are mainly detected under low-wind conditions (see, for instance [26,27,28]).

The ratio $\tau /E_K$ is shown in Figure 4 as a function of $R_B$ and compared with the Mauritsen and Svensson [21] empirical relation and the Zilitinkevich et al. [6] theoretical curve.

This ratio shows a relevant decrease with stability, in qualitative agreement with the literature, for both 30 min data and 2 min data. Note that the overall values of the ratio are larger for 2 min than for 30 min averages, and that our 30 min averages are smaller than those presented by Mauritsen and Svensson [21] which refer, for all but one dataset, to averaging times of 5 min or 10 min, supporting the view that TKE (to be more precise, the horizontal velocity variance) is more influenced than fluxes by the submeso motions.

Similar considerations hold for the normalized vertical heat flux $H=\langle w\theta \rangle /{\left[E_K\langle {\theta }^2\rangle \right]}^{1/2}$, shown in Figure 5. It is worth noticing that the normalized heat flux has a minimum for $0.01<R_B<0.1$, which is well documented by the observations [11,21] but not represented in the empirical curve by Mauritsen and Svensson [21] nor theoretically predicted by Zilitinkevich et al. [6].

The skewness of the vertical velocity $S_K=\langle w^3\rangle /{\langle w^2\rangle }^{3/2}$, averaged over the $2÷9$ m layer, for 30 min averaging time, is reported in Figure 6. It is positive under near neutral conditions, while under stable conditions, the median value approaches zero with a wide fluctuation range.

It is worth noticing that, under near-neutral conditions, the skewness is on median positive, basically due to the surface blocking effect. Consistently, the large variability (with zero median) under stable conditions corresponds to the reduction in the vertical scale length. Note also (right panel) that the skewness in near-neutral cases decreases with the vertical divergence of $\langle w^2\rangle$. From a flux–gradient perspective, skewness is expected to be zero for $I_w=0$. Actually, $\langle w^3\rangle /{\langle w^2\rangle }^{3/2}\approx 0$ at $I_w\approx 0.6$, again highlighting the importance of surface blocking. In the stable case, the skewness is in median near zero, consistent with the literature—as can be seen in Barberis [29], from the data of Higgins et al. [30]. It is worth noticing that the spread (25–75% percentiles) increases steadily with the index $I_{\langle w^2\rangle }$, consistent with the general consideration that the decoupling between the layers (and thus the vertical velocity distribution may be quite different from case to case) is the dominant factor characterizing such a surface layer.

4.5. The Behavior of the Turbulence Velocities

In this subsection, the behavior of $\sqrt{E_K}$ and of $u_{*}$, also called ‘turbulence velocities’, is analyzed as function of the wind velocity U, following the line of Sun et al. [3].

The plot of $u_{*}\left(z\right)$ vs. $U\left(z\right)$ is reported in Figure 7 and has the general features of the HOST paradigm:

• For near-neutral conditions:

MOST, i.e., Equation (9), is valid from approximately 4 m/s to approximately 8 m/s for $z=2$ m, to approximately 10 m/s for $z=9$ m; these threshold velocities are quite a bit larger than those observed in the CASES-99 experiment ([5], Figure 6); the roughness (derived from the fit of the log profile) is different for the two levels ($z_0\approx 2.5\times {10}^{-3}$ at 2 m, $z_0\approx {10}^{-3}$ at 9 m) suggesting the existence of an internal boundary layer;

The slope beyond the threshold velocity is $\approx 0.06$ (smaller than the $0.1$ value shown by Sun et al. ([5], Figure 6));

• For all the stability conditions, and for both the averaging times, $\tau \ne 0$ as $U\to 0$.

All these observations are the same for both 30 min and 2 min averaging intervals (not shown), quite consistent with the observation made before about the minor influence of the averaging time on the vertical momentum flux.

The relationship between $E_K^{1/2}$ and $U\left(z\right)$ is exploited in Figure 8, to be compared with Sun et al. ([3], Figure 1).

Under near neutral conditions, beyond the threshold velocity, the slope is approximately $0.17$ at $z=2$ m and $0.15$ at $z=9$ m for a 30 min averaging time (slightly less for 2 min averaging time at $z=9$ m), smaller than the value resulting from Sun et al. [3]. Under stable conditions, $E^{1/2}$ does not increase with the wind velocity for a 30 min averaging time. In general, the observations also support the existence of some level of TKE for the wind velocity decreasing to zero, associated with the low-frequency motions and not to the shear.

These results confirm the general HOST frame and at the same time stress the fact that the numerical values of the threshold velocity and of the slope are (at least) site dependent.

A final remark about turbulence velocities is suggested by the application of Equation (8), as reported in Figure 9.

We obtain the result that the data for near-neutral conditions approximately follow Equation (8), neglecting the values for small friction velocity, say for $u_{*}<0.1$, with a slope giving ${\alpha }_B\approx 0.22$ for the 2 min averaging time, in quite good agreement with the value suggested by Basu et al. [18]. For stable conditions, the TKE of 30 min averaged data is larger than that of the near-neutral ones, suggesting that part of the TKE is supplied by sources different from the shear production, while the 2 min averaged data show less TKE than the near-neutral ones, because the influence of the low-frequency motions is small or negligible.

5. Conclusions

The statistical analysis of the 5-year dataset of turbulence measurements at Jang Bogo Station confirmed some literature results and evidenced some important new features.

The effects of the averaging time were discussed—it was shown to be more relevant to the variances of the velocity components than on the vertical fluxes. This leads to recognize that the Obukhov length (as well as the Richardson number) cannot be used as the only parameter to model for instance the dependence of the TKE averaged for 30 min on $z/L$, because TKE may be influenced by low-frequency motions that only have a marginal effect on L and do not directly depend on stability, although they occur more frequently under stable conditions.

All the normalized second-order moments investigated herein exhibit a well-defined change in the range $0.1<R_B<1$, thus confirming literature results, although the absolute values can be site-dependent, and also time-averaging interval length dependent.

The vertical divergence of the second-order moments is shown to be a common feature of this boundary layer. In particular, the divergence of the vertical velocity variance is put into relation with the vertical velocity skewness, evidencing that in a relevant range of the parameter $I_w$, the flux–gradient relation is not satisfied. More generally, this observation is also supported by the behavior of the turbulent velocities for the mean velocity approaching zero.