2.1 Met Office Unified Model (MetUM)

The Met Office global simulations are calculated using the Met Office Unified Model (MetUM; Cullen, 1993), coupled to the Joint UK Land Environment Simulator (JULES) model for the land surface (Best et al., 2011; Clark et al., 2011). The MetUM solves the deep-atmosphere, non-hydrostatic, compressible equations of motion in spherical geometry using a terrain-following vertical coordinate. The global model was run in its operational version (at the time of the case study in 2018) using the Global Atmosphere configuration (GA6.1; Walters et al., 2017), which includes the dynamical core named “Even Newer Dynamics for the General Atmospheric Modelling of the Environment” (ENDGame; Wood et al., 2014). In the vertical there are 70 terrain-following levels, with a fixed model lid at 80 km. More detail on the model formulation can be found in Sanchez et al. (2020). The global high-resolution operational forecast during October 2018 ran with 10 km average horizontal grid spacing. In near real-time, a limited-area simulation with 4.4 km horizontal grid spacing was nested within the global operational forecast, using the tropical version of the Regional Atmosphere and Land 1 (RAL1) configuration summarised by Bush et al. (2020) with reduced air–sea drag at high wind speeds. The developers of the RAL1 model were aware that the model overestimated the air–sea drag at high wind speeds (∼50 m s⁻¹) and planned to address this in the next model version (RAL2) for greater consistency with available observations (see Bush et al., 2020, for further discussion of this point). In this study, we use the RAL1 model with the air–sea drag reduction implemented, in other words analogous to RAL1+. There are 80 vertical levels whose spacing increases quadratically with height, up to a fixed lid 38.5 km a.s.l. (above sea level). Both simulations were initialised at 12:00 UTC on 21 October 2018 and run out to 5 d. Subsequently, the global forecast was re-run at lower resolution (N768, 18 km average horizontal grid spacing) to enable a much greater number of multi-level fields to be output than in the operational forecast. The global re-run and limited-area simulations are examined in this study.

2.2 Semi-geotriptic (SGT) balance approximation tool

Geotriptic balance involves a three-way static balance between Ekman friction, Coriolis and pressure gradient forces. Therefore, on its own it cannot predict what will happen next from any given state. Semi-geotriptic (SGT) dynamics predicts what happens next by describing the way in which the system evolves through a sequence of balanced states. An essential part of the evolution is the component of the 3-D velocity that is not in geotriptic balance but can be deduced from the pressure field, similar to the way that vertical velocity in the semi-geostrophic model can be deduced using the pressure (or geopotential) field and the semi-geostrophic omega equation (e.g. Hoskins and Draghici, 1977). This component will be called the balanced component of ageotriptic wind. The SGT equations are an extension of semi-geostrophic dynamics accounting for Ekman friction in the atmospheric boundary layer as part of the balance. The SGT equations are a good approximation to the MetUM equations on scales larger than the deformation radius, which implies aspect ratios less than $f/N$. Thus, they may be applicable to shallow synoptic-scale disturbances near the Equator as studied here. Semi-geotriptic dynamics outside the boundary layer implies that both the zonal and meridional winds are close to geostrophic balance. In other types of tropical dynamics, only the zonal wind is close to geostrophic. The SGT tool is implemented on the sphere using deep-atmosphere equations and variable Coriolis acceleration terms, consistent with the formulation and numerical discretisation of the global MetUM (Cullen, 2018). The use of the MetUM, rather than a reanalysis dataset, is motivated by the consistency between the model and the SGT tool in terms of numerical discretisation and geometry, as well as the absence of lateral boundaries that would complicate the inversion process. These properties motivate the use of the MetUM to drive our diagnostic work with the SGT tool.

Taking N768 global MetUM data as input, the SGT balance tool of Cullen (2018) is used to partition the global flow into balanced and unbalanced components. For example, above the boundary layer the full meridional wind can be decomposed:

where v_g is the geostrophic wind calculated directly from the horizontal pressure gradient, v_a,SGT is the component of ageostrophic wind consistent with semi-geostrophic balance dynamics and obtained from the SGT tool (as described below), and v_r represents the unbalanced residual component of motion. The zonal and vertical wind components can be partitioned similarly, noting that geostrophic wind has no vertical component (w_g=0). Within the boundary layer, Ekman friction is included to calculate geotriptic balance in place of geostrophic balance. By manipulating the three components of the momentum equation and making the geotriptic momentum approximation (where the momentum advected is approximated by the geotriptic wind (u_g, v_g) but the advecting velocity is not approximated), it can be found that

