Applications of Non-Equilibrium Thermodynamics and Statistical Mechanics of Vortex Gases in Tornado Theory

Abstract

A recent numerical study of tornadogenesis suggests that processes and structures in the turbulent boundary layer below a supercell thunderstorm’s mesocyclone may play an important role in tornado formation. A theoretical explanation using concepts from turbulence theory would be desirable. This work puts into mathematical, statistical mechanics, and thermodynamical context the initial stages of the genesis of tornado-like vortices with the aim to be consistent with the current state of knowledge of the process of tornadogenesis. In particular, it discusses a mathematical foundation of the formation of coherent structures such as “cusps” and “hairpins” using variants of the nonlinear Schrödinger equation that arise via the Hasimoto transform of a vortex filament model. The behavior of such structures is then analyzed within a quasi-two-dimensional boundary layer model using the statistical mechanics of vortex gases to explain the rearrangement of cusps and other vertical vortex filaments into patches and possibly supercritical vortices. Non-equilibrium thermodynamics is used to obtain the entropic balance and the internal entropy production rate, and connect them to the turbulent heat flux. A formula for the non-equilibrium turbulent heat supply and formulas for the entropy supply and entropy production in the boundary layer are also provided. An example involving vorticity, potential vorticity, and the gradient of potential temperature is given in a thermodynamic non-equilibrium context with its implications for stretching and tilting of vorticity in the vertical direction. We conclude with some remarks on the equivalence of the Schrödinger and Gross-Pitaevskii equations in describing vortex filaments.

Share and Cite:

Bělík, P., Dokken, D.P. and Shvartsman, M.M. (2025) Applications of Non-Equilibrium Thermodynamics and Statistical Mechanics of Vortex Gases in Tornado Theory. Open Journal of Fluid Dynamics, 15, 175-205. doi: 10.4236/ojfd.2025.153011.

1. Introduction

A recent survey paper by Fischer et al. [1] has outlined the current state of knowledge about the process of tornadogenesis, reducing the process to a conceptual model with four stages. They are:

1) Mesocyclone formation, with the associated pressure perturbations driving upward acceleration and stretching near and below cloud base (critical for Stage 3 below).

2) Generation of vortex patches (amorphous surface-level vertical vorticity not yet with vortex characteristics), in simulations predominantly via reorientation of baroclinically generated vorticity in downdrafts.

3) Organization of one or multiple vortex patches to a single symmetric vortex and intensification via stretching.

4) Transition to a developed tornado with a vortex boundary layer and corner flow region in which horizontal vorticity is tilted into the vertical right above the surface.

Observational and numerical motivation for this conceptual model is discussed in the survey paper and additional problems that need further study or are unresolved are listed therein. Specifically, the organization and intensification of the pre-tornadic vorticity near the surface and its rapid development upward during Stage 3 are emphasized, together with the uncertainty of the 3D evolution of these vortex patches. The question of the role and importance of the streamwise vorticity current and the surface roughness in Stage 1 is raised, together with the role of turbulent eddies in the environmental boundary layer [2]. Other questions relate to the processes or mechanisms that determine the tornado size and intensity, both on the supercell scale (Stages 1 and 2) and on the vortex scale (Stages 3 and 4), as well as to the importance of supercell-external factors, such as storm mergers and mesoscale boundaries.

A recent computational paper by Parker [3] argues that two phenomena are important for the generation of tornado-like vortices. A careful study is performed and it is concluded that a sufficient amount of unstructured surface vertical vorticity, ζ , when coupled with a strong enough updraft, is sufficient for the formation of a tornado-like vortex. In particular, quasi-random ζ -fields are generated via the use of Perlin noise [4], and an updraft nudging technique [5] is used. It is concluded that regardless of the origin of the vertical vorticity, tornado-like vortices are produced when enough initial surface circulation and a large enough vertical velocity vertical gradient ( w/ z in the standard notation) near the ground are present. The former aspect would correspond to Stage 2 and the latter aspect to the intensification via stretching in Stage 3.

Relating to Stage 2, we will want to better understand the phenomena leading to the formation of vortex patches in the thin boundary layer near the ground. The case of three-dimensional turbulence in a shallow layer that can be approximated as two-dimensional is called quasi-two-dimensional turbulence. We use this term later in Section 3 to describe the situation studied there. In two-dimensional flows, energy E= 1 2 | v | 2 and enstrophy Ω= | ω | 2 are conserved quantities. In these formulas, v and ω are the velocity and vorticity fields, respectively, and denotes their space averages. The energy spectrum, E( k ) , and the enstrophy spectrum, Ω( k ) , are related via Ω( k )= k 2 E( k ) , where k is the wave number, which implies that in two dimensions the energy must cascade to larger scales and the enstrophy must cascade to smaller scales [6]. The extent to which this fails to hold in three dimensions is related to whether turbulence in thin surface layers is quasi-two-dimensional or not [7].

A simple laminar flow example of how surface friction can produce separation in the boundary layer is given in section 2.2 of [8]. Coherent structures in wall turbulence, specifically hairpin structures and their organization into packets and their dynamics in the boundary layer, are simulated and discussed in [9].

A different approach is studied in [10], in which a hybrid vortex filament numerical scheme is used to model a thin boundary layer within a flow past a thin horizontal plate and the subsequent development of turbulence.

The near-surface horizontal vortex filaments originally transverse to the flow develop perturbations referred to as furrows. This is shown in the top image of Figure 1 (see also Figs. 1 and 4 in [10]). The initially forward-tilting arch-like filaments over time turn into mushroom-like arrangements as seen in the bottom images of Figure 1 and in Figure 2. The rotational part of the vorticity is captured by visualizing the λ 2 -isosurfaces of rotation, shown as green “hairpin vortex legs” in Figure 1 and Figure 2. The filament segments create cusps between the hairpin legs, and the hairpin legs eventually pass through the lobes of the mushroom caps as seen in the top image of Figure 2. The vorticity within the hairpin legs is mostly aligned with the streamwise direction of the legs and is shown by the red and blue color coding in Figure 1 and Figure 2. The hairpin legs can occur as counter-rotating pairs as shown, or singly (not shown in this paper), and it is emphasized in [10] that there is much non-rotational vorticity engaged in the functioning of the vortical structures that produce hairpins and that neglecting the role of non-rotational vorticity has led to the illusion that hairpin vortices and packets are structures in their own right. These ideas are further elaborated and supported using direct numerical simulations in [11]. Consequently, cusp formation appears to be an essential mechanism leading to turbulence and possible eventual tornadogenesis, which is supported by computational evidence in [2] and shown in Figure 3. For a shallow surface layer, the behavior of the cusps can be approximated by a quasi-two-dimensional model.

Figure 1. Development of vortex filaments in the flow past a thin horizontal plate [10]. The onset of furrows and cusps, along with the green isosurfaces of rotation (top); filaments in the furrow viewed by a moving observer (bottom). Filaments closely aligned with the x-direction are colored red and blue.

Figure 2. Isosurfaces of rotation (green) and the mushroom-like cusps (top); filaments at a fixed time after mushroom-like cusps have developed and before the onset of a turbulent flow (bottom) [10].

Given the quasi-two-dimensional nature of the cusp and hairpin vortices in the thin boundary layer, they do not dissipate [7]. Instead, an inverse energy cascade forms to produce coherent (turbulent) structures that are later subject to vortex

Figure 3. Isosurfaces of the vertical velocity w=3 m∙s1 (gray) and near surface vortex lines passing through ( x,y )=( 9,12 ) at z=7.5 , 12.5, 17.5, 22.5, and 27.5 m [2]. A continuum of these vortex lines would form a corrugation of vertical vorticity as seen in figure 5 in [2] and furrows and cusps consistent with the results in [10].

stretching and may produce a tornado [12]. Due to the inverse energy cascade, the filaments concentrate into clusters of vorticity that satisfy the cyclostrophic balance conditions. The breakdown of the cyclostrophic balance between the pressure gradient and the centrifugal forces leads to the collapse of these vorticity clusters into narrow high-energy filaments.

If a quasi-two-dimensional model is used to model the behavior of cusps in a shallow boundary layer, the cusps will eventually be subject to stretching, which can then lead to storm intensification and tornado-like vortex development. In two dimensions, the concept of vortex stretching is nonexistent [13]. As stretching occurs, the cusps unravel (see Figure 3), and the quasi-two-dimensional turbulence model has to be replaced by a three-dimensional model.

The numerical simulations in [2] show that tornado development might be enhanced by turbulence and the dynamics of vortical structures in the boundary layer. Specifically, it is argued that boundary layer “coherent structures,” such as cusps and hairpins, play a significant role in tornadogenesis. The initial formation of cusps and hairpins occurs in a relatively thin boundary layer, whose thickness is estimated to be ≈6 meters [2]. When the height of the boundary layer is relatively small (see more on that in Section 3), we may view the filament model in [14] [15] and described in Section 2 as representing two-dimensional turbulence. We point out the similarity to the behavior of point vortices in two dimensions, where vorticity evolution can also be governed by other forms of the nonlinear Schrödinger equation. A recent work on verification of such analogies, particularly an analogy with the Gross-Pitaevskii equation, is presented in [16] and discussed further in Section 5. A deep mathematical analysis connecting fluid (Euler and Navier–Stokes) equations and various forms of the nonlinear Schrödinger equation is presented in [17].

