Analysis of Different Gray Zone Treatments in WRF-LES Real Case Simulations

When conducting meso-micro scale coupled atmospheric simulations, it is crucial to ensure an adequate treatment of gray zone or terra incognita resolutions in which a large portion of the kinetic energy is naturally produced by the momentum balance equations in the model, while the remaining part still needs to be parameterized. In this work, we conduct three multiday, real case, full-physics atmospheric simulations that are fully coupled from the meso to the micro scale and in which the only difference is the treatment of boundary layer physics at the gray zone domain. One simulation uses a well-established 5 parameterization, another uses its scale-aware version previously modified to accommodate gray zone resolutions, and a final one uses no parameterization at all and assumes that the gray zone domain can be run in large-eddy simulation (LES) mode. The simulated fields are cross-compared, and further compared to measurements collected during the Prince Edward Island Wind Energy Experiment. Use of LES in the gray zone domain influences the flow fields in a manner that is robust to temporal averaging. The best predictions of vertical wind shear were found for the simulations in which the gray zone is parameterized, 10 and the inclusion of a micro scale nest run in LES mode within the gray zone domains increased the model errors by producing overly homogeneous flow fields. The parameterized simulations also produced better agreement in terms of kinetic energy spectra at the two innermost simulation domains. In the gray zone domain, the energy decays as f−3 throughout most of the spectral range considered. In the micro scale domain, the same is only seen in the low-frequency end of the gray zone spectral range. In the high-frequency end, the energy decay follows a f−1 slope. Outside the gray zone spectral range, the micro scale 15 simulated spectra follow the expected f−5/3 slope and produce good agreement with measurements.

problem sufficiently for computational tractability and physical understanding (e.g., Schalkwijk et al. (2015); Shin and Hong (2013)).In this work, an important distinction is made between "idealized" and "real case" simulations based on the data being used as initial and boundary conditions.In idealized simulations, the model is often driven towards specified values of wind speed and temperature which are defined with the intent of representing canonical atmospheric conditions (e.g., Churchfield et al. (2012)).In real case simulations, initial and boundary conditions are specified to match a past event of interest and often a reanalysis data set is used to drive the model.
Another important distinction relevant for this work is that between "full-physics" and "fluid dynamics" models.In addition to focusing on idealized conditions, LES studies for wind farm aerodynamics have historically used fluid dynamics solvers that do not include the full range of atmospheric physics and dynamics that impact flow within wind farms (e.g., cloud physics and radiative processes).Thus far, full-physics real case LES focusing on atmospheric phenomena at micro scales have typically employed models that run un-coupled from the meso and macro scales (e.g., Kunz et al. (2000)) and for a time period on the order of hours (e.g., Liu et al. (2012)).For the purpose of this work, we distinguish between macro, meso, and micro atmospheric scales based on the duration (T ) and spatial dimensions (L) of the most energetic processes within them as given in Table 1.
Table 1.Scales of atmospheric motions as defined for the present work based on duration and spatial dimension of most energetic processes within them.

