Articles | Volume 17, issue 13
Model evaluation paper
05 Jul 2024
Model evaluation paper |  | 05 Jul 2024

Simulation of marine stratocumulus using the super-droplet method: numerical convergence and comparison to a double-moment bulk scheme using SCALE-SDM 5.2.6-2.3.1

Chongzhi Yin, Shin-ichiro Shima, Lulin Xue, and Chunsong Lu

The super-droplet method (SDM) is a Lagrangian particle-based numerical scheme for cloud microphysics. In this work, a series of simulations based on the DYCOMS-II (RF02) setup with different horizontal and vertical resolutions are conducted to explore the grid convergence of the SDM simulations of marine stratocumulus. The results are compared with the double-moment bulk scheme (SN14) and model intercomparison project (MIP) results. In general, all SDM and SN14 variables show a good agreement with the MIP results and have similar grid size dependencies. The stratocumulus simulation is more sensitive to the vertical resolution than to the horizontal resolution. The vertical grid length DZ  2.5 m is necessary for both SDM and SN14. The horizontal grid length DX < 12.5 m is necessary for the SDM simulations. DX  25 m is sufficient for SN14. We also assess the numerical convergence with respect to the super-droplet numbers. The simulations are well converged when the super-droplet number concentration (SDNC) is larger than 16 super-droplets per cell. Our results indicate that the super-droplet number per grid cell is more critical than that per unit volume at least for the stratocumulus case investigated here. Our comprehensive analysis not only offers guidance on numerical settings essential for accurate stratocumulus cloud simulation but also underscores significant differences in liquid water content and cloud macrostructure between SDM and SN14. These differences are attributed to the inherent modeling strategies of the two schemes. SDM's dynamic representation of aerosol size distribution through wet deposition markedly contrasts with SN14's static approach, influencing cloud structure and behavior over a 6 h simulation. Findings reveal sedimentation's crucial role in altering aerosol distributions near cloud tops, affecting the vertical profile of cloud fraction (CF). Additionally, the study briefly addresses numerical diffusion's potential effects, suggesting further investigation is needed. The results underscore the importance of accurate aerosol modeling and its interactions with cloud processes in marine stratocumulus simulations, pointing to future research directions for enhancing stratocumulus modeling accuracy and predictive capabilities.

1 Introduction

Marine stratocumulus clouds cover approximately one-quarter of the Earth's surface and play an important role in the planet's radiation budget (Wood, 2012; Matheou and Teixeira, 2019; Nowak et al., 2021). These clouds reflect the incident shortwave radiation and almost have no effect on the outgoing longwave radiation, resulting in a negative radiation flux (Wood, 2012). The temperature projection uncertainty in global warming simulation is mainly caused by the representation of marine low clouds in global climate models (Stephens, 2005; Bony and Dufresne, 2005; Bony et al., 2006; Boucher et al., 2013; Zelinka et al., 2020; Kawai and Shige, 2020); thus, the stratocumulus must be accurately represented in numerical models. The IPCC (AR6) states that aerosol–cloud-related processes introduce the greatest uncertainty among the radiative forcing assessment methods of major factors in the earth–atmosphere system. Therefore, we must understand the aerosol–cloud interaction of the stratocumulus.

In the cloud modeling community, two types of methods are commonly used to represent clouds in numerical models. The first type is to treat the cloud as a continuum in the Eulerian framework, namely Eulerian cloud models (ECMs). The second type is to treat the cloud as an ensemble of individual particles in the Lagrangian framework, that is, the Lagrangian particle-based cloud models (LCMs).

The bulk scheme (Kessler, 1969; Lin et al., 1983; Schoenberg Ferrier, 1994; Milbrandt and Yau, 2005) is one of the most widely used Eulerian microphysical schemes. It assumes a specific distribution (e.g., gamma distribution) to characterize the size distributions of aerosol and cloud particles; thus, only several predictors must be considered. This method is numerically efficient and saves computing resources. However, the cloud droplet size distribution (DSD) of the bulk scheme is a fixed and continuous function; thus, the calculation of microphysical processes depends on the set function properties, and the uncertainty of the cloud simulation results is high (Khain et al., 2015).

Another ECM category represents the cloud hydrometeors in discrete bins and is called the spectral bin microphysics scheme (Khain et al., 2000; Lynn et al., 2005; Morrison and Grabowski, 2010; Xue et al., 2010, 2012; Geresdi et al., 2017), which can explicitly predict the particle size or mass distribution but is computationally more costly. As a result, bin schemes suffer from the limitation of dimensionality (Shima, 2008; Shima et al., 2009; Grabowski et al., 2019). Most bin schemes are one-dimensional, which means they only predict the droplet size or mass distributions. The solute composition, mass, and soluble fraction within the cloud droplet all affect the droplet growth rate and determine the characteristics of particles remaining after the droplet has completely evaporated. In some cases, these factors are essential but are difficult to consider in the bin schemes (Shima, 2008; Shima et al., 2009; Grabowski et al., 2019; Dziekan et al., 2021). Another problem in bin microphysics comes from the limitation of the Smoluchowski equation used to represent the collision–condensation process (Smoluchowski, 1916). The Smoluchowski equation is deterministic, while the collision–coalescence of droplets is a stochastic process. Therefore, droplet collision, other than expected, can appear. The Smoluchowski equation also does not accurately predict even the mean behavior when the well-mixed volume is small, and the droplet discreteness is evident (see Alfonso and Raga, 2017; Dziekan and Pawlowska, 2017; Grabowski et al., 2019, and references therein). In addition, all ECMs are affected by numerical diffusion, which can lead to a simulated system that behaves differently from the expected physical system (Schoeffler, 1982). In bin microphysics, numerical diffusion results in the broadening of the unphysical droplet size distribution (Morrison et al., 2018; Grabowski et al., 2019). Due to the three abovementioned issues, ECMs still face difficulties in accurately simulating cloud microphysical processes. However, the recently developed Lagrangian particle-based method may be a viable solution for representing cloud and precipitation particle populations (Morrison et al., 2020).

Shima et al. (2009) proposed an LCM, called the super-droplet method (SDM), in which each super-droplet represents multiple numbers of aerosol/cloud/precipitation particles with the same attributes and position. The SDM has no numerical diffusion of liquid water and can provide more detailed microphysics information. Note that sub-grid-scale (SGS) diffusions are not represented in the original SDM, which may lead to under-diffused supersaturation and accelerate the mixing process (Grabowski and Abade, 2017; Abade et al., 2018; Hoffmann et al., 2019). The Monte Carlo collision–coalescence algorithm of the SDM is based on the stochastic process of collision–coalescence; hence, the SDM can be applied, even when the Smoluchowski equation is invalid (Dziekan and Pawlowska, 2017). Shima et al. (2009) theoretically estimated that when the number of attributes, which range from 2 to 4, becomes larger than a certain critical value, the SDM becomes computationally more efficient than the bin microphysics approach. With the increase in the supercomputer computing capacity, the number of studies using the SDM or other LCMs has increased in the past 10 years (e.g., Arabas and Shima, 2013; Naumann and Seifert, 2015; Dziekan and Pawlowska, 2017; Sato et al., 2017, 2018; Grabowski et al., 2018; Jaruga and Pawlowska, 2018; Schwenkel et al., 2018; Dziekan et al., 2019; Noh et al., 2018; Hoffmann and Feingold, 2019; Hoffmann et al., 2019; Seifert and Rasp, 2020; Unterstrasser et al., 2020; Shima et al., 2020; Dziekan et al., 2021; Richter et al., 2021; Chandrakar et al., 2022).

Various studies used ECMs to simulate marine stratocumulus. Some of them investigated the grid convergence characteristics during the simulation. In their work, Matheou et al. (2016) indicated that all flow statistics of stratocumulus simulated by the large eddy simulation (LES), except for those related to liquid water, converge for DX = DY = DZ < 2.5 m (i.e., DX and DY are the horizontal grid lengths, while DZ is the vertical grid length). A series of sensitivity experiments with seven numerical and physical parameters was conducted by Matheou and Teixeira (2019) to understand the source of difficulty in simulating stratocumulus by the LES. They not only used different grid spacings, but also changed the geostrophic wind, divergence, radiation parameterization, buoyancy formulation, surface fluxes, and the scalar advection numerical method. The grid convergence could merely be found at a very fine resolution. Moreover, the mean results of simulation of the finest resolution agree with the observations (Stevens et al., 2005). The entrainment rate and the mean profiles, except for the cloud liquid, were not sensitive to the grid resolution. The buoyancy perturbation run in the study of Matheou and Teixeira (2019) also suggested that the buoyancy reversal instability of the cloud top significantly enhances the entrainment rate. Some studies showed that a larger horizontal grid spacing leads to higher liquid water path (LWP) and cloud cover, whereas a larger vertical grid spacing has the opposite effect (Kurowski et al., 2009; Cheng et al., 2010; Yamaguchi and Randall, 2012; Pedersen et al., 2016).

