Downslope windstorm study in the Isthmus of Tehuantepec using WRF high-resolution simulations

Tehuantepecers or Tehuanos are extreme winds produced in the Isthmus of Tehuantepec, blowing south through the Chivela pass, the mountain gap across the isthmus, from the Gulf of Mexico into the Pacific Ocean. They are the result of the complex interaction between large-scale meteorological conditions and local orographic forcings around Chivela pass and occur mainly in winter months, due to cold air damming in the wake of cold fronts that reach as far south. They can generate localized 5 extreme events, such as downslope windstorms and hydraulic jumps, strong turbulent flows that have a direct effect on the Pacific side of the isthmus and the Gulf of Tehuantepec. This study focuses on investigating these phenomena using high horizontal and vertical resolution WRF (Weather Research and Forecasting) model simulations. In particular, we employ a four nested grid configuration, with up to 444 m horizontal spacing in the innermost domain and 70 hybrid-sigma vertical levels, 8 of which lie within the first 200 m above ground. We 10 select one 36 hour period in December 2013, when favorable conditions for a strong gap wind situation were observed. The high-resolution WRF experiment reveals a significant fine-scale structure in the strong Tehuano wind flow, beyond the well known gap jet. Depending on the Froude number upstream the topographic barrier, different downslope windstorm conditions and hydraulic jumps with rotor circulations develop simultaneously at different locations east of Chivela pass with varied crest height. A comparison with observations suggests that the model accurately represents the spatially heterogeneous intense 15 downslope windstorm and the formation of mountain wave clouds for several hours, with low errors in wind speed, wind direction, and temperature. Copyright statement. TEXT