Duration Spatial Dimension
Macro Scale 1 year < T < 1 day 10 4 km < L < 10 1 km Meso Scale 1 day < T < 1 hour 10 1 km < L < 10 0 km Micro Scale T < 1 hour L < 1 km Ongoing computational developments have only recently started to allow for full-physics, real case atmospheric simulations to be performed for wind energy applications.This can be achieved by coupling meso and micro scale models, and is a topic of active research for which two main approaches are being considered.One alternative is to perform meso scale simulations with an atmospheric model uncoupled from the micro scale and to subsequently use the meso scale flow to prescribe initial and boundary conditions to a secondary code that is then run on a smaller domain, a finer mesh, and in LES mode.Another alternative is to fully couple the meso and micro scales within a single code.As a consequence of these new developments, two important questions have emerged.Firstly, how best to prescribe initial and boundary conditions to the LES and whether to allow feedback between the scales when the coupling is performed within a single code.Secondly, how best to treat the transition from meso to micro scales within numerical models (Sanz Rodrigo et al., 2016).The concept of "gray zone" (GZ) resolutions or "terra incognita" was coined to describe the spatial scales O(10 2 − 10 3 ) m at which Navier Stokes is able to resolve a substantial fraction of the kinetic energy in the atmospheric boundary layer (ABL), while still needing to model the remaining part with physical parameterizations (Wyngaard, 2004).Under the first coupling approach, the GZ is avoided but turbulence generation methods must still be applied to the meso scale fields before they can be used as driving inflow to LES quantify differences between three simulations in which the GZ is treated differently by being run with a well-established ABLP, its scale-aware version, and no ABLP at all.The simulation output is compared to observational data collected during the Prince Edward Island Wind Energy Experiment (Barthelmie et al., 2016).This specific location was selected for our analysis because the terrain complexity and roughness changes found at the site warrant the use of LES in contrast to coarser, meso scale simulations which are unable to resolve the complex flow characteristics at the site.The analysis considers a 16-day period thus enabling an assessment of the simulations performance under a range of atmospheric conditions and several diurnal cycles.The methodology is discussed in Section 2. A cross-comparison of the model simulations is given in Section 3, and a comparison between the simulations and measurements in Section 4. A final discussion can be found in Section 5.

Boundary Layer Treatment at the Gray Zone
Parameterized phenomena in meso scale atmospheric models (such as WRF, when it is not run in LES mode) include radiative transfer, convection, and ABL physics.In ABLPs, a distinction is made between local vertical transport between adjacent grid levels and non-local vertical transport via strong updrafts that span the ABL depth (Frech and Mahrt, 1995).Non-local ABL schemes such as the well-established YSU (Hong et al., 2006) include a non-local term in addition to the local transport.The SH scheme (Shin and Hong, 2015) is an expansion of YSU in which the amount of SGS vertical transport that needs to be parameterized is regulated based on ∆ xy and on the strength of the non-local transport itself.
All simulations presented herein use the YSU ABLP for domains with a grid resolution coarser than the GZ resolution, and are run in LES mode for the innermost domain in which the grid resolution is finer than at the GZ.The WRF source code was modified to enable meso-micro scale coupling based on the potential temperature perturbation method of Muñoz-Esparza et al. (2014, 2015).Previous work has shown that this modification results in an improvement to simulations of wind speed and turbulence intensity under different meteorological regimes and terrain complexities (Montornès et al., 2016b, a).
Hereinafter, the three simulations conducted are referred to as YSU_LES, SH_LES, and LES_LES (Table 2) based on how they handle the ABL physics at the GZ and at the innermost domain.In the LES_LES simulation, the ABLP is switched off in the GZ domain and the SGS energy is modeled with a local closure in which a 1.5-order prognostic equation for turbulent kinetic energy (TKE) is used.This approach is recommended when the energy-containing eddies are substantially larger than the grid resolution (Shin and Hong, 2015) which is not necessarily the case at the GZ.Therefore the GZ treatment in our LES_LES simulation is not necessarily a physical choice, but can be justified by the existence of a nested domain which is appropriately run in LES mode.Because the differences between YSU and SH pertain to the strength of the vertical transport, we expect differences in vertical flow velocities and in vertical fluxes of heat and momentum between these two simulations.In idealized simulations, the system complexity is lower and substantial differences between YSU_LES and SH_LES are only expected during convective conditions (i.e.strong non-local transport) while under neutral and stable stratification, SH should default to its underlying YSU structure.On the other hand, the higher level of complexity in a real case simulation precludes the predetermination of the exact effects that SH will have on the results.This complexity further highlights the need for the present study in which these differences are investigated.Prior to the analysis, the only claim that can be made for the real case simulations conducted is that a distinction between YSU_LES and SH_LES is only expected during isolated periods throughout the simulation, when the differences between YSU and SH are noticeable.Instead, the largest difference is expected to be found between LES_LES and the other two simulations in which the ABL is parameterized.