The results of the numerical simulations in [2] show alternating streaks of wavelength 500 - 700 m of positive and negative vertical component of vorticity, ζ , on the warm side of the gust front (see Fig. 5 in [2]). In the positive- ζ streaks, ζ frequently reaches or exceeds 0.03 s1 at elevation 7.5 m. The ζ in the negative- ζ streaks is generally of lesser magnitude, rarely less than -0.02 s1. This compares favorably to the furrows and hairpin legs in [10] (see Figure 1 and Figure 2) and is reminiscent of the initialization of the surface vertical vorticity in [3] (see Fig. 5 in [3]). Fig. 7 in [2] shows streamwise oriented λ 2 -isosurfaces consistent with the hairpin legs in [10] and shown in Figure 1 and Figure 2. Similarly, the vortex lines shown in Figure 3 compare favorably to the solutions of the Schrödinger equation shown in Figure 4. We note that to generate the alternating vorticity streaks in [2], small random initial temperature perturbations were added and the boundary layer was given 12 hours to evolve to a quasi-steady state before storms were initiated via the introduction of a warm bubble. It may be interesting to explore the connection between this initialization process and the computational stage in [10] after the flow encounters the thin horizontal plate and before the onset of turbulence (see Figs. 1 and 2 of [10] for 0.5<x<0.8 ).

In the video “Mays’s Fury” [18], which documents the coverage of a violent tornado outbreak on May 3, 1999, a nearly horizontal vortex in or near a tornado, while the tornado crosses interstate I-35, is shown in the 41st minute. In addition, the horizontal vortex is visible in the region where the tornado appears to be the most intense. This is suggestive of a hairpin leg in or adjacent to the tornado and the associated unraveled narrow cusp (vertical filament) vortices, moving into the tornado, causing the tornado to intensify or maintain its intensity, and it is consistent with the computed vortex lines near ( x,y )=( 15,11 ) in Figure 3.

Investigation of atmospheric flows, especially those related to tornado-like events, requires a significant contribution from thermodynamic considerations, since a significant portion of a thunderstorm’s energy is generated by latent heat release associated with water vapor condensation. One possible mechanism of latent heat supplying energy to a tornado is proposed in [19]. The mechanism relies on water vapor from the cloud base being sucked into the tornado funnel due to a pressure gradient, condensing, and supplying additional buoyancy energy to the air surrounding the tornado via the released latent heat. The paper emphasizes an essential feature of this model, the central downdraft bringing the moisture, and thus energy, rich air from the cloud base into play and argues that this is the only way to bring water vapor into the tornado funnel. Citing [20]: There is a strong nexus with thermodynamics, because these thunderstorms are driven by the phase change of water vapor. There are lots and lots of things other than pure fluid dynamics in this field.

To bypass the complications of storm thermodynamics, researchers and forecasters utilize various indices, in particular, combined mechanical and thermodynamic parameters, such as the fixed layer significant tornado parameter (STP). The STP is a composite index that includes the 0 - 6 km bulk wind difference (6BWD), the 0 - 1 km storm-relative helicity (SRH1), the surface-based parcel convective available potential energy (sbCAPE), and the surface-based parcel lifting condensation level (sbLCL). It is defined as

STP= sbCAPE 1500 J kg 2000sbLCL 1000m SRH1 150 m 2 s 2 6BWD 20 m s .

The sbLCL is set to 1.0 when sbLCL < 1000 m, and it is set to 0.0 when sbLCL > 2000 m; the 6BWD term is capped at a value of 1.5 for 6BWD > 30 m/s, and set to 0.0 when 6BWD < 12.5 m/s [21].

A majority of significant tornadoes (EF2 or stronger on the Enhanced Fujita scale [22]) have been associated with STP values greater than 1, while most non-tornadic supercells have been associated with values less than 1 in a large sample of Rapid Refresh (RAP) analysis proximity soundings. The RAP is a continental-scale NOAA hourly-updated assimilation modeling system operational at the National Center for Environmental Prediction (NCEP) [23]. It covers North America and is comprised primarily of a numerical forecast model and an analysis assimilation system to initialize that model. The definitions of the updated indices are given on the NOAA website [24]. A detailed discussion of the combined and other tornado indices is given in [25].

In this paper we develop ideas from non-equilibrium thermodynamics and statistical mechanics of vortex gases in the context of a turbulent boundary layer, and we investigate how these ideas are related to the various stages of tornadogenesis outlined above and to the findings in [3]. We describe how the formation of coherent structures, such as cusps along an originally nearly straight horizontal vortex, can be modeled by the solutions of nonlinear Schrödinger equations. We will use the statistical mechanics of vortex gases to investigate the development of pre-tornadic vorticity that is first organizing and intensifying near the surface, and also the development of vortex patches and turbulent eddies in the environmental boundary layer. We use non-equilibrium thermodynamics to study the development of conditions in the boundary layer facilitating rapid vertical development of vorticity from the boundary layer to the cloud base. One of the aims of this work is to apply irreversible thermodynamics to incorporate non-equilibrium thermodynamics into the evolution of a tornado-like flow. We attempt to develop a better understanding of how non-equilibrium thermodynamics accounts for low-level boundary turbulence that leads to the generation of tornado-like vortices and their subsequent stretching. Unlike in [19], we will only consider one-constituent flows for the simplicity of developing these ideas. We plan to revisit the topic for two- or three-constituent flows later.

We consider as the main contribution of this paper the review of and bringing together ideas from previous research and applying them to better understand the various stages of tornadogenesis. It is our understanding that the connections we are making are novel and have not been made earlier. At this point, our contributions are primarily conceptual and qualitative, and we are just beginning to explore quantitative applications to real-life scenarios or tornado forecasting.

The paper is organized as follows. In Section 2 we discuss the nonlinear Schrödinger equations relevant for modeling the behavior of single vortex filaments that capture the self-induction and eventually also stretching consistent with the Navier–Stokes equations. We demonstrate that these models have the ability to model the creation of vortex cusps and are thus relevant for the early stages of tornadogenesis. In Section 3 we review the statistical mechanics of the relevant vortex gas models that allow us to better understand the quasi-two-dimensional turbulence in a shallow boundary layer. We connect these ideas to Stages 2 and 3 of the tornadogenesis process of [1]. Section 4 outlines the foundations of non-equilibrium thermodynamics, addresses the ideas of entropic balance in a tornado environment, and introduces the analysis of thermodynamic fluxes in the local entropic balance equation. We use the classical and extended irreversible thermodynamics to derive transport equations for the entropy density associated with non-equilibrium heat flux that, after the development of quasi-two-dimensional turbulence, combined with updrafts below the supercell wall cloud, will tilt and stretch the vortex lines and unravel their cusps leading to tornadogenesis. We also derive expressions for the heat and the entropy production in the boundary layer. In Section 5 we use the analogy between the filament equation and its counterpart in the nonlinear Schrödinger equation and develop it further to compare with the (mathematical but not modeling) structure of the special case of the Gross-Pitaevskii equation. We also provide references there on the equivalence of the Schrödinger and the Gross-Pitaevskii wave functions and fluid flow equations via the Madelung and Hasimoto transforms. The benefits in this case come from the fact that vortex equations use more fundamental, potential-like quantities, such as the vortex density and the curvature, rather than the original fluid flow fields (velocity and vorticity), to describe flow singularities.

Finally, Section 6 provides a summary and discussion of the main results and concludes with ideas for future investigation.

2. Asymptotic Equation for Vortex Filaments

In this section we address the behavior of single vortex filaments that could be in an originally nearly straight, horizontal position, similar to the vortex filaments prior to the onset of cusps and furrows in [10] and discussed in the introduction. We describe and highlight the relevant mathematics, focusing on the works [14] [15], with the intent to later tie the results to Stage 2 of tornadogenesis. One advantage of the approach in this section is that flow singularities can be described by focusing on properties of the vortex, in this case its curvature, rather than on the underlying velocity and vorticity fields that could suffer blow-ups that could be difficult to model [26]. While the models and equations in this section are purely mechanical, they are derived rigorously from the Navier–Stokes equations in the references provided. Based on the results arising from these models, we can conclude that vertical cusp-like structures can develop in the shallow boundary later, tying the results of this section to Stage 2 of tornadogenesis discussed in the introduction. Thermodynamical considerations will be used in later sections to address tilting, vertical stretching of vorticity, and vortex intensification.

Small-amplitude, short-wavelength perturbations of a single thin horizontal vortex filament are investigated in [14] [15] (see also [13]). This work relates to our understanding of the behavior of the initially nearly horizontal transverse vortices in [10] and captures the initial stages of the development of cusps and mushroom caps. Similar to [10], the background flow is a shear flow. In the following paragraphs, we will review the evolution of such a three-dimensional vortex filament using a model that does not allow for stretching and also show how cusps can form. Later in this section, we discuss the extension to a model that allows for stretching, also derived in [14] [15], and its implications.

Consider a narrow, three-dimensional vortex filament described by its center curve position function x( γ,t ) , where γ is the arc length curve parameter and t is the (rescaled) time variable. Let κ( γ,t ) be the curvature of the center curve and b( γ,t ) the binormal vector along the center curve of the filament (i.e., the cross product of the unit tangent and the unit normal vectors to the curve). The self-induction equation governing the evolution of the filament and consistent with the Navier-Stokes equations is then [14] [15]