A few studies employed LCMs for stratocumulus, but none of them investigated the grid convergence characteristics when using these LCMs for marine stratocumulus. Dziekan et al. (2019) studied the SDM sensitivity to the time steps of condensation and coalescence in two-dimensional simulations and compared the SGS turbulence models of different approaches in three-dimensional (3D) simulations using the setup of the second research flight of the second Dynamics and Chemistry of Marine Stratocumulus (DYCOMS-II (RF02)) study. They found that droplet condensational and collisional growth must be modeled with a 0.1 s time step. In addition, the simulation results using the Smagorinsky scheme and an algorithm for the SGS turbulent motion of computational particles were in the best agreement with the ECM results. They also tested various initial super-droplet number concentrations (SDNCs) ranging from 40 to 1000 super-droplets per grid cell and confirmed that 40 super-droplets per cell was sufficient in achieving the correct domain-averaged results for DYCOMS-II (RF02). However, they did not investigate the grid convergence. Hoffmann and Feingold (2019) applied a new modeling method L3 combining an LES, a linear eddy model, and an LCM to study stratocumulus. They found that the number of cloud holes (i.e., dry air parcels transported from the free atmosphere to the cloud layer) in the L3 simulation is higher and persists longer. Their simulations showed that reducing the number of cloud droplets during mixing results in larger remaining droplets. Their results also illustrated that inhomogeneous mixing does not increase the cloud droplet age because inhomogeneous mixing hinders the droplet evaporation at the cloud edge and makes the older droplets disappear from the cloud faster due to the faster sedimentation caused by diffusional growth. They did not assess the numerical convergence either but admitted that the vertical grid length of 35 m used in their study did not explicitly resolve all cloud holes. Another important factor they did not explicitly consider, which could affect the cloud hole persistence, is the impact of aerosols. They ignored the curvature and the solute effects in the condensational droplet growth and did not explicitly consider the activation–deactivation of aerosol particles. Chandrakar et al. (2022) studied the DSD evolution during the transition of closed cells to open cells through LES coupling with the SDM. They tracked the trajectories of some sample super-droplets and found that some droplets could rapidly grow to drizzle from the collision–coalescence process, mainly within downdraft. Their results showed that once the coalescence timescale becomes similar to the eddy turnover timescale, the coalescence growth could be enhanced, and then it increases a key driver of the closed-to-open cell transition. However, to save computational resources, they used a relatively coarse horizontal grid resolution of 100 m, which was too coarse for the precise cloud water simulation in a stratocumulus cloud (Matheou and Teixeira, 2019).

One of the aims of the present study is to determine how fine the spatial resolution must be for an accurate simulation of marine stratocumulus using the SDM. Several series of simulations based on the DYCOMS-II (RF02) setup with different horizontal and vertical resolutions are conducted in this work. The results are compared with those of the double-moment bulk scheme of Seiki and Nakajima (2014) (hereafter, SN14) and of the model intercomparison project (MIP) of Ackerman et al. (2009). The horizontal and vertical grid lengths ranged from 12.5 to 50 m and 2.5 to 10 m, respectively.

Considering the SDM accuracy, a large number of super-droplets could improve the simulation performance. The computational efficiency should be considered. Accordingly, the numerical convergence regarding the SDNC at different resolutions is discussed to find an optimized initial super-droplet number. Dziekan et al. (2019) suggested that 40 super-droplets per cell is sufficient in achieving the correct domain-averaged results for DYCOMS-II (RF02). However, they did not test smaller SDNCs, and it could be further reduced. Therefore, we choose eight different initial SDNCs ranging from 1 to 128 to investigate the super-droplet numerical convergence in the SDM.

We also compare the difference between a double-moment bulk scheme SN14 and the SDM using the same dynamical core. Some problems in cloud physics (e.g., entrainment-mixing mechanisms) have not been fully understood (Xu et al., 2022; Lu et al., 2023). Microphysics, thermodynamics, and turbulence simulations using high-resolution numerical models can help us understand these mechanisms in the absence of high-resolution observational instruments. Considering the abovementioned advantages of the SDM, our SDM results can provide a reference for the model setting of further studies on stratocumulus and improve our understanding of its macro- and microscopic properties. The time evolution of the aerosol number concentration and the size distributions through the aerosol–cloud interaction cannot be calculated by bulk models. This is difficult even when using bin models but can be accurately represented in particle-based models. Hence, this study on the aerosol–cloud interaction of stratocumulus clouds using particle-based models is important.

Section 2 introduces the basic information of the DYCOMS-II (RF02) simulation setup. Section 3.1 and 3.2 present the results of the SDM grid and the super-droplet number convergence, respectively. Section 3.3 shows the SN14 grid convergence results. Section 3.4 summarizes the SDM and SN14 numerical convergence characteristics. The SN14 and SDM results are compared with each other and with the DYCOMS MIP. Section 4 presents several sensitivity experiments conducted to investigate the mechanisms responsible for the differences. Section 5 summarizes the study findings and points out the shortcomings of our study and the future perspectives related to the numerical simulation of stratocumulus clouds and the aerosol–cloud interaction mechanisms.

2 Method

2.1 Model description

The numerical model used here was the Scalable Computing for Advanced Library and Environment (SCALE;, last access: 19 June 2024), which is a basic library for weather and climate models of the Earth and other planets (Nishizawa et al., 2015; Sato et al., 2015). For the cloud microphysics, the SDM (Shima et al., 2009) and the double-moment bulk scheme SN14 (Seiki and Nakajima, 2014) were used. We implemented SDM into SCALE version 5.2.6, so the model used in this study is referred to as SCALE-SDM 5.2.6-2.3.1.

Moist air fluid dynamics were solved by SCALE's dynamical core. We utilized a forward temporal integration scheme to solve the compressible Navier–Stokes equations for moist air using a finite volume method. The spatial discretization of Eulerian variables was performed on the Arakawa-C staggered grid (Arakawa and Lamb, 1977). The fourth-order central difference scheme was used for the dynamical variable advection. The third-order upwind scheme with Koren (1993) was utilized for the tracer advection. The second-order central difference scheme was employed for other spatial derivatives. We used the four-step Runge–Kutta scheme for the time integration of the dynamical variables and the Wicker and Skamarock (2002) three-step Runge–Kutta scheme for the time integration of tracers. For the SGS turbulence of moist air, unless otherwise stated, a Smagorinsky–Lilly-type scheme, including stratification effects (Lilly, 1962; Smagorinsky, 1963), was used. We added a fourth-order hyper-diffusion to stabilize the calculation. The nondimensional coefficient of the hyper-diffusion term defined in Eq. (A132) of Nishizawa et al. (2015) was set to 10−4.

In the SDM, the time evolution of the aerosol/cloud/precipitation particles is explicitly calculated by solving the elementary process equations of cloud microphysics. In this study, the considered cloud microphysics processes were advection and sedimentation; evaporation and condensation, including cloud condensation nuclei activation and deactivation; and collision–coalescence. To solve the condensation–evaporation process, the implicit Euler scheme was used to avoid stiffness. The Monte Carlo algorithm of Shima et al. (2009) was used for the collision–coalescence process. We employed the uniform sampling method to initialize the super-droplets (Sect. 5.3 of Shima et al., 2020). The SGS turbulence was not considered for the super-droplets. Please see the works of Shima et al. (2009) and Shima et al. (2020) for more details on the governing equations and numerical schemes.

SN14 is a double-moment bulk scheme in which the mixing ratio and the number concentration of the cloud and rain droplets are predicted in each grid, but not the aerosol number concentration. In this work, the microphysical processes of activation–deactivation, condensation–evaporation, and collision–coalescence were calculated at each time step. An aerosol nucleation scheme that estimates cloud condensation nuclei (CCN) activation based on traditional empirical formulas was adopted. This approach, following Twomey (1959), Rogers and Yau (1989), and Seifert and Beheng (2006), accounts for the supersaturation ratio and its influence on aerosol activity, with the maximum activated aerosol number concentration set at 1.5 times the CCN at a supersaturation ratio of 1 %. The scheme incorporates the effects of turbulence on nucleation as per Lohmann (2002), taking into account the effective vertical velocity and sub-grid turbulent kinetic energy. For the condensation process, an explicit condensation scheme rather than the saturation adjustment method was used. The collection processes were similar to those used by Seifert and Beheng (2001), Seifert and Beheng (2006), and Seifert (2008). The SGS turbulence affected the tracers in this bulk scheme. The specific calculation methods of the abovementioned processes were described in detail in the paper of Seiki and Nakajima (2014).

2.2 Numerical setup

We performed simulations of a drizzling marine stratocumulus case observed by the second research flight of DYCOMS-II (RF02) on 11 July 2001 off the coast of Southern California. This field campaign aimed to improve understanding on the stratocumulus characteristics (Stevens et al., 2003).

Table 1Setup of the sensitivity experiments.

a DT, DT_DYN, DT_PHY_SF, DT_PHY_TB, DT_PHY_MP, and DT_PHY_RD are the time steps of time integration and dynamical, surface, turbulence, microphysics, and radiation processes, respectively. b DT_cnd, DT_coa, and DT_adv are the time steps of the condensation, coalescence, and advection processes, respectively.

Download Print Version | Download XLSX

The initial vertical profiles of the wind, moisture air, and temperature followed that of Ackerman et al. (2009). The setup was for the model intercomparison project (hereafter, DYCOMS MIP) based on DYCOMS-II (RF02). DYCOMS MIP contained 14 different LES models with bulk or bin microphysics but no LCM. The domain area was 6 km× 6 km× 1.5 km. The simulation time was 6 h. However, due to constraints in computational resources, the simulation with a finer grid resolution of 12.5 m× 12.5 m× 2.5 m conducted using the SDM was limited to 5 h. The periodic boundary condition was imposed for the lateral boundaries. The simplified radiation model described in Ackerman et al. (2009) was used. Unlike in the work of Ackerman et al. (2009), the maximum supersaturation limited to 1 % in the first hour was used herein not only for droplet activation, but also for condensational growth. The initial aerosol number and the size distributions used for the SDM were the bimodal lognormal distribution specified in Ackerman et al. (2009). We reduced the initial aerosol number concentration from 100 to 70 cm−3 for the SN14 simulations to make the mean cloud droplet number concentration consistent with that of the DYCOMS MIP ( 55 cm−3). Constant latent and sensible heat from the surface was imposed. We also slightly decreased the constant surface latent heat flux from 93 to 86.7132 W m−2 to slightly reduce the predicted liquid water in the SN14 and SDM simulations, thereby avoiding overestimation. The momentum exchange between the super-droplets and the fluid was also considered in the SDM. The SN14 and SDM results were saved every minute.