Study Domain and Measurements
The simulations were conducted for a domain centered on the North Cape of Prince Edward Island for the period of the Prince

Simulations
We performed three simulations using the Advanced Research WRF core version 3.7.1 that differ only in their treatment of the GZ domain as described in Section 2.1.A total of five domains (with 70 x 70 grid points) were configured, centered at the meteorological mast location following a set of telescopic nests from 9 km to 111 m (Fig. 1b).Feedback between the scales is disabled so that the GZ fields can be analyzed without any influence from the turbulence that develops at the innermost domain which is run in LES mode across all simulations.Fifty vertical levels are defined throughout the domain, which are distributed every ∼ 10 m in the ABL and stretched in the free atmosphere.Initial and boundary conditions for all domains are taken from the Climate Forecast System Reanalysis (Saha et al., 2010) at a spatial resolution of ∼ 38 km.Topography for the highest resolution domain is taken from the Shuttle Radar Topography Mission (Jarvis et al., 2008) at a 90 m spatial resolution, and land use from GlobCover (Bontemps et al., 2011) at a ∼300 m spatial resolution.
The analysis of model output focuses on flow parameters in the ABL surface layer.The model was run for a 16-day period and time series were saved exclusively for the grid point closest to the meteorological mast location (Fig. 1a) for seven heights above ground between ∼ 10 m and ∼ 180 m.The model saves historic output with a temporal frequency of 4 Hz, which are immediately used to compute 10 minute-averaged diagnostics.

Variables and notation
All analyses are performed for a fixed location (i.e., the meteorological mast) and consider vertical profiles of 10-minute mean quantities.Hereinafter, the overbar ( ) is used to denote 10-minute averages when referring to (co)variances, but is omitted for all other variables to simplify the notation.We focus on flow variables of relevance to wind energy such as horizontal wind Differences in simulation output that arise from varying the model configuration are discussed in Section 3 without reference to the measurements.To avoid selecting one of the three cases as the reference benchmark, these differences are quantified by a deviation from each simulation to the ensemble average of the three simulations.The model performance relative to observations is quantified in a similar fashion by subtracting the modeled fields from the measured values.To facilitate a later comparison with measurements, the majority of the model cross-comparison section focuses on the measurement heights of 20, 40, and 60 m by linearly interpolating from the model levels to the target height.

Spectral Analysis
Spectral methods are used to compare the kinetic energy content of each data set in the frequency domain.Comparisons are presented between the three simulations, between the two innermost domains for each simulation, and between simulations and measurements.The power spectral density (PSD) values are computed for time series of 10-minute mean wind speeds using the Welch method.Previous work considering onshore and offshore sites has found that wind speed measurements collected at a frequency of 10 minutes are sufficient to accurately estimate the PSD in the atmospheric surface layer (Larsén et al., 2016).
Because the spectral analysis is based on relative comparisons between the data sets, it is not concerned with energy leakage at low frequencies that could be caused by non-stationarity in the time series, and a hanning window with a length of 1 day is applied when computing PSD values.With a temporal frequency of 10 minutes and a window length of 1 day, the spectra are limited to 1 day −1 < f < 20 min −1 which is sufficient when focusing on GZ phenomena.High frequency fluctuations are therefore not included in the spectra, but are instead contained in the analysis of turbulence intensity and TKE which are computed from the 4 Hz model output.

Cross-Comparison of Simulations
In this section the three simulations are cross-compared to quantify their sensitivity to different treatments of the ABL physics at the GZ scale.Section 3.1 focuses on horizontal and vertical wind speeds, Section 3.2 on vertical wind shear and turbulent kinetic energy, and Section 3.3 discusses the energy content of each simulation in the frequency domain.

