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

Prather, Michael1

Research facility: University of California, Irvine

Published May 15, 2024; Updated May 21, 2024 on Dryad. https://doi.org/10.5061/dryad.w9ghx3fxz

Data files

May 15, 2024 version files 1.75 MB

atmos_PTClds.dat
1.66 MB
ICA_nica640.out
52.48 KB
ICA-print1.out
26.61 KB
README.md
7.80 KB

Abstract

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.

README: Assessing Uncertainties and Approximations in Solar Heating of the Climate System

https://doi.org/10.5061/dryad.w9ghx3fxz

The primary goal of the computer code published here is to provide an algorithm for sorting an atmosphere of overlapping clouds into a complete and finite set of independent column atmospheres (ICAs), wherein each layer is horizontally homogenous (i.e., plane parallel) and either uniformly cloudy or clear. Many current solar heating codes can only estimate the distribution of ICAs with Monte Carlo sampling. This algorithm gets around that problem. It has been in use with UCI photochemistry codes, but it was so deeply embedded as to make it unusable by others. Here we extract that algorithm to enable ready use in other Earth system modeling codes.

This work presents a stand-alone FORTRAN 90 computer code Cloud-ICA that reads in a single column atmosphere (p at edges, T in layer) where each layer includes a cloud fraction (CldF = 0-1) plus the liquid- and/or ice-water path of the cloud in the layer (not the average over the clear and cloudy layer). The cloud data are used to calculate the visible optical depth of the cloud itself (CldTau), not the average optical depth of the layer. With these two quantities and an overlap algorithm, the code calculates a subset of independent column atmospheres (ICAs), each occupying a fraction of the original cell (wtICA) and each with layers that are either 100% cloudy or clear, thus enabling plane-parallel radiative transfer solutions.

Description of the data and file structure

There are two sample output files:

ICA_print1.out = detailed printout for the first cloudy atmosphere that has only 36 ICAs. This gives the values for the grouping indices and the profiles of the 36 ICAs in terms of cloud visible optical depth in leach model layer. See the rest of this readme file for a description of the variables.

ICA_nica640.out = a quick scan through 640 different tropical atmospherrs listing the number of ICAs in each (NICA) the resulting QCAs in terms of weights and column optical depth.

atmos_PTClds.dat = gives a standard atmospheric profile (p, T) and then 640 different profiles of cloud systems, sampled from ECMWF fields in the tropics.

Variables in the printout examples

>>>>>atmospheric data

L = model layer, must be from surface(L=1) to top-of-atmosphere (L=58)in this case. The inverse layering notation (L=1=toa)is not supported and may not work.

For physics, in our model L=1:LATM_

For clouds, it is limited to L=1:CLD_ (<LATM_)

z (km), geometric altitude at the bottom of the layer from integrating the hydrostatic equation

p (hPa) at the bottom edge of the layer

T (K) temperature of the layer.

CLF (dimensionless, 0.00 to 1.00) = cloud fraction in layer L = 1:LCLD, aka CLDF in some modules

LWP/IWP (g/m2) = liquid/ice water path in the cloud (not average in the layer)

WLC/WIC (g-water/kg-air) = cloud liquid/ice water content

CLTAU (optical depth)= approximate visible OD for cloud (not average in layer)

REFFL/REFFI (microns) = effective radius of liquid/ice water cloud used to calculate CLTAU

CLDCOR (0.00 to 1.00) decorrelation factor between adjacent maximally overlapped cloud groups. The decorrelation factor with a cloud group 2 above it is CLDCOR**2

>>>>>cloud overlap groups and quantization

FBIN (integer) number of quantized bins allowed for CLDF. Example, if FBIN = 10, then CLDF is quantized into fixed fractions: 0.0, 0.1, 0.2, … 0.9, 1.0. Note that if CLF < CLF_MIN = 0.02) then CLDF = 0.0 and the layer is cleared of clouds! The optical depth of the cloud (CLTAU) is scaled to conserve CLF*CLTAU. Quantized CLDF means any max-overlap group has at most 10 possible combinations.

NCLDF(1:LCLD_) = cloud fraction quantum number (1:FBIN) for layer. NCLDF = 0 means NO clouds

LNRG (= 0, 3, 6) index for Cloud-J to establish the number of semi-independent maximally overlapped cloud groups in an atmospheric column.