Table 1 summarizes the specific horizontal and vertical grid lengths, time steps, and SDNC.

  1. Grid resolution test (experiment groups A and B).

    Nine different grid resolution settings were used: DX (= DY) × DZ = 50 m× 2.5 m, 50 m× 5 m, 50 m× 10 m, 25 m× 2.5 m, 25 m× 5 m, 25 m× 10 m, 12.5 m× 2.5 m, 12.5 m× 5 m, and 12.5 m× 10 m. The aspect ratio (DX/DZ) was also considered. All grid cells in the SDM and SN14 runs were uniform and not stretched in space. The SDM and SN14 runs with all different grid spacings were categorized into groups A and B, respectively. In both groups, the runs with DX (= DY) × DZ = 50 m× 5 m were the benchmark runs also used in the DYCOMS MIP. The initial super-droplet number concentrations of the runs in Group A were all 64 per cell. To stabilize the numerical simulations, the time steps were reduced as the resolution became finer. The goal of groups A and B was to explore the grid convergence of SDM and SN14. By comparing these two groups, we expect to find the differences between SDM and SN14.

  2. Initial super-droplet number test (experiment groups C and D).

    Eight initial super-droplet number concentrations were set in Group C from smallest to largest: 1, 2, 4, 8, 16, 32, and 128 per cell. The same grid resolutions of the SDM runs in Group C (i.e., 50 m× 5 m) were compared with base run A8. We also investigated the impact of the grid resolution on the super-droplet number characteristics using Group D, which comprised a series of SDM simulations with super-droplet numbers 1, 4, 16, and 64 per cell at a 25 m× 2.5 m resolution. We expected groups C and D to help us understand the super-droplet convergence characteristic of the cloud droplet number concentration (CDNC).

Figure 1Time series of the (a) liquid water path (LWP), (b) cloud cover, (c) inversion height, (d) entrainment rate, (e) vertically integrated turbulent kinetic energy (TKE), (f) cloud droplet number concentration (CDNC), and (g) surface precipitation for the Group A runs (SDM). The ensemble range, interquartile range, and mean of the DYCOMS MIP results are denoted by the light and dark shading and solid lines, respectively. The ensemble mean from the simulations that included drizzle without sedimentation is denoted by the dashed lines (Ackerman et al., 2009).


3 Numerical convergence characteristics

In this section, we will compare the SCALE results with the DYCOMS MIP results. As specified in the work of Ackerman et al. (2009), the first 2 h was considered as the spin-up period. Moreover, the vertical profiles were averaged over the last 4 h. In all profiles, the y axis is defined as the height normalized by the inversion height zi, which is the mean height of the qt= 8 g kg−1 isosurface. The entrainment rate in the simulations was calculated as E=dzidt+Dzi, where D= 3.75 × 10−6s−1 is the uniform divergence of the large-scale horizontal winds. In the time series and vertical profiles, the ensemble range, interquartile range, and mean of the DYCOMS MIP results are denoted by the light and dark shading and solid lines, respectively. The ensemble mean from the simulations that included drizzle without sedimentation is denoted by the dashed lines (Ackerman et al., 2009).

Figure 2Evolution of the time average of the (a, b) LWP, (c, d) cloud cover, (e, f) inversion height, (g, h) entrainment rate, (i, j) vertically integrated TKE, (k, l) CDNC, and (m, n) surface precipitation for the Group A runs (SDM) with the grid resolution. The left and right columns represent the evolution of the variables with DZ and DX, respectively. Each point in these scatter plots represents the average of one variable from one SDM run. The error bars show the standard deviation of the detrended data.


3.1 SDM grid convergence

In this section, we will analyze the SDM results conducted in various grid resolutions (Group A) to assess the grid convergence characteristics.

3.1.1 Time series of the domain average

Figure 1 shows the time series of several domain-averaged quantities for the SDM results. The results were compared with the DYCOMS MIP results. Figure 2 depicts the statistics of the boundary layer and the cloud-related fields during the last 4 h (the last 3 h for run A1) versus the grid resolutions. The left and right columns represent the change of variables with DZ and DX, respectively. Each point in these plots represents a 4 h average of the variables for the corresponding SDM runs (a 3 h average for run A1). The error bars show the standard deviation of the time series.

The domain-averaged LWP decreased as DZ decreased when DZ  5 m (Figs. 1a and 2a) and was not very sensitive to DX (Figs. 1a and 2b). The LWP showed the trend of getting closer to the true solution as the grid resolution was being refined; however, the LWP changed rate in terms of DZ remained large, even when DZ = 2.5 m. This indicated that a DZ smaller than 2.5 m (DZ < 2.5 m) was needed to obtain a well-converged solution. We conclude from the LWP time series that DX less than or equal to 50 m (DX  50 m) was sufficient; however, in the subsequent paragraphs, this conclusion will be proven untrue for all fields. Unlike the DYCOMS MIP, our LWP increased with time in most of our simulations. The LWP values during the last 3 h in our simulations were all larger than the MIP average.

The cloud cover (CC) is the fraction of cloudy columns defined as columns with an LWP larger than 20 g m−2. Conversely, cloud holes are columns with LWP  20 g m−2. Figs. 1b and 2c showed that CC increased as DZ decreased. Its dependency on DX was relatively weak. CC exhibited the trend of getting closer to the true solution as the grid resolution was being refined but was still sensitive to DZ and weakly sensitive to DX, even when (DZ, DX) = (2.5 m, 12.5 m). This may indicate that DZ smaller than 2.5 m (DZ < 2.5 m) and DX smaller than 12.5 m (DX < 12.5 m) are necessary in achieving a converged solution. Our CC results were almost always higher than the MIP ensemble mean, except when DZ = 10 m. CC rapidly declined and deviated from the MIP results when DZ = 10 m.

The inversion height zi decreased as DZ decreased (Figs. 1c and 2e) and increased as DX decreased except run A1 (Figs. 1c and 2f). zi displayed the trend of also getting closer to the true solution as the grid resolution was being refined, but remaining strongly sensitive to DZ and weakly sensitive to DX, even when (DZ, DX) = (2.5 m, 12.5 m). In other words, DZ < 2.5 m and DX < 12.5 m are necessary in realizing a converged solution. The differences of zi among the Group A runs relative to zi were not big and were less than a few percent. zi of the SDM rapidly increased during the first 2 h and then flattened, whereas that of the MIP more rapidly and continuously increased. This suggests that the DYCOMS MIP has more instability near the cloud tops, leading to stronger upward wind that promotes cloud top growth, while the SDM has more stable cloud tops.

The entrainment rate (Fig. 1d) slowly decreased in time after the first hour then leveled off. By the definition (Eq. 1) adopted from the DYCOMS MIP, the entrainment rate was determined by the inversion height zi and its time derivative dzi/dt. Consequently, its dependency to the grid resolution was similar to that of the inversion height zi. The entrainment rate decreased as DZ decreased (Fig. 2g) and increased as DX decreased (Fig. 2h). It remained strongly sensitive to DZ and weakly sensitive to DX, indicating that DZ < 2.5 m and DX < 12.5 m are necessary for an accurate simulation. The entrainment rates of the SDM were positioned around the lower end of the DYCOMS MIP range.

The vertically integrated total turbulent kinetic energy (TKE), including the resolved TKE and the SGS TKE, rose rapidly during the first 40 min and then fell sharply and rose slowly after the first hour (Fig. 1e). It was insensitive to DX (Figs. 1e and 2j) and DZ (Figs. 1e and 2i). Therefore, DZ < 2.5 m and DX  50 m are necessary for an accurate simulation. For all the SDM runs, the vertically integrated TKE was almost always around the MIP ensemble average.

The CDNC (Fig. 1f) rapidly decreased in the first hour, with an average value of approximately 60 cm−3, which was slightly higher than the MIP ensemble average ( 55 cm−3). It was sensitive to DZ (Figs. 1f and 2k), but its dependence on it is unknown and less sensitive to DX (Figs. 1f and 2l). Therefore, DZ < 2.5 m and DX  50 m are necessary for an accurate simulation.

The surface precipitation (Fig. 1g) in all our simulations was much lower than that in the DYCOMS MIP. Although our SDM results greatly differed from bulk and bin microphysics of the DYCOMS MIP, they were consistent with those in the previous SDM study on this case by Dziekan et al. (2019). Figure 2m illustrates the surface precipitation increase with the decreasing DZ. However, its dependency on DX was not clear, whereas run A7 has the heaviest surface precipitation of all simulations (Fig. 2n).

Figure 3Vertical profiles of the (a) liquid water mixing ratio, (b) rain water mixing ratio, (c) total water mixing ratio, (d) liquid water potential temperature, (e) cloud fraction, (f) buoyancy production, (g) total TKE, (h) w variance (variance of the vertical wind speed), and (i) CDNC for the Group A runs (SDM). The meanings of light and dark shading, solid lines, and dashed lines are the same as in Fig. 1.


3.1.2 Horizontally averaged vertical profile

In addition to the time series in Figs. 1 and 2, we further investigated the vertical profiles to examine the grid convergence and the vertical structure of clouds. Figure 3 depicts the vertical profiles obtained by the horizontal average during the last 4 h (a 3 h average for run A1). The vertical axis in these plots represents the real height z scaled by the inversion height zi.

The liquid water potential temperature θl (Fig. 3d) and the total water mixing ratio qt (Fig. 3c) were not sensitive to the grid resolution. Consistent with the surface precipitation time series (Fig. 3g), the rain water mixing ratio qr profile (Fig. 3b) was almost 0.

The liquid water mixing ratio ql (Fig. 3a) increased as DZ decreased and was not sensitive to DX. It remained strongly sensitive to DZ; hence, (DZ  2.5 m DX  50 m) is necessary for an accurate simulation.

Figure 4Horizontal LWP distribution at 200 min in Series A (SDM). The size of the cloud holes decreased when the resolution became finer.

Figure 3e showed the cloud fraction (CF; fraction of the cloudy grid cells defined as the grid cells with CDNC > 20 cm−3). CF in the lower part of the cloud deck depicted the same grid resolution dependency as ql. The lower-part CF increased as DZ decreased. When DZ = 2.5 m, the lower-part CF increased with decreasing DX but was strongly sensitive to DZ. In other words, (DZ  2.5 m DX < 12.5 m) is needed for the lower-part CF. Conversely, when examining the CF profile around its peak within the cloud deck's midsection, the characteristics of grid convergence were ambiguous. The peak CF values were very similar across simulations with DZ  5 m. In short, (DZ  2.5 m DX < 12.5 m) is necessary for the accurate simulation of the CF. CF in all SDM runs was smaller in the upper part of the cloud deck but larger in the lower part of the cloud deck than that in the DYCOMS MIP. Figure 4 showed the horizontal LWP distribution at the end of the simulation (t= 200 min). Some cloud holes (areas with very low LWP) could be found in Group A. These cloud holes shrank as the grid resolution increased. Figure 3i presented the averaged CDNC profiles exclusively within cloudy cells, revealing an inverse correlation where the CDNC escalates with a decrease in DZ.