x t =κ( γ,t )b( γ,t ). (1)

The validity of the self-induction Equation (1) rests on the following four assumptions, in which Γ is the circulation around the filament, ν is the kinematic viscosity, ω( x,t ) is the vorticity, and r( γ,t )=κ ( γ,t ) 1 is the radius of curvature of the center curve:

  • the vorticity ω( x,t ) is nonzero only inside the filament of radius δ1 ;

  • the radius satisfies δ= ( Γ/ν ) 1/2 ;

  • the vorticity ω( x,t ) is large enough inside the tube so that Γ/ν as ν0 ;

  • the radius of curvature r( γ,t ) is bounded strictly away from zero.

We note that Equation (1) does not allow for changes of the length of the curve, and thus the dynamics described by (1) will not capture stretching effects [13].

Consider next the Hasimoto transform

ψ( γ,t )=κ( γ,t ) e iΦ( γ,t ) ,Φ( γ,t )= 0 γ τ( y,t )dy , (2)

where ψ( γ,t ) is the complex-valued filament function and τ( γ,t ) is the torsion of the filament curve [13] [14]. Recall that torsion can be viewed as a measure of the failure of a curve to be planar. A perturbation of (1) in the form

x t =κ( γ,t )b( γ,t )+ δ ˜ v ˜ ( γ,t ), (3)

where 0 δ ˜ 1 and v ˜ ( γ,t ) is a general perturbation of the velocity, together with the Hasimoto transform (2), then leads to the perturbed Schrödinger equation for the filament function ψ

1 i ψ t = 2 ψ γ 2 + 1 2 | ψ | 2 ψ δ ˜ ( ψ, x γ , 2 x γ 2 , v ˜ γ , 2 v ˜ γ 2 ), (4)

where is a smooth function whose expression is given by formula (7.21) in [13]. Since in this section we will be interested in the limiting case δ ˜ =0 , we won’t repeat it here. For the original self-induction case with δ ˜ =0 , Equation (4) turns into the nonlinear Schrödinger equation (with cubic nonlinearity)

1 i ψ t = 2 ψ γ 2 + 1 2 | ψ | 2 ψ. (5)

Note from (2) that the magnitude of the filament function appearing in (4) and (5) is the curvature of the filament, i.e.,

| ψ( γ,t ) |=κ( γ,t ).

The cubic Schrödinger Equation (5) has solutions (including explicit ones) that allow for large values of | ψ | and thus large curvatures. An example of an explicit solution exhibiting dependence on both γ and t is [27]

ψ( γ,t )=±2A e i( Bγ+( A 2 B 2 )t+ C 1 ) sech( Aγ2ABt+ C 2 ) (6)

with arbitrary real constants A , B , C 1 , and C 2 . For this solution we have max| ψ( γ,t ) |=2| A | , showing that the curvature can be made arbitrarily large. Note that (6) imitates a cusp-like and a soliton-like solution to (5), and the visualization of | ψ( γ,t ) | for fixed values of the constants and three different time values is shown in Figure 4.

Figure 4. A visualization of the exact solution (6) of the cubic nonlinear Schroödinger Equation (5) as a function of γ . The parameter values are A=B= C 1 =1 and C 2 =0 , and three different times, t=2 , 0, and 2, are shown left to right.

This solution indicates that the cusp structures that appear in the boundary layer and are observed in the simulations mentioned above may appear as a result of turbulence in the boundary layer before stretching occurs.

An asymptotic extension of the nonlinear Schrödinger filament Equation (5) is derived in [14] [15] to allow for stretching. An originally straight filament in the unit direction t 0 and of core size δ is replaced with a perturbed version

x( γ,t )=γ t 0 + ε 2 X( γ ε , t ε 2 )+O( ε 2 ),ε1,

where X( γ ε , t ε 2 ) is a perturbation term orthogonal to t 0 and rapidly decreasing in the first variable. For additional assumptions for this model see [14]. This ansatz then leads to an extension of the filament self-induction Equation (1) that has the form

x t =κ( γ,t )b( γ,t )+ ε 2 I[ X( γ,t ) ]× t 0 ,

which, recalling the Hasimoto transform (2),

ψ( γ,t )=κ( γ,t ) e iΦ( γ,t ) ,Φ( γ,t )= 0 γ τ( y,t )dy ,

and, assuming δε1 , leads to an asymptotic filament equation with self-stretching

1 i ψ t = 2 ψ γ 2 + ε 2 ( 1 2 | ψ | 2 ψI[ ψ ] ). (7)

In the above equations, the nonlocal linear operator I[ . ] is defined as

I[ ψ( γ,t ) ] 1 | h | 3 [ ψ( γ+h,t )ψ( γ,t )h ψ γ ( γ+h,t ) + 1 2 h 2 H( 1| h | ) 2 ψ γ 2 ( γ,t ) ]dh, (8)

where H is the Heaviside function limiting the influence of 2 ψ/ γ 2 to the nearby points. The operator I[ . ] is nonlocal in the sense that it includes information from all the points along the centerline of the filament. Note that in (7) the cubic term 1 2 | ψ | 2 ψ and the I[ ψ ] term compete on the same scale.

Equation (7) models the “birth” of strong cusps in the vortex filament from the interaction of the nonlocal operator (representing a strain flow) and the cubic nonlinearity, although the assumptions of the asymptotic theory become invalid as the cusp vortex is born [15]. Specifically, the assumptions of small-amplitude, short-wavelength perturbations of a thin straight vortex are violated just before cusps or any shapes resembling mushroom caps evolve. In the numerical calculations reported in [15] the vortex filament develops higher and much narrower peaks as time evolves when compared with the corresponding solutions of the cubic nonlinear Schrödinger Equation (5). The details of the evolution up to that time are fully discussed in [15]. The time sequence of the behavior of an initially helically perturbed straight filament was computed in [15] and is shown in Figure 5. The initially helical filament is shown through the middle of the figure, and three stages in the initial development of cusps along the filament are shown at times t=1.85 , 1.90, and 1.95, showing the formation of a cusp near where the arrows in Figure 5 are pointing.

At the end of Section 4, we give an example that incorporates the boundary layer shear that would induce hairpins and furrows of cusps, and the tilting of the transverse vortices with cusps, induced by a pressure deficit below a mesocyclone.

Figure 5. Time sequence of the filament curve showing the initial helical perturbation of a straight filament at time t=0 (thin black curve) and three later stages at t=1.85 (red dots), 1.90 (green dots), and 1.95 (yellow dots), showing the development of cusps along the filament (adapted from [15]).

Vortex filaments with vertical segments are important for the formation of vortex patches in Stage 2 of our understanding of tornadogenesis discussed in the introduction. Stage 3 starts with the organization of vortex patches into a single symmetric vortex that later intensifies via stretching. The following section reviews the ideas from statistical mechanics of vortex gases and ties them to the aforementioned organization and single vortex generation. It does not deal with any thermodynamical concepts, which will be discussed later in Section 4.

3. Statistical Mechanics of Vortex Gases

In statistical mechanics of gases, solving for the behavior of a large number of particles (gas molecules) is simplified by making assumptions about the mean behavior of the collection, rather than attempting to solve the intractable system of equations for the interactions of each particle with all the others. More specifically, a phase space containing information about the position and momentum of each particle is considered, and a probability distribution on the phase space is obtained. In the case of fluid dynamics, the phase space of vortex gases depends on the set up of the vortex gas, which could be point vortices in 2D (see, e.g., sections 6.2 and 7.3 of [13]), parallel or nearly parallel vortices in 3D (section 7.4 of [13]), or a single vortex filament [28]-[30]. Parallel vertical vortices are analogous to a corresponding set of point vortices in the plane in that they satisfy the same equations [13].

One has to be careful when working with concepts that have the same, or similar, nomenclature both in thermodynamics and in statistical mechanics. For example, the energy we use in this section is the kinetic energy of the ensemble in the usual sense, but the statistical mechanics entropy, to be used later in this section, and the thermodynamics entropy, used in Section 4, do not generally agree with each other [31]; similarly for temperature. This section will only deal with statistical mechanics concepts, so “entropy” will refer to the statistical mechanics entropy and “temperature” will refer to the statistical mechanics temperature. The corresponding thermodynamics concepts will only be used in Section 4.

In the case of an ideal gas, let the microstates be indexed by i . If we denote the entropy by S , the energy of microstate i by E i , and the total energy of the system by E , then the partition function corresponding to this scenario is

Z= i e β E i , (9)

where β= S/ E =1/ T s represents the inverse temperature (we will use T for the absolute thermodynamics temperature in Section 4). A similar approach, possibly with the sum in the partition function (9) replaced by an integral of e βE over all configurations of vortices under consideration, can be used in fluid dynamics.

From probabilistic considerations, when β=0 (or T s =± ), all energy levels are equally likely to occur, and the entropy S is at its highest point. For β>0 , configurations with lower energy are more likely and these correspond to more balled up vortices. Finally, for β<0 , configurations with higher energy are more likely and these correspond to straighter vortices. See [28]-[30] for details. Therefore, negative temperature vortices are hotter than positive temperature vortices, with vortices with T s =± ( β=0 ) temperatures in between them.

