Analysis of turbulent kinetic energy decay power law in atmospheric boundary layer models

The transitional periods that take place in the atmospheric boundary layer are challenging to model due to their non-stationary nature. Large-eddy simulation (LES) models have proven to be sufficiently accurate for modeling the evening transitions, and now a possibility to adapt turbulent kinetic energy (TKE) closure models for calculations in single-column models appears to be very plausible. In this study, evening transition modeling is analysed with emphasis on the pattern of a TKE decay which follows the power law. E(t) α t −αThe effects of different parameters on the results of the simulations are explored, along with the geostrophic wind effect on the model. It is shown that the model presented here behaves in a manner similar to that of a LES model, thus showing that the above adaptation is possible and worth being investigated further.


Introduction
The atmospheric boundary layer (ABL) is the lowest part of the atmosphere, where most of the transport processes take place. The ABL is directly influenced by the contact with the Earth's surface and responds to the surface forcing, such as frictional drag, heat transfer, pollutant emission, evaporation and transpiration, and terrain induced flow modifications. One of the key characteristics of the boundary layer is the diurnal variation taking place within it, for example, temperature variations, which is insignificant in the free atmosphere. Another one is turbulence, one of the most important transport processes, whose dynamics is affected by the heat flux from the ground. These processes are interdependent and strongly influence one another [1].
The diurnal cycle consists of the boundary layer changing its state between the stably stratified boundary layer (SBL) and the convective boundary layer (CBL). Two transitional periods are distinguished: the morning transition (from SBL to CBL) and the evening transition (from CBL to SBL). The morning transition is initiated when the surface heat flux becomes positive and a shallow entraining mixed layer grows into surface inversion. The evening transition starts when the surface heat flux becomes negative, and it consists of the decay of convective turbulence through dissipation. Therefore, when it comes to singling out transitional periods, one can use the following criteria taken from a BLLAST study [2]: The morning transition starts when the heat flux starts to grow and ends when it changes its sign (from negative to positive). The evening transition starts when the heat flux changes its sign (from positive to negative) and ends when the SBL is well-established. IOP Conf. Series: Earth and Environmental Science 611 (2020) 012014 IOP Publishing doi:10.1088/1755-1315/611/1/012014 2 Now that the temporal and spatial resolution in weather forecast and climate models has improved greatly throughout the recent decades, studying the dynamics that are influenced by the diurnal cycle, especially those that are taking place within transitional periods, is the next logical step towards further advancement of the models.
The evening transition is characterized by the decay of convective turbulence. The turbulence decay rate obeys the power law ,where is the normalized TKE and is the normalized time. The parameter can be found theoretically or empirically.
The value of is generally predicted to be for the decay of homogenous isotropic turbulence, which was obtained by using the hypothesis of the Loitsyansky invariant [3], supported by laboratory [4] and numerical [5] experiments. The assumption of a constant integral turbulence length scale, which reflects the strong turbulent kinetic energy concentration on the largest scales, results in for the decay of boundary layer turbulence [6,7], with the Reynolds number decreasing with time. These values are still used in modern climate and weather models, even though the assumption of homogenous isotropic turbulence results in simulations that differ marginally from the dynamics in the ABL. Meanwhile, there have been a number of large-eddy simulation (LES) experiments, with different LES models and setups, performed to study the evening decay and figure out its dynamics; therefore, it could be beneficial to redirect attention to the results obtained through those and see how they can be applied in the existing models.
Nieuwstadt & Brost pioneered the idea of studying the evening decay through a LES experiment [8]. First the CBL is formed through constant heating from the surface, with its convective velocity scale varying for different runs. Upon forming the CBL, the surface kinematic heat flux is changed to 0, simulating the period of evening transition, when the heat flux ceases to be positive. All experiments showed a tendency towards the same decay rate of . A significant advantage of the LES experiments over laboratory studies has been noted, in particular, in their ability to simulate some conditions found in the atmosphere, such as the influence of stratification.
Sorbian utilized this type of experiment in his study, comparing the experiment proposed in [8] with the one where the heat flux gradually decreases over time [9]. This has shown that the process of decay is governed by the relation of the external time scale to the convective time scale. In his experiments for an abrupt change of the heat flux to 0 and for a gradual change. Beare et al. extended this experiment to full transition, rather than only the part before sunset [10]. The grid used in those LES experiments was significantly finer than those in the aforementioned studies (a vertical grid length of 5-10 m against 30-50 m). This study confirmed the possibility of simulating accurate evening transitions in LES models, and also presented a problem of correct representation of the ageostrophic wind in the model due to its negative impact on the simulation result when misrepresented.
Pino et al. explored the influence of wind shear by adding geostrophic wind in his LES simulations ( m s -1 and setups) performed on a relatively fine grid with a vertical grid length of 16 m [11]. Those experiments were similar to the one in [8]: the surface was heated for 2 hours with the heat flux set at , then abruptly changed to , and then 2 more hours were simulated. As expected, the rate of the TKE decay was noticeably slower with the inclusion of geostrophic wind, with for these experiments, and for experiments with no influence of wind shear.
Some of the latest studies, in particular, the analysis of Lindenberg Inhomogeneous Terrain-Fluxes between Atmosphere and Surface (LITFASS-2003) observational data [12] and simulations of Cooperative Atmosphere-Surface Exchange Study-1999 (CASES-99) data [13], suggest much faster decay with . In case of [12], such rapid growth of the TKE decay rate is explained by the change of sign for buoyant production, where it no longer slows down the decay process, but contributes to it. Paper [13] points out the influence of stratification as the main factor for the increased decay rate. Paper [14] found through LES experiments that is between 1.6 and 1.7 in the main stage of the evening transition and the decay is slower at the beginning and end of the transition period.  [15,16], the RANS models are in dire need of improvement before they can be effectively used in the same way as the LES models are.