The TKE profiles (Fig. 3g) and the variance of the vertical velocity profiles (Fig. 3h) were almost within the ensemble range of the DYCOMS MIP in the cloud deck. Their grid dependency was similar to that of the vertically integrated TKE time series. The TKE profile increased as DZ decreased. It was less sensitive to DX, but stayed strongly sensitive to DZ. Hence, (DZ  2.5 m DX  50 m) is necessary for an accurate simulation.

The TKE buoyancy production profile (Fig. 3f) was relatively insensitive to the grid resolution in the upper part of the cloud layer. We could find a similar dependency to the TKE profiles (Fig. 3g) in the region around the cloud base (i.e., lower part of the cloud and sub-cloud layers). It increased as DZ decreased and was insensitive to DX. Considering the grid dependencies of the two regions, (DZ  2.5 m DX  50 m) is necessary for an accurate simulation. Note also that in the SDM, the buoyancy production below the cloud base is much bigger than that in the DYCOMS MIP.

3.1.3 Summary of the SDM grid convergence and interpretation

Based on the Group A results presented in Figs. 14, the results were qualitatively comparable with each other and got closer to the true solution as the grid resolution was being refined. However, the (DZ, DX) = (2.5 m, 12.5 m) resolution was not high enough yet. In particular, the result was still strongly sensitive to DZ.

Putting everything together, a much finer vertical resolution (DZ  2.5 m) is necessary for almost all quantities. A finer horizontal resolution (DX < 12.5 m) is necessary for CC, CF, zi, and entrainment rate.

We interpret the DZ dependency here. A more detailed cloud structure can be resolved when DZ is refined. This results in a higher CC (Figs. 1b and 4) and an enhanced cloud top cooling. The cooler boundary layer confirmed from the enlarged θl profile (not shown) led to a higher CF (Fig. 3e). Focusing on the lower part of the cloud layer (just above the cloud base) in Fig. 3e, we can observe a relatively large discrepancy of the CF between runs with DZ = 2.5 m and those with DZ = 10 m, explaining the difference of the ql profiles (Fig. 3a). The larger cloudy area near the cloud base resulted in a stronger TKE buoyancy production in the higher-resolution runs (Fig. 3f). In the marine stratocumulus case, the contribution of the buoyancy production dominated the TKE, consequently increasing the TKE (Fig. 3g) and developing the cloud structure that can be confirmed by the higher LWP in the time series (Fig. 1a) and thicker cloud layer (Fig. 3e).

3.2 Super-droplet number convergence

The super-droplet number convergence is discussed in this section. The SDM results at 50 m× 5 m with different initial super-droplet numbers ranging from 1 to 128 super-droplets per cell were compared (Group C). A similar comparison with a finer grid resolution of 25 m× 2.5 m (Group D) was also performed to determine the impact of the grid size on the super-droplet number convergence characteristics.

Figure 5Evolution of the time average of the (a) LWP, (b) cloud cover, (c) inversion height, (d) entrainment rate, (e) vertically integrated TKE, (f) CDNC, and (g) surface precipitation for the Group C and D runs and runs A8 and A4 with the SDNC. Each point in these scatter plots represents the average of one of those variables from one SDM run. Group C and D are marked in by a circle and x, respectively. The error bars show the standard deviation of the detrended data.


Figure 6The time series of (a) CDNC, (b) mean supersaturation within cloudy cells, and (c) maximum supersaturation for the SDM simulations with diverse initial super-droplet (SD) concentrations (Group C), conducted at a grid resolution of 50 m× 50 m× 5 m. Notably, the inset in panel (c) accentuates the first 20 min, underlining the inversely proportional relationship observed between initial super-droplet number concentrations and the magnitude of maximum supersaturation.


Figure 5 depicts the last 4 h average of the variables versus the initial super-droplet number concentration. In Group C, CDNC and precipitation decreased, and the entrainment rate and the inversion height increased with the increasing super-droplet number. CDNC converged to approximately 60cm−3. All the SDM results converged well in addition to the LWP and zi when the initial SDNC was greater than or equal to 16 per cell. The absolute relative errors between C5 (16 per cell) and C7 (128 per cell) in CC, zi, TKE, and CDNC were smaller than 5 %. We proposed herein an explanation for the CDNC decrease with the increasing SDNC (Figs. 5f and 6a). In a case with a quite small initial SDNC (e.g., 1 or 2 super-droplets per cell on average), the multiplicity of super-droplets may therefore be very high, and the phase relaxation time may be very long in some grids with almost no super-droplets and become extremely short in grids with relatively many super-droplets. Consequently, there is a greater potential for greater supersaturation, and thus more aerosols are activated to the cloud droplets. The time evolution of the supersaturation supports our point: the maximum supersaturation (Fig. 6c) significantly decreases as the SDNC increases during the beginning of the simulations. When the SDNC is small (e.g., 1 or 2), the average supersaturation in the cloudy grids (Fig. 6b) is also essentially at a relatively high value.

The variables in Group D (higher grid resolution) showed the same trend with the SDNC as those in Group C (lower grid resolution). The results in the high-resolution runs also converged well when CDNC  16 per cell. The absolute relative errors between D5 (16 per cell) and A4 (64 per cell) in the LWP, CC, zi, entrainment rate, CDNC, and surface precipitation were smaller than 5 %.

It is worth noting that for cases with large amounts of precipitation formation, a low super-droplet number per grid may not be sufficient. This may affect the precipitation formation rate and the spatial distribution of rain and cloud water. Similarly, for polluted cases with GCCN (giant CCN), a sufficient number of super-droplets may be needed to properly sample the aerosol size spectrum and capture the effect of GCCN on precipitation initiation. Dziekan et al. (2021) found that the addition of GCCN to the SDM simulation can significantly increase the precipitation of stratocumulus.

In general, the SDM results all converged well at the initial SDNC greater than or equal to 16 per cell. All, except precipitation, were within the ensemble range of the DYCOMS MIP. The comparison study on different grid resolution indicated that the super-droplet number per grid cell, not per unit volume, is essential for the super-droplet number convergence characteristic. Considering the balance of the computational cost and the simulation accuracy, the SDNC of 16 per cell is the optimal choice of the SDM for the stratocumulus simulations.

Figure 7Time series of the (a) LWP, (b) cloud cover, (c) inversion height, (d) entrainment rate, (e) vertically integrated TKE, (f) CDNC, and (g) surface precipitation for the Group B runs (SN14). The meanings of light and dark shading, solid lines, and dashed lines are the same as in Fig. 1.


Figure 8Evolution of the time average of the (a, b) LWP, (c, d) cloud cover, (e, f) inversion height, (g, h) entrainment rate, (i, j) vertically integrated TKE, and (k, l) CDNC and (m, n) surface precipitation for Group B runs (SN14) with the grid resolution. The left and right columns represent the evolution of the variables with DZ and DX, respectively. Each point in these scatter plots represents the average of one of the variables from one SN14 run. The error bars show the standard deviation of the detrended data.


Figure 9Vertical profiles of the (a) liquid water mixing ratio, (b) rain water mixing ratio, (c) total water mixing ratio, (d) liquid water potential temperature, (e) cloud fraction, (f) buoyancy production, (g) total TKE, (h) w variance, and (i) CDNC for the Group B runs (SN14). The meanings of light and dark shading, solid lines, and dashed lines are the same as in Fig. 1.


3.3 SN14 grid convergence

Similar to Sect. 3.1 for the SDM, this section discusses the conducted SN14 simulations in various grid resolutions (Group B) and the investigated grid convergence characteristics. The time series (Figs. 7 and 8) and the vertical profiles (Fig. 9) of Group B are shown.

The LWP (Figs. 7a and 8a), ql (Fig. 9a), and lower-part CF (Fig. 9e) depicted a similar dependency on the grid spacings; that is, they increased as DZ decreased. The dependency on DX was unclear (Fig. 8b). B1 (12.5 m× 2.5 m) and B4 (25 m× 2.5 m) were indistinguishable. In other words, (DZ  2.5 m DX  25 m) is necessary for these variables. The cloud deck of SN14 was thicker than that of the SDM. Consequently, it contained more liquid water. SN14 agreed with the MIP better than the SDM.

Figure 10Horizontal LWP distribution at 200 min in Series B (SN14).

Accordingly, CC (Figs. 7b and 8c) and maximum CF (Fig. 9e) showed similar dependencies to the grid spacings; that is, they increased as DZ decreased (Fig. 8d) but were not sensitive to the grid spacings if (DZ  5 m DX  50 m). CC and maximum CF of SN14 were nearly one, a reflection of minimal cloud holes as evidenced by the horizontal distribution of LWP observed in Figs. 4 and 10. This near-unity in CC and peak CF is indicative of a highly continuous cloud field.

Similar to the result in the SDM, the inversion height zi (Figs. 7c and 8e) decreased as DZ decreased (Fig. 8f). The DX dependency was unclear, but for all DZ, DX = 25 m and DX = 12.5 m agreed well, indicating that (DZ  2.5 m DX  25 m) is necessary for an accurate simulation. The entrainment rate (Figs. 7d and 8g, h) was determined by zi and dzi/dt. zi and the entrainment rate of SN14 were larger than those of the SDM due to the fast increase in zi during the spin-up time period. zi of SN14 was located around the upper end of the MIP, but the entrainment rate of SN14 during the last 4 h agreed well with the MIP result.