Using a different viewpoint, we can say that negative-temperature vortices have higher energy and are (more) straight or smooth, while positive-temperature vortices are (more) balled up and, from a physical point of view, dissipating. Therefore, energy is expected to flow in the direction of increasing β : straight, high energy vortices will become less straight and begin to fold up as β<0 increases towards 0. Such vortices will transfer their energy to larger scales. As β increases into positive values, the configurations become more balled up, starting the Kolmogorov energy cascade to smaller scales and dissipating [30] [32]. We note that for negative temperatures to occur, the phase space must be bounded, as well as the energies E i [33].

Let’s next consider a collection of nearly parallel vortices considered in [34] with their principal direction t 0 assumed to be vertical. This approach can be used to model the various nearly vertical vortices in a shallow layer. For future considerations, and, to use the theoretical results cited below, we assume that this shallow layer is situated above a bounded region on the ground and beneath a supercell. Similar to the description in Section 2, the general center curve of each of the filaments, x i , is assumed to have the asymptotic form

x i ( γ,t )=γ t 0 + ε 2 X i ( γ ε , t ε 4 )+O( ε 2 ),ε1,

where X i ( γ ε , t ε 4 ) is again a perturbation term orthogonal to t 0 . For additional assumptions for this model, similar to those in Section 2, see [34]. We point out that the separation distance between filaments has to be much larger than the filament core size. It is then possible to derive asymptotic limiting equations for N nearly parallel interacting vortex filaments

x i t =J[ α i Γ i 2 x i γ 2 ]+J[ ji N 2 Γ j x i x j | x i x j | 2 ],1iN, (10)

where J is the skew-symmetric matrix

J=( 0 1 1 0 ).

Here, Γ i is the circulation of the ith vortex and α i is a constant determined by the vortex core structure (see additional references in [34]). Using the identification ψ j ( γ,t )= x j ( γ,t )+i y j ( γ,t ) , where x j =( x j , y j ,γ ) , Equation (10) can be rewritten as coupled nonlinear Schrödinger equations. A special example illustrating their structure is for a pair of two interacting vortices with identical core parameters α 1 = α 2 =1 and circulations Γ 1 =1 and Γ 2 =Γ . In that case, making a change of variables ψ= ψ 1 ψ 2 and ϕ= ψ 1 + ψ 2 , system (10) becomes

1 i ϕ t = 1 2 ( 1+Γ ) 2 ϕ γ 2 + 1 2 ( 1Γ )[ 2 ψ γ 2 4 ψ | ψ | 2 ], 1 i ψ t = 1 2 ( 1Γ ) 2 ψ γ 2 1 2 ( 1+Γ )[ 2 ψ γ 2 4 ψ | ψ | 2 ], (11)

which interestingly becomes decoupled for the case of co-rotating filaments of equal strength ( Γ=1 ). Numerical evidence is presented in [34] for a finite-time collapse of a vortex pair with a negative circulation ratio Γ in the sense that the filaments would eventually come close to each other and violate the assumption on the separation distance between the vortices. On the other hand, numerical evidence is given that no finite-time collapse happens for a positive circulation ratio. Notice in Figure 1 and Figure 2 how the legs of the cusps, originally at a non-zero distance from each other at times t<0.7 appear to have merged after time t=0.8 , consistent with the first scenario. Also notice that the nonlinearity in (11) is not cubic, making these equations qualitatively different from those for a single vortex filament (Equations (5) and (7)).

The shallow boundary layer with cusps below the supercell has filaments in nearly vertical directions and therefore the nearly parallel vortex gas model should be applicable. A further simplification and insight come from assuming that the vortex filaments are actually vertical and parallel to each other. In that case, the second-derivative with respect to γ terms in (10) and (11) are equal to 0, and the equations become identical to those for point vortices in two dimensions (see section 7.3 of [13]). Hence, the term quasi-two-dimensional turbulence is used to describe three-dimensional turbulence in a shallow surface layer when the vortices can be approximated as vertical and straight. In view of the earlier discussion, such vortices, being perfectly straight, can be considered to have a negative statistical mechanics temperature.

We now review the results of the two-dimensional model and connect them to the formation of the vortex patches. The two-dimensional mean field theory was developed in [35] [36] and allows for a range of negative temperatures. This theory motivated a three-dimensional theory [37], which does not include negative temperatures; however, an extension of this work allows for a small range of negative temperatures provided the 3D vortex has a fractal cross section [38].

An N -vortex system in a bounded domain Λ in 2 with the associated canonical Gibbs measure at inverse temperature β ˜ is studied in [35] [36]. The vortices are all assumed to have the same intensity α (we are using the notation of [35], not to be confused with the 3D vortex core parameters α i earlier in this section). It is shown that as N , β ˜ /N β( 8π,+ ) , and αN1 , the stream function ψ associated with this problem satisfies the mean field equation

Δψ= e βψ Λ e βψ ,ψ=0onΛ. (12)

The limiting vorticity field then satisfies ω β =Δψ0 and Λ ω β =1 . For a disk domain Λ, the solutions of the mean field Equation (12) can be found and satisfy ω β δ x 0 as β8 π + , where x 0 is the center of the disk and δ x 0 ( x ) is the Dirac delta function. It can also be shown that this kind of vorticity concentration also exists in a simply connected domain sufficiently close to a disk and possibly also in other domains. However, there are domains for which this concentration does not occur (e.g., rectangles with large aspect ratios) [36]. Convergence to a smooth solution is also possible [35] [36]. If we were to consider the extension of this 2D model back to the case of co-rotating vertical parallel vortex filaments, the scenario in which vorticity tends to concentrate near a single point can be connected to the organization of vertical vorticity into a single symmetric intense vertical vortex that would later undergo stretching and intensify further.

Applying the previous considerations to the formation of furrow-like structures with cusps and adjacent hairpin legs present in the numerical simulations [2] [10], we let Λ be a narrow, bounded region with a smooth boundary below the furrow between the hairpin legs. In the presence of many furrows, each such vortex patch would be modeled by a corresponding region Λ. Figure 2 and Figure 3 showing the vertical vorticity cusps suggest that the total vorticity in the quasi-two-dimensional turbulence patches (regions Λ) should be zero, provided the positive and negative vorticity legs have the same strength. However, a 12-hour burn-in that included Coriolis effects (planetary vorticity) was used in the simulations of [2], which led to the positive vorticity leg having a larger vorticity value than the negative vorticity leg, or, in the language of [2], the positive- ζ anomalies in the boundary layer were stronger than the negative- ζ anomalies. Consequently, the total vertical vorticity induced by the cusps in the patches is expected to be positive.

If we think of the vertical vorticity cusps as straight vertical filaments and the turbulent boundary layer as made up of these filaments, we can use the two-dimensional turbulence model of [35] [36] discussed above. With the filaments being straight, we can assume that they have negative temperature and cluster into larger filaments with cyclonic and anti-cyclonic legs (the right image of Figure 2). As these larger clusters in quasi-two-dimensional turbulence are tugged up into the mesocyclone (in a streamwise sense), the negative side of the vertical vorticity cusp tilts 180 degrees, aligns with the positive side, and the vertical vorticity component becomes positive (see Figure 3). As this happens and the vortex stretches, a transition from a quasi-two-dimensional model to a three-dimensional model is needed.

Initially, the resulting vortex is in cyclostrophic balance, i.e., the pressure gradient force and the centrifugal force are in balance. Surface friction along with the stretching and rotation of the vortex causes the centrifugal force to weaken and allows the pressure gradient force to dominate. Consequently, the vortex collapses to a narrow filament with high energy density. This organization in the flow results in a decrease in entropy density, an increase in energy density, and a negative temperature of the vortex since ΔS/ ΔE <0 [30]. Such vortices are called supercritical or suction vortices and can be approximated by a delta function δ x 0 ( x ) supported on a vertical segment at the horizontal location x 0 . As discussed above, these supercritical vortices may then contribute to tornadogenesis. Further stretching of a negative temperature vortex would lead to its transferring energy to larger scales via folding, and eventual dissipation of energy via the Kolmogorov cascade, a direct energy cascade to smaller and smaller scales without a lower bound [30]. The processes described above may be occurring repeatedly in the region below the supercell and mesocyclone. The aforementioned example at the end of Section 4 will address the stretching leading to the transition to supercritical vortices.

4. Non-Equilibrium Thermodynamics and Entropy

In this section we outline the foundations of non-equilibrium thermodynamics and address entropic balance under highly non-equilibrium settings. Unlike in the previous section, entropy and temperature will now refer to their thermodynamic concepts. In subsection 4.1 we list expressions for the entropy flux, the heat flux, and the internal entropy production in classical irreversible thermodynamics and derive them for the case of extended irreversible thermodynamics. We also discuss entropic balance in the boundary layer that is instrumental in stretching and tilting of vorticity into the vertical direction, which is a part of Stages 3 and 4 of tornadogenesis described in the introduction. We provide the entropy balance equation (17) in the classical form and derive an extended form (25) that includes an internal entropy production term. In subsection 4.2 we utilize a connection between the potential temperature gradient and the local vorticity, and discuss non-equilibrium thermodynamics of the boundary layer including the derivation of the heat and the entropy production in the boundary layer.

We should note that in recent years studies of “standard” classical thermodynamics have developed in many branches influenced by non-equilibrium models in phenomenological and statistical formulations. In this paper, we are carefully specifying the thermodynamical context we use. For instance, this section applies “traditional” irreversible thermodynamics in subsection 4.1 and the recent developments of the much less traditional theory of entropic balance in subsection 4.2. In contrast, preceding Section 3 used statistical mechanics of vortex gases and these two approaches are independent of one another.