Experiment setup
Here we use a single-column standard two-equation -model implementation (which contains prognostic equations for the TKE and its dissipation rate). The standard -model is widely used in weather forecast and climate large-scale models and, while imperfect, does a better job at modelling the ABL than first-order models, such as the one used, for example, in an INM RAS climate model [17]. Thus, the idea of adapting two-equation models to parameterize turbulence mixing in large-scale models appears plausible and, as shown in recent studies (e.g. [18,19]), to be worth investigating.
The wind velocity and potential temperature (in this particular study, hereafter "temperature" means potential temperature) equations are averaged by the horizontal coordinates: is the production of the TKE by shear, is the buoyancy term, , , , , (for stable stratification) and (for unstable stratification) are the model constants [19]. All experiments were conducted on a one-dimensional uniform grid of 256 cells, with a vertical domain of 2314 m (thus, the vertical grid length was 9 m). The rest of the setup follows that of the CBL experiment setup in [7].
The scenario of the experiment follows that of [14]. During the 12-hour run, for the first 6 hours the surface kinematic heat flux is set to be K m s -1 to build up a CBL, then it changes abruptly to (the neutral boundary layer), K m s -1 (the SBL) and K m s -1 (the strong SBL) in runs 1, 2, and 3, respectively. The free atmosphere temperature gradient was set at K m -1 (simulation of the development of a shear-free CBL). To analyze the influence of the geostrophic wind, all 3 runs were simulated with and without it. The geostrophic wind was weakened with setting the parameters at m s -1 and m s -1 , and then introduced with the parameters m s -1 and m s -1 . The aerodynamic roughness was set at m.

