Super-Droplet Method (SDM) in Cloud Microphysics

• Super-Droplet Method (SDM) is a particle-based probabilistic framework that represents cloud and aerosol populations using weighted superdroplets.
• It employs Lagrangian tracking and Monte Carlo collision–coalescence to simulate condensation, evaporation, and precipitation processes with high fidelity.
• SDM bridges particle-resolved physics with Eulerian diagnostics, serving as a benchmark for reduced microphysical models and machine-learning emulators.

Searching arXiv for the cited SDM papers to ground the article in current literature. The Super-Droplet Method (SDM) is a Lagrangian, particle-based, probabilistic microphysics framework in which a finite set of computational particles, called superdroplets, statistically represents an underlying population of real droplets, aerosol particles, cloud droplets, drizzle, or raindrops. A superdroplet carries a multiplicity, namely the number of identical real particles it represents, together with continuous microphysical attributes such as size, composition, and position. In SDM, condensation or evaporation, sedimentation, and collision–coalescence are applied to these Lagrangian particles, while the carrier flow and thermodynamic fields are typically evolved on an Eulerian grid. The method was introduced as a particle-based and probabilistic alternative to mean-field formulations such as the Smoluchowski coagulation equation, and subsequent work has established SDM both as a high-fidelity research tool for cloud microphysics and as a benchmark for reduced parameterizations and machine-learning surrogates (Bartman et al., 2021, Arabas et al., 2012).

1. Conceptual basis and state representation

SDM replaces the explicit simulation of every physical droplet by a weighted ensemble of superdroplets. Each superdroplet represents a multiplicity $\xi$ of identical real droplets located in the same small region, and the evolution is Lagrangian: the algorithm follows each superdroplet’s properties along its trajectory rather than evolving a distribution function on an Eulerian grid. In this representation, the dispersed phase is described by superdroplet attributes such as multiplicity, mass or volume, composition, position, and related intensive properties. A convenient abstract form is an attribute vector

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$

with $\xi_j \in \mathbb{N}$ the multiplicity and $m_j$ a droplet mass or size variable (Sharma et al., 2024, Bartman et al., 2021).

A defining conceptual feature of SDM is that it need not impose separate prognostic categories for aerosol, cloud droplets, drizzle, and rain. In one LES application, the same superdroplet object can represent a dry aerosol, a CCN in equilibrium with subsaturated air, a cloud droplet, or a drizzle or rain drop, and diagnostic distinctions such as cloud water versus rain water are made a posteriori by thresholding radius. In that study, cloud and rain water mixing ratios were diagnosed by summing over particles with radius smaller and larger than $40\,\mu\text{m}$, respectively (Arabas et al., 2012).

This representation differs from bulk and spectral Eulerian schemes in two ways. First, SDM samples the droplet population directly rather than discretizing a number-density field on a size grid. Second, it avoids assuming a fixed analytic droplet size distribution during time evolution. In the McSnow-based warm-rain simulations used to train SuperdropNet, for example, the detailed droplet size distribution $f(x)$ over droplet mass $x$ is represented by the evolving superdroplet ensemble and can be sampled or reconstructed without imposing a gamma or lognormal form during the simulation (Sharma et al., 2024).

The same superdroplet formalism also provides a natural route from particle-resolved physics to Eulerian diagnostics. In warm-rain box simulations, the bulk state vector used for machine-learning emulation is

$y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$

where $L_c, L_r$ are cloud and rain water masses and $N_c, N_r$ are cloud and rain droplet numbers. With a cloud–rain threshold $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$0, these moments are computed from the superdroplet population as

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$1

with $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$2 the mass of superdroplet $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$3 and $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$4 its multiplicity (Sharma et al., 2024).

2. Microphysical processes and Monte Carlo collision–coalescence

In SDM, microphysical evolution is represented by explicit particle-level processes rather than bulk rate closures. The McSnow implementation discussed in the SuperdropNet study includes collision–coalescence, condensation or evaporation, and sedimentation, with collision–coalescence the primary focus (Sharma et al., 2024). In the RICO LES application, the superdroplets undergo gravitational settling, condensational growth or evaporation, and stochastic collision–coalescence; the resulting latent heating and vapor tendencies are then fed back to the Eulerian model (Arabas et al., 2012).