Horizontal and Vertical Flow
To compare the flow predictions across the three simulations, an ensemble average was first computed and then differences between each simulation and this reference ensemble were obtained at each time step.The values given in Fig. 2 show the mean and standard deviation of these differences at each hour of the day, including 16 days of simulation and 6 10-minute values in each hour.On average, the LES_LES simulation produced higher values of horizontal wind speed U and of vertical velocities w than the parameterized simulations.For both wind speed and vertical flow, there is no clear tendency of the mean differences with time of day.Rather, a consistent difference is seen that decreases with height for U and that remains of the same magnitude for w across the three heights.As expected, differences between YSU_LES and SH_LES are intermittent and average out when a long time series is considered.The standard deviation of these differences is also substantially larger for LES_LES than for the other two simulations, indicating that running the GZ domain in LES mode produces more fluctuations 10 in the solution than activating a traditional or scale-aware ABLP.Values are for the GZ domain.
When considering all model levels for which output was saved (Fig. 3), it is clear that LES_LES produces higher wind speed estimates only in the lower portion of the profile, and that the opposite is true above ∼ 60 m.Mean differences across the simulations for vertical flow decrease with height, while the standard deviation of these differences increases for the results suggest that running the GZ domain in LES mode has a large impact on the results, even when multi-day averages are considered.Mean horizontal wind speeds at some heights differ by more than 1 m s −1 on average, with standard deviations of the same magnitude which can be critical when using meso-micro coupled simulations for wind energy applications.The same is true for wind direction, with standard deviation values (as computed from the 10-minute periods) indicating that onethird of periods exhibit a difference in wind direction greater than 30 • .These differences cannot be neglected when performing WRF-LES simulations for wind plant optimization via layout design or wake-redirection control strategies, which is a potential application of these simulations as a robust meso-micro coupling strategy is developed.A better understanding of the physical reasons behind these differences is investigated in the next sections.

Shear and Turbulence
The probability density function of 10-minute mean values of vertical wind shear across a nominal turbine rotor plane (extending between ∼ 30−110 m) also differs substantially between LES_LES and the simulations in which the ABL is parameterized.
At the GZ and micro scale resolutions YSU_LES and SH_LES are almost indistinguishable (Fig. 4).However, the difference between these simulations and LES_LES is pronounced especially at the GZ domain where the shear distribution derived from LES_LES has a much lower mean and standard deviation, with all values below 0.03 s −1 .This is consistent with the results   T KE SGS is combined with T KE ∆ which is of the same magnitude for the three simulations, large differences result in the total TKE budget between LES_LES and the other two simulations.High TKE at ∆ xy = 333 m in LES_LES contributes to explaining the low vertical wind shear values seen in this simulation (Fig. 4).

9
In the micro scale nest, differences in TKE between LES_LES and the other two simulations are lower in magnitude but still clearly noticeable.These differences are similar to those seen in the GZ domain, thus suggesting propagation down the model chain.In the micro scale domain, these differences end up modulating the budget of SGS, modeled TKE and leading to differences in T KE SGS across the different simulations.While SH_LES and YSU_LES present median T KE SGS values that increase with wind speed up to 15 m s −1 , LES_LES shows a local TKE maximum between 10 and 12 m s −1 coinciding with the peak of subgrid TKE production in the GZ domain when run in LES mode.These results suggest that the subgrid TKE production at the GZ domain enabled by running it in LES mode inhibits the subgrid TKE production at the micro scale nest.
In other words, the use of ABLPs at the GZ results in lower TKE at the GZ domain but in higher TKE at the innermost nest.