4.1. Entropy in Irreversible Non-Equilibrium Thermodynamics

This subsection follows the standard development of the irreversible thermodynamic theory [39]. We begin with classical irreversible thermodynamics. The first step to modeling a system that is not in equilibrium is to employ the local equilibrium hypothesis, where all the relevant thermodynamic variables are defined locally and are in a quasi-equilibrium for a “small enough” system. In the fluid dynamics context, such a “small enough” system is an air parcel. We let s( x,t ) denote the specific thermodynamic entropy as a function of position, x , and time, t , and assume also that it is a function of the specific variables a j ( x,t ) , whose corresponding conjugate variables are denoted by F j ( x,t ) . The Gibbs’ equation

Ds( x,t )= j F j ( x,t )D a j ( x,t ) (13)

gives the relationship between these variables and leads to the differential equation

Ds( x,t ) Dt = j F j ( x,t ) D a j ( x,t ) Dt ,

where

D Dt = t +v

is the material derivative, i.e., the rate of the change of a quantity following the trajectory of the air parcel. The balance relation for the specific entropy is then

ρ  Ds Dt =B+σ, (14)

where ρ( x,t ) is the mass density, B( x,t ) is the entropy flux, and σ( x,t )0 is the entropy production rate per unit volume [39] [40]. The entropy production rate is assumed to be nonnegative for consistency with the second law of thermodynamics and, in general, incorporates the following sources [40]:

1) entropy production associated with heat transfer;

2) entropy production due to mass transfer;

3) entropy production because of viscous dissipation of fluid;

4) entropy production arising from chemical reactions.

For simplicity, we assume that all the thermodynamic functions involved are smooth (for a more general case with discontinuities see [41]) and our discussion will focus solely on the entropy production associated with heat transfer; however, we acknowledge that it would also be important to consider the entropy production resulting from mass transfer (moisture flux), which will be left for future work. Transport equations similar to (14) are valid for all the state variables a j . We use the basic thermodynamic model with two local state variables, the specific internal energy, u( x,t ) , and the specific volume, v( x,t )= ρ 1 ( x,t ) , so that the equation of state is s=s( u,v ) . For the purposes of this discussion, we assume that both state variables u and v are known functions of x and t . Consequently, we have two conjugate variables here,

s u = T 1 and s v =p T 1 , (15)

where T( x,t ) is the absolute thermodynamic temperature (later also referred to as the local equilibrium temperature) and p( x,t ) is the pressure.

In the case of a laminar motion of a one-constituent isotropic, incompressible ( v=0 ), and viscous fluid in the presence of a temperature gradient, without internal energy supply, we have

ρ  Du Dt =qP:Vandρ Dv Dt =0, (16)

where q is the heat flux, v is the velocity, V= 1 2 ( v+ ( v ) T ) , and P is the symmetric pressure tensor, all of these being functions of x and t . We recall that A:B=tr( A T B )= A ij B ij using the Einstein summation notation here. If we decompose the pressure tensor

P=( p+ p v )I+ P 0 ,

where p is the hydrostatic pressure and p v = 1 3 tr( PpI ) so that P 0 is traceless, then the balance Equation (16) can be substituted into

Ds Dt = s u Du Dt + s v Dv Dt

to obtain

ρ  Ds Dt =( T 1 q )+q T 1 T 1 ( P 0 : V 0 ), (17)

where V 0 =V 1 3 ( trV )I=V is the traceless part of V , coinciding with V in the case of an incompressible flow ( 1 3 trV=v=0 ). In the remainder of this section, we will therefore use V in place of V 0 . Comparing this expression with (14) and arguing as in [39] that since σ represents the rate of entropy production inside the body, its expression cannot contain a flux term like ( B T 1 q ) , we obtain the equations

B( x,t )=( T 1 ( x,t )q( x,t ) ) (18)

and

σ( x,t )=q( x,t ) T 1 ( x,t ) T 1 ( x,t )( P 0 ( x,t ):V( x,t ) ). (19)

Note that (18) can be interpreted as B= T 1 q+f , where f=0 . Since there is no constitutive formula for internal entropy production, the classical irreversible thermodynamics postulates the existence of a linear relationship between the thermodynamic fluxes in (19) and their respective conjugate “thermodynamic forces.” Specifically, we have the following phenomenological relationships [39]

q=λT=λ T 2 T 1 and P 0 =2ηV, (20)

where the first relationship is just the Fourier law with the heat conductivity λ and η is the dynamic shear viscosity. Using (20), we can rewrite (19) as

σ( x,t )=λ T 2 ( x,t ) | T( x,t ) | 2 +2η T 1 ( x,t )( V( x,t ):V( x,t ) ), (21)

where | | stands for a vector norm, and it follows that the entropy transport Equation (17) can be written, suppressing the explicit dependence on x and t , as

ρ  Ds Dt =( λ T 1 T )+λ T 2 | T | 2 +2η T 1 ( V:V ). (22)

To go beyond the local equilibrium hypothesis and to extend the relationships like the Fourier law, one can turn to extended irreversible thermodynamics. In a highly non-equilibrium system, heat flux may not be well approximated by the Fourier law, as a non-equilibrium system may have spatial and temporal non-local influence. In this approach, dissipative fluxes are elevated to the status of state variables. If we, for example, extend the specific entropy to depend on the heat flux, i.e., s=s( u,v,q ) , then

Ds Dt = s u Du Dt + s v Dv Dt + s q Dq Dt = T 1 Du Dt + s q Dq Dt , (23)

where we assume that the non-equilibrium temperature ( s u ) 1 is well approximated by the local equilibrium temperature T from (15) (see Chapter 7 in [39]), and where we again used the consequence of incompressibility that ρ  Dv Dt =0 as in (16). The remaining partial derivative is assumed to satisfy [39]

s q = T 1 vαq, (24)

where α is a function of u and v only, and thus we obtain a new transport equation for the extended irreversible case based on the fundamental Equation (14) in the form

ρ  Ds Dt =( T 1 q )+q( T 1 T 1 α  Dq Dt )+2η T 1 ( V:V ). (25)

Additionally, if we use the usual linearity assumption between the heat flux q and its corresponding thermodynamic force, corresponding to the term in parentheses in (25), we can write

T 1 T 1 α  Dq Dt =μq,

which, after a short manipulation, gives rise to a Cattaneo-type relation

τ R   Dq Dt =( q+λT ), (26)

where τ R =αλT denotes the heat flux relaxation time and μ= λ 1 T 2 . Equation (26) can be substituted into (25) to remove the dependence of the right-hand side on the material derivative of the heat flux to obtain

ρ  Ds Dt =( T 1 q )+ λ 1 T 2 | q | 2 +2η T 1 ( V:V ), (27)

where the first term in the right-hand side of (27) corresponds to the entropy flux and the last two terms to the entropy production σ0 . Note that if the left-hand side of (26) is negligible, then (27) will reduce to (22). In the present case, the dependence on the heat flux q cannot be suppressed because an explicit solution to (26) is not readily available.

If the heat flux is concentrated in the boundary layer, where the pressure changes (on average) are insignificant, we can use the Fourier law in (20) and the divergence theorem to estimate the heat supply Q to a region Ω (a thin 3D region with a piecewise smooth boundary) of the quasi-two-dimensional turbulent flow as

Q= Ω q( x,t )n( x )dxdt = Ω λΔT( x,t )dxdt , (28)

where n( x ) is the outer unit normal vector and we assume that λ is a constant. Using (18) and (20), the entropy supply to Ω, denoted by Δ S 1 , could be calculated as

Δ S 1 = Ω B( x,t )n( x )dxdt = Ω B( x,t )dxdt = Ω ( T 1 ( x,t )q( x,t ) )dxdt = Ω ( λ T 1 ( x,t )T( x,t ) )dxdt = Ω λ( ( T 1 )( x,t )T( x,t )+ T 1 ( x,t )ΔT( x,t ) )dxdt = Ω λ( T 2 ( x,t ) | T( x,t ) | 2 + T 1 ( x,t )ΔT( x,t ) )dxdt , (29)

where the first term in the last expression, which is clearly nonpositive, is due to the heat flux, and the second term is due to the heat diffusion. Finally, using (21), the entropy production in Ω, denoted by Δ S 2 , can be calculated as

Δ S 2 = Ω σ( x,t )dxdt = Ω λ T 2 ( x,t ) | T( x,t ) | 2 +2η T 1 ( x,t )( V( x,t ):V( x,t ) )dxdt , (30)

which is clearly nonnegative in agreement with the second law of thermodynamics. Notice that Equation (30) implies that the only flows with zero rate of production of entropy would be those that have constant temperature T and are a solid core rotation about an axis. We also note that expressions corresponding to (28)-(30) could also be obtained in the extended irreversible scenario, but they would need to be written in terms of the heat flux q .

The new theoretical expressions we obtained above, for the entropy flowing into the system Δ S 1 in (29) and the entropy produced inside the system Δ S 2 in (30), can be compared to experimental data from storm temperature and entropy profiles such as [42], which will be one of the directions of future work.

A detailed account of recent developments in the understanding of irreversible entropy production is given in [43]. Note that our new expression for the evolution of entropy far from the equilibrium (25) includes a contribution due to heat fluxes that, in turn, will affect vorticity development in Stage 3 described in the introduction. In the next section we will study the spatial gradient of the specific entropy, s , and its impact on the tilting, stretching, and amplification of vorticity associated with the quasi-two-dimensional turbulent motion leading to tornadogenesis.

