A new convective parameterization applied to Jupiter: implications for water abundance near the $24\degree$ N region

The RAS uses the same fundamental framework as all Arakawa-Schubert (AS) schemes. The general principle is to calculate all possible clouds within the vertical grid column, and determine the updraft profile, and subsequent changes to the environmental variables as a result of the updraft. All sub-grid scale clouds are assumed to have a base at the same level $z_{\rm B}$, but have varying cloud top altitudes $z_{\rm D}$ (the detrainment level), which, in the model, can start anywhere above the first layer up from $z_{\rm B}$. Therefore, several clouds (or ‘cloud types’) can co-exist within the same grid column in our model, and the $i$th cloud is defined as one where the $z_{\rm D}$ for cloud $i$ lies on the $i$th vertical model layer.

Convection is driven by the idea of buoyancy: pockets of air with lower density than the surrounding will freely rise, whereas those that are denser will sink. Therefore, convective cloud updrafts are driven by having positive buoyancy, and thus, a cloud that has a top at $z_{\rm D}$ is positively buoyant between $z_{\rm B}$ and $z_{\rm D}$ and has neutral buoyancy at $z_{\rm D}$ (the cloud must stop convecting freely when the positive buoyant forcing disappears). For each cloud invoked, the RAS scheme determines whether this cloud will be convective based on this buoyancy condition, and determines the updraft profile such that the buoyancy profile within the updraft matches this constraint. A second criteria for convection involves the conservation of energy within the grid cell. Convective upwellings must use the energy that is available for convection within the grid cell and thus, the ensemble of clouds within the grid cell must be limited by the total convective potential energy of all possible clouds within the column.

Following MS92, we calculate the buoyancy of a parcel within the updraft as,

where $g$ is the local acceleration due to gravity, $T_{\rm v}$ is the virtual temperature, $q_{\rm C}$ is the cloud specific humidity (ratio of cloud mass density to total density). Overbars denote atmospheric (grid scale) values, while variables without the overbars are in-cloud (sub-grid scale) values. The first term corresponds to the traditional buoyancy term (i.e. difference in density), while the second terms adds the contribution (i.e. reduction in buoyancy) from the weight of the in-cloud condensates.

Within the cloud, upwelling carries energy and mass to the upper levels, and thus the in-cloud temperature and cloud mixing ratios are determined from the vertical updraft profile. Generally, the vertical mass flux profile can be a complex function of the environmental variables but for simplicity, we follow AS74 and assume that the updraft entrains a constant fraction of mass flux $\lambda$, within a layer $\delta z$,

We can normalize the mass flux by the value at the cloud base, such that $M(z)=M_{\rm B}\eta(z)$, where $M_{\rm B}$ is the mass flux at the cloud base and $\eta(z)$ is the normalized mass flux profile. Therefore, the mass flux relation above simplifies to

In this fashion, $\lambda$ and $M_{\rm B}$ completely define a given updraft, where $\lambda$ defines the fraction of vertical mass flux that is entrained from the surrounding environment at any given level, and $M_{\rm B}$ defining the strength of the updraft. Therefore, we can solve for the normalized profile by integrating the differential equation. We integrate from the cloud base ($z_{\rm B}$, where the updraft begins), up to a level $z$, that is above the base,

where $z_{\rm B}$ is the altitude of the cloud base. This ensures that $\eta(\lambda,z_{\rm B})=1$, such that $M(\lambda,z_{\rm B})=M_{\rm B}$. However, in this form, solving for $\lambda$ is difficult. The RAS scheme instead approximates $\eta(\lambda,z)$ as,

where $\zeta=z-z_{\rm B}$ and $\xi=\frac{1}{2}\zeta^{2}$. This quadratic approximation reduces the computational cost of determining the updraft mass flux profile, while retaining the necessary accuracy.

With the given updraft parameterization, we solve for the vertical profiles of mass and energy with the vertical continuity equation. We define the dry ($s$), moist ($h$) and moist saturated ($h^{*}$) static energies, as,

