# Subgrid-scale variability of cloud ice in the ICON-AES
1.3.00

Sabine Doktorowski, Jan Kretzschmar, Johannes Quaas, Marc Salzmann, Odran Sourdeval

► To cite this version:

Sabine Doktorowski, Jan Kretzschmar, Johannes Quaas, Marc Salzmann, Odran Sourdeval. Subgrid-scale variability of cloud ice in the ICON-AES 1.3.00. Geosci. Model Dev., 2024, Geosci. Model Dev., 17, 10.5194/gmd-17-3099-2024. hal-04601719

HAL Id: hal-04601719

https://hal.univ-lille.fr/hal-04601719

Submitted on 5 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Subgrid-scale variability of cloud ice in the ICON-AES
1.3.00

Sabine Doktorowski, Jan Kretzschmar, Johannes Quaas, Marc Salzmann, Odran Sourdeval

► To cite this version:

Sabine Doktorowski, Jan Kretzschmar, Johannes Quaas, Marc Salzmann, Odran Sourdeval. Subgrid-scale variability of cloud ice in the ICON-AES 1.3.00. Geosci. Model Dev., 2024, Geosci. Model Dev., 17, 10.5194/gmd-17-3099-2024. hal-04601719

HAL Id: hal-04601719

https://lilloa.hal.science/hal-04601719v1

Submitted on 5 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Subgrid-scale variability of cloud ice in the ICON-AES 1.3.00

Sabine Doktorowski<fnref n="1" />, Jan Kretzschmar<fnref n="1" />, Johannes Quaas<fnref n="1" />, Marc Salzmann<fnref n="1" />, and Odran Sourdeval<fnref n="2" />

**Correspondence:** Sabine Doktorowski (sabine.hoernig@uni-leipzig.de)

Received: 7 February 2022 – Discussion started: 9 November 2023

Revised: 9 November 2023 – Accepted: 7 March 2024 – Published: 19 April 2024

**Abstract.** This paper presents a stochastic approach for the aggregation process rate in the ICOsahedral Nonhydrostatic general circulation model (ICON-AES), which takes subgrid-scale variability into account. This method creates a stochastic parameterization of the process rate by choosing a new specific cloud ice mass at random from a uniform distribution function. This distribution, which is consistent with the model's cloud cover scheme, is evaluated in terms of cloud ice mass variance with a combined satellite retrieval product (DARDAR) from the satellite cloud radar CloudSat and the Cloud-Aerosol Lidar and Infrared Pathfinder Observations (CALIPSO). The global patterns of simulated and observed cloud ice mixing ratio variance are in a good agreement, despite an underestimation in the tropical regions, especially at lower altitudes, and an overestimation in higher latitudes from the modeled variance. Due to this stochastic approach the yearly mean of cloud ice shows an overall decrease. As a result of the nonlinear nature of the aggregation process, the yearly mean of the process rates increases when taking subgrid-scale variability into account. An increased process rate leads to a stronger transformation of cloud ice into snow and therefore to a cloud ice loss. The yearly averaged global mean aggregation rate is more than 20% higher at selected pressure levels due to the stochastic approach. A strong interaction of aggregation and accretion, however, lowers the effect of cloud ice loss due to a higher aggregation rate. The new stochastic method presented lowers the bias of the aggregation rate.

## 1 Introduction

A correct representation of clouds and cloud-related microphysical process rates, which describes the time-dependent source and sink terms of cloud ice or liquid water, is one of the central challenges in global climate modeling. Since global climate models typically run at a rather coarse resolution (of the order of 100 km), it is important to look into the unresolved microphysical process rates. Most climate models use a cloud cover parameterization. Microphysical process rates are computed based on grid box mean in-cloud ice or liquid water mixing ratios. The in-cloud cloud ice mass mixing ratio is the cloud ice mass per cloudy area. Considering subgrid-scale variability of in-cloud variables reduces the biases of the nonlinear microphysical process rates (Pinус and Klein, 2000; Larson et al., 2001). For example, Weber and Quaas (2012) numerically integrated the process rate over the probability density function (PDF), which is a very accurate method but needs additional computational time. Another method, which works with no additional computational time, is a stochastic approach for the process rates by taking a randomly chosen value per time step and grid box (e.g., Palmer, 2001; Berner et al., 2017). With this method a randomly disturbed process rate is created in order to give a better representation of the state of the atmosphere. Another method is the Cloud Layers Unified By Binormals (CLUBB) (Golaz et al., 2002a, b), which works with a set of PDFs for all cloud types to avoid difficulties in coupling of stratocumulus and shallow convection parameterizations. Since CLUBB does not parameterize subgrid-scale variability of cloud ice, Thayer-Calder et al. (2015) include cloud ice in the CLUBB PDFs. The PDFs are sampled by the subgrid importance Latin hypercube sampler (SILHS), which connects CLUBB with the microphysics for stratiform and convective clouds.

<footnote n="1">¹Institute for Meteorology, Universität Leipzig, Leipzig, Germany</footnote>

<footnote n="2">²Laboratoire d'Optique Atmosphérique, Université de Lille, Lille, France</footnote>
<figure><img src="image_2.png" /></figure>

# Subgrid-scale variability of cloud ice in the ICON-AES 1.3.00

Sabine Doktorowski<fnref n="1" />, Jan Kretzschmar<fnref n="1" />, Johannes Quaas<fnref n="1" />, Marc Salzmann<fnref n="1" />, and Odran Sourdeval<fnref n="2" />

**Correspondence:** Sabine Doktorowski (sabine.hoernig@uni-leipzig.de)

Received: 7 February 2022 – Discussion started: 9 November 2023

Revised: 9 November 2023 – Accepted: 7 March 2024 – Published: 19 April 2024

**Abstract.** This paper presents a stochastic approach for the aggregation process rate in the ICOsahedral Nonhydrostatic general circulation model (ICON-AES), which takes subgrid-scale variability into account. This method creates a stochastic parameterization of the process rate by choosing a new specific cloud ice mass at random from a uniform distribution function. This distribution, which is consistent with the model's cloud cover scheme, is evaluated in terms of cloud ice mass variance with a combined satellite retrieval product (DARDAR) from the satellite cloud radar CloudSat and the Cloud-Aerosol Lidar and Infrared Pathfinder Observations (CALIPSO). The global patterns of simulated and observed cloud ice mixing ratio variance are in a good agreement, despite an underestimation in the tropical regions, especially at lower altitudes, and an overestimation in higher latitudes from the modeled variance. Due to this stochastic approach the yearly mean of cloud ice shows an overall decrease. As a result of the nonlinear nature of the aggregation process, the yearly mean of the process rates increases when taking subgrid-scale variability into account. An increased process rate leads to a stronger transformation of cloud ice into snow and therefore to a cloud ice loss. The yearly averaged global mean aggregation rate is more than 20% higher at selected pressure levels due to the stochastic approach. A strong interaction of aggregation and accretion, however, lowers the effect of cloud ice loss due to a higher aggregation rate. The new stochastic method presented lowers the bias of the aggregation rate.

## 1 Introduction

A correct representation of clouds and cloud-related microphysical process rates, which describes the time-dependent source and sink terms of cloud ice or liquid water, is one of the central challenges in global climate modeling. Since global climate models typically run at a rather coarse resolution (of the order of 100 km), it is important to look into the unresolved microphysical process rates. Most climate models use a cloud cover parameterization. Microphysical process rates are computed based on grid box mean in-cloud ice or liquid water mixing ratios. The in-cloud cloud ice mass mixing ratio is the cloud ice mass per cloudy area. Considering subgrid-scale variability of in-cloud variables reduces the biases of the nonlinear microphysical process rates (Pinус and Klein, 2000; Larson et al., 2001). For example, Weber and Quaas (2012) numerically integrated the process rate over the probability density function (PDF), which is a very accurate method but needs additional computational time. Another method, which works with no additional computational time, is a stochastic approach for the process rates by taking a randomly chosen value per time step and grid box (e.g., Palmer, 2001; Berner et al., 2017). With this method a randomly disturbed process rate is created in order to give a better representation of the state of the atmosphere. Another method is the Cloud Layers Unified By Binormals (CLUBB) (Golaz et al., 2002a, b), which works with a set of PDFs for all cloud types to avoid difficulties in coupling of stratocumulus and shallow convection parameterizations. Since CLUBB does not parameterize subgrid-scale variability of cloud ice, Thayer-Calder et al. (2015) include cloud ice in the CLUBB PDFs. The PDFs are sampled by the subgrid importance Latin hypercube sampler (SILHS), which connects CLUBB with the microphysics for stratiform and convective clouds.

<footnote n="1">¹Institute for Meteorology, Universität Leipzig, Leipzig, Germany</footnote>

<footnote n="2">²Laboratoire d'Optique Atmosphérique, Université de Lille, Lille, France</footnote>
They found improvements, e.g., of liquid water path (LWP), precipitable water, and shortwave cloud forcing, but also a degradation in precipitation. Including the subgrid-scale effect in the autoconversion and accretion rate in warm clouds reduces the bias significantly and leads to an enhancement of the process rate (Boutle et al., 2014; Lebsock et al., 2013). Since previous studies mainly focus on warm-rain formation processes (e.g., Morrison and Gettelman, 2008; Larson and Griffin, 2013; Boutle et al., 2014; Lebsock et al., 2013), it is also important to concentrate on snow formation effects by taking subgrid-scale effects into account.

Mülmenstädt et al. (2015) emphasized that most global rain is produced via ice-phase processes. Thus, we focus on cloud-ice-related process rates, especially on the precipitation processes initiated via the ice phase. In the ICOsahe- dral Nonhydrostatic general circulation model (ICON-AES) snow is formed by the aggregation process. Aggregation describes the process where cloud ice particles grow to snowflake sizes by sticking together. Therefore, in this study we implement a stochastic aggregation parameterization into the ICON-AES (Giorgetta et al., 2018) by taking subgrid-scale variability of cloud ice into account. Subgrid-scale vari- ability of the total water mixing ratio (sum of cloud liquid water, cloud ice, and water vapor) is already used for deter- mining the cloud cover (Sundqvist et al., 1989). Here, we use this uniform distribution approach to create a distribution of cloud ice within the cloudy part of the grid box. Instead of taking a grid box mean in-cloud ice mixing ratio for the non- linear aggregation parameterization, we feed the process rate with a cloud ice mass randomly chosen from the distribution of cloud ice mass assumed in the cloud scheme.

To evaluate the uniform distribution of cloud ice at a global scale, large-scale observations of cloud ice are necessary. The combined dataset of spaceborne radar from the CloudSat satellite (Stephens et al., 2002) and the lidar from the Cloud– Aerosol Lidar and Infrared Pathfinder Observations satellite (CALIPSO, Winker et al., 2010) allows retrieving a global dataset of cloud ice. Here, we use the lidar–radar (DAR- DAR; Delanoë and Hogan, 2008, 2010) dataset. Comparing observed cloud ice with modeled cloud ice is a challenge. Since the cloud ice from the ICON-AES does not include snow (precipitating cloud ice) and convective cloud ice, it is necessary to remove falling and convective cloud ice from the DARDAR dataset in order to give a meaningful compar- ison between model and observations. A flag method from Li et al. (2008, 2012) is used to estimate the precipitating and convective cloud ice part per grid box, which is removed from the dataset.

In this study, we investigate an important process rate which transforms ice to snow, the aggregation rate, and how it is treated in the ICON-AES. We include the stochastic ap- proach in the aggregation in order to quantify the influence of taking subgrid-scale variability into account. The selected distribution function of cloud ice is evaluated with the DAR- DAR dataset and the effect of the stochastic approach in the

ICON-AES is investigated. As an additional evaluation, we compare an unbiased process rate with the stochastic ap- proach in order to investigate how well this simple method performs.

## 2 Methods

For all simulations the ICOsahe- dral Nonhydrostatic general circulation model (ICON-AES-1.3.00, Giorgetta et al., 2018) in its global version is used. It includes the Max Planck Institute physics package based on the ECHAM6 physics (Stevens et al., 2013). All runs were performed for 6 years with prescribed sea surface temperature and sea ice bound- ary conditions for a period from 2004 to 2009 and with in- stantaneous diagnostics output every 6 h by using a model time step of 10 min, a horizontal resolution of 160 km, and 47 vertical hybrid sigma levels up to 80 km height (Crueger et al., 2018; Giorgetta et al., 2018). To avoid any effect of the model spin-up, we ignored the first year from the model results. Therefore all multiyear averages were done for the time period 2005–2009. Afterwards, the ICON-AES data are interpolated to selected pressure levels. The ICON-AES con- tains prognostic equations for water vapor and for cloud liq- uid water and cloud ice in stratiform clouds. Stratiform cloud cover is computed by using a diagnostic cloud cover scheme (Sundqvist et al., 1989). Stratiform cloud microphysics is pa- rameterized following Lohmann and Roeckner (1996). The prognostic equation for the grid box mean cloud ice mixing ratio ($\bar{q}_i$) including the different process rate terms ($Q$) is written as follows:

$$ \frac{d\bar{q}_i}{dt} = Q_{Ti} + Q_{sed} + Q_{dep} - Q_{mli} - Q_{sbi} + Q_{fr} - Q_{saci} - Q_{agg}, \quad (1) $$

including advection, parameterized turbulent diffusion, and convective detrainment of cloud ice $Q_{Ti}$; sedimentation of cloud ice $Q_{sed}$; deposition and sublimation $Q_{dep}$; melting of cloud ice $Q_{mli}$; instantaneous sublimation in the cloud-free part $Q_{sbi}$; homogeneous and heterogeneous freezing of liq- uid water $Q_{fr}$; accretion of cloud ice by snow $Q_{saci}$; and the aggregation of cloud ice $Q_{agg}$ (Giorgetta et al., 2013). Mi- crophysical processes are determined from in-cloud cloud ice and liquid water, respectively. The in-cloud values are calculated by dividing the grid box mean cloud ice mixing ratio and cloud liquid water mixing ratio by the fractional cloud cover ($C$). $C$ is calculated by a diagnostic cloud cover scheme by Sundqvist et al. (1989) (more details in Sect. 2.2). Cloud liquid water and cloud ice are considered to be well- mixed if they coexist. Therefore, there are no separate liq- uid or ice parts in the cloud. The ECHAM6 physics includes diagnostic rain and snow profiles in the columns. It is not transported by advection (Giorgetta et al., 2018).
They found improvements, e.g., of liquid water path (LWP), precipitable water, and shortwave cloud forcing, but also a degradation in precipitation. Including the subgrid-scale effect in the autoconversion and accretion rate in warm clouds reduces the bias significantly and leads to an enhancement of the process rate (Boutle et al., 2014; Lebsock et al., 2013). Since previous studies mainly focus on warm-rain formation processes (e.g., Morrison and Gettelman, 2008; Larson and Griffin, 2013; Boutle et al., 2014; Lebsock et al., 2013), it is also important to concentrate on snow formation effects by taking subgrid-scale effects into account.

Mülmenstädt et al. (2015) emphasized that most global rain is produced via ice-phase processes. Thus, we focus on cloud-ice-related process rates, especially on the precipitation processes initiated via the ice phase. In the ICOsahe- dral Nonhydrostatic general circulation model (ICON-AES) snow is formed by the aggregation process. Aggregation describes the process where cloud ice particles grow to snowflake sizes by sticking together. Therefore, in this study we implement a stochastic aggregation parameterization into the ICON-AES (Giorgetta et al., 2018) by taking subgrid-scale variability of cloud ice into account. Subgrid-scale vari- ability of the total water mixing ratio (sum of cloud liquid water, cloud ice, and water vapor) is already used for deter- mining the cloud cover (Sundqvist et al., 1989). Here, we use this uniform distribution approach to create a distribution of cloud ice within the cloudy part of the grid box. Instead of taking a grid box mean in-cloud ice mixing ratio for the non- linear aggregation parameterization, we feed the process rate with a cloud ice mass randomly chosen from the distribution of cloud ice mass assumed in the cloud scheme.

To evaluate the uniform distribution of cloud ice at a global scale, large-scale observations of cloud ice are necessary. The combined dataset of spaceborne radar from the CloudSat satellite (Stephens et al., 2002) and the lidar from the Cloud– Aerosol Lidar and Infrared Pathfinder Observations satellite (CALIPSO, Winker et al., 2010) allows retrieving a global dataset of cloud ice. Here, we use the lidar–radar (DAR- DAR; Delanoë and Hogan, 2008, 2010) dataset. Comparing observed cloud ice with modeled cloud ice is a challenge. Since the cloud ice from the ICON-AES does not include snow (precipitating cloud ice) and convective cloud ice, it is necessary to remove falling and convective cloud ice from the DARDAR dataset in order to give a meaningful compar- ison between model and observations. A flag method from Li et al. (2008, 2012) is used to estimate the precipitating and convective cloud ice part per grid box, which is removed from the dataset.

In this study, we investigate an important process rate which transforms ice to snow, the aggregation rate, and how it is treated in the ICON-AES. We include the stochastic ap- proach in the aggregation in order to quantify the influence of taking subgrid-scale variability into account. The selected distribution function of cloud ice is evaluated with the DAR- DAR dataset and the effect of the stochastic approach in the

ICON-AES is investigated. As an additional evaluation, we compare an unbiased process rate with the stochastic ap- proach in order to investigate how well this simple method performs.

## 2 Methods

For all simulations the ICOsahe- dral Nonhydrostatic general circulation model (ICON-AES-1.3.00, Giorgetta et al., 2018) in its global version is used. It includes the Max Planck Institute physics package based on the ECHAM6 physics (Stevens et al., 2013). All runs were performed for 6 years with prescribed sea surface temperature and sea ice bound- ary conditions for a period from 2004 to 2009 and with in- stantaneous diagnostics output every 6 h by using a model time step of 10 min, a horizontal resolution of 160 km, and 47 vertical hybrid sigma levels up to 80 km height (Crueger et al., 2018; Giorgetta et al., 2018). To avoid any effect of the model spin-up, we ignored the first year from the model results. Therefore all multiyear averages were done for the time period 2005–2009. Afterwards, the ICON-AES data are interpolated to selected pressure levels. The ICON-AES con- tains prognostic equations for water vapor and for cloud liq- uid water and cloud ice in stratiform clouds. Stratiform cloud cover is computed by using a diagnostic cloud cover scheme (Sundqvist et al., 1989). Stratiform cloud microphysics is pa- rameterized following Lohmann and Roeckner (1996). The prognostic equation for the grid box mean cloud ice mixing ratio ($\bar{q}_i$) including the different process rate terms ($Q$) is written as follows:

$$ \frac{d\bar{q}_i}{dt} = Q_{Ti} + Q_{sed} + Q_{dep} - Q_{mli} - Q_{sbi} + Q_{fr} - Q_{saci} - Q_{agg}, \quad (1) $$

including advection, parameterized turbulent diffusion, and convective detrainment of cloud ice $Q_{Ti}$; sedimentation of cloud ice $Q_{sed}$; deposition and sublimation $Q_{dep}$; melting of cloud ice $Q_{mli}$; instantaneous sublimation in the cloud-free part $Q_{sbi}$; homogeneous and heterogeneous freezing of liq- uid water $Q_{fr}$; accretion of cloud ice by snow $Q_{saci}$; and the aggregation of cloud ice $Q_{agg}$ (Giorgetta et al., 2013). Mi- crophysical processes are determined from in-cloud cloud ice and liquid water, respectively. The in-cloud values are calculated by dividing the grid box mean cloud ice mixing ratio and cloud liquid water mixing ratio by the fractional cloud cover ($C$). $C$ is calculated by a diagnostic cloud cover scheme by Sundqvist et al. (1989) (more details in Sect. 2.2). Cloud liquid water and cloud ice are considered to be well- mixed if they coexist. Therefore, there are no separate liq- uid or ice parts in the cloud. The ECHAM6 physics includes diagnostic rain and snow profiles in the columns. It is not transported by advection (Giorgetta et al., 2018).
## 2.1 Aggregation parameterization in the ICON-AES

In this study, we focus on the aggregation parameterization. In the current version of the ICON-AES, with ECHAM6 physics, the conversion rate of ice to snow by aggregation as a nonlinear process is given by Levkov et al. (1992) based on the work of Murakami (1990), with

$$Q_{\text{agg}} = \frac{q_i}{\Delta t_1}. \qquad (2)$$

$\Delta t_1$ is defined as the time needed for ice crystals to grow from the mean radius $\bar{R}_{vi}$ to the smallest radius of snow class particles $R_{S0}$:

$$\Delta t_1 = -\frac{2}{c_1} \log \left( \frac{R_{vi}}{R_{S0}} \right)^3, \qquad (3)$$

with

$$c_1 = \frac{q_i \rho a_I E_{ii} X}{\rho_i} \left(\frac{\rho_0}{\rho}\right)^{\frac{1}{3}}, \qquad (4)$$

where $a_I = 700 \text{ s}^{-1}$ is an empirical constant, $E_{ii} = 0.1$ the collection efficiency, $\rho_i = 500 \text{ kg m}^{-3}$ the density of cloud ice, and $\rho_0 = 1.3 \text{ kg m}^{-3}$ the reference density of air. Combining Eqs. (2)–(4) leads to the final aggregation rate equation:

$$Q_{\text{agg}} = C \gamma \frac{\rho \bar{q}_i^2 a_I E_{ii} X \left(\frac{\rho_0}{\rho}\right)^{\frac{1}{3}}}{-2 \rho_i \log \left(\frac{R_{vi}}{R_{S0}}\right)^3}. \qquad (5)$$

The process may be tuned with a tuning parameter $\gamma$, which is currently set to $\gamma = 95$. The parameterization uses the grid box mean in-cloud ice mixing ratio to calculate the process rate of aggregation from ice to snow, which introduces biases in the process rates, since the aggregation is a nonlinear process (Pincus and Klein, 2000). An unbiased aggregation rate is calculated by using an integral over a distribution function of subgrid-scale cloud ice mixing ratio ($\overline{Q}_{\text{agg}}(q_i)$).

## 2.2 Subgrid-scale variability of cloud ice in the aggregation process

The ICON-AES determines the fractional cloud cover according to Sundqvist et al. (1989). A uniform distribution function of the total water mixing ratio $q_t$ from $\bar{q}_t - \Delta q_t$ to $\bar{q}_t + \Delta q_t$ is considered (Fig. 1). The total water mixing ratio describes the sum of water vapor, cloud liquid water, and cloud ice mixing ratio. The saturation specific humidity $q_s$ is calculated from the grid box mean temperature considering the saturation with respect to ice at temperatures below 0°C and if $q_i$ is higher than a threshold value $\gamma_{thr} = 5 \times 10^{-7} \text{ kg kg}^{-1}$. The integral over the distribution from the saturation humidity $q_s$ up to the maximum of the distribution function $\bar{q}_t + \Delta q_t$ defines the fractional cloud

Cloud cover

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Schematic overview of cloud ice variability as part of the cloud cover parameterization: the cloud cover is defined as the part of the PDF which exceeds the saturation humidity $q_s$ up to $\bar{q}_t + \Delta q_t$. Within the cloudy part we define a distribution over the hydrometeors $q_i + q_l$. The half-width of the PDF of $\bar{q}_i$ is defined as $C$ multiplied with $\Delta q_t$ and the cloud ice fraction ($f_{\text{ice}} = \frac{\bar{q}_i}{\bar{q}_i + q_l}$).</figcaption></figure>

cover $C$, which depends on the calculation of the distribution width $\Delta q_t = \gamma q_s$. The parameter $\gamma$ varies with height from low values near the surface to larger ones in the free troposphere with a prescribed profile (Quaas, 2012). The final equation for $C$, which is used in the model, is given by

$$C = 1 - \sqrt{1 - \frac{r - r_0}{r_{\text{sat}} - r_0}}, \qquad (6)$$

where $r$ is the relative humidity, $r_{sat}$ is the saturation value ($r_{sat} = 1$), and $r_0$ is a function of pressure and depends on two different tuning parameters ($r_{top} = 0.8$ and $r_{surf} = 0.968$), which defines the condensation threshold.

To define a new subgrid-scale cloud ice mass, we use this fractional cloud cover approach (Fig. 1). The half-width of the cloud ice PDF ($\Delta q_i$) has to be re-scaled with $\Delta q_i = C \Delta q_t f_{\text{ice}}$, where $f_{\text{ice}}$ describes the cloud ice fraction ($\frac{\bar{q}_i}{\bar{q}_i + \bar{q}_l}$). Using the PDF over $q_t$ the all-sky cloud ice and cloud liquid water mixing ratio is defined as the integral over $q_t - q_s$ from $q_s$ to the maximum of the PDF($q_t$):

$$\bar{q}_i + \bar{q}_l = \int_{q_s}^{\bar{q}_t + \Delta q_t} (q_t - q_s) \text{PDF}(q_t) dq_t. \qquad (7)$$

Here, $q_i$ and $q_l$ are considered to be well-mixed within the cloudy part of the grid box, so $q_i$ and $q_l$ occur in the same volume. Solving the Eq. (7) yields

$$\bar{q}_i + \bar{q}_l = C^2 \Delta q_t. \qquad (8)$$

To get $\bar{q}_i$, Eq. (8) is multiplied with $f_{\text{ice}}$, which leads to

$$\bar{q}_i = (\bar{q}_i + \bar{q}_l) \cdot f_{\text{ice}} = C^2 \Delta q_t \cdot f_{\text{ice}}. \qquad (9)$$

As described above, the in-cloud cloud ice mixing ratio $q_i^c$ is defined as $q_i$ divided by $C$, and in combination with Eq. (8)
## 2.1 Aggregation parameterization in the ICON-AES

In this study, we focus on the aggregation parameterization. In the current version of the ICON-AES, with ECHAM6 physics, the conversion rate of ice to snow by aggregation as a nonlinear process is given by Levkov et al. (1992) based on the work of Murakami (1990), with