LNRG=0. Defines a break in the max-overlap cloud groups by the cloud optical depth quantum. NG_BRK=0 is default = clear layer, see FBIN and NCLDF below.

LNRG=3: Defines 3 separate max-overlap groups strictly by altitude, that is presently locked to our ECMWF L57 grid: GRP=1 is L=1:8 (surf to +1 km); GRP=2 is L=9 to uppermost liquid water cloud; GRP=3 is for all-ice clouds. Note that mixed-phase clouds count as liquid water since they dominate the liquid water cloud opacity. These groups may be merged if they are without clouds.

LNRG=6: Current recommendation is for correlated groups of max-overlap clouds. Uses altitude (actually calculated for the true model grid, not presumed to be ECMWF L57) to create 6 max-overlap groups, based on observed correlation lengths (see Cloud-J papers). GRP#1 is 0-1.5km, GRP#2 is 1.5-3.5km, GRP#3 is 3.5-6km, GRP#4 is 6-9km, GRP#5 is 9-13km, and GRP#6 is 13km+.

LCIRRUS = layer number of the top cirrus shield group (if there is one). New Group: if there is a CLF~1 cirrus shield, pull it off the top of the maximum GRP, and create a separate GRP#7. (Could be lower GRP# is some of GRP#1-6 have no clouds.

NRG (integer) final number of maximally overlapped groups, should be 1 to 7 in current code depending on numer of groups w/o clouds. Code is dimensioned to accommodate 9 total groups.

GMIN(1:NRG), GMAX(1:NRG) = layer number of the bottom and top of each max-ran group. All possible cloud layers (1:LCLD_) must be in a group. The numbers start with GMIN(1)=1 and end with GMAX(NRG) = LCLD_.

GNR(1:NRG) = no. of unique quantized cloud fractions in group, .le.FBIN+1

GBOT(1:NRG) = lowest model layer of max-overlap group (similar to GMIN)

GTOP(1:NRG) = uppermost model layer of max-overlap group (similar to GMAX)

GCMX(1:NRG) = maximum cloud fraction (quantum # 1:FBIN) of the max-overlap group

GLVL(1:NRG) = G6 group # (1:7,for LNRG=6 only) of the NRG group #. This keeps track of a missing (uncloudy) G6 group so that the decorrelation factor is calculated correctly, with extra factor for missing groups. Note that there can be 7 G6 groups if the upper cirrus shield is separated so as to be decorrelated with clouds below.

GFNR(1:NRG,1:GNR) = GFNR(N,G) = cloud fraction quantum no. (1 to CBIN_) for max-overlap group N’s cloud fraction # G.

>>>>>independent column atmospheres

ICA_ = dimensioned size for ICAs, currently 20,000. Can be larger. With LNRG=6 and FBIN=10, the maximum possible is about 5e6, but for LNRG=0, the maximum could be 2**(LCLD_/2) (i.e., alternating cloudy and clear layers).

NICA = total number of ICAs

WCOL(1:NICA) = fractional area weight of each ICA

OCOL(1:NICA) = column visible optical depth of each ICA

OPRF(1:LCLD_,1:NICA) = visible optical depth in each layer of each ICA

ISORT(1:NICA) = number of ICA (original sequence generated by groups) when sorted in order from smallest OCOL to largest OCOL.

>>>>>quadrature column atmospheres

N.B. These QCA variables are not returned in the main subroutine CLOUD-ICA because QCAs may not necessarily be invoked (need CLDFLAG = 7 or 8). In Cloud_j, the QCA calculation is called internally. These can be extracted if the QCAs are desired.

WTQCA(1:NQD_) = fractional area weight of QCA = 1:4, may be 0.

ODQCA(1:NQD_) = column cloud visible optical depth of QCA

QCOLN(1:LCLD_,1:NDQ_) = visible optical depth in each layer of each ICA

NQ1(1:NQD_) = starting ICA number in the QCA group (based on ISORT, sorted ICAs!)

NQ2(1:NQD_) = ending ICA number in the QCA group (from ISORT)

Code/Software

The core computer program and standalone driver is archived in the accompanying zenodo files. These are FORTRAN 90 files, written on a pc (cr/lf).

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 10⁴. 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 10⁶. 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 10⁴. 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 2^CLD_/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).