The TKE (Figs. 7e, 8i, and 9g), w variance (Fig. 9h), buoyancy production around the cloud base (Fig. 9f), and CDNC (Figs. 7f, 8k, and 9i) showed similar dependencies; that is, they increased as DZ decreased and decreased as DX decreased. In other words, (DZ  2.5 m DX  25 m) is necessary for an accurate simulation. No clear trend was observed in the grid dependency of the buoyancy production near the cloud top, but B1 (12.5 m× 2.5 m) and B4 (25 m× 2.5 m) were almost indistinguishable.

All variables characterizing turbulence (i.e., TKE, w variance, and buoyancy production) were larger in SN14 than in the SDM. In particular, the buoyancy production in the cloud layer was noticeably higher in SN14 than in the SDM. SN14 agreed better with the MIP.

The surface precipitation and qr were so small in all the Group B results that they did not affect the overall grid convergence characteristics. This was in agreement with the SDM results but much lower than the ensemble mean of the MIP.

SN14 showed a higher CDNC than the SDM in the first hour and a lower CDNC during the last 5 h (Figs. 1f and 7f). The higher peak of the CDNC in SN14 might be caused by the higher maximum supersaturation during the first hour. In that of Ackerman et al. (2009) and our SDM simulation, the maximum supersaturation was limited to 1 % during the convection spin-up to avoid precipitation suppression, which was not adapted in our SN14 simulation. We conducted an SN14 sensitivity test with the supersaturation limiter. The result (not shown) showed a lower CDNC during the first hour (maximum reduced from 120 to 90 cm−3), but it had little effect on the CDNC after the spin-up stage.

From Figs. 79, we conclude that the LWP, ql, CC, CF, TKE, w variance, and buoyancy production around the cloud base and CDNC increased, and the inversion height and the entrainment rate decreased with the decreasing DZ. However, the grid dependency on DX was relatively weak. The sensitivity of the variables to the grid resolution in SN14 was very similar to that in the SDM, but the SDM variables showed a stronger grid dependency on DZ and DX. In conclusion, the variation trend of the variables with the grid resolution in SN14 was more ambiguous than that in the SDM.

In summary, a much finer vertical resolution (DZ  2.5 m) was necessary for almost all quantities, except for CC and maximum CF. The horizontal resolution of DX = 25 m was sufficient for all quantities. The SN14 results agreed with the MIP better than the SDM results.

3.4 Summary of the numerical convergence characteristics

Section 3.1 and 3.3 revealed the grid convergence characteristics of SDM and SN14, respectively. The finest grid resolution tested was (DZ, DX) = (2.5 m, 12.5 m). Our analysis revealed that the grid convergence of both schemes has not yet been achieved. However, we observed a trend where the results approached the true solution. Overall, under the tested parameter range, both SDM and SN14 were strongly sensitive to the vertical resolution and relatively weakly sensitive to the horizontal resolution.

In conclusion, DZ  2.5 m was necessary for both schemes. Note, however, that the TKE and CDNC in the SDM and the CC and the maximum CF in SN14 were not any more sensitive to DZ.

In the SDM, a finer horizontal resolution (DX < 12.5 m) was necessary for CC, maximum CF, zi, and entrainment rate. In contrast, the horizontal resolution of DX  25 m was sufficient for SN14.

Figure 11Time series of the LWP, cloud cover, inversion height, entrainment rate, vertically integrated TKE, CDNC, and surface precipitation for runs A8, B8, and SDM without sedimentation and SDM with horizontal SGS turbulent mixing.


This grid convergence characteristic study showed that when the aspect ratio is unchanged, the LWP increases with the decreasing grid spacing, consistent with the previous studies (Pedersen et al., 2016; Mellado et al., 2018; Matheou and Teixeira, 2019) with isotopic grids (DX/DZ = 1). Some studies on the role of the grid resolution in the numerical simulations of the turbulent entrainment have also shown that the accurate entrainment rate simulation requires a vertical grid spacing no greater than the turbulent undulation scale, which can be 5–10 m for the inversion and turbulence levels typical of the subtropical marine stratocumulus (Stevens and Bretherton, 1999).

However, the stratocumulus simulation is notorious for such a slow grid convergence with respect to the vertical grid spacing DZ. The LES study of Matheou and Teixeira (2019) showed that the numerical convergence for the LWP is hard to achieve, even if the isotopic grid size of 1.25 m is used. However, the LES utilizing a 1.25 m grid resolution can reproduce a detailed cloud structure (e.g., elongated regions of low LWP, cloud holes, and pockets) (Matheou, 2018). Mellado et al. (2018) suggested that 2.5 m or less was necessary for their LES to approach the observation. Furthermore, the LWP was numerically converged when the Kolmogorov scale used in their direct numerical simulation (DNS), which was smaller than 0.7 m.

In Sect. 3.2, we also conducted a numerical convergence analysis on the super-droplet numbers at different grid resolutions. In conclusion, the initial super-droplet number concentration  16 per cell is sufficient for the tested stratocumulus case. The super-droplet number convergence characteristic was essentially determined by the super-droplet number per grid cell and not per unit volume. The LWP, CC, inversion height, entrainment rate, and TKE results also supported the finding that SDNC  16 per cell was good enough for an accurate stratocumulus simulation.

Note also that the SN14 results agreed with the MIP better than the SDM results. In the subsequent sections, we will focus on understanding the difference in the SDM, SN14, and MIP and conduct an in-depth analysis to elucidate the underlying mechanism.

4 Comparison of the SDM and SN14

In Sect. 3, we explored grid convergence characteristics of the SDM and SN14, noting their similar dependencies on grid size yet distinct outcomes in the cloud top height, entrainment rate, liquid water content, CF, CDNC, and turbulence characteristics.

Figure 12Vertical profiles of the (a) liquid water mixing ratio, (b) rain water mixing ratio, (c) total water mixing ratio, (d) liquid water potential temperature, (e) cloud fraction, (f) buoyancy production, (g) total TKE, (h) w variance, and (i) CDNC for the A8, B8, and SDM without sedimentation and SDM with horizontal SGS turbulent mixing.


The SDM simulations consistently yielded higher LWP values that exhibited an upward trend, unlike the SN14 simulations where LWP appeared lower and remained relatively constant past the initial hour (Fig. 11a). Furthermore, the SDM predicted enhanced ql (Fig. 12a) and a denser cloud base, evidenced by an increased cloud fraction near the cloud base, in contrast to the SN14 scheme which indicated a propensity towards more pronounced vertical cloud development, as demonstrated by a larger maximum cloud fraction near the cloud top (Fig. 12e).

The differences between the SDM and SN14 may stem from their distinct aerosol treatment, which caused the difference between SDM and SN14. Temporal changes of aerosol size distribution through wet deposition were explicitly considered in the SDM, while background aerosol particles were assumed to be uniform and unchanged in SN14. Considering the 6 h simulation time, this could modulate the characteristics of the stratocumulus deck.

Figure 13Vertical profiles of ql, qr, qt, θl, buoyancy production, total TKE, Nc, Np, qn, NSD, CF, and Nc for runs A8 and SDM without sedimentation. The solid and dashed lines in panels (a)(j) represent the profiles in the cloudy and hole columns, respectively.


To examine our hypothesis concerning aerosol wet deposition through cloud droplet sedimentation, an SDM simulation without sedimentation (SDM_no_sed) was conducted. In this simulation, all super-droplets are treated as passive tracers when updating their positions; i.e., their terminal velocities are always zero. Note also that their collision–coalescence by differential sedimentation is still considered in SDM_no_sed. The vertical profiles of the total particle number concentration Np, total particle number mixing ratio qn (ratio of Np and air density), and super-droplet number concentration NSD in the SDM simulations reveal local minima near the cloud top, consistent across scenarios with and without sedimentation (illustrated in Fig. 13h, i, and j, respectively). Notably, these minima were more pronounced in the simulations incorporating the sedimentation process, aligning with our conjecture and underscoring the sedimentation's significant impact on particle distribution near the cloud top. This aerosol particle depletion near the cloud top propagated to the hole volume (Fig. 13h, i, and j). Moreover, the aerosol particle depletion in the hole volume in the SDM would make the cloud hole persistent (Fig. 13k). Note also that Fig. 9b of Arabas et al. (2015) also showed the depletion of aerosol particles near the stratocumulus top and the associated hole volume. On the other hand, relatively high concentrations of particles near the cloud top with almost no cloud holes could be found in SDM_no_sed (Fig. 13h, i, j, and k). The above results demonstrated that the wet deposition of aerosols was responsible for the low CF near the cloud top and it also lowered the cloud tops of marine stratocumulus clouds and increased the cloud fraction near the cloud base. In our simulations, we assumed that the aerosol number concentration was initially uniform in space, including the free atmosphere. This should partially compensate for the aerosol particle reduction at the cloud top. In other words, the effect of the cloud top aerosol reduction on the cloud top volume reduction should be greater in the real world. Another Lagrangian cloud model, called UWLCM (Dziekan et al., 2019), presented a similar CF profile to ours in this stratocumulus case. Since all the super-droplets in SDM_no_sed do not sediment, larger droplets within the cloud are more prone to collide and coalesce. Consequently, cloud droplets are more likely to be collected in SDM_no_sed. This leads to an increase in qr (Fig. 13b) and a decrease in CDNC (Fig. 13i).

The underlying reasons for these discrepancies could also be attributed to the inherent modeling capabilities and approaches of the SDM and SN14. The SDM's design to circumvent numerical diffusion issues enables it to avoid the unphysical dispersion of moisture and cloud droplets across space. Consequently, cloud droplets were not artificially diffused into dry air to evaporate; instead, localized areas of higher supersaturation promote more water vapor to condense onto cloud droplets, leading to an increased liquid water content. This enhanced liquid water content further intensified radiative warming at the cloud base and radiative cooling at the cloud top. Radiative warming at the cloud base potentially made the airflows near the cloud base more active, thereby increasing TKE and facilitating more moisture transport into the cloud. As a result, in the SDM simulations, the cloud fraction, liquid water mixing ratio, buoyancy production, and CDNC were relatively higher near the cloud base (Fig. 12). Stronger radiative cooling at the cloud top reduced the temperature in this region (Fig. 12d), enhancing the stability of the cloud top, which in turn diminished cloud top entrainment mixing (Fig. 11d) and TKE (Fig. 11e). The localized subsidence that ensued acts to suppress the vertical development of the cloud layer, evidenced by a reduced cloud top height (Fig. 11c) and buoyancy production at the cloud top (Fig. 12f), leading to a decrease in CF near the cloud top (Fig. 12e).

