Dynamics of diabatically forced anticyclonic plumes in the stratosphere

A new class of vortices has been observed in the stratosphere following extreme wildfires (Canada 2017, Australia 2020) and volcanic eruptions (Raikoke 2019). These vortices are long‐lived mesoscale anticyclones (hundreds to 1,000 km in diameter) trapping plumes of aerosols and combustion/volcanic compounds. Owing to their unusual composition, these anticyclonically trapped plumes (ATPs) are associated with a significant radiative heating, which fuels their ascent through the stratosphere. This article investigates the dynamics of ATPs using two complementary approaches: analytically, in a potential vorticity (PV) perspective, and through idealised numerical simulations with the Weather Research and Forecasting (WRF) model. In both cases, we consider the vortical flow forced by a heating Lagrangian tracer. By reformulating the problem in the potential radius–potential temperature coordinate system introduced for tropical cyclones, we first clarify that ATP formation is concomitant with the injection of air into the stratosphere at extratropical latitudes. Then, we derive a set of simplified one‐dimensional equations describing the subsequent evolution of the flow after the injection. The equation obtained for the tracer is a variant of the classical Burgers' equation. In qualitative agreement with the three‐dimensional WRF simulations, this theoretical model predicts that ATPs develop an upper tracer front associated with sustained near‐zero anticyclonic PV, followed by a smooth tracer tail of cyclonic PV. Radiative relaxation of the temperature perturbations induced by the anticyclone and the presence of an initial PV anomaly tend to stabilise ATPs during their ascent. Finally, we note that the theory predicts a similar relationship between the plume and anticyclonic PV for cooled ATPs, which is supported by three‐dimensional simulations and may apply to the 2022 Hunga volcanic plume.


INTRODUCTION
Whereas the flow in the midlatitude lower stratosphere is generally characterised by slow diabatic vertical transport (Brewer, 1949;Butchart, 2014), extreme events such as intense pyroconvection (Fromm et al., 2010) or large volcanic eruptions (Carr et al., 2022;Holasek et al., 1996) can result in a rapid and direct injection of air from lower levels into the stratosphere.Such injections of tropospheric air charged with combustion products or volcanogenic compounds have the ability to alter stratospheric composition and general circulation for years (e.g., Brown et al., 2023).
They also have significant impacts on the global radiative budget (Schallock et al., 2023;Senf et al., 2023).For a few years, it has been recognised that stratospheric plumes may feature specific mesoscale dynamics.A typical example is the plume from the 2019-2020 Australian wildfires, which self-organised into several long-lived coherent smoke patches.These ellipsoidal smoke "bubbles", which were typically a few kilometres deep and hundreds of kilometres wide, ascended through the midlatitude stratosphere while maintaining a compact structure (Kablick et al., 2020;Khaykin et al., 2020;Pumphrey et al., 2021).Atmospheric reanalyses have revealed that the persistence of those structures despite stirring by the large-scale flow resulted from the confinement of the plume within anticyclonic vortices.The largest of those anticyclone-trapped plumes (ATPs) was associated with a 1,000 km diameter vortex and remained visible as an ascending patch of ash for 3 months, rising from 16 to 36 km.Similar smoke anticyclones have also been identified retrospectively in the aftermath of the 2017 wildfires in British Columbia (Lestrelin et al., 2021).For both the 2020 and 2017 cases, diabatic lofting of the smoke plume, reaching several degrees of potential temperature per day (e.g., Boers et al., 2010;de Laat et al., 2012;Khaykin et al., 2018;Khaykin et al., 2020;Yu et al., 2019), was attributed to the absorption of solar radiation by black carbon aerosols.An analogous phenomenon may have taken place in the ash-rich aerosol cloud from the 2019 eruption of the Raikoke volcano (Khaykin et al., 2022b), which coalesced into horizontally circular patches (Cai et al., 2022) and ascended from 18 to 25 km in the course of a few months (Chouza et al., 2020).
More recently, the eruption of the Hunga volcano generated the largest and highest stratospheric injection since at least the eruption of Mount Pinatubo in 1991, with the main detrainment located around 35 km altitude (e.g., Carr et al., 2022;Khaykin et al., 2022a;Proud et al., 2022).This time, the plume consisted of sulphur compounds and water vapour (Khaykin et al., 2022a;Millán et al., 2022), and the latter generated substantial in-plume diabatic cooling and a rapid descent (Sellitto et al., 2022).Although this feature was not captured by the European Centre for Medium-Range Weather Forecasts operational analyses, Legras et al. (2022) noted using a suite of space-borne instruments that the plume tended to roll up and was associated with a wind anomaly consistent with an anticyclonic rotation.This signature was observed to last for several weeks following the eruption.
The role of vorticity anomalies in isolating air masses and preserving them from dilution within the environment is well understood in the case of "adiabatic" geophysical flows (e.g., Garny & Randel, 2013;Schoeberl et al., 1992).However, the coupling between the dynamics and diabatic forcing within a mesoscale stratospheric plume has received less attention.Consequently, the processes responsible for the formation and maintenance of ATPs remain unclear in the literature.Doglioni et al. (2022) simulated smoke-charged vortices in the Goddard Earth Observing System Chemistry Climate Model global model by introducing sunlight-absorbing aerosols.Guimond et al. (2023) characterised the sensitivity of global simulations of stratospheric wildfire plumes to model resolution.Though both studies considered the instantaneous tendencies induced by the heating, they did not investigate the mechanisms governing the continuous evolution of the plume.Lestrelin et al. (2021) carefully analysed the structure of Australian and Canadian smoke vortices using space-borne cloud-aerosol lidar with orthogonal polarization (CALIOP) observations and the European Centre for Medium-Range Weather Forecasts Reanalysis v.5 (ERA5).They found that the ellipsoidal anticyclone-trapped smoke plumes had a vertical-to-horizontal aspect ratio of the order of f ∕N, where f and N are the Coriolis and Brunt-Väisälä frequencies respectively.Interestingly, at the centre of the anomaly, ERA5 indicated near-zero absolute vorticity and potential vorticity (PV).However, the degree of realism of the dynamical fields associated with the plume in ERA5 is unclear, since they are not generated directly by aerosol-induced heating, which is not included in the model, but appear as a consequence of data assimilation.Lestrelin et al. (2021) proposed to interpret low absolute-PV smoke anticyclones as resulting from the vertical transport of a smoke-charged tropospheric air mass, conserving both its low PV and its tropospheric composition during its diabatic ascent through the stratosphere.This conjecture was used by Khaykin et al. (2022b) to explain the confinement and anticyclonic motion of the Raikoke volcanic plume.Nevertheless, this theory is incomplete from a dynamical point of view: by considering Lagrangian conservation of PV, it overlooks the contribution of diabatic heating to its evolution (e.g., Haynes & McIntyre, 1987;Hoskins et al., 1985).
The goal of this article is to investigate the dynamics governing the formation and maintenance of heated (and cooled) ATPs and in particular their connection with the dynamics of PV.We propose a conceptual model based on theoretical considerations and confront it with observations and numerical simulations.The article is organised as follows.In Section 2, we revisit the observational evidence for plume vortices and highlight particular properties.In Section 3, we review the theory of PV dynamics in diabatically forced axisymmetric flows and derive a set of simplified equations describing the joint evolution of PV and a diabatically active tracer.Section 4 presents results from numerical simulations of a heated tracer plume and compares them with the theory.Section 5 contains a discussion.Conclusions and recommendations for future work are outlined in Section 6.

OBSERVATIONS OF THE 2020 SMOKE-CHARGED VORTEX KOOBOR
Observations of the 2020 Australian wildfire plume are revisited in order to highlight some special properties of the aerosol cloud and associated vortex that have not been emphasised in previous studies.We use the 532 nm total attenuated backscatter data from the space-borne Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP) onboard the Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observation (CALIPSO) satellite (Vaughan et al., 2004;Winker et al., 2010), specifically version 4.51 (Getzewich et al., 2018;Kar et al., 2018) of the level 1 product.Global Navigation Satellite System (GNSS) radio-occultation (RO) temperature soundings are obtained from "dry retrievals" provided by the Constellation Observing System for Meteorology Ionosphere and Climate (COSMIC) Data Analysis and Archive Center (CDAAC).Other temperature and PV fields are derived from the ERA5 (Hersbach et al., 2020).

Structure and ascent of the Koobor plume and vortex according to CALIOP and ERA5
CALIOP observations of attenuated backscatter ratio (SR) were instrumental in the discovery of stratospheric plume vortices (Kablick et al., 2020;Khaykin et al., 2020).The main plume, nicknamed Koobor in reference to an Aboriginal legend (Lestrelin et al., 2021), could be followed over 3 months with the lidar.Figure 1 (top row) displays selected vertical sections of CALIOP SR across Koobor between January 9 and 18, 2020, together with (middle row) vertical sections of temperature zonal anomaly and (bottom row) horizontal sections of PV, both from ERA5.The SR shown here is the ratio of the measured total attenuated backscatter profile over the theoretical molecular backscatter (Vaughan et al., 2004).It appears artificially reduced within and below thick plumes due to the attenuation of the laser beam.
In agreement with previous studies, CALIOP observations in Figure 1 evidence an ascent of Koobor over time.According to ERA5, the aerosol cloud is accompanied by a temperature anomaly dipole associated with the anticyclonic vortex, with the cold anomaly located above the plume.The latitude-longitude track of CALIPSO's orbit is displayed in the bottom row of Figure 1 and shows a good correspondence between anticyclonic PV in ERA5 and the smoke cloud.During this early phase of its life, Koobor was horizontally stagnating over the Pacific (bottom panels of Figure 1) but ascending rapidly from 19 to 22 km over 9 days.