Energy Spectra
The results presented in this section are based on the spectral analysis methods described in Section 2.5.The PSD for wind speeds considering the entire 16-day time series are shown for the three simulations in Fig. 6.As expected, YSU_LES and SH_LES exhibit higher PSD values when compared to LES_LES because the ABLP (enabled in YSU_LES and SH_LES) acts to generate turbulence at length scales higher than the grid size.Note that this analysis is different from the one presented in Section 3.2 which encodes turbulence information frequencies up to 4 Hz.An important result from this spectral analysis is that this difference between LES_LES and the other two simulations is seen throughout the spectrum, and is not concentrated in the spectral tail as might be expected.
Fig. 6 also indicates that the energy decay at the GZ domain follows a f −3 slope instead of the f −5/3 expected for meso scale and boundary layer turbulence at frequencies larger than 1 day −1 (Larsén et al., 2016).This result suggests that the GZ domain is lacking energy for the three simulations, and that a better treatment of turbulence at this scale is needed as will be further discussed in Section 4.3.At the micro scale domain, the spectra are substantially more energetic at higher frequencies and exhibit the f −1 decay seen in Muñoz-Esparza et al. ( 2017) for these frequencies, also exclusively in the nested LES domain.
Because simulation data are unavailable at temporal resolutions higher than 10 minutes, an in-depth analysis of the energy content at higher frequencies is not possible and we cannot determine whether the expected spectral behavior is achieved for f > 20 min −1 .

Comparison of Simulations and Measurements
In this section the three simulations are compared to measurements from sonic anemometers to quantify how much improvement can be gained from having a micro scale domain nested within a multi-scale simulation framework, and to determine which GZ treatment produces the overall best results in terms of mean flow velocities and turbulence.Section 4.1 focuses on horizontal wind speed, Section 4.2 on vertical wind shear and turbulent kinetic energy, and Section 4.3 compares the energy content of each data set in the frequency domain.Thus far the simulations have been cross-compared without reference to the measurements to illustrate the sensitivity of flow parameters to different GZ treatments in each simulation.In this section, we evaluate the simulated flow fields relative to the observational data set described in Section 2.2.Whenever not explicitly specified, the comparison considers simulation data from the micro scale domain.There is a systematic pattern of wind speed error distribution as a function of wind direction, across the three simulations.

Horizontal Wind Speed
During the few cases of winds between the northeastern and southeastern sectors (wdir < 165 • , flow over the island) the model mostly under-estimated the wind speed magnitude.Conversely, an over-estimation is seen for westerly and southerly winds (165 • ≤ wdir < 315 • ) which make up the majority of the data set and represent on and offshore flow (Fig. 1).In terms of differences between simulations, let us consider only the sectors with higher measurement density (wdir > 165 • ).We can see that LES_LES performed slightly better for southerly winds and the other two simulations for westerly winds.This result may be related to the onshore footprint at lower measurement heights (e.g., z = 20 m) at this site as discussed in Barthelmie et al. (2016).This onshore footprint would result in the development of an internal boundary layer brought on by changes in roughness and orography, which may fail to be reproduced in a simulation where the ABL is not parameterized thus explaining the slightly higher mean error values for LES_LES (when compared to YSU_LES and SH_LES) for westerly flow.
In terms of turbulence intensity, a clear pattern can also be seen.Overall, the median wind speed error is higher for higher levels of turbulence intensity where wind speeds tend to be lower.Regarding the differences across simulations, the model errors are larger for LES_LES at the low turbulence intensity bins (T I < 0.15) and lower for LES_LES for the bins 0.15 ≤ T I < 0.25.Note that the majority of the data lie within TI<0.15, leading to a temporal average of wind speed errors that is slightly higher for LES_LES than for YSU_LES and SH_LES.

Shear and Turbulence
Fig. 8 shows the distribution of shear across the measurement layer (z = 20 m to 60 m) during the entire period for which measurements are available after quality control.An important result from this analysis is that the best model performance is seen for the GZ domain when the ABL is parameterized.The simulation run in LES mode largely underestimates the observed shear across this layer, which is also true for the three simulations at the innermost domain (∆ xy = 111 m).This result indicates that the energy produced by the ABLP at the GZ is important to the simulation of vertical variations of wind speed in the surface layer, and that running a simulation in LES mode at the GZ compromises the model results in the presence of high shear (e.g., under stable atmospheric conditions which were commonly observed during the Prince Edward Island Wind Energy Experiment).These results also suggest that adding an extra domain with a higher horizontal resolution to the simulation with one-way coupling between the meso and micro scales (as is done here, with feedback being disabled between the domains) does not necessarily improve the simulated flow fields as might be expected.In other words, an adequate treatment of coarser domains is more important to obtaining accurate simulations of vertical wind shear than increasing the resolution of the grid by nesting the domains down to what can be assumed as an LES scale.