Although the avoidance of numerical diffusion is hypothesized to contribute to the enhanced liquid water content observed in SDM simulations, the precise impact of numerical diffusion on cloud dynamics necessitates further exploration. Subsequent sensitivity experiments are imperative to rigorously evaluate the extent of numerical diffusion's influence and to corroborate the mechanisms suggested by this study.

In addition to aerosol wet deposition and numerical diffusion effects, the treatment differences in CCN activation–deactivation between the SDM and SN14 may also have influenced the simulation outcomes. The nucleation scheme in SN14 potentially induces more frequent activation–deactivation events compared to the explicit calculations in SDM, possibly affecting latent heat exchange and cloud droplet formation (Hoffmann, 2016; Yang et al., 2023). This aspect of model divergence also requires additional investigation to substantiate its effects on buoyancy flux and cloud structure.

From the simulations presented herein, we could not conclude which mechanism is dominating the phenomenon. A detailed assessment of the proposed scenarios will be conducted in future studies.

5 Conclusions

In this study, we assessed the performance of the Lagrangian particle-based cloud microphysics scheme, called SDM, for the marine stratocumulus simulation. To do this, we conducted a series of numerical simulations based on DYCOMS-II (RF02). For comparison, we also tested a double-moment bulk scheme, called SN14, using the same dynamical core.

Our simulation results were compared with the results of the model intercomparison project DYCOMS-II (RF02) MIP (Ackerman et al., 2009). In general, all the SDM and SN14 variables mostly showed a reasonable agreement with the MIP results.

We also investigated the numerical convergence characteristics of both schemes. We first assessed their grid convergence and confirmed their similar grid size dependencies. The stratocumulus cloud simulation was more sensitive to the vertical resolution than the horizontal resolution. Both SDM and SN14 simulations showed that as DZ decreased, there was an increase in CC, ql, TKE, w variance, CDNC, CF, and buoyancy production near the cloud base. Conversely, the inversion height (zi) and the entrainment rate decreased with decreasing DZ. It is noteworthy that the LWP exhibited distinct sensitivities to DZ in the two schemes: in the SDM, the LWP decreased as DZ diminishes, particularly when DZ  5 m; conversely, in the SN14, the LWP increased with a decrease in DZ. Despite refining the grid resolution, these variables had not converged within the assessed grid spacing range, indicating that a DZ smaller than 2.5 m may be necessary for a well-converged solution. However, CC, CF, LWP, and ql appeared to be relatively insensitive to changes in DX. The previous studies suggested that an isotopic grid size of 2.5 m or less is needed for the LES (Mellado et al., 2018). Conclusively, DZ  2.5 m was necessary for both SDM and SN14, DX < 12.5 m was necessary for the SDM simulations, and DX  25 m was sufficient for SN14. Considering the huge computational resources required, we could not conduct in-depth and finer-resolution simulations to explore the grid size convergence properties of the SDM and SN14.

According to the super-droplet convergence results, the CDNC increased with the decreasing SDNC due to the longer phase relaxation time. The simulations numerically converged when the SDNC was larger than 16 super-droplets per cell, which is smaller than the 40 super-droplets per cell that Dziekan et al. (2019) confirmed. The entrainment rate and the inversion height also increased with the increasing SDNC. The super-droplet convergence study on different grid resolutions indicated that the super-droplet number per grid cell was more essential for the SDM simulation than that per unit volume. Considering the balance of the computational cost and the simulation accuracy, the SDNC of 16 per cell is the optimal choice of the SDM for this marine stratocumulus case.

In conclusion, our results suggest that due to the difficulty of the grid convergence of the cloud liquid, finer resolutions, especially vertical ones, are necessary for stratocumulus simulations using the SDM and the bulk scheme. Accordingly, to improve the computational efficiency, 16 super-droplets per grid cell should be enough for the SDM simulation of non-precipitating marine stratocumulus cases. When computational resources are limited, utilizing stretched grids can be beneficial. As mentioned previously, for simulating stratocumulus clouds, the liquid-water-related variables (e.g., LWP, CC, CF) exhibit a stronger dependence on vertical resolution than horizontal resolution. Therefore, it is possible to allocate a finer vertical resolution in the boundary layer while maintaining a coarser resolution in the free atmosphere to conserve computational resources. Additionally, some studies have highlighted the presence of a very thin inversion structure due to strong radiative cooling at the top of stratocumulus clouds, and turbulent entrainment through this thin layer can exert significant feedback effects on boundary layer and cloud properties (Mellado et al., 2018). Consequently, it is advisable to maintain a fine vertical resolution near the cloud top.

The main discrepancies between the SDM and SN14 simulations of marine stratocumulus are reflected in the content and distribution of liquid water as well as in the macrostructure of the cloud. The SDM simulations consistently yielded higher LWP values that increased over time, contrasting with SN14 simulations, where LWP remained lower and relatively stable past the initial hour. SDM predicted enhanced ql and a denser cloud base, evidenced by an increased CF near the cloud base. Conversely, SN14 tended towards more vertical cloud development, as shown by a larger maximum CF near the cloud top.

The differences observed between the SDM and SN14 simulations can be attributed to their distinct aerosol treatment. In the SDM approach, the dynamic process of aerosol size distribution through wet deposition is actively modeled, allowing for a realistic representation of aerosol behavior over time. Conversely, the SN14 scheme employs a static approach, treating background aerosol particles as constant and uniform throughout the simulation. This fundamental difference in aerosol treatment between the two models has the potential to significantly influence the behavior and structure of stratocumulus clouds over the course of a 6 h simulation period. The total particle number concentration Np, total particle number mixing ratio qn (ratio of Np and air density), and super-droplet number concentration NSD exhibit local minima near the cloud tops in simulations with and without the inclusion of sedimentation processes (SDM and SDM_no_sed, respectively). These minima are accentuated in simulations incorporating sedimentation, highlighting its significant role in altering particle distributions near the cloud top. Furthermore, aerosol depletion within the hole volume in SDM simulations promotes the persistence of these holes, contrasting with the SDM_no_sed scenario, where a higher particle concentration near cloud tops and a virtual absence of the hole volume are observed. This indicates wet deposition's pivotal role in reducing CF near the cloud top, lowering the cloud top, and enhancing CF near the cloud base in marine stratocumulus clouds.

The discrepancies noted could be attributed to the inherent modeling capabilities and methodologies of the SDM and SN14. Specifically, the SDM's strategy to evade numerical diffusion issues allows it to prevent the unphysical spread of moisture and cloud droplets. This characteristic facilitates the formation of localized regions of enhanced supersaturation, promoting the condensation of water vapor into cloud droplets. It concurrently minimizes the dispersion of cloud droplets into unsaturated areas, thereby mitigating evaporation. Consequently, this process contributes to an increase in the liquid water content within the cloud. Such an augmentation of liquid water content significantly amplifies radiative warming at the cloud base and enhances radiative cooling at the cloud top. This differential radiative effect further influences cloud dynamics by activating more vigorous airflows near the cloud base, which enhances TKE and facilitates the upward transport of moisture into the cloud. Additionally, it stabilizes the cloud top by reducing temperature, which in turn moderates cloud top entrainment-mixing processes and TKE, ultimately affecting the vertical development and structural integrity of the cloud. It is important to note that these proposed mechanisms represent initial hypotheses based on observed simulation differences and require further empirical validation through dedicated experimental studies.

While we also explored the potential impact of differences in CCN activation–deactivation treatments between the models, further evidence is required to support this hypothesis.

This study on numerical convergence can help researchers set up precise stratocumulus cloud simulations using the SDM and bulk schemes. Our comparison of the SDM and SN14 also indicates that accurate aerosol representation and its dynamic interaction with cloud processes play crucial roles in shaping marine stratocumulus characteristics as simulated by numerical models. Future studies are warranted for a more detailed examination of these mechanisms and their implications for cloud modeling and prediction.

Code and data availability

The source code of SCALE-SDM is available from (Yin, 2024a). All the model results used for this study can be reproduced by following the instructions included in the above repository. Due to the size of the model results, the data are deposited in local storage at the University of Hyogo in Kobe, Japan, and are available from the corresponding author upon request. DYCOMS MIP data are available on the website of the DYCOMS-II field campaign (, NASA Goddard Institute for Space Studies, 2018).

Video supplement

The video supplement related to this article is available online at (Yin, 2024b).

Author contributions

All the authors designed the experiments and CY carried them out. SS developed the model code and CY performed the simulations. CY prepared the manuscript with contributions from all co-authors.

Competing interests

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


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.


Chongzhi Yin and Shin-ichiro Shima would like to thank Yousuke Sato and Seiya Nishizawa for their generous support and informative discussions. We acknowledge the High-Performance Computing Center of the Nanjing University of Information Science and Technology for their support of this work. This study was supported by the National Key Scientific and Technological Infrastructure project “Earth System Science Numerical Simulator Facility” (EarthLab). This research partly used the computational resources of Kyushu University and Hokkaido University through the HPCI System Research Project (project IDs: hp160132, hp200078, hp210059, hp220062, and hp230166) and the computer facilities of the Center for Cooperative Work on Data science and Computational science, University of Hyogo. The SCALE library was developed by Team SCALE of RIKEN Center for Computational Sciences (, last access: 2 December 2022). The authors would like to thank Enago (, last access: 19 June 2024) for English language review of an earlier version of the manuscript.

Financial support

This research has been supported by the National Natural Science Foundation of China (grant no. 42325503); the Japan Society for the Promotion of Science (KAKENHI grant nos. 26286089, 20H00225, and 23H00149); the Ministry of Education, Culture, Sports, Science and Technology, Japan Society for the Promotion of Science (KAKENHI grant no. 18H04448); and the Japan Science and Technology Agency, Moonshot Research and Development Program (grant no. JPMJMS2286).