4.2. Mesoscopic Thermodynamics, Entropic Balance, and Vorticity

The variations in thermodynamic entropy s( x,t ) are responsible for the “mismatch” between macroscopic and microscopic descriptions, since they take place on different space and time scales. So, turbulent thermodynamic systems may be categorized as mesoscopic [39] [44]. Mismatches due to fluctuations in the velocity field are called turbulent fluxes [45]-[47]. Assuming that conduction is slower than convection and assuming that the atmosphere consists of dry air and water vapor, thus neglecting the liquid phase, implies that convection dominates the entropy production [40] [42] [48]. Tornado-like vortices, from radar observations and numerical simulations [2] [12] [30], exhibit exactly the behavior of a system far from equilibrium, so they can be treated as mesoscopic systems; in other words, they exhibit both macroscopic and microscopic features. It is important then to understand how far from equilibrium that system is.

The standard approach to connecting the equations of motion and thermodynamics uses buoyancy in the Boussinesq approximation [46]. However, the Boussinesq approximation imposes restrictions on horizontal density variations. In the general case, the connection between the dynamics of motion and the flow thermodynamics can be seen from the fact that the quantity 1 ρ θω must remain constant along the parcel trajectories for an inviscid isentropic flow, in which vorticity ω is tied to the potential temperature θ [49]. The relationship between different definitions of potential temperature and entropy is given in [48].

In [50], a theory of entropic balance is developed based on variational principles under the assumption of short time scales for phase changes. The governing equations of fluid flow are derived using a variational approach in the thermodynamic context, and a discussion of how entropic balance can be affected by a highly non-equilibrium storm environment is included. In particular, it is stated that the entropic balance must be complemented with a term accounting for non-equilibrium entropy fluctuations, which would be similar to the term described in (25). A comment is included on fluctuations not being averaged on all scales, thus leading to a mesoscopic model of the non-equilibrium system (see Equation (5.6b) in [50]).

More specifically, it is hypothesized in [50] that the “phase-change time scale is significantly shorter than the time scales of convective storms and tornadoes” and that “variations of the initial entropy levels are small enough and allow us to approximate them by their ensemble means”. The second of these is the local equilibrium hypothesis discussed in Section 4. Consequently, this entropic balance theory is a special detailed case of non-equilibrium thermodynamics and thus is fundamental in our study of the boundary layer. In [50], regions are identified as entropic sources and entropic sinks within the entire supercell storm. In our setting, we consider entropic sources and entropic sinks in the context of the boundary layer. The entropic sink is the evaporatively cooled boundary layer in the forward flank or rear flank of the supercell. The entropic source is the boundary layer, laden with moisture from the northward advection of water vapor from the Gulf of Mexico.

The region of the two indexed flows can be thought of as a box below the wall cloud, just on the cool side of the forward-flank downdraft. In that region the laws of non-equilibrium thermodynamics apply. We assume that the construction in [50] is applicable in the region of the forward-flank downdraft, not the entire supercell. The close proximity of the entropic source and sink determines the size of the gradient of the entropy. Using a variational argument, entropy considerations, the Clebsh transformation, and the definitions and discussion in [51] [52], it follows that

v= L m + L s s+ L X X,ω= L s ×s+ L X ×X, (31)

where L m , L s and L X are the Lagrange multipliers corresponding to the conservation of mass, entropy, and air parcel identifier X , respectively. We note that the decomposition in (31) is not unique, but it significantly complements the discussion of [49] on the nexus between vorticity and thermodynamics. We can rewrite the second equation in (31) as a vorticity decomposition in the form

ω= L s ×s+ L X ×X= ω bc + ω bt , (32)

where ω bc and ω bt are baroclinic and barotropic vorticities, respectively [52].

The boundary layer also undergoes strong wind shear, which leads to the development of a quasi-two-dimensional turbulent boundary layer. The vorticity in the moisture-rich boundary layer creates potential vorticity patches [49] [46]. The potential vorticity ω * is defined [51] [52] as

ω * = ωs ρ , (33)

and it is conserved for inviscid isentropic flows. Consequently, using (32), potential vorticity would be zero in the absence of barotropic vorticity ω bt , i.e., ω * =0 .1 For a hydrostatic flow and neglecting the Coriolis effect, we have [49]

ω * ζ θ   θ p , (34)

where ζ θ is the third component of vorticity in the isentropic coordinates ( x,y,θ ) , and θ is again the potential temperature.

The surface heat flux into the quasi-two-dimensional turbulent boundary layer, combined with buoyancy-generated pressure deficits, produces updrafts below the supercell wall cloud [3]. If tilting and stretching unravel the cusps generated by vortex lines, then such an environment can lead to tornadogenesis.

We now consider a near-surface (quasi-two-dimensional turbulence) example using the tilting terms in the vertical vorticity equation (see, e.g., Equation (12) on p. 341 of [49]). Assume a vertical shear of 20 m/s in 100 m in the westerly direction, leading to rotation about a south-to-north axis. The wind shear creates the vortex lines with cusps in the quasi-two-dimensional boundary layer. In addition, assume a change in vertical velocity from 0 m/s to 2 m/s in 100 m in the southerly direction. Then the corresponding rate of change of the vertical vorticity is 4 × 103 s2. At this rate of vorticity production, it would take approximately 80 s for a Rankine-combined tornado of 100-meter radius, vorticity of 0.32 s1, and velocity of 32 m/s to form. This is in general agreement with the results of Houser [53], where studies of seven cases of tornadogenesis show that tornadogenesis takes place in approximately 60 s. This example ignores vortex breakdown in the tornado formation. Taking 30 m/s as the velocity of the subcritical vortex, then, after the vortex breakdown, the supercritical vortex would have a maximum velocity of approximately v max =1.7×30m/s =51m/s [54]. The tilting of horizontal vorticity in the boundary layer creates vertical vorticity and potential vorticity. This, along with the potential vorticity patches in the quasi-two-dimensional turbulent boundary layer, contributes to the stretching in the vertical direction, leading to a vortex breakdown. Indeed, if the time scale is very small and thus the adiabatic assumption is reasonable, then (34) indicates that making θ/ p small and assuming potential vorticity constant, the vertical vorticity will be large.

We note that the ideas discussed in this paper apply to more general scenarios than the one discussed in the above example. For example, if the relative direction of the supercell storm motion and the surface flow vorticity production are in the opposite directions, as is the case of a storm moving east along an east-to-west oriented warm front, then the east-to-west oriented furrows of cusps would correspond to south-to-north vortex lines. As the updraft moves along a furrow, capturing successive segments of cusps, the associated vortex lines would unravel and go through vortex breakdown as in the previous paragraph, then be successively absorbed by the tornado and amplify its vorticity.

5. Considerations of the Gross-Pitaevskii Equation

We conclude the main ideas of this paper by considering the Gross-Pitaevskii equation, written in the dimensionless form [17] [55],

i  Ψ t =ΔΨ+β | Ψ( x,t ) | 2 Ψ+V( x )Ψ, (35)

where the unknown Ψ( x,t ) typically represents a wave function, V( x ) is an external potential, and β>0 . Notice the similarity to the cubic nonlinear Schrödinger Equation (7), particularly the cubic nonlinearity. Note also that the term I[ ψ( γ,t ) ] in (7) and defined in (8) is replaced by a simpler linear term V( x )Ψ . Equation (35) has several properties that may be useful in studying tornadogenesis. For example, it has been shown that its statistics include the −5/3 power law for the inverse energy cascade in two-dimensional turbulent flows [16]. Hence (35) can be used to model two-dimensional turbulence and we can try to use it to model quasi-two-dimensional turbulence in the boundary layer. Additionally, a particular modification of the Rankine-combined vortex is a solution to (35) [56]. Finally, we point out that we are not pursuing a quantum-mechanical analogy and do not claim any quantum mechanical effects. Instead, we focus only on the mechanical properties of (35) and its solutions.

If we rewrite the wave function Ψ in terms of its magnitude and phase ϕ ,

Ψ( x,t )=| Ψ( x,t ) | e iϕ( x,t ) ,

and use the identification

ρ( x,t )= | Ψ( x,t ) | 2 andv( x,t )=ϕ( x,t ),

then the Madelung transform of (35) leads to a continuity equation and a modified Euler equation with ρ and v playing the usual roles of mass density and velocity, respectively [56] [57]. In one spatial dimension, the Madelung transform agrees with the Hasimoto transform (2) [17]. In fact, the Madelung and Hasimoto transforms establish an equivalence of the Gross-Pitaevskii and the nonlinear Schrödinger equation classes to Newton’s and fluid flow equations including the Euler and the Navier–Stokes equations [17], and many features of turbulent vortices have been obtained using numerical simulations of the solutions of the Gross-Pitaevskii equation [56]. The hierarchy of the cubic nonlinear Schrödinger equations, including (35), is studied in [58]. As seen above, solutions of (35) can be shown, after an application of the Madelung transform, to satisfy the Euler equations [57], and, vice versa, the Euler equations can be turned into (35) using the inverse Madelung transform [59]. This equivalence between the Euler equations and the Gross-Pitaevskii equation indicates that (similarly to densely packed sources of rotation perpetrated by the Bose-Einstein statistical distribution) the turbulent interaction includes a “discrete-like” component [17] [56].