$$Q_{\text{agg}} = \frac{q_i}{\Delta t_1}. \qquad (2)$$

$\Delta t_1$ is defined as the time needed for ice crystals to grow from the mean radius $\bar{R}_{vi}$ to the smallest radius of snow class particles $R_{S0}$:

$$\Delta t_1 = -\frac{2}{c_1} \log \left( \frac{R_{vi}}{R_{S0}} \right)^3, \qquad (3)$$

with

$$c_1 = \frac{q_i \rho a_I E_{ii} X}{\rho_i} \left(\frac{\rho_0}{\rho}\right)^{\frac{1}{3}}, \qquad (4)$$

where $a_I = 700 \text{ s}^{-1}$ is an empirical constant, $E_{ii} = 0.1$ the collection efficiency, $\rho_i = 500 \text{ kg m}^{-3}$ the density of cloud ice, and $\rho_0 = 1.3 \text{ kg m}^{-3}$ the reference density of air. Combining Eqs. (2)–(4) leads to the final aggregation rate equation:

$$Q_{\text{agg}} = C \gamma \frac{\rho \bar{q}_i^2 a_I E_{ii} X \left(\frac{\rho_0}{\rho}\right)^{\frac{1}{3}}}{-2 \rho_i \log \left(\frac{R_{vi}}{R_{S0}}\right)^3}. \qquad (5)$$

The process may be tuned with a tuning parameter $\gamma$, which is currently set to $\gamma = 95$. The parameterization uses the grid box mean in-cloud ice mixing ratio to calculate the process rate of aggregation from ice to snow, which introduces biases in the process rates, since the aggregation is a nonlinear process (Pincus and Klein, 2000). An unbiased aggregation rate is calculated by using an integral over a distribution function of subgrid-scale cloud ice mixing ratio ($\overline{Q}_{\text{agg}}(q_i)$).

## 2.2 Subgrid-scale variability of cloud ice in the aggregation process

The ICON-AES determines the fractional cloud cover according to Sundqvist et al. (1989). A uniform distribution function of the total water mixing ratio $q_t$ from $\bar{q}_t - \Delta q_t$ to $\bar{q}_t + \Delta q_t$ is considered (Fig. 1). The total water mixing ratio describes the sum of water vapor, cloud liquid water, and cloud ice mixing ratio. The saturation specific humidity $q_s$ is calculated from the grid box mean temperature considering the saturation with respect to ice at temperatures below 0°C and if $q_i$ is higher than a threshold value $\gamma_{thr} = 5 \times 10^{-7} \text{ kg kg}^{-1}$. The integral over the distribution from the saturation humidity $q_s$ up to the maximum of the distribution function $\bar{q}_t + \Delta q_t$ defines the fractional cloud

Cloud cover

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Schematic overview of cloud ice variability as part of the cloud cover parameterization: the cloud cover is defined as the part of the PDF which exceeds the saturation humidity $q_s$ up to $\bar{q}_t + \Delta q_t$. Within the cloudy part we define a distribution over the hydrometeors $q_i + q_l$. The half-width of the PDF of $\bar{q}_i$ is defined as $C$ multiplied with $\Delta q_t$ and the cloud ice fraction ($f_{\text{ice}} = \frac{\bar{q}_i}{\bar{q}_i + q_l}$).</figcaption></figure>

cover $C$, which depends on the calculation of the distribution width $\Delta q_t = \gamma q_s$. The parameter $\gamma$ varies with height from low values near the surface to larger ones in the free troposphere with a prescribed profile (Quaas, 2012). The final equation for $C$, which is used in the model, is given by

$$C = 1 - \sqrt{1 - \frac{r - r_0}{r_{\text{sat}} - r_0}}, \qquad (6)$$

where $r$ is the relative humidity, $r_{sat}$ is the saturation value ($r_{sat} = 1$), and $r_0$ is a function of pressure and depends on two different tuning parameters ($r_{top} = 0.8$ and $r_{surf} = 0.968$), which defines the condensation threshold.

To define a new subgrid-scale cloud ice mass, we use this fractional cloud cover approach (Fig. 1). The half-width of the cloud ice PDF ($\Delta q_i$) has to be re-scaled with $\Delta q_i = C \Delta q_t f_{\text{ice}}$, where $f_{\text{ice}}$ describes the cloud ice fraction ($\frac{\bar{q}_i}{\bar{q}_i + \bar{q}_l}$). Using the PDF over $q_t$ the all-sky cloud ice and cloud liquid water mixing ratio is defined as the integral over $q_t - q_s$ from $q_s$ to the maximum of the PDF($q_t$):

$$\bar{q}_i + \bar{q}_l = \int_{q_s}^{\bar{q}_t + \Delta q_t} (q_t - q_s) \text{PDF}(q_t) dq_t. \qquad (7)$$

Here, $q_i$ and $q_l$ are considered to be well-mixed within the cloudy part of the grid box, so $q_i$ and $q_l$ occur in the same volume. Solving the Eq. (7) yields

$$\bar{q}_i + \bar{q}_l = C^2 \Delta q_t. \qquad (8)$$

To get $\bar{q}_i$, Eq. (8) is multiplied with $f_{\text{ice}}$, which leads to

$$\bar{q}_i = (\bar{q}_i + \bar{q}_l) \cdot f_{\text{ice}} = C^2 \Delta q_t \cdot f_{\text{ice}}. \qquad (9)$$

As described above, the in-cloud cloud ice mixing ratio $q_i^c$ is defined as $q_i$ divided by $C$, and in combination with Eq. (8)
it follows that

$$ \bar{q}_i^c = \frac{\bar{q}_i}{C} = C \Delta q_t \cdot f_{\text{ice}} = \Delta q_i. \qquad (10) $$

This implies that the width of the cloud ice distribution corresponds to the in-cloud grid box mean cloud ice mixing ratio $\bar{q}_i^c$. To obtain the subgrid-scale cloud ice mixing ratio, we use a Monte Carlo approach. Here, we choose an element from the cumulative distribution function (CDF) with the help of a random number $r \in [0, 1]$. This yields the following equation:

$$ q_{i, \text{new}}^c = \bar{q}_i^c + \Delta q_i(2r - 1). \qquad (11) $$

With Eq. (10),

$$ q_{i, \text{new}}^c = \bar{q}_i^c(2r). \qquad (12) $$

Finally, the subgrid-scale cloud ice mixing ratio only depends on the grid box mean cloud ice mass mixing ratio and the choice of the random number. This new specific cloud ice mass replaces the grid box mean cloud ice in Eq. (5). To compare the current, biased aggregation rate ($Q_{\text{agg}}(\bar{q}_i)$) with the unbiased process rate ($\bar{Q}_{\text{agg}}(q_i)$) we replace for each time step and grid box $\bar{q}_i$ mean in Eq. (5) with the integral over the entire distribution of $q_i$. This comparison will be shown in the last part of the Results section.

To evaluate the distribution function, the cloud ice variance is calculated for all-sky conditions because all output cloud variables are also all-sky variables. The variance can be written as follows:

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} (q_i - \bar{q}_i)^2 \text{PDF}(q_i) dq_i. \qquad (13) $$

As we focus on all-sky variance, this equation can be separated into two parts: the clear-sky part of the distribution ($q_i = 0$) and the cloudy part ($q_i > 0$). The cloud-free part can be described with Dirac's delta function, which yields the following equation (derivation is shown in the Appendix):

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} \left[ (1-C) \delta(q_i=0) + C \frac{1}{2\Delta q_i} \right] (q_i - \bar{q}_i)^2 dq_i \\ \sigma_{q_i}^2 = (q_i^c)^2 \cdot \left( \frac{4}{3} C - C^2 \right). \qquad (14) $$

## 2.3 Satellite retrievals of cloud ice water content

As mentioned above, we will focus on ice-phase processes and their impact on the global ice water content. To evaluate the results, a combined global ice water product of the CloudSat (Stephens et al., 2002) and Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) (Winker et al., 2010), the DARDAR-Nice dataset (Sourdeval

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Zonally averaged annual mean ice water content (mg kg⁻¹) from the DARDAR dataset. (a) Total ice mixing ratio ($q_{i,\text{total}}$), which includes cloud ice from any clouds, and precipitating ice. (b) Cloud ice mixing ratio ($q_i$), where precipitating ice and convective cloud ice are removed. (c) Difference between $q_{i,\text{total}}$ and $q_i$.</figcaption></figure>

et al., 2018), is used. Both satellites are part of the A-train constellation and are flying with a time interval of only 15 s between them. The W-band (94 GHz) cloud-profiling radar (CPR) on CloudSat provides vertical cloud profiles with a minimum detectable reflectivity of −28 dBz. CALIPSO contains the Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP), which measures co-aligned pulses of 532 and 1064 nm wavelength. Both satellite instruments provide global retrievals of clouds and their properties. A combined dataset of lidar and radar measurements provides global information about larger as well as smaller ice particles and therefore a more detailed overview of the cloud. Outside the temperature range of −40 to 0 °C, supercooled liquid water does not exist. The time period used from the dataset is January 2007 until December 2010, since the datasets are available for the entire years.

For an accurate comparison between satellite data and model data the two-dimensional DARDAR data have to be interpolated to the same horizontal grid as the ICON-AES. For each grid box, the mean of all variables is calculated and vertically interpolated to the same pressure levels as for the ICON-AES. The modeled $q_i$, which is used for the aggreg-
it follows that

$$ \bar{q}_i^c = \frac{\bar{q}_i}{C} = C \Delta q_t \cdot f_{\text{ice}} = \Delta q_i. \qquad (10) $$

This implies that the width of the cloud ice distribution corresponds to the in-cloud grid box mean cloud ice mixing ratio $\bar{q}_i^c$. To obtain the subgrid-scale cloud ice mixing ratio, we use a Monte Carlo approach. Here, we choose an element from the cumulative distribution function (CDF) with the help of a random number $r \in [0, 1]$. This yields the following equation:

$$ q_{i,\text{new}}^c = \bar{q}_i^c + \Delta q_i(2r - 1). \qquad (11) $$

With Eq. (10),

$$ q_{i,\text{new}}^c = \bar{q}_i^c(2r). \qquad (12) $$

Finally, the subgrid-scale cloud ice mixing ratio only depends on the grid box mean cloud ice mass mixing ratio and the choice of the random number. This new specific cloud ice mass replaces the grid box mean cloud ice in Eq. (5). To compare the current, biased aggregation rate ($Q_{\text{agg}}(\bar{q}_i)$) with the unbiased process rate ($\bar{Q}_{\text{agg}}(q_i)$) we replace for each time step and grid box $\bar{q}_i$ mean in Eq. (5) with the integral over the entire distribution of $q_i$. This comparison will be shown in the last part of the Results section.

To evaluate the distribution function, the cloud ice variance is calculated for all-sky conditions because all output cloud variables are also all-sky variables. The variance can be written as follows:

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} (q_i - \bar{q}_i)^2 \text{PDF}(q_i) dq_i. \qquad (13) $$

As we focus on all-sky variance, this equation can be separated into two parts: the clear-sky part of the distribution ($q_i = 0$) and the cloudy part ($q_i > 0$). The cloud-free part can be described with Dirac's delta function, which yields the following equation (derivation is shown in the Appendix):

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} \left[ (1-C) \delta(q_i=0) + C \frac{1}{2\Delta q_i} \right] (q_i - \bar{q}_i)^2 dq_i \\ \sigma_{q_i}^2 = (q_i^c)^2 \cdot \left( \frac{4}{3}C - C^2 \right). \qquad (14) $$

## 2.3 Satellite retrievals of cloud ice water content

As mentioned above, we will focus on ice-phase processes and their impact on the global ice water content. To evaluate the results, a combined global ice water product of the CloudSat (Stephens et al., 2002) and Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observations (CALIPSO) (Winker et al., 2010), the DARDAR-Nice dataset (Sourdeval

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Zonally averaged annual mean ice water content (mg kg⁻¹) from the DARDAR dataset. (a) Total ice mixing ratio ($q_{i,\text{total}}$), which includes cloud ice from any clouds, and precipitating ice. (b) Cloud ice mixing ratio ($q_i$), where precipitating ice and convective cloud ice are removed. (c) Difference between $q_{i,\text{total}}$ and $q_i$.</figcaption></figure>

et al., 2018), is used. Both satellites are part of the A-train constellation and are flying with a time interval of only 15 s between them. The W-band (94 GHz) cloud-profiling radar (CPR) on CloudSat provides vertical cloud profiles with a minimum detectable reflectivity of −28 dBz. CALIPSO contains the Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP), which measures co-aligned pulses of 532 and 1064 nm wavelength. Both satellite instruments provide global retrievals of clouds and their properties. A combined dataset of lidar and radar measurements provides global information about larger as well as smaller ice particles and therefore a more detailed overview of the cloud. Outside the temperature range of −40 to 0 °C, supercooled liquid water does not exist. The time period used from the dataset is January 2007 until December 2010, since the datasets are available for the entire years.

