Ice Particle Properties Inferred From Aggregation Modelling

We generated a large number 105,000 of aggregates composed of various monomer types and sizes using an aggregation model. Combined with hydrodynamic theory, we derived ice particle properties such as mass, projected area, and terminal velocity as a function of monomer number and size. This particle ensemble allows us to study the relation of particle properties with a high level of detail which is often not provided by in situ measurements. The ice particle properties change rather smoothly with monomer number. We find very little differences in all particle properties between monomers and aggregates at sizes below 1 mm which is in contrast to many microphysics schemes. The impact of the monomer type on the particle properties decreases with increasing monomer number. Whether, for example, the terminal velocity of an aggregate is larger or smaller than an equal‐size monomer depends mostly on the monomer type. We fitted commonly used power laws as well as Atlas‐type relations, which represent the saturation of the terminal velocity at large sizes (terminal velocity asymptotically approaching a limiting value) to the data set and tested the impact of incorporating different levels of complexity with idealized simulations using a 1D Lagrangian super particle model. These simulations indicate that it is sufficient to represent the monomer number dependency of ice particle properties with only two categories (monomers and aggregates). The incorporation of the saturation velocity at larger sizes is found to be important to avoid an overestimation of self‐aggregation of larger snowflakes.


Introduction
The terminal velocity v term of ice monomers and aggregated ice particles and its relation to size has manifold impacts on precipitation and radiative effects of ice containing clouds. For example, Morales et al. (2019) show that parameters describing v term of aggregates have the largest impact on the precipitation of simulated orographic clouds. Experiments with global climate simulations revealed that also radiative fluxes are very sensitive to changes in v term (Jakob, 2002). Sanderson et al. (2008) found that v term of ice is the second most influential parameter for the climate sensitivity in their multimember perturbed physics General circulation model ensemble. Constraining v term of cloud ice and aggregated ice particles can reduce the degrees of freedom in model tuning (e.g., to improve top of atmosphere radiative fluxes Schmidt et al., 2017) and improve the physical consistency in atmospheric models. dependence on size which are implemented in widely used two-moment schemes are shown in Figure 1b. When comparing these relations with observations (Figure 1a), we miss the saturation behavior of v term for larger sizes in most relations. This discrepancy is expected as most schemes use power laws, which are unable to represent a saturation behavior. Alternative "Atlas-type" three-parameter fits have been suggested (Seifert et al., 2014), but so far, they have not been tested thoroughly. The recent Predicted Particle Properties (P3) scheme (Morrison & Milbrandt, 2015) uses only one ice category and a look-up table approach for v term , which better matches the saturation at large sizes. At the smaller size range, the snow category is found for all schemes to fall significantly faster than the ice category with the same size. Considering that v term depends strongly on m and A of the particle, it might sound plausible, that for example, an aggregate of a few plates should fall faster than a single plate of the same size. Unfortunately, most observations do not provide sufficiently detailed information about monomer number and type which would be needed to answer the question of whether there exists a "jump" in v term for the number of monomers exceeding a certain threshold. Fairly direct observations of the particles' m and A are only available from manual, particle-based observations (e.g., Locatelli & Hobbs, 1974).
An interesting new tool to better understand the underlying principles of aggregation and its effects on particle properties are aggregation models (Hashino & Tripoli, 2011;Leinonen & Moisseev, 2015;Ori et al., 2014;Przybylo et al., 2019;Westbrook et al., 2004a). Those models use idealized monomer shapes (e.g., dendrites, needles, plates, and columns) with particle properties matched to in situ observations. Aggregates simulated with the model by Westbrook et al. (2004a) helped to better understand theoretical scaling relations associated to aggregation such as the increase of aggregate mass with size by a power of two (Westbrook et al., 2004b), which was known from several previous in situ observations. This model has been extended by Leinonen and Moisseev (2015) providing a large number of monomer shapes and also provides an option to rime the aggregate (Leinonen & Szyrmer, 2015). This allowed to better understand the evolution of size and mass of a large number of aggregates which were increasingly rimed (Seifert et al., 2019).
To infer v term from modeled ice particles or aggregates, computational fluid dynamics is an accurate but also computational costly method. It has been recently applied to idealized ice particle shapes (Bürgesser et al., 2019;Hashino et al., 2016;Nettesheim & Wang, 2018), and more computations with more complex shapes can be expected shortly. Hydrodynamic theory is a computational cheaper alternative to calculate v term based on a number of bulk particle characteristic, rather than the complex 3D shape (e.g., Böhm, 1992;Heymsfield & Westbrook, 2010;. The accuracy of hydrodynamic theories has recently been evaluated by ice particle analogs falling in an oil tank (Westbrook & Sephton, 2017). The experimental results show deviations smaller than 20% for the Heymsfield and Westbrook (2010) theory. A problematic aspect of these theories is still the formulation of the scaling toward higher Reynolds number (i.e., large particles) and the simulation of more complex particle shapes (Westbrook & Sephton, 2017).
Aggregation models in combination with hydrodynamic theory have recently been used to study v term of aggregates (Hashino & Tripoli, 2011;Schmitt et al., 2019). Hashino and Tripoli (2011) identified a dependency of the aggregation rate and aggregate mass on the mean size and type of the monomers. Schmitt et al. (2019) analyzed v term and its variability of simulated aggregates composed of hexagonal prisms taken from a monodisperse monomer size distribution. They found that the variability of v term is caused by the variability of the number of monomers N mono and the monomers' aspect ratio.
In this study, we aim to study the dependency of m, A, and v term on size, monomer number, and type. For this, we create a large number of aggregates with various monomer types including also mixtures of different monomer types. The monomer size is sampled from a size distribution rather than a constant size to better represent real ensembles of aggregates. Central questions of this study are, how important is the monomer number and type information for parameterizing aggregate properties and how well can they be parameterized by different functional relations?
To answer these questions, we describe in section 2 the aggregation model and the created data set of unrimed aggregates as well as the hydrodynamic theory to calculate v term based on m and A of these particles. The simulated particle properties are compared to in situ observations in section 3. Section 4presents several parameterizations of the particle properties. Finally, in section 5, we use a 1D Lagrangian particle model to test the impact of including different complexity of particle properties for aggregation 2. Methods

Aggregation Model
We use the aggregation model developed by Leinonen and Moisseev (2015) which includes a large number of realistic monomers (hexagonal plates, dendrites, columns, and needle). Originally, the aggregation model was designed to produce realistic snow particle structures which can then be used to calculate their scattering properties (Leinonen & Moisseev, 2015;Leinonen et al., 2018). The model has also been used to systematically investigate microphysical processes, such as riming (Seifert et al., 2019).
The shape characteristics (length, thickness, etc.) of the monomers are predefined by geometric relations based on in situ observations (Leinonen & Moisseev, 2015). The aggregation process starts with generating N mono monomers with sizes following a predefined inverse exponential probability density function p d (D max ), (1) where λ −1 is the size parameter of the monomer distribution and D max is the maximum size of the monomer. The higher λ −1 , the larger are the sizes of the monomers.
The monomers sizes are sampled from the monomer distribution and assembled until an aggregate consisting of N mono monomers is build up. In each aggregation step, pairs of particles are selected according to a simplified gravitational collection kernel. The probability distribution of collision among each possible particle pair is calculated as being proportional to the particle geometric cross sections and differential fall speed (Westbrook et al., 2004a). The two colliding particles form an aggregate which then becomes one of the candidates for the next aggregation step. This process includes the collision between aggregates. The aggregation code is publicly available at https://github.com/jleinonen/aggregation, and more details on the implementation can be found in Leinonen and Moisseev (2015). During the aggregation process, the collecting particles are partially aligned with the principal axis in the x-y plane. Rotations around the principal axis are performed randomly with a standard deviation of 40°. The collected particles are randomly aligned, which mimics the complex flow in the vicinity of other particles (Leinonen & Moisseev, 2015).
The aggregation simulations performed in this study differ from previous studies in two main aspects. The first aspect is the resolution of the particle structure. The particle is internally represented by a three-dimensional lattice with a predefined distance of the volume elements of typically 40 μm. This distance was found to be sufficiently small for scattering computations, while being coarse enough in order to keep the numerical costs for the scattering computations in a reasonable range. However, we discovered that for small particle sizes, the theoretical relations for certain particle properties (see Fig. 1 in Leinonen and Moisseev (2015)) are not exactly matched by the discretized particle. This discrepancy can be easily explained when considering for example that plate monomers with D max < 3.03 mm consist of only one layer of volume elements if the default resolution of 40 μm is used. This does not necessarily affect the aggregate properties of those monomers as shown in Leinonen and Moisseev (2015); however, in our study, the focus is to investigate the transition from small to larger sizes particles. Hence, we need to refine the resolution especially for small particles.
As a compromise between computational feasibility and having fine enough resolved particles, aggregates with N mono ≤ 100 are simulated with a resolution of 5 μm, while aggregates with N mono ≥ 100 are simulated with 10-μm resolution. With a resolution of 5 μm (10 μm), a plate monomer with D max = 3 mm has a thickness of 4 (8) volume element layers. It should be noted that the sensitivity to resolution is smaller for monomer types with less extreme aspect ratios (e.g., columns).
The second major difference to previous aggregation studies using the model by Leinonen and Moisseev (2015) is that we extended the code in a way that we can also generate aggregates composed of monomers with different habits. The motivation for this new feature was based on observations that larger snowflakes often consist of a mixture of dendrites and needles (Lawson et al., 1998). The modified code extends Equation 1 to be the joint distribution of multiple mono-dispersed distributions. Each monomer distribution is defined by its own settings (e.g., monomer type, mean size, and truncation). The joint distribution is defined by the relative weights of each mono-dispersed distribution. These modifications have been merged to the main aggregation code and are also publicly available.
In order to account for a large variability of naturally observed particle shapes (Bailey & Hallett, 2009), we simulated a large suite of aggregates consisting of plates, columns, dendrites, needles, and mixtures of dendrites and columns. The m−D max and A−D max relations for the monomers are given in Table 1. Two sets of aggregates with mixed monomer types were created. For the first mixture, the selection of the monomer type is random with the same probability density function for both monomer types ("Mix1"). This would represent a scenario, where dendrites and needles coexist with similar PSD and likelihood of aggregation. For the second mixture, the monomers with D max < 1 mm are columns, while dendrites are taken for larger monomers ("Mix2"). This choice is motivated by the fact that at temperatures below −20°C, the particle shape is less distinct but mostly described by polycrystals, while at temperatures between −20°C and −10°C, one finds more planar and dendritic crystals (Bailey & Hallett, 2009). Considering a thick cloud, we could assume that the small polycrystal or columnar crystals forming in the upper part of the cloud begin to form the first aggregate and then further grow by collection of larger dendrites at lower layers. Of course, both scenarios are quite ad hoc, and more detailed studies are needed to better understand the real properties of mixed-monomer aggregates. Our mixtures are thus rather intended to qualitatively analyze the differences of mixed monomer aggregates compared to single-monomer type aggregates (as done in another recent study by Dunnavan et al., 2019).
The aggregation process strongly depends on the number concentration of particles and their relative terminal velocity differences. In conditions which are less favorable for aggregation (e.g., low number concentration), the particles can grow by depositional growth to relatively large sizes before aggregation becomes the dominant process. It is therefore possible that aggregation involves very different monomer sizes. In order to account for this variability, we vary λ −1 in a large range from 50 μm to 10 mm with 500 different values of λ −1 , spaced evenly in the logarithmic space. The monomer distribution is limited to sizes of 100 μm up to 3 mm following Leinonen and Moisseev (2015) in order to be consistent with the typical size range of observed ice particles. Due to this truncation of the inverse exponential distribution, the mean monomer size differs from λ −1 and ranges from 150 μm to 1.48 mm.
The spacing of the monomer number (Table 2) is finer at low N mono and becomes more coarse at larger numbers. In this way, we can investigate the changes at small monomer numbers with greater detail. In fact, we expect the largest changes in snow properties at the transition from single monomers to aggregates composed of few pristine crystals as shown in earlier studies (Dunnavan et al., 2019;Schmitt & Heymsfield, 2010). The coarser spacing of N mono also limits computational costs. With our settings, we obtain maximum aggregates sizes ranging from 3 to 5 cm which means that we include also the typically observed large snowflakes during intense snowfall on the ground (Lawson et al., 1998).
In Figure 2, several examples of similar sized aggregates simulated with different combinations of λ −1 , N mono , and monomer types are shown. In total, 105,000 particles were simulated. Apart from the visual differences of shapes and structure, also the particle properties such as mass, area, or terminal velocity show a wide range of values although all aggregates have maximum sizes ranging between 3 and 5 mm.

Hydrodynamic Models
Hydrodynamic models are needed in order to derive the terminal velocity v term from the particle's mass m, projected area A, and maximum size D max . The most commonly used hydrodynamic models are Böhm (1992, hereafter B92), Khvorostyanov and Curry (2005, hereafter KC05), and Heymsfield and Westbrook (2010, hereafter HW10). All models are based on particle boundary layer theory and rely on the Best number (X) approach (Abraham, 1970). v term is calculated via, Note. All monomers have a grid resolution of 5 μm. The shapes are predefined in the aggregation model and mostly based on Pruppacher and Klett (1998) (see Fig. 1 in Leinonen and Moisseev, 2015).
where η is the dynamic viscosity, Re the Reynolds number (parameterized as a function of X), and ρ a is the air density. X is defined as where C d is the drag coefficient. The proportionality of X to the particle properties is given by for B92.
For this study, we decided to use B92 because it best represents the saturation of v term for our simulated particles at larger aggregate sizes ( Figure A2) in accordance with observations ( Figure 1). B92 includes an Note. D max denotes the maximum size range of the generated aggregates in the data set. empirical correction of X due to wake turbulence which increases the drag of large particles. X depends on the aspect ratio α, which is larger than one for prolate and smaller than one for oblate particles. For this study, we set α to 1.0, because aggregates with small values of N mono are not easily classifiable as either prolate or oblate and show in general a large variability of α .
To be able to interpret the dependency of v term on N mono in section 4.3, we sketch here how v term scales with D max in the simplified case of Re≪1 (Stokes drag) and Re ≫ 1 (Newtonian drag). For Re≪1, C D is approximately proportional to 1/Re. Inserting this approximation and Equations 3 and 4 into Equation 2 yields v term ∼ mD −0:5 max A −0:25 : In both extreme cases of Re, v term increases the faster with size the higher b m −0.25b A is, and we expect this also to be in between these cases where Re transitions from Re∼X to Re∼X 0.5 . This has certain implications for the dependency of v term on N mono (section 4.3).
The differences between the three hydrodynamic models as well as an analysis of the potential impact of changing to different hydrodynamic models is discussed in the Appendix A2.

Mass-and Area-Size Relations
Particle properties, such as m, A, and D max , are used in hydrodynamic models to calculate v term (section 2.2). We compare our relations of these particle properties and v term with frequently used relations that are based on in situ measurements from Locatelli 306 and Hobbs (1974, LH74) and Mitchell (1996, M96). LH74 defined an equivalent diameter that is equal to "the diameter of the smallest circle into which the aggregate as photographed will fit without changing its density." M96 collected observations as a function of D max without specifying the exact definition. The definitions of particle size used in these studies are limited by the observation equipment used, and the conversion from one to the other is not trivial. In our simulation study, we can access the full 3D structure of the particles and use the true maximum size (i.e., the maximum distance between any two points of the particle) as size definition.
Except for the aggregates of dendrites, which have a considerably lower density than LH74 aggregates of dendrites, the absolute value of m of the simulated aggregates is similar to the observations, where the same monomer type is available ( Figure 3). The slope of the m−D max relation from this study is comparable to the slope from M96, while LH74 report lower slopes for the aggregates of dendrites. The m−D max relation of the mixed aggregates ("aggregates of unrimed radiating assemblages of plates, side planes, bullets, and columns," LH74 mix), however, has a similar slope to the simulated Mix2 aggregates. The mixS3 and sideplane aggregates from M96 are similar to many simulated aggregates (composed of different monomers).
M96 derived A−D max relations for "assemblages of planar polycrystals in cirrus clouds" (M96 polycrystal in Figure 3) based on observations in a relatively small size range and applied them to other aggregate types. This A−D max relation is also used in several microphysics schemes (Brdar & Seifert, 2018;Morrison & Milbrandt, 2015). The absolute value of A given in M96 is slightly higher than A of the simulated particles from this study (except for the aggregates of plates). The slope of the A−D max relations is slightly higher (b A = 1.88) in M96 observations compared to the relations from this study (1.79 < b A < 1.88). Observations of aggregates composed of the same monomer types than the one used in these studies are not available.

Terminal Velocity-Size Relations
Observations of v term versus size have been reported using several different definitions of the diameter (Szyrmer & Zawadzki, 2010). To facilitate a consistent comparison between the observations from the PIP instrument (which are described in Section A0.1) and v term of the simulated aggregates, we use the same bin sizes as the PIP instrument to derive the median v term . Moreover, we derive the maximum dimension from a side projection of the modeled particle in the same way as in the observations from the PIP instrument described by ( At small sizes (D max < 1 mm), v term of the simulated aggregates of dendrites is close to v term of the monomers from Kajikawa (1972, K72, Figure 3c). The plate monomers in K72 are reported with a similar v term as the aggregates of plates, needles, and Mix1 (which all have similar values). Note that v term of plates and dendrites from K72 and v term of all aggregates simulated in this study (except for the aggregates of columns and "Mix2") are considerably smaller than v term of the aggregates from the PIP-CARE data set and LH74. The . Particle properties of simulated aggregates from this study (green and black), from previous studies (Kajikawa, 1972;Locatelli & Hobbs, 1974;Mitchell, 1996) (M96, LH74, and K72) and measurements of ice particle observed by PIP at the CARE site (brown, see text). (a) m versus D max ; (b) A versus D max ; (c) median (and 25th and 75th percentile for PIP CARE) of v term versus side projected maximum dimension D max,side for data from this study and versus the size definition of the respective study (v term is directly observed in K72 and LH74 and calculated with B92 from the m−D max and A−D max relations of M96) (d) same as (c) but for larger sizes. Note that K72 observations are for single monomers.
observations from LH74 are within the 25th and 75th percentile of the PIP-CARE data set. The median of v term of the simulated aggregates of this study increases faster with size compared to the in situ observations at sizes of several millimeters ( Figure 3d). Only v term of the mixture of small columns and large dendrites ("Mix2") has a comparably low slope. Potential reasons for this mismatch are limitations of the observations at these sizes (Brandes et al., 2008), turbulence affecting the observations (Garrett & Yuter, 2014), missing processes in the aggregation model (e.g., depositional growth on aggregates), imperfect parameterizations in the hydrodynamic model, or the dominance of monomer type mixtures in the aggregates.
Figures 3c and 3d also show v term calculated with B92 and the m−D max and A−D max relations from M96 (which did not measure v term directly). The simulated slope of v term from M96 observed aggregates is similar to the one simulated in this study, while the absolute value is slightly higher.
At sizes larger than about 5 mm, the simulated and the observed v term reach a saturation value close to 1 m s −1 . The median of v term of most simulated aggregates lies within the 25th and 75th percentile in the subcentimeter range, except the aggregates with the most extreme density (aggregate of dendrites and aggregates of columns). Thus, based on this comparison, these aggregates can be considered most representative for many aggregates found in the atmosphere.

Parameterization of Particle Properties
The relationships between hydrometeor properties such as mass, size, projected area, and velocity are key components in any ice microphysics scheme, and they strongly influence various microphysical processes (e.g., sedimentation, depositional growth, aggregation, or riming). Different microphysics schemes require a more or less simplified parameterization of particle properties. To address these different needs, we derive in this section fits for m and A as a function of D max and N mono that can be used in microphysics schemes, which can predict m and N mono given a certain D max (section 4.2). Of course, most bulk schemes require less detailed fits, and hence, we also derive fits of m, A, and v term as a function of D max or the mass-equivalent diameter D eq . This also allows us to assess the potential error of the less detailed fits (section 4.5), while their impact on modeled processes is studied later in section 5.

Fitting Approach for Monomer Number-Dependent Particle Properties
The particle properties of the monomers are defined a priori in the aggregation model and based on well-established observations. In contrast, the aggregate properties are determined by the aggregation process and change with increasing N mono . As we are particularly interested in quantifying how key particle properties of aggregates differ from the properties of the same-sized monomers, we normalize the aggregate properties by the property of a monomer with the same D max p represents the particle properties (mass or area), p(D max ,N mono =1) is the property of single monomers (given in Table 1), and f p is the normalizing function. A normalizing function which is larger (smaller) than 1 indicates that the aggregate properties are larger (smaller) than its composing monomer with the same size ( Figure 4).
To fit f p to various monomer types, we parameterize f p by a power law and express the coefficients by rational functions to fit the dependency on N mono similar to the approach presented in Frick et al. (2013).
The coefficients of f p for all monomer types can be found in Table 3. Note that we excluded the mixture of monomer types from the monomer-dependent analysis because our normalization approach cannot be applied to monomer mixtures.

Dependence of Aggregate Mass and Area on Monomer Number
Motivated by the common classification of unrimed ice hydrometeors in cloud ice and snow in many bulk schemes, we will investigate in this section how mass and area change when building up an aggregate with an increasing number of monomers. In particular, we want to explore whether the properties change smoothly with monomer number or whether they show any sharp transition at certain monomer numbers.
When we compare the mass of an aggregate with the mass of its monomer of the same size, we find in some conditions the aggregate to be heavier or lighter than the monomer. The relevant mechanisms which explain this behavior are illustrated in Figure 4 for aggregates of plates. Note that we assume for simplicity a monodisperse monomer distribution in Figure 4. When we consider pure depositional growth, we obtain a specific m−D max relation for each monomer type (Table 1; black line in Figure 4). One extreme aggregation scenario, which leads to the maximal size of an aggregate with a given number of monomers (which in this simplified case of a monodisperse distribution also determines its mass), would be if we assume that all monomers align along their maximum dimension. Clearly, the resulting aggregate would have a smaller m than a monomer of the same size. Of course, this maximal elongated assemblage of monomers is rather unlikely, and thus, the aggregate will have a more compact structure. If we imagine rearranging the monomers inside the aggregate in a progressively more packed configuration (indicated by the horizontal arrow in Figure 4), we might be able to reach the point where the size of the aggregate equals the one of the equal-mass monomer. At this point, it might be even possible to pack the monomers in a way that their size is smaller than an equal-mass monomer. A simple example of such an extreme packing would be to stack a number of plates on top of each other, that is, along their smallest axis. Whether an aggregate can be smaller than an equal-mass monomer is of course also dependent on how close the monomer m−D max relation is to the theoretical maximum packing of an equal-mass sphere.
The dependency of A on N mono can be understood analogously. Also, for A, the maximal elongated assemblage of the monomers leads to a lower A of the aggregate compared to the monomer of the same size, but in reality, the monomers will assemble in a more compact way. In addition, we have to consider that A is not simply additive as it is the case for m. Overlap (in the horizontally projected plane) and nonhorizontal alignment of the constituting monomers lead to a smaller A than the sum of A of the constituting monomers. Based on these simplified considerations, it becomes clear that the dependency of m and A on N mono is determined by the exponent of the monomer power laws and the overall "compactness" of the aggregates.
When considering the monomer dependence of all simulated aggregates, we find the most different behavior for plate and needle aggregates. For plate aggregates, m and A steadily decrease for a given D max with an increasing number of monomers (Figures 5b and 5d). From the principal considerations discussed in Figure 4, this behavior can be well understood. The plate monomers have the largest exponent (b m,1 = 2.48) of all monomers (Table 1), while the monomers itself show relatively loose connections within the aggregate (Figures 2a-2c). Interestingly, the aggregate mass for very small N mono can be slightly larger than the

Journal of Advances in Modeling Earth Systems
equal-size monomer, while A is immediately decreasing for N mono > 1. This effect can be easily understood when considering, for example, two plates that connect in a 90°angle of their major axes.
An opposite behavior is found for needle aggregates (Figures 6b and 6d). With increasing N mono , both m and A of the aggregates become larger than the equal-size monomers. In contrast to plates, the needle monomers have the lowest exponents for the m and A power laws (Table 1). The aggregates of the more one-dimensional needles also show a more compact packing.
The deviation of the particle properties of the individual simulated particles from the fit is characterized by the mean absolute error (Table A2), which is smallest for plates (0.1190 for f m and 0.0816 for f A ) and largest for needles (0.3737 for f m and 0.3926 for f A ). The mean absolute error also shows that the monomer number-dependent fit is superior to the more simple power law fit (section 4.4) when there is a substantial dependence of the particle property on N mono .
Dendrite and column aggregates have been analyzed similarly (according figures can be found in Supporting Information S1). The dendrites are similar to plates, while the columns are similar to needles. However, for all aggregate types, we find on average a relatively smooth transition of m and A when changing from single monomers to aggregates. For these two particle properties, we are unable to identify a "jump" due to the onset of aggregation. The next sections will show whether this behavior will change when deriving terminal velocity from m and A.

Dependence of Terminal Velocity on Monomer Number
The terminal velocity for all aggregates was calculated with the hydrodynamic model of B92 (section 2.2). In Figure 7a, v term is shown as a function of D max for plate aggregates. Note that the fits have been derived by applying B92 to the m−D max and A−D max fits (Table 3) rather than fitting them directly to the cloud of individual v term . In this way, we are consistent with the way how v term relations are usually connected to m−D max in bulk schemes. The terminal velocity of plate aggregates steadily decreases with increasing N mono . This dependency is much less pronounced at small D max as compared to the largest sizes. However, it should be noted that the fits for very small monomer numbers are probably unrealistic for large D max as we do not

10.1029/2020MS002066
Journal of Advances in Modeling Earth Systems expect aggregates of centimeter sizes to be composed of only a few large plates. In fact, the here used geometrical relations for the plate monomers are only valid up to a size of 3 mm (Pruppacher & Klett, 1998).
We find a similar decreasing v term with increasing N mono for dendrites (see supporting information S1). As we might expect from the different change of m and A with N mono seen in Figure 7a, also the behavior of v term with increasing N mono is different for needles (Figure 7). Needle aggregates seem to fall slightly faster when their monomer number increases. Interestingly, all aggregates reveal a very low dependence of v term on monomer number at small sizes which is in contrast to assumptions in some microphysics schemes that distinguish between monomers and aggregates (e.g., Seifert & Beheng, 2006;Tsai & Chen, 2020). Besides, all aggregates reveal a saturation of v term at large (centimeter) sizes which is in good agreement with   (Table 3) fits rather than fitting them directly to the cloud of individual v term .

Journal of Advances in Modeling Earth Systems
observations ( Figure 1). However, the absolute value of the saturation vterm ranges from 0.8 to 1.6 m s −1 depending on the monomer type.
Because v term of monomers and aggregates is converging toward the same value at small sizes (Figure 7), we can use the previously derived scaling relation (Equations 6 and 7) to relate the dependency of v term on N mono to the exponents b m and b A of the monomers (b m,1 and b A,1 ) and aggregates (b m,agg and b A,agg ) in the m−D max relation. Starting from a similar value of v term at small sizes, v term of an average aggregate increases slower than v term of a monomer if s monodep =b m,agg −b m,1 −0.25(b A,agg −b A,1 )<0 (cf. Equations 6 and 7). As a result, at larger sizes, v term of the aggregate is lower than v term of the monomer. In an analog way, v term of an aggregate is larger than v term of the monomer if s monodep >0. As b m,agg and b A,agg are similar for all aggregates (Table 4), the sign of v term with increasing N mono depends mainly on b m,1 and b A,1 . For plates and needles, s monodep equals −0.21 and 0.12, respectively.
How the particle properties change with increasing N mono as well as the absolute values of calculated v term depends on the choice of the hydrodynamic model. Finding the optimal formulation of hydrodynamic models for ice and snow particles is still an active field of research (Nettesheim & Wang, 2018;Westbrook & Sephton, 2017) and outside the scope of this study. In Appendix A2, we tested the sensitivity of the results to the choice of the hydrodynamic model for plate aggregates. HW10 seems to yield overall similar results to

Journal of Advances in Modeling Earth Systems
B92 except for the saturation at large diameters. Due to the absence of the turbulence correction in HW10, v term increases also at large diameters. For KC05, the monomer dependence is much weaker. However, all hydrodynamic models show an overall small monomer dependence at small particle sizes.
It has also been observed (e.g., Garrett & Yuter, 2014) that tumbling of particles caused for example by turbulence might decrease the effective projected area and therefore increase v term . We also tested the sensitivity of our results to different degrees of tumbling (section A0.2.2). As expected, the effect of tumbling is largest for single crystals (due to their more extreme aspect ratio) but strongly decreases for aggregates. Certainly, for aggregates, the choice of the hydrodynamic model has a larger effect of v term than different assumptions on particle tumbling.

Mean Particle Properties of Monomers and Aggregates of Different Monomer Types
The relatively continuous change of particle properties with N mono found in the last section justifies a simplified fit, which is also necessary for implementing the results into common bulk microphysics schemes. These schemes often only contain two classes for unrimed ice particles, usually denoted as cloud ice (monomers) and snow (aggregates). Table 1. Similar fits of m and v term to aggregates of any monomer number larger than 1 are shown in Figures 8c and 8d; the fit coefficients can be found in Table 4.

Figures 8a and 8b show the derived m−D max relations for single monomers (N mono =1) and the derived v term based on the m−D max and A−D max relations summarized in
The m−D max relations for monomers show a larger spread especially for larger sizes as compared to the aggregates. This is expected considering that the exponents for monomers range between 1.89 and 2.48 (Table 1), while the exponents for aggregates are between 1.95 and 2.22 (Table 4). The values for aggregates agree well with theoretical aggregation studies (Westbrook et al., 2004b) as well as in situ observations (section 3.1). Despite the similar exponent, the effective density of the aggregates varies considerably (compare m at a given size in Figure 8c), which is in agreement with previous studies (Dunnavan et al., 2019;Hashino & Tripoli, 2011), even though their approaches to simulate aggregates are very different from the approach used in this study. Aggregates of columns exhibit the highest density, while aggregates of dendrites show the lowest density.
The differences in the m−D max relation are linked to the resulting v term −D max relation (Figures 8c and 8d). At D max = 5 mm, the v term of different monomers spread nearly 1 m s −1 . The differences are in general smaller for aggregates. Interestingly, most aggregate types reveal very similar v term . The main exceptions are dendrite aggregates with the slowest and column aggregates with the fastest v term . v term of the Mix2 aggregates increases slower with increasing D max compared to the other aggregates.
Similar to the previous monomer number dependent fits, also the "two-category" fits show similar v term at small sizes. The monomer type appears to have in general a much larger impact on v term than the classification into certain N mono regimes.

Power-Law and Atlas-Type Fits for Terminal Velocity
Power-law fits for m, A, and v term are commonly used in bulk schemes. Especially for v term , the power law introduces inconsistencies with observations because a saturation value for v term as observed for raindrops or snowflakes cannot be represented. Instead of using standard power laws in the form, and the two fit parameters a vDmax and b vD max , Atlas et al. (1973) proposed a three-parameter (α Deq , β Deq , and γ Deq ) formulation

Journal of Advances in Modeling Earth Systems
Formulating this "Atlas-type" fit with the mass equivalent diameter D eq instead of D max has been found to provide optimal fit quality for snow aggregates (Seifert et al., 2014). For small (large) values of D eq , v term approaches α Deq − β Deq (α Deq ). With increasing values of γ, the transition from small to larger values of v term is shifted toward larger values of D eq . Approximations for bulk collision rates based on Atlas-type fits can be found in Seifert et al. (2014) which makes them usable in bulk microphysics schemes without the necessity of look-up tables.
Power-law and Atlas-type relations have been applied to the various aggregates and the fit coefficients are summarized in Table 5. For the fitting, we did not use v term of the individual particles but directly applied to fit to v term derived with B92 and the existing m−D max and A−D max relations.
In Figure 9, the different fits are compared for plate monomers and their aggregates. Note that the saturation region (D max > 1 cm) has been excluded for the power-law fits. It can be seen in Figure 9b that the Atlas-type fit is very close to the theoretical line calculated with B92 and the m−D max and A−D max relations. The power-law fits (Figure 9a) provide only a close fit to the theoretical values at the smaller size range. Between 300 μm and 4 mm, they cause a slight underestimation, while at larger sizes, they increasingly overestimate v term . Similar fits have been derived for all aggregate types (Table 5, figures for other monomer types similar to Figure 9 can be found in the supporting information S1).
When we compare the calculated v term with some widely used microphysics schemes (Figure 9c), we find most schemes to overestimate v term at small sizes (except of the cloud ice category in Morrison et al., 2005).
The absolute values for v term at small sizes are strongly dependent on monomer type, and hence, additional constraints should be provided by additional observations. However, the aggregation model shows independent on the monomer type that at submillimeter sizes, there should be no strong "jump" in v term between ice particles and small aggregates. Also, in the centimeter-size range, models using a power-law formulation are strongly overestimating v term for all aggregate types.

Application and Sensitivity Tests in the Lagrangian Particle Model McSnow
In this section, we will test the possible impact of implementing particle properties with different amount of complexity (monomer number dependence) or different fitting functions (power law vs. Atlas type) on the simulation of sedimentation, aggregation, and depositional growth. For this, we use a one-dimensional setup of the Lagrangian particle model McSnow (Brdar & Seifert, 2018), which provides the flexibility to implement the different particle property formulations.
For simplicity, only sedimentation, depositional growth, and aggregation are considered in our simulations. Aggregation is calculated with a Monte-Carlo algorithm following Shima et al. (2009), and the sticking efficiency of Connolly et al. (2012) is used. McSnow is based on the Lagrangian super particle approach (Shima et al., 2009), which allows deriving not only the particle mass and its multiplicity X mult but also predicts the number of monomers the particle is composed of. This information is key to test the N mono dependent particle relations. The thermodynamic profiles and the overall setup is similar to previous simulation studies Power-law (dashed-dotted, a) and Atlas-type fits (dashed-dotted, b) have been applied to the directly calculated v term (solid lines) rather than the individual points. (c) v term used in microphysics schemes (same as in Figure 1b).

Journal of Advances in Modeling Earth Systems
with McSnow in Brdar and Seifert (2018) and Seifert et al. (2019). Particles are initialized at the upper boundary of the 5-km thick domain with a mass flux of F m =2·10 −5 kg s −1 and a mean mass of the particle size distribution of m mean = 2·10 −10 kg. The initial ice particles follow a generalized gamma distribution of particle mass with a shape parameter of 0 and a dispersion parameter of 1/3 (following Eq. 9 in Khain et al., 2015). The temperature decreases linearly from 273.1 K at z = 0 km to 242.2 K at z = 5 km. The supersaturation over ice is held constant at 5% with respect to ice in the whole column and is not consumed by the growth of the particle. The simulations are performed with 250 vertical levels, which result in a vertical resolution of 20 m. The model time step is set to 5 s, and the initial multiplicity is chosen to be 1,000. The simulations are run for 10 hr, from which the last 5 hr are averaged in 10-min intervals to reduce noise in the analyzed profiles.  Table 3), the separation between N mono =1 and N mono > 1 ("binary"; Tables 1 and 4) and single relation (the one fitted to all aggregates) for all N mono ("constant"; Table 4) for plates. For each individual super particle, B92 is used directly to calculate v term . Shown are height profiles of (a and b) number flux F N , (c and d) mass flux F m , and (e and f) mean mass m mean . The particles are categorized into N mono = 1 (left) and N mono > 1 (right).
In the following, we will focus the comparison on particle number flux (F N ), mass flux (F M ), and mean mass m mean (which is the ratio between the integrated mass density q m and the integrated number density q N ).
In the first simulation experiment shown in Figure 10, we include particle properties for which the full N mono dependence is taken into account (Table 6). In this setup, we call hereafter the control simulation ("CTRL"). Profiles are separated into single monomers (N mono = 1) and aggregates (N mono > 1) to better distinguish the effects on what we might define as "cloud ice" and "snow" category in a bulk scheme. This separation might be important considering that there can be cases of weak aggregation. With weak aggregation, most of the particles will remain monomers, and thus, it is especially important to match profiles of these particles accurately.
In general, aggregation decreases the number flux (F N ), while the increase in the mass flux (F m ) is due to depositional growth. The mass flux of aggregates increases also due to conversion from monomers to aggregates by aggregation. The combination of both processes is causing m mean to continuously increase toward the surface. Aggregation rates in McSnow are proportional to the gravitational collection kernel (Eq. 21 in Brdar & Seifert, 2018). Thus, the probability of collision for two particles is high if they have strongly different v term and if the sum of their cross-sectional areas is large. F N of the monomers (N mono = 1) decreases monotonically with decreasing height because the monomers are converted into aggregates (N mono > 1) by the aggregation process and there is no source of monomers like nucleation considered (Figure 10a). This decrease of F N (and increase of m mean ) is especially strong at heights between 2 and 3 km. This region of enhanced aggregation is found at heights where the temperature is close to −15°C where the sticking efficiency has a local maximum. As a result of the conversion of monomers to aggregates, F N of the aggregates increases at heights higher than about 3.5 km (Figure 10b). At lower heights, the number of aggregate-aggregate collisions outweighs the number of monomer-monomer collisions, and thus, F N of the aggregates decreases.
The signature of the conversion from monomers to aggregates is also seen in F m of the monomers (Figure 10c). Especially in the region of enhanced aggregation, this leads to a strong decrease of F m . In the heights above this region, depositional growth outweighs the loss of mass of the monomers to the aggregates, and thus, there is an increase of F m with decreasing height. F m of the aggregates increases monotonously due to both depositional growth of the aggregates and conversion from monomers to Note. The second column specifies the monomer type from which the m−D max and A−D max (and v term −D max for the Atlas and power law run) fit is taken. The third column denotes how the N mono dependency is represented. f p (N mono ,D max ) is the normalizing function with full N mono dependence (section 4.1), f p (N mono =1;N mono >1,D max ) denotes only a binary seperation in N mono =1 and N mono >1, and f p =f (N mono >1,D max ) indicates that the fit for all aggregates N mono >1 is taken for all particles (section 4.4). The fourth column indicates whether v term is calculated using B92 or with a parameterized relation of v term −D max (section 4.5). The fifth column shows the precipitation rate (F m (z = 0 m)) and in brackets its deviation from the "CTRL" run. The last column lists the mean mass m mean at the surface and the ratio of m mean between the simulation and its "CTRL" run (in brackets).
aggregates (Figure 10d). In this setup, the change of F m with height is governed by v term and q N at a given height. For example, a combination of low v term and high q N at upper layers leads to a large increase in F m . Simply speaking, a large number of slow falling ice crystals can grow efficiently by deposition which increases F m .

Sensitivity to Representation of Monomer Number Dependency
The "CTRL" simulation is now compared to simulations with a binary separation into single-monomer particles and aggregates of any monomer number larger than 1 (binary). An additional simulation is performed with no monomer number dependence (constant). Here, the particle properties that were fitted to the mean of all aggregates are used for all particles. All simulations are done for plate and needle monomers and aggregates because we found the monomer dependence to be most pronounced for these monomer types. For the other monomer types, the effect of N mono can be expected to be smaller.
The most apparent difference between the simulations with different representations of the N mono dependencies for plate monomers and aggregates of plates is the faster decrease of F N and F m and faster increase of m mean of the monomers (N mono = 1) in the "constant" simulation ( Figure 10). A slightly

10.1029/2020MS002066
Journal of Advances in Modeling Earth Systems faster decrease of F N (faster increase of m mean ) for aggregates (N mono > 1) with decreasing height can be seen for both the "binary" and the "CTRL" simulation. However, all of the simulations show very similar profiles. Figure 11 shows the same experiment as Figure 10 but using the parameterizations for needles instead of plates. Also, for needles, the most remarkable difference between the simulations is the difference between the "constant" and the "CTRL" run (Figures 11a and 11e). Also, aggregate-aggregate collections are less effective in the "CTRL" run (Figures 11b and 11f). Note that all monomers have been depleted by aggregation at a height of about 1,000 m, and thereby, m mean is not defined below.
Overall, the differences of m mean at the ground of the total ice particle population are small (factor of 1.2 and 1.4 higher m mean for the "binary" and "constant" simulation for plates and factor of 0.8 lower m mean for the "binary" and "constant" simulation for needles, Table 6).  Tables 1 and 4) and power law (without ("powerlaw") and with imposing an upper limit on v term ("powerlawlimit"), which is consistent with the saturation value of the Atlas-type relation) and Atlas-type v term −D max relations for plate monomers and aggregates of plates (see Table 5). Overlayed is the CTRL/monodep simulation in gray (see also Figure 10). Shown are height profiles of (a and b) number flux F N , (c and d) mass flux F m , and (e and f) mean mass m mean . The particles are categorized into N mono = 1 (left) and N mono > 1 (right).
Also, the differences in the precipitation rates (F m ) are small (less than 5%; see Table 6). These small differences are due to the small difference of the absolute value of v term at small sizes ( Figure 7) and q N at upper heights, which lead to a similar mass uptake ( Figure 10). However, the precipitation rate between the "Plate CTRL" simulation and the "Needle CTRL" simulation is relatively large (Table 6), which might be due to the strongly different v term of the monomers.
The N mono dependency is even weaker for aggregates composed of other monomer types (sections 4.2 and 4.3). In summary, the simulation experiments with different monomer dependency indicate that a binary separation between single monomers and aggregates performs similarly well as relations which take into account a more detailed monomer dependency. Some but still small differences are found if no monomer dependency is taken into account, that is, a single ice class for all monomer numbers is assumed. In our simulation, particles with low N mono are only prevalent at cold temperatures, where aggregation is less important due to the low sticking efficiency. Additional simulations (shown in the supporting information S1) with lower F m and therefore weaker aggregation show that the "binary" simulations stay very close to the "CTRL" simulation, while the "monodep" simulations show considerably larger deviations. Hence, we find that the classical separation in cloud ice (monomers) and snow (aggregates) is sufficient for the aspects of monomer number-dependent particle properties.