elevation of the Chivela pass is 224 m, whereas mountain peaks in the side sierras reach 2000 m, creating ideal conditions to generate a powerful wind corridor (Romero-Centeno et al., 2003). In winter, cold high-pressure systems originating in North America move over the Gulf of Mexico in the wake of south-reaching cold fronts, and large pressure differences develop across the isthmus between the bay of Campeche and the Gulf of Tehuantepec, on the Pacific side. This pressure gradient results in a northerly wind situation, in which the flow is accelerated southward by cold air damming, traveling through Chivela Pass to 5 finally blow violently outward into the Pacific Ocean. In the Gulf of Tehuantepec, the strong sea surface wind stress generates intense upwelling and vertical mixing in the upper ocean (Hong et al., 2018). These powerful mountain gap winds are called Tehuantepecers or Tehuanos, and have been the focus of several previous studies (Brennan et al., 2010;McCreary et al., 1989;Jaramillo and Borja, 2004) detailing the general setting, drivers and dynamics (Steenburgh et al., 1998) of strong wind situations in the isthmus. The largest number of these events tend to occur in December, with a mean duration of 48h (Brennan 10 et al., 2010). Little knowledge exists, however, on the fine-scale structure of the Tehuantepecer flow, which is prone to result in downslope windstorms (DSWS hereafter) potentially producing severe turbulent phenomena such as rotors and hydraulic jumps (Durran, 1986a;Sheridan and Vosper, 2006) on the Pacific side of the isthmus, to the lee of local orography, particularly during the cold season. There is evidence from observations and earlier numerical studies (Steenburgh et al., 1998) that the low elevation topography of Chivela pass can excite mountain waves, and also that these are not only restricted to the pass itself 15 but extend into the much higher mountain crests to the west and especially to the east, as the cold air pool is often thick enough to surpass them. The ability to understand and forecast these events is very relevant, since the Isthmus has been an important development site for wind farms since the 2000's (Coldwell et al., 2017). Currently this region allocates 76.8% of the wind power capacity installed in Mexico, with approximately 2360 MW (Baxter et al., 2017), which is expected to double to 5076 MW by 2020 (AMDEE ,2018). In addition, several accidents related to the strong winds are reported by the Oaxaca´s Civil 20 Protection Commission (Santiago, 2018;Hernández, 2018;Televisa, 2018;Rodríguez, 2018) every year during some Tehuano occurrences.
The main goal of the present work is to study the variability of flow behavior in Tehuano wind episodes across the isthmus of Tehuantepec, depending on topographic barrier height and thermodynamic conditions of the air mass, using high-resolution simulations with the Weather Research and Forecast (WRF) model. Many studies have successfully employed WRF to analyze 25 this kind of mountain-flow events in other parts of the world. In the US, for example, B. Pokharel et al. Pokharel et al. (2017b) study a DSWS and hydraulic jumps to the lee of the Medicine Bow Mountains in southeast Wyoming. Another DSWS in this same area is also investigated by Grubišić et al. (2015) with WRF. Other studies such as Cao and Fovell (2016), Pokharel et al. (2017a), Prtenjak andBelusic (2009), Jung-Hoon andChung (2006), Ágústsson and Ólafsson (2014), Priestley et al. (2017) also use WRF at high resolution to analyze mountain wave flows in different locations. To the best of our knowledge, there is 30 not any previous work that studies in detail lee wave phenomena in the Mexican state of Oaxaca.
The high-resolution WRF simulations employed in the present study allow us to obtain a more complete knowledge of these events, from the synoptic scale to the microscale, focusing on the downslope winds and the hydraulic jumps that develop along the Chivela pass and the neighboring mountain ranges. The article is organized as follows: in section II the methodology is explained in detail, from the climatology of the region to the model configuration. In section III, the primary results obtained are shown, divided by the synoptic-mesoscale situation, the upstream-downstream structure of the phenomena and the microscale situation. Finally, in section IV the conclusions reached are discussed.
1.1 Mountain wave phenomena, hydraulic analog and relation with the Froude number Downslope windstorms result from the intense flow acceleration occurring on the lee slope of a mountain under certain crossbarrier wind conditions. They resemble the behavior of water flow in an open channel when encountering an obstacle, such 5 as when a relatively slow river increases its speed when flowing in a thin layer over a mill's dam or over a rock. As in the case of the river, downslope windstorms often end abruptly with a return to the state upstream of the obstacle through a turbulent hydraulic jump somewhere downstream. The similarity in both fluid's behavior suggests that the physical processes behind are also alike, and that shallow-water theory could be applied somehow to the atmosphere (Long, 1953). However, the complexity of the unbounded atmospheric flow, without a free surface, makes it difficult to make the analog so simple, because 10 gravity waves in the atmosphere propagate vertically in addition to horizontally as in shallow water, and non-linear effects are important. In shallow water theory, the Froude number (Fr), which is the ratio of the mean speed to the intrinsic gravity wave phase speed in the fluid, determines its behavior when encountering the obstacle, depending on whether gravity is balanced mostly by acceleration (Fr>1) or pressure gradient forces (Fr<1). Hydraulic jumps and significant flow acceleration to the lee of the obstacle occur when the fluid transitions from a subcritical (Fr <1) to supercritical (Fr>1) regime at the top of the barrier. 15 In this case the fluid thins and speeds up while approaching the crest like in subcritical conditions, but instead of returning to the prior state, it keeps thinning and accelerating on the way down as a supercritical flow (Durran, 1986b). Fr must, thus, be around 1 from the start, or in other words, the relation among flow speed, flow depth, which determines gravity wave speed, and obstacle height should be about right.
Considerable observational and numerical experiments have been developed to elucidate the dynamics of atmospheric moun-20 tain lee flows (see Durran (1990) for a review). The current view is that in the atmosphere, downslope windstorms are observed when stable air at low levels flows, similarly to shallow water, as like having a free surface somewhere above the obstacle preventing vertical energy dissipation. This occurs when gravity waves do not propagate vertically and break because of the presence of a critical level  or due to overturning for having too large amplitude (Clark and Peltier, 1977), or when without wave breaking, there is an interface separating highly stable lower layers from less stable air above 25 (Durran, 1986b;Klemp and Durran, 1987;Bacmeister and Pierrehumbert, 1988). Wave breaking creates a well-mixed region to the lee of the obstacle that induces flow separation; the generation of a dividing streamline between undisturbed flow above and trapped energy and flow analogous to hydraulics in the lower surficial branch (Smith, 1985b). A highly stable lower layer topped by less stable air results in reflecting and decaying waves aloft, enhanced non-linear effects, and the atmosphere also flowing like a two-layer fluid and behaving like shallow water when encountering the barrier. The height of the aforementioned 30 dividing streamline or the thickness of the highly stable layer at low levels play a similar role to that of the depth of shallow water for atmospheric flows over an obstacle. In either case, strong flow acceleration on the lee slope and hydraulic jumps are produced when there is a transition from a subcritical like regime on the windward side of the mountain to supercritical like conditions on the lee, as in hydraulic theory. Increasing speeds upslope are due to pressure gradient forces in the standing gravity wave generated by the mountain dominating over gravity (subcritical flow). If the barrier is high enough, velocities at the crest may become so large that lee side gravity-wave induced pressure gradient forces are unable to decelerate the flow, gravity dominates (supercritical flow) and speeds further increase on the way down attaining very high values (Durran, 1990).
It is not straightforward to define a Froude number to determine the regime of atmospheric flows as it is for shallow water, because there is no clear analog to flow depth controlling gravity wave speed and pressure perturbations. Key factors determin-5 ing wave properties and flow characteristics in this case are stratification, wind speed and barrier height, or rather, the relative values among them. A quantity that combines these three variables and is referred to as Froude number for shallow atmospheric flows over mountains by many authors, is the following (Smith, 1989): where U is the flow speed, N , the Brunt-Väisäla frequency (Equation (2)), and H, the mountain height.
with g the acceleration of gravity, dθ/dt, the potential temperature gradient in the stable layer, and θ i , the potential temperature at the base of this layer. To avoid confusion with the classical Froude number, the inverse of Fr defined as above is often used instead and called the non-dimensional mountain height. Fr is a measure of flow deceleration and stagnation upwind from the mountain (Baines, 1987). It can also be regarded as an estimate of nonlinearity. When Fr«1 there are significant non- 15 linear effects and blocking, whereas for Fr»1 the opposite occurs (Smolarkiewicz and Rotunno, 1989). Fr around 1 indicates a transitional regime between the two states and favorable conditions for the formation of DSWS and hydraulic jumps. The flow dynamics related to a DSWSs and HJs in the atmosphere is complex (Clark and Peltier, 1984), and there are different Fr calculations that take into account further factors (Sheridan and Vosper, 2006;Smith, 1989Smith, , 1985a. However, due to the flow type for the period studied, a steady wind due south, and a fairly constant stratification with height in the lower layers, Fr in 20 (Eq3) can be considered a good indicator of the flow characteristics in the considered area, especially in HJ cases.

WRF Configuration
We use the Advanced Research WRF (ARW) model (Skamarock et al., 2008) version 3.9 (WRFV3.9) to perform the simulations. Based on a fully compressible and non-hydrostatic dynamic core, WRFV3.9 is a limited-area mesoscale and microscale 25 model, with a terrain-following hydrostatic-pressure vertical coordinate, designed for operational forecasting, as well as research. For the experiments, we employ a nested domain configuration, in order to achieve sufficiently high resolution in the innermost grids to capture the small scale structure of the flow, while reproducing the synoptic phenomenology conducive to local DSWS in the parent one ( Figure 1a). and gradually transitions to constant pressure at higher levels. The benefit of this vertical coordinate system is a numerical noise reduction in the upper-layers over mountains (Powers et al., 2017). We maintain this fine vertical grid spacing in all the domains, to capture as wide a range of motions as possible over the depth of the boundary layer. other domains, terrain and land use data are from the WRF global standard database, both at 30" resolution for d03 and 2' for d02 and d01.
We simulate a 36-hour period, from 2013-12-23 12:00 to 2013-12-25 00:00 UTC, which registers DSWS conditions in the observational data. Regarding the main physics options, the simulations use the tropical suite configuration (Table1), introduced in WRF version 3.9., except for the planetary boundary layer, which is parametrized by the Shin-Hong scale-aware scheme 5 (S-H) (Shin and Hong, 2015). The next table summarizes this physics configuration used.

Microphysics
Hong and Lim (Hong and Lim, 2006) Cumulus Zhang and Wang. (Zhang et al., 2011) * disabled in d04 and d05 Long wave radiation RRTMG (Iacono et al., 2008) Short wave radiation RRTMG (Iacono et al., 2008) Planet boundary layer Shin and Hong (Shin and Hong, 2015) Surface layer option Revised MM5 surface layer (Jiménez et al., 2012) Land-surface physics Noah land-surface (Tewari et al., 2004)  The S-H planetary boundary layer option is more suitable for the high resolution of the innermost domain (444 m) because it helps to mitigate a double counting effect of the small-scale processes in gray-zone resolutions. Apart from this, this scheme provides a turbulent kinetic energy (TKE) diagnostic variable useful for our analyses.
2.2 Global model and real data 10 Global Forecast System (GFS) analysis data from the National Centers for Environmental Prediction (NCEP) is used as initial and boundary conditions for the WRF model, with a 3-h update interval. The horizontal resolution of this dataset for all variables is 0.5 x 0.5 deg, with 32 vertical levels ranging from 1000 to 10 hPa. The observational data used in this work is provided by the Mexican National Laboratory of remote sensors (https://clima.inifap.gob.mx/lnmysr/Estaciones/MapaEstaciones), collected every 15 minutes at two meteorological stations whose location is presented in Table 2  Model data are extrapolated from their native sigma levels to the height of the meteorological station using equation (3), which relates wind speed with friction wind speed, which is also a diagnostic variable in the model.
Where ws * in the friction velocity, K is the von Kármán constant, z is the height and z 0 is the rugosity.   inversion associated to the subsidence within the high pressure system aloft is also clearly apparent just below 500hPa. Above this level, winds are weak and veer from being southeasterly to southwesterly in the upper troposphere. Below 500 hPa, winds back from an easterly to a northeasterly direction at about 800hPa, and more strongly in lower levels, becoming westerly at the surface, indicating intense cold air advection. There is a pronounced reverse wind shear in the lower troposphere.
As mentioned earlier, the only available observations in the area are from stations MET1 and MET2, whose position is occurrence. At about 0 UTC on the 29 th , the episode decays and winds go back to local breeze regimes. The extent of the simulated period covers the first 36h of the event (shaded in Fig.4), corresponding with its highest intensity in both locations.
These observations away from Chivela pass show evidence, as the simulations suggest, that the neighboring mountains may also induce strong winds, even though certainly not as far reaching as the mountain gap does. The tight isentropes on the windward side of the mountains (to the right in Fig. 5a, d, g and j), indicate a quite similar stable 10 stratification in the whole lower tropospheric column in both cases; stronger in the layers below about 2.5 km and weaker above. The temperature profile in Fig. 3d evidences that the latter level corresponds with the depth of the cold air. Winds are northerly in these lower layers to the north of the mountains, with somewhat higher speeds above 20m/s upstream from MET1 than further east, north of MET2. In both cases, winds above the cool pool are much weaker, and back to an easterly component at about 4000m and above. The latter height is therefore a critical level, since the cross-mountain flow component 15 becomes null; thus, triggered gravity waves will break and dissipate when they approach it. Markowski and Richardson (2010) outline seven conditions conducive to DSWS, albeit not all of them absolutely necessary. These conditions are: a mountain with steeper lee slope (1) crossed by strong winds (> 15m/s) (2) mostly normal to the barrier (3). A stable layer above the top and less stable above that (4) with cold air advection and large-scale subsidence to maintain the stability (5). Apart from this, reverse wind shear above (6) and no cool pool in the lee (7), is also desirable. These conditions are all perfectly met for both 20 locations analyzed, as discussed previously, and indeed intense downslope windstorms occur in both cases.
The stably stratified cross barrier flow displays wave activity from early on, and wave breaking, as the vertical isentropes suggest ( Fig. 5a and d), enhances turbulent mixing ( Fig. 5c and f) and yields a region of weak stability and reverse flow immediately downwind from the mountain crests (Fig 5g and j). In both cases, isentropes on the windward side sink sharply under these layers of low stability on the lee side, much more pronouncedly for the tallest mountain (Fig. 5j). Encompassing the 25 well mixed region to the lee, a split streamline develops (Smith, 1985b), and below its lower branch there is flow thinning and a significant increase in wind speed. The particular features existent on the lee side differ, however, depending of the height of the topographic obstacle. The strong accelerated flow bounded by intense turbulence extends for many kilometers downwind from the lowest mountain (Fig 5i), while ending with a hydraulic jump and rotors at the foot when the barrier is higher (Fig   5j). The formation of either of these lee wave events is related to the Froude number upstream (Smith, 1989).  wavelength and small amplitude and a modest hydraulic jump that propagates downstream can be seen at the initial stages of the episode at 15 UTC December 23 (Fig. 5a). 12 hours later, during nighttime, the aforementioned strong jet extending for tens of kilometers downwind is fully formed, with values above 35 m/s at about 750m above ground and strong turbulence at the surface and in the layers above the jet, where stability is much reduced. The temperature profile upwind (at 3 UTC December 24 and location UP1 in Fig 5) shows stable conditions for the lowest troposphere, specially marked at crest level, Results from the innermost nested grid d05 with the finest resolution (Fig 7) suggest that trapped lee waves develop in the MET1 case within the high stability layer where the strongest winds are found, just below the low stability region and the critical level aloft that prevent their vertical propagation. This wavelike pattern is a common feature in DSWS periods (Pokharel et al., 2017b;Hertenstein and Kuettner, 2005) and fully formed 12 h later (Fig. 7c), extends for more than 100 km downstream 20 aligned with the general orientation in the northwest-southeast direction of the Sierra in the region. Waves are absent further east in the MET2 cross section, where the topographic barrier is higher and a stationary HJ forms instead, as discussed above.
Our results suggest that during Tehuano wind events, the Pacific side of the Isthmus of Tehuantepec east of Chivela pass is very prone to host extreme winds' phenomena. The formation of DSWSs in the area increases the impacts of the already strong mountain gap winds.

Validation
Finally, we contrast our simulation results with the very few data available for validation at meteorological stations MET1 and MET2. The two plots in Figure 8 compare the simulated wind speed (in d04 and d05) with observations from both stations.
Wind speed results from d04 and d05 at location MET1 are quite similar, and fare well with respect to observations ( Figure   8a), slightly better so those from d05, with a mean absolute error (MAE) of 1.55 m/s (Table2). The similarity in the low mean 30 error between both simulations and their high correlation throughout the period are due to the nature of the event in that area, an intense and mostly steady jet that the d04 domain resolution (1.3 km) is already capable of resolving accurately. However, results in MET2, which registers the HJ situations, present more differences between d04 and d05 ( Figure 8b) and there is a significant improvement in d05 with respect to its parent domain d04. Wind speeds in d04 are overestimated (Mean error ME = 2.76 m/s), and present a daily cycle that is absent or very subtle in the observations. The complexity and fast variability of HJs formation in this area are better resolved in the higher resolution grid, which perhaps reproduces more accurately the stagnant flow and rotor formations downstream from the HJ. With regard to wind direction, errors are small for MET1 and significantly higher for MET2, due to the same reasons. As in the case of wind speeds, wind direction results from the finer grid d05 are also better than in d04. Temperature errors are equal or below 1 K in both locations and domains, hence the surface thermal 5 evolution is well captured. Lee waves can promote orographic cloud formation at different scales, depending on the amplitude of the wave and elevation (Armi and Mayr, 2011;Szmyd, 2016). Model results in d05 suggest that lenticular clouds form at the crests of the trapped lee waves depicted in Figure 7. Figure Table 3. Wind speed (WS), wind direction (WD), and temperature (T) mean errors (ME) and mean absolute errors (MAE) at MET1 and MET2 locations during the simulated period and for the two higher resolution domains (d04 and d05).
were indeed observed in the area. The actual existence of these lenticular clouds with the same location and pattern as in the simulation further validates the model results.

Conclusions
In the present work, we studied lee wave phenomena occurring during Tehuano events on the Pacific side of the Isthmus of Tehuantepec using WRF high-resolution simulations. Orographic forcings at different scales result in the well-known gap wind 5 jet off Chivela pass, but also in downslope windstorms and hydraulic jumps in the neighboring mountains. We analyzed these phenomena in an episode in December 2013 having the typical genesis of Tehuantepecer wind events. An Arctic air mass in North America pushed as far south as the bay of Campeche due to cold air damming east of the Rockies continuing to the east of the Sierra Madre Oriental range in Mexico. The displacement of the associated high-pressure system on the wake of the cold front created large pressure differences across the Isthmus of Tehuantepec, ultimately producing the strong mountain gap  winds through the low elevation of Chivela pass. The model simulates these intense winds, blowing with speeds at the surface of more than 25 m/s that extend for many kilometers from the mountain gap, fanning out well into the Gulf of Tehuantepec, as it is commonly observed in Tehuano events (Steenburgh et al., 1998). The depth of the cold surge on the coastal plains of the Gulf of Mexico side of the Isthmus is about 2500m, therefore thick enough to surmount the lower elevations to the west, and especially to the east of Chivela pass. The flow over these mountains 5 results in intense downslope wind storms to their lee, with the generation of intense turbulence, hydraulic jumps and rotors, depending on the particular height of the topography. We focus on two locations of different barrier elevation from where there are surface observations downwind: one of 963m closer to Chivela pass and another further east, with height increasing to 1736m The thermodynamic characteristics of the air mass are rather uniform upwind of both mountains, with strong stability within 10 the cool pool and weaker above, and intense northerly winds that back and weaken aloft to a more easterly component parallel to the barrier. The critical level where the cross-mountain wind component becomes zero, inhibiting wave propagation, is about 4000m. Mountain waves are generated in both cases, with smaller amplitudes for the lower mountain, where the Brunt-Väisälä frequency at crest height is between 0.020 and 0.025 s −1 , than for the higher mountain, where the Brunt-Väisälä frequency is about half. Wave breaking produces mixing and generates a region of low stability to the lee of the mountains, which is deeper 15 where the waves have higher amplitude. The Froude number is around 2.5 at crest height in the lower barrier and the flow presents a supercritical behavior. The region of low stability to the lee of the mountain lies above about 1500m and leads to a packing of the isentropes and streamlines underneath, resulting in strong stability and flow acceleration. An intense jet develops with wind speeds of 35 m/s at about 750m above ground extending for tens of kilometers downwind from the mountains. Wind speeds are reduced closer to the surface due to intense turbulence. Trapped lee waves form at about 1500m, just below the well 20 mixed layer aloft that prevents their vertical propagation. The Froude number decreases to about 1 further east as elevation rises and the flow presents a critical regime. Isentropes on the windward side of the mountain sink much more pronouncedly under the wider mixed layer generated by wave breaking to the lee, and are tightly packed in a shallow layer above the surface. This generates an intense wind storm on the lee slope of the mountain, with surface wind speeds up to 35 m/s. The accelerated flow down the mountain ends abruptly at its foot, where the flow turns to subcritical state and a marked stationary hydraulic jump forms, with vertical velocities of 6 m/s. A rotor circulation develops further downstream from the jump.
Only limited observations are available to validate our model results. Errors in surface wind speeds, directions and temperature are small at the only two stations available on the Pacific coastal plain downwind from the mountains. In addition, 5 lenticular clouds similar in location and pattern to those produced by the model, are apparent in satellite imagery of the day of the event, providing valuable indication that the mountain lee wave phenomena simulated indeed correspond to a real scenario.
Our model results suggest that when the cold air mass intruding from the north on the Gulf of Mexico side of the Isthmus is thick enough to surmount the mountains, extreme wind events develop in the area during Tehuano events beyond the gap wind jet. These include downslope wind storms and hydraulic jumps, which are intense and highly turbulent flows that can 10 have a substantial impact on the existent wind farm industry in the region. It is likely that the depth of the cold pool and how it compares with topographic barrier height, is a key factor determining the extent, location and intensity of these lee wave phenomena and whether they take place at all.