Peculiar properties of Koobor
Besides the ascent of Koobor and the associated vortex, we further note that the attenuated SR in Figure 1 (top panels) is not uniform within the plume.Instead, it is systematically larger near the top of the smoke cloud than in its lower part.This vertical gradient could arise from actual variations in aerosol properties, but also from the partial attenuation of CALIOP's laser beam.To disentangle the two contributions, we retrieved profiles of attenuation-corrected SR at the centre of the smoke plume for a subset of the sections in Figure 1.The lidar ratio (assumed constant) needed for the inversion is constrained for each profile using the attenuation of the molecular backscatter signal in the clear-sky stratosphere below the smoke cloud, following the method described in Fernald et al. (1972).The values obtained lie between 50 and 75, on the lower side of the range of lidar ratios reported by Ohneiser et al. (2020) for the 2020 Australian smoke plume.
The resulting profiles of attenuation-corrected SR are presented in Figure 2.They reveal large variations along the vertical, which demonstrates that the vertical gradients noted in Figure 1 are not only due to the attenuation of CALIOP's laser beam but also stem from actual heterogeneities in aerosol properties within the plume.Assuming that the optical and microphysical properties of the aerosol particles are homogeneous throughout the plume, the SR is proportional to the aerosol mixing ratio (in either mass or number).Here, SR profiles exhibit an asymmetrical, sawtooth-like pattern, with a sharp front at the top of the cloud and a smoothly sloping (and sometimes  perturbed) tail trailing behind.Furthermore, the bottom of the cloud seems to rise more slowly than its top, and the depth of the plume to increase with time.
Figure 2b shows profiles of temperature anomalies with respect to the zonal mean at the location of the SR profiles in Figure 2a, according to ERA5.The minimum of the temperature anomaly remains located at the front of the aerosol sawtooth, whereas the maximum gets broader with time as the aerosol tail develops.In this oceanic area void of radiosonde stations, ERA5 primarily relies on the assimilation of satellite observations of the temperature field.The capacity of the assimilation system to follow the observed temperature profiles is illustrated in Figure 3 for an encounter between Koobor and one of its siblings, denoted as vortex V2, which occurred on January 15, 2020.Koobor originated from an intense event of pyroconvection on December 29-31, 2019, whereas V2 formed after a similar case on January 4, 2020 (Peterson et al., 2021).Figure 3a shows the two clouds in a single cross-section of CALIOP, and Figure 3c the corresponding track of CALIPSO, which intersects the two vortices.In Figure 3b, ERA5 temperature profiles at the location of the cores of the two ATPs are compared with three nearby GPS RO (GPS-RO) soundings (their exact location is indicated in Figure 3c).The superposition of two anticyclones is associated with a double temperature dipole, which ERA5 follows with remarkable fidelity.Though the reanalysis accurately reproduces the temperature structure of the vortices, other meteorological fields, such as PV, may be less reliable.In particular, ERA5 reconstructs a single-signed, quasi-ellipsoidal PV anomaly matching the size of the aerosol cloud, but neither the exact shape of the vortices nor the distribution of the PV are directly constrained by the assimilated observations.Nevertheless, we note that, over its lifetime, ERA5 consistently indicates a near-zero PV at the centre of Koobor (Figure 3c, also see Lestrelin et al., 2021).This low absolute value of PV is very uncommon in the midlatitude stratosphere, but we will see in the following sections that it is predicted both by the theory and by the numerical simulations.

Potential vorticity budget of an axisymmetric flow
Given their typical scales, ATPs in rotating stratified fluids can be studied using tools from PV theory (Hoskins et al., 1985).As a first approximation, these structures can be considered axisymmetric, an assumption supported by the observations (see Figure 3).Isolated axisymmetric vortices are a long-studied topic in fluid mechanics.Regarding atmospheric applications, Eliassen (1952) introduced the model of a quasi-stationary axisymmetric vortex, in which the azimuthal wind field is in gradient wind balance with a radial temperature anomaly.This equilibrium state of the flow is sometimes referred to as the Eliassen balanced vortex model (Schubert & Hack, 1983).Eliassen's seminal work was later developed to investigate tropical cyclone dynamics from the point of view of PV (e.g., see, Schubert, 2018, and references therein).In the following, we will rely on a framework introduced for tropical cyclones by Schubert and Alworth (1987) and Schubert (2018), adapting it to the case of ATPs.

PV in geometric coordinates
The starting point is the dynamics of Ertel's PV (Ertel, 1942), noted P, which is expressed in Cartesian coordinates as where  is the density,  the potential temperature, and  a the absolute vorticity vector.In an inviscid (i.e., neglecting non-conservative mechanical forcing) adiabatic flow, PV is materially conserved and behaves like an inert tracer.On the contrary, in the case of vortices subject to diabatic heating or mechanical forcing, Lagrangian conservation of P is no longer guaranteed.In Cartesian coordinates, the PV tendency equation reads (e.g., Andrews et al., 1987) with Q = . = D∕Dt the diabatic heating rate and F the mechanical forcing.
In this article, we restrict ourselves to the case of an axisymmetric flow on the f -plane, where f is the Coriolis parameter.With r the radial coordinate and u and v the components of the horizontal wind in the radial and azimuthal directions respectively, PV may be expressed as (e.g., Schubert, 2018;Shutts, 1991) where we use the notation for the Jacobian determinant and R is the potential radius (e.g., Schubert & Hack, 1983), which is related to the absolute angular momentum and defined by Lagrangian changes in R arise solely due to drag in the azimuthal direction F  ; namely: . (5) Hence, R has the interesting property of Lagrangian conservation in an axisymmetric flow in the absence of mechanical forcing (F = 0).

PV in a transformed coordinate system
To simplify the equations further, a convenient change of coordinates consists in switching from r to R as horizontal coordinate (e.g., Shutts & Thorpe, 1978;Thorpe, 1985) and from z to  as vertical coordinate (e.g., Schubert & Alworth, 1987).Such coordinate transformation requires that there exists a bijective mapping from (r, z) to (R, ), which is guaranteed as long as (fP) > 0 (or simply P > 0 in the Northern Hemisphere); that is, as long as the flow is stable with respect to symmetric and convective instabilities.When P = 0 in a region of the fluid (symmetrically neutral flow), the area with strictly zero PV in (r, z) cross-sections collapses to a line in (R, ) coordinates, and the coordinate transformation is no longer bijective.This special case will be briefly discussed in the next subsection, but let us assume fP > 0 for the remainder of this subsection.
Recalling that the Jacobian determinant represents the change in differential volume under coordinate transformation, the infinitesimal volume element is given by and the corresponding mass is This relation allows us to define a potential pseudo-density (Schubert & Alworth, 1987), which is the density in the (R, ) coordinate system.From the expression of P in Equation (3), we have Hence, using the continuity equation in (R, ) coordinates, we obtain the PV tendency equation in this system: Another derivation of this equation, starting directly from Equation (2), is provided in Schubert and Alworth (1987).Note that the form of Equation ( 11) a priori guarantees thatP remains single-signed.

3.1.3
Extended PV impermeability theorem in (R, ) coordinates It is enlightening to put the equivalence between PV and potential pseudo-density , Equation ( 9), in the perspective of the PV impermeability theorem of Haynes and McIntyre (1987).Haynes and McIntyre (1987) showed that the (Eulerian) tendency equation of the quantity P =  a ⋅  may be recast into a flux-form budget equation in geometric coordinates (here, (r, z)): with J a flux vector.On this ground, they proposed an analogy between PV and chemical substances; namely, that the quantity P be seen as the concentration of a PV substance (PVS) of mixing ratio P. Haynes and McIntyre (1987) then proved that, contrary to fluid mass or chemical substances, there is no net flux J of PVS across isentropic surfaces even in the presence of diabatic heating and mechanical forcing.This result, named the impermeability theorem, follows solely from the definition of PV, and does not rely on the equations of motion.
For our axisymmetric problem in (R, ) coordinates, the analogous quantity to P is just P = f .Hence, the impermeability theorem is expressed trivially as This relation shows that J vanishes entirely in (R, ) space; or, in other words, there is no net PVS flux across either R or  surfaces.The mass of PVS contained within a volume delimited by iso-R and iso- surfaces remains constant in time.
3.1.4PV anomaly induced by an injection of air into the stratosphere According to Equations 9 and 13, injecting a mass of fluid within a volume delimited by iso-R and iso- surfaces is formally equivalent to creating anticyclonic (i.e., low) PV, whereas removing fluid produces cyclonic PV.Such a relationship between PV and mass redistribution in rotating stratified fluids is not new.Considering the Lagrangian conservation of  and R, Gill (1981) proposed to create anticyclones by injecting a mass of fluid at its level of neutral buoyancy in a rotating stratified flow, a method applied in tank experiments by Griffiths and Linden (1981) and a number of studies since then (e.g., Aubert et al., 2012).A similar line of reasoning was followed by Shutts et al. (1988), who used it to investigate the balanced state resulting from the inviscid and adiabatic adjustment of the atmosphere to penetrative convection.Here, we recall that this concept also applies to injections originating from volcanic eruptions-as discussed by Baines and Sparks (2005), albeit from a different point of view-and from pyroconvection.
More quantitatively, if we consider an injection of fluid in a background flow of potential pseudo-density  in (R, ) coordinates, the potential pseudo-density after the injection,  i , is increased due to the intrusion augmenting the mass of fluid contained within a given range of R and .The mass fraction of injected air relative to the total mass of air is equal to  i = ( i − )∕ i ; that is, to the relative density anomaly with respect to the background.Then, based on the correspondence between  and PV expressed in Equation ( 9),  i can also be written as a function of PV: where P() = f ∕ is the PV of the background state and P i the actual PV after the injection.Equation ( 14) provides a relationship between the initial mixing ratio of a plume tracer (black carbon aerosols, for instance), which is proportional to  i , and the relative PV anomaly with respect to P: that is,  i = −P an (R, , t = 0).

PV and tracer evolution in an inviscid flow
In the previous subsection, we proposed a possible mechanism for the initial formation of ATPs, namely the adiabatic and inviscid adjustment of an injection of air from lower levels into the stratosphere.We now consider the subsequent evolution of the plume, and wish to understand the ascent of the anticyclonic vortex together with plume aerosols and trace gas anomalies.Since our work focuses on the response to the diabatic heating induced by the plume, we will neglect mechanical forcing in the remainder of this section (i.e., .R = 0).