where the first term is a vector where the three components are the ageotriptic wind, ${\mathit{u}}_{\mathrm{a}}=(u_{\mathrm{a}},v_{\mathrm{a}},w)$, weighted by the product of a constant matrix B and matrix Q^′ which depends on gradients in geostrophic variables u_g, v_g and θ_V, as given in Eqs. (17) and (21) of Cullen (2018). The second term is the time tendency of the pressure gradient (where $\mathit{\pi }=(p/p_{\mathrm{0}}{)}^{\mathit{\kappa }}$ is the Exner pressure, c_p is the specific heat of dry air at constant pressure, $\mathit{\kappa }=R_{\mathrm{d}}/c_{\mathrm{p}}$, R_d is the gas constant for dry air and θ_V is the virtual potential temperature). H^′ is a vector forcing term including separate terms for the effect of diabatic heating, friction and geostrophic forcing of ageotriptic motion (see Sanchez et al., 2020, for the details of this term). Since the matrix B multiplies both the wind term on the left and the forcing on the right, the balanced ageotriptic flow obtained by inverting the equation can also be partitioned into that forced by large-scale geotriptic motion and diabatic heating:

This process is analogous to inverting the quasi-geostrophic omega equation to obtain vertical motion (e.g. Hoskins and James, 2014) except that the solution obtains the 3-D vector ageotriptic wind, rather than just vertical velocity. The SGT tool employs numerical discretisation of the governing equations consistent with the global MetUM, but with the geostrophic momentum approximation, namely that the momentum advected by the full wind is approximated by its geostrophic value (Hoskins, 1975). One key purpose of the approach is that the tool can calculate the response of the balanced flow to forcing associated with diabatic and frictional processes.

As well as the ageotriptic wind, the SGT tool solves for the tendency of the Exner pressure associated with balanced motion, given only the pressure field, diabatic heating and frictional forcing as input. This property means that in principal the pressure field can be stepped forwards in time and then the updated pressure used to diagnose the next time-step values of u_g, v_g, θ_V and the balanced ageotriptic wind. This sequence enables calculation of updated matrix Q^′, and the process can be repeated. Therefore, this SGT system is a complete model of the balance dynamics with only one prognostic variable – here formulated in terms of pressure. Those familiar with quasi-geostrophic or semi-geostrophic dynamics may expect to use PV as the single prognostic variable of the balance model because of its conservation property, but in semi-geostrophic theory PV is only materially conserved if the Coriolis parameter is regarded as constant. No such approximation is made here, which is vital for a global solution. The SGT tool enables investigation of the degree to which the vortex can be understood in terms of balanced dynamics and the role of diabatic heating and friction in the structure, intensification and maintenance of the vortex (in a way that cannot be achieved without the concept of balance).¹

The pressure input to the SGT tool must always satisfy the condition that the equation solved in the SGT diagnostic is elliptic, which implies that PV must be positive. This condition will not always be satisfied by higher-resolution MetUM data, particularly near the Equator. Therefore, the MetUM data passed to the SGT tool are interpolated to a coarser horizontal resolution and smoothed. In addition, a zonal filter is applied within a latitude band either side of the Equator, both to the input data and to geostrophic winds as part of the solution process. The values are filtered towards a zonal mean close to the Equator, in order to ensure that the wind and temperature fields do not vary in the zonal direction on the Equator, which would violate SGT balance. This zonal filter represents a linear blend of the zonal mean and the full 3-D field, in which the coefficient of blending varies with latitude. Heavy smoothing is applied between the Equator and a defined half-width (ϕ₁=4^∘), with lighter smoothing between ϕ₁ and a defined total width (ϕ0=8^∘). An additional 2-D smoothing filter, defined by a Gaussian convolution function applied isotropically in latitude and longitude between the Equator and ϕ₀, ensures that the pressure and geostrophic wind are smooth at scales appropriate to the balance approximation. The kernel used in the convolution function has a Gaussian half-width of 3.5 grid cells (each grid cell in the tropics is approximately 130 km in length). This additional filter removes unrealistic vorticity anomalies along the boundaries of the zonal filter region at 8^∘ either side of the Equator.

The degree to which the meridional wind from the MetUM is approximated by both the full ($v_{\mathrm{g}}+v_{\mathrm{a},\mathrm{SGT}}$) and the geostrophic meridional wind from the SGT tool for the case study analysed in this paper is shown in Fig. 1. The synoptic-scale pattern in the MetUM simulation is characterised by a wavelike pattern in the easterly flow across the South China Sea, with east–southeasterlies around 105^∘ E and two regions of east–northeasterlies either side (Fig. 1a). Although the flow pattern over Peninsular Malaysia is different, this wavelike structure is largely represented by both the full wind from the SGT tool (Fig. 1c) and the smoothed geostrophic wind (Fig. 1e). Longitude–height cross-section plots of meridional wind at 6^∘ N indicate that the SGT tool mostly captures the large-scale flow as represented by the MetUM (Fig. 1b, d and f). We will return to this case later in more detail.