where $c_{\rm p}$ is the specific heat capacity at constant pressure, $T$ is the gas temperature, $\Phi$ is the geopotential, $L$ is the latent heat per unit mass, and $q$ is the vapor specific humidity of the moisture species in the air. $h^{*}$ is the saturated moist static energy, which is similar to $h$, but defined from $q^{*}$, which is the saturation vapor specific humidity at a given pressure and temperature, i.e., the specific humidity at which the atmosphere is fully saturated. During adiabatic motion, $s$ is a conserved quantity, while during moist ascent, $h$ is conserved.

With these definitions, we can write the cloud and atmospheric virtual temperatures $T_{\rm v}$ and $\overline{T}_{\rm v}$ in terms of the moist static energy. Therefore, the first condition (neutral buoyancy at the cloud top), then becomes, using the definition of the buoyancy from Equation 1,

where $\overline{h}^{**}$ is the grid-scale virtual moist saturated static energy, and $\tilde{L}$ is latent heat term in Equation 7, accounting for virtual effects. Effectively, in the limit of low cloud condensate density, the cloud top buoyancy condition is matched when the updraft moist static energy profile ($h$) matches the grid-scale virtual moist saturated static energy profile ($\overline{h}^{**}$). Figure 2 shows a schematic of how this cloud top buoyancy condition is matched at different heights. Since $\overline{h}$ is usually smaller than $\overline{h}^{**}$ for a subsaturated atmosphere, increasing the entrainment from the surrounding generally decreases $h$ within the updraft, shifting the curve to the left. The value of $\lambda$ is the slope of $h$ that intersects $\overline{h}^{**}$ at $z_{\rm D}$.

With this parameterization of $\eta$, we can determine the cloud top buoyancy condition ($B(z_{\rm D})=0$), which gives,

where $a$, $b$ and $c$ are functions of only the grid-scale variables (i.e. model values). This gives two solutions for the entrainment parameter $\lambda$, where we choose the highest positive value that is less than $10^{-4}$ m^-1 (i.e. an entrainment of 10% mass flux per kilometer of updraft). This is an arbitrary value currently chosen to maintain numerical stability in our model. This upper limit of $\lambda$ essentially prohibits very shallow updrafts, as shown by the lower grey region in Figure 2. This parameter can be tweaked to optimize the false positive rate of updraft triggers for shallow convection. Generally, for an atmospheric profile similar to the one shown in Figure 2, where $\overline{h}$ decreases with height until near the tropopause, it is impossible to get two positive solutions, since $h$ will generally decrease monotonically with height, making both $a$ and $b$ negative. The derivation of $a$, $b$ and $c$ are given in Moorthi & Suarez (1999), and in the Appendix, along with the vertical tendencies for the moisture and thermodynamic variables.

Prior to calculating the updraft tendencies, the cloud base mass flux, $M_{\rm B}$ must be determined. We need to determine $M_{\rm B}$ which describes the strength of the convective upwelling, through a closure (given by energy conservation) for our system of equations. Convection is triggered by the increase in convective available potential energy (CAPE) or the decrease in convective inhibition (CIN), through external forcing (i.e., large-scale dynamics). CAPE defines the total potential energy gained by a parcel when it rises through a region of positive buoyancy,

while CIN defines the total energy required by the parcel to rise through a region of negative buoyancy to get the base of the upwelling layer,

Therefore, the trigger function for convection within a column is tied to change in buoyancy for clouds of different heights. We implement this using the idea of large-scale changes in atmospheric CAPE, similar to Zhang (2002). The change in CAPE for cloud type $i$ due to cumulus effects is obtained by integrating the change in buoyancy with height,

where $F$ defines the fractional change in CAPE for cloud type $i$ due to a unit mass flux, $F$, and the cloud top height ($z_{\rm top}$) in Equation 11 is replaced with the detrainment level for cloud type $i$. We determine the large scale effects on the CAPE for cloud type $i$ using finite difference,

Therefore, we can then determine the cloud base mass flux, as,

The derivation of the CAPE formalism used in this scheme, and the equation for $F$ are given in the Appendix.