Sensitivity to the Parameterization of Terminal Velocity
In this section, we test the sensitivity of the simulations to different implementations of the v term −D max formulation. In Figure 12, v term of plate monomers and aggregates is parameterized either as power-law or Atlas-type fit.
As we saw in Figure 9, the power-law and Atlas-type fits match very closely at small particle sizes. This explains the very close matching of the three simulations in the upper part of the simulated profiles ( Figure 12) where the PSD is dominated by small particles. As soon as the aggregation becomes stronger (below ca. 3 km), F N in the simulations using the power law ( Figure 12b) is much lower than for Atlas-type. The decreasing number flux of aggregates with lower height (Figure 12b) also indicates that especially the self-collection of aggregates is stronger than for Atlas type. In the same height region, the mean mass of the aggregates (Figure 12f) is strongly increased for the power law (factor of 5). Instead of using an Atlas-type fit to consider the saturation of the terminal velocity, one can also think of imposing an upper limit on v term in the power law relation. In the "Powerlawlimit" simulation, we chose the saturation value of the Atlas-type fit (α Deq ) as an upper limit. This limit does not only affect the sedimentation but also all processes which depend on v term (e.g., aggregation). In this way, the overestimation of m mean , caused by an unlimited increase of v term , can indeed be prevented, but the height profile of F N and m mean is not as well matched as with the Atlas-type approximation. As expected, the continuously increasing v term in the unlimited power law leads to much stronger growth of aggregates as compared to relations which include the saturation velocity at large particle sizes. This is an interesting finding and could be one reason for the overestimation of radar reflectivities found at lower layers in ice clouds simulated with the Seifert-Beheng scheme (Heinze et al., 2017).
Although m mean of the aggregates is much larger for the power law, the difference to the Atlas type in precipitation rates is very small (smaller than 5%; Figure 12d and Table 6). Note that in more realistic cases, as for example in presence of stronger sublimation layers, the difference in m mean can induce larger differences in the precipitation rate because larger particles can fall through a thicker layer of subsaturated air before they sublimate completely.