Figure 1(a) Meridional wind (shaded; m s⁻¹) and horizontal wind (m s⁻¹; reference vector = 10 m s⁻¹) at 7 km, from the global MetUM simulation initialised at 12:00 UTC on 21 October 2018, valid at 12:00 UTC on 23 October 2018 (T+48); (b) longitude–height cross-section of meridional wind (shaded; m s⁻¹) along 6^∘ N. Overlain are the 850 hPa vorticity centre identified by the tracking algorithm (black star) and the height of the x–y section in (a) (dashed black line at 7 km); (c) and (d), as in (a) and (b) but for the full balanced meridional wind from the SGT tool, $v_{\mathrm{g}}+v_{\mathrm{a},\mathrm{SGT}}$, with diabatic forcing; (e) and (f), as in (c) and (d) but for the geostrophic component of the meridional wind, v_g, after smoothing and the application of the tropical filter.

Example output from the SGT tool and the driving global MetUM are compared across a much larger domain in Fig. 2 for another event in December 2018, characterised by twin tropical cyclones either side of the Equator over the Indian Ocean. The extratropical horizontal flow (poleward of 20^∘ N) is dominated by the SGT balance flow (cf. Fig. 2a and c) including a Rossby wave packet extending eastwards from China and the Korean Peninsula. Figure 2 demonstrates the ability of the SGT tool to partition the 3-D balanced ageostrophic flow into that forced by diabatic heating and geostrophic forcing. Vertical velocity in this northern part of the domain is largely driven by geostrophic forcing (Fig. 2b and f). A notable exception is the ascent ahead of a mid-latitude trough on the west side of Japan which is amplified by heating (Fig. 2d). The unbalanced residual component, v_r, which represents the difference between the wind from the MetUM simulation and that output by the SGT tool, highlights the flow features in the MetUM that are not captured by SGT balance flow (Fig. 2e). This term has relatively small amplitude poleward of 20^∘ N, except the area near Japan next to the heating already mentioned. However, the meridional component of the unbalanced residual wind, v_r, is larger in the tropics and changes sign at about 10^∘ S (converging from the north and south at this latitude). Given the time of year, this pattern is likely a signature of the low-level convergence of the Hadley circulation. Note that there are also large zonal flows in the unbalanced residual converging on Sumatra at this time (Fig. 2e). Zonal variation in wind along the Equator cannot be described by SGT balance, including equatorial wave dynamics. The strong circulation around the twin tropical cyclones is only partially represented by the SGT balance, with a smaller unbalanced residual on their poleward flanks and a large residual in between (since the equatorial westerlies between them cannot be part of the SGT balance). Also, the balance is not good nearer the middle of the tropical cyclones because semi-geostrophic balance is inaccurate when trajectory curvature is much tighter than the Rossby radius of deformation. However, later Fig. 10 shows that semi-geostrophic balance can describe qualitatively the rotational flow around the Borneo vortex (except equatorward of about 2–3^∘ N) and the vertical motion in balance with latent heating on the synoptic scale.

Figure 2(a) Meridional wind (shaded; m s⁻¹) at 2 km from the N768 MetUM simulation initialised at 12:00 UTC on 11 December 2018, valid at 12:00 UTC on 13 December 2018 (T+48). (c) Meridional component of the full balanced wind, $v_{\mathrm{g}}+v_{\mathrm{a},\mathrm{SGT}}$, and (e) unbalanced residual meridional wind, v_r, at 2 km from the SGT tool at 12:00 UTC on 13 December 2018; (b), (d) and (f) Vertical velocity (shaded; cm s⁻¹) at 2 km inverted from the SGT tool at 12:00 UTC on 13 December 2018 with (b) full forcing, (d) diabatic forcing only, and (f) geostrophic forcing only. The wind vectors represent the full horizontal wind in (a)–(c), the forced ageostrophic component in (d) and (f), and the unbalanced residual component in (e), with the reference vector in the upper left corner for all panels.

2.3 Verification data