The distinction between the AS and RAS scheme is that the RAS relaxes the sounding towards a quasi-equilibrium over several timesteps, rather than within a single one. This is due to the fact that solving directly for equilibrium is generally too strict for convective parameterization (Yano & Plant, 2020). Consequently, we only allow a fraction $\alpha$ of the mass flux to affect the sounding. This acts similarly to, for example, under-relaxed Jacobi’s iteration method for linear system solvers, where the true solution to the system of equation is obtained after several iterations.

In the case of moist convection, this means that the quasi-equilibrium will be achieved after $n$ timesteps, or over a timescale of $\tau=n\Delta t$, where each timestep steps the model forward by $\Delta t$. $\tau$ is then the timescale over which inter-cloud interactions self-stabilize the large-scale forcing. Therefore, with this principle, we can define the relaxation parameter, $\alpha$, as,

where $\tau\sim 10^{3}-10^{4}$ s on Earth (AS74). On Jupiter, this is currently unknown from observations. However, modeling attempts show that these convective storms stabilize on the timescale of ${\sim}10$ hours (Hueso & Sánchez-Lavega, 2001). Therefore, we test different values of $\tau$ to see the sensitivity of the parameter to jovian convection. On Earth, increasing $\tau$ linearly increased the amount of iterations required to relax the model to the critical value of CAPE (MS92), so we expect the process to work similarly on Jupiter.

The model is discretized vertically as shown in Figure 3. EPIC features a vertically staggered grid, with thermodynamic and dynamic variables stored on the grid interface (i.e. top/bottom of the cell) rather than the center, and so that index $k$ denotes the $k$th edge from the top (Figure 3). $K$ defines the cell edge corresponding to the bottom of the cloud, i.e., $z_{K}=z_{\rm B}$, and $i$ denotes the grid interface corresponding to the cloud top. In the original RAS scheme, values are stored in the cell center, and the cloud starts at the cell center $K$ and detrains at the center $i$. However, in EPIC, we retain the numbering scheme but note that integer values correspond to edges, rather than centers. Correspondingly, the mass fluxes are defined on the cell centers in EPIC, (i.e., at half-integer indices), while the tendencies (mass and heat) are located on the cell edges (integer indices) so as to be able to update the thermodynamic variables. $\eta_{K-\frac{1}{2}}=1$ since that is the base of the updraft and $\eta_{i-\frac{1}{2}}=0$ because the updraft stops at $i$. Thermodynamic values at half-integer indices are interpolated similarly to Konor & Arakawa (1997). The vertical discretization is shown in Figure 3.

Fig 9 shows the cloud phenomonology for the different abundance cases 100 and 150 days into the simulation. The water cloud densities are shown in blue and the ammonia clouds are shown in orange. Cases II.1-4 show large scale cloud formation. The clouds form predominantly near the equator, south of the 24$\degree$ N jet and near $30\degree$ S, and they evolve between the two times. For all cases, the cloud formations are similar in terms of size and location. The thickest, and tallest water clouds reach the base of the ammonia layer, and are driven by moist convective motion rather than large scale dynamics. The water storms are also accompanied by upper level ammonia clouds. The largest clouds are in the equatorial zone and in the belt just south of the $20\degree$ N latitude. There are some storms that form at around $20\degree$ S.

Cases II.5-7 show little cloud activity throughout the model. Both cases show short clouds dispersed throughout the region which persist at both times without much difference. At day 150, Case II.7 shows a storm forming near the equator. However, this storm is still driven by large scale motion, as indicated by the wind vectors in its vicinity. Therefore, even here, this is moreso a dynamical response to the instability in the atmosphere (i.e dry convection), rather than one driven by moist convection (which would not require large scale vertical motion).