Results
As compared to the original LES experiment [14], both similarities and differences can be observed in the results of the RANS series of the experiments (Figure 2).   The similarities are best observed in the results obtained for Run 1, where the heat flux is changed to . The value of in the experiment without geostrophic wind closely resembles that obtained in the aforementioned LES experiments (Figure 2a). The same could be said for the   (Figure 2b) [11]. Although there is no particularly clear following of the power law, slowing down of the decay rate for Run 1 is apparent.
As for the differences, they show up first and foremost in the dynamics of the TKE decay in two ways. First, in the LES experiments there was a slight delay before the main stage of the decay commences (where it is the fastest); at the very beginning of the transition, for a short but noticeable period of time, the decay is much slower before it picks up speed. In our case the delay is also present, however, it is barely noticeable. Second, by the end of the 3rd hour of the experiment the decay rate gradually becomes weaker, with the value of approaching zero. This does not take place in our RANS experiment: the decay rate appears to be relatively uniform after it is established.
Another difference is present in the decay patterns for Runs 2 and 3. In the LES experiments with no geostrophic wind they follow the pattern of Run 1 very closely, while in these RANS experiments they clearly diverge, with the decay rate being noticeably faster than for Run 1 (Figure 2a).
Changes in the vertical distribution of the temperature and the TKE for the experiments with geostrophic wind are shown in Figure 3. A slight cooling of the layer near the surface can be observed closer to the end of the transition. The height of this "cool" layer is larger in Run 2, but the temperature of the air is lower in Run 3. In the absence of geostrophic wind the results for all 3 runs are similar to Run 1 in case of the temperature and to Run 3 in case of the TKE.
This calls for looking into how different closure parameters influence the outcome of the experiment. One of such parameters, the turbulent Schmidt number , is of a particular interest due to its influence on how fast the length scale reaches the state of equilibrium [19]. We can see that changing the value of has an influence on the overall decay rate, but most noticeably on the end of the evening transition and how it reaches the supposed stationary state (Fig. 4a). We can see that the larger is the value of (weaker diffusion), the slower is the decay rate. Upon the introduction of geostrophic wind (Figure 4b) the speed of decay inevitably reaches the state of equilibrium, and the larger the value of , the sooner it happens.
Changes in the value of ( Figure 5b) result in a noteworthy increase of the aforementioned delay before the main part of the decay. The delay also increases with the change of the (Figure 5a), but not as prominently. . For , a slight change can be observed in the decay rate, however, it barely seems to pose any significance to the experiment.
Another variable of interest in this study not shown here was the subsidence rate . In the original setup, its value is at 0 through the entire experiment. Then the value of was set to change to a non-zero value (in the range from -0.0001 to -0.01) at hour 1 of the experiment. It turns out that, in fact, the difference between the presence and absence of the subsidence rate in the model, as well as between different values it can take, is practically negligible, regardless of geostrophic wind.
The TKE balance at certain points in time can be observed in Figure 6. Here the vertical profiles for the TKE balance and its components are shown for time points during the first hour, 30 minutes apart. It can be seen that the strongest TKE dynamics is during the first half an hour of the evening decay simulationat 6:15, there is a lot of processes taking place, and they are actually uneven all the way to the top. By 6:45 the TKE balance is much weaker, nearing 0 along the entire ABL. Thus, it is fair to conclude that the majority of the decay takes place within the first hour. The vertical profiles for Runs 2 and 3 are similar to the one for Run 1 shown here.

Conclusions
In this study, a number of single-column experiments were conducted in order to understand the influence of the above-chosen parametrization on the experiment. It has been shown that there is a difference between LES and RANS experiment results when the setups are nearly identical, but also a number of similarities in the experiment results, such as the case of the influence of geostrophic wind closely resembling results from [11]. There are parameters which have little to no influence on the dynamics (such as and ), and there are parameters that change the dynamics of the experiment when tuned properly (such as , and ). This opens the possibility of further simulation improvement by choosing correctly the parameters.
The choice of a proper parametrization of the transitional periods is known to be a particularly challenging task, since the existing models do not represent strong stratification and non-stationary processes (which include the transitional periods) reasonably well. It is worth investigating further how and why the different parameters influence the model's behavior in RANS experiments, compared to LES ones, given that there are several ways in which their results are already similar, and also some evidence of how they could be further improved.