To verify the representation of rainfall within the vortex by the MetUM, Integrated Multi-satellitE Retrievals for Global Precipitation Measurement (GPM-IMERG; Huffman et al., 2019) satellite data are analysed. The product used is Level 3 Final Run Precipitation, which combines precipitation estimates from GPM satellites (see https://gpm.nasa.gov/missions/GPM/constellation, last access: 12 October 2023) and Global Precipitation Climatology Centre precipitation rain gauges and has a horizontal grid spacing of 0.1^∘ and an output interval of 30 min. The product combines infrared geostationary satellite data with passive microwave radiometer data to produce the best observational precipitation data available in each 30 min interval.

The ERA5 reanalysis dataset (Hersbach and Coauthors, 2020) is used to further verify the ability of the MetUM to represent the 3-D structure of the vortex. The ERA5 dataset has an output interval of 1 h, a spectral horizontal resolution of TL639 (approximately 31 km grid spacing at the Equator), and 137 hybrid sigma–pressure levels in the vertical between the surface and 1 Pa.

2.4 Vortex tracking

The vortex event chosen for analysis herein was identified by tracking 850 hPa relative vorticity at T63 resolution at 6 h intervals in ERA-Interim reanalysis data (Dee and Coauthors, 2011). The dataset was produced as part of the Newton Fund project under the auspices of the WCSSP Southeast Asia project by Dr. Kevin Hodges of the National Centre for Atmospheric Science and Department of Meteorology, University of Reading. The algorithm used to track vortices in this study has been used previously on a range of synoptic-scale and mesoscale vortices including midlatitude cyclones (e.g. Priestley et al., 2020), tropical cyclones (e.g. Manganello and Coauthors, 2016) and polar lows (e.g. Zappa et al., 2014), as well as the Borneo vortex (Liang et al., 2021). In the case of Borneo vortices in this study, disturbances with relative vorticity greater than $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{6}}$ s⁻¹ that last for more than 1 d during the period October to March are retained.

5.1 Three-dimensional Borneo vortex structure

As discussed in Sect. 2.2, a key purpose of the SGT balance approximation tool is that it enables investigation of the degree to which the vortex can be understood in terms of balanced dynamics and the role of diabatic heating in the structure, intensification and maintenance of the vortex.

When considering 3-D vortex structure and movement, it is important to characterise the large-scale flow. The vertical and horizontal shear in the background zonal flow are crucial ingredients, affecting system tilt and possible growth mechanisms. The mean flow during the case study (12:00 UTC on 21 October to 12:00 UTC on 26 October 2018) is generally easterly, when averaged over a longitude band covering the region of interest (Fig. 9a). Easterlies strengthen with height up to around 7 km, with a weaker gradient above this height. At 7 km, there is a broad maximum in easterly flow across 8–11^∘ N. There is a strong cyclonic shear on the southern flank of the easterlies across 5–7^∘ N. This easterly flow pattern is typical of the northeast winter monsoon (e.g. Chang et al., 2005a).

Figure 9(a) Mean zonal wind from the N768 MetUM simulation initialised at 12:00 UTC on 21 October 2018, averaged over time (12:00 UTC on 21 October to 12:00 UTC on 26 October 2018) and longitude (95 to 120^∘ E). (b) The quantity rN² (shaded; 10⁻⁷ s⁻²), averaged across 5.5–7^∘ N, combines specific humidity and static stability to identify changes between moist air with high static stability (large rN²) and drier air (smaller rN²), valid at 00:00 UTC on 22 October 2018 (T+12). Blue contours indicate the diabatically forced component of the balanced vertical velocity from the SGT tool (−2, 2, 6, 10 cm s⁻¹).

Download

A second key ingredient of the basic state is stability. Colour shading in Fig. 9b shows rN², the specific humidity multiplied by the static stability, $N^{\mathrm{2}}=(g/{\mathit{\theta }}_{\mathrm{0}})\partial \mathit{\theta }/\partial z$. It the next section it will be shown that the vertical gradient in this quantity can define a moist stability interface along which “diabatic Rossby waves” (e.g. Parker and Thorpe, 1995; Moore and Montgomery, 2004; Boettcher and Wernli, 2013) can propagate. The pale-blue contours show the diabatically forced vertical motion, w_diab, consistent with balanced motion deduced using the SGT tool. There is an obvious wave at the level of the strong gradient in rN² (4.5–5 km) and above. It is important that the moist stability interface occurs within the region of large-scale vertical shear, given that this vertical shear enables baroclinic interaction between the low-level vortex and mid-level wave.