Fig 10 shows the distribution of the storm frequency as a function of latitude for the different abundances. For cases II.2-4, the storm formation is generally distributed throughout the domain. However, there are also specific regions where storm formation is noticeably higher. These include the region around $30\degree$ S, between $5\degree-20\degree$ N latitude and near the peak of the $24\degree$ N jet. These regions are more prominent in Cases II.3 and II.4 compared to Case II.2, with the regions around $30\degree$ S and near $24\degree$ N showing the most convective activity. For cases II.1, 5, 6, 7, there are no noticeable regions that display strong convective potential. As expected, the moist convective events follow regions of dynamical instabilities, as a lot of convective peaks occur near the jet peaks (which are on the boundary between cyclonic and anticyclonic shear, by definition), since sub grid scale moist convection is forced by large scale instabilities. Fig 11 shows the mean distribution of storms for different water abundances. This clearly shows the drop-off in convective activity above $5-6\times$ solar water abundance is evident here, as a result of the ‘mass-loading’ effect of the heavy condensibles.

The distribution of moist convective activity for cases II.2-4 compare favorably with the detection of lightning using the Juno Microwave Radiometer (MWR) data (Brown et al., 2018), which shows peaks in lightning detection around $20\degree$ S, and between $10\degree-25\degree$ N. Our results, particularly match the equatorial peak at ${\sim}8\degree$ N latitude. We also observe less convective activity near the equatorial region for all cases. However, there are differences between the observed convective profile our model results, which we attribute to the following:

Brown et al. (2018) suggest that the discrepancy in belt/zone moist convection might be due to variation in the water abundance in those regions. Indeed, Li et al. (2017) shows that this is true for the distribution of ammonia, which is highly variable, even below $1$ bar, at different latitudes. In our model, we assume that the deep water abundance is constant at different latitudes, which would explain the nearly constant meridional distribution.

Furthermore, the equatorial region contains a significantly higher concentration of ammonia vapor (Li et al., 2017), compared to the rest of the planet. In this region, the increased ammonia concentration would lead to higher atmospheric molar mass, compared to other regions. This would consequently lead to a lower value for CAPE, since the parcel density would be closer to the atmospheric density here, compared to one that is drier.

We assume a constant forcing as a function of latitude, which is likely true, if we account only for radiation. However, moist convective activity is also driven meridional circulation cells, vorticity, etc., which cannot be simulated effectively in this 2-dimensional setting. Brown et al. (2018) also suggest the increase in moist convective activity above $50\degree$ latitude is suggestive of preferentially poleward distribution of the internal heat flux.

The large scale dynamics also act as a response to the increasing upper level instability. Therefore, the response to the forcing is split partially into the response from water convection and, partly from the changes in the atmospheric structure from large scale dynamics. The increased convective activity at $30\degree$ S in our model is particularly striking, compared to what was determined by the MWR data. The region around $30\degree$ S shows moist convective activity across all the cases, implying that this is indeed a feature of the model that seemingly does not agree with the natural conditions on Jupiter. Fig 12 shows the region from Juno’s JunoCam, with the latitude values annotated on the right side. The region between $30\degree$ and $40\degree$ S are filled with turbulent cyclonic features, suggesting that dynamical instability are a common feature in this latitudinal band. However, the lack of lightning detections suggests that moist convective triggers are minimal. Therefore, it is likely that the instabilities manifest quickly as cyclonic vortices, thereby reducing, or completely negating the need for a moist convective outburst. This is indeed similar to the findings of Showman (2007), who found that adding energy to the atmosphere at higher latitudes resulted in the emergence of vortices, rather than jets. He found that the critical latitude for the transition between jets and vortices was around $45\degree$, which is comparable to the latitude of this region.

Alternatively, storm formation in the region near the $30\degree$ S latitude might not be strong enough to generate lightning, or not have the proper atmospheric structure to generate detectible lightning storms. This region might promote shallow convection, due to the response from the dynamics reducing the instability in the upper layer. It should be noted, however, that there are likely detections of lightning below a pressure level of 5 bars (Vasavada & Showman, 2005), where the temperatures are too warm for ice to form. The liquid clouds are fairly inefficient in generating sufficient charge separation for lightning to form (Gibbard et al., 1995). Therefore, these observations are currently unexplained.

