## Abstract

Organized rainstorms and their associated overturning circulations can self-emerge over an ocean surface with uniform temperature in cloud-resolving simulations. This phenomenon is referred to as convective self-aggregation. Convective self-aggregation is argued to be an important building block for tropical weather systems and may help regulate tropical atmospheric humidity and thereby tropical climate stability. Here the author presents a boundary layer theory for the horizontal scale *λ* of 2D (*x*, *z*) convective self-aggregation by considering both the momentum and energy constraints for steady circulations. This theory suggests that *λ* scales with the product of the boundary layer height *h* and the square root of the amplitude of density variation between aggregated moist and dry regions in the boundary layer, and that this density variation mainly arises from the moisture variation due to the virtual effect of water vapor. This theory predicts the following: 1) the order of magnitude of *λ* is ~2000 km, 2) the aspect ratio of the boundary layer *λ*/*h* increases with surface warming, and 3) *λ* decreases when the virtual effect of water vapor is disabled. These predictions are confirmed using a suite of cloud-resolving simulations spanning a wide range of climates.

## 1. Introduction

Large-scale tropical circulations, including the Hadley circulation and the Walker circulation, arise owing to either meridionally or zonally asymmetric external forcings (Schneider 2006; Gill 1980). However, steady overturning circulations can develop spontaneously with uniform sea surface temperatures (SSTs) and boundary conditions in cloud-resolving models (CRMs) without rotation (Bretherton et al. 2005; Muller and Held 2012; Muller and Bony 2015; Holloway and Woolnough 2016). This phenomenon is referred to as convective self-aggregation. The atmosphere is anomalously humid in the ascending patch, where deep convection is ubiquitous. In contrast, the atmosphere is anomalously dry in the descending patch, where deep convection rarely occurs.

Understanding convective self-aggregation may have broad implications on understanding the tropical weather and climate (Bretherton et al. 2005; Bretherton and Khairoutdinov 2015; Arnold et al. 2015; Emanuel et al. 2014). A number of studies have proposed that tropical cyclones and the Madden–Julian oscillation (MJO) can be considered as forms of convective aggregation with the Coriolis effect and the beta effect, respectively. This is partly because both tropical cyclones (Boos et al. 2016; Held and Zhao 2008) and the MJO (Pritchard and Yang 2016; Arnold and Randall 2015; Yang and Ingersoll 2013, 2014) can be simulated under uniform boundary conditions. Convective self-aggregation may also help shape the mean climate state. Emanuel et al. (2014) suggested warmer climates favor aggregated convection, which significantly dries the environmental atmosphere and in turn increases the overall outgoing longwave radiation (OLR). By means of this negative climate feedback loop, convective self-aggregation may help stabilize the tropical climate (Hohenegger and Stevens 2016).

The research focus has been on understanding what leads to convective self-aggregation by developing simple instability theories (Craig and Mack 2013; Bretherton et al. 2005; Emanuel et al. 2014; Beucler and Cronin 2016) and by diagnosing these instability processes in CRMs (Wing and Emanuel 2014; Muller and Bony 2015). Although details are still under debate, it was suggested by independent modeling studies that water vapor–radiation feedback and moisture–precipitation feedback might be important for simulating convective self-aggregation.

However, it is unclear what sets the basic features of the equilibrium state of convective self-aggregation, for example, its horizontal scale. Wing and Cronin (2016) proposed that the horizontal scale of convective self-aggregation is the distance over which an air parcel from the dry patch is remoistened by the surface moisture flux. This distance depends on the boundary layer height and the moisture exchange coefficient at the surface. However, by varying model parameters, the authors found that the remoistening length does not agree with the size of convective self-aggregation. This might result from neglecting other processes that affect the specific humidity within the boundary layer, for example, cold pools and shallow circulations. These processes are dynamical and are difficult to account for from their thermodynamic perspective. This may suggest that only considering thermodynamic constraints is not sufficient to address this question. Grabowski and Moncrieff (2004) proposed that the horizontal scale of convective self-aggregation could be a distance over which moisture–convection feedback is effective. This distance was then the product of the large-scale horizontal wind speed (~5 m s^{−1}) and the temporal scale of moisture–convection feedback (~10 days). The horizontal scale was, therefore, about 5000 km, which is a reasonable horizontal scale of convective self-aggregation. Because it is unclear about what sets the horizontal wind speed, this theory is not closed and does not offer further predictability.

In the rest of the paper, we present a new theory for the horizontal scale of convective self-aggregation in 2D (*x*, *z*). In contrast to previous studies focusing on thermodynamics, here we consider both dynamic and thermodynamic constraints on a steady circulation—the equilibrium state of self-aggregation—from a boundary layer perspective. This theory suggests that the length scale of convective self-aggregation is proportional to the product of the boundary layer height and the square root of the amplitude of the density variation between aggregated moist and dry regions in the boundary layer. We show that the density variation mainly arises from the moisture variation due to the virtual effect of water vapor. We further show that the horizontal scale of convective self-aggregation becomes smaller if the virtual effect is disabled in CRM simulations.

## 2. Theory

Figure 1 illustrates the large-scale circulation associated with convective self-aggregation. In the dry zone, the atmosphere cools by radiation, which requires compensating adiabatic heating from subsiding motions. Because of continuity, the atmosphere flow diverges horizontally within the boundary layer, and this horizontal motion requires a horizontal pressure gradient to balance the momentum damping. In this section, we develop a linear boundary layer model for convective self-aggregation and employ the Boussinesq approximation, which assumes the density variation is much smaller than the reference density .

### a. Derivation

The steady state, linear, Boussinesq equations around the domain-mean equilibrium state are

where

*u*,*w*the perturbation horizontal and vertical velocities, respectively (m s^{−1}),*p*the perturbation pressure (Pa),τ the momentum damping time scale (days),

the static stability (K km

^{−1}),*g*the gravity acceleration (m s^{−2}),*c*_{p}the specific heat of dry air at constant pressure (J kg^{−1}K^{−1}),the domain-averaged temperature profile (K), and

*Q*the net heating-rate anomaly (K day^{−1}), which includes both the radiative cooling and convective heating.

We use Rayleigh friction to parameterize momentum damping due to vertical momentum transport and surface drag. This is a widely used parameterization in toy models (Gill 1980; Lindzen and Nigam 1987; Neelin 1989) and has been recently revisited against CRM simulation results (Kuang 2012). Equations (1) and (3) will be integrated over the entire boundary layer, but (2) will only be evaluated at the top of the boundary layer (*z* = *h*), which is defined from the control volume perspective as opposed to a material surface. Here we assume that *h* varies little with *x*, which is tested and discussed in the next section. The horizontal boundary condition is periodic, and the vertical velocity *w* = 0 at surface.

The net diabatic forcing *Q* results from the dynamics and feeds back to it. This interaction could be particularly important for the onset of self-aggregation (e.g., Emanuel et al. 2014). Since (1)–(3) describe a system that has already reached a steady state, we, therefore, do not account for this feedback and assume that *Q*(*x*, *z*) has the same spatial structure as that of the homogeneous part of the equations.

We treat the aggregated circulations as linear waves:

where *k* is the horizontal wavenumber. This implicitly assumes that , which is consistent with spectral analysis results in self-aggregation simulations (Wing and Cronin 2016). We substitute (4) into (1)–(3), integrate (1) and (3) over *z* from 0 to *h*, and solve for the wavenumber *k*. The wavelength of convective self-aggregation is then given by

We define and . So represents the amplitude of the vertically averaged horizontal pressure anomaly, and represents the diabatic forcing anomaly. Without loss of generality, we choose that *δp* is positive and *δQ* is negative.

The horizontal pressure variation within the boundary layer only depends on the density variation within the boundary layer. This is because gravity waves can efficiently smooth out pressure anomalies in the free troposphere (Sobel et al. 2001). For the large-scale steady circulation, the dominant balance in the vertical momentum equation is the hydrostatic balance. The pressure variation within the boundary layer is then given by

where measures the amplitude of density anomaly averaged over the boundary layer. We substitute (6) into (5) and get

where is nondimensional. In the rest of this paper, we propose that the variation of *λ* can be explained by the variations of *h* and , assuming *β* as a constant scaling factor.

Equation (7) is our key result and deserves further analysis. According to the ideal gas law for moist air, the density variation is given by

where *T*_{υ} is the virtual temperature. We compare the two rhs terms and find that in the boundary layer of our simulations. The density variation in the boundary layer is then given by

where and is the specific humidity. Equation (9a) implies that moist air is more buoyant, and this effect is referred to as the virtual effect of water vapor, which arises because the mean molecular weight of water vapor is lower than that of the dry air. Here we assume that the moisture variation dominates over the temperature variation and approximate (9a) by

We then perform the first order expansion of rhs and get , where RH is the relative humidity and is the saturation specific humidity. We find that the second term is much smaller than the first term: . We can simplify (9b) to

### b. Independent parameters and physical interpretation

where is the aspect ratio of the boundary layer and is also a nondimensional measure of the wavelength and . Equation (10) is a nondimensional form of (7), with two independent parameters: and in (10a); and in (10c).

Equation (10) suggests that the aspect ratio will increase with the density variation if remains constant. The aspect ratio of the boundary layer is determined by the ratio between *u* and *w*, which are, in turn, determined by the horizontal pressure gradient and the diabatic heating rate at steady state. The aspect ratio then self-adjusts to the horizontal variations of pressure and diabatic heating. If the diabatic heating rate and its associated vertical mass flux are constant ( is constant), the aspect ratio only depends on , which determines the horizontal pressure variation.

We will test (10a)–(10c) individually, and it would not be a surprise if the quantitative accuracy decreases from (10a) to (10c) owing to further approximations. However, (10c) is more useful to make qualitative predictions. If *γ* remains constant, the only unknown input parameter in (10c) would be the boundary layer temperature, which increases with surface warming. Equation (10c) then explicitly predicts that the aspect ratio increases with SST.

### c. Predictions

Order of magnitude of

*λ*: The typical parameter values in the boundary layer are , , and . Equation (7) then gives .Increasing aspect ratio with SST: this prediction is based on (10c) with constant

*γ.*The virtual effect and

*λ*: this theory predicts that the wavelength*λ*becomes smaller and varies less with SST when the virtual effect is disabled according to (7)–(10).

In the next section, we will use CRM simulations to test the three predictions from our theory.

## 3. Testing theory in CRM simulations

### a. Experiment setup

We study 2D (*x*, *z*) convective self-aggregation with the System for Atmosphere Modeling (SAM, version 6.10.8) (Khairoutdinov and Randall 2003). SAM is an anelastic model. The radiation scheme is identical to that of the National Center for Atmospheric Research (NCAR) Community Atmosphere Model, version 3 (CAM3) (Collins et al. 2006). We turn off the diurnal cycle for simplicity and fix the incoming solar radiation at 413.9 W m^{−2} to match the annual-mean insolation on the equator. The microphysics is the SAM one-moment parameterization.

The horizontal domain size is 12 288 km, and the model top is at 63.6 km. The horizontal resolution is 3 km. The vertical resolution is 50 m for the lowest 1 km and increases to 600 m above 3 km. A sponge layer is employed for the upper 6 km of the model domain. The model is coupled to an ocean surface with a fixed SST, and we perform simulations across a wide range of SSTs from 290 to 325 K. Each simulation is integrated for 120 days, and the results from the last 20 days are presented. This paper presents two groups of simulations: the control group and the other group without the virtual effect of water vapor. Here we switch off the virtual effect by removing the water vapor dependence in the buoyancy equation.

### b. Prediction I: The order of magnitude of λ

Figure 2 shows the Hovmöller diagram of column-integrated precipitable water (PW) in the corresponding simulations. Convective self-aggregation is successfully simulated over all SSTs. Distinct moist and dry patches coexist, and the sizes of moist and dry patches are comparable across this set of simulations.

We first construct composites for convective self-aggregation and then diagnose the wavelength as the horizontal scale of the composite. Here we apply the standard technique of constructing composites for convectively coupled overturning circulations, for example, the mock Walker circulation and MJO (Kuang 2012; Arnold and Randall 2015). First, we filter out short-lived signals by performing moving average in time for PW with the window width of 24 h. The results are robust: different choices of the window width do not significantly affect the moving time-averaged PW. Second, we filter out small-scale (e.g., wavenumber >20 in our simulations) signals by performing Fourier transform. The new PW field slowly varies both in time and space, so that the aggregation pattern can easily stand out. Third, we identify the center of dry patches as a local minimum of PW that deviates from the domain-averaged PW by 0.5 standard deviations. We then align all the dry centers together and average them to make the composite. This final step is repeated for all fields using the same dry-center index from the PW field. The wavelength is then diagnosed as the distance between neighboring maximum boundary layer convergence.

### c. Prediction II: Increasing aspect ratio with SST

Testing this prediction requires testing the proposed scaling in (10), which then requires diagnosing *h* and .

We diagnose *h* in the dry and moist patches separately by calculating the altitude where exceeds −0.2 km^{−1} above the level of the minimum gradient (Fig. 4a). In a given simulation, we find that the difference in *h* between the moist and dry patches is much smaller than the reference value. This allows us to use a fixed boundary layer height in our theory. We also find that the horizontal pressure difference at *h* is much smaller than that in the boundary layer (Fig. 4b). This confirms that the horizontal pressure difference in the boundary layer is mainly determined by the density difference in the boundary layer [see (6)].

Figure 4c plots *h* against SST. The boundary layer height almost monotonically decreases with SST, which is consistent with previous simulation results (Takahashi 2009; Wing and Cronin 2016). The decrease of *h* might be due to decreasing buoyancy flux with surface warming, which is mainly associated with the decrease of surface sensible heat flux. The ratio between the sensible heat flux and latent heat flux is defined as the Bowen ratio, which decreases with warming (Mitchell et al. 1987). The latent heat flux slowly increases with warming owing to the energy constraint (Held and Soden 2006). As a result, the sensible heat flux decreases with warming, leading to weaker buoyancy flux and thereby shallower boundary layer in warmer climates (Tan et al. 2017).

We define as the amplitude of the density difference between the dry and moist centers and calculate using two methods. In the first method, we vertically integrate the equation of state [see (9)] to get . In the second method, we diagnose and use the hydrostatic balance [see (6)] to get . Because the second method takes into account of the small pressure perturbation at *z* = *h*, comparing results from the two estimates will help further test the weak pressure gradient approximation at the top of the boundary layer. The boundary layer height is a free parameter in calculating . To achieve robust results, we calculate by using *h* in the moist and dry areas and their mean values, respectively. The two estimates agree very well and both increase with surface warming (Fig. 5a). This good agreement suggests that the pressure perturbation at *z* = *h* is negligible so that the boundary layer pressure variation is mainly set by the density variation. This increase with warming may be due to increasing moisture variation. Figure 5b shows that the moisture variation dominates the density variation, especially in the warm climates. Because the moisture variation increases with warming, the density variation also increases with warming. This is consistent with our approximation in deriving (9b). These results are robust over all choices of *h* (Figs. 5a,b).

We further analyze and show that the increasing moisture variation is mainly due to the increasing saturation specific humidity in the boundary layer. Figure 5c plots , , and against SST, where 0.25 is the typical in our simulations and is calculated using temperature and pressure in the boundary layer according to the Clausius–Clapeyron relation. We find that the three estimates are consistent and that accurately reproduces the moisture variation, except for the warmest case.

Figure 6 tests our final scaling in (10). Because we have verified all main assumptions and approximations, we expect good agreement between the proposed scaling theory and simulation results. Figure 6a plots against , where , and has been varied by over an order of magnitude. All data points fall along a one-to-one line, suggesting that the density variation determines the aspect ratio of the boundary layer. Figures 6b and 6c test the approximate theory [see (10b) and (10c)]. Although some data points deviate from the one-to-one line owing to further approximations, the approximate theory still provides reasonable estimation for the aspect ratio.

The aspect ratio increases with SST in the 305–325-K simulations (Figs. 6a–c). In this SST range, the density variation can be accurately reproduced by the moisture variation (Fig. 5), which increases with SST mainly because of increasing saturation specific humidity . This increase of then leads to the increasing aspect ratio with SST according to (10c). In the 290–305-K simulations, both the density variation and aspect ratio remain fairly constant with SST (Fig. 6a), suggesting that *h* dominates the variation of *λ*.

Our results also explain why the wavelength varies nonmonotonically with SST (Fig. 3). For the 290–310-K simulations, *h* decreases linearly at a relatively fast rate, and increases moderately since temperature variation competes with moisture variation, as illustrated by the distance between the dashed line and the markers in Fig. 5b. For the 310–325-K simulations, *h* remains almost constant, and the moisture variation dominates the density variation, which increases exponentially. As a result, *λ* first decreases and then increases with surface warming.

### d. Prediction III: The virtual effect and λ

The role of the virtual effect is already illustrated by Fig. 5. Here, to provide causal evidence, we perform mechanism-denial experiments. We remove the buoyancy dependence on water vapor (i.e., switching off the virtual effect) and then repeat all our simulations. Figure 7 shows the Hovmöller diagram of normalized PW of the new simulations. The atmosphere still self-organizes into moist and dry patches, but the scale does not vary with SST as much as in the control simulations.

Figure 3 plots and compares the two sets of simulations. The new simulations show a smaller wavelength at a given SST, and the variation of wavelength is within a factor of 2, which is significantly reduced. This is because both and *h*, as well as their variations over the 35-K SST range, are reduced when the virtual effect is disabled. The decrease of is due to removing the virtual effect, and this effect is dominant in the warm simulations. The decrease of *h* is likely due to the reduced buoyancy flux, which could also result from removing the virtual effect. The decrease of *h* is dominant in the cold simulations, in which the boundary layer is deeper (Fig. 4a). These results provide causal evidence for the hypothesis that the virtual effect helps set the horizontal scale of convective self-aggregation.

## 4. Discussion

This paper presents a linear boundary layer theory for the horizontal scale of convective self-aggregation *λ*, a nondimensional form of which is the aspect ratio of the boundary layer. The mass balance of the boundary layer suggests that the aspect ratio is determined by the horizontal mass flux and the vertical mass flux. Because the vertical mass flux varies weakly with SST , the aspect ratio is then only a function of the horizontal mass flux, which is determined by horizontal pressure or density variation. This gives , where is significantly influenced by the virtual effect of water vapor. We then further simplify this relation and get , in which the boundary layer temperature is the only unknown parameter. This theory does not only have explanatory power but also successfully predicts the following: 1) *λ* ~ *O*(2000) km, 2) *λ*/*h* increases with SST (in the 305–325-K simulations), and 3) *λ* and its variation with SST both decrease when the virtual effect is switched off.