Energy Spectra
Fig. 11 shows the spectral energy content in the 10-minute mean wind speeds at each measurement height, computed using the methods described in Section 2.5 and considering only the days in which the measurement data set was complete.Once again, the difference between YSU_LES and SH_LES can hardly be discerned.The spectral interval characterized by the frequencies f ∈ 4 × 10 −5 , 5 × 10 −4 Hz (time scales between ∼ 7 hours and ∼ 30 minutes) reveals a clear and consistent underestimation of the energy content by all of the simulations and at all three heights, clearly marking the terra incognita described in Section 1.The analysis for this period and site suggest that running the GZ domain with an ABLP produces better results in terms of variance of the streamwise velocity than configuring this domain in LES mode.This analysis complements the discussion presented in Section 3.3 by demonstrating that the three simulations perform well below and above GZ scales, where the energy decay follows a f −5/3 slope which was not evident from the model cross-comparison alone.While the observations show a f −5/3 relationship throughout the frequencies considered, the simulations present a f −1 slope at the GZ scales as previously discussed.
16 for observations (black), LES_LES (gray), YSU_LES (blue), and SH_LES (red).Top shows modeled fields at GZ domain, and bottom at innermost domain.Spectra consider all 10-minute mean time stamps in which measurements and simulations were available and an entire 24hour period was complete.Gray shading marks GZ spectral region.Top figures are for the GZ domain, and bottom figures for the micro-scale domain.