Condensational growth is represented by a classical diffusion- and heat-conduction-limited growth law. In the trade-wind cumulus study, the radius evolution is summarized as

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$5

or, more generally,

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$6

where $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$7 is supersaturation with respect to liquid water and $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$8 is the equilibrium supersaturation from the Köhler curve. Because SDM integrates individual particle radii directly, that work describes the condensation scheme as “diffusive-error free” in radius space (Arabas et al., 2012).

The collision–coalescence step is the characteristic probabilistic element of SDM. Within each collision volume or LES grid cell, superdroplets are randomly paired, a collision kernel is evaluated for each pair, and a Monte Carlo draw determines whether a collision occurs and how many real droplets coalesce. In the Shima-type formulation described in the PySDM case study, with $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$9 superdroplets in a cell and timestep $\xi_j \in \mathbb{N}$0, the pairwise collision probability is

$\xi_j \in \mathbb{N}$1

where $\xi_j \in \mathbb{N}$2 is the collision kernel and $\xi_j \in \mathbb{N}$3 the collision volume. An integer number of collision events is then sampled via

$\xi_j \in \mathbb{N}$4

with $\xi_j \in \mathbb{N}$5. This construction permits multiple collisions between the same superdroplet pair in one time step when $\xi_j \in \mathbb{N}$6 (Bartman et al., 2021, Ware et al., 5 Sep 2025).

After an accepted collision, multiplicities and extensive attributes are updated so that water mass is conserved and the number of superdroplets is kept constant or non-increasing. In one notation, if $\xi_j \in \mathbb{N}$7 is the collector and $\xi_j \in \mathbb{N}$8 the donor, the realizable number of events is capped by

$\xi_j \in \mathbb{N}$9

If donor multiplicity remains positive,

$m_j$0

where $m_j$1 denotes an extensive attribute. If the donor multiplicity is exhausted, the collector is split to preserve the superdroplet count. These pairwise rules conserve extensive quantities such as water mass exactly at the update level (Ware et al., 5 Sep 2025, Bartman et al., 2021).

The physical interpretation of traditional warm-rain process labels changes under SDM. Autoconversion, accretion, and self-collection are not separate parameterized tendencies but emergent outcomes of the same collision–coalescence mechanism operating across different parts of the size spectrum. In warm-rain SDM, autoconversion occurs when cloud-sized droplets collide and a product crosses the cloud–rain threshold $m_j$2; accretion arises from rain–cloud collisions; self-collection occurs when collisions do not change category. This unified treatment is repeatedly contrasted with bulk schemes in which these rates are given by analytic formulas (Sharma et al., 2024, Arabas et al., 2012).

3. Coupling to Eulerian flow models and observational diagnostics

SDM is commonly embedded in an Eulerian–Lagrangian framework. The fluid variables—velocity, temperature, water vapor, and related thermodynamic fields—are evolved on an Eulerian grid, while the dispersed phase is represented by Lagrangian superdroplets. In the CReSS-coupled RICO configuration, the two components are sequentially advanced: during a given time step, each component sees the other’s state as constant. The LES provides interpolated velocity and thermodynamic fields to the superdroplets; SDM returns latent heat release and water vapor tendencies as Eulerian source terms (Arabas et al., 2012).

That study exemplifies the operational workflow of an LES–SDM system. The RICO case used a $m_j$3 km quarter domain with coarse, middle, and high resolutions of $m_j$4, $m_j$5, and $m_j$6, corresponding to $m_j$7 and $m_j$8 m and $m_j$9 and $40\,\mu\text{m}$0 m. The SDM time stepping was split into condensational growth at $40\,\mu\text{m}$1 s, collisional growth at $40\,\mu\text{m}$2 s, and particle motion at $40\,\mu\text{m}$3 s or $40\,\mu\text{m}$4 s depending on resolution. Superdroplet sampling densities ranged from $40\,\mu\text{m}$5 to $40\,\mu\text{m}$6 per grid cell in the coarse configuration, $40\,\mu\text{m}$7 to $40\,\mu\text{m}$8 per cell in the middle configuration, and $40\,\mu\text{m}$9 per cell in the high-resolution run (Arabas et al., 2012).