Prediction 2 of our theory depends on the constraint that or remains constant with SST. The variation of over the entire SST range is on the order of 50%. Assuming is, therefore, appropriate, although it may be responsible for small deviations from the one-to-one line in the cold simulations (Fig. 6a). The slow variation of is partly because the heating anomaly at the boundary layer top varies much slower with SST than water vapor and its associated contribution to the density variation. Because boundary layer clouds are not properly represented in CRM simulations, it would be worthwhile to examine if the heating anomaly and still vary slowly in simulations with proper representations of boundary layer clouds (e.g., with higher resolutions or with a shallow cloud parameterization). This would help link our theory to more realistic convective organizations.

Prediction 3 would suggest that *λ* decreases more in warm simulations when the virtual effect is switched off. However, *λ* decreases significantly in the cold simulations as well. This is because *h* decreases significantly in cold simulations, which is likely due to the reduced buoyancy flux. This is intuitive but has not been included in our theory. A more complete theory then requires quantitative understanding of what determines the boundary layer height.

Although highly idealized, this theory could be a good starting point for understanding convective organizations in 3D atmospheres and may have two potential implications. One is to help explain why 3D convective self-aggregation might be easier to be simulated around 300-K SST and with a sufficiently large horizontal domain (Muller and Held 2012; Emanuel et al. 2014). If the 3D aggregation has similar temperature dependence to 2D aggregation, the natural horizontal scale for convective self-aggregation varies with SST and reaches its minimum value around 300 K. A successful simulation may require a sufficiently large CRM domain—at least comparable to *λ*, and the domain-size requirement should vary with SST accordingly. For simulations with SST around 300 K, the domain size *D* can be small owning to small *λ*. For simulations with SST either much higher or much lower than 300 K, a successful simulation requires a much larger domain size than the typical domain size in 3D CRM simulations (about a few hundred kilometers). This might be the reason why convective self-aggregation did not occur over either very high or very low SSTs in Wing and Emanuel (2014). However, a recent study suggested that convective self-aggregation can be simulated in a snowball Earth configuration with interactive SST (Abbot 2014). This may suggest that other physical processes could also be important. For example, may have different values in 2D and 3D simulations.

