Cloud_ICA: A deterministic cloud-overlap algorithm for generating a complete set of independent column atmospheres

In calculating solar radiation, climate models make many simplifications, in part to reduce computational cost and enable climate modeling, and in part from lack of understanding of critical atmospheric information. Whether known errors or unknown errors, the community's concern is how these could impact the modeled climate. The simplifications are well known and most have published studies evaluating them, but with individual studies it is difficult to compare. Here, we collect a wide range of such simplifications in either radiative transfer modeling or atmospheric conditions and assess potential errors within a consistent framework on climate‐relevant scales. We build benchmarking capability around a solar heating code (Solar‐J) that doubles as a photolysis code for chemistry and can be readily adapted to consider other errors and uncertainties. The broad classes here include: use of broad wavelength bands to integrate over spectral features; scattering approximations that alter phase function and optical depths for clouds and gases; uncertainty in ice‐cloud optics; treatment of fractional cloud cover including overlap; and variability of ocean surface albedo. We geographically map the errors in W m−2 using a full climate re‐creation for January 2015 from a weather forecasting model. For many approximations assessed here, mean errors are ∼2 W m−2 with greater latitudinal biases and are likely to affect a model's ability to match the current climate state. Combining this work with previous studies, we make priority recommendations for fixing these simplifications based on both the magnitude of error and the ease or computational cost of the fix.

Methods
A unique and useful advantage of Cloud-ICA is that with low computational cost it readily generates a complete set of ICAs such that the sum of wtICAs equals one.  Cloud-ICA can then calculate a set of up-to-four quadrature column atmospheres (QCAs) that can be used to approximate the integral over the full set of ICAs.  The method of parsing partially overlapping clouds within a single column atmosphere (SCA) into a large number of ICAs and then into a set of only four QCAs originated in Neu, Prather & Penner (2007) and was fully developed with observation-based overlap models in the chemistry photolysis code Cloud-J v7.3c in Prather (2015).  Cloud-J has continued to develop (Prather & Hsu, 2019;  2019d; Hsu & Prather, 2021) and the current Cloud-ICA code is taken from Cloud-J v8.0c (Prather, 2023). 

Note added (2024/05/20):  The program as written limits the number of ICAs to the dimension size ICA_ = 20,000 = 2 x 104.  If using the preferred 6 max-overlap vertical blocks plus cirrus shield (LNRG=6) and the quantization of cloud fraction in units of 1/10, then the maximum number of ICAs is limited to about 5 x 106.   While this maximum number has not been found by us with realistic atmospheres, it is possible. We did not want to dimension and scan/sort such large numbers of ICAs when running multi-year CTM simulations and thus truncated it at 2 x 104.  This truncation is a soft landing and just combines upper level decorrelated blocks into a single max-overlap block until the ICA_ limit is met.  The invocation of this limit does invoke a 'write' statement that may disrupt parallel optimization, and it can be commented away.

NICA greater than ICA_       22680       20000        7560           7           6

Beware that if you run LNRG=0, then alternating cloudy and clear layers can generate 2CLD_/2 ICAs.

The primary Cloud-ICA references are:

Neu, J.L., M.J. Prather, J.E. Penner (2007) Global atmospheric chemistry:  integrating over fractional cloud cover, J. Geophys. Res., 112, D11306, doi:10.1029/2006JD008007

Prather, M.J. (2015) Photolysis rates in correlated overlapping cloud fields: Cloud-J 7.3c, Geosci. Model Dev., 8, 2587-2595, doi:10.5194/gmd-8-2587-2015

Prather, M.J. and J.C. Hsu (2019) A round Earth for climate models, Proc Natl Acad Sci, 116 (39) 19330-19335; https://doi.org/10.1073/pnas.1908198116.

Prather, Michael and Hsu, Juno (2019d), Solar-J and Cloud-J models version 7.6c, v2, UC Irvine Dash, Dataset, doi.org/10.7280/D1096P

Hsu, J.C. and M.J. Prather (2021), Assessing uncertainties and approximations in solar heating of the climate system. Journal of Advances in Modeling Earth Systems, 13, e2020MS002131. doi: 10.1029/2020MS002131