A major methodological advantage of SDM is the direct comparability of model output with particle-counting observations. Because the method predicts individual particle properties, one can construct synthetic diagnostics matched to specific instrument size ranges and sampling conventions. In the RICO comparison, Fast-FSSP-like cloud microstructure diagnostics included cloud droplet number concentration, effective radius

$f(x)$0

liquid water content, the parameter

$f(x)$1

and the radius standard deviation $f(x)$2, restricted to particles in the $f(x)$3–$f(x)$4 range and in-cloud cells with CDNC $f(x)$5. Precipitation spectra were compared with OAP-2DS observations using 61 bins spanning the instrument’s size range (Arabas et al., 2012).

The reported results show both the strengths of the method and its sensitivity to resolution. The simulations reproduced height-resolved structure in effective radius, liquid water content, and rain spectra that compared favorably with previously published RICO aircraft observations. They also reproduced the observed lack of hydrometeors smaller than $f(x)$6 in diameter within rain shafts. At the same time, coarse-grid simulations inhibited convection, and low superdroplet densities produced nearly monodisperse spectra, with $f(x)$7 and $f(x)$8. In that case study, approximately $f(x)$9 superdroplets per grid cell were found sufficient for convergence of $x$0 and $x$1, while larger increases yielded diminishing returns for those metrics (Arabas et al., 2012).

This suggests a broader methodological role for SDM. Because it preserves physically interpretable particle-level information while remaining bi-directionally coupled to the resolved dynamics, SDM can serve simultaneously as a process model, a synthetic-observation generator, and a reference against which simpler microphysical closures are evaluated.

4. Bulk moments, reduced models, and machine-learning emulation

Although SDM evolves a particle ensemble, many atmospheric models require a bulk interface. The SuperdropNet work formalizes this interface by defining a time-stepping map $x$2 on the four bulk moments $x$3, sampled every $x$4, and seeking

$x$5

for trajectories generated by superdroplet simulations. The emulator uses an explicit Euler-like update

$x$6

where the primary inputs are the current bulk moments and the additional inputs $x$7 include $x$8, $x$9, $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$0, $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$1, and $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$2 (Sharma et al., 2024).

This design choice is significant because it differs from earlier machine-learning approaches that predicted process rates such as autoconversion and accretion and then recovered moment tendencies through a bulk scheme’s analytic formulas. SuperdropNet instead learns the net tendency

$y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$3

thereby avoiding the inherited assumptions of the Seifert–Beheng 2001 two-moment scheme. In the same work, the reference bulk scheme SB2001 evolves the same four moments but parameterizes tendencies with analytic formulas,

$y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$4

where $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$5 is autoconversion, $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$6 accretion, and $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$7 self-collection. SDM does not require such fixed analytic closures, because the corresponding transitions are emergent from the particle collisions themselves (Sharma et al., 2024).

The training data for SuperdropNet were generated with McSnow in a zero-dimensional box of volume $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$8, multiplicity $y_t = [L_c(t),\, L_r(t),\, N_c(t),\, N_r(t)],$9, and warm-rain microphysics including collision–coalescence, condensation or evaporation, and sedimentation. Initial size distributions were gamma distributions with

$L_c, L_r$0

giving $L_c, L_r$1 distinct initial conditions. Because SDM collision–coalescence is stochastic, the authors ran $L_c, L_r$2 realizations per initial condition and trained on trajectories averaged over those realizations, explicitly noting that training on all individual runs led to overfitting and worse performance (Sharma et al., 2024).

The resulting neural network is a fully connected feedforward architecture with $L_c, L_r$3 hidden layers, $L_c, L_r$4 neurons each, ReLU activations, and approximately $L_c, L_r$5 trainable parameters. Its training procedure emphasizes autoregressive stability. Rather than minimizing only the one-step loss