Our model grid for the 3D cases is the same as the one in S21 and covers the region from 15$\degree$-$33.5\degree$ N with 80 grid points in the meridional direction. Zonally, the model covers $0\degree$-$120\degree$ longitude, with 256 grid points. This produces a resolution of $0.47\degree$ per grid point zonally and $0.23\degree$ per grid point meridionally. The angular resolution is higher in the meridional direction to compensate for the geometric effect at the mid-latitudes, ensuring that the physical extents of the grid points are roughly equal (about $200$ km per side, in this case). The zonal wind profile used is from Voyager data from Limaye (1986) and is shown in Figure 13 for the model region. We use the same vertical wind shear assumption as in our 2D test cases, and our vertical temperature structure is the same as the 2D cases described above, and is shown in Figure 8.

One of the main goals of this study is to interpret the effects of different water abundances on the resulting atmospheric and cloud structure, and on the convective capability of the atmosphere. To compare the new convective parameterization with previous EPIC simulations of this region without the RAS scheme, we use the same initialization as in S21, and test $0.5\times$, $1\times$ and $2\times$ solar for water. Given the lack of difference between the ammonia abundance on both the convective potential and the dynamics of the wave, we maintain the deep ammonia abundance at $2\times$ the solar [N/H] ratio, which is consistent with Li et al. (2017).

Similarly to S21, we also allow the model to evolve for 200 days without any cloud processes enabled, in order to let the model adjust to the initial parameters. This is to ensure that there are no unphysical cloud formation in the first few timsteps, which affect model results. We also use the same perturbation profile to remove the zonal symmetry set up in our model as a result of this ‘spin up’ phase, whereby we induce 100 small vortices (up to 10 m/s in speed) randomly throughout the model between a pressure of 500-5000 hPa. We use three different set of these random vortices to ensure that our model results is not dependent on the profile of the perturbation. We test three different profiles for each deep abundance value of water. The list of test cases is shown in Table 2.

For the RAS scheme, we tested simulations with $\tau=30$ min, $\tau=1$ hour and $\tau=2$ hours for all cases. We noticed that the simulations with $\tau=30$ mins produced excessive horizontal and vertical cloud coverage that was inconsistent with observations. We did not observe this effect in the simulations with $\tau=1$ and $\tau=2$ hours. Therefore, for the following analysis, we will only use $\tau=1$ hour for the $0.5\times$ and $1\times$ cases and $\tau=2$ hours for the $2\times$ cases. The different relaxation times lead to the same result towards the end of the simulation, but the higher relaxation times for the $2\times$ case allows for better stability in the beginning.

With the moist convective scheme active, the atmosphere quickly evolves into a similar wave-driven structure seen in the stratiform case (S21). The dynamics of the motion of the wave results in up- and downwelling that leads to distinct cloud features. There are arc shaped ammonia clouds and deep water clouds that are extremely localized due to the same circulation described in S21. In the sections below, we look at the differences between the stratiform cases and the addition of the new convective module, to investigate the effect of moist convection in this region. We also draw new conclusions on the convective potential of this region, and implications for the water abundance.

Figure 14 shows the change in the zonal wavenumber of the wave with time. The trend for the different cases is downward with time (i.e. the wave stretches zonally with time), similar to the stratiform case detailed in S21. The $0.5\times$ and $1\times$ stabilize at a wavenumber of $\sim 10$, while the $2\times$ case drops to a wavenumber of $\sim 8$ after about 40 days. The stratiform case reaches this state a few days later.

Overall, the moist convective upwellings do not affect the wave, and the dissipation (i.e., lengthening of the wave over time) happens at roughly the same pace as before. Understandably, this is likely due to the moist convective plumes reaching a peak well below the location of the wave, and thus is a negligible change to the thermodynamic structure of the upper atmosphere.