Prather, Michael (2023d). An updated cloud-overlap photolysis module for atmospheric chemistry models, UCI Cloud-J v8.0, with near-UV H2O absorption [Dataset]. Dryad. https://doi.org/10.7280/D1Q398

A description of the cloud overlap algorithm and it basis in observations is take directly from Prather (2015), see references within.

“Our recommended cloud-overlap model uses the information on vertical correlations (Pincus et al., 2005; Naud and DelGenio, 2006; Kato et al., 2010; Oreopoulis et al., 2012), which shows cloud decorrelation lengths on the order of 1.5 km in the lower atmosphere increasing to 3 km or more in the upper troposphere. Since a true COR model scales as 2**NL and becomes rapidly impractical for high-resolution models, we define vertical groups of cloud layers globally according to the decorrelation lengths: 0–1.5 km altitude, 1.5–3.5, 3.5–6, 6–9, 9–13, and >13 km. We assume that the cloud layers within a decorrelation length are highly correlated with one another and thus form a MAX group. When such MAX groups are adjacent they have a mean separation of one decorrelation length, and we choose a cloud correlation factor of cc=0.33, similar to 1 e-fold. When there is a clear-sky gap between a pair of G6 layers, the MAX groups are separated by more than one decorrelation length; thus, we reduce the factor cc with successive multiples (i.e., with two missing G6 MAX groups between two cloudy layers,the effective cc=0.333**2 =0.036). This model is denoted G6/.33. Two other G6 models were tested: cc=0.00 corresponds to randomly overlapped adjacent groups (MAXRAN, G6/.00); and cc=0.99 is almost maximally overlapped (MAX, G6/.99).

In looking at how this model aligned the clouds for realistic FCAs, we found that extensive cirrus fractions in the uppermost layers prevented the expected overlap of small fraction cumulus below. Thus, a seventh MAX group is added if there was a cirrus shield (defined from top down as adjacent ice-only clouds with f >0.5). Because of the cloud fraction binning into 10% intervals, the number of ICAs is bounded by 5x106 (including the cirrus shield). This limit is resolution independent and was never reached in any FCAs examined here (highest number of ICAs for one FCA was 3500). The major computational cost comes with the Fast-J computation, and the methods for approximating the average of J values over all ICAs (Sect. 3) use at most four Fast-J calculations no matter how many ICAs.  Two other cloud-overlap models tested here are the MAXRAN groupings G0 and G3 (Feng et al., 2004; Neu et al., 2007). Model G0 assumes that all vertically adjacent cloudy layers are a MAX group (maximally overlapped), and all such groups separated by a clear layer are RAN overlapped.  This model seems logical but has difficulty finding a clear layer when the FCA has been averaged over several hours or taken from a parameterized cloud-resolving model.  It our tests, using meteorological data with NLD36, the maximum number of G0 ICAs was 375. Model G3 has at most three MAX-RAN groups demarcated by atmospheric regimes: a fixed altitude (1.5 km, stratus top) and temperature (the liquid-to-ice cloud transition). The maximum possible number of ICAs per FCA for G3 is 103, and in our tests we found 288.

Our recommended cloud-overlap model is G6/.33 since it is based on the observed–modeled cloud decorrelation lengths. For a given FCA, we treat the J values calculated by summing Fast-J over all the ICAs generated by G6/.33 as the correct value.  We calculate errors for the other cloud-overlap models (here) or various ICA-approximation models using the G6/.33 model (Sect. 3).

We use a high-resolution snapshot from the European Center for Medium-rangeWeather Forecasts, similar to what is used (at lower resolution) in the UC Irvine and University of Oslo chemistry-transport models (Søvde et al., 2012; Hsu and Prather, 2014). The 640 FCAs are a 3 h average of a single longitudinal belt just above the Equator (T319L60 Cycle 36) and have clouds only in the lowermost 36 layers. Profiles of temperature and ozone are taken from tropical mean observations; the Rayleigh-scattering optical depth at 600 nm is about 0.12, and a mix of aerosol layers has a total optical depth of 0.23. J value errors are calculated separately for each FCA and then averaged. The number of ICAs per FCA averages 169 for model G6, 21 for model G3, and 19 for model G0; see Fig. 2 for the probability distribution of ICA numbers.

Further discussion about the deterministic cloud cover generator is found in Hsu and Prather (2021)