$L_c, L_r$6

the model is trained with a multi-step loss

$L_c, L_r$7

using the “pushforward trick,” with $L_c, L_r$8 and $L_c, L_r$9 for $N_c, N_r$0, and a curriculum that increases rollout length from $N_c, N_r$1 to $N_c, N_r$2 (Sharma et al., 2024).

Physical constraints are built into this emulator to mirror warm-rain SDM behavior. Total liquid mass is conserved by predicting only the update to $N_c, N_r$3 and enforcing

$N_c, N_r$4

At inference time, positive updates to cloud mass or cloud number are set to zero, and negative values of $N_c, N_r$5 are clipped to zero with consistent mass adjustments. In the reported results, SuperdropNet predicted hydrometeor states and cloud-to-rain transition times more accurately than previous ML emulators and matched or outperformed bulk moment schemes in many cases, with especially strong performance for number-based transition timing (Sharma et al., 2024).

A separate line of work extends machine learning within LES–SDM to condensational growth rather than collision–coalescence. Using DNS of entrainment and mixing, a 2024 study derives the physically correct superdroplet growth law

$N_c, N_r$6

with effective supersaturation

$N_c, N_r$7

where $N_c, N_r$8 denotes the DNS droplets represented by superdroplet $N_c, N_r$9. For multiplicity $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$00, filtered supersaturation at the superdroplet location explained only $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$01 of the variance in $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$02, whereas an MLP trained on filtered DNS variables achieved $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$03. This shows that machine learning can also be used to correct a specific SDM closure—the well-mixed supersaturation assumption—without abandoning the superdroplet representation itself (Divyaprakash et al., 2024).

5. Computational complexity, stochasticity, and algorithmic corrections

A central attraction of SDM is that its collision algorithm can be made linear in the number of superdroplets. If $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$04 superdroplets occupy a cell, a naive all-pairs method would scale as $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$05. The Shima-style SDM instead samples only $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$06 random non-overlapping pairs per step and rescales the pairwise probability accordingly, yielding approximately linear cost per time step. This linear scaling, together with the fact that the state vector length remains constant or non-increasing, underlies the method’s suitability for CPU and GPU parallelization (Bartman et al., 2021, Arabas et al., 2012).

Nevertheless, the method remains computationally expensive in realistic atmospheric settings. The SuperdropNet study emphasizes that droplet-based Lagrangian schemes are “underutilized due to their large computational overhead,” noting that cost scales at least linearly and, for some algorithms or naive pairing, effectively quadratically with the number of superdroplets. In the McSnow warm-rain box model, reducing stochastic noise required a large box volume of $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$07 together with multiplicity $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$08, and the authors explicitly identified SDM as a computationally expensive “gold standard” (Sharma et al., 2024).

The stochasticity of the Monte Carlo collision step has both physical and numerical consequences. Repeated SDM simulations from the same initial condition produce different trajectories because of random collision sampling; part of that variability reflects physical stochasticity, and part is an artifact of finite superdroplet number, box size, and multiplicity. In the SuperdropNet dataset design, larger box volume and averaging over $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$09 repeated realizations were used to suppress sampling noise. The reported analyses also showed that larger shape parameter $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$10, corresponding to narrower initial size distributions, was associated with more variability in bulk moments across repeated runs (Sharma et al., 2024).

More recent work has identified an additional systematic error source in the collision algorithm, namely the collision deficit. For a pair $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$11 with predicted collision count $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$12 and realizable count $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$13, the deficit is defined as

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$14

the number of real donor droplets that should have collided according to the stochastic rate but cannot be represented because of the multiplicity constraint. This is a one-sided bias: it systematically underestimates collision activity rather than contributing symmetric Monte Carlo spread (Ware et al., 5 Sep 2025).

The 2025 adaptive-timestep analysis shows that the collision deficit increases with timestep and superdroplet count in the classical Safranov–Golovin test case. It also depends on initialization strategy: constant-multiplicity sampling tends to produce larger deficits than schemes with a dynamic range of multiplicities, such as uniform-in-$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$15 or uniform-in-$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$16. This result is notable because increasing $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$17, which ordinarily improves spectral resolution, can in a deficit-dominated regime increase the RMSE relative to the analytic solution when fixed timesteps are used (Ware et al., 5 Sep 2025).