Figure 15 shows the difference between the stratiform and convective cloud morphology for the $2\times$ solar cases (Cases III.3, 6 and 9) at day 45. The blue and red contour show the location of the wave as determined by calculating the eddy potential vorticity (i.e., the departure of the Ertel potential vorticity from the zonal mean value). The red contour corresponds to regions that have a cyclonic (positive) PV anomaly and the blue corresponds to anticyclonic (negative) PV anomaly. There is very little difference in the large scale structure of the ammonia clouds within the jet. The arc shape is still present, and the dynamics of the ammonia cloud match that of the stratiform cases. Understandably, this is due to the fact that ammonia does not form convective clouds in the model, and the effect of moist convection on ammonia would be through the upwelling of vapor from the depth, with the water-induced storms.

Indeed, the moist convective cases show an abundance of small scale outbursts, highlighted by the black arrows in Figure 15. These outbursts result in thick ammonia clouds as the moist convection brings ammonia rich vapor to the upper atmosphere, and are primarily concentrated to the south of the jet. They do disrupt the structure of the arc shapes within the wave, but are very short-lived, lasting on the order of less than a day. The number of such small scale moist convective updrafts increases with water abundances, as shown in Figure 16. This is expected, as the driver for convection is the water.

Figure 17 shows the water cloud density for perturbation 1 at day 30 for different water abundances, along with the convective mass flux in black. The mass flux, and the cloud density increase with the abundance of water, as expected. Furthermore, the storm formation is generally constrained to specific zonal bands, and is the highest near the jet peak. The storms within the jet lead the cyclonic packet of the wave (i.e. are on the eastern side of the cyclonic packets). The location of the wave, and the effect on the cloud density is consistent between the stratiform and convective cases, and also between the different water abundances.

Figure 18 shows the zonal slice across the jet peak for perturbation 1, for different water abundances, 30 days into the convective simulations. We see the same cloud morphology as the stratiform cases, with the thicker clouds existing below the cyclonic part of the wave, beginning to the east of the upwelling, while the anticyclonic part of the wave is still generally clear, and has downwelling. The cloud density corresponds directly to the amount of water, with the $2\times$ case having the thickest clouds.

The upwelling drives water clouds up to the 1-bar level, in the $2\times$ cases, which changes the habit for water ice, both in terms of the air density and ambient temperature, compared to the water cloud base. This changes the microphysical profile of the water ice particles within these updrafted water clouds, which is the key difference between simulations with, and without, the RAS scheme.

Figure 19 shows the cloud density for water and ammonia clouds in the top panel and cloud particle size in the bottom panel 30 days into the simulation with the first perturbation profile. Particles falling with terminal velocity greater than 50 cm/s are shown in black. As expected, particle sizes are largest within the upwelling, and thus have the largest terminal velocities. The red contours show the location of snow particles, which occurs where particle sizes exceed $500~{}\mu$m. Snow forms at the top of the upwelling and fall with the downdraft. Therefore, both the precipitation and the downwelling deplete the updrafted cloud particles on the eastern edge of the wave, which returns the water mass back to the deeper levels.

Within the updraft, the particle sizes are generally between 300-500 $\mu$m for the ice clouds. The fastest falling particles correspond to those just at the critical diameter of $500~{}\mu$m and fall with a terminal velocity of 70 cm/s, at a pressure level of 800-900 hPa. As these particles fall, their terminal velocity drops to about 58 cm/s near the cloud base at 4 bars. At an average velocity of 65 cm/s, this journey takes about 16 hours from the top of the storm (900 hPa) to the cloud base – a total height of about 37 km.

While snow particles are larger, their larger surface area decreases the terminal velocity compared to the bullet-shaped cloud ice particles. Therefore, while snow particles reach diameters up to 1 mm, near the cloud base, the terminal velocity of snow within a grid cell reaches at most 85-90 cm/s, which is not significantly higher than that of cloud ice.

Horizontally, the distance from the eastern edge of the cyclonic packet to the western edge of the next one is about $10\degree$ or about 11000 km at this latitude. At the 1-4 bar level, the horizontal wind speed is about 160 m/s, and thus advecting the water cloud particles from a region of downwelling to the next upwelling takes about 19 hours, which is a little longer than the vertical fall speed of the particles. Therefore, even with the moist convection lofting particles to the upper atmosphere, precipitation still is a viable method of removing volatiles from the convective plume, resulting in a rarefied atmosphere below the anticylonic packet.