Using the usual approach of the separation of variables, stationary solutions of (35) of the form

Ψ( x,t )=Z( x ) e iμt (36)

can be sought. Substituting (36) into (35), we obtain the problem

μZ( x )=ΔZ( x )+β | Z( x ) | 2 Z( x )+V( x )Z( x ),

which is to be solved for the eigenvalue μ and the eigenfunction Z( x ) , which then fully determine Ψ( x,t ) [55].

In light of the above connections between the Schrödinger equations and the Gross-Pitaevskii equation, it appears plausible that applications of the variants of (35) could possibly lead to broadening the type of vortex patches in our consideration, i.e., the narrow furrows of the pre-turbulent flow, to larger patches in a more fully developed quasi-two-dimensional turbulence [16] [56] [60].

6. Summary

In this paper, we integrate in a novel and original way the mathematical, statistical mechanics, and thermodynamical concepts that allow us to better understand the development of a three-dimensional vortex from an originally mostly non-turbulent flow. These concepts come from a variety of different disciplines. From numerical simulations and turbulence theory, one observes that surface friction and vertical wind shear lead to the development of horizontal vortex filaments, which can later develop into more vertically coherent structures such as cusps and furrows, exhibiting apparent hairpin legs [2] [10] [11].

Since these coherent structures occur in a shallow region near the surface (see Figure 2), we suggest that their behavior can be modeled and understood using quasi-two-dimensional models. In two-dimensional and therefore also quasi-two-dimensional models, energy cascades to larger scales, further supporting the development of structures or patches discussed in [1]. Then such structures are pulled into the updraft of the supercell flow [3], get further tilted, unraveled, and stretched in the vertical direction [2] (see Figure 3), and a three-dimensional model is required. When the cyclostrophic balance between the pressure gradient force and the centrifugal force is disrupted due to the surface friction, these stretched filaments collapse into narrow, very intense filaments. In Section 3 we interpreted these vortices as having a negative temperature, and such filaments, further stretched by the updraft, then transfer energy to the surrounding flow, possibly leading to or supporting tornadogenesis [28]-[30]. This paper therefore combines ideas from tornado theory, turbulence, and two- and three-dimensional statistical mechanics of vortex gases.

We integrate mathematical fluid dynamics into the discussion of the quasi-two-dimensional boundary layer. Mathematical fluid dynamics provides an explanation of the development of cusp-like features from horizontal vortex filaments via nonlinear Schrödinger equations such as (5) and (7). In particular, the nonlocal interaction term I[ ψ ] in (7) is shown in [15] to be responsible for the development of cusps in an initially straight, helically perturbed, vortex filament. These models are rigorously derived in [14] [15] under appropriate and reasonable mathematical assumptions. The models’ assumptions cease to be satisfied as the cusps are formed; however, when such cusps are formed, their behavior can be approximated by using models appropriate for collections of vortices, such as a two-dimensional vortex gas model or a quasi-two-dimensional model for parallel vertical vortices. In Section 3 we use these models and the results of [35] [36] to model vortex patches in the quasi-two-dimensional boundary layer.

The results of the numerical simulations in [2] appear to be consistent with those in [10] in that λ 2 -isosurfaces in the boundary layer of the environmental air mass align with the direction of the flow, and transverse vortex filaments develop cusps and potentially more complicated shapes between these “hairpin legs.” It is important to point out that the horizontal hairpin legs do not get tilted into the vertical. Instead, the transverse vortex filaments with cusps appear to be tilted by a combination of potential vorticity in the boundary layer and the pressure drop in the mesocyclone. These filaments then get pulled into the updraft and appear to contribute significantly to tornadogenesis. The work in this paper therefore integrates concepts from a variety of different areas to contribute to the explanation of tornadogenesis.

We have derived the entropy transport equation in a non-equilibrium context based on irreversible thermodynamics and various ideas on entropy redistribution based on calculus of variations studied in the works of Davies-Jones and Sasaki [50] [52], in which regions of entropic sources and sinks are identified and studied in the context of the entire supercell storm. The interaction of these sources and sinks leads to a thermodynamic non-equilibrium in the supercell. In the context of tornadogenesis, this theory can be applied to the boundary layer. The entropic sink is the evaporatively cooled boundary layer in the forward flank or rear flank of the supercell. The entropic source is the boundary layer, laden with moisture from northward advection of water vapor from the Gulf of Mexico and the absorption of water vapor from crops undergoing evapotranspiration. The close proximity of the entropic source and the entropic sink determines the size of the entropy gradient, and their interaction, together with other factors such as strong wind shear or surface heat flux, leads to the (potential) vorticity in this quasi-two-dimensional layer that appears to be essential for tornadogenesis. Combined with updrafts below the supercell wall cloud, this mechanism leads to tilting, unraveling, and stretching of the cusps along vortex lines in this region, contributing to tornadogenesis.

Finally, we list some questions to consider for future investigation. What are the thresholds for surface winds to produce quasi-two-dimensional turbulence with robust vortex patches? Prior to tornadogenesis, how does the near-storm thermodynamic environment affect how close the vortex patches get to the triple point, i.e., the point where the surface air from the forward-flank and the rear-flank downdrafts join with the surface environmental air? What roles do the ground flux and the moisture flux in the non-equilibrium thermodynamics play in tornadogenesis? Is there a way to continuously perturb the linear local term into a nonlocal term in the cubic nonlinear Schrödinger equation in such a way that the evolution of the solutions reflects the increasing turbulence? If there is a homotopy of evolution from the Schrödinger to the Gross-Pitaevskii equation mentioned above, can we then find entropy based purely on a time-dependent wave function without using a quantum analogy? How would the results of numerical simulations change if the classical or extended irreversible thermodynamics equations derived above were added to the model? How do these equations change for multiphase flows and what additional insight do they provide, i.e., what roles do, for example, the combination of heat flux, dry air flux, water vapor flux, and moisture flux in non-equilibrium thermodynamics play in tornadogenesis?

NOTES

