Abstract
A supercell convective storm is simulated by using a cloud-resolving model. Numerical experiments have been performed in 3D by using the same domain size, with a different spatial and temporal resolution of the model. High-resolution cloud model has been shown to represent convective processing quite well. Running the model in a high-resolution mode gives a more realistic view of the life cycle of convective storm, internal structure, and storm behavior. The storm structure and evolutionary properties are evaluated by comparing the modeled radar reflectivity to the observed radar reflectivity. The comparative analysis between physical parameters shows good agreement among both model runs and compares well with observations, especially using a fine spatial resolution. The lack of measurements of these species in the convective outflow region does not allow us to evaluate the model results with observations. A three-dimensional simulation using higher grid resolution mode exhibits interesting features which include a double vortex circulation, cell splitting, and secondary cell formation.
1. Introduction
Convective clouds and storms represent one of the most important and challenging problems for forecasters. The severe local storms and deep convective clouds are characterized by the enhanced transport of heat and moisture in the upper layers, very strong self-organized flow fields, very complex microphysical transformations and stratospheric penetrations, and rapid evolution and dissipation processes. The precipitation processes are activated in very limited time interval and space, and their intensities are manifested by large natural variability. Supercell storms are perhaps the most violent of all storm types and are capable of producing damaging winds, large hail, and weak-to-violent tornadoes. They are most common during the spring across the mid-latitudes when moderate-to-strong atmospheric wind fields, vertical wind shear, and instability are present. The degree and vertical distribution of moisture, instability, lift, and especially wind shear have a profound influence on convective storm type. It is generally recognized that the environmental buoyancy and vertical wind shear have important effect on the characteristics of convective storms. Much of our understanding of the sensitivity of convective storms to these environment parameters has been derived from modeling studies that tested a variety of, but often idealized, environmental conditions. Numerical cloud models have contributed substantially to our understanding of supercell storms.
A number of three-dimensional cloud models have been developed to simulate the structure, intensity, and movement of convective clouds [1–9]. Many previous studies using high resolution cloud-resolving models (or convective cloud models) have shown that case-specific simulations are able to represent the storm structure, upward transport of air and movement, radar reflectivity, wind speed and direction, and outflow heights [10]. Convective-scale model or cloud resolving models can be used to obtain general characteristics of these subgrid processes.
The main objective of numerical experiments performed here is to simulate the supercell storm structural and evolutionary properties by using the same initialization and different spatial resolution and time step chosen to satisfy the CFD conditions. Accurate storm cell identification and tracking that reflects splitting and merging cells is a major challenge in this study. The model is initialized on upper-air sounding representing the initial vertical profiles of meteorological data. Three-dimensional (3D) numerical experiments have been carefully set up in order to simulate storm dynamics, microphysics, and rainfall processes. The storm structure is evaluated by comparing the modeled and simulated radar reflectivity through examination of its horizontal and vertical cross-sections. In Section 2, we briefly describe the convective cloud model, numerical technique, and boundary conditions. The numerical experiment, the initial conditions, and the experimental setup are represented in Section 3. Then we focus on the results of sensitivity experiments with constructive discussion about the thermodynamic conditions and the physical properties of clouds, as well as the radar reflectivity comparison. Then modeled physical data are compared with observations. Finally, results are discussed and summarized in the last section.
2. The Cloud Model
The convective cloud model is a three-dimensional, nonhydrostatic, time-dependent, compressible system using the dynamic scheme from Klemp and Wilhelmson [4], Lin et al. [11] microphysics, and Orville and Kopp [12] thermodynamics. It includes ten prognostic equations: three momentum equations, the pressure and thermodynamic equations, four continuity equations for the water substances, and a subgrid kinetic energy equation. The equations are specified in the Cartesian coordinate system.
2.1. Dynamics and Thermodynamics
The dynamical part of the model is based on the pressure equation and nonhydrostatic compressible equations of motion introduced by Klemp and Wilhelmson in [4]. The equations of motion are derived from Navier-Stokes equations using Boussinesq approximation, and take account of advection, turbulent transport, buoyancy (either due to warming or loading hydrometeors), and pressure gradient force. The pressure equation is derived by combining the compressible continuity equation and the thermodynamic equation.
The subgrid transport and variance terms that required parameterization are derived by performing Reynolds averaging on each of the prognostic variables. The subgrid turbulence equations, obtained in that way, are simplified, and only the prognostic equation for subgrid turbulence kinetic energy is retained. This equation depends on local buoyancy, and shear and dissipation and is used to specify the eddy mixing coefficients. First-order closer is applied for the calculation of turbulence terms in equations for microphysics and thermodynamics.
The thermodynamic energy equation is based on Orville and Kopp [12] with the effects of the snow field added. This equation takes account of the advection, turbulent mixing, and the heating effect when freezing the hydrometeors, or cooling effect when melting, or sublimational cooling effect when hail and snow are outside a cloudy environment.
2.2. Cloud Microphysics
Bulk water parameterization is used for simulation of microphysical processes. Six categories of water substance are included: water vapour, cloud water, cloud ice, rain, snow, and graupel or hail. Cloud water and cloud ice are assumed to be monodisperse, with zero terminal velocities. Rain, hail, and snow have the Marshal-Palmer type size distributions with fixed intercept parameters. Details can be found in Lin et al. [11]. The source reference for the scheme to allow coexistence of cloud water and cloud ice in the temperature region of (40°C to 0°C) is Hsie et al. [13]. Instead of using the hail spectrum from zero to infinity (idealized spectrum), Curic and Janc [14, 15] proposed considering the hail size spectrum which includes only hail-sized particles (larger than 0.5 cm in diameter; hereafter called realistic hail spectrum). Four prognostic conservation equations for the exchanges of water substances are considered in the model. One of the prognostic variables is the sum mixing ratios for water vapour, cloud water, and cloud ice. Other prognostic variables are the mixing ratios of rain, graupel or hail, and snow. It takes into account 6 water variables (water vapour, cloud droplets, ice crystals, rain, snow, and graupel). The graupel hydrometeor class is represented as hail with a density of 0.9 g cm-3. Natural cloud ice is normally initiated by using a Fletcher-type equation for the ice nuclei number concentration. The change of water vapour takes place due to condensation of cloud water, evaporation of cloud water and rain, sublimation of graupel, and sublimation and depositional growth of snow. Cloud water is produced by condensation of water vapour and melting of cloud ice. The processes which transform cloud water to other substances are: evaporation, the depositional growth of cloud ice and snow, the homogeneous freezing of cloud water, the autoconversion of cloud water to form rain, the accretion of cloud water by rain, snow (°K) and hail, and the wet growth of graupel.
Cloud ice changes take place due to melting, the depositional growth from cloud water, the homogeneous freezing of cloud water, the autoconversion to snow, the accretion by snow, rain and graupel, and the transfer to snow by Bergeron process (deposition and riming). Cloud ice is initiated by using a Fletcher-type equation for the ice nuclei number concentration.
Snow may be produced by the following processes: the autoconversion of snow ice, the Bergeron’s growth of cloud ice and transfer of cloud water, the accretion of cloud ice and rain by snow, the accretion of cloud ice by rain, and sublimation and depositional growth of snow. Transformations of snow to the other substances are caused by melting, the autoconversion of snow to form graupel, the accretion of snow by graupel and rain, and precipitation flux.
Rain is produced by the autoconversion of cloud water, the accretion of cloud water by rain and snow , melting of the snow and hail, and shedding during the wet growth of hail. The processes: evaporation, the accretion of rain by cloud ice, snow and graupel, the probabilistic freezing of rain to form graupel, and the precipitation flux are responsible for the reduction of rain.
Hail is produced by the autoconversion of snow, the accretion of snow, cloud ice and rain by graupel, the probabilistic freezing or rain, sublimation, and wet growth. The processes which reduce hail are evaporation, melting, and the precipitation flux.
The graupel hydrometeor class is represented as hail with a density of 0.9 g cm-3. Natural cloud ice is normally initiated by using a Fletcher-type equation for the ice nuclei number concentration. In this version of the model, cloud ice may be also produced by the Hallett-Mossop ice multiplication. The continental cloud droplet is taken to be 0.05 m in diameter. The corresponding value for tropical clouds is 0.15 m. The equivalent radar reflectivity factors for hail and rain are computed by using equations from Smith et al. [16], and the empirical equation for snow is taken from Sekhon and Srivistava [17]. More detailed information regarding the hydrodynamic equations, microphysics equations, turbulent closure, and numerical methods could be found in Telenta and Aleksic [18] and Spiridonov and Curic [19, 20].
2.3. Numerical Techniques
Model uses a nonmonotonic numerical scheme. Model equations are solved on a semistaggered grid, C-grid. For the dissipation we have applied a four-order accuracy filter. All velocity components are defined at one-half grid interval , while scalar variables are defined at the midpoint of each grid. The horizontal and vertical advection terms are calculated by the centered fourth- and second-order differences, respectively. Since the model equations represent a compressible fluid, a time splitting procedure is applied to achieve numerical efficiency. The scalar prognostic equations, except that for pressure, are stepped from to by a single leapfrog step. The terms which are not responsible for sound wave generation in the equations of motion and the pressure equation are evaluated at the central time level . In grid points adjacent to lateral boundaries, the normal horizontal advection terms are approximated using second-order differences instead of the fourth-order ones used elsewhere. At lateral boundaries, the normal derivatives for all prognostic variables are calculated with first-order accuracy, through one-sided differences lagged at time to provide stability. While the size of the model domain was the same, configured to a km3, the resolution of the model was different. The first numerical experiment is performed with resolution km3 with a temporal resolution of 10 s for long time. The second run of the model was set up at a very high horizontal resolution of km3, using a smaller time step of s. Time-splitting procedure is applied in both model runs by using a smaller time step of 2 s. for solving sound waves. At the top of the model, a rigid lid is used; a damping layer at the top of the domain was not included.
2.4. Boundary Conditions
Boundary conditions are defined so that the normal component of the velocity approaches zero along the top and bottom of the model domain. To ensure that a rigid top boundary assumption does not cause vertical oscillations in the numerical simulation, the authors have upgraded the model with a radiative upper boundary condition, as suggested by Klemp and Durran [21]. The lateral boundaries are opened and time dependent, so those disturbances can pass through with minimal reflection by Durran [22]. When the component of velocity normal to the boundary is directed toward the domain (inflow boundary), normal derivatives are set to zero. At outflow boundaries, for all variables except pressure, normal derivatives are calculated by the upstream difference, with a time lag of a large time step in order to ensure stability. Pressure boundary conditions are calculated from other boundary values to maintain consistency.
2.5. Initial Conditions and Initialization
The model is initialized on upper-air sounding representing initial vertical profiles of meteorological data for the 9th of August 2008 near Skopje (Figure 1). Upper-air sounding indicates unstable atmospheric conditions favorable for convection. The main characteristic of the upper-air sounding is weak wind veering in the surface layer and strong wind shear at the middle and upper part of the atmosphere. Vertical profile indicates moisture deficit at the 500 hPa pressure level, and increased moisture content at the 700 and the 300 hPa level.
For initialization we use artificial initiation of convection. It means that we use a standard methodology of perturbation with an ellipsoidal form of thermal bubble where temperature decreases exponentially from the bubble centre towards boundaries in order to produce a singe-cell cloud based on observations which show that isolated convective clouds form over ground where a thermal bubble with the similar characteristics is generated. Thus, the initial impulse for convection is an ellipsoidal warm bubble of the form where Here, the subscript refers to the location of the center of the perturbation, while , , and are radial dimensions of the bubble. Values used for these simulations are , and the temperature perturbation is maximum in the bubble center, distinguishing (°C) and exponentially decreasing to zero on the bubble boundary. The initial perturbation of water vapor mixing ratio, caused by this initial temperature perturbation, is calculated with the assumption that relative humidity has the same value as it had before the perturbation. Of course the initialization of the model with random perturbation including ground surface that heats via solar radiation and evaporates water will give a more realistic storm initiation, with the horizontal scale being less influenced by the modeler.
The horizontally homogeneous initial fields of temperature, humidity, pressure, and horizontal wind are specified from sounding observed in the vicinity of the cloud system under investigation.
3. Results
3.1. Characteristic of the Case Study of the 9th of August 2008
One class of severe thunderstorm that produces some of the most damaging weather is the supercell. These storms contain a cyclonically (anticyclone) rotating updraft and have a propagation component to the right (left) of the vertical wind shear vector. For vertical wind profiles typical of severe weather environments in Macedonia, a cyclonically rotating updraft (mesocyclone) also moves to the right of the mean wind. An anticyclonically rotating updraft (meso-anticyclone) moves to the left of the mean wind. Statistical proceeding reports that right-moving mesocyclones are more prevalent than left-moving mesoanticyclones; left-moving supercells do exist and can produce severe weather.
A severe left-moving thunderstorm occurred on the 9th of August 2008 in the northeastern part of Macedonia between the cities of St. Nikola and Kochani, along the river valley of Bregalnica. Over 2-hour lifetime, the storm was responsible for reports of severe showers and hail, identified as a weak tornado. An upper-air sounding data was taken 12 hours prior to observe the horizontal west wind direction, updraft wind 35 m/sec at 7500 m.a.s.l., and great energy of atmospheric instability.
The atmosphere over Macedonia on the 9th of August, 2008 contained the necessary ingredients for the development of severe weather, namely, moisture, instability, and a lifting mechanism. It is evident from radar observations (made by upgraded WSR 74 S) that the supercell due to topography and anticyclone rotating is splitting the thunderstorm in two single-cloud cells: a left-moving one and a stationary one. The thunderstorm that formed along the river loop and the mountain cusp produced a left mover that traveled eastward along the initial thunderstorm’s low-level outflow boundary. This long-lived left mover contained a meso-anticyclone and was responsible for numerous severe weather reports.
By 1553 UTC, thunderstorm that was formed was producing 53-dBZ reflectivity at 10 km distance of the radar location. At 1612 UTC, the super cell is splitting inside. The leftward moving single cell is continuing moving toward the east with maximum radar reflectivity over 50 dBz. At 1636 UTC, the convective cells are separated at a distance of 30 km from radar site.
Over the next 10 minutes, however, the stationary rightward storm cell rapidly reintensified to the 50-dBZ range, intensity equal to that of the initial storm, which at this point could be considered as the right mover. Still, the left mover was considerably smaller than the right mover in area extent. This size discrepancy remained throughout the left mover’s lifetime.
3.2. A Three-Dimensional Simulation of Supercell Storm
The structure and evolution of the supercell storm, observed over northeastern part of Macedonia on the 9th of August, 2008, is described through a combined observational radar analysis and a numerical modeling study. This convective cloud system was long-lived and exhibited characteristics similar to those of classic supercells, including a cell splitting. The development and evolution of the supercell storm was simulated using a cloud resolving model with upgrade version of bulk parameterization microphysics scheme. The convection was initiated with a warm bubble (1.3°C perturbation) oriented in a WSW to ENE line according to the main convective mass movement. Differences related to storm physical processes which reflect structural and evolutionary properties are mainly due to the different spatial and temporal resolution used in the model simulations. Showing how each model run responds to the same initiation is valuable in itself. Simulations were integrated for a 1.5-hour period. This was a very specific situation to simulate. The main characteristics of convective storm, structural and evolutionary properties, are examined by analysis of the basic dynamical, microphysical, and radar reflectivity parameters.
Three-dimensional simulations of the 9th of August 2008 supercell storm case indicate that the results are sensitive to the initializations. It means that a small temperature perturbation value in the centre of the thermal bubble is needed to initiate convection. Figure 2 shows three-dimensional views of the supercell storm life cycle at 15, 30, 45, 60, 75, and 90 min of the simulation time using a coarser spatial resolution of km3 and time step of s. The general supercell storm appearance is shown through distribution of mixing ratio of cloud water, cloud ice, hail, snow, and rainwater during the simulation time. Initial cloud water has occurred within 15 minutes after the initiation. Hail is formed in 25 minutes while cloud ice, rainwater, and snow have occurred in 30 minutes, respectively. The modeled cloud penetrates the stable layer and then experiences an intensive growth, developing into a vigorous supercell storm with formation of large amount of ice crystals. Numerical simulation of supercell storm in 25, 35, 45, 60, 75, and 90 minutes using finer spatial grid resolution of ( km3), and a smaller time step of s is depicted on Figure 3. It is clearly illustrated that the model run with high resolution mode shows a more realistic view of the supercell storm structure and evolution. Initial cloud water in this model run has occurred 12 minutes after initiation. The numerical experiment using finer resolution illustrates an early formation of hail, snow, and cloud ice relative to the model run with coarser grid resolution. It is also evident that the supercell storm exhibits cell splitting especially early in its lifetime. Even though this storm started small (left), it had no problem dividing itself in two cells, as the result of surrounding winds that supported both leftward movers which tend to spin clockwise (anticyclonic) and rightward moving cell which turns counterclockwise (cyclonic). The reasons for supercell storm splitting involve concepts of fluid dynamics and treatment of subgrid scale processes, both in the original storm and its environment.
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
3.3. Comparison with Observations
3.3.1. Cloud Top History
According to the observations, the best estimate of the cloud top history is that the top was at about 180 mb (59°C) at the time of the first penetration; then it rose to about 110 mb (63.5°C) or 16 km height in its developing stage and sank gradually in mature stage. The model cloud, on the other hand, consistently has a lower cloud top. At 25 minutes, the time we identified as corresponding to the first penetration, and the cloud top was at about 218 mb (46.3°C), rising up to the 173 mb (56°C) or 13.1 km height m.a.s.l. and then gradually sinking back. Thus, the model underestimates cloud top for about 63 mb.
3.3.2. Updraft Velocity
Both cases show a rapid increase in peak updraft velocity at the beginning of the simulation. The maximum updraft velocity of 27.7 m s-1 in the first numerical experiment is calculated in 60 minutes of the simulation time in the cloud mature stage. The model run using finer grid resolution shows an increased updraft velocity of 31.3 m s-1 in the early stage of supercell storm evolution. The height of the peak updraft reaches 7.5 km m.s.l., which is similar but somewhat higher than observations. Model runs with finer and coarser resolutions maintain peak updrafts during the remainder of the simulation for about 15.3 m s-1 and 11.4 m s-1, respectively.
3.3.3. Liquid Water Contents
Radar reflectivity information is often displayed in two dimensions, making it difficult to extract the structural characteristics of convective storms. The maximum radar reflectivity and the vertical profile of liquid water distribution in a vertical column of a convective cell are used to determine the structural and intensity classification of the cell. Data set provided for the same convective case gives the pass-average values of the liquid water contents. We have calculated average LWC values for the entire horizontal cross-sections of the cloud at different vertical levels. Timing of model calculated averaged values of liquid water contents is consistent with radar observations. Comparison between the modeled and measured LWC values in Table 1 shows a relatively good agreement, with the slight underestimate in heavy precipitation period which is attributed to the fallout of the precipitation in the model. These systematic differences are more evident in numerical simulation using a coarser spatial and temporal resolution, as the result of the increase in modeled precipitation.
3.3.4. Radar Reflectivity History
The storm structure can be evaluated by comparing the modeled radar reflectivity to the observed radar reflectivity. In order to achieve that we have compared horizontal and vertical cross-sections of radar reflectivity calculated in different simulation times with observed parameters. The first radar reflectivity echo at 15:37 UTC shown in (Figure 4), viewed on the 10-sm radar reflectivity maps, indicates existence of an isolated convective core with 10 km diameter, slowly moving in a west-southwest line with an anvil spreading to the east-northeast. In the next 60 minutes, the air mass thunderstorm is successively extended over the affected area, separating in two cores (Figure 7). The frontal core is illustrating increased radar reflectivity patterns compared to the backward core. At 16:38 UTC, a multicellular convective system has two separate radar patterns with maximum reflectivity echoes greater than 60 dBz. Vertical cross-section of the simulated reflectivity clearly illustrates 2 cells which reach only 11.5 km, a.m.s.l. (see Figures 5 and 6). Modeled radar reflectivity (dBz) at km m.s.l. after 40 minutes of simulation has shown a similar pattern with a slight increase of magnitude of the reflectivity compared to observations (Figures 8 and 9). The observations show the reflectivity top to be 14.5 to 16.5 km a.m.s.l. However, during the mature stage of the storm, 2 to 4 convective cells were observed. After 1 hour of simulation, the results from the model have 2 to 3 convective cores oriented west-northwest-northeast which is in line with the observations. The magnitude of the reflectivity is similar with observations. The only slight difference is due to (1) treatment of graupel or hail, (2) model resolution, and (3) single-moment versus multimoment microphysics parameterizations. The width of the anvil varies among models. The observed reflectivity has an anvil width of 32–40 km at 16:12 UTC, while model results range from 12.5 km to 45 km. Seifert and Weisman [23] noted that double-moment microphysics parameterizations tend to produce broader anvils than single-moment microphysics parameterizations do. The results from our study do not distinctly show this correlation. Other factors contributing to the anvil width are the graupel or hail characteristics used (which influence the particle’s fall speed), the dynamics formulation, the vertical or horizontal resolution, and the number of bubbles used to initiate the convection.
(a)
(b)
(c)
(d)
(e)
(a)
(b)
(c)
(d)
(e)
(a)
(b)
(c)
(d)
(e)
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
Both cases show approximately the same reflectivity magnitudes somewhat extensive (60 dBz) in cloud developing stage and slight decrease of about 10 to 15 dBz in the cloud mature stage. The maximum height of the modeled reflectivity varies among runs. The reflectivity for first model run reaches 11.5 km and 12.5 km m.s.l., respectively. Early formation of precipitation, a total accumulated amount of 51.6 mm, with relatively short (20 minutes) heavy precipitation period appeared during the midlatitude simulation.
4. Conclusions
Three-dimensional numerical simulations have been performed running by a model with the same initializations and using only different spatial and temporal resolution. Numerical simulations of the cloud system duplicate the general observational features, including horizontal and vertical dimensions, cyclic behavior, and convective core and anvil characteristics.
The comparison of the horizontal and vertical cross-sections of radar reflectivity echoes in different stages of the multicell storm evolution in general shows a relatively good agreement with observations. This difference in rainfall efficiency is probably attributed to differences in the interaction of cloud dynamics and microphysics, and precipitation flux processes. Both model runs have reproduced the observed convection with radar reflectivity reaching 50 dBz. Both numerical experiments simulated the development of supercell storm structure. However the model run with finer grid resolutions has been more realistic and accurate in simulation of storm splitting. A three-dimensional simulation using higher grid resolution mode exhibits interesting features which include a double vortex circulation, cell splitting, and secondary cell formation. We found that a very small change in the grid resolution of the model can produce very different behaviours of storms after their splitting. In spite of the assumptions mentioned in the explanation of the model, the results of the simulations are self-consistent, have a clear physical explanation, and in general show an encouraging degree of realism. Fine resolution model runs even with the uniform bottom boundary conditions show a promising result in respect to better treatment of the subscale processes of convective transport in the convective storms caused by the complex turbulent motion caused by the effects of buoyancy, wind shearing, Coriolis force, and viscous energy dissipation.
One of the important aspects of this research might be possibilities in opening new ways in treatment of the model initialization by replacing artificial initiation of convection by using a random perturbation and including ground surface with heat via solar radiation. That may give a more realistic storm initiation with the horizontal scale less influenced by the modeler. This model can improve our understanding of the changes within the atmosphere during convective storms.
Acknowledgments
The authors would kindly like to acknowledge the Hydrometeorological Meteorological Service of Macedonia for the provision of initial data, radar images, and rainfall data for the convective case experiment. The anonymous referees have made many substantial remarks that increased the completeness and clarity of this paper. The authors would like to express their special thanks to Dr. George Kalos for his important contribution giving them constructive comments during the final revision of the paper.