“In a manner similar to Hogan and Bozzo's (2018) deterministic cloud-cover generator that goes from MAXRAN to EXP-RAN, Cloud-J developed a deterministic ICA generator for MAX-RAN and then adapted it to use vertical decorrelation lengths in its MAX-COR algorithm (M. J. Prather, 2015). Chemistry models need the selection of ICAs for any overlap method to be deterministic because many critical applications require perturbation-control pairs without stochastic noise (e.g., M. J. Prather & Hsu, 2010). Thus Solar-J cannot use a stochastic cloud generator (e.g., Räisänen et al., 2004), and this drove the structure of our cloud overlap algorithm. MAX-COR was designed to be (i) deterministic, (ii) linear in cost with increasing numbers of layers, and (iii) robust when cloud data are averaged in time or space, because such averaging tends to eliminate cloud-free layers and revert to MAX overlap. Based on observations of decorrelation length (Kato et al., 2010; Naud et al., 2008; Oreopoulos et al., 2012; Pincus et al., 2005), MAX-COR defines 6-layer groupings by altitude range. Because decorrelation is small across the vertical range of each group, we assume MAX overlap within each group and a decorrelation of the overlap of each MAX group with its neighbor. Adopting terminology of climate community, MAX-COR is effectively a MAX-EXP algorithm. By quantizing the cloud fraction to the nearest 10% and allowing an independent cirrus shield at the top, the absolute maximum number of ICAs under MAX-EXP is <5 × 106 and thus ICAs can be rapidly defined and binned with low computational overhead. Deterministic EXP-EXP or EXP-RAN models in our code would have to enumerate up to 233 ICAs for our model that has potentially 33 cloudy layers, which is truly prohibitive and not linearly scalable with resolution. We believe that a MAX-COR or MAX-EXP algorithm is likely the most stable and scalable deterministic ICA generator for vertical cloud decorrelation algorithms. The RRTMG v4.0 code available at the time of this study uses primarily MAX-RAN cloud overlap, but the new v5.0 code includes an EXP-RAN option. Thus, our comparisons of cloud-overlap results with the RRTMG code are limited to MAX-RAN. Within Solar-J we can run both MAX-RAN (SJ/RAN) and the standard MAX-COR (SJ) and thus compare with J. K. P. Shonk and Hogan (2010), as discussed below.

Let us accept that ICAs generated by cloud overlap algorithms can be solved with 1D RT as horizontally homogeneous plane parallel layers, then the next step is how to solve the RT problem for all ICAs and average the results. The number of ICAs are often numerous enough that no practical climate RT code can solve them all, and most codes do not even count them all (Räisänen et al., 2004). RRTMG randomly selects an ICA for each wavelength bin in the RT solution, a method designated Monte Carlo ICA McICA, Pincus et al., 2003). McICA has errors at each time step by mixing ICAs across wavelengths and by not accurately sampling the average of ICAs (e.g., average cloud optical depth) in that time step. McICA is intended to deliver the correct mean when averaged enough times over the same cloud system, but it has hourly grid-cell rms errors of 40 W m−2 (H. W. Barker et al., 2008; Pincus et al., 2003). A key underlying premise is that solar heating errors propagate symmetrically and linearly in the climate system and average out, as was found for simple forecast models. Assessing net bias errors caused by noisy heating rates would need to examine nonlinear processes in hydrology, cloud systems, ecosystem productivity, and air quality in Earth system models (e.g., Pincus & Stevens, 2013).

With a deterministic ICA generator, we can calculate an "exact" non-stochastic answer as was done for limited test cases in M. J. Prather (2015), but we could not afford to do this for our January climate metric.  Solar-J identifies and sorts all ICAs by cloud optical depth and then selects up to four representative quadrature column atmospheres (QCAs) each with a fractional area to represent the distribution of ICAs. The full-wavelength RT solutions are completed for each QCA (Neu et al., 2007). See Figure S3 for a global picture of the average frequency of occurrence of the 4 QCA bins for January 2015. Cloud quadrature does a very good job of averaging over the ICAs with net bias errors of ∼1% in solar intensity and rms errors of 2%–4%. To reach equivalent accuracy for a single time step using random selection would require about 50 ICAs each with full wavelength calculation (not as in McICA) versus an average of 2.8 QCAs (many grid cells have less than 4 QCAs).