Forced, Balanced, Axisymmetric Shallow Water Model for Understanding Short-Term Tropical Cyclone Intensity and Wind Structure Changes

: A minimal modeling system for understanding tropical cyclone intensity and wind structure changes is introduced: Shallow Water Axisymmetric Model for Intensity (SWAMI). The forced, balanced, axisymmetric shallow water equations are reduced to a canonical potential vorticity (PV) production and inversion problem, whereby PV is produced through a mass sink (related to the diabatic heating) and inverted through a PV/absolute–angular–momentum invertibility principle. Because the invertibility principle is nonlinear, a Newton–Krylov method is used to iteratively obtain a numerical solution to the discrete problem. Two versions of the model are described: a physical radius version which neglects radial PV advection (SWAMI-r) and a potential radius version that naturally includes the advection in the quasi-Lagrangian coordinate (SWAMI-R). In idealized numerical simulations, SWAMI-R produces a thinner and more intense PV ring than SWAMI-r, demonstrating the role of axisymmetric radial PV advection in eyewall evolution. SWAMI-R always has lower intensiﬁcation rates than SWAMI-r because the reduction in PV footprint effect dominates the peak magnitude increase effect. SWAMI-r is next demonstrated as a potentially useful short-term wind structure forecasting tool using the newly added FLIGHT+ Dataset azimuthal means for initialization and forcing on three example cases: a slowly intensifying event, a rapid intensiﬁcation event, and a secondary wind maximum formation event. Then, SWAMI-r is evaluated using 63 intensifying cases. Even though the model is minimal, it is shown to have some skill in short-term intensity prediction, highlighting the known critical roles of the relationship between the radial structures of the vortex inertial stability and diabatic heating rate. Because of the simplicity of the models, SWAMI simulations are completed in seconds. Therefore, they may be of some use for hurricane nowcasting to short-term (less than 24 h) intensity and structure forecasting. Due to its favorable assumptions for tropical cyclone intensiﬁcation, a potential use of SWAMI is a reasonable short-term upper-bound intensity forecast if the storm intensiﬁes.