In Fig. 10, the vertical velocity (w, shown at 7 km) and the horizontal wind vectors (shown at 2 km) are partitioned into the balanced flow response to diabatic heating ($\dot{\mathit{\theta }}$) and geostrophic forcing using the SGT tool, valid at 12:00 UTC on 23 October (T+48 in the driving global MetUM simulation) when the vortex circulation peaks in ERA5 (Fig. 5b). The qualitative similarity between w from the SGT tool and the 1 h accumulated heating from the driving global MetUM (Fig. 10a) shows that the regions of ascent identified by the SGT tool are representative of the flow in the MetUM, given that on scales larger than several hundred kilometres the primary balance in the thermodynamic equation in the tropics is between diabatic heating and vertical motion multiplied by static stability (as expressed in Eq. 2). The coherent region of heating to the south of Vietnam on the northwestern flank of the cyclonic circulation at 2 km is collocated with the maximum in total precipitable water and the area of organised precipitation associated with the vortex at this time (cf. Figs. 4c and 10a). An organised region of ascent is located in a similar region relative to the 2 km cyclonic circulation in the balanced flow from the SGT tool (Fig. 10b). This qualitative similarity indicates that the vortex is described by the balanced flow. Examining the partition of 3-D ageostrophic motion into that forced by $\dot{\mathit{\theta }}$ (Fig. 10c) and large-scale geostrophic forcing (Fig. 10d) reveals that upward motion within the vortex is primarily forced by $\dot{\mathit{\theta }}$. The 2 km ageostrophic wind vectors in Fig. 10c converge from the southwest and northeast underneath the heating in an asymmetric pattern, with the main inflow from the equatorward side. In contrast, the large-scale geostrophic forcing produces an anticyclonic circulation that opposes the geostrophic wind (Fig. 10d), although the difference in vector scaling should be noted. This result implies that the full flow is sub-geostrophic.

Figure 10(a) The 1 h accumulated heating at 7 km (shaded; K) and horizontal wind at 2 km (vectors; m s⁻¹) from N768 MetUM simulation initialised at 12:00 UTC on 21 October 2018, valid at 12:00 UTC on 23 October 2018 (T+48). (b–d) Balanced vertical velocity at 7 km (shaded; cm s⁻¹) and horizontal wind at 2 km (vectors; m s⁻¹) inverted using the SGT tool at 12:00 UTC on 23 October 2018 with (b) full forcing, (c) diabatic forcing only, and (d) geostrophic forcing only. In (a) and (b), the wind vectors represent the full horizontal wind. In (c) and (d), the wind vectors represent the ageostrophic wind only. Note that the scale of the vectors is different in all panels, with the key in the upper left corner.

The time sequence of wave propagation is shown in Fig. 11 over a 2 d interval. The wave grows to significant amplitude and therefore develops into pronounced meridional undulation in the eastward flow. The southward–northward velocity dipole moves westward together with the 7 km PV and w_diab maxima, the locations of which are marked onto the figure panels using symbols. The easterlies at 7 km also extend westwards with the wave. The vertical cross-sections are calculated along 6^∘ N. The asterisk just above the boundary layer (at 1.4 km) shows the diagnosed location of the vortex centre, and this centre is also marked on the 7 km maps. The PV maximum is close to the zero node of v, as expected for Rossby wave propagation, and at most times the maximum in w_diab is slightly shifted to the west. The upper-level wave propagates faster to the west than the lower-level vortex – this property will be explained in Sect. 5.3. Also, note how v is partitioned into two distinct layers, split approximately at 3–4 km. The vortex circulation is confined below this level and is approximately untilted up to 3 km a.s.l.

Figure 11(a) Meridional wind (shaded; m s⁻¹) and horizontal wind (m s⁻¹; reference vector = 10 m s⁻¹) at 7 km, from the global MetUM simulation initialised at 12:00 UTC on 21 October 2018, valid at 12:00 UTC on 22 October 2018 (T+24); (b) longitude–height cross-section of meridional wind (shaded; m s⁻¹) along 6^∘ N. Overlain are the 850 hPa vorticity centre identified by the tracking algorithm (black star), maximum PV at 7 km from the MetUM simulation (purple triangle), maximum vertical velocity forced by diabatic heating from the SGT tool at 7 km (blue diamond), and the height of the x–y section in (a) (dashed black line at 7 km); (c) and (d), as in (a) and (b) but valid at 12:00 UTC on 23 October 2018 (T+48); (e) and (f), as in (a) and (b) but valid at 12:00 UTC on 24 October 2018 (T+72). For reference, panels (c) and (d) are equivalent to Fig. 1a and b.

Figure 12c and d show the horizontal flow at 7 km in the earliest stages of the low-level vortex development (cf. Fig. 5a); the asterisk marks the vorticity centre diagnosed using the tracking algorithm. Note how the centre is situated on the southern flank of the strong easterlies of a subtropical anticyclone and associated with an easterly cold surge event. Within the shear zone between 4–7^∘ N, a well-defined wave has developed in both the v and balanced w components, as well as in Ertel PV. This structure is even more apparent in vertical cross-sections along 6^∘ N in Fig. 12a and b. At this time, the signature in v is strongest in the layer 4–9 km, upwards from the moist stability interface at about 5 km, but with some signature extending downwards too. The signature in balanced w_diab is stronger at higher levels than in v, reflecting the extension of the latent heating in deep convection into the upper troposphere. As characteristic of a Rossby wave, relative vorticity has to be in quadrature with v and so are the PV anomalies at 7 km (this analysis level is chosen between the height of the moist stability interface, at 5 km, and maximum w_diab, at 9 km). The positive w_diab maxima are shifted westwards relative to the positive vorticity maxima into the southerly sectors of the wave. This structure means that typically the ascent leads the positive PV anomalies in the propagation of the wave westwards. As explained in the next section, this structure is a key signature of diabatic Rossby wave propagation in the mid-troposphere.