1It is noted in that (32) may not allow a full study of the three-dimensional evolution of vorticity, as it would require the additional Clebsh potentials to match the number of degrees of freedom.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1] Fischer, J., Dahl, J.M.L., Coffer, B.E., Houser, J.L., Markowski, P.M., Parker, M.D., et al. (2024) Supercell Tornadogenesis: Recent Progress in Our State of Understanding. Bulletin of the American Meteorological Society, 105, E1084-E1097.[CrossRef]
[2] Markowski, P.M. (2024) A New Pathway for Tornadogenesis Exposed by Numerical Simulations of Supercells in Turbulent Environments. Journal of the Atmospheric Sciences, 81, 481-518.[CrossRef]
[3] Parker, M.D. (2023) How Well Must Surface Vorticity Be Organized for Tornadogenesis? Journal of the Atmospheric Sciences, 80, 1433-1448.[CrossRef]
[4] Perlin, K. (1985) An Image Synthesizer. In: Cole, P., Heilman, R. and Barsky, B.A., Eds., Proceedings of the 12th Annual Conference on Computer Graphics and Interactive Techniques, Association for Computing Machinery, 287-296.[CrossRef]
[5] Naylor, J. and Gilmore, M.S. (2012) Convective Initiation in an Idealized Cloud Model Using an Updraft Nudging Technique. Monthly Weather Review, 140, 3699-3705.[CrossRef]
[6] Kraichnan, R.H. (1967) Inertial Ranges in Two-Dimensional Turbulence. The Physics of Fluids, 10, 1417-1423.[CrossRef]
[7] Alexakis, A. (2023) Quasi-Two-Dimensional Turbulence. Reviews of Modern Plasma Physics, 7, Article No. 31.[CrossRef]
[8] Chorin, A.J. and Marsden, J.E. (1992) A Mathematical Introduction to Fluid Mechanics. Springer-Verlag.[CrossRef]
[9] Adrian, R.J. (2007) Hairpin Vortex Organization in Wall Turbulence. Physics of Fluids, 19, Article 041301.[CrossRef]
[10] Bernard, P.S. (2011) The Hairpin Vortex Illusion. Journal of Physics: Conference Series, 318, Article 062004.[CrossRef]
[11] Bernard, P.S. (2019) On the Inherent Bias of Swirling Strength in Defining Vortical Structure. Physics of Fluids, 31, Article 035107.[CrossRef]
[12] Orf, L., Wilhelmson, R., Lee, B., Finley, C. and Houston, A. (2017) Evolution of a Long-Track Violent Tornado within a Simulated Supercell. Bulletin of the American Meteorological Society, 98, 45-68.[CrossRef]
[13] Majda, A.J. and Bertozzi, A.L. (2002) Vorticity and Incompressible Flow. Cambridge University Press.[CrossRef]
[14] Klein, R. and Majda, A.J. (1991) Self-Stretching of a Perturbed Vortex Filament I. the Asymptotic Equation for Deviations from a Straight Line. Physica D: Nonlinear Phenomena, 49, 323-352.[CrossRef]
[15] Klein, R. and Majda, A.J. (1991) Self-Stretching of Perturbed Vortex Filaments: II. The Structure of Solutions. Physica D: Nonlinear Phenomena, 53, 267-294.[CrossRef]
[16] Müller, N.P. and Krstulovic, G. (2024) Exploring the Equivalence between Two-Dimensional Classical and Quantum Turbulence through Velocity Circulation Statistics. Physical Review Letters, 132, Article 094002.[CrossRef] [PubMed]
[17] Khesin, B., Misiolek, G. and Modin, K. (2018) Geometric Hydrodynamics via Madelung transform. Proceedings of the National Academy of Sciences, 115, 6165-6170.[CrossRef] [PubMed]
[18] KFOR-TV, Documentary (2024) Rapid Refresh Products.
https://www.youtube.com/watch?v=E5Dw3tXcIvA
[19] Ben-Amots, N. (2016) Dynamics and Thermodynamics of a Tornado: Rotation Effects. Atmospheric Research, 178, 320-328.[CrossRef]
[20] Cutlip, K. (2015) Interview with Richard Rotunno, a Senior Scientist at UCAR. Weatherwise, No. 8, 56-57.
[21] NOAA (2024) Significant Tornado Parameter (Fixed Layer).
https://www.spc.noaa.gov/exper/mesoanalysis/help/help_stor.html
[22] NOAA (2024) EF Scale.
https://www.weather.gov/oun/efscale
[23] NOAA (2024) Rapid Refresh Products.
http://www.nco.ncep.noaa.gov/pmb/products/rap
[24] NOAA (2024) Composite Indices.
https://www.spc.noaa.gov/exper/mesoanalysis/
[25] Doswell, C.A. and Schultz, D.M. (2021) On the Use of Indices and Parameters in Forecasting Severe Storms. E-Journal of Severe Storms Meteorology, 1, 1-14.[CrossRef]
[26] Elgindi, T.M. (2021) Finite-time Singularity Formation for C1,α Solutions to the Incompressible Euler Equations on R3. Annals of Mathematics, 194, 647-727.[CrossRef]
[27] Novikov, S.P., Manakov, S.V, Pitaevskii, L.B. and Zakharov, V.E. (1984) Theory of Solitons. The Inverse Scattering Method. Plenum Press.
[28] Bělík, P., Dokken, D.P., Potvin, C.K., Scholz, K. and Shvartsman, M.M. (2017) Applications of Vortex Gas Models to Tornadogenesis and Maintenance. Open Journal of Fluid Dynamics, 7, 596-622.[CrossRef]
[29] Bělík, P., Dahl, B., Dokken, D., K. Potvin, C., Scholz, K. and Mikhail Shvartsman, (2018) Possible Implications of Self-Similarity for Tornadogenesis and Maintenance. AIMS Mathematics, 3, 365-390.[CrossRef]
[30] lík, P., Dokken, D.P., Shvartsman, M.M., Bibelnieks, E., Laskowski, R. and Lukanen, A. (2023) Equilibrium Energy and Entropy of Vortex Filaments in the Context of Tornadogenesis and Tornadic Flows. Open Journal of Fluid Dynamics, 13, 144-176.[CrossRef]
[31] Paglietti, A. (2023) Why Thermodynamic Entropy and Statistical Entropy Are Two Different Physical Quantities. Current Physical Chemistry, 13, 233-245.[CrossRef]
[32] Chorin, A.J. (1994) Vorticity and Turbulence. Springer-Verlag.[CrossRef]
[33] Newton, P.K. (2001) N-Vortex Problem: Analytical Techniques. Springer.[CrossRef]
[34] Klein, R., Majda, A.J. and Damodaran, K. (1995) Simplified Equations for the Interaction of Nearly Parallel Vortex Filaments. Journal of Fluid Mechanics, 288, 201-248.[CrossRef]
[35] Caglioti, E., Lions, P.L., Marchioro, C. and Pulvirenti, M. (1992) A Special Class of Stationary Flows for Two-Dimensional Euler Equations: A Statistical Mechanics Description. Communications in Mathematical Physics, 143, 501-525.[CrossRef]
[36] Caglioti, E., Lions, P.L., Marchioro, C. and Pulvirenti, M. (1995) A Special Class of Stationary Flows for Two-Dimensional Euler Equations: A Statistical Mechanics Description. Part II. Communications in Mathematical Physics, 174, 229-260.[CrossRef]
[37] Lions, P. and Majda, A. (2000) Equilibrium Statistical Theory for Nearly Parallel Vortex Filaments. Communications on Pure and Applied Mathematics, 53, 76-142.[CrossRef]
[38] Flandoli, F. and Gubinelli, M. (2002) The Gibbs Ensemble of a Vortex Filament. Probability Theory and Related Fields, 122, 317-340.[CrossRef]
[39] Lebon, G., Jou, D. and Casas-Vazquez, J. (2008) Understanding Non-Equilibrium Thermo-Dynamics. Springer-Verlag.[CrossRef]
[40] Demirel, Y. and Gerbaud, V. (2019) Nonequilibrium Thermodynamics Approaches. In: Demirel, Y. and Gerbaud, V., Eds., Nonequilibrium Thermodynamics, Elsevier, 663-709.[CrossRef]
[41] Zavolzhenskii, M.V. (1980) Hydro-and Thermodynamics of Tornadoes and Oceanic Waterspouts. Journal of Applied Mechanics and Technical Physics, 21, 538-548.[CrossRef]
[42] Liu, Y. and Liu, C. (2008) Entropy Flow and the Evolution of a Storm. Entropy, 10, 430-440.[CrossRef]
[43] Landi, G.T. and Paternostro, M. (2021) Irreversible Entropy Production: From Classical to Quantum. Reviews of Modern Physics, 93, Article 035008.[CrossRef]
[44] Bertini, L., De Sole, A., Gabrielli, D., Jona-Lasinio, G. and Landim, C. (2015) Macroscopic Fluctuation Theory. Reviews of Modern Physics, 87, 593-636.[CrossRef]
[45] Wang, A., Pan, Y. and Markowski, P.M. (2020) The Influence of Turbulence Memory on Idealized Tornado Simulations. Monthly Weather Review, 148, 4875-4892.[CrossRef]
[46] Holton, J.R. and Hakim, G.J. (2019) An Introduction to Dynamic Meteorology. 6th Edition, Academic Pres.
[47] Dokken, D. and Shvartsman, M. (2006) Time Averaging, Hierarchy of the Governing Equations, and the Balance of Turbulent Kinetic Energy. Taylor & Francis Group, LLC, 155-164.
[48] Hauf, T. and Höller, H. (1987) Entropy and Potential Temperature. Journal of the Atmospheric Sciences, 44, 2887-2901.[CrossRef]
[49] Dutton, J.A. (1976) The Ceaseless Wind: An introduction to the Theory of Atmospheric Motion. McGraw-Hill.
[50] Sasaki, Y.K. (2014) Entropic Balance Theory and Variational Field Lagrangian Formalism: Tornadogenesis. Journal of the Atmospheric Sciences, 71, 2104-2113.[CrossRef]
[51] Salmon, R. (1998) Introduction to Geophysical Fluid Dynamics. Oxford University Press.[CrossRef]
[52] Davies-Jones, R. (2015) Formulas for Parcel Velocity and Vorticity in a Rotating Cartesian Coordinate System. Journal of the Atmospheric Sciences, 72, 3908-3922.[CrossRef]
[53] Houser, J.L., Bluestein, H.B., Thiem, K., Snyder, J., Reif, D. and Wienhoff, Z. (2022) Additional Evaluation of the Spatiotemporal Evolution of Rotation during Tornadogenesis Using Rapid-Scan Mobile Radar Observations. Monthly Weather Review, 150, 1639-1666.[CrossRef]
[54] Fiedler, B.H. and Rotunno, R. (1986) A Theory for the Maximum Windspeeds in Tornado-Like Vortices. Journal of the Atmospheric Sciences, 43, 2328-2340.[CrossRef]
[55] Bao, W. and Cai, Y. (2013) Mathematical Theory and Numerical Methods for Bose-Einstein Condensation. Kinetic & Related Models, 6, 1-135.[CrossRef]
[56] Barenghi, C.F. (2024) Tangled Vortex Lines: Dynamics, Geometry and Topology of Quantum Turbulence. In: Ricca, R.L. and Liu, X., Eds., Lecture Notes in Mathematics, Springer Nature Switzerland, 243-279.[CrossRef]
[57] Roitberg, A. and Ricca, R.L. (2021) Hydrodynamic Derivation of the Gross-Pitaevskii Equation in General Riemannian Metric. Journal of Physics A: Mathematical and Theoretical, 54, Article 315201.[CrossRef]
[58] Chen, T. and Pavlović, N. (2014) Derivation of the Cubic NLS and Gross-Pitaevskii Hierarchy from Manybody Dynamics in d = 3 Based on Spacetime Norms. Annales Henri Poincaré, 15, 543-588.[CrossRef]
[59] Carr, L.D., Miller, R.R., Bolton, D.R. and Strong, S.A. (2012) Nonlinear Scattering of a Bose-Einstein Condensate on a Rectangular Barrier. Physical Review A, 86, Article 023621.[CrossRef]
[60] Davis, M.J., Ballagh, R.J. and Burnett, K. (2001) Dynamics of Thermal Bose Fields in the Classical Limit. Journal of Physics B: Atomic, Molecular and Optical Physics, 34, 4487-4512.[CrossRef]

Copyright © 2026 by authors and Scientific Research Publishing Inc.

Creative Commons License

This work and the related PDF file are licensed under a Creative Commons Attribution 4.0 International License.