For an accurate comparison between satellite data and model data the two-dimensional DARDAR data have to be interpolated to the same horizontal grid as the ICON-AES. For each grid box, the mean of all variables is calculated and vertically interpolated to the same pressure levels as for the ICON-AES. The modeled $q_i$, which is used for the aggreg-
gation parameterization, does not include precipitating and convective cloud ice. Therefore, an estimate of the convective and precipitating cloud ice has to be removed from the satellite dataset. To provide a dataset that is comparable with the ICON-AES, a flag method after Li et al. (2008, 2012) is used. All columns that are flagged as "precipitating at the surface" are not considered to ensure that larger falling ice particles are removed from the dataset. This is because the model distinguishes between precipitating ice and cloud ice in two bulk classes. In contrast, it is challenging for the satellite retrievals to separate falling ice and cloud ice. Additionally, cloud ice, which is flagged as originating from "deep convection" or "cumulus" (2B-CLDCLASS dataset), is also removed. In order to compare modeled cloud ice with observations, alternative methods are possible (e.g., removing ICON nonzero surface precipitation points). Since we want to evaluate the cloud ice distribution that is used for the aggregation, we had to adjust the DARDAR data in order to find the most consistent way.

Cloud variables simulated by the ICON-AES are output for grid box mean all-sky conditions. Therefore, grid box mean cloud ice and cloud ice variance are calculated for all-sky conditions in the observational data as well.

Figure 2 shows a multiyear mean of zonally averaged total ice mixing ratio ($q_{i,\text{total}}$) and cloud ice mixing ratio ($q_i$) obtained from the DARDAR product. $q_{i,\text{total}}$ shows the highest values over the Equator at lower pressure fields and over the midlatitudes at higher pressure fields. Most of the removed cloud ice is detected due to the precipitation flag. As a result, the remaining cloud ice water content (Fig. 2b) is much lower, especially in the midlatitudes at lower regions and in the tropics. The maxima of cloud ice are shifted upward, since most of the precipitating cloud ice is below or in the lower regions of the cloud. The plot highlights the importance of removing precipitation from the DARDAR dataset in order to give a more reliable comparison with the ICON-AES. This is due to the fact that the remaining $q_i$ is much less than $q_{i,\text{total}}$. But we have to keep in mind that the result of $q_i$ strongly depends on the accuracy of the different flags because it selects which cloud ice is detected as convective or precipitating cloud ice and therefore which cloud ice is removed from the dataset. Additionally, the initial satellite data are measured on a 1D curtain, while the model uses 3D grid boxes. Hill et al. (2015) calculated a measure of the difference in standard deviation considering a 2D cloud field compared to a 1D cross-section. They estimated a 30% larger standard deviation in 2D fields compared to the 1D track. Therefore, we should consider that this also has an effect on our cloud ice variance calculation. However, there are limitations in availability of satellite data. Therefore we tried to be as consistent as possible in the comparison between simulations and observations.

# 3 Results

## 3.1 Employing a distribution of cloud ice in the aggregation parameterization of the ICON-AES

As described above, for the evaluation of the cloud ice distribution the all-sky cloud ice variance is calculated for each grid box. Equation (14) was used for the modeled cloud ice variance, while for the DARDAR data the spatial variance within each GCM grid box was calculated from all footprints within a grid box. But before comparing the cloud ice variance, we have to make sure that the cloud ice mixing ratios from DARDAR data and ICON-AES show the same order of magnitude. Figure 3 shows the annually averaged global distribution of cloud ice from DARDAR and ICON as well as the differences between them at different pressure levels. ICON-AES shows higher values in the midlatitudes at higher pressure fields. At 400 hPa the main difference is around the Equator, where DARDAR shows higher values than ICON-AES. ICON-AES still has maximum values in the midlatitudes. However, there is good agreement between the pattern of cloud ice mixing ratio from ICON and DARDAR.

The cloud ice water path (CIWP), which is the column-integrated ice water content, is given in Fig. 4. ICON overestimates the CIWP in the middle and higher latitudes, while it underestimates the CIWP in the tropics. As is already visible in Fig. 3, DARDAR shows larger cloud ice values down the lower altitudes over the tropics, which leads to an increased CIWP compared to the modeled CIWP. Especially in the midlatitudes at higher pressure fields, the model tends to overestimate the cloud ice. One should keep in mind that the way the DARDAR data are filtered to get the CIWP or $q_i$, which is comparable with the model data, is not perfect.

Figure 5 shows the comparison of cloud ice variance between ICON-AES and DARDAR data for the same pressure levels chosen in Fig. 3. In general, the pattern of cloud ice variance follows the pattern of the global cloud ice distribution. Higher values are visible over the tropics at higher altitude and over the storm tracks in the midlatitudes at lower altitudes. Besides the good agreement in the distribution pattern, there are some regional differences. At 300 hPa, the maxima over the northern part of North America and central Africa are more intense in the DARDAR data than in the ICON data. At 400 hPa ICON underestimates the cloud ice variance compared to DARDAR. At lower altitudes the cloud ice variance is higher in the ICON data in most of the midlatitudes compared to the DARDAR data. ICON underestimates the cloud ice variance in the tropics, especially in lower regions, since ICON shows less cloud ice in the tropical region at the same pressure levels. However, we should keep in mind that the measured variance contains uncertainties regarding the method of filtering precipitation and convection. Additionally the modeled cloud ice variance makes an assumption of a distribution, while the DARDAR data show the variance from discrete measurements. Despite these dis-
gation parameterization, does not include precipitating and convective cloud ice. Therefore, an estimate of the convective and precipitating cloud ice has to be removed from the satellite dataset. To provide a dataset that is comparable with the ICON-AES, a flag method after Li et al. (2008, 2012) is used. All columns that are flagged as "precipitating at the surface" are not considered to ensure that larger falling ice particles are removed from the dataset. This is because the model distinguishes between precipitating ice and cloud ice in two bulk classes. In contrast, it is challenging for the satellite retrievals to separate falling ice and cloud ice. Additionally, cloud ice, which is flagged as originating from "deep convection" or "cumulus" (2B-CLDCLASS dataset), is also removed. In order to compare modeled cloud ice with observations, alternative methods are possible (e.g., removing ICON nonzero surface precipitation points). Since we want to evaluate the cloud ice distribution that is used for the aggregation, we had to adjust the DARDAR data in order to find the most consistent way.

Cloud variables simulated by the ICON-AES are output for grid box mean all-sky conditions. Therefore, grid box mean cloud ice and cloud ice variance are calculated for all-sky conditions in the observational data as well.

Figure 2 shows a multiyear mean of zonally averaged total ice mixing ratio ($q_{i,\text{total}}$) and cloud ice mixing ratio ($q_i$) obtained from the DARDAR product. $q_{i,\text{total}}$ shows the highest values over the Equator at lower pressure fields and over the midlatitudes at higher pressure fields. Most of the removed cloud ice is detected due to the precipitation flag. As a result, the remaining cloud ice water content (Fig. 2b) is much lower, especially in the midlatitudes at lower regions and in the tropics. The maxima of cloud ice are shifted upward, since most of the precipitating cloud ice is below or in the lower regions of the cloud. The plot highlights the importance of removing precipitation from the DARDAR dataset in order to give a more reliable comparison with the ICON-AES. This is due to the fact that the remaining $q_i$ is much less than $q_{i,\text{total}}$. But we have to keep in mind that the result of $q_i$ strongly depends on the accuracy of the different flags because it selects which cloud ice is detected as convective or precipitating cloud ice and therefore which cloud ice is removed from the dataset. Additionally, the initial satellite data are measured on a 1D curtain, while the model uses 3D grid boxes. Hill et al. (2015) calculated a measure of the difference in standard deviation considering a 2D cloud field compared to a 1D cross-section. They estimated a 30% larger standard deviation in 2D fields compared to the 1D track. Therefore, we should consider that this also has an effect on our cloud ice variance calculation. However, there are limitations in availability of satellite data. Therefore we tried to be as consistent as possible in the comparison between simulations and observations.

# 3 Results

## 3.1 Employing a distribution of cloud ice in the aggregation parameterization of the ICON-AES

As described above, for the evaluation of the cloud ice distribution the all-sky cloud ice variance is calculated for each grid box. Equation (14) was used for the modeled cloud ice variance, while for the DARDAR data the spatial variance within each GCM grid box was calculated from all footprints within a grid box. But before comparing the cloud ice variance, we have to make sure that the cloud ice mixing ratios from DARDAR data and ICON-AES show the same order of magnitude. Figure 3 shows the annually averaged global distribution of cloud ice from DARDAR and ICON as well as the differences between them at different pressure levels. ICON-AES shows higher values in the midlatitudes at higher pressure fields. At 400 hPa the main difference is around the Equator, where DARDAR shows higher values than ICON-AES. ICON-AES still has maximum values in the midlatitudes. However, there is good agreement between the pattern of cloud ice mixing ratio from ICON and DARDAR.

The cloud ice water path (CIWP), which is the column-integrated ice water content, is given in Fig. 4. ICON overestimates the CIWP in the middle and higher latitudes, while it underestimates the CIWP in the tropics. As is already visible in Fig. 3, DARDAR shows larger cloud ice values down the lower altitudes over the tropics, which leads to an increased CIWP compared to the modeled CIWP. Especially in the midlatitudes at higher pressure fields, the model tends to overestimate the cloud ice. One should keep in mind that the way the DARDAR data are filtered to get the CIWP or $q_i$, which is comparable with the model data, is not perfect.

Figure 5 shows the comparison of cloud ice variance between ICON-AES and DARDAR data for the same pressure levels chosen in Fig. 3. In general, the pattern of cloud ice variance follows the pattern of the global cloud ice distribution. Higher values are visible over the tropics at higher altitude and over the storm tracks in the midlatitudes at lower altitudes. Besides the good agreement in the distribution pattern, there are some regional differences. At 300 hPa, the maxima over the northern part of North America and central Africa are more intense in the DARDAR data than in the ICON data. At 400 hPa ICON underestimates the cloud ice variance compared to DARDAR. At lower altitudes the cloud ice variance is higher in the ICON data in most of the midlatitudes compared to the DARDAR data. ICON underestimates the cloud ice variance in the tropics, especially in lower regions, since ICON shows less cloud ice in the tropical region at the same pressure levels. However, we should keep in mind that the measured variance contains uncertainties regarding the method of filtering precipitation and convection. Additionally the modeled cloud ice variance makes an assumption of a distribution, while the DARDAR data show the variance from discrete measurements. Despite these dis-
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Multiyear mean of the cloud ice mixing ratio (kg kg⁻¹) at four different pressure levels calculated for the DARDAR data (a, d, g, j), the ICON data (b, e, h, k), and the difference between ICON and DARDAR (c, f, i, j).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Multiyear cloud ice water path (CIWP) for the DARDAR data and ICON data as well as the difference between DARDAR and ICON.</figcaption></figure>

crepancies, we conclude that the simple uniform distribution of cloud ice as written in Eq. (14) basically captures the measured distribution of cloud ice in each grid box. Therefore it is suitable for use in the aggregation parameterization.

As described in Sect. 2.2, with the help of this cloud ice distribution a new specific cloud ice mass is used for the aggregation parameterization. Figure 6 shows its influence on the zonally averaged cloud ice adapted from Trömel et al. (2021), where the same tuning parameters were used. The control run of the ICON model (CTRL) gives a maximum of cloud ice at levels of lower pressure over the tropics. There are also two maxima over the midlatitudes between 600 and 400 hPa, with a more pronounced maximum in the Southern Hemisphere than in the Northern Hemisphere. An overall re-

duction of cloud ice is visible due to the stochastic aggregation scheme (AGGstoch).

Since aggregation is a nonlinear process, it becomes stronger due to taking a distribution of cloud ice into account. A stronger aggregation process leads to a higher conversion rate from cloud ice to snow. Therefore, more cloud ice is removed due to the process. Figure 7 shows the effect of using the stochastic approach instead of the default parameterization directly for the aggregation rate in a joint histogram and the corresponding multiyear zonally averaged aggregation rate. While the CTRL simulation produces an intense maximum aggregation rate around ca. 30 mg kg⁻¹ h⁻¹, the AG-Gstoch run reveals a larger spread around a less intense maximum. The higher spread of the maximum stems from the
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Multiyear mean of the cloud ice mixing ratio (kg kg⁻¹) at four different pressure levels calculated for the DARDAR data (a, d, g, j), the ICON data (b, e, h, k), and the difference between ICON and DARDAR (c, f, i, j).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Multiyear cloud ice water path (CIWP) for the DARDAR data and ICON data as well as the difference between DARDAR and ICON.</figcaption></figure>

crepancies, we conclude that the simple uniform distribution of cloud ice as written in Eq. (14) basically captures the measured distribution of cloud ice in each grid box. Therefore it is suitable for use in the aggregation parameterization.

As described in Sect. 2.2, with the help of this cloud ice distribution a new specific cloud ice mass is used for the aggregation parameterization. Figure 6 shows its influence on the zonally averaged cloud ice adapted from Trömel et al. (2021), where the same tuning parameters were used. The control run of the ICON model (CTRL) gives a maximum of cloud ice at levels of lower pressure over the tropics. There are also two maxima over the midlatitudes between 600 and 400 hPa, with a more pronounced maximum in the Southern Hemisphere than in the Northern Hemisphere. An overall re-