Review statement

This paper was edited by Nina Crnivec and reviewed by two anonymous referees.


Abade, G. C., Grabowski, W. W., and Pawlowska, H.: Broadening of cloud droplet spectra through eddy hopping: Turbulent entraining parcel simulations, J. Atmos. Sci., 75, 3365–3379,, 2018. 

Ackerman, A. S., VanZanten, M. C., Stevens, B., Savic-Jovcic, V., Bretherton, C. S., Chlond, A., Golaz, J.-C., Jiang, H., Khairoutdinov, M., and Krueger, S. K.: Large-eddy simulations of a drizzling, stratocumulus-topped marine boundary layer, Mon. Weather Rev., 137, 1083–1110,, 2009. 

Alfonso, L. and Raga, G. B.: The impact of fluctuations and correlations in droplet growth by collision–coalescence revisited – Part 1: Numerical calculation of post-gel droplet size distribution, Atmos. Chem. Phys., 17, 6895–6905,, 2017. 

Arabas, S. and Shima, S.-I.: Large-eddy simulations of trade wind cumuli using particle-based microphysics with Monte Carlo coalescence, J. Atmos. Sci., 70, 2768–2777,, 2013. 

Arabas, S., Jaruga, A., Pawlowska, H., and Grabowski, W. W.: libcloudph++ 1.0: a single-moment bulk, double-moment bulk, and particle-based warm-rain microphysics library in C++, Geosci. Model Dev., 8, 1677–1707,, 2015. 

Arakawa, A. and Lamb, V. R.: Computational design of the basic dynamical processes of the UCLA general circulation model, Methods in Computational Physics: Advances in Research and Applications, 17, 173–265,, 1977. 

Bony, S. and Dufresne, J. L.: Marine boundary layer clouds at the heart of tropical cloud feedback uncertainties in climate models, Geophys. Res. Lett., 32, L20806,, 2005. 

Bony, S., Colman, R., Kattsov, V. M., Allan, R. P., Bretherton, C. S., Dufresne, J.-L., Hall, A., Hallegatte, S., Holland, M. M., and Ingram, W.: How well do we understand and evaluate climate change feedback processes?, J. Climate, 19, 3445–3482,, 2006. 

Boucher, O., Randall, D., Artaxo, P., Bretherton, C., Feingold, G., Forster, P., Kerminen, V.-M., Kondo, Y., Liao, H., and Lohmann, U.: Clouds and aerosols, in: Climate change 2013: the physical science basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA, 571–657,, 2013. 

Chandrakar, K. K., Morrison, H., and Witte, M.: Evolution of Droplet Size Distributions During the Transition of an Ultraclean Stratocumulus Cloud System to Open Cell Structure: An LES Investigation Using Lagrangian Microphysics, Geophys. Res. Lett., 49, e2022GL100511,, 2022. 

Cheng, A., Xu, K. M., and Stevens, B.: Effects of resolution on the simulation of boundary-layer clouds and the partition of kinetic energy to subgrid scales, J. Adv. Model. Earth Sy., 2, 3,, 2010. 

Dziekan, P. and Pawlowska, H.: Stochastic coalescence in Lagrangian cloud microphysics, Atmos. Chem. Phys., 17, 13509–13520,, 2017. 

Dziekan, P., Waruszewski, M., and Pawlowska, H.: University of Warsaw Lagrangian Cloud Model (UWLCM) 1.0: a modern large-eddy simulation tool for warm cloud modeling with Lagrangian microphysics, Geosci. Model Dev., 12, 2587–2606,, 2019. 

Dziekan, P., Jensen, J. B., Grabowski, W. W., and Pawlowska, H.: Impact of Giant Sea Salt Aerosol Particles on Precipitation in Marine Cumuli and Stratocumuli: Lagrangian Cloud Model Simulations. J. Atmos. Sci., 78, 4127–4142,, 2021. 

Geresdi, I., Xue, L., and Rasmussen, R.: Evaluation of orographic cloud seeding using a bin microphysics scheme: Two-dimensional approach, J. Appl. Meteorol. Clim., 56, 1443–1462,, 2017. 

Grabowski, W. W. and Abade, G. C.: Broadening of cloud droplet spectra through eddy hopping: Turbulent adiabatic parcel simulations, J. Atmos. Sci., 74, 1485–1493,, 2017. 

Grabowski, W. W., Dziekan, P., and Pawlowska, H.: Lagrangian condensation microphysics with Twomey CCN activation, Geosci. Model Dev., 11, 103–120,, 2018. 

Grabowski, W. W., Morrison, H., Shima, S.-I., Abade, G. C., Dziekan, P., and Pawlowska, H.: Modeling of cloud microphysics: Can we do better?, B. Am. Meteorol. Soc., 100, 655–672,, 2019. 

Hoffmann, F.: The Effect of Spurious Cloud Edge Supersaturations in Lagrangian Cloud Models: An Analytical and Numerical Study, Mon. Weather Rev., 144, 107–118,, 2016. 

Hoffmann, F. and Feingold, G.: Entrainment and Mixing in Stratocumulus: Effects of a New Explicit Subgrid-Scale Scheme for Large-Eddy Simulations with Particle-Based Microphysics, J. Atmos. Sci., 76, 1955–1973,, 2019. 

Hoffmann, F., Yamaguchi, T., and Feingold, G.: Inhomogeneous mixing in Lagrangian cloud models: Effects on the production of precipitation embryos, J. Atmos. Sci., 76, 113–133,, 2019. 

Jaruga, A. and Pawlowska, H.: libcloudph++ 2.0: aqueous-phase chemistry extension of the particle-based cloud microphysics scheme, Geosci. Model Dev., 11, 3623–3645,, 2018. 

Kawai, H. and Shige, S.: Marine Low Clouds and Their Parameterization in Climate Models, J. Meteorol. Soc. Jpn. Ser. II, 98, 1097–1127,, 2020. 

Kessler, E.: On the distribution and continuity of water substance in atmospheric circulations, in: On the distribution and continuity of water substance in atmospheric circulations, Springer, Boston, MA., 1–84,, 1969. 

Khain, A., Ovtchinnikov, M., Pinsky, M., Pokrovsky, A., and Krugliak, H.: Notes on the state-of-the-art numerical modeling of cloud microphysics, Atmos. Res., 55, 159–224,, 2000. 

Khain, A., Beheng, K., Heymsfield, A., Korolev, A., Krichak, S., Levin, Z., Pinsky, M., Phillips, V., Prabhakaran, T., and Teller, A.: Representation of microphysical processes in cloud-resolving models: Spectral (bin) microphysics versus bulk parameterization, Rev. Geophys., 53, 247–322,, 2015. 

Koren, B.: A robust upwind discretization method for advection, diffusion and source terms, Centrum voor Wiskunde en Informatica, Amsterdam, 1993. 

Kurowski, M. J., P. Malinowski, S., and W. Grabowski, W.: A numerical investigation of entrainment and transport within a stratocumulus-topped boundary layer, Q. J. Roy. Meteor. Soc., 135, 77–92,, 2009. 

Lilly, D. K.: On the numerical simulation of buoyant convection, Tellus, 14, 148–172,, 1962. 

Lin, Y.-L., Farley, R. D., and Orville, H. D.: Bulk parameterization of the snow field in a cloud model, J. Appl. Meteorol. Clim., 22, 1065–1092,<1065:BPOTSF>2.0.CO;2, 1983. 

Lohmann, U.: Possible aerosol effects on ice clouds via contact nuclei, J. Atmos. Sci., 59, 647–656, 2002. 

Lu, C., Zhu, L., Liu, Y., Mei, F., Fast, J. D., Pekour, M. S., Luo, S., Xu, X., He, X., Li, J., and Gao, S.: Observational study of relationships between entrainment rate, homogeneity of mixing, and cloud droplet relative dispersion, Atmos. Res., 293, 106900,, 2023. 

Lynn, B. H., Khain, A. P., Dudhia, J., Rosenfeld, D., Pokrovsky, A., and Seifert, A.: Spectral (bin) microphysics coupled with a mesoscale model (MM5). Part I: Model description and first results, Mon. Weather Rev., 133, 44–58,, 2005. 

Matheou, G.: Turbulence structure in a stratocumulus cloud, Atmosphere, 9, 392,, 2018. 

Matheou, G. and Teixeira, J.: Sensitivity to physical and numerical aspects of large-eddy simulation of stratocumulus, Mon. Weather Rev., 147, 2621–2639,, 2019. 

Matheou, G., Chung, D., and Teixeira, J.: On the synergy between numerics and subgrid scale modeling in LES of stratified flows: Grid convergence of a stratocumulus-topped boundary layer, International Symposium on Stratified Flows, San Diego, CA, 31 August 2016, (last access: 19 June 2024), 2016. 

Mellado, J.-P., Bretherton, C., Stevens, B., and Wyant, M.: DNS and LES for simulating stratocumulus: Better together, J. Adv. Model. Earth Sy., 10, 1421–1438,, 2018. 

Milbrandt, J. and Yau, M.: A multimoment bulk microphysics parameterization. Part I: Analysis of the role of the spectral shape parameter, J. Atmos. Sci., 62, 3051–3064,, 2005. 

Morrison, H. and Grabowski, W. W.: An improved representation of rimed snow and conversion to graupel in a multicomponent bin microphysics scheme, J. Atmos. Sci., 67, 1337–1360,, 2010. 

Morrison, H., Witte, M., Bryan, G. H., Harrington, J. Y., and Lebo, Z. J.: Broadening of modeled cloud droplet spectra using bin microphysics in an Eulerian spatial domain, J. Atmos. Sci., 75, 4005–4030,, 2018. 

Morrison, H., van Lier-Walqui, M., Fridlind, A. M., Grabowski, W. W., Harrington, J. Y., Hoose, C., Korolev, A., Kumjian, M. R., Milbrandt, J. A., Pawlowska, H., Posselt, D. J., Prat, O. P., Reimel, K. J., Shima, S.-I., van Diedenhoven, B., and Xue, L.: Confronting the Challenge of Modeling Cloud and Precipitation Microphysics, J. Adv. Model. Earth Sy., 12, e2019MS001689,, 2020. 