General considerations on the zero-PV plume
To reformulate the problem in a simple framework, we model the radiative effect of the plume by introducing a chemically inert Lagrangian tracer, which will be termed A. At this point, A might represent any plume tracer: black carbon aerosols, ash, water vapour, carbon monoxide, and so on.The evolution of the mixing ratio  of A is governed by an advection-diffusion equation: where  is a diffusion operator (linear in ).
Ignoring tracer diffusion as a first-order approximation, the right-hand sides of Equations 2 and 16 differ in general, due to the diabatic forcing term in the PV tendency equation, which is absent from the tracer tendency equation.In fact, a joint diabatic ascent of P and  anomalies can only take place if, everywhere where Q ≠ 0 and  ≠ 0, the PV cancels (in this case, this must be realised through  a = 0).Then, the irrotational frictionless flow remains irrotational, and PV remains zero even in the presence of diabatic forcing; see Equation (2).The Lagrangian tendencies of PV and tracer mixing ratio hence coincide.Although this P = 0 condition is consistent with the value of PV found at the core of Koobor according to ERA5 (Section 2), it is unreasonable to conceive an isolated, strictly zero-PV tracer patch persisting over time.Viscosity and tracer diffusion will impose continuity of the PV and A fields between the patch and the environment.In the following, we investigate the evolution of an area of non-zero (albeit small) PV (i.e., a finite injection of fluid mass) associated with a tracer plume (i.e., an injection of tracer) in the (R-) framework.

3.2.2
Tracer-PV relationship in the case of a (non-zero) low-PV intrusion Assuming .R = 0 and neglecting tracer diffusion across iso-R surfaces, the general solution to Equations 11 and 16 may be obtained by solving, separately for each R, a couple of one-dimensional (1D) partial differential equations (with respect to ).This set of 1D solutions concatenated along R is then the solution of the two-dimensional (2D) problem.In each 1D domain, the Lagrangian tendencies of PV and  are given by where Q  = (Q∕) R and it is now understood that D∕Dt = ∕t + Q∕.As noted by Schubert and Alworth (1987), the solution to Equation (17a) along material (Lagrangian) air trajectories (which are known as characteristic curves in the terminology of nonlinear waves) reads and highlights exponential growth/decay of the PV with time  for constant Q  (with P 0 the initial PV value along the trajectory).
Before considering a more specific form for Q and relating it to , we note that, when neglecting tracer diffusion ( = 0 in Equation ( 16)), where   = (∕) R .This relation is valid for any Q and any perfectly Lagrangian (inert and non-diffused) quantity.By integrating Equation ( 19) in time  along characteristic curves, we obtain that is, that PV evolves as the inverse of the tracer gradient.Equation ( 20) shows that PV decreases following the flow when the tracer gradient becomes steeper in  coordinates, and vice versa.

Flow response to a heating (or cooling) tracer
Let us now consider the case in which the tracer A heats (or cools) the air.From this point on, A stands for a radiatively active species; for instance, black-carbon aerosols (which absorb sunlight) or water vapour (which increases infrared emission and cooling).Our goal is to determine how the tracer field evolves, and how this affects the PV.To the best of our knowledge, this configuration has never been considered.

Heating rate formulation
We will assume that the heating rate expressed as a temperature tendency  is proportional to the mixing ratio  of A; that is: where  is a constant factor.For the wildfire plumes, A consists of black carbon particles heating the air and  > 0. Equation ( 21) yields the following expression for the potential temperature heating rate Q: where strictly speaking the temperature field T depends on time, R, and , but may in practice be replaced by its background profile T().
An alternative to Equation ( 22) would be to replace  by Q directly in Equation ( 21).Though this choice is slightly more convenient to work with, Equations 21 and 22 constitute a better approximation encapsulating the linear relationship between  and tracer mixing ratio in the case of constant and uniform radiative properties of the tracer.
Note that, in Equation ( 21), any attenuation of the radiative fluxes, which could arise in thick plumes, is neglected.In the case of the plumes from the 2020 Australian fires, this appears as a reasonable approximation, since CALIOP shows limited attenuation of the lidar beam at 532 nm (Section 2).Furthermore, Sellitto et al. (2023) recently reported that the radiative properties of the aerosols remained stable over time-scales ranging from weeks to months, which justifies using a constant .Finally, it should be mentioned that, besides simplifying the heating induced by the absorbing/emitting layer, Equation ( 21) also ignores the radiative damping of the temperature anomaly associated with the anticyclone (Lestrelin et al., 2021).This component is left out for now, but it will be included in Section 4.

Reformulation of the tracer equation into Burgers' equation
With this form of the heating rate, and after a change of the vertical coordinate to Lagrangian conservation of the tracer mixing ratio in the (R, Θ) coordinate system is expressed by a form of Burgers' equation: where on the right-hand side we have introduced a simple parametrisation of  with K Θ ≥ 0 a constant "diffusivity" in Θ coordinates.Burgers' equation is a canonical example of a nonlinear wave equation, and it has been extensively studied (e.g., Cole, 1951;Hopf, 1950).For details, we refer the reader to standard textbooks on partial differential equations (e.g., Whitham, 1999).In Appendix A, we present briefly the general properties of the solution to Burgers' equation initialised with a localised patch of tracer, sometimes referred to as the "single hump" solution.With K Θ = 0, a singularity occurs at a finite time  s , indicating an instability in the physical system leading to potential mixing.For small K Θ > 0, the solution is regularised and converges asymptotically towards a "triangular wave": a sharp front, which translates itself and decreases in amplitude over time, and is followed by a linear tail.Appendix A also contains the analytical solution for the corresponding PV profile, governed by Equation (17a).The latter aspect is specific to the present study.

3.3.3
Behaviour of the PV-tracer system We now propose to illustrate graphically the dynamics of the system.To this end, we computed the solution of the PV and tracer equations (Equations 17a and 23).The background flow is assumed at rest, so that the reference density is simply  = ∕(∕z) r ; we assume an isothermal profile at T = 215 K, representative of the extratropical stratosphere (see Figure 3).The initial distribution of the tracer  i is taken equal to the fraction of the air originating from the injection-see Equation ( 14)-and is designed to be smoothly varying in R and  using the cosine function: in unit of kg ⋅ kg −1 (kg of plume air per kg of air).P i follows from Equation ( 14).Motivated by the main 2020 Australian vortex Koobor, we take guidance from ERA5 on January 9, 2020, to set the magnitude and size of the injection in (R, ) coordinate.This suggests R min = 25 km, R max = 350 km,  min = 360 K, and  max = 460 K (not shown).The amplitude, (1 − )∕4, with  ≪ 1, is such that the maximum value of the tracer distribution (reached for R ∈ [0, R min ] and  = ( min +  max )∕2) is  i = 1 −  and corresponds to a nearly vanishing PV (P i = P(1 −  i ) = P). represents the mass fraction of environmental air mixed with the plume at the core of the intrusion.When the pseudo-density of the background  is low enough at the detrainment level,  may be tiny even for a limited total amount of injected fluid.This point will be discussed further in Section 5.In practice, we use  = 5 × 10 −2 .

T
Figure 4 displays the time evolution of the tracer and PV profile at the centre of the ATP (R = 0).Over 37 days of simulation, the symmetrical, smooth initial bulge of the tracer profile gradually transforms itself into a strongly asymmetrical triangular wave.At the top of the cloud, the tracer slope initially steepens until it stabilises as a sharp tracer front.With this parametrisation of , the instantaneous ascent speed of the front after its formation at t =  s is  max ∕2 in the Θ coordinate ((∕T) max ∕2 in the  coordinate).Owing to the decay of the tracer jump  max as t −1∕2 (Whitham, 1999), this ascent speed decreases gradually with time and the front asymptotically rises as √ t in the Θ coordinate.At the bottom of the cloud, the tracer gradient is continuously decreasing, leading to the formation and deepening of a long tail (Figure 4a).The resulting sawtooth pattern resembles the observed aerosol backscatter profiles shown in Figure 2.
There are two effects contributing to the formation and displacement of PV anomalies.The first is the upward advection of low PV in the background PV profile, which always results in a negative (anticyclonic) PV anomaly, since PV increases rapidly with  in the stratosphere (Figure 4b).The second contribution is the Lagrangian diabatic tendency, Equation (17a).This tendency leads to a decrease of PV right at the front, where it is brought to zero exponentially with time (linearly in the absence of diffusion), and to an increase of PV elsewhere in the tracer cloud, as explained in Appendix A. Consequently, a near-zero absolute PV minimum (recall that P ≥ 0) ascends together with the tracer front.The minimum PV value decreases in magnitude over time, both in relative terms (since P increases with , P − P decreases) and absolute terms (in agreement with Equation ( 19), P decreases).In the tail, diabatic and advective PV tendencies are of opposite signs.A region of anticyclonic anomaly, originating in the vertical advection of both the low-PV initial condition and the background gradient, is initially maintained just below the front.It gradually shrinks, and the diabatic tendency ultimately prevails, leading to a cyclonic anomaly in the tracer tail.
Figure 5 shows a 2D representation of the set of tracer and PV solutions for a range of R and selected times.In each vertical column (fixed R), the tracer and PV profiles exhibit an evolution qualitatively similar to that seen in the central column and described earlier herein.This evolution entails the ascent of a tracer front accompanied by a low-PV anomaly, and a cyclonic tracer tail behind.Quantitatively, the ascent speed of the front decreases with increasing R, due to the decrease with R of the maximum tracer concentration  max in the initial condition.
Overall, given its degree of idealisation, it is noteworthy that this analysis predicts features that are consistent with observations of ATPs; in particular: • the joint ascent of a heating tracer and anticyclonic PV anomaly; • the formation of near-zero PV and its maintenance within the tracer cloud; • the formation of a tracer front at the top of the tracer cloud.
However, there is no observational evidence for the cyclonic tail, and the slowing down of the ascent rate of the plume over time appears faster than observed.Also note that the sharpness of the front is likely enhanced in  coordinates compared with altitude coordinates, since the absolute PV (and hence ∕z) is lower in that region.This effect cannot be assessed in the current framework.Moreover, some of the assumptions underlying the present theoretical development deserve further scrutiny.Although supported by the observations, the preservation of the initial axisymmetry of the flow is not a priori guaranteed, since the non-diffusive tracer equation has the intrinsic tendency to produce singularities.Another weakness of the present approach is its lack of flexibility; for instance, to incorporate a more realistic diabatic forcing.Indeed, it is expected that the heating rate in ATPs includes a significant contribution from radiative relaxation (Lestrelin et al., 2021).Estimating this component within the (R-) framework would require inverting the temperature anomaly from the PV field (e.g., Thorpe, 1985), and this would not address the question of flow stability.As an intermediate framework between analytical approaches and real-case simulations, we therefore turn to fully three-dimensional (3D) idealised numerical simulations, in which radiative relaxation can be readily included and the validity of our assumptions directly tested.