Conclusions
Due to ongoing advances in computational power, meso-micro scale coupling for wind energy applications has recently become a topic of active research.We are currently moving away from the need to choose between (i) idealized simulations at high spatial and temporal resolutions or (ii) real case simulations limited to the meso scale.In this study, we conduct multi-day, real case, full-physics simulations of the atmosphere down to a spatial resolution of 111 m, and consider different ways of treating 5 the terra incognita or GZ which characterizes the transition from meso to micro scales.At this transition, a large portion of the kinetic energy is naturally produced by the momentum balance equation in the model, while the remaining part still needs to be parameterized.We conducted three simulations that only differ in their treatment of the ABL in the GZ domain.One of them uses the traditional YSU parameterization, another uses its scale-aware version SH, and a final one uses no parameterization at all and runs the GZ domain in LES mode.We seek to quantify how different GZ treatments affect flow simulations at the GZ domain and at the micro scale domain nested within it, and to determine which simulation performs best relative to measurements.
We find that running the GZ domain in LES mode produces more fluctuations in the horizontal and vertical wind fields, and that it has a large impact on the simulated flow, with temporally averaged differences in the horizontal wind speed and direction as large as 1 m s −1 and 10 • between LES_LES and an ensemble mean of the three simulations.The distribution of wind speed errors showed clear trends with atmospheric conditions, but the tendency was uniform across the three simulations.
Differences in mean wind speed errors between the simulations are small, with SH_LES and YSU_LES performing slightly better for low turbulence intensity and for offshore flow than LES_LES.
The LES_LES simulation produced higher vertical velocities and TKE at the gray zone domain, which led to more homogeneous flows and lower levels of shear when compared to SH_LES and YSU_LES.The best shear predictions were those by the parameterized runs at the GZ domain, suggesting that the treatment of the ABL at GZ resolutions is more critical to accurate simulations of shear than increasing the spatial resolution with an extra nested domain.These results are likely to differ if feedback is enabled between the domains, which was not considered in the present study.
Differences between YSU_LES and SH_LES were found to be negligible when multi-day temporal averages are considered, and much lower in magnitude than the model error relative to measurements.However, differences between these two simulations are discernible over short intervals not only at the GZ domain but also at the inner nest, which uses data from the GZ domain as boundary conditions.These differences are most pronounced for kinematic heat flux, where SH_LES produces peaks of slightly lower magnitude than YSU_LES during daytime.None of the simulations capture the diurnal cycle of the heat flux near the surface, indicating that this systematic error might be related to surface layer physics instead of ABL physics.
The spectral analysis indicates that including a micro scale domain within the model chain leads to a substantial recovery in the decaying tails of the spectral energy content, and clearly reveals the existence of a GZ or terra incognita in the frequency range between 7 hours −1 and 30 minutes −1 .The PSD computed for wind speeds simulated in the GZ domain follow a f −3 slope throughout most of the spectral range considered.In the micro scale domain, the same is only seen in the low-frequency end of the GZ spectral range.In the high-frequency end, the energy decay follows a f −1 slope.Outside the GZ spectral range, the micro scale simulated spectra follow the expected f −5/3 slope and produce good agreement with measurements.Overall, the parameterized simulations SH_LES and YSU_LES produce streamwise velocity variance values that more closely match the observations at both resolutions.
Overall, we found that for multi-day simulations the present question is whether to run the GZ in LES mode or with a ABLP at all, and not whether to consider scale-aware parameterizations which are still in early development phases.For the period considered, very small differences were seen between YSU and its scale-aware version, SH.For shorter time periods, the differences between SH_LES and YSU_LES can be substantial and the GZ treatment should be considered more carefully.With the current data set, we cannot yet determine whether SH_LES could result in a substantial improvement of the simulated fields under non-idealized simulations.We generally found that parameterizing the GZ domain results in more accurate predictions of shear, TKE, and spectral energy content.Further work is needed to generalize the development of scale-aware parameter- Site satellite imagery and elevation.(b) WRF simulation domains.

Figure 1 .
Figure 1.Measurement site terrain elevation, satellite imagery from Environmental Systems Research Institute, and location of meteorological mast (triangle) (a).Simulation domains from ∆xy = 9 km to ∆xy = 111 m (b).
Wind Energ.Sci.Discuss., https://doi.org/10.5194/wes-2017-61Manuscript under review for journal Wind Energ.Sci. Discussion started: 15 January 2018 c Author(s) 2018.CC BY 4.0 License.speed, wind direction, vertical wind shear, and turbulence intensity.Turbulence intensity is only calculated for time stamps in which the horizontal wind speed is ≥ 3 m s −1 which are most relevant to wind energy applications.While the focus of this work is on variables of interest to wind energy, we also choose to include vertical flow velocities and vertical fluxes of heat and momentum which are crucial to characterizing turbulence and atmospheric stability.Perturbations are calculated by subtracting the 10-minute mean from the 4 Hz instantaneous values as the simulation output is produced, and then used to obtain turbulence diagnostics at a 10-minute frequency.The turbulence diagnostics computed are kinematic heat flux and turbulent kinetic energy T KE = T KE ∆ + T KE SGS .T KE ∆ represents the amount of TKE that is explicitly resolved by the model at a grid resolution ∆ xy including the ABLP budget, and T KE SGS is the amount that is produced by the turbulence closure model when no ABLP is used and the simulation is run in LES mode.The SGS partition is diagnosed by the model directly, and the resolved partition is calculated as Wind Energ.Sci.Discuss., https://doi.org/10.5194/wes-2017-61Manuscript under review for journal Wind Energ.Sci. Discussion started: 15 January 2018 c Author(s) 2018.CC BY 4.0 License.