duction of cloud ice is visible due to the stochastic aggregation scheme (AGGstoch).

Since aggregation is a nonlinear process, it becomes stronger due to taking a distribution of cloud ice into account. A stronger aggregation process leads to a higher conversion rate from cloud ice to snow. Therefore, more cloud ice is removed due to the process. Figure 7 shows the effect of using the stochastic approach instead of the default parameterization directly for the aggregation rate in a joint histogram and the corresponding multiyear zonally averaged aggregation rate. While the CTRL simulation produces an intense maximum aggregation rate around ca. 30 mg kg⁻¹ h⁻¹, the AG-Gstoch run reveals a larger spread around a less intense maximum. The higher spread of the maximum stems from the
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Same as Fig. 3 but for cloud ice variance (kg² kg⁻²).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Zonally averaged multiyear mean cloud ice (mg kg⁻¹). (a) Control run of the ICON-AES. (b) Difference between the new stochastic aggregation parameterization and the control run (adapted from Trömel et al., 2021).</figcaption></figure>

randomly picked cloud ice mass for the aggregation process. Therefore, there are values which are much higher and much lower than the original grid box mean cloud ice mass. Instead of having this intense maximum aggregation rate, higher and lower aggregation rate cases are visible. Therefore, the an-

nual average over all of these aggregation rate cases yields a higher value for the AGGstoch run compared to CTRL run, caused by the nonlinearity of the process, which is also visible in Fig. 7f.

The microphysical processes in clouds are strongly interlinked with each other and react whenever one process rate changes. Figure 8 shows how the different cloud-ice-related process rates change using the stochastic aggregation parameterization. The different vertical profiles of the global mean process rates are depicted. Negative process rates indicate a cloud ice loss, while positive process rates lead to a cloud ice gain. So, a more negative process rate is linked with a stronger cloud ice loss. It is visible that the deposition rate is the most effective cloud-ice-gaining process. It is triggered by the Wegener–Bergeron–Findeisen (WBF) (Wegener, 1911; Bergeron, 1935; Findeisen et al., 2015) process. This process takes place if the temperature is between –35 and 0 °C and the $q_i$ exceeds a specific threshold (Giorgetta et al., 2013). The strongest negative process rate is the accretion rate. The accretion rate (Aci) describes the growth of snow by collecting the surrounding ice crystals. Therefore, the accretion rate is strongly linked with aggregation rate, since the aggregation process produces the snow first. Aggregation and accretion describe the formation and growing of snow by coagulation of cloud ice, and thus they lead to a cloud ice reduction. The riming process (not shown here) is described like the linear accretion process and leads to an increase in snowfall, but with the difference that the equa-
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Same as Fig. 3 but for cloud ice variance (kg² kg⁻²).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Zonally averaged multiyear mean cloud ice (mg kg⁻¹). (a) Control run of the ICON-AES. (b) Difference between the new stochastic aggregation parameterization and the control run (adapted from Trömel et al., 2021).</figcaption></figure>

randomly picked cloud ice mass for the aggregation process. Therefore, there are values which are much higher and much lower than the original grid box mean cloud ice mass. Instead of having this intense maximum aggregation rate, higher and lower aggregation rate cases are visible. Therefore, the an-

nual average over all of these aggregation rate cases yields a higher value for the AGGstoch run compared to CTRL run, caused by the nonlinearity of the process, which is also visible in Fig. 7f.

The microphysical processes in clouds are strongly interlinked with each other and react whenever one process rate changes. Figure 8 shows how the different cloud-ice-related process rates change using the stochastic aggregation parameterization. The different vertical profiles of the global mean process rates are depicted. Negative process rates indicate a cloud ice loss, while positive process rates lead to a cloud ice gain. So, a more negative process rate is linked with a stronger cloud ice loss. It is visible that the deposition rate is the most effective cloud-ice-gaining process. It is triggered by the Wegener–Bergeron–Findeisen (WBF) (Wegener, 1911; Bergeron, 1935; Findeisen et al., 2015) process. This process takes place if the temperature is between -35 and 0 °C and the $q_i$ exceeds a specific threshold (Giorgetta et al., 2013). The strongest negative process rate is the accretion rate. The accretion rate (Aci) describes the growth of snow by collecting the surrounding ice crystals. Therefore, the accretion rate is strongly linked with aggregation rate, since the aggregation process produces the snow first. Aggregation and accretion describe the formation and growing of snow by coagulation of cloud ice, and thus they lead to a cloud ice reduction. The riming process (not shown here) is described like the linear accretion process and leads to an increase in snowfall, but with the difference that the equa-
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Joint histograms of aggregation rate (mg kg⁻¹ h⁻¹) and cloud ice (mg kg⁻¹) for CTRL (a), AGGstoch (b), and the difference between the two (c). Zonally averaged aggregation rate of the multiyear mean for CTRL (d), AGGstoch (e), and the difference (f).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Vertical profiles of the global mean process rates which are related to cloud ice microphysical processes for CTRL (solid lines) and AGGstoch (dashed lines): (Aci), freezing rate (Frz), deposition rate (Dep), melting rate (Mlt), evaporation rate (Evp), sedimentation rate (Sed), aggregation rate (Agg), and accretion rate (transport terms are not included).</figcaption></figure>

tion depends on cloud water. However, the accretion rate is the more efficient process compared to aggregation. Figure 8 confirms the higher cloud ice loss in the AGGstoch compared to the CTRL simulation due to the more negative aggregation rate. Moreover, due to the stronger aggregation rate at levels of high pressure, less cloud ice may be transported into lower regions. Hence, the sedimentation process slightly decreases. Due to the change in the aggregation process rate, the accretion rate becomes less intense and leads to more cloud ice. The reduced accretion rate may be explained by the cloud ice loss due the stronger aggregation

rate. Less cloud ice is collectable for growing snowflakes. No significant changes are visible for freezing, deposition, melting, and evaporation rate. Due to the stochastic method, the global mean aggregation rate is intensified by more than 20% in the middle and upper troposphere. In contrast, the maximum change in global mean accretion rate is less than 8.9%. Since the accretion rate is much higher, this smaller change leads to a compensation of the aggregation rate increase as described above. However, the increase in aggregation rate by more than 20% is significant. Due to the small change in the microphysical properties, no important change in multiyear global mean shortwave (+0.161 W m⁻²), long-wave (−0.195 W m⁻²), and net radiation (+0.03 W m⁻²) at the top of the atmosphere is visible (not shown here). Overall, we found a cloud ice loss in the AGGstoch run of up to 5%, but the reduction of cloud ice is compensated for by a less intense accretion rate. The effect on radiation could be increased if this stochastic approach was implemented in other processes, since there has to be a change of a factor of 2 or more in microphysical process rates to see an effect on radiation (Michibata et al., 2020; Imura and Michibata, 2022). Using the new approach just for one single process makes it easier in the beginning to see the effect of changing one process, since all processes are connected. Besides that, the aggregation process is the only nonlinear cloud ice process rate in the ICON-AES. Since we focused on cloud-ice-related processes we just implemented the subgrid-scale approach in the aggregation parameterization. In future studies one can think about additionally including the subgrid-scale approach in cloud-water-related processes in order to see a
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Joint histograms of aggregation rate (mg kg⁻¹ h⁻¹) and cloud ice (mg kg⁻¹) for CTRL (a), AGGstoch (b), and the difference between the two (c). Zonally averaged aggregation rate of the multiyear mean for CTRL (d), AGGstoch (e), and the difference (f).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Vertical profiles of the global mean process rates which are related to cloud ice microphysical processes for CTRL (solid lines) and AGGstoch (dashed lines): (Aci), freezing rate (Frz), deposition rate (Dep), melting rate (Mlt), evaporation rate (Evp), sedimentation rate (Sed), aggregation rate (Agg), and accretion rate (transport terms are not included).</figcaption></figure>

tion depends on cloud water. However, the accretion rate is the more efficient process compared to aggregation. Figure 8 confirms the higher cloud ice loss in the AGGstoch compared to the CTRL simulation due to the more negative aggregation rate. Moreover, due to the stronger aggregation rate at levels of high pressure, less cloud ice may be transported into lower regions. Hence, the sedimentation process slightly decreases. Due to the change in the aggregation process rate, the accretion rate becomes less intense and leads to more cloud ice. The reduced accretion rate may be explained by the cloud ice loss due the stronger aggregation

rate. Less cloud ice is collectable for growing snowflakes. No significant changes are visible for freezing, deposition, melting, and evaporation rate. Due to the stochastic method, the global mean aggregation rate is intensified by more than 20% in the middle and upper troposphere. In contrast, the maximum change in global mean accretion rate is less than 8.9%. Since the accretion rate is much higher, this smaller change leads to a compensation of the aggregation rate increase as described above. However, the increase in aggregation rate by more than 20% is significant. Due to the small change in the microphysical properties, no important change in multiyear global mean shortwave (+0.161 W m⁻²), longwave (−0.195 W m⁻²), and net radiation (+0.03 W m⁻²) at the top of the atmosphere is visible (not shown here). Overall, we found a cloud ice loss in the AGGstoch run of up to 5%, but the reduction of cloud ice is compensated for by a less intense accretion rate. The effect on radiation could be increased if this stochastic approach was implemented in other processes, since there has to be a change of a factor of 2 or more in microphysical process rates to see an effect on radiation (Michibata et al., 2020; Imura and Michibata, 2022). Using the new approach just for one single process makes it easier in the beginning to see the effect of changing one process, since all processes are connected. Besides that, the aggregation process is the only nonlinear cloud ice process rate in the ICON-AES. Since we focused on cloud-ice-related processes we just implemented the subgrid-scale approach in the aggregation parameterization. In future studies one can think about additionally including the subgrid-scale approach in cloud-water-related processes in order to see a
stronger effect on radiation, cloudiness, and precipitation as well.

As already mentioned, using a stochastic approach in non-linear process rates lowers the bias of the process rate. To create an unbiased aggregation rate (AGGsample), as described above, we make use of the entire cloud ice distribution function. Figure 9 shows a comparison of the averaged aggregation process rate at different $q_i$ bins and the 2D histograms of CTRL, AGGstoch, and AGGsample. AGGstoch shows good agreement with AGGsample, while CTRL produces much lower process rates. Both AGGstoch and AGGsample produce a higher aggregation rate compared to CTRL, especially for higher $q_i$ values. Figure 9b gives the corresponding joint histogram of CTRL and (c)–(d) the difference histogram to CTRL of cloud ice and the aggregation rate. AGGstoch shows a higher spread in the distribution, as already mentioned, while in AGGsample the distribution is shifted towards higher aggregation rates. However, both methods produce higher aggregation rates on average for different $q_i$ sections. The main difference between these methods is that AGGsample needs an additional integration over $q_i$ for the calculation of the aggregation rate, while AGGstoch needs no extra computational time. Therefore, the described stochastic method, which is used in AGGstoch, allows an improvement of the process rate bias against CTRL.

## 4 Summary

We introduce a stochastic approach of the aggregation parameterization in the ICON-AES by considering, consistently with the models’ cloud scheme, a uniform distribution function of cloud ice and randomly choose a cloud ice water content from this distribution function for each grid box and time step. The chosen distribution function of cloud ice is evaluated with the help of a combined lidar–radar dataset (DARDAR). An estimate of the precipitating and convective cloud ice mass is removed from the dataset in order to allow a more consistent comparison to the cloud ice in the model. The comparison of the simulated and the observed cloud ice variance shows good agreement in pattern, with just a few regional differences.

Overall we show that the uniform distribution function of cloud ice is usable for the microphysical process rates, e.g., the aggregation. From this uniform distribution function a randomly chosen cloud ice value is implemented into the nonlinear aggregation rate in order to represent subgrid-scale variability. Due to this stochastic method, the aggregation rate is intensified on average, since aggregation is a non-linear process. As a result, more cloud ice is transformed to snow, which leads to a cloud ice loss. However, the decrease in the accretion rate that results from using the stochastic aggregation scheme acts against the more intense aggregation rate. Therefore, the cloud ice loss is not as strong as expected. This indicates that changing only one of the two

process rates of snow formation does not lead to a change as large as one might have initially expected because aggregation and accretion interact very strongly with each other. However, the effect of taking subgrid-scale variability into account for process rates has a significant impact on the microphysical process (20% stronger global averaged process rate).

An unbiased process rate is calculated by integrating the aggregation rate over the entire cloud ice distribution function. The new stochastic method shows better agreement of the aggregation rates with the unbiased method compared to the control run. It follows that the new method lowers the bias of the aggregation rate, which does not need additional computational time. Therefore, this study suggests that using a stochastic approach for microphysical process rates helps to improve the representation of clouds and precipitation processes in global climate models.

## Appendix A: Variance of cloud ice

Here, the calculation of cloud ice variance, which is used in the ICON-AES, is given step by step. The general equation is written as follows:

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} (q_i - \bar{q}_i)^2 \text{PDF}(q_i) dq_i. \quad (A1) $$