Figure 12(a) Longitude–height cross-sections of meridional wind (shaded; m s⁻¹), and (b) Ertel PV (shaded; PVU) along 6^∘ N, from the global MetUM simulation initialised at 12:00 UTC on 21 October 2018, valid at 00:00 UTC on 22 October 2018 (T+12). (c) Meridional wind (shaded; m s⁻¹) and (d) Ertel PV (shaded; PVU) with horizontal wind vectors overlaid (m s⁻¹; reference vector = 10 m s⁻¹) at 7 km. In all panels, balanced vertical velocity forced by diabatic heating is overlain (blue contours; −2, 2, 6 and 10 cm s⁻¹).

5.2 Theory for diabatic Rossby wave propagation

Diagnosis of the MetUM simulations above has shown that a westward-propagating wave-like disturbance exists with a coherent signature in v, as well as Ertel PV and w in the mid-troposphere, and θ near the lower boundary. Here, following Bretherton (1966) and Hoskins et al. (1985), we will describe any such wave in PV as a “Rossby wave”. This formulation includes waves that exist on a basic state (here sector zonal average) meridional PV gradient or θ gradient at the lower boundary. The common feature is that the propagation mechanism relies on meridional advection of PV (or lower boundary θ) and the meridional flows that the chain of PV (or θ) anomalies induces. Hoskins et al. (1985) argue that the propagation direction of a wave and the nature of its interaction with a second wave at different level (mutual growth or decay) can be anticipated without explicit calculation of the flows induced by each wave using a specific form of balance or PV inversion relation. All that is required from balance is that v induced by one wave is in quadrature with PV. Then, the sense of propagation and interaction depends only on the phase difference between variables and also between a pair of waves. Heifetz et al. (2004) derived evolution equations for the amplitude and phase of such counter-propagating Rossby waves (CRWs) that can provide a mechanistic explanation for baroclinic instability. De Vries et al. (2010) extended the framework to moist dynamics by parameterising $\dot{\mathit{\theta }}$ in terms of w and then using the omega equation to relate induced w to PV in the waves. Following other authors (e.g. Moore and Montgomery, 2004; Boettcher and Wernli, 2013) the term “diabatic Rossby wave” was used to describe PV waves where the zonal propagation is dependent on the coupling with w and latent heat release. The argument for propagation is similar but depends on the phase of induced w rather than v (or a combination of both). Here we will use this CRW approach to analyse the output from the MetUM and understand whether the zonal phase speeds of the waves are consistent with the propagation mechanism and the role of latent heating in the waves.

Begin with the Ertel PV equation linearised about a basic state zonal flow, including the effects of $\dot{\mathit{\theta }}$:

where q is perturbation Ertel PV and ζ is the vertical component of basic state absolute vorticity. U(y,z), Q(y,z) and ρ_r(z) are the basic state zonal flow, PV and density, respectively. Note that the PV is materially conserved following the full flow in the absence of heating, including vertical advection, a major difference from quasi-geostrophic theory where only horizontal advection by the geostrophic flow is considered. Here we will take the approach that, although the balance approximation is not known precisely, Eq. (4) describes the evolution of the flow, where v and w are understood as the anomalous winds induced by the PV waves. The SGT tool has also been used to show that the signature of the balanced semi-geostrophic wind is indeed similar to that of the full wind in the MetUM and is the justification for the analysis. This result enables application of PV thinking and the CRW approach, where the propagation rates of the waves can be estimated and also attributed to the effects of $\dot{\mathit{\theta }}$ and the conditions for growth or decay through wave interaction can be deduced simply from phase differences between waves.

De Vries et al. (2010) investigate two types of closure in the parameterisation of heating: “large-scale rain” and wave-CISK (Conditional Instability of the Second Kind), which differ in the assumed relationship between $\dot{\mathit{\theta }}$ and w. Here, the large-scale rain approach is used, in which the heating rate is proportional to the product of basic state specific humidity, r, and w at each location (no meridional variation is assumed for simplicity):