Therefore, as water clouds travel below the wave, they are lofted up to the east of the upwelling below the cyclonic packet, until they reach the eastern edge, where convection stops, and they sediment to the lower levels. Sedimentation occurs over the anticylonic packet, until they reach the next cyclonic packet, where they undergo moist convective upwelling again, and the process repeats. This cycle is shown in Figure 20, with the red arrow denoting the packet that demonstrates this effect. As this parcel reaches the eastern edge, the moist convective intensity (shown with the black contours) and the amount of snowfall both decrease, and is restored as it crosses the western edge of the next cyclonic packet. Note how the water cloud dissipates near day 22 as it is pushed downwards, and thus decreases the relative humidity, and consequently, the cloud density.

For the $1\times$ solar case, the microphysical structure of the atmosphere is much more simple (Figure 21). The water clouds are shorter, reaching below the 1-2 bar level, and do not form the persistent convective towers in the eastern edge of the cyclonic wave packet. Indeed, both the snow formation and precipitation occur near the cloud base, rather than at 1 bar. Therefore, the water clouds in the $1\times$ cases are much more transient, compared to the $2\times$ case, even with the convective module. This is due to the relative lack of convection in the $1\times$ cases, compared to the $2\times$ cases.

Figure 22 shows the corresponding mass flux profile for the simulations in Figure 18. Here, we can see the increase in convective mass flux directly result in increased cloud density. Only the $2\times$ cases show deep convective storms, while the $1\times$ shows shallow convection. The $0.5\times$ show very little convection.

For the $2\times$ case, the storms are tied to the location of the convective packet, and more specifically, to the east of the updraft. The convective fluxes stop just to the west of the downwelling. Interestingly, the upwelling itself does not directly result in convection. Indeed, the lack of convection at the western edge of the cyclonic packet actually shows that the upwelling inhibits sub-grid scale convection. Instead, the updrafts adds moisture to the east of its location, which increases the CAPE below the bulk of the cyclonic wave packet, which in turn triggers convection.

Figure 23 shows the average water cloud densities for all nine simulations, and the corresponding location of the wavepackets. The black contours show the locations where storm frequencies exceed 1 every two days. As expected, the $0.5\times$ cases show very little storm frequency and convective mass flux, but the $1\times$ and $2\times$ cases show very strongly localized regions of convection. For the $1\times$ cases, the structure of the wave correlates well with the location of convective events, and the wave itself is maintained well. In the $2\times$ cases, the storms are more disorganized, and bleed into the belt to the south. The wave itself loses its structure from the upwelling (especially in the case of perturbation 2).

Figure 24 shows the zonal and temporal mean of the storm frequency and intensity as a function of latitude and abundance. The shaded region denotes the jet, and the vertical black line corresponds to the jet peak. As stated before, the jet itself is the region of highest convective activity across all the cases, and interestingly, does not show uniform distribution of convective storms. Most of the storm formation is to the south, with a decrease in storm formation just to the north of the peak, ending with more storms on the northern edge. Elsewhere, the belt on the southern boundary shows some convective activity, as does the belt to the north. Convection seems to be primarily concentrated near the peaks of the jets at 17.5$\degree$ N and $23.7\degree$ N. To the north of $25\degree$ latitude, the convection is distributed throughout, without any distinguishable structure.

Interestingly, the $2\times$ case show a pretty uniform distribution of storm intensity across the region, even though there is a distribution in storm frequency. Therefore, for the $2\times$ case, the moist convective flux is nearly uniformly everywhere, with differences in the atmospheric structure between the regions resulting in either recurring, low intensity storms, or infrequent, high intensity storms. Near regions of dynamical instability (e.g., at jet peaks), the former is true, while elsewhere, it is the latter. The uniformity is likely a function of the initial atmospheric structure of this region, in terms of both the distribution of water and the temperature at depth.