Separating the integral in the cloud-free part, which is expressed as a Dirac function, and the cloudy parts yields

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} \left[ \underbrace{(1-C)\delta(q_i=0)}_{\text{cloud-free part}} + \underbrace{C \frac{1}{2\Delta q_i}}_{\text{cloudy part}} \right] (q_i - \bar{q}_i)^2 dq_i. \quad (A2) $$

Solving the Dirac function of the cloud-free part yields

$$ \sigma_{q_i}^2 = (1-C) \bar{q}_i^2 + \frac{C}{2\Delta q_i} \int_0^{2\Delta q_i} (q_i - \bar{q}_i)^2 dq_i. \quad (A3) $$

Solving the integral with $\bar{q}_i = C \Delta q_i$ yields the following (see Eq. 8).

$$ \begin{aligned} \sigma_{q_i}^2 &= (1-C) \bar{q}_i^2 + \frac{C}{6\Delta q_i} \left[ (2\Delta q_i - \bar{q}_i)^3 + \bar{q}_i^3 \right] \\ \sigma_{q_i}^2 &= (1-C) (C \Delta q_i)^2 + \frac{C}{6\Delta q_i} \left[ (2\Delta q_i - C \Delta q_i)^3 + (C \Delta q_i)^3 \right] \\ \sigma_{q_i}^2 &= (1-C) (C \Delta q_i)^2 + \frac{C}{6\Delta q_i} \left[ 8\Delta q_i^3 - 12C \Delta q_i^3 + 6C^2 \Delta q_i^3 \right] \\ \sigma_{q_i}^2 &= (\Delta q_i)^2 \cdot \left( \frac{4}{3} C - C^2 \right) \end{aligned} \quad (A4) $$
stronger effect on radiation, cloudiness, and precipitation as well.

As already mentioned, using a stochastic approach in non-linear process rates lowers the bias of the process rate. To create an unbiased aggregation rate (AGGsample), as described above, we make use of the entire cloud ice distribution function. Figure 9 shows a comparison of the averaged aggregation process rate at different $q_i$ bins and the 2D histograms of CTRL, AGGstoch, and AGGsample. AGGstoch shows good agreement with AGGsample, while CTRL produces much lower process rates. Both AGGstoch and AGGsample produce a higher aggregation rate compared to CTRL, especially for higher $q_i$ values. Figure 9b gives the corresponding joint histogram of CTRL and (c)–(d) the difference histogram to CTRL of cloud ice and the aggregation rate. AGGstoch shows a higher spread in the distribution, as already mentioned, while in AGGsample the distribution is shifted towards higher aggregation rates. However, both methods produce higher aggregation rates on average for different $q_i$ sections. The main difference between these methods is that AGGsample needs an additional integration over $q_i$ for the calculation of the aggregation rate, while AGGstoch needs no extra computational time. Therefore, the described stochastic method, which is used in AGGstoch, allows an improvement of the process rate bias against CTRL.

## 4 Summary

We introduce a stochastic approach of the aggregation parameterization in the ICON-AES by considering, consistently with the models’ cloud scheme, a uniform distribution function of cloud ice and randomly choose a cloud ice water content from this distribution function for each grid box and time step. The chosen distribution function of cloud ice is evaluated with the help of a combined lidar–radar dataset (DARDAR). An estimate of the precipitating and convective cloud ice mass is removed from the dataset in order to allow a more consistent comparison to the cloud ice in the model. The comparison of the simulated and the observed cloud ice variance shows good agreement in pattern, with just a few regional differences.

Overall we show that the uniform distribution function of cloud ice is usable for the microphysical process rates, e.g., the aggregation. From this uniform distribution function a randomly chosen cloud ice value is implemented into the nonlinear aggregation rate in order to represent subgrid-scale variability. Due to this stochastic method, the aggregation rate is intensified on average, since aggregation is a non-linear process. As a result, more cloud ice is transformed to snow, which leads to a cloud ice loss. However, the decrease in the accretion rate that results from using the stochastic aggregation scheme acts against the more intense aggregation rate. Therefore, the cloud ice loss is not as strong as expected. This indicates that changing only one of the two

process rates of snow formation does not lead to a change as large as one might have initially expected because aggregation and accretion interact very strongly with each other. However, the effect of taking subgrid-scale variability into account for process rates has a significant impact on the microphysical process (20% stronger global averaged process rate).

An unbiased process rate is calculated by integrating the aggregation rate over the entire cloud ice distribution function. The new stochastic method shows better agreement of the aggregation rates with the unbiased method compared to the control run. It follows that the new method lowers the bias of the aggregation rate, which does not need additional computational time. Therefore, this study suggests that using a stochastic approach for microphysical process rates helps to improve the representation of clouds and precipitation processes in global climate models.

## Appendix A: Variance of cloud ice

Here, the calculation of cloud ice variance, which is used in the ICON-AES, is given step by step. The general equation is written as follows:

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} (q_i - \bar{q}_i)^2 \text{PDF}(q_i) dq_i. \quad (A1) $$

Separating the integral in the cloud-free part, which is expressed as a Dirac function, and the cloudy parts yields

$$ \sigma_{q_i}^2 = \int_0^{2\Delta q_i} \left[ \underbrace{(1-C) \delta(q_i=0)}_{\text{cloud-free part}} + \underbrace{C \frac{1}{2\Delta q_i}}_{\text{cloudy part}} \right] (q_i - \bar{q}_i)^2 dq_i. \quad (A2) $$

Solving the Dirac function of the cloud-free part yields

$$ \sigma_{q_i}^2 = (1-C) \bar{q}_i^2 + \frac{C}{2\Delta q_i} \int_0^{2\Delta q_i} (q_i - \bar{q}_i)^2 dq_i. \quad (A3) $$

Solving the integral with $\bar{q}_i = C \Delta q_i$ yields the following (see Eq. 8).

$$ \begin{aligned} \sigma_{q_i}^2 &= (1-C) \bar{q}_i^2 + \frac{C}{6\Delta q_i} \left[ (2\Delta q_i - \bar{q}_i)^3 + \bar{q}_i^3 \right] \\ \sigma_{q_i}^2 &= (1-C) (C \Delta q_i)^2 + \frac{C}{6\Delta q_i} \left[ (2\Delta q_i - C \Delta q_i)^3 + (C \Delta q_i)^3 \right] \\ \sigma_{q_i}^2 &= (1-C) (C \Delta q_i)^2 + \frac{C}{6\Delta q_i} \left[ 8\Delta q_i^3 - 12C \Delta q_i^3 + 6C^2 \Delta q_i^3 \right] \\ \sigma_{q_i}^2 &= (\Delta q_i)^2 \cdot \left( \frac{4}{3} C - C^2 \right) \end{aligned} \quad (A4) $$
<figure><img src="image_11.png" /><figcaption>Figure 9. Averaged aggregation rate at different cloud ice mixing ratio bins and the corresponding 2D histograms for CTRL, AGGstoch, and AGGsample.</figcaption></figure>

Replacing $\Delta q_i$ with Eq. (10) yields

$$ \sigma_{\bar{q}_i}^2 = (\bar{q}_i^c)^2 \cdot \left( \frac{4}{3} C - C^2 \right). \qquad (A5) $$

In the end, the variance of cloud ice just depends on the cloud cover and grid box mean of in-cloud ice. So, the variance, which is calculated for the aggregation process, can be checked directly for the specific amount of cloud ice available for the aggregation.

**Code availability.** The source code of this study is available at https://doi.org/10.5281/zenodo.10980022 (Doktorowski, 2024). Updates of the code can be found at https://github.com/shoernig/code_subgrid_scale (last access: 16 April 2024). Information about the license and the respective README files are also included. Access to ICON source code can be obtained at https://doi.org/10.35089/WDCC/IconRelease01 (ICON partnership, 2024).