This approximation is motivated by the ascent of saturated air masses where condensation results in latent heat release. Note that the typical balance of the thermodynamic equation in the tropics is between large-scale w and $\dot{\mathit{\theta }}$, in which case ϵr≈1. However, when there is also geostrophic forcing of w, this parameterisation is only valid for ϵr<1 at all levels. It is typical to take ϵ to be a constant. Linearisation of the dynamics requires that there is a symmetric cooling where there is descent; this requirement is partly justified by arguing that there is enhanced long-wave cooling from clear air regions. If it can be assumed that rN² varies more quickly in the vertical than w does at a level that we will call a “moist stability interface” (see Fig. 9b), so that w can be taken outside the heating gradient $\partial \left(\mathit{ϵ}r{N}^{\mathrm{2}}w\right)/\partial z$, then Eq. (4) becomes

where $G=-\mathit{ϵ}(\mathit{\zeta }/{\mathit{\rho }}_{\mathrm{r}})({\mathit{\theta }}_{\mathrm{0}}/g)\partial \left(r{N}^{\mathrm{2}}\right)/\partial z$ encapsulates the basic state moist stability gradient that the diabatic PV wave propagates along as a result of the gradient in heating coupled with w. Together with the equation for the evolution of θ on the lower boundary, Eq. (4) completely specifies the balanced evolution because the PV can be inverted to obtain the balanced wind field at any instant. Appendix A shows how the CRW theory can be used to deduce expressions for the zonal phase speed of a PV wave at upper levels (labelled wave-2 with induced velocity amplitude v₂, w₂) and a wave in θ along the lower boundary (wave-1 with velocity v₁, w₁):

The zonal wavelength is $\mathrm{2}\mathit{\pi }/k$, A₂ represents the upper wave amplitude in terms of Ertel PV and A₁ the lower wave amplitude in θ. The phase speed involves advection by the basic state zonal flow at the level where the wave activity is concentrated for each CRW. The upper CRW propagates relative to the zonal flow, U₂, through a mechanism associated with the effect of w and associated heating in the presence of the moist stability gradient (seen in Fig. 9b between 4 and 5 km) as well as meridional advection of the basic state meridional PV gradient, Q_y, and vertical advection of the vertical PV gradient, Q_z. In the case examined, G is positive and therefore contributes to propagation of the wave to the west, adding to the advection which is also towards the west (U₂<0). $c_{\mathrm{2}}^{\mathrm{1}}$ represents the modification of the phase speed of CRW-2 associated with interaction with v₁ and w₁ induced by CRW-1 at the level of CRW-2.

The lower CRW propagates against the zonal flow, U₁, through meridional advection of the basic state meridional θ gradient, Θ_y, and the interaction with the upper CRW, represented by $c_{\mathrm{1}}^{\mathrm{2}}$. The condition w=0 at the lower boundary is used here. Heifetz et al. (2004) have shown that in a growing normal-mode phase-locked configuration the definition of the CRW structures is such that $c_{\mathrm{2}}^{\mathrm{1}}=-c_{\mathrm{1}}^{\mathrm{2}}$, although in general initial value problems this relation will not hold. These interaction terms are expected to be smaller than the self-propagation rates.

The self-propagation rate of the diabatic Rossby wave along the moist stability interface (zone of sharp vertical gradient in rN²) can be calculated as

where r(z) and N(z) are the background specific humidity and static stability profiles, evaluated in an upper layer (U) and lower layer (L) separated by an interface zone with a characteristic depth, Δz, and θ separation, Δθ. This factor represents the strength of the basic state interface in terms of the contrast in moist stratification that the wave propagates along. The stronger the contrast, the faster the propagation. The basic state absolute vorticity is calculated from f−U_y. The final component is the amplitude of the wave in the non-dimensionalised diabatic heating rate, ${\dot{\mathit{\theta }}}_{\mathrm{2}}/\mathrm{\Delta }\mathit{\theta }$, which is assumed proportional to the vertical velocity, w₂, in balance with the diabatic Rossby wave motion. This amplitude is estimated from the heating rate represented on scales resolved by the global MetUM; the SGT tool shows that the balanced response to heating dominates the balanced w in the tropics, as expected from scale analysis of the thermodynamic equation. Note that for the lower boundary thermal wave, the amplitude A₁ is measured in terms of θ anomalies, and the model level at approximately 1.4 km is used to represent the anomalies and basic state meridional gradient (just above the boundary layer).

5.3 Estimate of zonal propagation rates and baroclinic interaction

