the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Statistical properties of a stochastic model of eddy hopping
Takeshi Watanabe
Toshiyuki Gotoh
Statistical properties are investigated for the stochastic model of eddy hopping, which is a novel cloud microphysical model that accounts for the effect of the supersaturation fluctuation at unresolved scales on the growth of cloud droplets and on spectral broadening. Two versions of the model, the original version by Grabowski and Abade (2017) and the second version by Abade et al. (2018), are considered and validated against the reference data taken from direct numerical simulations and largeeddy simulations (LESs). It is shown that the original version fails to reproduce a proper scaling for a certain range of parameters, resulting in a deviation of the model prediction from the reference data, while the second version successfully reproduces the proper scaling. In addition, a possible simplification of the model is discussed, which reduces the number of model variables while keeping the statistical properties almost unchanged in the typical parameter range for the model implementation in the LES Lagrangian cloud model.
The purpose of the present paper is to investigate the statistical properties of the stochastic model of eddy hopping proposed by Grabowski and Abade (2017). This stochastic model, referred to hereinafter as the eddyhopping model, was developed in order to account for the effect of the supersaturation fluctuation at unresolved (subgrid) scales on the growth of cloud droplets by the condensation process. In a turbulent cloud, cloud droplets arriving at a given location follow different trajectories and thus experience different growth histories, which leads to significant spectral broadening. This mechanism, referred to as the stochastic condensation theory, has been investigated since the early 1960s by a number of researchers (mostly Russian; see Sedunov, 1974; Clark and Hall, 1979; Korolev and Mazin, 2003), but the importance of this mechanism was later reinforced by Cooper (1989) and LasherTrapp et al. (2005). For this mechanism, Grabowski and Wang (2013) emphasized the importance of largescale eddies (turbulent eddies with scales not much smaller than the cloud itself) and proposed the concept of largeeddy hopping. Grabowski and Abade (2017) formulated this concept and developed the eddyhopping model. Abade et al. (2018) extended the model by introducing a term accounting for the relaxation of supersaturation fluctuations due to turbulent mixing. For clarity, we hereinafter refer to the model by Grabowski and Abade (2017) as the original version, and the model by Abade et al. (2018) as the second version. For the following study using the eddyhopping model, readers are referred to Thomas et al. (2020).
It should be noted that the turbulent entrainment mixing is another important mechanism for the supersaturation fluctuation generation other than the stochastic condensation and that the effects of the turbulent entrainment mixing are not included in the eddyhopping model considered in the present study. Abade et al. (2018) investigated the effects of the turbulent entrainment mixing and entrained CCN activation by using the entraining parcel model.
In the present paper, we take a rather theoretical approach to obtain various statistical properties of the eddyhopping model, such as the variance, covariance, and autocorrelation function of the supersaturation fluctuation. These statistical properties are used to validate the model against the reference data taken from direct numerical simulations (DNSs) and largeeddy simulations (LESs). We show that the original version of the eddyhopping model fails to reproduce a proper scaling for a certain range of parameters, resulting in the deviation of the model prediction from the reference data, while the second version successfully reproduces the proper scaling. We show how the relaxation term introduced by Abade et al. (2018) leads to this improvement. We also discuss the possibility of simplification of the model, which reduces the number of model variables while keeping the statistical properties almost unchanged in the typical parameter range for the model implementation in the LES Lagrangian cloud model.
The remainder of the present paper is organized as follows. Section 2 describes the governing equations of the original version. Section 3 presents a theoretical analysis and numerical experiments and demonstrates the improper scaling in the model prediction by the original version. Section 4 describes the second version. Finally, Sect. 5 discusses the possibility of simplification of the model.
The original version of the eddyhopping model proposed by Grabowski and Abade (2017) consists of the following two evolution equations. First, the fluctuation of the vertical velocity of turbulent flow at the droplet position, w^{′}(t), is modeled by the Ornstein–Uhlenbeck process:
where δt is the time increment, ψ is a Gaussian random number with zero mean and unit variance drawn every time step, ${\mathit{\sigma}}_{{w}^{\prime}}$ is the standard deviation of w^{′}, and τ is the integral time, or the largeeddy turnover time of the turbulent flow. Here, ${\mathit{\sigma}}_{{w}^{\prime}}$ and τ are used as the model parameters. Second, the supersaturation fluctuation at the droplet position, S^{′}(t), is governed by
Here, the first term on the righthand side represents the effect of adiabatic cooling and warming due to air parcel ascent and descent caused by the vertical velocity w^{′}(t). The parameter a_{1} has the unit of a scalar gradient. The second term on the righthand side represents the effect of condensation or evaporation of droplets. The timescale τ_{relax} is referred to as the phase relaxation time and is inversely proportional to the average of the number density and radius of the droplets (Politovich and Cooper, 1988; Korolev and Mazin, 2003; Kostinski, 2009; Devenish et al., 2012).
Equation (1) can also be written as the following derivative form (Pope, 2000):
Here, the term ${F}_{{w}^{\prime}}\left(t\right)$ is statistically independent of S^{′} and obeys the Gaussian random process that has zero mean and twotime covariance defined by
where the angle brackets indicate an ensemble average and δ( ) is the Dirac delta function. In the following theoretical analysis, Eqs. (3) and (2) are used as the governing equations of the original version.
We now obtain the analytical expression for the standard deviation of the supersaturation fluctuation, ${\mathit{\sigma}}_{{S}^{\prime}}$, in a statistically steady state. Starting from Eqs. (3) and (2), the result is provided in Eq. (13).
First, multiplying Eq. (3) by S^{′} and taking an ensemble average, we obtain
because of statistical independence $\left(\langle {S}^{\prime}{F}_{{w}^{\prime}}\rangle =\mathrm{0}\right)$. Second, multiplying Eq. (2) by w^{′} and taking an ensemble average, we obtain
Summing Eqs. (5) and (6), we obtain
Next, we consider a statistically steady state. Since an ensembleaveraged variable does not change in time ($\mathrm{d}\u2329\circ \u232a/\mathrm{d}t=\mathrm{0}$) and $\langle {{w}^{\prime}}^{\mathrm{2}}\rangle ={\mathit{\sigma}}_{{w}^{\prime}}^{\mathrm{2}}$ in the statistically steady state, we obtain the flux of the supersaturation in the vertical direction as follows:
where Da is the Damköhler number (Shaw, 2003) defined as
Next, multiplying Eq. (2) by S^{′} and taking an ensemble average, we obtain
In the statistically steady state, we have
Combining Eqs. (8) and (11), we obtain
or equivalently,
Here, ${\mathit{\sigma}}_{{S}^{\prime}}$ in Eq. (13) has two important asymptotic forms, as shown below:

Large scale limit.
For τ→∞ (or equivalently, Da→∞, L→∞, where $L={\mathit{\sigma}}_{{w}^{\prime}}\mathit{\tau}$ is the integral scale), ${\mathit{\sigma}}_{{S}^{\prime}}$ in Eq. (13) is approximated as$$\begin{array}{ll}{\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}& {\displaystyle}\approx {a}_{\mathrm{1}}D{a}^{\mathrm{1}/\mathrm{2}}{\mathit{\tau}}_{\mathrm{relax}}^{\mathrm{1}/\mathrm{2}}{\mathit{\tau}}^{\mathrm{1}/\mathrm{2}}{\mathit{\sigma}}_{{w}^{\prime}}\\ \text{(14)}& {\displaystyle}& {\displaystyle}={a}_{\mathrm{1}}{\mathit{\tau}}_{\mathrm{relax}}{\mathit{\sigma}}_{{w}^{\prime}}.\end{array}$$For the case of a constant dissipation rate of turbulent kinetic energy ε, ${\mathit{\sigma}}_{{w}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}$ (see Appendix B), and we have the following scaling:
$$\begin{array}{}\text{(15)}& {\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}.\end{array}$$ 
Small scale limit.
For τ→0 (or equivalently, Da→0, L→0), ${\mathit{\sigma}}_{{S}^{\prime}}$ in Eq. (13) is approximated as$$\begin{array}{}\text{(16)}& {\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}\approx {a}_{\mathrm{1}}{\mathit{\tau}}_{\mathrm{relax}}^{\mathrm{1}/\mathrm{2}}{\mathit{\tau}}^{\mathrm{1}/\mathrm{2}}{\mathit{\sigma}}_{{w}^{\prime}}.\end{array}$$For the case of a constant dissipation rate of turbulent kinetic energy ε, ${\mathit{\sigma}}_{{w}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}$ and $\mathit{\tau}\sim {L}^{\mathrm{2}/\mathrm{3}}$ (see Appendix B), and we have the following scaling:
$$\begin{array}{}\text{(17)}& {\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{2}/\mathrm{3}}.\end{array}$$
The above asymptotic forms of ${\mathit{\sigma}}_{{S}^{\prime}}$ in the two limits can be validated through comparison with the result of the scaling argument by Lanotte et al. (2009). From their argument, we should have ${\mathit{\sigma}}_{{S}^{\prime}}\sim {a}_{\mathrm{1}}{\mathit{\tau}}_{\mathrm{relax}}{\mathit{\sigma}}_{{w}^{\prime}}$ for the large scale limit, which is consistent with Eq. (14). On the other hand, we should have ${\mathit{\sigma}}_{{S}^{\prime}}\sim {a}_{\mathrm{1}}\mathit{\tau}{\mathit{\sigma}}_{{w}^{\prime}}$ for the small scale limit, which is inconsistent with Eq. (16). Therefore, the original version given by Eqs. (3) and (2) does not reproduce the proper scaling for the small scale limit.
Figure 1 compares the scale dependence of ${\mathit{\sigma}}_{{S}^{\prime}}$ for the analytical expression given by Eq. (13) (orange curve) with the results of the numerical integration of the original version given by Eqs. (1) and (2) (blue squares). Here, numerical integration is conducted in the same manner as that by Thomas et al. (2020) (Sect. 5 of their study), except that the integration time is increased from 6τ to 10τ (see Appendix A for details). After the integration time of 10τ, all of the experimental results achieved a statistically steady state and agreed with the theoretical curve (compare the orange curve and the blue squares). As expected based on the analysis, the theoretical curve shows the scaling ${\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}$ for large scales (approximately L>10^{1}m) and the improper scaling ${\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{2}/\mathrm{3}}$ for small scales (approximately $L<{\mathrm{10}}^{\mathrm{1}}$m). These results are contrary to the results of DNSs and LESs (scaledup DNSs) conducted by Thomas et al. (2020) (black dots in Fig. 1), which show the proper scalings both for large and small scales (${\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}$ and ∼L^{1}, respectively).
Note that Fig. 1 also shows the results of the numerical integration of the original version reported by Thomas et al. (2020) (red triangles), and their results disagree with the results of the present study. A possible reason for this discrepancy might be that their results did not achieve a statistically steady state. For details, see Appendix C.
The original version of the eddyhopping model given by Eqs. (1) and (2) shows the improper scaling for small scales because of the assumption made in the formulation of the model. Originally, Eq. (2) (corresponding to Eq. 8 in Grabowski and Abade, 2017) was formulated under the assumption of large scales (or Da≫1), since this assumption usually holds for typical situations in atmospheric clouds. Thus, it is reasonable that the original version does not reproduce the proper scaling for small scales.
We next consider the second version of the eddyhopping model by Abade et al. (2018), which is written as follows:
Note that, for subsequent use, we write the governing equations in a slightly generalized form by introducing two parameters c_{1} and c_{2}. The second version by Abade et al. (2018) has ${c}_{\mathrm{1}}={c}_{\mathrm{2}}=\mathrm{1}$.
The important change introduced by Abade et al. (2018) into the original version is the term proportional to ${S}^{\prime}/\mathit{\tau}$ in Eq. (19). Physically, this term represents the damping effect on S^{′} due to turbulent mixing (eddy diffusivity). This type of term is commonly included in stochastic models used in cloud turbulence research (Sardina et al., 2015, 2018; Chandrakar et al., 2016; Siewert et al., 2017; Saito et al., 2019a). The timescale of the damping effect due to turbulent mixing is characterized by the integral time τ, whereas that due to condensation or evaporation of cloud droplets is characterized by the phase relaxation time τ_{relax}. The relative importance of these two effects is characterized by the Damköhler number ($Da=\mathit{\tau}/{\mathit{\tau}}_{\mathrm{relax}}$), where the damping effect due to turbulent mixing is dominant for Da≪1 (corresponding to small scales). Below we show that the term ${S}^{\prime}/\mathit{\tau}$ plays an essential role in reproducing the proper scaling.
Note that it is possible to further extend the second version by additionally introducing the Wiener process term representing smallscale fluctuations or mixing into Eq. (19) for the supersaturation fluctuation. Readers are referred to Paoli and Shariff (2009) and Sardina et al. (2015) for the Langevin model including such terms, and also to Chandrakar et al. (2021), who implemented such a Langevin model into the LES Lagrangian cloud model and investigated the impact of entrainment mixing and turbulent fluctuations on droplet size distributions in a cumulus cloud. In the present study, however, we focus on statistical properties of the second version with Eqs. (18) and (19).
Applying the analytical procedure described in Sect. 3 to the second version given by Eqs. (18) and (19), we first obtain
instead of Eq. (8). Next, instead of Eq. (11), we have
Finally, the analytical expression corresponding to Eq. (13) is
Asymptotic forms of ${\mathit{\sigma}}_{{S}^{\prime}}$ in Eq. (22) for the large and small scale limits are, respectively, given as follows:

Large scale limit.
For τ→∞ (or equivalently, Da→∞, L→∞), ${\mathit{\sigma}}_{{S}^{\prime}}$ in Eq. (22) is approximated as$$\begin{array}{ll}{\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}& {\displaystyle}\approx {c}_{\mathrm{2}}{a}_{\mathrm{1}}D{a}^{\mathrm{1}}\mathit{\tau}{\mathit{\sigma}}_{{w}^{\prime}}\\ \text{(23)}& {\displaystyle}& {\displaystyle}={c}_{\mathrm{2}}{a}_{\mathrm{1}}{\mathit{\tau}}_{\mathrm{relax}}{\mathit{\sigma}}_{{w}^{\prime}}.\end{array}$$For the case of a constant dissipation rate of turbulent kinetic energy ε, we have
$$\begin{array}{}\text{(24)}& {\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}.\end{array}$$ 
Small scale limit.
For τ→0 (or equivalently, Da→0, L→0), ${\mathit{\sigma}}_{{S}^{\prime}}$ in Eq. (22) is approximated as$$\begin{array}{}\text{(25)}& {\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}\approx {\mathrm{2}}^{\mathrm{1}/\mathrm{2}}{c}_{\mathrm{1}}{a}_{\mathrm{1}}\mathit{\tau}{\mathit{\sigma}}_{{w}^{\prime}},\end{array}$$which indicates that ${\mathit{\tau}}_{\mathrm{relax}}^{\mathrm{1}/\mathrm{2}}$ in Eq. (16) has been replaced by ${\mathit{\tau}}^{\mathrm{1}/\mathrm{2}}$ by introducing the term ${S}^{\prime}/\mathit{\tau}$ in Eq. (19). For the case of a constant dissipation rate of turbulent kinetic energy ε, we have
$$\begin{array}{}\text{(26)}& {\displaystyle}{\mathit{\sigma}}_{{S}^{\prime}}\sim L.\end{array}$$
Therefore, the second version successfully reproduces asymptotic forms ${\mathit{\sigma}}_{{S}^{\prime}}\sim {a}_{\mathrm{1}}{\mathit{\tau}}_{\mathrm{relax}}{\mathit{\sigma}}_{{w}^{\prime}}$ and $\sim {a}_{\mathrm{1}}\mathit{\tau}{\mathit{\sigma}}_{{w}^{\prime}}$ for the large and small scale limits, respectively, which are both consistent with the result of the scaling argument by Lanotte et al. (2009).
Figure 2 (orange curve) shows the theoretical curve given by Eq. (22) for the second version (${c}_{\mathrm{1}}={c}_{\mathrm{2}}=\mathrm{1}$). The second version reproduces the proper scalings both for large and small scales (${\mathit{\sigma}}_{{S}^{\prime}}\sim {L}^{\mathrm{1}/\mathrm{3}}$ and ∼L^{1}, respectively) and demonstrates better agreement with the reference data (black dots) than the original version for L<10^{0} m. Figure 2 (green diamonds) also shows the results of the numerical integration of the second version, which agree with the theoretical curve (orange curve), as expected. Here, the numerical integration was conducted in the same manner as in the previous section (see Appendix A for details).
Although improved, the second version still slightly over and underestimates the supersaturation fluctuations for $L<\mathrm{3}\times {\mathrm{10}}^{\mathrm{1}}$ and $L>\mathrm{2}\times {\mathrm{10}}^{\mathrm{0}}$ m, respectively, as shown in Fig. 2 (compare the orange curve with black dots). This deviation from the reference data can be further reduced by adjusting two parameters c_{1} and c_{2} in Eqs. (18) and (19). The analytical expression (22) and its asymptotic forms (23) and (25) show how c_{1} and c_{2} work. These types of parameters are not new. For example, a parameter corresponding to c_{1} is commonly used in the Langevin stochastic equation in turbulence research (Sawford, 1991; Marcq and Naert, 1998). Formally, the inverse of c_{1} is referred to as the drift coefficient, and the coefficients for the velocity and scalar equations should be distinguished. However, we treat these coefficients as the same parameter in Eqs. (18) and (19) for simplicity. On the other hand, the importance of a parameter corresponding to c_{2} has been demonstrated in a recent study on turbulence modulation by particles (Saito et al., 2019b).
Here, we do not consider any physical meaning for c_{1} and c_{2} and use them just as tuning parameters. Two parameters c_{1} and c_{2} are determined by comparing the theoretical curve given by Eq. (22) with the reference data taken from DNSs and LESs in Thomas et al. (2020). The best fit is given by c_{1}=0.746 and c_{2}=1.28. Figure 3 (solid curve) shows the theoretical curve given by Eq. (22) with these values of c_{1} and c_{2}, which agrees almost perfectly with the reference data (black dots). Although the improvement from the second version with ${c}_{\mathrm{1}}={c}_{\mathrm{2}}=\mathrm{1}$ is slight, this result shows that the eddyhopping model can be easily tuned to reproduce the reference data almost perfectly.
Finally, we discuss the possibility of simplification of the eddyhopping model. Here, our discussion is based on the second version given by Eqs. (18) and (19), but the same argument also applies to the original version given by Eqs. (1) and (2).
The eddyhopping model consists of two evolution equations for the supersaturation and vertical velocity fluctuations, S^{′} and w^{′} respectively, and these two variables fluctuate randomly according to the Ornstein–Uhlenbeck process. However, if we have S^{′} that fluctuates with a proper amplitude and autocorrelation function, then we do not need the evolution equation for w^{′}, because only S^{′} is used in the growth equation of the droplet size. As described in Sect. 4, we obtained an analytical expression for ${\mathit{\sigma}}_{{S}^{\prime}}$, i.e., the standard deviation of the supersaturation fluctuation in the statistically steady state given by Eq. (22). On the other hand, the autocorrelation function for S^{′} in a statistically steady state can also be obtained analytically. The derivation is described in Appendix D. The result is given in Eq. (D14) and is as follows:
where τ_{1} and τ_{2} are, respectively, defined as
We can also obtain the autocorrelation time for S^{′} by time integration of A(t) (see Eq. D16 in Appendix D), which is given as
The autocorrelation function A(t) in Eq. (28) and the autocorrelation time τ_{0} in Eq. (30), with τ_{1} and τ_{2} defined by (29), have important asymptotic forms in two limits. First, for the large scale limit (Da→∞), the asymptotic forms of A(t) and τ_{0} are given by
respectively. Second, for the small scale limit (Da→0), the asymptotic forms of A(t) and τ_{0} are given by
respectively.
Based on analytical expressions for the fluctuation amplitude and the autocorrelation function for S^{′} (Eqs. 22 and 30, respectively), a simplified version of the eddyhopping model is defined as follows:
where ${\mathit{\sigma}}_{{S}^{\prime}}$ and τ_{0} are given by Eqs. (22) and (30), respectively. Note that the simplified model given by Eq. (33) is a singleequation model, as compared to the twoequation model given by Eqs. (18) and (19) before the simplification. The autocorrelation function for S^{′} in the simplified model given by Eq. (33) is as follows:
which has the following two asymptotic forms. First, for the large scale limit (Da→∞),
which agrees with the corresponding asymptotic form given by Eq. (31) for the second version. Second, for the small scale limit (Da→0),
which disagrees with the corresponding asymptotic form given by Eq. (32) for the second version.
Figures 4a through e compare the autocorrelation function for the simplified model (B(t) in Eq. 35: blue dashed curve) and that for the second version (A(t) in Eq. 28: red solid curve) for five cases ranging from Da≪1 to Da≫1. Here, ${c}_{\mathrm{1}}={c}_{\mathrm{2}}=\mathrm{1}$. Note that the time t is normalized by the autocorrelation time τ_{0} for each case. Although B(t) and A(t) share the same autocorrelation time, B(t) deviates from A(t) for cases with Da of order unity or smaller, as shown in Fig. 4a through c. On the other hand, for Da≫1, B(t) agrees with A(t) very well, as shown in Fig. 4d and e.
The simplified model has the desirable convergence property. The autocorrelation function for the simplified model (B(t) in Eq. 35) converges to that for the second version in the largescale limit (Da→∞), as shown in Eq. (36). As confirmed in Figs. 4d and e, the two autocorrelation functions are almost identical for an integral length L greater than 10 m (or Da≥10). In the implementation of the eddyhopping model to the LES Lagrangian cloud model, the integral length L is supposed to roughly correspond to the grid size, which is often greater than several meters to several tens of meters. Therefore, the assumption of large scales (or Da≫1) usually holds, in which case the statistical properties of the simplified model are expected to be almost unchanged after the simplification.
Figure 5 compares the time evolutions of the supersaturation fluctuations obtained from the numerical integration of the original version (dashed red curve), the second version (dotted blue curve), and the simplified model (solid green curve). Here, the numerical integration was conducted in the same manner as described in Appendix A. All results were obtained by using the same random number series. The results are shown in the time range $\mathrm{10}\mathit{\tau}\le t\le \mathrm{20}\mathit{\tau}$, where all cases are already in a statistically steady state.
For small scales, the simplified model produces qualitatively different trajectories of S^{′} from the second version, as shown in Fig. 5a ($L={\mathrm{10}}^{\mathrm{2}}$ m; compare the solid green and dotted blue curves), even though these two models share the same fluctuation amplitude ${\mathit{\sigma}}_{{S}^{\prime}}$ and the autocorrelation time τ_{0} in the statistically steady state. The difference is smaller for intermediate scales (Fig. 5b, L=10^{0} m). For sufficiently large scales (Fig. 5c, L=10^{2} m), the simplified model and the second version produce almost identical results.
The purpose of the present paper was to obtain various statistical properties of the eddyhopping model, a novel cloud microphysical model, which accounts for the effect of the supersaturation fluctuation at unresolved scales on the growth of cloud droplets and on spectral broadening. Two versions of the model are considered: the original version by Grabowski and Abade (2017) and the second version by Abade et al. (2018). Based on derived statistical properties, we first showed in Sect. 3 that the original version fails to reproduce a proper scaling for smaller Damköhler numbers (corresponding to small scales), resulting in a deviation of the model prediction from the reference data taken from DNSs and LESs, as shown in Fig. 1. In Sect. 4, we showed that the second version successfully reproduces the proper scaling and agrees better with the reference data than the original version for small scales (L<10^{0} m in Fig. 2). We also showed that, by adjusting two parameters c_{1} and c_{2}, the second version can almost perfectly reproduce the reference data. In Sect. 5, we discussed the possibility of simplification of the model. The simplified model consists of a single stochastic equation for the supersaturation fluctuation, as in Eq. (33), with amplitude and time parameters given by the corresponding analytical expressions for the model before the simplification. We showed that, for larger Damköhler numbers (corresponding to large scales), the autocorrelation function of the supersaturation fluctuation for the simplified model converges to that for the model before the simplification. This convergence property is desirable because the assumption of large scales usually holds in the typical parameter range for the model implementation in the LES Lagrangian cloud model.
The results of the numerical integration of the original version (Eqs. 1 and 2) and the second version (Eqs. 18 and 19) are shown in Figs. 1 (blue squares) and 2 (green diamonds), respectively. For these experiments, we used the same setting as that in Sect. 5 in Thomas et al. (2020), except that the integration time was increased from 6τ to 10τ. We set ${a}_{\mathrm{1}}=\mathrm{4.753}\times {\mathrm{10}}^{\mathrm{4}}$ m^{−1}, ε=10 cm^{2} s^{−3}, and τ_{relax}=3.513 s, and the integral time τ as
As described in Appendix B, for the case of a constant dissipation rate of turbulent kinetic energy ε, ${\mathit{\sigma}}_{{w}^{\prime}}$ is given as a function of L. We time integrated the governing equations of the model using 12 values of L: L=0.0128, 0.0256, 0.064, 0.128, 0.256, 0.512, 1.024, 2.56, 6.4, 12.8, 25.6, and 64.0 m. The time step δt is set as $\mathrm{1}/\mathrm{1000}$ of τ, and the integration time is 10τ. The numerical scheme is the forward Euler method. The initial condition is such that ${w}^{\prime}\left(\mathrm{0}\right)={\mathit{\sigma}}_{{w}^{\prime}}\mathit{\psi}$ and ${S}^{\prime}\left(\mathrm{0}\right)=\mathrm{0}$. Each result in Figs. 1 (blue squares) and 2 (green diamonds) is obtained by averaging the results for 1000 ensembles with different seeds of random numbers.
We consider classical homogeneous isotropic turbulence, in which energy is mainly injected into the system at large scales, cascaded to smaller scales by nonlinear interaction, and finally dissipated by the molecular viscosity in the smallest scales. In a statistically steady state, the dissipation rate of turbulent kinetic energy is defined as ε. If ε is fixed and the integral scale L is changed, then the kinetic energy E scales as follows (Thomas et al., 2020):
The black dots in Fig. B1 show the relation between L and E in the reference data taken from DNSs and LESs by Thomas et al. (2020) (Table 2 of their study). In their simulation, the dissipation rate was fixed to ε=10 cm^{2}s^{−3}. The orange curve in Fig. B1 indicates the function $E=\mathit{\alpha}{\mathit{\epsilon}}^{\mathrm{2}/\mathrm{3}}{L}^{\mathrm{2}/\mathrm{3}}$, where α is the fitting parameter. The best fit is given by α=0.475. The rootmeansquare turbulent velocity is calculated as a function of L by ${u}_{\mathrm{rms}}={\mathit{\sigma}}_{{w}^{\prime}}=\sqrt{(\mathrm{2}E/\mathrm{3})}$, and ${\mathit{\sigma}}_{{w}^{\prime}}$ is used as the parameter in the eddyhopping model. Note that Thomas et al. (2020) used the same type of largescale forcing as that used by Kumar et al. (2012), where the integral length L is set to be equal to the box length L_{box}.
We confirm that all of the results of the numerical integration of the eddyhopping model in the present study achieved statistically steady states. For this purpose, we first derive the analytical expression for the time evolutions of the variance and covariance of the variables in the model and then compare these analytical expressions with the results of the numerical integration.
The governing equations given by Eqs. (3) and (2) can be rewritten in generalized forms as
where τ_{1} and τ_{2} are the relaxation times for w^{′} and S^{′}, respectively, and the forcing term ${F}_{{w}^{\prime}}\left(t\right)$ satisfies Eq. (4). Evolution equations for the variance and covariance of the variables are derived as follows:
where ${V}_{{w}^{\prime}}\left(t\right)$, C(t), and ${V}_{{S}^{\prime}}\left(t\right)$ are, respectively, defined as
For the numerical integration of the eddyhopping model by Thomas et al. (2020), τ_{1}=τ and τ_{2}=τ_{relax}. Since the initial conditions for w^{′}(t) and S^{′}(t) are set to ${w}^{\prime}\left(\mathrm{0}\right)={\mathit{\sigma}}_{{w}^{\prime}}\mathit{\psi}$ and ${S}^{\prime}\left(\mathrm{0}\right)=\mathrm{0}$ in Thomas et al. (2020), the corresponding initial conditions for the variance and covariance are given by
Solving Eqs. (C3) through (C5) with the initial conditions given by Eqs. (C9) through (C11), we obtain
where τ_{3} and τ_{4} are, respectively, defined as
Figure C1 compares the analytical expression given by Eq. (C14) (cyan dots) with the results of the numerical integration of the original version given by Eqs. (1) and (2) (black crosses). The setting for the numerical experiment is the same as that used in Fig. 1, except that the integration time is 0.6τ in Fig. C1a and 6τ in Fig. C1b. The results of the numerical integration (black crosses) agree well with the analytical expression (cyan dots), and both approach the theoretical curve for the statistically steady state (orange curve in each panel) as the integration time increases. Figure C1a also indicates that the results of the numerical integration of the eddyhopping model by Thomas et al. (2020) (red triangles) are fairly close to our results at 0.6τ. Thus, it might be possible that the integration time of their numerical experiment was not long enough to achieve a statistically steady state.
We derive the analytical expression for the autocorrelation function of the supersaturation fluctuation S^{′}(t) in the eddyhopping model. As in Appendix C, we start from the generalized form of the eddyhopping model as follows:
where τ_{1} and τ_{2} are the relaxation times for w^{′} and S^{′}, respectively, and the forcing term ${F}_{{w}^{\prime}}\left(t\right)$ satisfies Eq. (4). We consider that the system is in a statistically steady state.
First, multiplying Eq. (D2) by ${\mathrm{e}}^{t/{\mathit{\tau}}_{\mathrm{2}}}$ and applying the product rule of differentiation, we obtain
Integrating Eq. (D3) from t=0 to t, we obtain
(Note that we chose the integration range [0,t] for simplicity of notation. Since we consider a statistically steady state, the following discussion is unchanged if the integration range is $[{t}_{\mathrm{0}},{t}_{\mathrm{0}}+t]$.) Applying a similar procedure to that above to Eq. (D1) with the integration range $t:\mathrm{0}\to \mathit{\xi}$, we obtain
Substituting Eq. (D5) into (D4) and calculating some of the integrations, we obtain
Multiplying Eq. (D8) by S^{′}(0) and taking an ensemble average, we obtain
because of the statistical independence $\left(\u2329{F}_{{w}^{\prime}}\left(\mathit{\zeta}\right){S}^{\prime}\left(\mathrm{0}\right)\u232a=\mathrm{0}\right)$. Next, as in the derivation of Eq. (11), we multiply Eq. (D2) by S^{′} and consider the statistically steady state. We obtain
Substituting Eq. (D10) into (D9), we have
Therefore, the autocorrelation function of the supersaturation fluctuation S^{′}(t) for the eddyhopping model given by Eqs. (D1) and (D2) in the statistically steady state is written as follows:
The autocorrelation time τ_{0} is obtained by timeintegrating A(t) as
For the original version of the eddyhopping model given by Eqs. (1) and (2), we have
For the second version given by Eqs. (18) and (19), we have
Programs for numerical integration of the eddyhopping model are written in Fortran 90 and are available upon request.
The data in this study are available upon request.
The authors declare that they have no conflicts of interest.
Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
We are grateful to Kei Nakajima for his technical support. The present study was supported by MEXT KAKENHI grant nos. 20H00225 and 20H02066, by JSPS KAKENHI grant no. 18K03925, by the Naito Foundation, by the HPCI System Research Project (project ID: hp200072, hp210056), by the NIFS Collaboration Research Program (NIFS20KNSS143), by the Japan High Performance Computing and Networking plus Largescale Data Analyzing and Information Systems (jh200006, jh210014), and by High Performance Computing (HPC 2020, HPC2021) at Nagoya University.
This paper was edited by Markus Petters and reviewed by three 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, https://doi.org/10.1175/JASD180078.1, 2018. a, b, c, d, e, f, g, h, i
Chandrakar, K. K., Cantrell, W., Chang, K., Ciochetto, D., Niedermeier, D., Ovchinnikov, M., Shaw, R. A., and Yang, F.: Aerosol indirect effect from turbulenceinduced broadening of clouddroplet size distributions, P. Natl. Acad. Sci. USA, 113, 14243–14248, https://doi.org/10.1073/pnas.1612686113, 2016. a
Chandrakar, K. K., Grabowski, W. W., Morrison, H., and Bryan, G. H.: Impact of entrainment–mixing and turbulent fluctuations on droplet size distributions in a cumulus cloud: An investigation using Lagrangian microphysics with a sub–grid–scale model, J. Atmos. Sci., https://doi.org/10.1175/JASD200281.1, 2021. a
Clark, T. L. and Hall, W. D.: A numerical experiment on stochastic condensation theory, J. Atmos. Sci., 36, 470–483, https://doi.org/10.1175/15200469(1979)036<0470:ANEOSC>2.0.CO;2, 1979. a
Cooper, W. A.: Effects of variable droplet growth histories on Droplet size distributions. Part I: Theory, J. Atmos. Sci., 46, 1301–1311, https://doi.org/10.1175/15200469(1989)046<1301:EOVDGH>2.0.CO;2, 1989. a
Devenish, B. J., Bartello, P., Brenguier, J.L., Collins, L. R., Grabowski, W. W., IJzermans, R. H. A., Malinowski, S. P., Reeks, M. W., Vassilicos, J. C., Wang, L.P., and Warhaft, Z.: Droplet growth in warm turbulent clouds, Q. J. Roy. Meteor. Soc., 138, 1401–1429, https://doi.org/10.1002/qj.1897, 2012. a
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, https://doi.org/10.1175/JASD170043.1, 2017. a, b, c, d, e, f, g
Grabowski, W. W. and Wang, L.P.: Growth of cloud droplets in a turbulent environment, Annu. Rev. Fluid Mech., 45, 293–324, https://doi.org/10.1146/annurevfluid011212140750, 2013. a
Korolev, A. V. and Mazin, I. P.: Supersaturation of water vapor in clouds, J. Atmos. Sci., 60, 2957–2974, 2003. a, b
Kostinski, A. B.: Simple approximations for condensational growth, Environ. Res. Lett., 4, 015005, https://doi.org/10.1088/17489326/4/1/015005, 2009. a
Kumar, B., Janetzko, F., Schumacher, J., and Shaw, R. A.: Extreme responses of a coupled scalar–particle system during turbulent mixing, New J. Phys., 14, 115020, https://doi.org/10.1088/13672630/14/11/115020, 2012. a
Lanotte, A. S., Seminara, A., and Toschi, F.: Cloud droplet growth by condensation in homogeneous isotropic turbulence, J. Atmos. Sci., 66, 1685–1697, https://doi.org/10.1175/2008JAS2864.1, 2009. a, b
LasherTrapp, S. G., Copper, W. A., and Blyth, A. M.: Broadening of droplet size distributions from entrainment and mixing in a cumulus cloud, Q. J. Roy. Meteor. Soc., 131, 195–220, https://doi.org/10.1256/qj.03.199, 2005. a
Marcq, P. and Naert, A.: A Langevin equation for the energy cascade in fully developed turbulence, Physica D, 124, 368–381, https://doi.org/10.1016/S01672789(98)002371, 1998. a
Paoli, R. and Shariff, K.: Turbulent condensation of droplets: Direct simulation and a stochastic model, J. Atmos. Sci., 66, 723–740, https://doi.org/10.1175/2008JAS2734.1, 2009. a
Politovich, M. K. and Cooper, W. A.: Variability of the supersaturation in cumulus clouds, J. Atmos. Sci., 45, 2064–2086, https://doi.org/10.1175/15200469(1988)045<1651:VOTSIC>2.0.CO;2, 1988. a
Pope, S. B.: Turbulent flows, Cambridge University Press, Cambridge, 2000. a
Saito, I., Gotoh, T., and Watanabe, T.: Broadening of cloud droplet size distributions by condensation in turbulence, J. Meteor. Soc. Japan, 97, 867–891, https://doi.org/10.2151/jmsj.2019049, 2019a. a
Saito, I., Watanabe, T., and Gotoh, T.: A new time scale for turbulence modulation by particles, J. Fluid Mech., 880, R6, https://doi.org/10.1017/jfm.2019.775, 2019b. a
Sardina, G., Picano, F., Brandt, L., and Caballero, R.: Continuous growth of droplet size variance due to condensation in turbulent clouds, Phys. Rev. Lett., 115, 1–5, https://doi.org/10.1103/PhysRevLett.115.184501, 2015. a, b
Sardina, G., Picano, F., Brandt, L., and Caballero, R.: Broadening of cloud droplet size spectra by stochastic condensation: effects of mean updraft velocity and CCN activation, J. Atmos. Sci., 75, 451–467, https://doi.org/10.1175/JASD170241.1, 2018. a
Sawford, B. L.: Reynolds number effects in Lagrangian stochastic models of turbulent dispersion, Phys. Fluids AFluid, 3, 1577–1586, https://doi.org/10.1063/1.857937, 1991. a
Sedunov, Y. S.: Physics of Drop Formation in the Atmosphere, John Wiley and Sons, New York, 1974. a
Shaw, R. A.: Particleturbulence interactions in atmospheric clouds, Annu. Rev. Fluid Mech., 35, 183–227, https://doi.org/10.1146/annurev.fluid.35.101101.161125, 2003. a
Siewert, R., Bec, J., and Krstulvic, G.: Statistical steady state in turbulent droplet condensation, J. Fluid Mech., 810, 254–280, https://doi.org/10.1017/jfm.2016.712, 2017. a
Thomas, L., Grabowski, W. W., and Kumar, B.: Diffusional growth of cloud droplets in homogeneous isotropic turbulence: DNS, scaledup DNS, and stochastic model, Atmos. Chem. Phys., 20, 9087–9100, https://doi.org/10.5194/acp2090872020, 2020. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p
 Abstract
 Introduction
 Governing equations
 Statistical properties of the original version
 Statistical properties of the second version
 Possibility of simplification of the model
 Summary and conclusions
 Appendix A: Numerical integration of the eddyhopping model
 Appendix B: Scalings for the case of a constant dissipation rate of turbulent kinetic energy
 Appendix C: Achievement of a statistically steady state
 Appendix D: Derivation of autocorrelation function
 Code availability
 Data availability
 Author contributions
 Competing interests
 Disclaimer
 Acknowledgements
 Review statement
 References
 Abstract
 Introduction
 Governing equations
 Statistical properties of the original version
 Statistical properties of the second version
 Possibility of simplification of the model
 Summary and conclusions
 Appendix A: Numerical integration of the eddyhopping model
 Appendix B: Scalings for the case of a constant dissipation rate of turbulent kinetic energy
 Appendix C: Achievement of a statistically steady state
 Appendix D: Derivation of autocorrelation function
 Code availability
 Data availability
 Author contributions
 Competing interests
 Disclaimer
 Acknowledgements
 Review statement
 References