3D NUMERICAL SIMULATIONS
In this section, we study the response of the flow to an injection of the active tracer A, as in Section 3, but including more complete dynamics.The validity of our theoretical model is assessed using a realistic numerical model which, although still set in an idealised configuration, solves a more exact set of governing equations for the atmosphere.This allows us to bridge our theoretical results with the observations.The set-up may be summarised as follows: we consider a background atmosphere initially at rest on the f -plane.The chosen Coriolis parameter (f = 10 −4 s −1 ) and background stratification profile mimic the typical midlatitude troposphere-stratosphere transition, with a tropopause located at 12 km altitude and associated with a stratification jump in background Brunt-Väisälä frequency N from N = 10 −2 s −1 to N = 2 × 10 −2 s −1 .At the initial time, a low-PV "intrusion" is placed in the stratosphere, collocated with a patch of the tracer A. As in the previous section, we investigate the midterm (tens of days) evolution of the flow starting from this configuration.

Model choice and set-up of the numerical experiments
The 3D numerical simulations have been performed with the Advanced Research Weather Research and Forecasting (WRF) model (Skamarock et al., 2008), Version 4.2.2.The model is a community model, designed for limited-area simulations of mesoscale flows.It is used operationally for weather forecasting, but it also supports idealised configurations that can be flexibly used for research purposes (e.g., Bui et al., 2019;Foussard et al., 2019).Two key advantages of this model are its realism and broad range of applications, and hence of validation.It is important to note, however, that there is a significant distance between The fraction of total mass of the tracer A contained within regions where its mixing ratio is larger than 0.5 kg⋅kg −1 (F >0.5 ), right after the initial adjustment and after 18 days of simulation.Thus, F >0.5 is a proxy of the dilution of the tracer plume.
the highly idealised situation considered in the previous section and the flows that are simulated with WRF, which retains the full complexity of the compressible stratified atmosphere.Choices need to be made for the boundary conditions in the simulations.In addition, as expected in any numerical simulation, a limited resolution will induce numerical diffusion, with possibly significant implications for derived quantities such as PV.Hence, it is not the purpose of these simulations to reproduce exactly the solutions of the idealised model studied in Section 3, nor the observations.Rather, they aim at identifying which features of the theoretical analysis are robust and shared with the 3D simulations, and which are not.The domain extends from the surface to about 55 km altitude and spans about 3,800 km in the zonal and meridional directions.The horizontal grid is regular in Cartesian coordinates, with the consequence that the initial axisymmetry of the flow cannot be strictly preserved.Hence, even in the absence of azimuthal perturbations at the initial stage, non-axisymmetric instabilities have the potential to develop, and the simulations may be used to assess the stability of the axisymmetric vortex.The horizontal resolution is 15 km, whereas the (161) vertical levels are approximately regularly spaced in altitude, corresponding to a vertical resolution of 300 m between 10 and 38 km (the range of altitude of interest).Our idealised set-up has open lateral boundary conditions, a free-slip lower boundary condition, and a 15 km deep sponge layer in the upper part of the 55 km deep domain.Since we are interested in the dry stratospheric flow, no microphysics parametrisation is included.
To represent radiative processes acting during the maintenance stage of the plume, we use a simple approach and consider the heating tracer A, but also include radiative damping in the form of Newtonian relaxation.Again, the mixing ratio  of A is conserved following the flow, except for numerical and explicit diffusion, and follows Equation (17b), but with this time () standing for an (anisotropic) 3D diffusion operator, which includes both parametrised and numerical contributions.In the default model configuration, a horizontal Smagorinsky first-order closure (Smagorinski, 1963) handles diffusion in the horizontal, whereas a second-order diffusion scheme with constant diffusivity equal to 0.1 m 2 ⋅s −1 is employed in the vertical.Newtonian cooling dampens deviations of the temperature from its reference profile T(z).The full diabatic heating rate  = .T = (T∕)Q is given by where  is the heating (or cooling) rate per unit of tracer mixing ratio, and  is the temperature relaxation rate.
In the case of the smoke vortices, the first term () is positive and represents the short-wave (solar radiation) effect, whereas the second term ((T − T(z))) represents the radiative relaxation resulting from the perturbed balance between long-wave (infrared) absorption and emission in the presence of a temperature anomaly.Notice that, in our simulations,  is uniform and does not depend on the tracer concentration.Table 1 summarises the simulations presented in this article; all of them started from an axisymmetric initial state.Different values of  were tested.In the smoke vortex case (UP-REF simulation), it is chosen to roughly match the observed ascent rate of the Koobor plume-as seen, for example, in Khaykin et al. (2020)-at initialisation.The radiative damping rate is estimated from ERA5 to be of the order of  = 1∕7 day −1 (Lestrelin et al., 2021).

Model initialisation
Ideally, we would like the initial state to be identical to the one used in Section 3 and given in (R-) coordinates by Equation ( 24), but the initialisation is complicated by the transformation from (R-) to model coordinates.In In all panels, white contours are isopleths of tracer A (0.1, 0.5, and 0.9 kg⋅kg −1 ).
practice, similar initial conditions are achieved in model coordinates by starting from an atmosphere at rest (i.e., u = v = 0).In this case, r = R and  is equal to the isentropic density in r- coordinates (i.e. = z∕).Profiles of pressure, geopotential, and potential temperature profiles can be derived from (r, ) by integrating the hydrostatic balance equation numerically in each column of the model.We set  = ∕(1 −  i ) to be consistent with Equation ( 14), with  i given by Equation ( 24), noting again that R = r.
Instead of starting from a state of rest and letting the flow adjust to reach the cyclo-geostrophically balanced state corresponding to the initial geopotential and temperature perturbations, we could have directly initialised the simulation with a balanced vortex structure in (r, z), determined from the PV distribution in (R, ) coordinates using an inversion approach similar to that of Thorpe (1985).Besides its simplicity, our procedure has practical advantages; notably, it is consistent with the discretised equations of the model, and it can be easily applied to more complex background flows, such as sheared flows.This flexibility comes at some additional computational cost, as a small model time step (30 s) is required for numerical stability during the adjustment process.
In the simulations, the adjustment triggers a substantial emission of gravity waves.After a spin-up of about 2 days, the balanced state that is left consists of an upper level anticyclone, in agreement with the theoretical study of extratropical penetrative convection by Shutts et al. (1988).All times in the following are defined with respect to the end of the 2-day spin-up.
Figure 6 depicts cross-sections of the dynamical fields on day 1, at the end of the adjustment.In agreement with expectations for an elliptical anticyclonic vortex in gradient wind balance, the flow consists of an anticyclonic wind torus of magnitude ∼15 m⋅s −1 sandwiched within a temperature anomaly dipole of amplitude ±5-6 K.The size of this vortex is comparable to that of Koobor (Figure 1), with a vertical depth H ≃ 4.3 km and a diameter L ≃ 1, 000 km.The aspect ratio H∕L ≃ 4.3 × 10 −3 is close to the value predicted by Charney's quasi-geostrophic (QG) scaling (f ∕N = 5 × 10 −3 ), although the vortex, with its vanishing PV, is very far from geostrophy.Hassanzadeh et al. (2012) and Aubert et al. (2012) proposed another formula for the aspect ratio of an elliptical vortex embedded in a Boussinesq fluid at rest with background Brunt-Väisälä frequency N and in gradient wind balance with the environment.It is given by where  v and N v are respectively the vertical component of the absolute vorticity and the Brunt-Väisälä frequency inside the vortex.For a low absolute PV vortex ( v ≃ 0), Equation ( 26) can be rewritten as with T v ≃ 215 K the (unperturbed) temperature at the middle of the vortical structure, and Horizontal cross-sections of potential vorticity anomaly (colours) and tracer mixing ratio (white contours: 0.1, 0.3, 0.5, 0.7, 0.9 kg⋅kg −1 ) at the level of (top row) the tracer front Z front and (bottom row) the barycentre of the tracer field Z centre , as follows from the UP-REF simulation. ( is the lapse rate anomaly within the vortex.The value estimated from Equation ( 27), namely H∕L = 4 × 10 −3 , is in good agreement with the model's result.
The diabatic PV tendency induced by the heating field at the start of the simulation is shown in Figure 6c.It features a vertical dipole with a negative PV tendency above and a positive PV tendency below the heated area, which is consistent with theoretical expectations-see Equation (11) and Hoskins et al. (2003)-and with ERA5 data assimilation tendencies around Koobor (Lestrelin et al., 2021).This pattern initially leads to an upward translation of the PV anomaly (Figure 6b) that accompanies the ascent of the tracer.However, we cannot prefigure the later stages solely on this basis: the flow is essentially nonlinear and, since the ascent takes place against stratification, it involves a sort of convection.We will now investigate the actual evolution of the structure.