Introduction
Tropical cyclones (TCs) intensify through latent heating in deep convective clouds near the center of circulation [1,2]. This upward convective mass flux is fueled by sensible and latent heat fluxes from the underlying warm ocean, creating high moist entropy air in the boundary layer which sustains the convection. The surface fluxes, boundary layer, and moist convection interact with the vortex in complex and nonlinear ways. TC intensity and structure prediction remain great challenges because of these internal interactions as well as external factors such as vertical wind shear and dry air intrusions [3][4][5][6].
Although the boundary layer (particularly in the eyewall region) and outflow layer in TCs are often unbalanced, the layer in between has been observed to be in an approximate state of gradient wind [7] and hydrostatic balance (away from the eyewall). While the large inventory of observations analyzed in [7] suggest this is true in average sense, caution should always be taken for individual storms and periods of the storm life cycle; strict balance depends on how dominant the rotational flow is in relation to the divergent flow [8]. Although the flow field in TCs is three-dimensional and complex, the azimuthal wavenumber-zero (axisymmetric) component is usually the most significant, followed by the wavenumber-one asymmetry. Thus, axisymmetric dynamically-and energeticallybased theoretical frameworks have been developed over the years to understand the evolution of the TC vortex [2,[9][10][11]. One of the simplest frameworks to understand vortex evolution is Eliassen's axisymmetric balanced vortex model [12]. In Eliassen's model, a second-order linear transverse circulation equation was derived to understand the slow vortex response to heat and momentum sources under the assumptions of gradient and hydrostatic balance. Significant advances have been made through the application of this model to TC genesis, intensity, and structure change [13][14][15][16][17][18][19][20][21][22][23]. The foundational study for the present work is found in [13]. The authors obtained analytic solutions to the Eliassen balanced vortex model under some simplifying assumptions, and found that heating in the high inertial stability region near the storm center is far more efficient in producing tangential velocity accelerations and mid-level warming than heating outside this region. Recently, the work in [24] has provided observational support for this theory using airborne Doppler wind velocity composites. The authors of [13] illustrated that perhaps the two most critical factors governing hurricane structure and intensity change are the radius-height structures of inertial stability and diabatic heating.
Related to understanding TC evolution through balanced dynamics is the concept of potential vorticity (PV) thinking [25]. Any balanced equation set can be reduced to an equation governing the evolution of PV and an associated invertibility principle to recover the other quasi-balanced fields [26,27]. The PV framework is thus both a simpler and more elegant way to understand balanced dynamics than the full balanced equations. While PV is materially conserved in the absence of friction or diabatic heating, PV is generally not materially conserved in TCs because of strong latent heating in the eyewall and rain bands and near-surface frictional effects. PV production in the TC eyewall can be particularly rapid because the dot product of the absolute vorticity vector and gradient of diabatic heating is very large there [28]. Because of this rapid PV production, very large values of PV have been simulated [29,30] and observed [31,32] in the hurricane inner-core region. As TCs have highly curved flow, the appropriate horizontal balance is nonlinear balance. The nonlinearity of TC invertibility principles creates some challenges to obtaining wellbehaved numerical solutions. A recent example of a TC model that uses PV prediction and inversion is the work in [33], the authors of which obtained analytic solutions to the shallow water equations using the wave-vortex approximation to understand the intensification of TCs.
The purpose of the present work is to develop a minimal model for understanding short-term (lead times of less than 24 h) intensity and wind structure changes in real TCs using PV production and inversion. The model is based upon the forced, balanced, axisymmetric equations in a barotropic framework. Motivated by the works in [13,33], the evolution of the TC vortex in the minimal model is dependent on the radial structure of tangential velocity (or equivalently, inertial stability) and the diabatic heating rate. After describing the model, we will first demonstrate its usefulness in understanding TC evolution in an idealized framework. Then, we will demonstrate how the model can be initialized and forced with real data and be used to predict short-term changes in azimuthal mean wind structure and intensity. In Section 2, the model equations and simplified models are described. In Section 3, we present results of the models for an ideal scenario. The real-case initialization procedures and some limitations of the minimal model are given in Section 4. An evaluation of the model for real-case prediction is presented in Section 5. The summary and conclusions are given in Section 6.

Divergent Barotropic Model
The dynamical model is based upon the forced divergent barotropic (shallow water) equations in polar coordinates on an f -plane. The radial momentum, tangential momentum, and continuity equations are ∂Φ ∂t where r is the radius from the origin, h is the fluid depth, Φ = gh is the geopotential, φ is the azimuthal angle, u is the radial velocity, v is the tangential velocity, Q is the mass sink, and f is the Coriolis parameter. By taking the curl of (1) and (2) and combining with (3), while eliminating the divergence, we obtain the potential vorticity principle where P = ( f + ζ)/h is the potential vorticity, ζ = ∂(rv)/r∂r − ∂u/r∂φ is the relative vorticity, and D/Dt = (∂/∂t) + u(∂/∂r) + v(∂/r∂φ) is the material derivative. Equation (4) indicates that P is not materially conserved because of the mass sink term PQ. Equations (1)-(3) are the complete forced shallow water equations, valid for both axisymmetric and asymmetric flows and forcings, and containing slow (vortex Rossby) and fast (gravity) mode waves. We now describe two reduced models with two key simplifying assumptions on (1)-(3): (i) axisymmetry and (ii) balanced flow. The first model uses a physical radius coordinate and the second model uses a potential radius coordinate. The models are designed to describe the TC evolution in the lower troposphere.

Reduced Models Using PV Production and Inversion
As our ultimate intent will be to develop simplified models from (1)-(3) that use PV production and inversion, Equation (4) will be a sufficient starting point to develop both models. Under the assumption of axisymmetry (∂/∂φ = 0) in the material derivative, the PV production equation is ∂P ∂t Because of the radial PV advection term, an analytic solution to (5) is not possible. However, with neglecting of this term (see Section 4.4 for discussion of this limitation), the analytic solution to (5) is By removal of the radial advection term in (4), the mass sink directly forces local PV production, rather than contributing to both radial advection and local PV production. Furthermore, as this model is balanced, no energy from the mass sink can go into gravity waves, which can occur in numerical simulations of the full shallow water Equations (1)- (3). In this minimal model, even outside the region of high inertial stability, the mass sink will directly contribute to the tangential velocity spin-up there. The minimal model evolves PV according to (6) for a period of time, and then at any instant the geopotential and tangential velocity may be obtained by an invertibility principle, subject to gradient wind balance, where M = rv + (1/2) f r 2 is the absolute angular momentum per unit mass. Three invertibility principles are given as follows: (i) potential vorticity and geopotential (P-Φ), (ii) potential vorticity and absolute angular momentum (P-M), and (iii) potential vorticity and tangential velocity (P-v). The P-Φ invertibility principle is derived by eliminating the ∂M/∂r term in the numerator of P = (g∂M/r∂r)/Φ using the gradient wind balance Equation (7). The P-M invertibility principle is derived by eliminating the ∂Φ/∂r and Φ terms in (∂/∂r)[PΦ − g∂M/r∂r] = 0 using the gradient wind balance equation and the P definition. The P-v invertibility principle is obtained by eliminating the ∂Φ/∂r term in the gradient wind balance equation using the P definition. The P-Φ invertibility principle is Equation (8) is solved for unknown Φ given known P, which is produced by (6). The boundary conditions are Φ r→∞ = Φ ref (Dirichlet) and ∂Φ/∂r r=0 = 0 (Neumann). Once Φ is known, the balanced tangential velocity v may be obtained by solving the gradient wind balance equation. The P-M invertibility principle is Equation (9) is solved for unknown M given known P which is produced by (6). The Dirichlet boundary conditions are M r=0 = 0 and M r→∞ = (1/2) f r 2 r→∞ . Once v is known, Φ is then obtained by radially integrating the gradient wind balance equation. The P-v invertibility principle is Equation (10) is solved for unknown v given known P which is produced by (6). The Dirichlet boundary conditions are v r=0 = 0 and v r→∞ = 0. Once M is known, the tangential velocity v is obtained by M = rv + (1/2) f r 2 , and Φ is then obtained by solving the gradient wind balance equation. Note that the no-vortex solution of v = 0 satisfies (10) exactly, the no-vortex solution of M = (1/2) f r 2 and P = f g/Φ satisfies (9) exactly, and the no-vortex solution of P = f g/Φ satisfies (8) exactly.
An advantage of potential radius coordinates is that the radial advection term is absorbed into the coordinate allowing for a natural analytic solution to the PV production equation. The potential radius is defined as the radius a parcel would need to brought to in order to reduce the tangential velocity to zero [34], i.e., (1/2) f R 2 = rv + (1/2) f r 2 . Transforming from (r, t) to (R, T) coordinates, where ∂/∂T = ∂/∂t +Ṙ∂/∂r, T = t anḋ R = 0 (parcels remain at a fixed R), The same P-Φ and P-M invertibility principles can be used for the potential radius model. However, the solution procedure for this model is different than the physical radius model because R changes in time due to changes in v. Hereafter, the minimal modeling system will be referred to as the Shallow Water Axisymmetric Model for Intensity (SWAMI). The physical radius version without radial PV advection will be called SWAMI-r and the potential radius version will be called SWAMI-R. In both models, we have elected to use the P-M invertibility principle because the numerical solution procedures are simpler than the P-Φ principle. Numerical solution procedures to the P-v invertibility principle are similar to the P-M invertibility principle.
In order to halt unbridled growth of PV and intensification, a brake is needed when the TC intensity approaches the maximum potential intensity (MPI). One way to accomplish this would be to add a brake to the PV production Equations (5) and (11) from frictional drag, however such a brake would then prevent obtaining an analytic solution to the PV production equation and violate the balance assumption in the model. We have elected to construct the braking mechanism through a logistic limiter on the PV production equation, similar to that in [35]. An MPI vortex can be constructed using an energetically-based V mpi [10] which is valid at the radius of maximum winds, using assumptions on the radial vortex structure (e.g., a Rankine vortex). Once the radial structure of the MPI vortex tangential velocity is determined, the MPI vortex geopotential can be found using the gradient wind balance equation. With the known geopotential and tangential velocity of MPI vortex as a function of radius, the MPI vortex PV P mpi (r) can be computed using the PV definition. With the MPI limiter P mpi (r), the PV production equation becomes a logistic growth equation. The logistic growth equation and analytic solution in physical radius coordinates are Equation (14) indicates that PV production is reduced as P(r, t) approaches P mpi (r). A similar logistically-limited equation set also exists in (R, T) coordinates, The numerical solution procedures for SWAMI-r and SWAMI-R are given in Appendix A.

Ideal Case Setup and Results
We first describe how the models are initialized and forced in ideal cases, and present the results. The models require the initial PV and mass sink: P(r, 0) and Q(r, t) for SWAMIr and P(R, 0) and Q(R, T) for SWAMI-R. In order to illustrate the utility of the model for ideal cases, we use simple functions for both quantities. The initial vortex is a Gaussian vortex [36], where Γ = 16, 000, 000 m 2 s −1 and b = 50, 000 m are constants that determine the intensity and radial structure, respectively. To compute P(r, 0), the gradient wind balance equation is solved using v(r, 0) to determine Φ(r, 0). Then, both quantities are substituted into the P definition. The heating function is a Gaussian ring that is invariant in time defined by where µ R = 200,000 m, σ R = 50, 000 m, and Q t = 30 K. For SWAMI-R, Q(R) can change each outer loop iteration (Appendix A) as R changes due to changes in v. For SWAMI-r, Q(R) is remapped to Q(r) on the first iteration using R(r), and remains invariant in space. A constant radial grid spacing ∆r = 1000 m is used with N r = 500, for a total domain size of r f = 500 km. The discrete problem radial boundary conditions are M r=0 = 0 and M r f = rv(r f ) + (1/2) f r 2 f . SWAMI-R uses an outer loop (see Appendix A) number of iterations N m = 10 and inner loop maximum iterations N nk = 20. SWAMI-r uses an inner loop with maximum iterations N nk = 200. Convergence is based on the absolute value of the residuals of the nonlinear problem becoming less than the tolerance of 10 −12 . The rate of convergence depends most strongly on the complexity of the radial structures of the mass sink and tangential velocity. For the ideal cases, residuals less than 10 −12 are rapidly obtained well before N m N nk = 200 for the R-solver and before N nk = 200 for the r-solver.
In Figures 1-3, the ideal case results results are shown forθ = 0, 3, 6 K h −1 , respectively (see discussion in Section 4.3 on the conversion of the diabatic heating rate to an equivalent mass sink). PV is grown analytically from t = 0-12 h. Forθ = 0 K h −1 (Figure 1), SWAMI-r and SWAMI-R forecasts are identical because there is no forcing. This result demonstrates that SWAMI-r and SWAMI-R can recover the all the quasi-balanced fields solely from an axisymmetric PV profile along with boundary conditions. Withθ = 3 K h −1 (Figure 2), both simulations evolve the Gaussian PV profile to a hollow profile, and the mass sink and PV are thinner in SWAMI-R than SWAMI-r. The PV ring thinning effect is amplified significantly withθ = 6 K h −1 in Figure 3. For all positive heating rates, SWAMI-r has a larger intensification rate than SWAMI-R because of the lack of radial PV advection and no PV ring thinning. This result is consistent with numerical simulations of the forced shallow water Equations (1)-(3) (discretized on a Cartesian grid) from in [37]. In that study, radial PV advection erodes the outer edge of the PV ring under forcing from the mass sink, thinning the ring, and yielding a lower intensity because of a reduced PV footprint. SWAMI-R captures this effect in a simple way using the R-coordinate. In both SWAMI-r and SWAMI-R, higher heating rates cause higher peak vortex intensities at t = 12 h.
To demonstrate the MPI limiter on PV production, theθ = 6 K h −1 case is simulated again using an MPI vortex that is constructed from (17) with Γ = Γ mpi = 24,000,000 m 2 s −1 . The MPI PV P mpi is obtained from a balanced vortex assumption on the tangential velocity. In Figure 4, the SWAMI-r and SWAMI-R forecasts are shown with the MPI limiters [ (14) and (16) with invertibility principle (9)]. Both models show little intensification because of the MPI limiter on PV production. Another test was done with Γ = Γ mpi = 16,000,000 m 2 s −1 , so that the MPI vortex is exactly the same as the initial vortex. In this case, the results are identical to theθ = 0 K h −1 case ( Figure 1). These results demonstrate that both SWAMI-r and SWAMI-R can be used with the MPI limiter when a tropical cyclone is close to its MPI to prevent unbridled intensification past the MPI.

Real Case Initialization
The SWAMI initialization and forcing procedures for real cases are more complicated than for ideal cases and are described in the next three subsections. First, we describe the aircraft reconnaissance data used, the data processing procedures, and the real-case model initialization and forcing procedures. Then, we provide a discussion of the limitations of SWAMI for real hurricane intensity and wind structure prediction.

Observational Data
For real cases, the initial PV and mass sink are computed from the Extended Flight Level Dataset for Tropical Cyclones (FLIGHT+) [38]. The FLIGHT+ Dataset consists of flight-level (FL) data for a number of meteorological parameters and surface radial FL data from the Stepped Frequency Microwave Radiometer (SFMR) instrument from all aircraft reconnaissance missions into North Atlantic and eastern North Pacific TCs from 1997 to present. The research-grade dataset standardizes FL data collected from the National Oceanic and Atmospheric Administration (NOAA) Hurricane Hunters and U.S. Air Force Reserve into convenient forms for use in process studies, composite studies, and other scientific and industry endeavors. These include (i) an earth-relative frame, in which the data are given in latitude and longitude coordinates; (ii) a storm-relative frame, in which data can be referenced by the zonal and meridional distance from the storm center; and (iii) a frame moving with the storm center and referenced in radius space. Additionally, FL data have undergone visual and automated quality control checks and have been parsed by radial leg and interpolated to a standardized radial grid.
An automated azimuthal mean computational capability was recently added to the FLIGHT+ Dataset. The code computes azimuthal means of 12 different quantities using Level-3 radial flight leg data from the FLIGHT+ Dataset. The current quantities are FL radial velocity, FL tangential velocity, FL vertical velocity, SFMR surface rain rate, SFMR surface wind speed, FL pressure, sea level pressure, FL temperature, FL inflow angle, FL dew point temperature, FL absolute vertical vorticity, and FL inertial stability. The L3 data consist of storm-relative radial flight legs from aircraft reconnaissance missions into a given storm. The procedure for computing azimuthal means is as follows. First, an azimuthal mean time window is specified. We are currently using a time window size of 3 h. The purpose of a smaller window is to capture radial flights legs that are close to each other in time so that the azimuthal mean is more representative of an instantaneous value. The window of 3 h works well for most flights, capturing the azimuthal mean structure after at least 4 radial legs are completed during a typical aircraft reconnaissance mission. Next, a requirement is enforced that data must exist in each of the four quadrants at a given radius. This requirement ensures that the azimuthal mean has each quadrant represented, and is not erroneously representative of one part of the storm. Finally, depending on how many radial legs passed through each quadrant, the azimuthal mean is computed using weights at each radius. Azimuthal means are computed at each radius, and therefore some radial segments may have a valid azimuthal mean while others may not have a valid azimuthal mean. As an example, consider a case where a flight had 6 radial legs through the storm in the 3 h time window in the following quadrants: 0-90 • : legs 1-2, 90-180 • : legs 3-4, 180-270 • : leg 5, and 270-360 • : leg 6. The azimuthal mean of quantity A(r i ) (where r i is the radius of the L3 data at radial grid index i), In both SWAMI-R and SWAMI-r, the FLIGHT+ Dataset azimuthal mean FL tangential velocityv and azimuthal mean SFMR rain rate are used for initialization and forcing, respectively. The UTC time of the azimuthal mean group is the average time of all radial leg center passes in the azimuthal mean group, and does not necessarily correspond to a 0000, 0600, 1200, or 1800 UTC synoptic time. A new Level-4 netCDF output file is created with the azimuthal means of the 12 quantities.
The azimuthal mean processing is illustrated in Figure 5a for azimuthal mean group 14 of Hurricane Dorian (2019). In this case, nine valid radial legs were found in a 3 h period. The tangential velocity for those legs are shown in red, and the azimuthal mean is shown in blue. The first 14 valid azimuthal mean tangential velocities are then shown in Figure 5b. The azimuthal mean profiles illustrate the intensification of Hurricane Dorian from an initially flat profile to an intense vortex, with a peak azimuthal mean tangential velocity of 60 m s −1 and a radius of maximum azimuthal mean tangential velocity of approximately 20 km.

Data Processing
The FLIGHT+ Dataset azimuthal mean quantities are relative noisy in space, particularly the SFMR rain rate. Additionally, the valid azimuthal means typically do not extend much farther than 100-150 km from the vortex center because reconnaissance missions are typically focused on the inner core. Therefore, both extrapolation and smoothing procedures are necessary for use of these data in the PV inversion models which require a larger radial domain. Smoothing is particularly important because nonlinear elliptic problems such as (9) can be very sensitive to high frequency spatial variations.
First, the tangential velocity and SFMR rain rate are interpolated from the FLIGHT+ radial grid to the model radial grid. The FLIGHT+ radial grid has a grid spacing of 100 m, and the model grid spacing is typically ∆r = 500-1000 m. For the tangential velocity extrapolation, a modified Rankine vortex assumption is made for the outer wind field decay. The modified Rankine decay parameter is α = ln(v last /v m )/ ln(r m /r last ), where v m is the peak FLIGHT+ azimuthal mean tangential velocity at r = r m , and v last is the last valid azimuthal mean point in the outer wind field at r = r last . The decay parameter α could also be computed using the National Hurricane Center (NHC) wind radii data, such as the maximum sustained 1-min wind (VMAX) and radius of maximum winds (RMW) with the 34 kt wind (V34) and the radius of 34 kt winds (R34). For cases where the FLIGHT+ data end abruptly after the RMW, the NHC values can be useful to construct a reasonable outer wind field. However, some caution should be taken for use in SWAMI since the NHC values are not azimuthal mean values. The azimuthal mean SFMR rain rate typically approaches zero with increasing radius, so this zero value is extrapolated outward. The extrapolations yield values at all radii on the model grid. Smoothing of the tangential velocity and SFMR rain rate are done using a Savitsky-Golay filter. The filter works by specifying a radial window and polynomial order to fit a polynomial to the data within the window. For the tangential velocity, a window of 401 radial grid points is used with a polynomial order 3. For the SFMR rain rate, a window of 201 points with a polynomial order of 3 is used. As an example of the extrapolation and smoothing procedures, the FLIGHT+ Dataset azimuthal means along with the smoothed and extrapolated values to be used in the models are given in Figure 6 for azimuthal mean group 4 of Hurricane Dorian (2019).

Model Initialization and Forcing
Having described the FLIGHT+ Dataset azimuthal means and data processing, we now describe how SWAMI is initialized and forced with these real data. First, the FLIGHT+ Dataset azimuthal mean groups are filtered to only include those with average FL pressures within r < 150 km between 600 and 900 hPa. Second, the gradient wind balance Equation (7) is used to obtain the axisymmetric geopotential Φ from the azimuthal mean FLIGHT+ Dataset tangential velocityv. Then, the initial azimuthal mean PVP(r, 0) is obtained by substituting the geopotential and tangential velocity in the PV definition (where the axisymmetric absolute vorticity is obtained via finite differencing on the azimuthal mean tangential velocity). While we have elected to use a balance assumption in computing the initial PV, PV can also be directly computed from radial flight legs without a balance assumption (Appendix B). The mass sink is obtained via a conversion using the SFMR rain rate. This accomplished by first converting the SFMR rain rate to an equivalent diabatic heating rate, and then converting the diabatic heating rate to an equivalent mass sink for the shallow water models. The Coriolis parameter f is set to the average latitude of the azimuthal mean group.
The conversion of the SFMR rain rate to an estimated mid-level diabatic heating rate is based upon Figure 4 in [39]. This figure depicts a calibration curve between the diabatic heating rate and precipitation rate for numerical simulations of a squall line. Based on this curve, the conversion is 1.53 K h −1 = 1 mm h −1 . We acknowledge that there is a some amount of uncertainty in attempting to ascertain a mid-tropospheric peak diabatic heating rate from the near-surface SFMR rain rate. Additionally, our model assumes that the SFMR rain rate is collocated radially with the mid-tropospheric diabatic heating maximum. However, there is some uncertainty in this assumption because of slantwise convection in real hurricanes.
The conversion of the diabatic heating rate to the mass sink is based upon the works in [21,37], and is described below. First, an analogy is made between the shallow water continuity Equation (3) and the continuously stratified continuity equation using an isentropic vertical coordinate, ∂σ ∂t where σ = −(1/g)(∂p/∂θ) is the pseudo-density andθ is the diabatic heating rate. When σ varies slowly in the vertical (which often occurs at midlevels), the right hand side of Equation (19) can be approximated by −σ(∂θ/∂θ). In an analogy between (3) and (19), the axisymmetric shallow water geopotential Φ(r, t) is analogous to σ(r, θ, t). Next, we assume a two-layer version of the θ-coordinate model is used, with a lower tropospheric layer defined by 300 < θ < 330 K and a peak diabatic heating rate at θ = 330 K ofθ max (r, t).
With these assumptions, ∂θ/∂θ can be approximated byθ max (r, t)/(30 K). In summary, the SFMR rain rate is converted to an equivalent mass sink by where Q is the mass sink in Equation (3),q r is the near-surface SFMR rain rate, and the constant c is SWAMI-R is initialized and forced in the same manner as SWAMI-r. The only difference is thatP(r, 0),Q(r, t) are mapped toP(R, 0),Q(R, T) using R(r) at each outer loop iteration m (Appendix A). Then, Equation (9) is also solved every outer loop iteration. For SWAMI-R to produce a thinning PV ring in time for real cases similar to the ideal cases, at each outer loop iteration m,Q(R, T) must be fit to an analytic function of R. Considering the complex SFMR rain rate structures, high-order polynomials must be used to capture all the radial variations, and these polynomials may not be able to reproduce the actual SFMR rain rate azimuthal mean radial structure. For this reason, we currently use SWAMI-r for the real cases.

Limitations of the Reduced Models
SWAMI-R and SWAMI-r are reduced models, and as such they cannot possibly capture the wealth of factors in real TCs that are responsible for intensity variability. Because their simplicity, the models would be expected to have poorer prediction skill in comparison to "full physics" dynamical models. Below, we list how some key assumptions in SWAMI affect the veracity of the predicted intensity change. First, the models are axisymmetric balance models. The models cannot capture any asymmetric intensification mechanisms nor unbalanced intensification mechanisms associated with the planetary boundary layer. It is well known that asymmetric and boundary layer processes are important for TC intensity change [40]. Second, the heat source is specified, rather than being allowed to develop naturally through microphysical equations a result of surface fluxes, boundary layer convergence, and other environmental factors. This is different from the windinduced surface heat exchange (WISHE) [10] mechanism which includes the positive feedback loop between the surface fluxes and intensity. Third, all of the mass sink energy directly forces vortex spin-up. In real hurricanes, energy can be lost to the environment via gravity waves, particularly for convection outside the high inertial stability region [14]. Fourth, no negative environmental factors are included. The reason for not including these factors is that they are often asymmetric and cannot be properly accounted for in an one-layer axisymmetric model (e.g., vertical wind shear and dry air intrusions). Fifth, related to the first point, real TCs have significant azimuthal variability in wind fields and convection. By initializing and forcing with azimuthal mean quantities, the possible importance of this variability is not accounted for. As will be shown, the individual SFMR rain rate and tangential velocity radial flight legs exhibit significant variability about the the azimuthal mean value. Sixth, the SFMR rain rate is assumed to be linearly related to mid-level peak diabatic heating rate, and radial profiles of the diabatic heating rate and mass sink are assumed to directly correspond to the radial profile of the SFMR rain rate. Because of slantwise convection, is it possible that there is a radial offset between the mid-level diabatic heating and near-surface SFMR rain rate. Seventh, the convective forcing based upon the SFMR rain rate is assumed constant in time for the entire intensity forecast. In real hurricanes, the convective structure can fluctuate on short time scales. For this reason, SWAMI should primarily be used for nowcasting to short-term (less than 24 h) intensity and structure forecasting. Eighth, SWAMI-r neglects the radial advection of PV, which is an important physical process. The neglecting of this term will tend to yield larger balanced intensification rates than with this term included in SWAMI-R (see Figures 2 and 3). Additionally, the radial movement of the radius of maximum winds cannot properly be simulated using SWAMI-r. The assumptions of axisymmetry, no gravity wave response to convection, and no factors causing weakening will lead to SWAMI producing a higher intensification rates than observed on average. In this context, the most suitable use for SWAMI is for a reasonable upper bound on the intensification rate.

Real-Case Studies
In this section, we present the results of SWAMI-r (hereafter SWAMI) on a set of real cases. The list of real cases is shown in Table 1. There are 63 cases total for North Atlantic and eastern North Pacific tropical cyclones from 2015 to 2019. These cases are all intensifying events, defined by a positive change in the maximum FLIGHT+ Dataset azimuthal mean tangential velocity from the first date-time-group of the azimuthal mean group to the subsequent (second) date-time-group. SWAMI is initialized at the first datetime-group in the second column and verified at the next date-time-group in the third column. The verifying date-time-group is always defined as the next valid azimuthal mean group in the FLIGHT+ Dataset azimuthal means Level-4 netCDF file. SWAMI forecasts for real cases are on a radial grid of length 500 km with a radial grid spacing ∆r = 1 km. Before describing the results of the entire sample, we first describe some example cases: a slow intensifying event of Hurricane Joaquin (2015), formation of a secondary wind maxima in Hurricane Irma (2017), and a rapid intensification event of Hurricane Ignacio (2015). The MPI limiter on SWAMI is not used for these real cases because of the short term forecasts and no storm is close to its MPI.
In Figure 7, SWAMI forecasts of a slow intensifying event of Hurricane Joaquin (2015) are given. By the FLIGHT+ Dataset azimuthal means, Joaquin intensified from approximately 49 m s −1 to 53 m s −1 in a period of approximately 12 h. Figure 7 shows that the large mass sink inferred from the azimuthal mean SFMR rain rates within r = 20-50 km produces large PV there, spinning up the vortex. The verification (fourth panel of Figure 7) shows that the SWAMI forecast produces an accurate peak azimuthal mean tangential velocity and reasonably accurate wind structure forecast for r < 130 km.
In Figure 8, we demonstrate the ability of SWAMI to forecast a secondary wind maximum in Hurricane Irma (2017). The peak values of the mass sink are between r = 50-100 km. PV production primarily occurs between r = 30-50 km and r = 50-100 km because of the mass sink structure. The production between r = 50-100 km enhances the existing weak secondary PV ring there. In verifying using the next azimuthal mean group (fourth panel of Figure 8), SWAMI forecasts a distinct secondary wind maximum at r = 80 km. The FLIGHT+ Dataset azimuthal mean observations have a secondary wind maximum at a slightly smaller radius of r = 60 km. These results show the potential use of SWAMI to provide short-term forecasts of the secondary wind maxima associated with the formation of the secondary eyewall.
In Figure 9, the performance of SWAMI for a rapid intensification event of Hurricane Ignacio (2015) is shown. PV is produced rapidly between r = 30-50 km at the intersecting region of the large PV and large mass sink, leading to a SWAMI forecast of an increase in the peak azimuthal mean tangential velocity of approximately 10 m s −1 in 14 h. In verifying with the next azimuthal mean group (fourth panel of Figure 9), SWAMI is not able to capture the rapid intensification that occurred. In Figure 10, we double the heating (or equivalently, the mass sink forcing), and in this case SWAMI is able to capture the rapid intensification. One of the most uncertain parameters in SWAMI is Q 0 , or the relationship between the SFMR rain rate with the mid-level diabatic heating rate and subsequently the mass sink. Figure 10 illustrates that larger values of Q 0 could help capture rapid intensification, however those larger values may not yield as accurate of forecasts for other weaker intensification events. Because of the balance assumption in SWAMI, it may not be as well suited to capture extreme RI events which have significant unbalanced influences. In any event, Q 0 can be better calibrated within the range of uncertainty in SWAMI through testing of numerous cases; it is indeed the parameter with the most uncertainty in the model.
The performance of SWAMI for all 63 intensifying events is given in Figure 11, with an average lead time of 13.6 h. In order to establish a reasonable baseline comparison for SWAMI, a persistence forecast of the initial (current) intensity is also shown. There is some utility of SWAMI in short-term forecasting of the maximum azimuthal mean tangential velocity, as shown by the linear correlation. The SWAMI mean error (defined using model forecast minus the observations) of the peak azimuthal mean tangential velocity is 1.9 m s −1 and the mean absolute error is 5.8 m s −1 . The persistence forecast mean error of the peak azimuthal mean tangential velocity is −5.7 m s −1 and the mean absolute error is 5.7 m s −1 . Thus, SWAMI has an improved bias over the persistence baseline and similar errors. The SWAMI mean error of the tangential velocity structure for r < 150 km is 1.8 m s −1 and the mean absolute error is 4.3 m s −1 . The mean errors indicate that on average, SWAMI has positive biases in the peak azimuthal tangential velocity forecast and outer wind field structure forecasts. As discussed previously, the assumptions of the model and inclusion of no negative environmental factors contribute to the positive biases. The mean absolute errors are somewhat large for short-term forecasts (average lead time of 13.6 h). However, these errors should be interpreted in the context that SWAMI is a minimal model using only two parameters from the FLIGHT+ Dataset azimuthal means to forecast TC intensity change. State-of-the-art TC forecasting models are vastly more complicated that SWAMI and include all physical processes contributing to TC intensity change. Thus, the errors of those models will be lower than SWAMI. As an example, average intensity forecast errors of the operational intensity consensus IVCN are approximately 4.1 m s −1 at the 24-h lead time [41].     Table 1 (63 cases). The mean forecast time of the SWAMI forecasts is 13.6 h, and the linear correlation coefficient is R = 0.86.
We also examined the correlations of FLIGHT+ Dataset azimuthal mean quantities relevant to this study. A scatterplot of the relationship between the maximum azimuthal mean intensification rate by change in the tangential velocity to the subsequent azimuthal mean group (including the 63 intensifying cases as well as the neutral and weakening cases) versus the maximum azimuthal mean SFMR rain rate is given in Figure 12a. A scatterplot of the relationship between the maximum azimuthal mean tangential velocity versus the maximum azimuthal mean SFMR rain rate is given in Figure 12b. There is no correlation between the intensification rate and the peak SFMR rain rate, while there is a positive correlation between the peak SFMR rain rate and the peak tangential velocity. Figure 12b is consistent with the findings in [42], whose authors found significant correlations between total condensate and the peak tangential velocity in tropical cyclones using CloudSat data. Thus, they argued that total condensate could be used as a metric for the intensity of the storm. There is no such correlation with intensity change (Figure 12a). The results herein demonstrate that the SFMR rain rate can be used to predict intensity change, but only when used as a proxy for heating in an appropriate dynamical model.

Conclusions
We have introduced a minimal modeling system, referred to as the Shallow Water Axisymmetric Model for Intensity (SWAMI), for understanding short-term (lead times less than 24 h) TC intensity and wind structure changes. The forced, balanced, axisymmetric shallow water equations are cast in a canonical PV production and inversion problem in order to isolate the fundamental way that deep convection in TCs interacts with the potential vorticity monopole or ring. Two versions of this modeling system are described: The first version is a physical radius version that neglects the radial PV advection term (SWAMI-r), and the second version is a potential radius version in which radial PV advection is naturally included in the coordinate (SWAMI-R). In both versions, PV is produced through the term PQ, and thus intersecting regions of large values of the mass sink and potential vorticity can cause exponential growth of PV, leading to large intensification rates. A simple MPI limiter is proposed to prevent unbridled PV growth when a storm is at or near its MPI. The diabatically-produced PV is inverted through a PV-absolute angular momentum invertibility principle in order to obtain the other balanced fields.
SWAMI is tested using both ideal and real case studies. In the ideal studies, the critical role of radial PV advection in causing an intensifying an thinning PV ring in time is demonstrated using the potential radius coordinate. Although the peak PV is larger in SWAMI-R in comparison to SWAMI-r, the thinner PV ring always yields a lower intensification rate. Additionally, the critical relationship between the radial structures of the mass sink and initial PV in governing the subsequent intensification is demonstrated.
The novel aspect of this study is the application of SWAMI to real cases. Because SWAMI requires estimates of azimuthal mean PV and mass sink radial profiles in real storms, these data were obtained by adding an automated azimuthal mean computational capability to the FLIGHT+ Dataset. The code produces a new Level-4 netCDF output file containing azimuthal means of 12 quantities. To initialize SWAMI, the azimuthal mean FL tangential velocity is used to estimate an azimuthal mean PV using a gradient wind balance assumption. SWAMI is forced by converting the azimuthal mean SFMR rain rate to an equivalent mass sink under the assumption that the SFMR rain rate is proxy for the mid-level diabatic heating rate. We first demonstrated the usefulness of SWAMI in capturing a slowly intensifying event, a rapid intensifying event, and a secondary wind maximum formation case. Through the real case examples, SWAMI was shown to have some short-term intensity and wind structure prediction skill. However, the model did not capture the observed RI of Hurricane Ignacio (2015), indicating that the balanced dynamics assumption may not be adequate for some RI events. The model was tested on a larger sample of 63 intensifying cases, and was demonstrated to be able to predict the 13.6 h intensity change of those cases reasonably well. This clearly illustrates that the inferred azimuthal mean radial structures of heating and vortex inertial stability are indeed important for TC intensity change, as illustrated by past studies [13,24]. SWAMI captures the interaction of the TC vortex with heating in a fundamental way, which has been known for many years to be critical for hurricane intensification. SWAMI can be useful for understanding short-term balanced responses of the vortex azimuthal mean tangential velocity to diabatic heating.
While we have demonstrated the potential usefulness of SWAMI in a hindcast mode, SWAMI could also be run in real-time right after an aircraft reconnaissance mission, provided Level-3 FLIGHT+ Dataset files are produced shortly after the flight and radial flight legs are executed into all four quadrants with an estimated center. Because SWAMI simulations are completed in seconds, they could potentially be useful for nowcasting and short term (less than 24 h) intensity and structure forecasting. Since the favorable assumptions in SWAMI always produce intensification, a potential use of the model is a reasonable upper bound on short-term intensification. Additionally, the current FLIGHT+ Dataset azimuthal mean structures of the 12 variables could potentially be useful to forecasters, and in particular the relationship between the radial structures of the SFMR rain rate and tangential velocity.
Due to of its simplicity and lack of inclusion of many factors governing hurricane intensity change (Section 4.4), SWAMI is unlikely to have the skill of full-physics dynamical models or sophisticated statistical-dynamical models at 0-24 hour lead times during intensifying events. Therefore, we recommend that SWAMI be merged in an intelligent way with other skillful prediction tools in order to be of potential use for real hurricane intensity forecasting. Additionally, as a step in complexity above SWAMI, the two-dimensional shallow water Equations (1)-(3) could be directly forced with the SFMR rain rate estimated heating. Integration of these equations would allow for more realistic effects such as barotropic instability of the PV ring, PV mixing between the eyewall and eye, and gravity wave responses to heating.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Abbreviations
The following abbreviations are used in this manuscript:

Appendix A. Numerical Solution Methods for SWAMI-r and SWAMI-R
The physical radius model SWAMI-r is discretized on a radial grid with constant grid spacing ∆r, N r radial grid points, and grid point index i. Time is discretized with a time interval ∆t, N t time points, and time index j. With Q(r, t) = Q(r) constant in time, Equation (5) becomes P i,j = P i,0 exp Q(r i )t j .
If an MPI limiter is chosen on PV production, Using finite differencing, the second-order-accurate discrete form of (9) is where A i = 1, B i = −(1/r i + (1/P i,j ))(∂P/∂r) i , C i = −P i,j /(gr 2 i ), and D i = P i,j f 2 r 2 i /(4g). Equation (A2) results in a system of N r nonlinear equations which can be compactly written in vector form as F(M) = 0, where F is the vector of nonlinear operators and M = (M 0 , . . . , M i−1 , M i , M i+1 , . . . , M N r ) is the solution vector. A Newton-Krylov method [43] is used to iteratively solve the nonlinear problem (with N nk iterations), and a residual of less than 10 −12 is required for the discrete problem to be considered converged. The Krylov solver used is the Loose Generalized Minimal Residual Method (LGMRES). As discussed in [44], LGMRES can be used to accelerate the convergence of the GMRES algorithms.
The numerical solution procedure for the potential radius model SWAMI-R is similar to the physical radius model, but requires and outer loop to update coefficients and then solve the invertibility principle multiple times with the updated coefficients. One outer loop iteration m solves the following sequential equations: If an MPI limiter is not chosen on PV production, , and the loop repeats. The outer loop is executed to N m iterations, and the Newton-Krylov inner loop (N nk iterations) is run at each outer loop iteration m. The same convergence criterion is used for the potential radius model. In SWAMI-R, R is the independent variable and r is the dependent variable. Thus, the model captures changes in the radial structure of Q(r), P(r) in physical space in response to Q(R). As shown in the main paper, the potential radius model captures effects such as a thinning PV ring which do not occur in the physical radius model with the neglecting of the PV advection term.
Using potential temperature θ as the vertical coordinate, the dry PV in isentropic coordinates (i.e. Ertel's PV) for the quasi-balanced axisymmetric equations is Assuming an aircraft is able to fly along a constant potential temperature (isentropic) surface θ f l ≈ 315 K (near p = 700 hPa), Ertel's PV can be computed as below. The flight level f l and level of PV estimation l are similar to above, however the subscripts now denote isentropic surfaces. The PV estimation level l can be considered here to be approximately half way between the flight level and surface, or at θ l = 307 K.
The above PV computations can either be performed for one radial leg, or the average of multiple radial legs. It is recommended that the average of multiple radial legs are used, so that the PV will be a representation of the azimuthal mean structure, and not reflect an asymmetry. For more realism, moist potential vorticity could also be computed, provided virtual temperature is known from aircraft observations of the moisture content at both the flight level and the surface.
The above PV values can be converted to an equivalent shallow water PV by normalizing the denominator by adding a reference environmental value to the numerator. In that case, the azimuthal mean PV could be directly used to initialize SWAMI-r or SWAMI-R.