**Data availability.** The 2B-CLDCLASS (Sassen and Wang, 2008) CloudSat data and DARDAR-Nice dataset (https://doi.org/10.25326/09, Sourdeval et al., 2018) are downloaded from the ICARE Data and Services Center Archive (https://www.icare.univ-lille.fr/asd-content/archive/, last access: 16 April 2024), which requires a user account.

**Author contributions.** SD and JQ conceived this study. JK provided expertise on how to run the model. OS gave advice about the satellite data. JK, JQ, and MS helped by checking the theoretical approach of this study. SD prepared the article with contributions from all co-authors.

**Competing interests.** The contact author has declared that none of the authors has any competing interests.

**Disclaimer.** Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.

**Acknowledgements.** This work used resources of the German Climate Computing Centre (Deutsches Klimarechenzentrum, DKRZ) granted by its Scientific Steering Committee (WLA) under project ID 1143. We thank the AERIS/ICARE Data and Services Center for providing access to the DARDAR and 2B-CLDCLASS data used in this study. We thank Takuro Michibata and an anonymous reviewer for comments that led to improvements in the analysis and paper.

**Financial support.** This research received funding from the German Research Foundation (Deutsche Forschungsgemeinschaft, DFG) via its priority program "Fusion of Radar Polarimetry and Atmospheric Modelling" (SPP-2115, PROM), specifically the project PARA (FKZ QU 311/21-1). Johannes Quaas further acknowledges funding by the EU Horizon 2020 project "FORCES" (GA 821205).

**Review statement.** This paper was edited by Sylwester Arabas and reviewed by Takuro Michibata and one anonymous referee.
<figure><img src="image_11.png" /><figcaption>Figure 9. Averaged aggregation rate at different cloud ice mixing ratio bins and the corresponding 2D histograms for CTRL, AGGstoch, and AGGsample.</figcaption></figure>

Replacing $\Delta q_i$ with Eq. (10) yields

$$ \sigma_{\bar{q}_i}^2 = (\bar{q}_i^c)^2 \cdot \left(\frac{4}{3}C - C^2\right). \qquad (A5) $$

In the end, the variance of cloud ice just depends on the cloud cover and grid box mean of in-cloud ice. So, the variance, which is calculated for the aggregation process, can be checked directly for the specific amount of cloud ice available for the aggregation.

**Code availability.** The source code of this study is available at https://doi.org/10.5281/zenodo.10980022 (Doktorowski, 2024). Updates of the code can be found at https://github.com/shoernig/code_subgrid_scale (last access: 16 April 2024). Information about the license and the respective README files are also included. Access to ICON source code can be obtained at https://doi.org/10.35089/WDCC/IconRelease01 (ICON partnership, 2024).

**Data availability.** The 2B-CLDCLASS (Sassen and Wang, 2008) CloudSat data and DARDAR-Nice dataset (https://doi.org/10.25326/09, Sourdeval et al., 2018) are downloaded from the ICARE Data and Services Center Archive (https://www.icare.univ-lille.fr/asd-content/archive/, last access: 16 April 2024), which requires a user account.

**Author contributions.** SD and JQ conceived this study. JK provided expertise on how to run the model. OS gave advice about the satellite data. JK, JQ, and MS helped by checking the theoretical approach of this study. SD prepared the article with contributions from all co-authors.

**Competing interests.** The contact author has declared that none of the authors has any competing interests.

**Disclaimer.** Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.

**Acknowledgements.** This work used resources of the German Climate Computing Centre (Deutsches Klimarechenzentrum, DKRZ) granted by its Scientific Steering Committee (WLA) under project ID 1143. We thank the AERIS/ICARE Data and Services Center for providing access to the DARDAR and 2B-CLDCLASS data used in this study. We thank Takuro Michibata and an anonymous reviewer for comments that led to improvements in the analysis and paper.

**Financial support.** This research received funding from the German Research Foundation (Deutsche Forschungsgemeinschaft, DFG) via its priority program "Fusion of Radar Polarimetry and Atmospheric Modelling" (SPP-2115, PROM), specifically the project PARA (FKZ QU 311/21-1). Johannes Quaas further acknowledges funding by the EU Horizon 2020 project "FORCES" (GA 821205).

**Review statement.** This paper was edited by Sylwester Arabas and reviewed by Takuro Michibata and one anonymous referee.
## References

Bergeron, T.: On the physics of clouds and precipitation, Proc. 5th Assembly UGGI, Lisbon, Portugal, 156–180, https://worldcat.org/oclc/31921934 (last access: 16 April 2024), 1935.

Berner, J., Achatz, U., Batté, L., Bengtsson, L., de la Cámara, A., Christensen, H. M., Colangeli, M., Coleman, D. R. B., Crommelin, D., Dolaptchiev, S. I., Franzke, C. L. E., Friederichs, P., Imkeller, P., Järvinen, H., Juricke, S., Kitsios, V., Lott, F., Lucarini, V., Mahajan, S., Palmer, T. N., Penland, C., Sakradzija, M., von Storch, J.-S., Weisheimer, A., Weniger, M., and Williams, P. D.: Stochastic parameterization: Toward a new view of weather and climate models, B. Am. Math. Soc., 98, 565–588, https://doi.org/10.1175/BAMS-D-15-00268.1, 2017.

Boutle, I. A., Abel, S. J., Hill, P. G., and Morcrette, C. J.: Spatial variability of liquid cloud and rain: observations and microphysical effects, Q. J. Roy. Meteor. Soc., 140, 583–594, https://doi.org/10.1002/qj.2140, 2014.

Crueger, T., Giorgetta, M. A., Brokopf, R., Esch, M., Fiedler, S., Hohenegger, C., Kornblueh, L., Mauritsen, T., Nam, C., Naumann, A. K., Peters, K., Rast, S., Roeckner, E., Sakradzija, M., Schmidt, H., Vial, J., Vogel, R., and Stevens, B.: ICON-A, The Atmosphere Component of the ICON Earth System Model: II. Model Evaluation, J. Adv. Model Earth Sy., 10, 1638–1662, https://doi.org/10.1029/2017MS001233, 2018.

Delanoë, J. and Hogan, R. J.: A variational scheme for retrieving ice cloud properties from combined radar, lidar, and infrared radiometer, J. Geophys. Res., 113, D07204, https://doi.org/10.1029/2007JD009000, 2008.

Delanoë, J. and Hogan, R. J.: Combined CloudSat-CALIPSO-MODIS retrievals of the properties of ice clouds, J. Geophys. Res., 115, D00H29, https://doi.org/10.1029/2009JD012346, 2010.

Doktorowski, S.: Subgrid-scale variability of cloud ice in the ICON-AES 1.3.00, Zenodo [code], https://doi.org/10.5281/zenodo.10980022, 2024.

Findeisen, W., Volken, E., Giesche, A. M., and Brönnimann, S.: Colloidal meteorological processes in the formation of precipitation, Meteorol. Z., 24, 443–454, https://doi.org/10.1127/metz/2015/0675, 2015.

Giorgetta, M. A., Roeckner, E., Mauritsen, T., Bader, J., Crueger, T., Esch, M., Rast, S., Kornblueh, L., Schmidt, H., Kinne, S., Hohenegger, C., Möbis, B., Krismer, T., Wieners, K.-H., and Stevens, B.: The atmospheric general circulation model ECHAM6-model description, Berichte zur Erdsystemforschung/Max-Planck-Institut für Meteorologie, 135, 172 pp., https://doi.org/10.17617/2.1810480, 2013.

Giorgetta, M. A., Brokopf, R., Crueger, T., Esch, M., Fiedler, S., Helmert, J., Hohenegger, C., Kornblueh, L., Köhler, M., Manzini, E., Mauritsen, T., Nam, C., Raddatz, T., Rast, S., Reinhert, D., Sakradzija, M., Schmidt, H., Schneck, R., Schnur, R., Silvers, L., Wan, H., Zängl, G., and Stevens, B.: ICON-A, the Atmosphere Component of the ICON Earth System Model: I. Model Description, J. Adv. Model Earth Sy., 10, 1613–1637, https://doi.org/10.1029/2017MS001242, 2018.

Golaz, J.-C., Larson, V. E., and Cotton, W. R.: A PDF-based model for boundary layer clouds. Part I: Method and model description, J. Atmos. Sci., 59, 3540–3551, https://doi.org/10.1175/1520-0469(2002)059<3540:APBMFB>2.0.CO;2, 2002a.

Golaz, J.-C., Larson, V. E., and Cotton, W. R.: A PDF-based model for boundary layer clouds. Part II: Model results, J. Atmos. Sci., 59, 3552–3571, https://doi.org/10.1175/1520-0469(2002)059<3552:APBMFB>2.0.CO;2, 2002b.

Hill, P. G., Morcrette, C. J., and Boutle, I. A.: A regime-dependent parametrization of subgrid-scale cloud water content variability, Q. J. Roy. Meteor. Soc., 141, 1975–1986, https://doi.org/10.1002/qj.2506, 2015.

ICON partnership (DWD; MPI-M; DKRZ; KIT; C2SM): ICON release 2024.01, World Data Center for Climate (WDCC) at DKRZ [code], https://doi.org/10.35089/WDCC/IconRelease01, 2024.

Imura, Y. and Michibata, T.: Too Frequent and Too Light Arctic Snowfall With Incorrect Precipitation Phase Partitioning in the MIROC6 GCM, J. Adv. Model. Earth Sy., 14, e2022MS003046, https://doi.org/10.1029/2022MS003046, 2022.

Larson, V. E. and Griffin, B. M.: Analytic upscaling of a local microphysics scheme. Part I: Derivation, Q. J. Roy. Meteor. Soc., 139, 46–57, https://doi.org/10.1002/qj.1967, 2013.

Larson, V. E., Wood, R., Field, P. R., Golaz, J.-C., Haar, T. H. V., and Cotton, W. R.: Systematic biases in the microphysics and thermodynamics of numerical models that ignore subgrid-scale variability, J. Atmos. Sci., 58, 1117–1128, https://doi.org/10.1175/1520-0469(2001)058<1117:SBITMA>2.0.CO;2, 2001.

Lebsock, M., Morrison, H., and Gettelman, A.: Microphysical implications of cloud-precipitation covariance derived from satellite remote sensing, J. Geophys. Res.-Atmos., 118, 6521–6533, https://doi.org/10.1002/jgrd.50347, 2013.

Levkov, L., Rockel, B., Kapitza, H., and Raschke, E.: 3D mesoscale numerical studies of cirrus and stratus clouds by their time and space evolution, Beitr. Phys. Atmosph., 65, 35–58, 1992.

Li, J.-L., Waliser, D., Chen, W.-T., Guan, B., Kubar, T., Stephens, G., Ma, H.-Y., Deng, M., Donner, L., Seman, C., and Horowitz, L.: An observationally based evaluation of cloud ice water in CMIP3 and CMIP5 GCMs and contemporary reanalyses using contemporary satellite data, J. Geophys. Res., 117, D16105, https://doi.org/10.1029/2012JD017640, 2012.

Li, J.-L. F., Waliser, D., Woods, C., Teixeira, J., Bacmeister, J., Chern, J., Shen, B.-W., Tompkins, A., Tao, W.-K., and Köhler, M.: Comparisons of satellites liquid water estimates to ECMWF and GMAO analyses, 20th century IPCC AR4 climate simulations, and GCM simulations, Geophys. Res. Lett., 35, L19710, https://doi.org/10.1029/2008GL035427, 2008.

Lohmann, U. and Roeckner, E.: Design and performance of a new cloud microphysics scheme developed for the ECHAM general circulation model, Clim. Dynam., 12, 557–572, https://doi.org/10.1007/BF00207939, 1996.

Michibata, T., Suzuki, K., and Takemura, T.: Snow-induced buffering in aerosol-cloud interactions, Atmos. Chem. Phys., 20, 13771–13780, https://doi.org/10.5194/acp-20-13771-2020, 2020.

Morrison, H. and Gettelman, A.: A New Two-Moment Bulk Stratiform Cloud Microphysics Scheme in the Community Atmosphere Model, Version 3 (CAM3). Part I: Description and Numerical Tests, J. Climate, 21, 3642–3659, https://doi.org/10.1175/2008JCLI2105.1, 2008.

Mülmenstädt, J., Sourdeval, O., Delanoë, J., and Quaas, J.: Frequency of occurrence of rain from liquid-, mixed-, and ice-phase clouds derived from A-Train satellite retrievals, Geophys. Res.
## References

Bergeron, T.: On the physics of clouds and precipitation, Proc. 5th Assembly UGGI, Lisbon, Portugal, 156–180, https://worldcat.org/oclc/31921934 (last access: 16 April 2024), 1935.

Berner, J., Achatz, U., Batté, L., Bengtsson, L., de la Cámara, A., Christensen, H. M., Colangeli, M., Coleman, D. R. B., Crommelin, D., Dolaptchiev, S. I., Franzke, C. L. E., Friederichs, P., Imkeller, P., Järvinen, H., Juricke, S., Kitsios, V., Lott, F., Lucarini, V., Mahajan, S., Palmer, T. N., Penland, C., Sakradzija, M., von Storch, J.-S., Weisheimer, A., Weniger, M., and Williams, P. D.: Stochastic parameterization: Toward a new view of weather and climate models, B. Am. Math. Soc., 98, 565–588, https://doi.org/10.1175/BAMS-D-15-00268.1, 2017.

Boutle, I. A., Abel, S. J., Hill, P. G., and Morcrette, C. J.: Spatial variability of liquid cloud and rain: observations and microphysical effects, Q. J. Roy. Meteor. Soc., 140, 583–594, https://doi.org/10.1002/qj.2140, 2014.

Crueger, T., Giorgetta, M. A., Brokopf, R., Esch, M., Fiedler, S., Hohenegger, C., Kornblueh, L., Mauritsen, T., Nam, C., Naumann, A. K., Peters, K., Rast, S., Roeckner, E., Sakradzija, M., Schmidt, H., Vial, J., Vogel, R., and Stevens, B.: ICON-A, The Atmosphere Component of the ICON Earth System Model: II. Model Evaluation, J. Adv. Model Earth Sy., 10, 1638–1662, https://doi.org/10.1029/2017MS001233, 2018.

Delanoë, J. and Hogan, R. J.: A variational scheme for retrieving ice cloud properties from combined radar, lidar, and infrared radiometer, J. Geophys. Res., 113, D07204, https://doi.org/10.1029/2007JD009000, 2008.

Delanoë, J. and Hogan, R. J.: Combined CloudSat-CALIPSO-MODIS retrievals of the properties of ice clouds, J. Geophys. Res., 115, D00H29, https://doi.org/10.1029/2009JD012346, 2010.

Doktorowski, S.: Subgrid-scale variability of cloud ice in the ICON-AES 1.3.00, Zenodo [code], https://doi.org/10.5281/zenodo.10980022, 2024.

Findeisen, W., Volken, E., Giesche, A. M., and Brönnimann, S.: Colloidal meteorological processes in the formation of precipitation, Meteorol. Z., 24, 443–454, https://doi.org/10.1127/metz/2015/0675, 2015.

Giorgetta, M. A., Roeckner, E., Mauritsen, T., Bader, J., Crueger, T., Esch, M., Rast, S., Kornblueh, L., Schmidt, H., Kinne, S., Hohenegger, C., Möbis, B., Krismer, T., Wieners, K.-H., and Stevens, B.: The atmospheric general circulation model ECHAM6-model description, Berichte zur Erdsystemforschung/Max-Planck-Institut für Meteorologie, 135, 172 pp., https://doi.org/10.17617/2.1810480, 2013.

Giorgetta, M. A., Brokopf, R., Crueger, T., Esch, M., Fiedler, S., Helmert, J., Hohenegger, C., Kornblueh, L., Köhler, M., Manzini, E., Mauritsen, T., Nam, C., Raddatz, T., Rast, S., Reinhert, D., Sakradzija, M., Schmidt, H., Schneck, R., Schnur, R., Silvers, L., Wan, H., Zängl, G., and Stevens, B.: ICON-A, the Atmosphere Component of the ICON Earth System Model: I. Model Description, J. Adv. Model Earth Sy., 10, 1613–1637, https://doi.org/10.1029/2017MS001242, 2018.

Golaz, J.-C., Larson, V. E., and Cotton, W. R.: A PDF-based model for boundary layer clouds. Part I: Method and model description, J. Atmos. Sci., 59, 3540–3551, https://doi.org/10.1175/1520-0469(2002)059<3540:APBMFB>2.0.CO;2, 2002a.

Golaz, J.-C., Larson, V. E., and Cotton, W. R.: A PDF-based model for boundary layer clouds. Part II: Model results, J. Atmos. Sci., 59, 3552–3571, https://doi.org/10.1175/1520-0469(2002)059<3552:APBMFB>2.0.CO;2, 2002b.

Hill, P. G., Morcrette, C. J., and Boutle, I. A.: A regime-dependent parametrization of subgrid-scale cloud water content variability, Q. J. Roy. Meteor. Soc., 141, 1975–1986, https://doi.org/10.1002/qj.2506, 2015.

ICON partnership (DWD; MPI-M; DKRZ; KIT; C2SM): ICON release 2024.01, World Data Center for Climate (WDCC) at DKRZ [code], https://doi.org/10.35089/WDCC/IconRelease01, 2024.

Imura, Y. and Michibata, T.: Too Frequent and Too Light Arctic Snowfall With Incorrect Precipitation Phase Partitioning in the MIROC6 GCM, J. Adv. Model. Earth Sy., 14, e2022MS003046, https://doi.org/10.1029/2022MS003046, 2022.

Larson, V. E. and Griffin, B. M.: Analytic upscaling of a local microphysics scheme. Part I: Derivation, Q. J. Roy. Meteor. Soc., 139, 46–57, https://doi.org/10.1002/qj.1967, 2013.

Larson, V. E., Wood, R., Field, P. R., Golaz, J.-C., Haar, T. H. V., and Cotton, W. R.: Systematic biases in the microphysics and thermodynamics of numerical models that ignore subgrid-scale variability, J. Atmos. Sci., 58, 1117–1128, https://doi.org/10.1175/1520-0469(2001)058<1117:SBITMA>2.0.CO;2, 2001.

Lebsock, M., Morrison, H., and Gettelman, A.: Microphysical implications of cloud-precipitation covariance derived from satellite remote sensing, J. Geophys. Res.-Atmos., 118, 6521–6533, https://doi.org/10.1002/jgrd.50347, 2013.

Levkov, L., Rockel, B., Kapitza, H., and Raschke, E.: 3D mesoscale numerical studies of cirrus and stratus clouds by their time and space evolution, Beitr. Phys. Atmosph., 65, 35–58, 1992.

Li, J.-L., Waliser, D., Chen, W.-T., Guan, B., Kubar, T., Stephens, G., Ma, H.-Y., Deng, M., Donner, L., Seman, C., and Horowitz, L.: An observationally based evaluation of cloud ice water in CMIP3 and CMIP5 GCMs and contemporary reanalyses using contemporary satellite data, J. Geophys. Res., 117, D16105, https://doi.org/10.1029/2012JD017640, 2012.

Li, J.-L. F., Waliser, D., Woods, C., Teixeira, J., Bacmeister, J., Chern, J., Shen, B.-W., Tompkins, A., Tao, W.-K., and Köhler, M.: Comparisons of satellites liquid water estimates to ECMWF and GMAO analyses, 20th century IPCC AR4 climate simulations, and GCM simulations, Geophys. Res. Lett., 35, L19710, https://doi.org/10.1029/2008GL035427, 2008.

Lohmann, U. and Roeckner, E.: Design and performance of a new cloud microphysics scheme developed for the ECHAM general circulation model, Clim. Dynam., 12, 557–572, https://doi.org/10.1007/BF00207939, 1996.

Michibata, T., Suzuki, K., and Takemura, T.: Snow-induced buffering in aerosol-cloud interactions, Atmos. Chem. Phys., 20, 13771–13780, https://doi.org/10.5194/acp-20-13771-2020, 2020.

Morrison, H. and Gettelman, A.: A New Two-Moment Bulk Stratiform Cloud Microphysics Scheme in the Community Atmosphere Model, Version 3 (CAM3). Part I: Description and Numerical Tests, J. Climate, 21, 3642–3659, https://doi.org/10.1175/2008JCLI2105.1, 2008.

Mülmenstädt, J., Sourdeval, O., Delanoë, J., and Quaas, J.: Frequency of occurrence of rain from liquid-, mixed-, and ice-phase clouds derived from A-Train satellite retrievals, Geophys. Res.
Lett., 42, 6502–6509, https://doi.org/10.1002/2015GL064604, 2015.

Murakami, M.: Numerical modeling of dynamical and microphysical evolution of an isolated convective cloud, J. Meteorol. Soc. Jpn., 68, 107–128, https://doi.org/10.2151/jmsj1965.68.2_107, 1990.

Palmer, T. N.: A nonlinear dynamical perspective on model error: A proposal for non-local stochastic-dynamic parametrization in weather and climate prediction models, Q. J. Roy. Meteor. Soc., 127, 279–304, https://doi.org/10.1002/qj.49712757202, 2001.

Pincus, R. and Klein, S. A.: Unresolved spatial variability and microphysical process rates in large-scale models, J. Geophys. Res.-Atmos., 105, 27059–27065, https://doi.org/10.1029/2000JD900504, 2000.

Quaas, J.: Evaluating the “critical relative humidity” as a measure of subgrid-scale variability of humidity in general circulation model cloud cover parameterizations using satellite data, J. Geophys. Res.-Atmos., 117, D09208, https://doi.org/10.1029/2012JD017495, 2012.

Sassen, K. and Wang, Z.: Classifying clouds around the globe with the CloudSat radar: 1-year of results, Geophys. Res. Lett., 35, L04805, https://doi.org/10.1029/2007GL032591, 2008.

Sourdeval, O., Gryspeerdt, E., Krämer, M., Goren, T., Delanoë, J., Afchine, A., Hemmer, F., and Quaas, J.: Ice crystal number concentration from satellite lidar-radar observations (DARDARNice), AERIS [data set], https://doi.org/10.25326/09, 2018.

Stephens, G. L., Vane, D. G., Boain, R. J., Mace, G. G., Sassen, K., Wang, Z., Illingworth, A. J., O’connor, E. J., Rossow, W. B., Durden, S. L., Miller, S. D., Austin, R. T., Benedetti, A., Mitrescu, C., and the CloudSat Science Team: THE CLOUDSAT MISSION AND THE A-TRAIN: A New Dimension of Space-Based Observations of Clouds and Precipitation, B. Am. Meteorol. Soc., 83, 1771–1790, https://doi.org/10.1175/BAMS-83-12-1771, 2002.

Stevens, B., Giorgetta, M., Esch, M., Mauritsen, T., Crueger, T., Rast, S., Salzmann, M., Schmidt, H., Bader, J., Block, K., Brokopf, R., Fast, I., Kinne, S., Kornblueh, L., Lohmann, U., Pincus, R., Reichler, R., and Roeckner, E.: Atmospheric component of the MPI-M Earth System Model: ECHAM6, J. Adv. Model Earth Sy., 5, 146–172, https://doi.org/10.1002/jame.20015, 2013.

Sundqvist, H., Berge, E., and Kristjánsson, J. E.: Condensation and cloud parameterization studies with a mesoscale numerical weather prediction model, Mon. Weather. Rev., 117, 1641–1657, https://doi.org/10.1175/1520-0493(1989)117<1641:CACPSW>2.0.CO;2, 1989.

Thayer-Calder, K., Gettelman, A., Craig, C., Goldhaber, S., Bogen-schutz, P. A., Chen, C.-C., Morrison, H., Höft, J., Raut, E., Grif-fin, B. M., Weber, J. K., Larson, V. E., Wyant, M. C., Wang, M., Guo, Z., and Ghan, S. J.: A unified parameterization of clouds and turbulence using CLUBB and subcolumns in the Community Atmosphere Model, Geosci. Model Dev., 8, 3801–3821, https://doi.org/10.5194/gmd-8-3801-2015, 2015.

Trömel, S., Simmer, C., Blahak, U., Blanke, A., Doktorowski, S., Ewald, F., Frech, M., Gergely, M., Hagen, M., Janjic, T., Kalesse-Los, H., Kneifel, S., Knote, C., Mendrok, J., Moser, M., Köcher, G., Mühlbauer, K., Myagkov, A., Pejcic, V., Seifert, P., Shrestha, P., Teisseire, A., von Terzi, L., Tetoni, E., Vogl, T., Voigt, C., Zeng, Y., Zinner, T., and Quaas, J.: Overview: Fusion of radar polarimetry and numerical atmospheric modelling towards an improved understanding of cloud and precipitation processes, Atmos. Chem. Phys., 21, 17291–17314, https://doi.org/10.5194/acp-21-17291-2021, 2021.

Weber, T. and Quaas, J.: Incorporating the subgrid-scale variability of clouds in the autoconversion parameterization using a PDF-scheme, J. Adv. Model Earth Sy., 4, M11003, https://doi.org/10.1029/2012MS000156, 2012.

Wegener, A.: Thermodynamik der atmosphäre, JA Barth, http://worldcat.org/oclc/39667532 (last access: 16 April 2024), 1911.

Winker, D. M., Pelon, J., Coakley, J. A., Ackerman, S. A., Charlson, R. J., Colarco, P. R., Flamant, P., Fu, Q., Hoff, R. M., Kit-taka, C., Kubar, T. L., Treut, H. L., Mccormick, M. P., Mégie, G., Poole, L., Powell, K., Trepte, C., Vaughan, M. A., and Wielicki, B. A.: The CALIPSO Mission: A Global 3D View of Aerosols and Clouds, B. Am. Meteorol. Soc., 91, 1211–1230, https://doi.org/10.1175/2010BAMS3009.1, 2010.
Lett., 42, 6502–6509, https://doi.org/10.1002/2015GL064604, 2015.

Murakami, M.: Numerical modeling of dynamical and microphysical evolution of an isolated convective cloud, J. Meteorol. Soc. Jpn., 68, 107–128, https://doi.org/10.2151/jmsj1965.68.2_107, 1990.

Palmer, T. N.: A nonlinear dynamical perspective on model error: A proposal for non-local stochastic-dynamic parametrization in weather and climate prediction models, Q. J. Roy. Meteor. Soc., 127, 279–304, https://doi.org/10.1002/qj.49712757202, 2001.

Pincus, R. and Klein, S. A.: Unresolved spatial variability and microphysical process rates in large-scale models, J. Geophys. Res.-Atmos., 105, 27059–27065, https://doi.org/10.1029/2000JD900504, 2000.

Quaas, J.: Evaluating the “critical relative humidity” as a measure of subgrid-scale variability of humidity in general circulation model cloud cover parameterizations using satellite data, J. Geophys. Res.-Atmos., 117, D09208, https://doi.org/10.1029/2012JD017495, 2012.

Sassen, K. and Wang, Z.: Classifying clouds around the globe with the CloudSat radar: 1-year of results, Geophys. Res. Lett., 35, L04805, https://doi.org/10.1029/2007GL032591, 2008.

Sourdeval, O., Gryspeerdt, E., Krämer, M., Goren, T., Delanoë, J., Afchine, A., Hemmer, F., and Quaas, J.: Ice crystal number concentration from satellite lidar-radar observations (DARDARNice), AERIS [data set], https://doi.org/10.25326/09, 2018.

Stephens, G. L., Vane, D. G., Boain, R. J., Mace, G. G., Sassen, K., Wang, Z., Illingworth, A. J., O’connor, E. J., Rossow, W. B., Durden, S. L., Miller, S. D., Austin, R. T., Benedetti, A., Mitrescu, C., and the CloudSat Science Team: THE CLOUDSAT MISSION AND THE A-TRAIN: A New Dimension of Space-Based Observations of Clouds and Precipitation, B. Am. Meteorol. Soc., 83, 1771–1790, https://doi.org/10.1175/BAMS-83-12-1771, 2002.

Stevens, B., Giorgetta, M., Esch, M., Mauritsen, T., Crueger, T., Rast, S., Salzmann, M., Schmidt, H., Bader, J., Block, K., Brokopf, R., Fast, I., Kinne, S., Kornblueh, L., Lohmann, U., Pincus, R., Reichler, R., and Roeckner, E.: Atmospheric component of the MPI-M Earth System Model: ECHAM6, J. Adv. Model Earth Sy., 5, 146–172, https://doi.org/10.1002/jame.20015, 2013.

Sundqvist, H., Berge, E., and Kristjánsson, J. E.: Condensation and cloud parameterization studies with a mesoscale numerical weather prediction model, Mon. Weather. Rev., 117, 1641–1657, https://doi.org/10.1175/1520-0493(1989)117<1641:CACPSW>2.0.CO;2, 1989.

Thayer-Calder, K., Gettelman, A., Craig, C., Goldhaber, S., Bogen-schutz, P. A., Chen, C.-C., Morrison, H., Höft, J., Raut, E., Grif-fin, B. M., Weber, J. K., Larson, V. E., Wyant, M. C., Wang, M., Guo, Z., and Ghan, S. J.: A unified parameterization of clouds and turbulence using CLUBB and subcolumns in the Community Atmosphere Model, Geosci. Model Dev., 8, 3801–3821, https://doi.org/10.5194/gmd-8-3801-2015, 2015.

Trömel, S., Simmer, C., Blahak, U., Blanke, A., Doktorowski, S., Ewald, F., Frech, M., Gergely, M., Hagen, M., Janjic, T., Kalesse-Los, H., Kneifel, S., Knote, C., Mendrok, J., Moser, M., Köcher, G., Mühlbauer, K., Myagkov, A., Pejcic, V., Seifert, P., Shrestha, P., Teisseire, A., von Terzi, L., Tetoni, E., Vogl, T., Voigt, C., Zeng, Y., Zinner, T., and Quaas, J.: Overview: Fusion of radar polarimetry and numerical atmospheric modelling towards an improved understanding of cloud and precipitation processes, Atmos. Chem. Phys., 21, 17291–17314, https://doi.org/10.5194/acp-21-17291-2021, 2021.

Weber, T. and Quaas, J.: Incorporating the subgrid-scale variability of clouds in the autoconversion parameterization using a PDF-scheme, J. Adv. Model Earth Sy., 4, M11003, https://doi.org/10.1029/2012MS000156, 2012.

Wegener, A.: Thermodynamik der atmosphäre, JA Barth, http://worldcat.org/oclc/39667532 (last access: 16 April 2024), 1911.

Winker, D. M., Pelon, J., Coakley, J. A., Ackerman, S. A., Charlson, R. J., Colarco, P. R., Flamant, P., Fu, Q., Hoff, R. M., Kit-taka, C., Kubar, T. L., Treut, H. L., Mccormick, M. P., Mégie, G., Poole, L., Powell, K., Trepte, C., Vaughan, M. A., and Wielicki, B. A.: The CALIPSO Mission: A Global 3D View of Aerosols and Clouds, B. Am. Meteorol. Soc., 91, 1211–1230, https://doi.org/10.1175/2010BAMS3009.1, 2010.