Tracking the cyclonic vorticity centres at the upper and lower levels (Fig. 11), the propagation rate of the upper wave is 5–6^∘ d⁻¹ (6–8 m s⁻¹) westward. The lower wave and Borneo vortex propagate more slowly: at about 4–5^∘ d⁻¹ (5–6 m s⁻¹) from 22 to 23 October and 2^∘ d⁻¹ (2.6 m s⁻¹) from 23 to 24 October. The theory shows that the movement of the waves is a combination of advection by the basic state zonal flow, self-propagation and a change in propagation rate through interaction between the waves (Eq. 7). The basic state is estimated by averaging over a domain along the strip of wave activity: 5.5–7.0^∘ N and 104–115^∘ E, and also averaging over the 5 d time window between 12:00 UTC on 21 and 26 October 2018. This calculation gives basic state zonal flows $U_{\mathrm{2}}\approx -\mathrm{7.1}m$s⁻¹ and $U_{\mathrm{1}}\approx -\mathrm{0.9}m$s⁻¹, estimated at 4.5–5.0 km (the level of the moist stability interface in Fig. 9) and 1.4 km respectively. Changing the latitude extent by 1^∘ or the time average or length of longitude strip (e.g. 104–120^∘ E) gives a variation in these estimates of about ±0.5 m s⁻¹ for U₂ and ±0.3 m s⁻¹ for U₁.

The self-propagation rate of the lower θ wave relative to the basic state shear can be estimated approximately using Eq. (7). The chief difficulty is estimating the meridional velocity v₁ attributable to the lower wave. However, since the flow induced by the upper wave is expected to be much weaker than the flow induced by the lower wave near the lower boundary (Heifetz et al., 2004), the amplitude of the induced flow is estimated from the standard deviation of v over the strip (defined as above) at 1.4 km. Similarly, the wave amplitude A₁ is estimated from the standard deviation of the θ perturbation at the same level. The zonal wavelength is approximately 900 km, and the basic state meridional θ gradient in this region is a positive, but weak, $+\mathrm{0.7}\times {\mathrm{10}}^{-\mathrm{6}}$ K km⁻¹. This calculation yields a propagation rate of 1.5–1.6 m s⁻¹ westwards, adding to advection by U₁ to give a phase speed of 2.4–2.7 m s⁻¹ westwards, consistent with the wave movement observed in the later period.

The self-propagation rate of the upper diabatic Rossby wave is harder to estimate. Equation (8) can be used to estimate the magnitude of the propagation rate associated with heating in a sheared environment on the moist stability interface. The basic state moist stability term in Eq. (8) is calculated using 2–4 km to define the lower layer and 6–10 km for the upper layer specific humidity and static stability. The quantities $\stackrel{\mathrm{‾}}{r}$ and ${\stackrel{\mathrm{‾}}{N}}^{\mathrm{2}}$ are given by the average of upper and lower layer values. The moist stability interface depth is Δz=2 km and Δθ=27 K. The zonal wavelength is 900 km, as for the lower wave. The upper wave amplitude in Ertel PV, A₂, is estimated from the standard deviation of PV in the strip domain at 7 km. The wave amplitude in heating rate is estimated from the standard deviation of the diabatic temperature tendency (both convection and large-scale rain parameterisation in the MetUM) also at 7 km. Bringing all these terms together the self-propagation rate is 1.3–1.7 m s⁻¹ westwards. Averaged over the same strip at 7 km, the meridional PV gradient is negative (being in the region between the strong easterlies of the cold surge and the cyclonic vorticity strip on its southern edge) with estimated $Q_y\approx -\mathrm{0.57}\times {\mathrm{10}}^{-\mathrm{13}}$ in SI units and v₂≈3.4 m s⁻¹, yielding an eastwards contribution to propagation of 0.9 m s⁻¹. The vertical PV gradient term involving Q_z is not very coherent and has not been estimated. In contrast, G and Θ_y are positive over the whole domain of interest, meaning that propagation rate associated with both of these terms is negative definite (westwards), while the Q_y term is weak and eastwards. Taking all terms together yields a total phase speed of 7.2–8.4 m s⁻¹ westwards, slightly faster than the observed movement. Note that the cyclonic vorticity strip and upper wave exist 1 d or so before the lower wave and Borneo vortex, so it is possible that the upper wave arises first through barotropic instability of this vorticity strip.

In summary, the self-propagation rates of both waves are considerably weaker than the shear in the background zonal flow, and the propagation through interaction ($c_{\mathrm{2}}^{\mathrm{1}}$) must be smaller still. So even if the configuration of the waves were favourable, they would not be able to stay phase-locked and must shear apart. However, over the early stages of development the observed difference in wave phase speeds is only 1–3 m s⁻¹, and therefore the estimates of propagation rates from theory are comparable. These self-propagation speeds (1.3–1.7 m s⁻¹) provide an upper bound on the strength of baroclinic interaction between the waves (Heifetz et al., 2004), and $\mathrm{1}/\left(ck\right)$ gives a growth timescale of 1–1.5 d. Therefore, although the strength of baroclinic interaction is relatively weak compared with the shear rate, baroclinic interaction is a plausible explanation for the emergence of the lower wave and Borneo vortex.