4.3
Flow response to a heating tracer: simulation of ascending smoke-charged vortices

Evolution in physical space
The UP-REF simulation is run for 40 days following the spin-up.The first question we wish to investigate concerns the maintenance of axisymmetry.In Figure 7, we present successive horizontal cross-sections of the tracer and relative PV anomaly P an , defined by Equation ( 15), at the level of the plume core (centroid of the  field) and the level of the plume front (maximum of |d∕dz| at the centre of the domain).The figure shows that the flow remains axisymmetric throughout the entire simulation.This property justifies the use of vertical cross-sections to provide a full account of the dynamics.
Figure 8 shows such vertical cross-sections for the azimuthal wind and the temperature anomaly T ′ = T − T(z) (top row), and for the relative PV anomaly (bottom row), at various times spanning the duration of the simulation.Let us first consider the first 13 days (first two columns in Figure 8).On day 13, the low PV and tracer patches have risen together by ∼3 km.The tracer distribution, which is initially vertically symmetrical, has gained skewness with a steepening of the vertical tracer gradients near the top of the plume and a flattening in its tail.At this stage, the pattern of the tracer mixing ratio already resembles the SR cross-sections through Koobor shown in Figure 1.The formation of a tracer front and tail is qualitatively consistent with the theory presented in the previous section.As expected, the evolution of the PV field is already significantly different from that of the tracer.Though the PV anomaly has mostly remained of a single sign, and the initial anticyclonic PV patch, like the tracer plume, has both ascended and deepened, we also note the creation of slightly cyclonic PV anomalies, and the enhancement of the anticyclonic PV anomaly near the tracer front.The wind and temperature anomalies related to this PV structure have been both translated upward and elongated vertically.The anticyclonic wind ring has retained a roughly constant magnitude.Regarding temperature, the cold anomaly is strengthened, whereas the underlying warm anomaly has become deeper and weaker with time.Such short-term evolution of the temperature field is generally compatible with that seen in ERA5 near the Koobor vortex (Figure 1).
At later times (last columns in Figure 8), the rising tracer cloud maintains a roughly constant horizontal Azimuthal wind (colours) and temperature anomaly with respect to the background (black contours).Bottom row: Relative potential vorticity (PV) anomaly from the background P an (colours) and potential radius R (light green contours).The area of negative absolute PV (P an < −100%) is dotted in red.White contours are isopleths of tracer A (from 0.1 to 0.9 kg⋅kg −1 every 0.2 kg⋅kg −1 ).
extent, with only a slight lateral expansion (in the order of 14% between day 1 and day 37).Along the vertical direction, the tracer front undergoes a significant ascent (from ∼16 to 26 km) whereas the base remains steady, such that the plume depth has increased by a factor of 3 between day 1 and day 37. Whereas the bulk of the tracer remains located within an anticyclone over the entire duration of the simulation, a smaller and weaker cyclonic PV anomaly develops in the lower part of the tracer cloud, and gradually deepens and intensifies.At the same time, the overlying patch of anticyclonic PV is not simply being translated upward, but also changes shape, evolving from a vertically symmetric ellipsoid to a strongly asymmetric pattern.Thus, on day 37, at the end of the simulation, anticyclonic PV is found at the top front of the tracer cloud in the central part of the domain, but also further down around the edges of the plume (about 6 km below the central tracer front), laterally surrounding the cyclonic PV anomaly (Figure 8).The overall PV field at the end of the simulation looks like a "hollow tooth", with a cyclonic PV core and an anticyclonic shell above and around it.As a result of this peculiar PV pattern, the cyclonic anomaly is partly shielded by the lateral anticyclonic PV anomaly, which tends to reduce signature on the wind and temperature fields.Hence, although this PV structure is strikingly different from the anticyclonic PV monopole diagnosed by ERA5 (Khaykin et al., 2020;Lestrelin et al., 2021), the discrepancy appears less pronounced for the meteorological fields accessible to the observation and assimilated in the reanalysis.Furthermore, we expect that missing dynamical ingredients, such as shear of the background wind, would promote the evacuation of the cyclonic vorticity in the tail.
Despite the presence of this cyclonic anomaly, it seems visually that most of the tracer remains within the anticyclone, avoiding dilution within the environment.In order to quantify the dilution of the plume, we computed the fraction of the total mass of tracer contained in regions with mixing ratio of the tracer A larger than 0.5 (kg⋅kg −1 ), denoted as F >0.5 .F >0.5 is a proxy for plume dilution.At the beginning of the simulation, we have F >0.5 = 78%; after 18 days it is reduced to F >0.5 = 45% (see Table 1).This relatively weak dilution is related to the deflection of the air surrounding the anticyclone during the ascent of the latter.This is illustrated in Figure 8, where the green contours represent isolines of the potential radius R (i.e., isolines of angular momentum), which approximately coincide with actual streamlines (not shown).
Overall, although the flow shown by the WRF simulation is more complex, we note qualitative similarities with Section 3 regarding the evolution of the tracer and PV fields, with the formation of a vertically skewed tracer cloud at the vortex centre and an accompanying asymmetric vertical PV dipole, with an anticyclone lying above a small cyclonic tail.Though part of the flow develops a negative absolute P < 0 (at the edge of the hollow tooth in particular, dotted in red in Figure 8), most of the domain has a single-sign absolute PV P > 0. This property and the preservation of axisymmetry on the horizontal make it possible to analyse the PV and A fields in (R-) coordinates.Together with the qualitative agreement with Section 3, this motivates us to compare the theory and simulation in more detail using the transformed coordinates.

4.3.2
Evolution in (R, ) space Figure 9 displays the same snapshots of the UP-REF simulation as Figure 8, but with the fields first averaged along the azimuthal direction and then transformed to (R-) coordinates (following the constant R surfaces depicted in Figure 8).If radiative relaxation and drag ( .R) were completely negligible, we would expect the evolution in (R-) to be entirely consistent with the theory of Section 3. Striking qualitative similarities are indeed found between the UP-REF simulation (Figure 9) and the theoretical model (Figure 5), showing the relevance of the latter to this more realistic case.
At any given R, the tracer profile along  (Figure 9) evolves into a shape similar to the triangular wave seen in Section 3(Figure 5).In agreement with the theory, which predicts an ascent speed proportional to  max , the varying ascent speed of the tracer front at different potential radii R reflects the initial value of the tracer maximum  max (R) in each column.Having a closer look at the PV, Figure 9 shows that, in the short term (up to day 13), the PV anomaly field tends to remain negative everywhere and is characterised by an ascending anticyclonic PV patch with little cyclonic counterpart, as expected from the maintenance property of an initially irrotational flow (characterised by P ≃ 0).Later on, the PV profiles become vertically asymmetrical.PV exhibits a pronounced vertical minimum at the tracer front (near −100% relative PV anomaly) located within a broader area of anticyclonic PV lying over a positive PV anomaly in the tail.Overall, the anticyclonic front and cyclonic tail are qualitatively consistent with Figure 5. Switching from r to R coordinate furthermore reveals that the lateral ring of anticyclonic PV in r-z coordinates (Figure 8) is actually always located above the cyclone R- coordinates (Figure 9), which fits within the theoretical picture.
There are also noticeable differences between UP-REF and the theoretical model.Some are related to changes in the initial WRF model fields compared with the initialisation of the theoretical model, Equation ( 24), and are already present on day 1 (Figure 9).The changes include a homogenisation of  max (R) along R for R ≤ 100 km and the generation of two local maxima along  for 100 ≤ R ≤ 250 km.Both are most likely due to the unavoidable mixing and friction (i.e., diffusion of momentum) occurring during the adjustment.At later times, the evolution of these initial differences can nevertheless still be understood in the framework of the theory.In particular, the two local maxima in the tracer mixing ratio generate the transient formation of two tracer fronts, each of which is developing its own low-PV filament (around day 13).As the lower front is characterised by a higher  max and a faster ascent speed, it closes in on the upper front, and the two entities eventually merge.
Other differences between UP-REF and the theoretical model integration cannot be attributed to the initialisation and are not as easily reconciled with the theory.First, the decrease in the tracer maximum  max is faster, and hence the ascent rate (proportional to  max ) slower in the theoretical model than in UP-REF.Then, although both Figures 5 and 9 highlight a PV dipole, the anticyclonic PV anomaly gets quickly (15 days) confined to the front in the case of the Burgers' solution but remains spread out around the front in the 3D simulation.The respective roles of radiative damping and viscous drag in creating those differences are discussed in the following subsection.