Figure 2 .
Figure 2. Differences between each simulation and their ensemble average for horizontal wind speed (top) and vertical flow velocity (bottom) for the heights 20 m (left), 40 m (middle) and 60 m (right).Error bars give the mean ± standard deviation of differences at each hour of day.
Wind Energ.Sci.Discuss., https://doi.org/10.5194/wes-2017-61Manuscript under review for journal Wind Energ.Sci. Discussion started: 15 January 2018 c Author(s) 2018.CC BY 4.0 License.three simulations.Finally, a consistent difference is seen for wind direction, where LES_LES predicts flow that is veered (i.e., rotated in the clockwise direction) by ∼ 10 • relative to YSU_LES and SH_LES throughout the vertical layer considered.These

Figure 3 .
Figure 3. Vertical profile of differences between each simulation and the ensemble average for horizontal wind speed (left), vertical flow velocity (middle) and wind direction (right) for the GZ domain.Errorbars give mean ± standard deviation considering the entire time series.
Wind Energ.Sci.Discuss., https://doi.org/10.5194/wes-2017-61Manuscript under review for journal Wind Energ.Sci. Discussion started: 15 January 2018 c Author(s) 2018.CC BY 4.0 License.discussed in Section 3.1 which indicate higher vertical velocities for the LES_LES simulation, resulting in more homogeneous flow and less vertical wind shear in the surface layer.This result is also evident in the micro scale domain, where the mean and standard deviation of shear are still lower for LES_LES than for the other two simulations, evidence of the importance of the ABL treatment at the GZ when conducting simulations with multiple nested domains.

Figure 4 .
Figure 4. Curve-fitted histograms for vertical wind shear between ∼ 30 − 110 m, considering the 16 days of simulation.Values shown for GZ and micro scale domains.

Fig. 5
Fig.5shows the TKE budget from the three simulations for the highest measurement height (z = 60 m) and allows us to examine how the treatment of the ABL physics at the GZ resolution affects the total TKE budget within the two innermost domains.At the GZ domain, differences between LES_LES and the two other simulations are obvious, and result from the fact that LES_LES produces large amounts of subgrid TKE for wind speeds > 8 m s −1 , while the other two produce none.When

10Figure 5 .
Figure 5. Median (markers) TKE at 60 m [m 2 s −2 ] as a function of wind speed [m s −1 ] discretized in 1 m s −1 bins in the GZ (top) and micro scale (bottom) domains.Solid lines delimiting shaded area give 25 th and 75 th percentile values.

Figure 6 .
Figure 6.Power spectral density as a function of frequency for the 16-day time series of 10-minute mean horizontal wind speeds at z = 60 m for the three simulations at the GZ and micro scale domains.

Fig. 7
Fig.7shows the distribution and magnitude of horizontal wind speed errors as a function of different atmospheric conditions in terms of wind speed, wind direction, and turbulence intensity.An obvious and important result from this analysis is that the YSU_LES and SH_LES error distribution is almost indistinguishable regardless of the atmospheric conditions considered, indicating that the differences between these two simulations discussed in Section 3 are lower in magnitude than the model error.The data distribution histogram above the leftmost subplot shows that most of the observed wind speeds are between 2 and 14 m s −1 , and that wind speed errors exhibit a decreasing tendency with increasing wind speed.Although time stamps in which the meteorological mast was subjected to a wind turbine wake were filtered out, the model still over-estimates the wind

Figure 7 .
Figure 7. Wind speed errors [m s −1 ] at z = 60 m averaged over different measured atmospheric conditions characterized at the same height: wind speed (left), wind direction (middle), and turbulence intensity (right).Markers represent the median error, and whiskers give 25 th and 75 th percentile values for each bin.Bars indicate how much of the total data considered belong in each bin.

Figure 8 .
Figure 8. Histogram for vertical wind shear ∆U/∆z [s −1 ] between 20 and 60 m, for all 10-minute periods for which coincident model output and observations are available.

Table 2 .
Treatment of atmospheric boundary layer physics at each simulation domain for the three simulations conducted in the experiment.