This paper presents a case study from an intensive observing period (IOP) during the Plains Elevated Convection at Night (PECAN) field experiment that was focused on a bore generated by nocturnal convection. Observations from PECAN IOP 25 on 11 July 2015 are used to evaluate the performance of high-resolution Weather Research and Forecasting Model forecasts, initialized using the Gridpoint Statistical Interpolation (GSI)-based ensemble Kalman filter. The focus is on understanding model errors and sensitivities in order to guide forecast improvements for bores associated with nocturnal convection. Model simulations of the bore amplitude are compared against eight retrieved vertical cross sections through the bore during the IOP. Sensitivities of forecasts to microphysics and planetary boundary layer (PBL) parameterizations are also investigated. Forecasts initialized before the bore pulls away from the convection show a more realistic bore than forecasts initialized later from analyses of the bore itself, in part due to the smoothing of the existing bore in the ensemble mean. Experiments show that the different microphysics schemes impact the quality of the simulations with unrealistically weak cold pools and bores with the Thompson and Morrison microphysics schemes, cold pools too strong with the WDM6 and more accurate with the WSM6 schemes. Most PBL schemes produced a realistic bore response to the cold pool, with the exception of the Mellor–Yamada–Nakanishi–Niino (MYNN) scheme, which creates too much turbulent mixing atop the bore. A new method of objectively estimating the depth of the near-surface stable layer corresponding to a simple two-layer model is also introduced, and the impacts of turbulent mixing on this estimate are discussed.
Thunderstorm frequency over the Great Plains during the summer occurs near midnight with a maximum in rainfall later in the night (Wallace 1975), yet nocturnal convection is predicted particularly poorly in NWP models (e.g., Davis et al. 2003; Surcel et al. 2010). Atmospheric bores play a critical role in both the initiation and maintenance of such nocturnal convection in the Great Plains (e.g., Koch et al. 1991; Karyampudi et al. 1995; Knupp 2006; Parker 2008; Coleman and Knupp 2011; Parsons et al. 2018). A bore is a quasi-permanent increase in the depth of a stable layer, often with superimposed undulations. Bores are generated when a stable air mass flows toward an obstacle relative to the obstacle’s movement, such as a convectively generated density current moving through the nocturnal boundary layer, resulting in blocked or partially blocked flow (Rottman and Simpson 1989). Bores are a characteristic feature of the nocturnal convective environment in the Great Plains (Haghi et al. 2017) due to the frequent occurrence of blocked or partially blocked flows and the presence of nocturnal low-level jets, which act to duct the bore energy. The common occurrence of bores and the relatively poor performance of NWP models in representing these nocturnal systems motivate our aim to investigate and quantify the uncertainties and errors in bore simulations of a well-observed case study. The purpose of this effort is to guide future work to improve NWP forecasts of bores and their impacts on nocturnal convection.
An accurate bore simulation is dependent on an accurate representation of both the evolution of the nocturnal boundary layer and low-level jet and a correct forecast of the formation, location, temperature, and depth of the density current. These features can be separated into three key ingredients that are necessary for successful NWP forecasts of bores. First, the mesoscale environment in which the bore occurs must be accurately described. Some important parameters of the environment that must be predicted include the Froude number of the undisturbed fluid,1 the depth of the obstacle, and the depth of the undisturbed stable fluid (Baines and Davies 1980; Baines 1984; Rottman and Simpson 1989). These parameters can be combined to predict the bore amplitude, as reviewed in section 2 below. Once a bore is generated, a criterion for its continued maintenance is having vertical profiles of buoyancy and wind shear that are conducive for gravity wave ducting (Crook 1986; Koch and Clark 1999). The second main ingredient needed for successful NWP forecasts of bores is an accurate forecast of the parent convection that generates the density current, which in turn triggers the bore. This ingredient is satisfied in this study by initializing the simulations from analyses that include the parent convection. This is done using data assimilation that incorporates radar, surface, and upper-air observations. However, it will also be shown that the model physics play an important role in the bore simulations through determining the characteristics of the density current. The third key ingredient is the requirement to have sufficiently small model grid spacing in order to resolve the bore. Since the horizontal wavelength of bores is typically around 10 km, a horizontal grid spacing of no more than ~1–2 km is needed to resolve such features (Skamarock 2004). This was confirmed by the Koch et al. (2008a) study, which demonstrated successful bore simulations at 2-km grid spacing, although some details of the bore that were resolved with 700-m grid spacing were not fully resolved with 2-km grid spacing. Koch et al. (2008a) used the simulations to better understand the observed bore processes and turbulent mixing in particular. Simulations with 1-km grid spacing have also resolved bores in the studies of Hartung et al. (2010) and Blake et al. (2017). Hartung et al. (2010) used the simulation to connect the observed bore dissipation to the weakening of the wave duct in the near-surface stable layer. Blake et al. (2017) used the simulation to investigate the role of a bore in the maintenance of a nocturnal MCS.
The basic theory of bores as hydraulic jumps in the depth of the stable layer, maintained by the ducting of vertically propagating wave energy, has been confirmed by numerous laboratory experiments and by analytically using a simple two-layer model of the atmosphere (e.g., Maxworthy 1980; Simpson 1982; Smith et al. 1982; Crook and Miller 1985; Rottman and Simpson 1989). However, the continuous thermodynamic stratification and vertical wind shear of the real atmosphere can affect density currents and bores differently than the two-layer models of early studies (Lindzen and Tung 1976; Liu and Moncrieff 2000). In particular, continuous stratification allows for vertical propagation of wave energy, such that gravity wave features are not perfectly trapped by a duct (Lindzen and Tung 1976). In an atmosphere that is stably stratified above the nocturnal boundary layer, the characteristics of the low-level jet determine the depth of the wave duct. A bore is therefore not simply a phenomenon of the nocturnal boundary layer (e.g., Koch and Clark 1999; Koch et al. 2008b; Haghi et al. 2017; Toms et al. 2017). Furthermore, the thermodynamic stability and vertical wind shear profiles vary throughout the night as the boundary layer cools and deepens (e.g., LeMone et al. 2014) and the low-level jet evolves (e.g., Shapiro et al. 2016). Therefore, the feasibility of using explicit model simulations, rather than theory-based parameters to predict the presence and characteristics of a bore, is an important aspect of the present study. Furthermore, objectively calculating a bore height from observations can be nontrivial, and a new method of doing so is presented herein.
Improving nocturnal convection forecasts was one of the motivations for the Plains Elevated Convection at Night (PECAN; Geerts et al. 2017) field project. Improving the prediction of bores themselves, in addition to their impacts on subsequent convection, is also a specific focus of the PECAN goals. Following the PECAN focus on the understanding and prediction of nocturnal convection, the present study aims to quantify the ability to predict the amplitude of a bore observed during a PECAN intensive observing period (IOP) using the multiscale Gridpoint Statistical Interpolation (GSI)-based data assimilation and ensemble forecasting system (Johnson et al. 2015; Wang and Wang 2017). Correctly predicting the bore amplitude is a prerequisite for predicting how the bore will modify its environment for subsequent convection. Correctly predicting the bore speed is also important for predicting exactly where the bore will be at a given time.
Given the small scale of bore features, it is hypothesized that the resolution of the data assimilation system, not just the forecast model, will play a role in the quality of bore predictions at short lead times (i.e., < ~6 h). It will be shown that the 1-km-resolution ensemble Kalman filter (EnKF) used in this study, described in section 3, is sufficient to initialize both the parent convection and the mesoscale environment, satisfying the first two ingredients described above. However, it is not clear if even the finest (1 km) resolution considered in this study will be sufficient to analyze the bore itself. Therefore, simulations initialized from analyses already containing the bore will first be compared to simulations initialized from analyses that contain the parent convection of the bore as a result of the radar data assimilation, relying on the model dynamics to generate the bore from that convection during the forecast period. The evolution of the density current, which strongly influences bore amplitude, is also expected to be sensitive to the microphysical parameterization of the model (e.g., Dawson et al. 2010). Furthermore, the bore speed and morphology are expected to be sensitive to the parameterized mixing processes (Christie 1989; Koch et al. 2008a). An investigation of the sensitivity of the bore predictions to model parameterizations of microphysics and the planetary boundary layer is therefore also undertaken.
The structure of this paper is as follows. Section 2 overviews the two-layer theory of bores from past studies, and section 3 describes the model and data assimilation system configuration. Section 4 presents an overview of the case selected for study and the observations obtained during the PECAN IOP. Model results are presented in section 5, while section 6 contains a summary and conclusions.
2. Two-layer bore theory
Rottman and Simpson (1989) review the hydraulic theory describing bores resulting from a density current impinging on a low-level temperature inversion and then conduct laboratory experiments in a simple two-layer fluid. In a fluid comprising two layers of different density and no mixing, a reduced gravity can be defined as , where and are the densities of the denser and less dense layers, respectively, and g is the gravitational constant. Hydraulic theory then predicts that the existence of the blocked or partially blocked flow necessary for bore formation and the height of the bore are determined by two nondimensional parameters, the Froude number F0 and a nondimensional height H, defined as
In Eqs. (1) and (2), U is the component of the wind speed that is normal relative to the density current. The wind speed is considered beneath the inversion in the ambient environment. In these equations, h0 is the stable layer depth in the ambient environment, and d0 is the depth of the density current. In reality, the atmosphere is a continuously stratified fluid, rather than a simple two-layer fluid. However, numerous studies of bores generated by convection over the Great Plains in the presence of nocturnal stable boundary layers have shown that the two-layer theory provides a relatively accurate representation of the bore characteristics (Koch et al. 2008a,b; Tanamachi et al. 2008; Marsham et al. 2011; Haghi et al. 2017; Blake et al. 2017; Toms et al. 2017).
Rottman and Simpson (1989) derive the bore height for a partially blocked flow (their Fig. 2). Once the two parameters H and F0 are known, the bore amplitude h1/h0 can be determined by iteratively solving the following equations [Rottman and Simpson 1989; their Eqs. (2.4)–(2.6)]:
Here, h1 is the depth of the stable layer after bore passage, C is the bore speed, g′ is reduced gravity, and u1 is the bore-relative flow (different from the density current relative flow U). It has been shown (Klemp et al. 1997; Kingsmill and Crook 2003) that the following equation improves upon the derivation for C in Eq. (3) and Rottman and Simpson (1989) by including the effects of energy dissipation by undulations and turbulence and the effects of the finite depth of the troposphere:
where S = h1/h0 is the bore amplitude. Equation (6) shows that, in theory, knowing the bore amplitude together with g′ and h0 is equivalent to knowing the speed of the bore. In the atmosphere, this relationship may not hold because of the effects of the ambient mean wind on the bore speed, which changes in space, height, and time, especially at night when the low-level jet is present. However, for the purpose of this study, bore speed is sufficiently dependent on bore amplitude so that we focus hereafter only on bore amplitude. Bore amplitude is also the more important characteristic for determining how the bore will affect the nocturnal environment through mean lifting and turbulent mixing.
Once it is determined that the hydraulic response of the environment to the density current is to generate a bore of a certain amplitude, the ability of the environment to duct, or trap, the gravity wave energy at low levels determines whether the bore response will be maintained or will be quickly radiated away as internal gravity waves. As originally derived by Scorer (1949), the presence of a ducting layer can be diagnosed using the Scorer parameter l2, defined as
where is the moist Brunt–Väisälä frequency, U is the horizontal wind speed normal to the bore, Cb is the speed of bore propagation, and k and m are the horizontal and vertical wavenumbers, respectively, of gravity waves that may or may not be trapped. In Eq. (7), k is often neglected due to the relatively small horizontal wavenumber of the waves of interest with wavelengths ~10 km (Hartung et al. 2010). Thus, a negative (or small, relative to k2) Scorer parameter l2 is interpreted as a layer that cannot support vertical propagation of waves with a real vertical wavelength. Such a layer thus ducts the wave energy beneath it, provided that the layer beneath it has a positive Scorer parameter and can therefore support the wave. Physically, this occurs when the environment is only weakly stable and/or the wind speed curvature (in the bore-normal direction) is very large. In the Great Plains, the Scorer parameter calculation is typically dominated by the wind curvature term (Koch and Clark 1999; Toms et al. 2017).
Toms et al. (2017) compared the observed bore speed and strength to the values predicted by the above two-layer model equations based on surface and upper-air soundings during bore passage. The observed values were in agreement with expectations from this two-layer hydraulic theory to within ~15%. Such deviations may be contributed to by the neglecting of wind shear and turbulence in development of the theory in Toms et al. (2017), as well as by stable stratification aloft, in contrast to the neutral stratification assumed by theory (Crook 1986). These deviations are partially accounted for in the theory described above since Eq. (6) accounts for energy dissipation by turbulence, but not the effects of wind shear.
Since an advantage of NWP at convection-permitting resolution is the ability to explicitly resolve bore characteristics and, therefore, their impact on subsequent convection, the present study focuses on the fidelity of the explicitly simulated bore amplitude in the numerical model to the corresponding observations described below in section 4. However, since this two-layer model has been shown to describe bore characteristics in the continuously stratified atmosphere with reasonable accuracy, it also provides a conceptual model for understanding sensitivities of the explicitly simulated bore characteristics in our numerical experiments.
3. Model and DA configuration
The experiments in this study use the multiscale GSI-based EnKF and ensemble forecast system, where both the in situ and convective-scale radar data are assimilated (Johnson et al. 2015; Johnson and Wang 2017; Johnson et al. 2017). The forecast model is the Advanced Research version of the Weather Research and Forecasting (WRF-ARW; Skamarock et al. 2005) Model version 3.6.1. The outermost model domain has 12-km grid spacing (Fig. 1), with a second domain with 4-km grid spacing (black box in Fig. 1) and an innermost domain with 1-km grid spacing (red box in Fig. 1). The data assimilation ensemble physics configuration (Table 1) is the same as used in Johnson et al. (2017), except that the Mellor–Yamada–Nakanishi–Niino (MYNN) (Nakanishi and Niino 2009) planetary boundary layer (PBL) parameterization is used instead of quasi-normal scale elimination (QNSE; Sukoriansky et al. 2005). Initial tests with this case showed better performance when making this change, in terms of providing the best preconvective and prebore environment.
In this study, assimilation of surface and upper-air observations from the operational National Centers for Environmental Prediction (NCEP) data stream on the 12-km domain is conducted every 3 h from 0300 UTC 10 July through 0000 UTC 11 July. NEXRAD observations, together with the North American Mesoscale Forecast System (NAM) model data assimilation system (NDAS) observations, are then assimilated every 10 min on the 1-km domain from 0000 through 0600 UTC 11 July. The PECAN IOP observations are not assimilated for these experiments. Instead, they are left out of the data assimilation system to be used as independent validation data for the purpose of understanding current model performance and sensitivities for bore prediction. The impact of assimilating PECAN observations will be included in another forthcoming paper. Deterministic forecasts, summarized in Table 2, are then initialized from the ensemble mean analysis at different times during this period and using different physics configurations for the forecast model.
The preprocessing of the NEXRAD observations is different than our previous work, necessitated by the nocturnal convection that is directly assimilated in this study. After performing the same preprocessing as described in Johnson et al. (2015) with Warning Decision Support System–Integrated Information (WDSS-II; Lakshmanan et al. 2007a,b; 2010), large areas of low-level nonmeteorological echoes remained. These echoes were likely biological (i.e., birds or insects; Martin and Shapiro 2007), had reflectivity values up to 20–25 dBZ, and became more pronounced during the nocturnal period. These echoes were causing spurious convection to be analyzed by the data assimilation system. Therefore, while echoes below 5 dBZ are still assimilated as nonprecipitation observations, only echoes greater than 30 dBZ are also assimilated. Radar observations where the reflectivity is greater than 5 dBZ and less than 30 dBZ are simply omitted from consideration since we do not know a priori whether they are representative of clear air echoes or precipitating hydrometeors. This approach is similar to the method used by Jones et al. (2016), except that the minimum assimilated reflectivity is 30 dBZ instead of 10 dBZ, due to the enhanced nonmeteorological echoes in these nocturnal observations. In the future, dual-polarization radar variables should be added to the quality control (QC) process to better remove nonmeteorological reflectivity (e.g., Mueller and Larkin 1985). Here, we simply omit the questionable data as described above since we find this method to maintain both the assimilation of robust convection and the suppression of spurious convection where observed reflectivity is less than 5 dBZ, while also eliminating the spurious convection induced by nonmeteorological echoes near the radar locations. Furthermore, a simple superob method is applied to the reflectivity observations to reduce the noisiness of the observations and the large number of observations to be assimilated. The superob method consists of dividing the domain into a grid with 0.03° (i.e., ~3 km) horizontal grid spacing and 500-m vertical grid spacing, then averaging the reflectivity observations within each grid box before assimilating them. We also note that radial velocity observations were only assimilated where reflected exceeded the higher 30-dBZ threshold.
4. Overview of case and observational data
a. Synoptic environment
The large-scale environment at 0000 UTC 11 July 2015 is summarized using the outer-domain ensemble mean analysis from the data assimilation system described in section 3 (Fig. 2). At 500 hPa (Fig. 2a), the region of Kansas and Nebraska was generally characterized by westerly flow of ~10–15 m s−1. There was also a mesoscale maximum in the 500-hPa wind speed in southwest Nebraska (Fig. 2a) that may have played a role in convection initiation around this time by increasing cyclonic vorticity advection aloft. At 700 hPa, there was a maximum in water vapor mixing ratio greater than ~9 g kg−1 in southwestern Nebraska and western Kansas (Fig. 2b). The wind barbs indicated advection of the moist midlevel air into the region of southern Nebraska and northern Kansas (Fig. 2b) that the mesoscale convective system (MCS) and bore would be observed in over the several hours following this analysis. At the surface, generally southerly flow in northern Kansas and far southern Nebraska gave way to a stronger easterly component of the flow in central Nebraska (Fig. 2c), indicating that the synoptic-scale warm front had already lifted into central Nebraska.
b. Reflectivity and IOP observations
The evolution of the MCS and subsequent bore are qualitatively summarized in Fig. 3 by a composite reflectivity mosaic, including the NEXRADs at KTWX, KUEX, and KGLD. At 0300 UTC (i.e., 2200 LST) 11 July (Fig. 3a), there were three discrete supercells in south-central Nebraska, which had grown upscale into a small MCS by 0500 UTC (Fig. 3c). At 0610 UTC (Fig. 3d), a reflectivity fine-line feature could be seen ahead of the southeastward MCS motion (parallels red curve on Figs. 3d–i). This feature, which will be shown to be associated with a bore, was transected by both the DC8 and UWKA aircraft during the 10-min window centered at 0610 UTC (Fig. 3e). Aircraft transects of the bore also occurred at 0630, 0700, 0730, and 0740 UTC (Figs. 3e–h). Although new convection initiates northwest of the original MCS (Figs. 3f–h), the original MCS began to dissipate after about 0740 UTC, as indicated by the decreasing area and intensity of reflectivity near the bore front, especially after 0740 (Figs. 3f–i). By 0900 UTC, the clear air radar fine lines associated with the bore had moved nearly 100 km away from the parent convection and into the array of fixed and mobile PECAN Integrated Sounding Array (PISA) instruments (Fig. 3i). Although not seen by the NEXRADs, the continued existence of multiple parallel fine lines in the PISA array was confirmed using the S-band dual-polarization (SPOL) radar located near Hays, Kansas, for the PECAN experiment (not shown). The locations of the SPARC, FP3, and CLAMPS AERI and MP2 Doppler wind lidar (DWL) emphasized in this study are indicated in Fig. 3. The aircraft flight segments during the 10-min period centered on each radar valid time are also shown in Fig. 3 as blue (green) dots for the DC8 (UWKA) aircraft.
The three-dimensional structure of the bore and its parent convective outflow were captured by an array of observations spaced along the extent of the bore front during the PECAN IOP. Of particular interest for this study are the potential temperature and water vapor retrievals from the three Atmospheric Emitted Radiance Interferometer (AERI; Turner and Löhnert 2014; Turner 2016, 2017; Wagner et al. 2016) instruments (Fig. 3), as well as the water vapor retrievals from the differential infrared absorption lidar (DIAL) instrument onboard the NASA DC8 aircraft (Ferrare et al. 2016) and the Raman lidar onboard the University of Wyoming King Air (UWKA) aircraft (Wang et al. 2016) (locations indicated in Fig. 3). These observations provided eight samples from which we independently infer the bore amplitude (Table 3). The UWKA transects at 0610 and 0630 UTC are omitted due to postbore cloud contamination (not shown) and redundancy with nearly collocated DC8 transects (Fig. 3). These observations are used instead of rawinsonde and surface observations due to the inability of rawinsonde and surface observations to provide full vertical cross sections of the observed bore.
The vertically pointing DWLs at FP3 (Hanesiak and Turner 2016) and MP2 (Knupp and Wade 2016) qualitatively showed the bore passage at these locations at about 1000 and 0815 UTC, respectively (Figs. 4a,b). The bore undulations were indicated by the alternating pattern of positive and negative vertical velocity after bore passage, with a periodicity at the surface of ~30 min. The vertical velocity averaged over the first two full wavelengths of the bore undulations (i.e., the period between vertical black lines in Figs. 4a and 4b) is shown in Fig. 4c. The mean positive vertical velocity at low levels (with the exception of the lowest ~200 m) was consistent with the mean lifting expected with bore passage and confirmed that the undulations were indeed related to the passage of an atmospheric bore. While these data qualitatively confirm the presence of a bore, their limited vertical extent and partial cloud blocking requires the additional thermodynamic retrievals to fully characterize the three-dimensional bore structure.
Figure 5 shows the eight time–height cross sections of the bore structure in terms of water vapor, as well as potential temperature for the AERI retrievals in Figs. 5a–c. The DIAL instrument did not retrieve temperature, so water vapor is used here for the flight retrievals. To maintain finescale bore structure, no vertical or temporal averaging was applied to the raw aircraft retrievals. The potential temperature retrievals from AERI revealed a subjectively clear picture of the bore structure (Figs. 5a–c). In particular, there was a mean increase in the height of the potential temperature contours at approximately 1000, 0930, and 0745 UTC in Figs. 5a–c, respectively. This abrupt cooling aloft was more pronounced than the gradual upward slope of potential temperature contours at low levels resulting from nocturnal radiative cooling. Undulations in the potential temperature contours were observed after the cooling aloft, consistent with the passage of the undular bore. Both the undulations and lifting of potential temperature contours were consistent with the vertical velocity observations described above. In the DC8 DIAL and UWKA Raman lidar retrievals (Figs. 5d–h), abrupt moistening aloft was also consistent with the mean layer lifting resulting from bore passage and was therefore also used to quantify the bore amplitude, as described below.
For data from each of the eight transects with the UWKA, a time period is subjectively chosen to represent both the prebore and immediate postbore profiles of potential temperature or water vapor as indicated in Fig. 5. To quantify the bore amplitude, we use the fact that water vapor and potential temperature are both conserved in the atmosphere in the absence of mixing or condensation or evaporation. Cloud formation is ruled out from the aircraft data because the presence of cloud would attenuate the signal to the point of unusable data (e.g., ~0634 UTC in Fig. 5e and ~0656 UTC in Fig. 5g). Cloud contamination is ruled out from the AERI data in Fig. 5 using the accompanying cloud base height and liquid water path data (not shown). Significant mixing can occur in bores based on both theory and model simulations (Koch et al. 2008a,b) and is also implied by observations of bores (e.g., Koch et al. 2008a,b). Therefore, subjective evaluation of the calculated bore amplitude is also presented in order to contribute to the understanding of mixing processes in bores.
The bore amplitude is defined as the ratio h1/h0, where h0 is the depth of the prebore stable fluid and h1 is the depth of the postbore stable fluid, in the simple two-layer model (Rottman and Simpson 1989). Previous studies have identified h0 subjectively by inspecting a vertical profile for a strong temperature inversion (e.g., Koch et al. 1991, 2008b). However, since the stability both within and above the inversion can vary greatly from case to case, it is not clear how to objectively calculate h0, such that the same method can be applied to many cases or many model simulations. One possibility is to use a threshold of the vertical potential temperature gradient. However, we find that different thresholds correspond to the subjectively apparent inversion height in different cases (not shown). Here, we propose to define the prebore stable fluid depth in the continuously stratified atmosphere for a given transect by first calculating the prebore and postbore heights for every possible prebore height in that transect. If our above assumptions approximately hold, then the lifting inferred by the height difference should be approximately proportional to the mean vertical velocity during bore passage. We further assume that the vertical velocity during bore passage is maximized near the top of the stable fluid depth and take the prebore height that gives the maximum inferred lifting to be h0. The validity of this assumption is investigated and discussed below in section 5 using model data. The height at which the average water vapor (potential temperature in the case of AERI data) at h0 in the prebore column was observed in the postbore column is taken to be h1. This method is applied to each transect in Fig. 5, and the results are summarized in Fig. 6. We note that the one exception to taking the height of maximum inferred lifting to be h0 is the 0740 UTC DC8 data (Fig. 6d, dot–dash line). There were two distinct maxima of inferred lifting, but the overall maximum at ~550 m was clearly above the very shallow moist layer in Fig. 5f. We believe that this moist layer corresponded to the stable boundary layer at this time and location, and therefore we take the lower maximum at 240 m AGL (~1500 MSL) to have been h0 for this transect. The validity of assuming that the vertical velocity maximum occurs at or near h0 is further investigated below in section 5a using model simulations.
The bore amplitudes calculated using the method described above are summarized in Table 3. The calculated prebore stable fluid depths (column 1 of Table 3) are different largely due to the different times and places that the observations are made. For example, the AERI instruments (Fig. 6b) showed the shallowest prebore stable fluid depth for SPARC, which was the earliest AERI to experience bore passage (~0745 UTC), and the deepest prebore stable fluid depth for FP3, which was the latest AERI to experience bore passage (~1015 UTC). Despite these differences, the FP3, CLAMPS, and SPARC bore amplitudes calculated with this method are all quite similar, ranging from 1.39 to 1.43.2 It is also worth noting that Fig. 6a shows an approximately exponential decrease of h1/h0 (horizontal axis) with increasing h0 (vertical axis). The aircraft-retrieved bore amplitudes (Figs. 6c,e) also showed a similar trend, except with abrupt increases in amplitude at certain heights (e.g., DC8_0740 at ~500 m in Fig. 6c). This pattern of approximately exponential decay of magnitude with height, punctuated by abrupt increases, was not as clear in the AERI data, as the aircraft retrievals due to decreases in the vertical resolution of the AERI data with height, although there were levels where the decreasing amplitude with increasing h0 was halted or reversed in Fig. 6a. This pattern emphasizes that the stratification of the real atmosphere is neither perfectly consistent with the two-layer model nor perfectly continuous.
The bore amplitudes determined from aircraft-retrieved water vapor transects show similar values to each other of 1.85–2.16 that are much higher than the first three rows of Table 3. Subjective examination of Fig. 5 suggests that the difference may be due to the influence of turbulent mixing, which violates the assumed conservation of water vapor and potential temperature that is used to estimate the bore lifting. The impact is particularly pronounced in Figs. 5g and 5h, where subjectively following the strongest vertical water vapor gradient would indicate less mean lifting than is estimated by assuming nonturbulent conservation of water vapor during the lifting process. The difference is caused by the moistening above the level of strong vertical water vapor gradient near h0, which is consistent with the turbulent vertical transport of water vapor that is expected to be occurring in this region of the bore (Koch et al. 2008a). While this effect is most clearly seen in Figs. 5g and 5h, it is also apparent in the other water vapor–based bore amplitude estimates (Figs. 5d–f). Therefore, the three AERI-retrieved potential temperature-based estimates (Figs. 5a–c) are averaged together for quantitative comparisons to model simulated bore amplitudes in section 5. We interpret the results of Fig. 6 and Table 3 as indicating that the observed bore amplitude is about 1.41 (average of the first three rows of Table 3), while the impact of vertical mixing creates an apparent bore amplitude of up to about 2.16. The impact of vertical mixing may have more impact on the water vapor–based estimate because the enhanced vertical gradient of water vapor near h0 is much stronger than the enhanced potential temperature gradient. Furthermore, we note that the AERI vertical resolution is insufficient to resolve the strong vertical gradient of water vapor.
To properly compare the model-simulated bores to the bore observed at different times and locations in section 2, seven cross sections are drawn through each model-simulated bore as shown in Fig. 7. For qualitative purposes, the central cross section (thicker black line in Fig. 7) will be shown, while the mean and standard deviation of bore amplitude from the seven cross sections will be used for quantitative purposes. The cross sections are chosen to focus on the southward-propagating section of the bore that was sampled by the observations presented above in section 4.
a. Bores in forecasts with different lead times
A first look at the bores and cold pools generated by forecasts at different lead times is provided by the vertical velocity at 1 km AGL and 2-m temperature, respectively, in Fig. 7. For all initialization times, the surface temperatures (Figs. 7e–h) show slight warming as a result of vertical mixing, behind the leading edge of the wave train seen in the vertical velocity (Figs. 7a–d), consistent with the expected bore structure (Koch et al. 2008a,b; Haghi et al. 2017). The location of the bore at 0700 UTC is very similar among the 4-, 3-, 2-, and 1-h forecasts in Figs. 7a–d, respectively. There are variations in the bore structure among the forecasts from all four initialization times. These differences are likely caused by differences in the evolution of the parent MCS and its associated cold pool (not shown). In particular, the vertical velocity magnitude and number of undulations in the bore are both reduced in the forecast initialized at 0600 UTC, compared to the others (Fig. 7d). The different bore structure for the 0600 UTC initialization is seen even more clearly in the vertical cross section (Fig. 8) through the thick black line segments on Fig. 7. While there is a hydraulic jump in the depth of the nocturnal stable layer, evidenced by increased height of the potential temperature contours in Fig. 8d, the undulations that characterize the observed bore and earlier-initialized simulated bores are greatly reduced. This difference is likely because the bore starts to pull away from the parent convection and cold pool after ~0500 UTC in the observations (Fig. 3). The data assimilation system, even cycling at 1-km grid spacing, is not providing an adequate analysis to fully resolve the bore itself in the forecast initial conditions. In contrast, the forecasts initialized earlier have an adequate analysis of the parent convection to spin up the bore using the model dynamics during the early forecast period. Therefore, in contrast to the typical expectation that a shorter forecast is a better forecast, properly simulating the qualitative bore structure in this case with this data assimilation system requires a forecast initialization time before the bore pulls away from the parent convection as a distinct entity.
Similar to the observed bore amplitudes presented in section 4, the model bore amplitudes are quantified by first calculating the inferred lifting as a function of height. The inferred lifting profile is then averaged over the seven cross sections to remove some of the noise in the profiles of inferred lifting and gives a more representative value of h0 over the southward-propagating section of the bore (Figs. 9a,d). The value of h1 is then calculated for each cross section using potential temperature, as described in section 4 for the AERI-retrieved potential temperature in the observed bore. In section 4, the assumption was made that the height of the inferred lifting maximum is an approximation of the maximum vertical velocity in the bore hydraulic jump. The vertical profile of simulated vertical velocity is also averaged over the three grid points centered at the vertical velocity maximum in the leading bore wave, along each of the seven cross sections in each experiment (Figs. 9b,e). The heights of the vertical velocity maximum in the model are, in general, similar to the heights of the inferred lifting maximum in the model, although the height of vertical velocity maximum is systematically ~200 m lower than the height of inferred lifting maximum. This discrepancy is likely a result of vertical mixing of potential temperature and is consistent with the similar, but more pronounced, impact of vertical mixing of water vapor discussed above in section 4. Although the inferred lifting maximum is a slightly biased estimate of the height of the vertical velocity maximum, it is used as the estimate for h0 to be consistent with the method applied to the observations (where full profiles of vertical velocity are not available) and because this small bias does not appear to significantly affect the conclusions of this study. We also note that the wave ducting layer, as indicated by the level where the Scorer parameter transitions from a positive to a negative value, is only ~400–500 m in the model simulations (Figs. 9c,f). This suggests that the height of maximum vertical velocity, at least in these model simulations, is determined more by the hydraulic jump of h in the bore than by the top of the wave duct that traps the associated gravity waves.
The means and standard deviations of the seven bore amplitudes (h1/h0) from each experiment are shown in Table 4. The amplitudes of the simulated bores with WRF single-moment 6-class (WSM6; Hong and Lim 2006) microphysics and MYNN PBL scheme (i.e., the first four rows of Table 3), which have been found to work relatively well for predicting nocturnal convection (e.g., Johnson and Wang 2017), show a consistently positive bias (1.59–1.78), compared to the observed amplitude of 1.41, which is in part related to model physics errors. The sensitivity of the results to the model microphysics and planetary boundary layer parameterization schemes are therefore evaluated in the following subsection.
b. Bore predictions with different model physics
The model-simulated bores are expected to be sensitive to the microphysics parameterization, which affects the parent density current, and the boundary layer parameterization, which affects the vertical mixing of surface radiative cooling and therefore the low-level stability. Therefore, the forecast initialized at 0300 UTC with WSM6 microphysics and MYNN PBL scheme was also run with the same 0300 UTC initial conditions but different forecast physics (Table 2 and Figs. 10 and 11). The microphysics schemes considered include Thompson et al. (2008), Morrison et al. (2009), and WRF double-moment 6-class (WDM6; Lim and Hong 2010). As expected, the cold pool structure is very sensitive to the microphysics parameterization (Figs. 11a–d). As a result, there is much more diversity in the simulated bores with different microphysics schemes (Figs. 10a–d) than with different forecast initialization time (Fig. 8). The different simulated bore structures are even more pronounced in the cross sections in Figs. 12a–d, ranging from a much stronger bore with WDM6 microphysics (Fig. 12b) to no bore at all with Thompson microphysics (Fig. 12d). A bore amplitude could not be calculated for the Thompson experiment because there is not an apparent hydraulic jump or density current in the area at the time (Fig. 10d).
Changing the microphysics scheme results in different bore amplitudes as a result of different cold pool characteristics. The general tendencies toward weaker and warmer cold pools with Thompson and Morrison, and a stronger and cooler cold pool with WDM6, are seen in Figs. 11a–d. As a result, the bore amplitude for WDM6 (MORR) is stronger (weaker) than for 0300init (Table 4). The reasons for the cold pool differences are further illuminated by Fig. 13. The two main differences among 1-h forecasts (i.e., at 0400 UTC, when the cold pool that later impacts the bore is being generated) with the different microphysics schemes, averaged over the rectangular box in Fig. 11a, are 1) the partitioning of the hydrometeor species among rain, snow, and graupel and 2) the rate of low-level rain evaporation inferred by the decrease in mixing ratio with decreasing height at low levels (i.e., below cloud base). In particular, Thompson and Morrison have especially high snow mixing ratios at 8–10 km AGL at the expense of rain and graupel in Thompson and at the expense of rain in Morrison (Figs. 13b–d). Below about 3 km, the rain mixing ratio decreases sharply with decreasing height in WDM6 (Fig. 13b), consistent with the cooler low-level temperatures for WDM6 (Fig. 11d). In contrast, Thompson shows nearly uniform rain mixing ratio below ~2.5 km, which may reflect a less active evaporation process in this scheme for this case.
The experiments with different PBL schemes, including MYNN, QNSE, MYJ (Janjić 1994), Yonsei University (YSU; Hong et al. 2006), and the asymmetric convective model, version 2 (ACM2; Pleim 2007), show similar cold pool structures to each other (Figs. 11a,e–h) and similar strength and number of undulations in the bore (Figs. 10a,e–h). However, unlike the MYNN experiment in Fig. 12a, the cross sections in Figs. 12e–h all show, to some degree, the expected undular bore structure (e.g., Koch et al. 2008a). The bore simulated with MYNN is thus markedly less consistent with the expected bore structure than the other PBL schemes. Thus, while MYNN during data assimilation provided the most realistic prebore environment, it also provided the least realistic simulation of the bore itself. The bore amplitudes (Table 4) also support this conclusion, since the bore amplitude for MYJ, QNSE, YSU, and ACM2 of 1.46, 1.46, 1.44, and 1.54, respectively, are much closer to the observed value of 1.41 from Figs. 4a–c. It is also noted that the bore structure qualitatively varies in the along-bore direction (Fig. 14). The six MYJ cross sections from the east and west of the central cross section shown in Fig. 12e reveal differences in the amplitude ordering of the bore undulations. In particular, the cross sections to the west of the central cross section show the potential temperature contours in the second undulation, reaching a higher height than in the first undulation (Figs. 14a–c), while farther east, the heights of the first two undulations are more similar (Figs. 14d–f). Although the observation network was still not dense enough to confirm that this variability is realistic, along-bore variability of model bore structure emphasizes the need to consider multiple cross sections when quantifying bore amplitude.
The unrealistic bore structure with MYNN is explained by comparing cross sections of turbulent kinetic energy (TKE) among the local mixing schemes (Fig. 15). TKE is a predicted quantity that can be used to diagnose the amount of vertical mixing in the local schemes of MYNN, MYJ, and QNSE. MYNN is different from MYJ and QNSE in terms of the magnitude of TKE near the top of the cold pool–bore system (Fig. 15). The larger values of TKE explain the fact that the bore undulations have been completely mixed out in the MYNN forecast. This mixing should only occur for bore amplitude greater than 4 (Rottman and Simpson 1989), which is much larger than the bore shown here (Table 4).
6. Summary and discussion
Bores are a characteristic part of the nocturnal convective environment over the Great Plains. They influence nocturnal convective initiation and maintenance through mean upward layer lifting and turbulent mixing over mesoscale regions (Koch et al. 2008a; Haghi et al. 2017; Parsons et al. 2018). Past studies have shown that NWP models are often able to resolve bores that appear similar to observed bores. However, there has not yet been, to the authors’ knowledge, an investigation into the skill and uncertainties of NWP-based bore predictions in the Great Plains nocturnal environment and, in particular, the dependence on model physical parameterizations. The present paper begins to address this need by validating bore forecasts during a PECAN IOP against the special observations collected during the IOP. Convection-permitting model forecasts are evaluated for different initialization times and with different model physics configurations. Ongoing future work will investigate improvements to bore forecasts over many bore-focused IOPs that are obtained by assimilating the special IOP observations and by further increases to model horizontal and vertical resolutions.
For this study, the IOP observations were withheld from the data assimilation system in order to be used as completely independent model validation data. The comprehensive observations of this bore during the PECAN IOP were first used to quantify the bore amplitude, which determines how the bore affects the nocturnal convective environment for potential subsequent convection. A new method of objectively calculating bore amplitude consistently in both observations and model simulations is proposed in this study. The method consists of determining the inferred lifting for parcels starting at different altitudes and using the parcel with the maximum inferred lifting to identify h0. This method assumes that such height approximately corresponds to the level of maximum vertical velocity. This assumption is shown to be a reasonable approximation in the model simulations, although with a positive bias that is hypothesized to be due to the impacts of vertical mixing. Koch et al. (2008a) have also demonstrated that in addition to the lifting of low-level layers of air that may be potentially unstable (Schultz et al. 2000), bore-induced mixing also alters the temperature and water vapor profiles of the postbore environment. Here, we speculate that this effect impacted calculations of observed bore amplitude based on water vapor retrievals. The bore amplitude obtained by potential temperature retrievals was less affected by the mixing and therefore used for model validation. Bore speed is related in theory to bore amplitude and therefore not emphasized in this study. However, bore speed can also be affected by model processes independent of the bore itself, such as the low-level jet. Future work will therefore also directly evaluate bore speed in addition to amplitude.
For this case, forecasts initialized before 0600 UTC show a more accurately simulated bore structure than the 0600 UTC initialization. This difference is because the earlier forecasts spin up the bore from its parent convection during the forecast period, while the later forecast is not able to properly maintain the slightly smoothed bore already present in the ensemble mean analysis. Differences exist between simulations with different initialization times before bore formation, but to a lesser degree than compared to an initialization after bore formation. Recent studies have begun to use dual-resolution ensemble data assimilation techniques to provide a higher-resolution single analysis than the resolution of the ensemble members (Wang and Wang 2017; Lu et al. 2017). Short-term forecasts initialized after bore formation may require such a data assimilation system in order to efficiently provide a higher-resolution analysis without the extra smoothing caused by ensemble averaging.
Bore forecasts are primarily sensitive to the microphysics parameterization scheme because of the impact of the microphysics on the convectively generated cold pool. While the WSM6 cold pool agrees well with observations for this case, the Morrison and (especially) Thompson cold pools are much too weak, and the WDM6 cold pool is much too strong. These differences correspond to the partitioning between rain and snow/graupel hydrometeors and the evaporation rate of rain hydrometeors below cloud base inferred by their decrease of mixing ratio with decreasing height. Specifically, the Thompson and Morrison schemes generate more snow and ice hydrometeors at the expense of generating rainwater. This may reduce the amount of precipitation at the surface because of the lower fall velocity of snow and small ice particles in comparison to rain drops. Furthermore, the WDM6 scheme shows a particularly rapid decrease of rainwater mixing ratio with distance below cloud base, likely indicating enhanced evaporation rates, compared to the other microphysics schemes.
The bore behavior with different PBL schemes was generally similar, with the pronounced exception of the MYNN scheme. The MYNN scheme generated very high TKE at the bore/cold pool top, which completely mixed out the undulations in the bore. While this behavior is expected for bore amplitudes greater than 4, the bore amplitude here is less than 2. Therefore, it appears that the MYNN PBL scheme is not optimally suited for explicit bore prediction. Compared to the similar MYJ scheme, MYNN uses new closure constants, a diagnostic equation for the turbulent length scale that increases with decreasing stability, and a term in the TKE prognostic equation that considers “buoyancy effects on pressure-strain and pressure-temperature-gradient covariances” (Nakanishi and Niino 2009). Further research is needed on how these changes respond to the strong stability of the nocturnal inversion, especially in the presence of vertical velocity and buoyancy perturbations associated with undular bores. Although MYNN poorly simulated the undular bore structure, in practice it must also be considered that MYNN has been found to forecast the locations of nocturnal convective rainfall better in general (Johnson and Wang 2017; Blake et al. 2017). Meanwhile, the QNSE scheme has been shown to provide the best representations of the nocturnal low-level jet (Smith et al. 2016) and also provided quite good simulations of the explicitly simulated bore in this study. These results suggest that with the exception of the MYNN scheme, the model representation of the nocturnal boundary layer and its evolution in time is a secondary source of forecast uncertainty in this case, compared to the more pronounced sensitivity to the microphysics parameterization.
This initial study, focused on a single case study, provides a better understanding of the issues relating to practical predictability of a bore that is hypothesized to be representative of typical bore behavior in the Great Plains nocturnal convective environment. Although much caution is warranted when generalizing from a single case study, it appears that the general NWP forecasting of nocturnal convection might be improved by a focus on typical cold pool structures of the microphysics parameterization and the mixing processes within and near nocturnal convective systems that are at least partially parameterized by the PBL scheme. These results also suggest that although uncertainties exist in the model treatment of the stable PBL (Holtslag et al. 2013), microphysics parameterization is a more important source of uncertainty in NWP predictions of nocturnal convection and its interaction with bores. It is also noted that all of the simulations considered in this study failed to resolve the eastward-propagating section of the bore that was seen in observations. Therefore, none of the model configurations evaluated in this study were truly optimal. Work is ongoing to determine if further improvements to the model resolution and observation network for data assimilation (i.e., assimilation of the IOP observations) can lead to improved bore forecasts and therefore improved nocturnal convection forecasts. There were approximately 10 official and unofficial bore-focused IOPs during the PECAN field campaign. These data will allow us to build on the present initial study by systematically determining a baseline level of NWP-based bore prediction skill and quantify the improvements to skill provided by assimilation of the experimental PECAN observations and enhanced horizontal and vertical model resolutions.
The work is primarily supported by NSF Award AGS-1359703. The third and fourth authors also acknowledge NSF Award AGS-1237404. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation Grant ACI-1053575. Data were provided by NCAR/EOL under the sponsorship of the National Science Foundation. The various participants and agencies that made the PECAN project possible and contributed to IOP data collection are especially appreciated. The authors would like to specifically acknowledge the Atmospheric Radiation Measurement program for their support of the FP3 AERI instrument that collected data used in this study. The manuscript was greatly improved by the feedback from three anonymous reviewers.
This article is included in the Plains Elevated Convection At Night (PECAN) Special Collection.
In this study, the undisturbed stable fluid is the nocturnal boundary layer. The former terminology is used for continuity with earlier studies using the more general, though simplified, framework of a two-layer model.
Linear interpolation of potential temperature between vertical levels is used to reduce the impact of the decreasing vertical resolution of AERI retrievals with height.