NASA Goddard Institute for Space Studies: GCSS LES Intercomparison Based on Flight RF02, in: GCSS-DIME, DYCOMS-II Field Campaign, NASA Goddard Institute for Space Studies [data set], (last access: 19 June 2024), 2018. 

Naumann, A. K. and Seifert, A.: A Lagrangian drop model to study warm rain microphysical processes in shallow cumulus, J. Adv. Model. Earth Sy., 7, 1136–1154,, 2015. 

Nishizawa, S., Yashiro, H., Sato, Y., Miyamoto, Y., and Tomita, H.: Influence of grid aspect ratio on planetary boundary layer turbulence in large-eddy simulations, Geosci. Model Dev., 8, 3393–3419,, 2015. 

Noh, Y., Oh, D., Hoffmann, F., and Raasch, S.: A Cloud Microphysics Parameterization for Shallow Cumulus Clouds Based on Lagrangian Cloud Model Simulations, J. Atmos. Sci., 75, 4031–4047,, 2018. 

Nowak, J. L., Siebert, H., Szodry, K.-E., and Malinowski, S. P.: Coupled and decoupled stratocumulus-topped boundary layers: turbulence properties, Atmos. Chem. Phys., 21, 10965–10991,, 2021. 

Pedersen, J. G., Malinowski, S. P., and Grabowski, W. W.: Resolution and domain-size sensitivity in implicit large-eddy simulation of the stratocumulus-topped boundary layer, J. Adv. Model. Earth Sy., 8, 885–903, 2016. 

Richter, D. H., MacMillan, T., and Wainwright, C.: A Lagrangian Cloud Model for the Study of Marine Fog, Bound.-Lay. Meteorol., 181, 523–542,, 2021. 

Rogers, R. R. and Yau, M. K.: A Short Course in Cloud Physics, 3rd edn., Butterworth-Heinemann, Burlington, MA, ISBN 9780750632157, 1989. 

Sato, Y., Nishizawa, S., Yashiro, H., Miyamoto, Y., Kajikawa, Y., and Tomita, H.: Impacts of cloud microphysics on trade wind cumulus: which cloud microphysics processes contribute to the diversity in a large eddy simulation?, Progress in Earth and Planetary Science, 2, 1–16,, 2015. 

Sato, Y., Shima, S. I., and Tomita, H.: A grid refinement study of trade wind cumuli simulated by a Lagrangian cloud microphysical model: the super-droplet method, Atmos. Sci. Lett., 18, 359–365,, 2017. 

Sato, Y., Shima, S. I., and Tomita, H.: Numerical Convergence of Shallow Convection Cloud Field Simulations: Comparison Between Double-Moment Eulerian and Particle-Based Lagrangian Microphysics Coupled to the Same Dynamical Core, J. Adv. Model. Earth Sy., 10, 1495–1512,, 2018. 

Schoeffler, P.: Dissipation, dispersion and stability of numerical schemes for advection and diffusion, SMHI, Report No. 34, Swedish Meteorological and Hydrological Institute, Norrköping, Sweden, ISSN 0347-2116., 1982. 

Schoenberg Ferrier, B.: A double-moment multiple-phase four-class bulk ice scheme. Part I: Description, J. Atmos. Sci., 51, 249–280,;2, 1994. 

Schwenkel, J., Hoffmann, F., and Raasch, S.: Improving collisional growth in Lagrangian cloud models: development and verification of a new splitting algorithm, Geosci. Model Dev., 11, 3929–3944,, 2018. 

Seifert, A.: On the parameterization of evaporation of raindrops as simulated by a one-dimensional rainshaft model, J. Atmos. Sci., 65, 3608–3619,, 2008. 

Seifert, A. and Beheng, K. D.: A double-moment parameterization for simulating autoconversion, accretion and selfcollection, Atmos. Res., 59, 265–281, 2001. 

Seifert, A. and Beheng, K. D.: A two-moment cloud microphysics parameterization for mixed-phase clouds. Part 1: Model description, Meteorol. Atmos. Phys., 92, 45–66,, 2006. 

Seifert, A. and Rasp, S.: Potential and Limitations of Machine Learning for Modeling Warm-Rain Cloud Microphysical Processes, J. Adv. Model. Earth Sy., 12, e2020MS002301,, 2020. 

Seiki, T. and Nakajima, T.: Aerosol effects of the condensation process on a convective cloud simulation, J. Atmos. Sci., 71, 833–853,, 2014. 

Shima, S., Sato, Y., Hashimoto, A., and Misumi, R.: Predicting the morphology of ice particles in deep convection using the super-droplet method: development and evaluation of SCALE-SDM 0.2.5-2.2.0, -2.2.1, and -2.2.2, Geosci. Model Dev., 13, 4107–4157,, 2020. 

Shima, S.-I.: Estimation of the Computational Cost of Super-Droplet Method (Fast Algorithms in Computational Fluids: theory and applications), Kyoto University Research Information Repository, RIMS Kokyuroku, 1606, 110–118, 2008. 

Shima, S.-I., Kusano, K., Kawano, A., Sugiyama, T., and Kawahara, S.: The super-droplet method for the numerical simulation of clouds and precipitation: A particle-based and probabilistic microphysics model coupled with a non-hydrostatic model, Q. J. Roy. Meteor. Soc., 135, 1307–1320,, 2009. 

Smagorinsky, J.: General circulation experiments with the primitive equations: I. The basic experiment, Mon. Weather Rev., 91, 99–164,;2, 1963. 

Smoluchowski, M.: Drei Vorträge über Diffusion, Brownsche Molekularbewegung und Koagulation von Kolloidteilchen, Z. Phys., 17, 557–585, 1916. 

Stephens, G. L.: Cloud feedbacks in the climate system: A critical review, J. Climate, 18, 237–273,, 2005. 

Stevens, B., Lenschow, D. H., Vali, G., Gerber, H., Bandy, A., Blomquist, B., Brenguier, J.-L., Bretherton, C., Burnet, F., and Campos, T.: Dynamics and chemistry of marine stratocumulus—DYCOMS-II, B. Am. Meteorol. Soc., 84, 579–594,, 2003.  

Stevens, B., Moeng, C.-H., Ackerman, A. S., Bretherton, C. S., Chlond, A., de Roode, S., Edwards, J., Golaz, J.-C., Jiang, H., and Khairoutdinov, M.: Evaluation of large-eddy simulations via observations of nocturnal marine stratocumulus, Mon. Weather Rev., 133, 1443–1462,, 2005. 

Stevens, D. and Bretherton, C.: Effects of resolution on the simulation of stratocumulus entrainment, Q. J. Roy. Meteor. Soc., 125, 425–439,, 1999. 

Twomey, S.: The nuclei of natural cloud formation part II: The supersaturation in natural clouds and the variation of cloud droplet concentration, Geofisica Pural e Applicata, 43, 243–249,, 1959. 

Unterstrasser, S., Hoffmann, F., and Lerch, M.: Collisional growth in a particle-based cloud microphysical model: insights from column model simulations using LCM1D (v1.0), Geosci. Model Dev., 13, 5119–5145,, 2020. 

Wicker, L. J. and Skamarock, W. C.: Time-splitting methods for elastic models using forward time schemes, Mon. Weather Rev., 130, 2088–2097,;2, 2002. 

Wood, R.: Stratocumulus clouds, Mon. Weather Rev., 140, 2373–2423,, 2012. 

Xu, X., Lu, C., Liu, Y., Luo, S., Zhou, X., Endo, S., Zhu, L., and Wang, Y.: Influences of an entrainment–mixing parameterization on numerical simulations of cumulus and stratocumulus clouds, Atmos. Chem. Phys., 22, 5459–5475,, 2022. 

Xue, L., Teller, A., Rasmussen, R., Geresdi, I., and Pan, Z.: Effects of aerosol solubility and regeneration on warm-phase orographic clouds and precipitation simulated by a detailed bin microphysical scheme, J. Atmos. Sci., 67, 3336–3354,, 2010. 

Xue, L., Teller, A., Rasmussen, R., Geresdi, I., Pan, Z., and Liu, X.: Effects of aerosol solubility and regeneration on mixed-phase orographic clouds and precipitation, J. Atmos. Sci., 69, 1994–2010,, 2012. 

Yamaguchi, T. and Randall, D. A.: Cooling of entrained parcels in a large-eddy simulation, J. Atmos. Sci., 69, 1118–1136,, 2012. 

Yang, F., Hoffmann, F., Shaw, R. A., Ovchinnikov, M., and Vogelmann, A. M.: An Intercomparison of Large-Eddy Simulations of a Convection Cloud Chamber Using Haze-Capable Bin and Lagrangian Cloud Microphysics Schemes, J. Adv. Model. Earth Sy., 15, e2022MS003270,, 2023. 

Yin, C.: wangyouyue/GMD_2024_code_Yin: SCALE-SDM_Yin2024_v1.0.4, Zenodo [code],, 2024a. 

Yin, C.: Supplement to the article “Simulation of marine stratocumulus using the super-droplet method: Numerical convergence and comparison to a double-moment bulk scheme”, Zenodo [video],, 2024b. 

Zelinka, M. D., Myers, T. A., McCoy, D. T., Po-Chedley, S., Caldwell, P. M., Ceppi, P., Klein, S. A., and Taylor, K. E.: Causes of higher climate sensitivity in CMIP6 models, Geophys. Res. Lett., 47, e2019GL085782,, 2020. 

Short summary
We investigate numerical convergence properties of a particle-based numerical cloud microphysics model (SDM) and a double-moment bulk scheme for simulating a marine stratocumulus case, compare their results with model intercomparison project results, and present possible explanations for the different results of the SDM and the bulk scheme. Aerosol processes can be accurately simulated using SDM, and this may be an important factor affecting the behavior and morphology of marine stratocumulus.