The distribution of CAPE, as shown in Figure 25 reveals a similar profile. Here we show the CAPE determined from the cloud base, following the convective updraft profile up to the top of that cloud type, using Equation 11. The $0.5\times$ solar cases show very little CAPE, except at the jet peak. In the $1\times$ and $2\times$ solar cases, the CAPE is concentrated near the jet peak and to the south, which corresponds closely with the distribution of storm formation in Figure 24.

The value of CAPE in our model average between 500-1000 J/kg for parcels reaching the 1 bar level, with the higher values corresponding to the larger water abundances. For shallow convection, CAPE was only about 30 J/kg for the $1\times$ cases and about $100$ J/kg for the $2\times$ cases. While this is inconsistent with our stratiform cases, note that CAPE defined here follows the parcels ascent profile such that the convecting cloud is at most neutrally buoyant at the cloud top. Therefore, CAPE defined here also accounts for the decrease in buoyancy from the entrainment of dry material, making it much more consistent with an actual updrafting parcel, and also much smaller than the CAPE determined for the stratiform cases.

Since CAPE corresponds to the total potential energy gain from convection, we can equate this to the kinetic energy of the parcel at the top of the updraft, and determine an approximate mean vertical velocities within these convective plumes, as,

which means that for shallow convection, the mean updraft vertical velocity is under 10 m/s, travelling only about 5-10km. These updrafts thus take less than half an hour to reach the cloud top. With a mean mass flux of about $10^{-3}$ kg/m²/s (Figure 24), which is approximately equal to the base mass flux, since these clouds only travel 1-2 layers vertically, this means that these shallow updrafts only move a few kilograms of air for every square meter that they cover.

For the deeper convective cells, reaching the 1 bar level, and a CAPE of about 1000 J/kg, we get a peak updraft vertical velocity of about 45 m/s, and mean velocities of around 20 m/s, which is consistent with moist convective modeling jovian storms (Hueso & Sánchez-Lavega, 2001). Given a vertical height of about 30 km from the cloud base to the top, the updrafts should take under 20 mins, making these storms very efficient at quickly carrying energy upwards. While the cloud base mass flux for these systems are comparable to those for shallow convection, the increase in mass flux, as the cloud entrains, results in a net flux that is an order of magnitude, or more, larger. Note that there are multiple storm cells that reach the layers between 800-2000 hPa, each with a mean cloud base flux of $\sim 10^{-3}$ kg/m²/s. Therefore, these deep convective cells transport a correspondingly larger amount of mass and energy vertically, as expected.

The net energy flux is more complicated to determine. The convection happens as a response to large scale forcing, and so the cloud base flux is created to reduce the increase in CAPE. In our simulations, the mean rate of change of CAPE in the deeper levels is only about $10^{-3}$ W/kg, and an order of magnitude larger for the deep convection for the $2\times$ cases. Curiously, the $1\times$ cases show increased CAPE in the 800-1000 hPa region, but very little convection occurs in this region. Indeed, while the CAPE is large, the rate of change of CAPE is nearly negligible, meaning that there is not sufficient forcing in the atmosphere to initiate convection up to the $1000$ hPa level in the $1\times$ cases.

Intuitively, we can determine the net energy flux observed as,

where $P_{c}$ is the convective energy flux (i.e., rate of change of CAPE), $\rho$ is the atmospheric density and $F$ is the emitted power. Note that this is a trivial, order-of-magnitude estimation, rather than a rigorous calculation.

Figure 26 shows the net energy flux from convective activity as a function of latitude, for different water abundances. Here, we see that the general trends as described above, with the distribution of convection both vertically and meridionally being similar to what was described above. The key difference, however, is in the stark difference in the convective energy flux between the water abundances: for the $0.5\times$ cases, there is very little convective flux, which has a mean value $0.1$ W/m², the $1\times$ has a mean value of $6.3$ W/m², while the $2\times$ is around $13$ W/m². These mean values are represented well in the northern half of the model domain for both the $1\times$ and $2\times$ cases, but near the jet, there is an extremely large increase in the convective energy flux for the $2\times$ case, with a peak of around $80$ W/m².