Radiative relaxation
Compared with the UP-REF simulation, the theoretical solution presented in Section 3 overestimates the extent and magnitude of cyclonic PV anomalies.It also shows a more localised anticyclonic PV, as well as excessive tracer dilution.We expect that these discrepancies partly result from the radiative relaxation introduced in the UP-REF simulation and ignored so far in the theoretical derivation.The potential mechanism is illustrated in Figure 6d.By redistributing heating over a deeper layer and decreasing the heating gradient in the area of the strongest tracer gradient, radiative relaxation tends to produce a negative PV tendency above and below the plume, beyond the area where the tracer lies, thereby increasing the spatial extent of the localised negative PV tendency induced by the tracer heating only.By smoothing and spreading out the heating field, radiative relaxation also tends to slow down the sharpening of the tracer front, thus partly shielding the tracer cloud and limiting its dilution during the ascent.
In order to quantitatively assess the impact of radiative relaxation, we first attempted to include this term in the heating Q used in the reduced theoretical model of Section 3.This new 1D model, solving Equations 17a and 17b with a more complete Q, is described in Appendix B. A significant addition is the estimation of radiative relaxation, which requires inverting the temperature anomaly distribution from the PV field.For internal consistency, this is achieved without any information from the UP-REF simulation, but instead using a QG ansatz and assuming a typical horizontal scale  h ( h = 250 km) for the vortex (see Appendix B).To enable a quantitative comparison, the updated 1D model is initialised from the axial profile of tracer and PV of the UP-REF simulation on day 1 (i.e., after the 2-day spin-up).
Figure 10 compares the evolution of the tracer, PV, and temperature profiles at the centre of the domain in the UP-REF simulation with integrations of the 1D model for an adjusted set of parameters ( = 5 K⋅day −1 ⋅kg −1 ⋅kg), 1∕ = 5 days, K  = 17 K 2 ⋅day −1 ,  h = 250 km; see Appendix B).An excellent agreement is achieved regarding the shape of the tracer profile.Compared with the 1D simulation in Section 3 (Figure 4), including radiative relaxation effectively slows down the decay of tracer maximum, leading to an increased ascent rate of the triangular wave, improving quantitatively the comparison with UP-REF.The distribution of PV also depicts a good quantitative agreement with UP-REF, in particular regarding the sharpness of the PV minimum.Nevertheless, we note that the PV jump below the front is more pronounced in the 1D model.The temperature anomalies (Figure 10c,f) are well aligned at the initial time.However, the amplitude of the negative anomaly grows over time in the UP-REF simulation whereas it is decaying in the 1D model (by a factor of 2 over the course of the integration).This is likely related to the assumption of a constant horizontal scale  h in the QG inversion, whereas the UP-REF simulation actually shows a lateral expansion of the anticyclone with time (Figure 8).Despite those remaining disagreements, this reduced 1D model simulation is closer to the WRF UP-REF simulation.Along with more extensive sensitivity calculations with the 1D model, it suggests that increasing the temperature relaxation rate  slows down the decay of the maximum tracer mixing ratio  max near the front, although without preventing its formation.This is in turn responsible for a higher vertical ascent rate of the front and the associated PV minimum (approximately proportional to  max ∕2).
In order to confirm the results of the 1D model, we also performed another numerical experiment with the WRF model in which radiative relaxation is switched off (simulation NO-RELAX).The early evolution of the flow is qualitatively similar to UP-REF.However, over time, the area of anticyclonic PV becomes significantly shallower and the positive PV significantly larger in absolute value (not shown).In agreement with the 1D model results, the dilution of the tracer in NO-RELAX is increased compared with UP-REF: after 18 days, F >0.5 is reduced to 36% in NO-RELAX (vs.45% in the presence of relaxation; Table 1).As already discussed, this may be explained by the shielding of the plume by the relaxation-induced part of the Q field in UP-REF, which is lacking in NO-RELAX.At later stages, this smaller tracer mixing ratio leads to a slower ascent of the plume, as in the 1D model.Overall, both the 1D and 3D numerical experiments point to an important role of the long wave radiative transfer in the evolution of anticyclone-trapped plumes, which may preserve them from dilution as they ascend.Here, we only included a simplistic parametrisation of radiative relaxation with a uniform damping parameter , which only accounts for the temperature anomaly.The impact of the complex mixture of compounds radiatively active in the infrared contained in wildfire plumes, including both black carbon aerosols and greenhouse gases (water vapour, etc.), or of a depletion in other species (ozone), remains to be assessed.

Friction
Besides radiative relaxation, some differences between the theory and the UP-REF simulation must be a consequence of the diffusion of momentum.Indeed, the reduced 1D model relies on the assumption of an inviscid flow ( .R = 0), which is clearly violated in WRF, both due to explicit and implicit (numerical) viscosity.
There are specific features of the 3D simulation that may be attributed to the diffusion of momentum generating a cross-R flow.In particular, the top of the tracer front exhibits a protrusion towards large R (Figure 9).This protrusion corresponds to the tracer ring above the main plume in Figure 8.Such a feature cannot occur from an inviscid evolution of the initial condition.We interpret this protrusion as resulting from a drag occurring at the top of the vortical flow where the vertical derivative of the balanced azimuthal wind is the largest, leading to positive .R = (r∕R)F  .As a proxy for potential instabilities and drag, we examined the distribution of areas of negative absolute PV (dotted in red in Figure 8) and of low Richardson number (not shown), which indeed emphasise both the shallow front and the lateral tail of the tracer cloud as regions with likely significant diffusion/mixing.
The impact of friction may be grasped quantitatively by considering the integral of Equation ( 10) along : only mechanical forcing can result in mass fluxes across iso-R surfaces ( .R ≠ 0).Over the course of the simulation, the mass of air enclosed within R < 100 km and 360 <  < 800 K is reduced by 30%, demonstrating significant viscous drag effects that are not captured in the 1D framework.This effect is even larger for the flux of A, emphasising either an enhanced exchange of mass through in regions of high tracer content, or the existence of cross-R diffusion of A. Overall, this suggests a significant role of parametrised and numerical (e.g., associated with the Cartesian geometry of the grid) drag in the evolution of anticyclonic plumes, an effect that cannot be accounted for in a purely 1D framework.This process may be associated with dynamic instabilities, but they do not generate a breaking of the axial symmetry at the resolution of the simulations.

Isolation of the plume from environmental air: Impact of the initial PV anomaly
A special property of the ascending ATPs is their sustained isolation from their environment.Potential radiative instabilities caused by the heating may have been expected to result in significant mixing between the plume and the environment.We have seen that this is not the case, for two reasons.First, the presence of a low PV anomaly is associated with streamlines avoiding the tracer cloud, thus preventing its dilution within background stratospheric air (Section 4.3).Second, the additional heating rate provided by the long wave radiative relaxation of the temperature perturbation associated with the PV anomaly shields the tracer cloud and reduces mixing with outside air.Without an initial PV anomaly, there is no temperature anomaly, and this shielding is absent at the beginning of the simulations.
Thus, we expect that the co-injection of a heating tracer and of air, resulting in a plume initially trapped within a low absolute PV anticyclone, plays a key role in preventing the dilution of the plume along its ascent, through the effects of both dynamical isolation in a low PV vortex and radiative relaxation.However, a few previous studies (e.g., Doglioni et al., 2022;Guimond et al., 2023) have reported on simulations of ATPs that did not include an initial PV anomaly.In order to assess the actual impact of this parameter on the evolution of the plume, we performed a simulation including the same initial tracer field as UP-REF but not the initial PV anomaly; that is,  i =  (simulation UP-NOPV).As shown in Figure 11, an anticyclone is still forming at the level of the tracer front in this simulation.However, it only encompasses a small fraction of the tracer plume after 2 weeks.As a result, the dilution of the tracer is larger in this simulation: after 18 days, F >0.5 = 24% (compared with 45% in the UP-REF simulation, 36% in the NO-RELAX simulation).In contrast with UP-REF, the initial axisymmetry of the plume is not maintained in UP-NOPV.On the contrary, we observe a destabilisation of the initially axisymmetric field.Whereas there seems to be a relaxation towards axisymmetry at later times, the anticyclonic plume trapping is still weak on day 19 of this simulation.
The difference between cases with and without initial mass injection (i.e.,  i = ∕(1 −  i ) and  i = ) is reduced when we consider an absorbing tracer injection below the tropopause, between 10 and 12 km, as used by Doglioni et al. (2022) (not shown).Although a breaking of the axisymmetry of the flow still occurs, large anticyclonic PV anomalies are created by the diabatic tropopause crossing, and they tend to merge more rapidly (within less than 2 weeks) than when the tracer is directly introduced in the stratosphere.
In realistic dynamical settings, including vertical shear and large-scale deformation, we expect that the initial presence of low PV anomalies collocated with the tracer will prevent its rapid dilution both below and above the tropopause and enable the rapid formation of a compact plume.Indeed, even if distributed among a few patches, we speculate that such low PV anomalies favour the emergence of a single coherent plume through vortex merger (e.g., Dritschel, 2002;Reinaud & Dritschel, 2002).In the case of the Australian fires, several (about 10) large pyrocumulonimbii reached the stratosphere (Peterson et al., 2021) and probably later self-organised into a small number of ATPs, including Koobor.