To remove this bias, the adaptive scheme constrains the local collision substep so that, for every tested pair,

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$18

At the cell level, the substep is chosen as

$$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$19

and repeated until the base step is completed. In the reported results, this made the collision deficit identically zero by construction, restored monotonic accuracy improvement with increasing $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$20, and, in 2-D flow-coupled simulations, removed a systematic delay in precipitation onset as quantified by the $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$21 metric (Ware et al., 5 Sep 2025).

At the implementation level, PySDM demonstrates how SDM can be organized for heterogeneous hardware. The package defines a common API for CPU and GPU backends, with data structures such as Storage, IndexedStorage, Index, and PairIndicator. The CPU backend uses Numba with NumPy arrays and multi-threading, whereas the GPU backend uses ThrustRTC and CURandRTC. In the reported benchmark, the GPU backend achieved roughly a factor of $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$22 speedup over the CPU backend for large state vectors on commodity hardware, while retaining the same algorithmic workflow and collision logic (Bartman et al., 2021).

6. Scope, limitations, and current directions

Across the cited studies, SDM is consistently presented as a high-fidelity framework for warm-cloud and warm-rain microphysics, with extensions toward broader multiphase applications. Its strengths include direct particle-level treatment of aerosol activation, condensational growth, collision–coalescence, and sedimentation; minimal reliance on bulk parameterizations; exact conservation of extensive quantities at the pairwise update level; and straightforward generation of synthetic particle-counting diagnostics (Arabas et al., 2012, Bartman et al., 2021).

The method also has clearly documented limitations. In the RICO LES application, missing processes included turbulence enhancement of collision efficiency, drop breakup, aerosol sources, and subgrid mixing below LES resolution. In the SuperdropNet work, the SDM-generated training data covered only warm rain in a zero-dimensional box, with no vertical motion, no supersaturation dynamics, and no coupled dynamics–microphysics feedbacks. In the DNS-informed condensational-growth study, the learned correction was trained for a specific regime—forced homogeneous isotropic turbulence with $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$23, filter width $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$24, multiplicities $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$25, and purely condensational physics. These constraints delimit the regimes in which current results can be interpreted directly (Arabas et al., 2012, Sharma et al., 2024, Divyaprakash et al., 2024).

Current research directions fall into two broad classes. One extends SDM itself. The adaptive-timestep analysis suggests that algorithmic details of the collision step remain consequential even for well-established SDM formulations, especially when precipitation onset is sensitive to a small subset of strongly collisional regions. The same work also argues that the adaptive strategy should transfer naturally to SDM extensions involving collisional breakup, ice, and mixed-phase processes (Ware et al., 5 Sep 2025).

The other direction uses SDM as a training or calibration target for reduced models. SuperdropNet treats SDM as a reference microphysical dynamics generator and learns a bulk-interface emulator designed as a drop-in replacement for bulk warm-rain schemes, while the DNS-informed LES–SDM study learns a local correction to the superdroplet growth law by replacing filtered supersaturation $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$26 with an ML-predicted effective supersaturation $$\mathbf{A}_j = \big( \xi_j,\; m_j,\; \text{other attributes}\big),$$27. Taken together, these results suggest two complementary strategies: emulate SDM bulk evolution directly, or retain the SDM structure and correct one of its constituent closures (Sharma et al., 2024, Divyaprakash et al., 2024).

In that sense, SDM now occupies a dual position in atmospheric microphysics. It remains a process-resolving Lagrangian Monte Carlo method in its own right, but it also functions as a physically detailed reference system for benchmarking bulk schemes, constructing synthetic observations, identifying algorithmic biases, and training reduced-order surrogates. A plausible implication is that the future role of SDM will continue to be shaped both by advances in particle-based numerics and by hybrid workflows in which SDM-generated data support more economical but SDM-consistent parameterizations.