Summary and Conclusions
In this study, we generated a large ensemble of ice aggregates (ca. 105,000 particles) using an aggregation model and hydrodynamic theory to study the change of particle properties such as mass m, projected area A, and terminal velocity v term as a function of monomer number N mono and size. The aggregates were composed of various monomers types (plates, dendrites, needles, and columns), monomer sizes, and monomer numbers. In order to test the impact of habit mixtures, we also included in our analysis two different mixtures of dendrites and columns. The choice of mixing specifically dendrites and columns was motivated by in situ observations of the composing monomers in large aggregates sampled on the ground (Lawson et al., 1998).
When comparing our aggregate properties with in situ observations, we find m and A to be very similar to the results presented in Mitchell (1996) but the slope of our m−D max relations is larger than the slope given in Locatelli and Hobbs (1974). A better agreement with Locatelli and Hobbs (1974) and also with theoretical considerations in Westbrook et al. (2004b) is reached for mixtures of small columns and larger dendrites (Mix2). Interestingly, this monomer mixture also achieves the best agreement with observed v term −D max relations. Considering the large spread in the observations (Figure 3), we can overall conclude that our aggregate ensemble matches the observed range of variability and does not show any substantial bias.
Our synthetic aggregate ensemble allowed us to investigate the transition of particle properties from single crystals to aggregates with increasing number of monomers in a level of detail which is currently unavailable from in situ observations. For m and A as a function of size, we find the relations to change rather smoothly with increasing N mono . The differences introduced by the choice of the monomer type are found to be overall larger than due to the number of monomers. We find the exponents in the A−D max and m−D max relations of the monomers to be closely connected to the resulting change with N mono .
The derived A−D max and m−D max relations including the monomer type and number dependence were then used to calculate v term −D max relations. Again, we find a rather smooth transition from single crystals to aggregates rather than a "jump" as found in several microphysics schemes (Figure 1b). For small sizes below a few millimeters, our results suggests that the "ice" and "snow" category of microphysics schemes should have similar properties. At larger sizes, the aggregates v term are found to deviate more from the monomers. Again, the monomer type is found to have a larger impact than the monomer number. Aggregates of plates tend to be faster, while aggregates of needles are slower than the equal-size monomer. In accordance to in situ observations, our simulations reveal for all aggregate types a saturation of v term at centimeter sizes. However, the saturation value varies for the different aggregate types from 0.8 to 1.6 m s −1 .
In order to potentially implement our results in microphysics schemes, we derived two-parameter powerlaw fits and three-parameter Atlas-type fits for single monomers (N mono = 1) and aggregates (N mono > 1) representing the commonly used ice and snow classes in models. The new power-law fits match the small sizes well and avoid unrealistic "jumps" found in current schemes. However, the power laws are unable to represent the saturation of v term at larger sizes. The Atlas-type fits are found to match the entire size range well and should thus be considered to be implemented in ice microphysics schemes as they do not substantially increase the computational costs while strongly improving the realism of the relations.
We finally tested the impact of implementing monomer dependence, habit type, and velocity fitting method on idealized aggregation simulations. For this, we used a new 1D Lagrangian Monte Carlo model which allowed us to implement the derived relations with different degree of complexity. The simulations experiments revealed that there is only a very small impact of using a relation of only two monomer categories (single particle and aggregate) as compared to a continuous monomer number dependence. A single category which does not take any monomer number into account shows slightly larger deviations, but the variability due to monomer type is in general larger than the impact of monomer number.
In a second simulation experiment, we investigated the impact of using a power law or an Atlas-type fit for v term . The simulations show very small differences in the upper part of the cloud where the profiles are dominated by small particles which are fitted similarly well with the two relations. Once aggregation becomes more dominant and the spread of particles sizes shifts to larger sizes, the simulations using the power law lead to a much stronger aggregation and in particular stronger self-aggregation of particles as compared to the Atlas-type fit. The impact of the widely used power-law relations for v term should thus be further studied for bulk schemes as it seems to be likely that they might cause an overestimation of aggregation and snow particle sizes.
We also shortly investigated the sensitivity of our derived relations to particle tumbling and the choice of the hydrodynamic theory. While tumbling can significantly affect the properties of single monomers, it has a surprisingly small effect on our results for the aggregates. The choice of the hydrodynamic theory is a larger source of uncertainty which should be further investigated in future studies. It seems to be important in the future to better constraint the composition of aggregates regarding the monomer type. This question could be approached by improved in situ techniques but also with detailed models that allow to predict the particle habit such as presented in ( area ratio using different parametrizations of the hydrodynamic theory (Bohm, 1989;Heymsfield & Westbrook, 2010;Mitchell & Heymsfield, 2005). For each snowfall event, each of these parameterizations is calculated and the one which minimizes the error in the estimate of the liquid water equivalent precipitation with respect to the precipitation gauge is selected for that event. This procedure and the related uncertainties are described more in detail in von Lerber et al. (2017). Additionally, observations during 5-min intervals, where the mean horizontal wind speed exceeds 4 m s −1 , are excluded to reduce turbulence effects (similar to Brandes et al., 2008).
After applying these filters, the data set, which covers the winters from 2014 to 2017 with 48 snowfall events, contains about 4.3 million ice particles. It should be noted that PIP is providing a measurement of the ensemble of particles and no particle by particle-based classification is performed. Hence, the measurement volume includes mixtures of different habits.

A.2. Sensitivity of the Terminal Velocity to the Hydrodynamic Model and Tumbling A.2.1. Hydrodynamic Models
As mentioned in section 2.2, the hydrodynamic models of B92, KC05, and HW10 differ in several aspects. The Re(X) relation requires assumptions about particle surface roughness, which are differently implemented in the models. Also, the definition of X is different (Table A1). While in B92 X is proportional to mD 0:5 max A −0:25 , X is proportional to mD max A −0.5 in HW10 and mD 2 max A −1 in KC05. As a result in B92 and HW10, v term increases slower with decreasing area ratio (A r =4Aπ −1 D −2 ) than in the formulation of KC05. The empirical correction of X due to wake turbulence is also applied in KC05 but not in HW10.
These differences affect the behavior of v term at large sizes and the monomer number dependency (which we quantify by s monodep ). Without the empirical correction of X (which considers wake turbulence), v term only saturates if v term,Re ≫ 1 ∼D 0 . For example with HW10, the saturation would be reached for b m −0.5b A −1=0 (Table A1). This is, for example, not the case for aggregates of plates simulated in this study, and therefore, HW10 does not predict a saturation of v term at larger sizes ( Figure A1a). Also, the sign and the strength of the increase/decrease of v term with increasing N mono depends on the formulation of X. In section 4.3, we introduced s monodep as a measure for this monomer number dependency. Applying this measure to the aggregates of plates yields s monodep = −0.21 for HW10 and s monodep = −0.06 for KC05. Both HW10 and KC05 show the decrease of v term with increasing N mono which we saw when using B92, but this decrease is very weak for KC05.

A.2.2. Tumbling
To investigate the effect of the tumbling of the aggregates (as reported, e.g., by Garrett & Yuter, 2014) on the projected area A and v term , the particles are tilted with a standard deviation of 0°, 20°, 40°, and 60°, around the principal axis ( Figure A2). This is done only after the final aggregate is assembled and thereby does not influence the structure of the aggregates. This rotation reduces A, and in turn, v term increases.
The monomers (top panel in Figure A2) are stronger effected by tumbling (especially at large D max ) due to their lower aspect ratio (not shown). The largest increase in v term with increasing tumbling is found for KC05 due to the largest increase in the Best number with decreasing A (see section 2.2). B92 shows the least influence of tumbling, which increases v term at maximum by about 0.1 m s −1 and has a negligible effect on v term for the aggregates. Figure A2. v term based on m/A−D fits (Tables 1 and 4) and different hydrodynamic models. The particles are horizontally aligned ("no tumbling") rotated by 20°or 40°around the principal axis to mimic different strength of tumbling. (a) plate monomers; (b) aggregates of plates.

A.3. Mean Absolute Error of the Mass and Area-Size Relations
In sections 4.2 and 4.4, we provided fit relations for mass and area with and without taking into account the monomer number dependency of the simulated aggregates. The mean absolute error of the fits shown in Tables 3 and 4 (normalized by the properties of the monomers; e.g., shown for plates with the green dotted lines in Figure 5) is shown in Table A2.