Environmental conditions for anticyclonic plume formation: Altitude and latitude of the injection
In Section 3 and 5.1, we argued that the (cross-isentropic) injection of an air mass near its level of neutral buoyancy was responsible for the initial formation of an ATP.From a dynamical point of view, neither the vertical gradient of PV nor the difference between tropospheric and stratospheric PV play a direct role, since Lagrangian conservation of PV does not apply.However, for a similar amount of air being injected, the smaller the background isentropic density  is, the larger the background PV and the relative PV anomaly are.Given the background profile of PV (see Figure 4), it comes with no surprise that anticyclonic plumes are more likely to occur higher up in the atmosphere, where the fluid is lighter and more stratified.This likely explains why ATPs have not been reported after tropospheric injections.
A second factor regarding the initial formation is the latitude of the injection.This important parameter has been ignored so far in our investigations adopting the F I G U R E 11 Horizontal cross-sections from the UP-NOPV simulation of potential vorticity anomaly P an (colours) and tracer (white contours: 0.1, 0.3, 0.5, 0.7, 0.9 kg⋅kg −1 ), at the level of (top row) the tracer front and (bottom row) the tracer maximum.
f -plane approximation.We expect that the efficiency of the anticyclonic confinement depends on the ratio of the relative vorticity (equal to −f if P ≃ 0) to the change in planetary vorticity between the centre and the edge of the structure.In polar and midlatitude regions, encountering quasi-zero PV is extremely rare, and the P ≃ 0 vortex provides strong isolation.In contrast, in tropical regions, ambient PV is closer to zero and the isolation provided by low PV is limited.
The non-dimensional number comparing the magnitude of latitudinal vorticity variations to the local vorticity anomaly (≃ −f in a zero-PV vortex) is (L)∕|f |, where  = (2Ω∕a) cos() (a = 6, 400 km is the Earth radius,  the latitude, and Ω ≃ 7.27 × 10 −5 s −1 the Earth rotation rate) and L is the horizontal scale of the vortex.Assuming that the typical horizontal size of the structure L is of the order of the Rossby radius of deformation R d , as suggested by its aspect ratio, we have For an injection of depth H ≃ 1 km in the stratosphere (N ≃ 2 × 10 −2 s −1 ), the condition (L)∕|f | < 0.2 ≪ 1 is satisfied poleward of about 20 • .Closer to the Equator (L∕|f | > 1), the efficiency of the trapping is reduced as the latitude shift that is necessary to reach P ≃ 0 on an isentropic surface diminishes.Hence, confinement within an anticyclone is weaker in tropical regions-provided it occurs.The dynamical response to heating (or is affected by the equatorial beta effect, as suggested for the evolution of the Hunga plume by Schoeberl et al. (2023).We leave investigations of this interesting aspect for future dedicated studies.

5.3
The case of a cooling tracer (descending water-rich/ozone-poor plumes) Up to now, we have been focusing on the case of a heating tracer, which was motivated by observations of stratospheric wildfire plumes.A volcanic plume reaching deep in the stratosphere can detrain large amounts of water vapour, as anticipated by, for example, Glaze et al. (1997), and recently observed in the aftermath of the January 2022 Hunga eruption (Khaykin et al., 2022a;Millán et al., 2022).Following the latter event, Sellitto et al. (2022) suggested that radiative cooling from the injected water vapour is likely to strongly dominate the diabatic forcing.Indeed, evidence of compact descending patches of Hunga volcanic compounds was presented in Legras et al. (2022).Motivated by this example, we explore the response of the stratospheric flow to a plume of tracer cooling the fluid.Contrary to the heating tracer case, however, we do not intend here to explain a specific event.Instead, we consider the hypothetical case of a cooling tracer injected in the midlatitude stratosphere at 35 km altitude.The cooling rate per unit tracer is now  = −5 K⋅day −1 ⋅kg −1 ⋅kg; that is, opposite to the heated case.The resulting cooling rate remains within the range of plausible values; for instance, I G U R E 12 Successive cross-sections from the DOWN-REF simulation of (top row) wind (colours) and temperature anomaly (black contours) and (bottom row) potential vorticity anomaly from the background (colours) and momentum radius R (light green contours).White contours are isopleths of A (from 0.1 to 1 every 0.1 kg⋅kg −1 ).Sellitto et al. (2022) estimated the cooling rate in the Hunga plume to be of the order of .T = −10 K⋅day −1 at the day +3 stage.
Figure 12 shows the time evolution of the flow in (r, z) coordinates.Similar to the ascending case, the tracer profile develops a front and a tail, although their vertical orientation is reversed.A fundamental difference with the ascending case is the orientation of the motion with respect to the background PV (and density) gradient.This results in the faster dilution of the tracer in denser air in (R, ) space, Equation ( 16), and its contraction in physical space.Quantitatively, F >0.5 decreases from 70% after the adjustment to 32% after 18 days of simulation, slightly faster than in the ascending case (Table 1).Owing to its downward advection into a lower PV environment, the magnitude of the anticyclonic PV anomaly reduces over time, whereas that of the positive PV anomaly grows.The effect of PV advection is here opposite to the ascending case and limits the development of the anticyclonic anomaly.Nevertheless, this experiment shows that a cooled vortex may also persist during its descent and form a long-lived ATP.

SUMMARY AND CONCLUSIONS
Diabatically forced ATPs were discovered in the extratropical stratosphere following the extreme 2019-2020 Australian wildfires.The largest of these Australian ATPs, a 1,000 km diameter, radiatively heated vortex nicknamed Koobor, was observed to ascend from 18 to 35 km.More recently, the plume of the 2022 Hunga eruption underwent a spectacular diabatic descent due to the cooling induced by an extreme water vapour anomaly.In that case, satellite observations suggested the presence of a vorticity anomaly collocated with the tracer plume.By hindering the dilution of the plumes, vortical structures increase their longevity and modify their impact on radiation, chemistry, and dynamics.From the point of view of geophysical fluid dynamics, these ATPs exhibit several peculiarities, such as a virtually zero PV maintained at the centre of the 2020 Koobor vortex according to ERA5.
In this study, we adapted a theoretical framework developed for tropical cyclones (Schubert & Alworth, 1987) to investigate the dynamics of potential vorticity in diabatically forced axisymmetric anticyclones on the f -plane.Based on an extension of the impermeability theorem of Haynes and McIntyre (1987), we argued that the creation of anticyclonic PV in ATPs occurs concomitantly with the injection of a mass of air directly into the stratosphere by pyroconvection or volcanic eruptions.The formation of large ATPs may additionally involve the merging of a number of distinct initial anticyclones, through vortex merger.
The subsequent evolution of an ATP over time-scales ranging from weeks to months is controlled by radiative heating or cooling within the plume.At first order, radiatively active species, such as sunlight-absorbing aerosols or infrared-emitting water vapour, can be represented by a chemically inert tracer heating or cooling the air at a rate proportional to its mixing ratio.If the active tracer is initially located within an isolated patch of exactly zero absolute PV, the Lagrangian tendencies of PV and tracer become the same, and we expect their coincidence to be maintained over time.Considering the more realistic case of a low, but non-zero initial PV collocated with the tracer, the problem is greatly simplified by switching from cylindrical coordinates to angular momentum-potential temperature coordinates.Then, neglecting radiative relaxation, the diabatic motion of the tracer reduces to a variant of Burgers' equation, whereas PV passively follows the tracer evolution.The dynamics of Burgers' equation leads to the formation of a vertically asymmetrical tracer cloud, with a sharp front ahead of the plume (i.e., at its top for heated vortices) and a smoothly sloping tail at its rear.In the case of a heating tracer, this pattern is reminiscent of CALIOP SR cross-sections through Koobor.Near the tracer front, anticyclonic PV is sustained, whereas the tail is associated with cyclonic PV.Though the joint ascent of the anticyclonic PV and tracer plume is consistent with the observations, the cyclonic tail predicted by this theory is not present in the case of Koobor according to ERA5.
To assess the validity of the theory, we performed more complete numerical simulations with the WRF model, including Newtonian relaxation of the temperature anomalies in the heating rates.The simulated flow remained axisymmetric, and showed qualitative agreement with the theoretical picture regarding the evolution of both tracer and PV.In particular, despite the larger extent of the anticyclonic PV, which encompasses most of the tracer plume, a cyclonic PV anomaly still develops in the tail of the plume.This cyclone is smaller and weaker than the overlying anticyclone, and appears partly shielded by the anticyclonic PV anomaly, which may partly explain why it is not detected in the observations.Furthermore, a sensitivity simulation demonstrated that the inclusion of Newtonian cooling (radiative relaxation) tends to weaken this cyclonic anomaly.Radiative relaxation of temperature perturbations also stabilises the anticyclonic structure and reduces the mixing of the plume with environmental air.Consistently, introducing a low PV anomaly together with the tracer at the beginning of the simulations was found to favour the persistence of a compact plume and to preserve it from dilution as it ascends through the stratosphere.Compared with the inviscid theory, the WRF simulations also emphasise the role of diffusion of momentum, which becomes significant for the long-term evolution of ATPs.The effect of mechanical forcing can be directly quantified in (R, ), and this aspect deserves further dedicated investigations.
Several processes were ignored in our set-up.Notably, vertical shear of the horizontal wind may contribute to the evacuation of the cyclonic tail and the maintenance of the ellipsoidal shape of the anticyclone.Variations of the Coriolis parameter influence the latitudinal drift of the anticyclone and the evolution of tropical plumes.Other crude simplifications in our investigations concern the radiation and, for absorbing particles, aerosol microphysics.The parametrised radiative heating/cooling rate did not vary in time, and depended linearly on the tracer mixing ratio.Although we do not expect diurnal variations of the heating rates to have a strong impact on the (overall-balanced) plume, nonlinearities in the heating function (e.g., saturation above a given tracer concentration) could affect the dynamics.Regarding microphysics, the vertical motion of the tracer was assumed purely Lagrangian, which is well adapted to gases but fails to describe the sedimentation of aerosol particles over the course of a few months.The respective roles of these different processes should be assessed in future work.
Extreme wildfires are expected to become more frequent in the coming decades due to climate change.Meanwhile, there are important stakes around understanding volcanic perturbations of the stratospheric aerosol layer due to the need to assess the potential effects and likely risks of stratospheric geoengineering.In this context, it is crucial to better constrain the dynamical behaviour of stratospheric volcanic and wildfire plumes as a prerequisite to reliably quantifying their long-term impact on ozone chemistry, stratospheric circulation, and global climate.
(multivalued tracer profile) occurs as characteristic curves intersect and the PV reaches zero.We note that, since this collapse of the tracer cloud in R- coordinates corresponds to a PV P ≲ 0, a marginally unstable flow with respect to symmetric and/or convective instability develops.Characterising the exact evolution of this potentially unstable state can only be done in r-z coordinates and is beyond the scope of our study, as it would require high-resolution simulations.Instead, the focus here is on the general behaviour of the system and on stabilisation mechanisms (e.g., radiative relaxation).
To regularise the tracer equation and prevent the formation of the singularity in practice, we adopt the commonly used approach to introduce a non-zero diffusivity K Θ > 0, albeit small (compared with ∫  dΘ), yielding Equation ( 23).This procedure may be seen as a relaxation towards a stable state, or as a parametrisation of the mixing expected from instability development.
Formally, for K Θ → 0, it is still possible to construct physically acceptable ("vanishing viscosity") discontinuous solutions to Burgers' equation (Whitham, 1999).If the initial condition is a localized patch of tracer, the tracer profile asymptotically converges towards the "triangular wave" solution, which exhibits a linear tail with a steady base and an upper discontinuity.The magnitude of the discontinuity decreases as and the singularity initially located at Θ 0 f ascends along Θ at a rate dΘ f ∕dt = ( max ∕2); that is, as Note that the relations in Equations A.1 and A.2 remain valid everywhere except right at the singularity, and in particular on the lower side of the tracer cloud (the fraction of the cloud upstream of the jump and where (  ) 0 > 0).

A.2 Frontal PV dynamics in the presence of tracer diffusion
For a non-vanishing but small and constant diffusivity K Θ > 0, Equation ( 23) is the viscous Burgers' equation and the equivalent of the triangular wave is known as the "single hump" solution (Whitham, 1999).The discontinuity is replaced by a tracer front of width ΔΘ ≃ K Θ ∕( max ).The front location and magnitude correspond to that of the jump Θ f expected in the K Θ → 0 limit-that is, in the case of the triangular wave solution -and do not depend on the exact value of K Θ .Furthermore, following Equation ( 18), the PV decreases exponentially in the area of the front (Q  < 0).An illustration of this behaviour is provided in Section 3.3.3using numerical solutions of the PV and tracer equations.Here, we investigate further the case of a steadily propagating tracer front, which differs from the single hump solution but is more convenient, as both the tracer and PV equations admit explicit analytical solutions.
The steadily propagating front solution to the viscous Burgers' equation is given by the function In Equation A.4,  max is the tracer maximum (height of the front), and the integration constant Θ 0 f corresponds to the position of the front at t = 0. Here, the front has width K Θ ∕( max ) and it is travelling at constant speed  max ∕2 in  coordinates.Assuming the tracer profile given by Equation (A.4), the hyperbolic PV Equation (17a) may be solved along its characteristic curves (here, material trajectories) for each initial value  0 of  = where we have further assumed that the background profile is isothermal at a constant temperature T; that is,  =  ref exp(Θ∕T).As shown in Figure A1, the characteristic curves gather in the frontal region of the tracer distribution due to the Characteristic curves of the potential vorticity equation for  = 5 K⋅day −1 ⋅kg −1 ⋅kg, K Θ = 5 K 2 ⋅day −1 ,  max = 0.5 kg⋅kg −1 and a front initially at  0 f = 500 K.(b) Evolution of the potential vorticity ratio P∕P 0 along the characteristics for these conditions.advective effect of the heating, and PV decays exponentially to zero there, as demonstrated by Equation (A.7).

APPENDIX B. COUPLED TRACER-PV 1D MODEL INCLUDING RADIATIVE RELAXATION B.1 Model equations
In order to assess the role of radiative relaxation, we built a 1D model solving the coupled system made of Equation (17a) for the PV and Equation (17b) for the tracer and including the full heating rate Q.The tracer equation includes a diffusion term accounting for the background density gradient and is given by the following: where we use K  = 17 K 2 ⋅day −1 .The sensitivity to K  is limited for the range of values tested; in particular, it does not affect the decay of the tracer maximum, or only very marginally.The value of K  for the tracer was adjusted to reproduce the width of the front observed in the WRF simulation.
In practice, the evolution of the PV field is computed using the flux-form budget equation of the potential pseudo-density  = f ∕P, which reads

B.2 Estimation of the temperature anomaly
Contrary to the case examined in Section 3, it is no longer strictly valid to consider the model columns independent.Indeed, knowledge of the radial PV distribution is required in order to invert the local temperature profile and estimate the contribution of radiative relaxation to the heating rate, (T − T).Here, for simplicity, we assume a constant horizontal scale  h for the structure, and use a QG ansatz to derive a direct formula relating the temperature anomaly T ′ = (T − T) to the perturbation PV profile.Under the QG approximation (Berrisford et al., 1990), the full PV is approximated as where  = −1∕g(p∕) x,y is a mean vertical profile that depends only on , and q is the QG potential vorticity, which is defined as Based on this, we can derive the temperature anomaly from the anomaly P ′ = P − P of the PV along the axis for which q = P ′ .The assumed horizontal scale  h is used to replace the horizontal Laplacian in Equation (B.3); that is, ∇ 2  M = −(1∕ 2 h )M.For the sake of simplicity, we also assume a uniform background temperature T = 215 K in the vertical domain of the vortex.The explicit formula for the perturbation temperature on the axis is then where

B.3 Choice of the parameters and sensitivity study
Note again that the temperature anomaly obtained from this QG inversion only takes the central PV profile into account and is not expected to reproduce exactly the actual evolution of the temperature field in the WRF simulation.The value of  h (250 km) was chosen in order to reproduce the depth of the temperature anomaly at the initial time of the simulation.The slight lateral expansion of the vortex as it ascends is not represented in this constant- h framework, and the agreement between UP-REF and the 1D model temperature anomaly deteriorates over time, since an underestimated  h results in an underestimated amplitude of the temperature anomaly.This decrease of the temperature anomaly ahead of the front in the 1D model likely explains the need for an increased value of  (1∕ ≃ 5 days) to compensate this effect and improve the agreement with UP-REF (which has 1∕ = 7 days).
Regarding the sensitivity to model parameters, increasing the temperature relaxation rate  slows down the steepening of the front and the reduction in the value ofthe tracer maximum  max .As a consequence, the ascent of the front is accelerated for larger .Decreasing the horizontal scale parameter  h reduces the amplitude of the temperature anomaly and shallows its vertical scale.Though we restricted ourselves to tracer diffusivities K  small enough such that diffusion only marginally affects the results (in particular the decay of the tracer maximum or its propagation speed remain insensitive to this parameter), this parameter still controls the width of the tracer front.

F
Top row: Selected sections of cloud-aerosol lidar with orthogonal polarization (CALIOP) 532 nm attenuated backscatter ratio.The crosses indicate the altitude and latitude of the associated anticyclone in European Centre for Medium-Range Weather Forecasts Reanalysis v.5 (ERA5.Middle row: Corresponding height-longitude sections of ERA5 temperature anomaly at the latitude of the vortex centre, defined with respect to the zonal mean.The crosses indicate the altitude and longitude of the vortex centre.Bottom row: Horizontal charts of potential vorticity at the potential temperature level of the vortex centre (whose horizontal location is shown by the cross).The purple curves represent the track of CALIOP corresponding to the sections in the top row.F I G U R E 2 (a) Vertical profiles of cloud-aerosol lidar with orthogonal polarization (CALIOP) 532 nm (attenuation-corrected) backscatter ratio across Koobor at the latitude of the associated vortex centre for four night orbits between January 9 and 16, 2020, chosen for the proximity of CALIPSO's orbit to the vortex.The full cross-sections of attenuated backscatter ratio and exact times of the orbits are provided in Figure1.(b) Collocated ERA5 temperature zonal anomaly profiles.(c) As (a), but with potential temperature (inferred from the European Centre for Medium-Range Weather Forecasts Reanalysis v.5 [ERA5]) as the vertical coordinate.

F
I G U R E 3 (a) Section of cloud-aerosol lidar with orthogonal polarization (CALIOP) 532 nm total attenuated backscatter ratio on January 15, 2020, at 2103 UTC (time J).The two crosses indicate the latitude and altitude of Koobor (K) and the second vortex (V2) at time J in European Centre for Medium-Range Weather Forecasts Reanalysis v.5 (ERA5).(b) Vertical profiles of temperature at the location of the two vortices according to ERA5 (black curves) and to three nearby GPS soundings on January 15, denoted as GPS1 (KOMPSAT5, 1429 UTC), GPS2 (COSMIC2, 1522 UTC) and GPS3 (TSX, 1232 UTC).The selected profiles are located within 600 km of both vortices.The horizontal lines with crosses indicate the altitude of the two vortices (20.1 km for Koobor, 15.9 km for V2).(c) ERA5 potential vorticity maps of Koobor (colour) and V2 (contours) at the isentropic level of their centre (503 K for Koobor, 408 K for V2).The black lines indicate the trajectories of the two vortices in the time interval [J − 18 hr, J + 18 hr].The large crosses indicate the location of the two vortices at time J.The locations of the GPS radio-occultation soundings are shown as crosses (colours as in (b)), and the locations of the vortices at the time of CALIOP's section are depicted by diamonds.The pink curve represents the track of CALIOP corresponding to the section in (a).
Numerical solution of a set of Equations 23 and 17a concatenated along the R-dimension, with  = 5 K⋅day −1 ∕(kg ⋅ kg −1 ), K Θ = 5 K 2 ⋅day −1 , starting from the initial condition given in Equation (24) and a background PV profile typical of the midlatitudes.The top row shows the mixing ratio of tracer A and the bottom row shows the relative PV anomaly with respect to the background, P an .Similar conditions are used for initialisation and background in the reference (UP-REF) Weather Research and Forecasting model simulation.

TA B L E 1 Fa
List of idealized WRF simulations presented in this article and values of the parameters varied across the simulations a .The altitude of the tracer injection and the presence or absence of a potential vorticity (PV) anomaly associated with the tracer anomaly., the heating rate coefficient per unit tracer mixing ratio;  r , the temperature relaxation time-scale, Equation (25).b Simulations UP-REF and DOWN-REF are presented in Section 4. The sensitivity simulations NO-RELAX and UP-NOPV are presented in Section 5. c Dynamical fields from the UP-REF simulation on day 1.Cross-sections of (a) temperature anomaly and wind, (b) Ertel's potential vorticity (PV) through the adjusted vortex following the injection, (c) Lagrangian PV tendency due to diabatic heating, (d) Lagrangian PV tendency due to the radiative relaxation (corresponding to the long-wave [LW] component of the heating).
Successive cross-sections through the centre of the domain of different fields from the UP-REF simulation.Top row:

F
I G U R E 9 Successive cross-sections in R- coordinates of (top row) tracer mixing ratio  and (bottom row) relative potential vorticity anomaly from the background P an , as follows from the UP-REF simulation.

F
I G U R E 10 Profiles of (a, d) tracer mixing ratio, (b, e) potential vorticity, and (c, f) temperature anomaly.Panels (a)-(c) correspond to numerical integration of Equations B.2 and B.1 with the diabatic heating defined by Equation (25) and the temperature inverted using the formula provided in Appendix B with  h = 250 km.Panels (d)-(f) are average profiles for 0 ≤ R ≤ 60 km in the UP-REF simulation.
C p T ′ +  ′ is the anomaly in the Montgomery potential with respect to the basic stratification and  the background density profile.q satisfies