The other potential implication is to provide an upper limit to the size of mesoscale convective systems (MCSs). It is well established that the sustainability—the area that can sustain convection—limits the size of MCSs (Houze 2004). The moist patch favors convection, and therefore its size sets the maximum scale of MCSs in our simulations. Because moist and dry patches have similar scales, the sizes of moist patches in the 295- and 300-K simulations are ~1000 km. This is in theory the maximum size of MCSs in the simulations at SSTs similar to the current tropical SST. This scale could offer an upper bound to the largest MCSs observed in the tropical atmosphere (Houze 2004). A more accurate estimation of this upper bound would require detailed diagnosis of the boundary layer properties of MCSs. This hypothesis is based on our 2D simulation results, and it would be worthwhile to further explore this idea using 3D CRM simulations.

## Acknowledgments

This work was supported by Laboratory Directed Research and Development (LDRD) funding from Berkeley Lab, provided by the Director, Office of Science, of the U. S. Department of Energy under Contract DE-AC02-05CH11231. The author was also supported by a Miller Research Fellowship at University of California, Berkeley. The author thanks N. Arnold, J. Edman, M. Herman, A. Ingersoll, N. Jeevanjee, J. Seeley, Z. Tan, and anonymous reviewers for constructive comments on an earlier draft of this paper. The author thanks C. Bretherton, K. Emanuel, Z. Kuang, C. Muller and D. Romps for helpful discussions.

## REFERENCES

^{2}climate sensitivity and model dependence of results

## Footnotes

© 2018 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).