Super-Droplet Method for the Numerical Simulation of Clouds and Precipitation: a Particle-Based Microphysics Model Coupled with Non-hydrostatic Model

0.54 MB·English

by  Shin-ichiro Shima

Super-Droplet Method for the Numerical Simulation of Clouds and Precipitation: a
Particle-Based Microphysics Model Coupled with Non-hydrostatic Model
Shin-ichiro Shima,∗ Kanya Kusano, Akio Kawano, Tooru Sugiyama, and Shintaro Kawahara
The Earth Simulator Center, Japan Agency for Marine-Earth Science and Technology, Yokohama, Japan
A novel simulation model of cloud microphysics is developed, which is named Super-Droplet
Method (SDM). SDM enables accurate calculation of cloud microphysics with reasonable cost
in computation. A simple SDM for warm rain, which incorporates sedimentation, condensa-
7 tion/evaporation, stochastic coalescence, is developed. The methodology to couple SDM and a
0 non-hydrostaticmodelisalso developed. It isconfirmed thattheresult of ourMonteCarlo scheme
0 for the coalescence of super-droplets agrees fairly well with the solution of stochastic coalescence
2 equation. A preliminary simulation of a shallow maritime cumulus formation initiated by a warm
bubbleispresentedtodemonstratethepracticalityofSDM.Furtherdiscussionsaredevotedforthe
n
extensionandthecomputationalefficiencyofSDMtoincorporatevariouspropertiesofclouds,such
a
as, several types of ice crystals, several sorts of soluble/insoluble CCNs, their chemical reactions,
J
electrification, and the breakup of droplets. It is suggested that the computational cost of SDM
9
becomes lower than spectral (bin) method when the number of attributes d becomes larger than

some critical value, which may be 2∼4.
]
h
p
- I. INTRODUCTION lightning.
o
SDM can be regarded as a sort of Direct Simula-
a
tion Monte Carlo (DSMC) method, which was initially
.
s Althoughcloudsplayacrucialroleinatmosphericphe-
proposed to simulate the Boltzmann equation for pre-
c
nomena, the numerical modeling of cloud is not well es-
i dicting rarefied gas flows (Bird 1994). Particularly,
s tablished up to the present. The macro-scale processes
SDM has similarity to the Extended version of the No-
y
such as the fluid motion of moist air associated with
h Time-Counter (ENTC) method, which was developed
clouds is called ”cloud dynamics” and the micro-scale
p by Schmidt and Rutland (2000) for the simulation of
[ processes such as condensation/evaporation, stochastic spray flows. Though the basic idea of using a compu-
coalescence, and sedimentation of water droplets, are
tational particle with varying multiplicity is common in
1
called ”cloud microphysics”. These two processes mu-
both SDM and ENTC, there is a difference in the proce-
v
tually affect each other, and thus we can say that
3 dure of the Monte Carlo scheme.
cloud formation and precipitation development are typi-
0 The main purpose of the present paper is to de-
1 calmultiscale-multiphysicsphenomena. Clouddynamics velop a simple SDM for warm rain, which incorporates
1 model to describe the fluid motion in the atmosphere
sedimentation, condensation/evaporation, and stochas-
0 has been well developed (Jacobson 2005). However,it is
tic coalescence, and discuss the theoretical foundations
7 still difficult to perform an accurate simulation of cloud
of SDM. The methodology to couple SDM and the non-
0
microphysics though several simulation methods, such
/ hydrostatic model is also developed. Further discussions
s as bulk parameterization method (Kessler 1969; Ziegler
on the extensions and the computational efficiency of
c
1985; Murakami 1990; Ferrier 1994; Meyers et al. 1997),
i SDM are also carried out, even though briefly.
s spectral (bin) method (Berry 1967; Soong 1974; Bott
The organizationofthe presentpaper is as the follow-
y
1998, 2000), and the exact Monte Carlo method (Gille-
h ing. InSec.IIweintroducetheprinciplemodel,whichisa
spie 1975; Seeßelberg et al. 1996), have been proposed.
p veryprimitivemicrophysics-macrodynmicscoupledcloud
Numericalmethodstoaccuratelysimulatebothprocesses
: model and a good starting point for our discussions. In
v andtheirinteractionsarerequiredtounderstandandpre-
Sec. III we review the traditional methods, i.e., the ex-
i
X dict cloud-related phenomena. actMonteCarlomethod,bulkparameterizationmethod,
r We have developed a novel, particle-based simula- and spectral (bin) method and see how they are derived
a
tion model of cloud microphysics, named Super-Droplet fromtheprinciplemodelandwhataretheproblemsthey
Method (SDM), which enables accurate numerical simu- involve. In Sec. IV, we develop a simple SDM for warm
lationofcloudmicrophysicswithreasonablecostincom- rain, which can be regarded as a coarse-grained model
putation. Though several extensions and validations are of the principle model, and discuss the theoretical foun-
stillnecessary,weexpectthatSDMprovidesusanewap- dationandcharacteristicsofSDM.Numericalsimulation
proach to the cloud-related open problems, such as the scheme is also proposed for each process. Especially, a
cloud and aerosol interactions, the cloud-related radia- novel Monte Carlo scheme for the stochastic coalescence
tive processes,andthe mechanismof thunderstorms and processis developedandvalidated. InSec.Vaverypre-
liminary result of a shallow maritime cumulus formation
initiated by a warm bubble is presented to demonstrate
thepracticalityofSDM.InSec.VIwebrieflydiscussthe
∗Electronicaddress: s [email protected] possibilitytoextendSDMtoincorporatevariousproper-
2
tiesofclouds,suchas,severaltypesoficecrystals,several 2. motion of real-droplets
sortsofsoluble/insolubleCCNs,theirchemicalreactions,
electrification, and the breakup of droplets. The esti- Ingeneral,the motionequationofareal-dropletisde-
mation of the computational cost of SDM suggests that termined by the drag force from the ambient air and the
SDMwillbe computationallymoreefficientcomparedto gravitationalforce:
spectral (bin) method when the number of attributes d
becomes larger than some critical value, which may be m dvi =m g+F (v ,U(x ),R ), dxi =v , (1)
i i D i i i i
2 4. Finally, a summary and concluding remarks are dt dt
∼
presented in Sec. VII.
where m := (4π/3)R3ρ is the mass of the droplet i
i i liq
with ρ = 1.0 g cm−3, g is the gravitational constant,
liq
F is the atmospheric drag, and U(x ) is the wind ve-
D i
locity at the droplet position x (t). Assuming that all
i
II. THE PRINCIPLE MODEL the droplets reach their terminal velocity immediately
enough, v (t) is evaluated as v (t) = U(x ) zˆv (R ),
i i i ∞ i
−
which is the steady solution of (1). We give v (R )
Thissectionisdevotedfortheintroductionoftheprin- ∞ i
according to the semi-empirical formulas developed by
ciple model. The principle model is a very primitive
Beard (1976). Consequently, v (t) is no longer an inde-
microphysics-macrodynmics coupled cloud model. The i
pendent attribute of real-droplet in our model.
characteristic point of this model is the very precise de-
scriptions of cloud microphyscs processes, which is re-
alized by following all the water droplets in the atmo-
3. condensation/evaporation of real-droplets
sphere. Though being not so useful because of the high
costincomputation,thismodelgivesusaclearoverview
We consider the mass change of real-droplets through
of the relationship between each methods to be intro-
the condensation/evaporation process according to
duced later, i.e., the exact Monte Carlo method, bulk
K¨oheler’s theory, which takes into account the solution
parameterization method, spectral (bin) method, and
and curvature effects on the droplet’s equilibrium vapor
super-droplet method.
pressure (Ko¨heler 1936; Rogers and Yau 1989). The
growth equation of the radius R is derived as the fol-
i
lowing.
A. cloud microphysics for the principle model R dRi = (S−1)− Rai + Rb3i , (2)
i
dt F +F
k d
1. definition of real-droplets L Lρliq ρliqRvT
F =( 1) , F = .
k d
R T − KT De (T)
v s
Warmcloudsandrainconsistoflargenumberofwater Here,S istheambientsaturationratio;F representsthe
k
droplets with broadrangeof radius from micrometers to thermodynamic term associated with heat conduction;
several millimeters. Several sorts of aerosol particles are F is the term associated with vapor diffusion; the term
d
also floating in the atmosphere,which can be a seed of a a/R represents the curvature effect which expresses the
i
cloud droplet. These particles are called Cloud Conden- increase in saturation ratio over a droplet as compared
sationNuclei(CCN).Weregardeachoftheseparticlesin- to a plane surface; the term b/R3 shows the reduction
i
cluding both water droplets and CCNs as a real-droplet. in vapor pressure due to the presence of a dissolvedsub-
Let Nr(t) be the number of real-droplets floating in stance and b depends on the mass of solute Mi dissolved
the atmosphere at time t. Each real-droplet is located in the droplet. Numerically, a 3.3 10−5/T (cm) and
≃ ×
at position xi(t), i = 1,2,...,Nr(t), with velocity vi(t) b ≃ 4.3iMi/ms (cm3), where T (K) is the temperature,
and radius Ri(t). Soluble CCN is dissolved in each real- i≃2isthedegreeofionicdissociation,ms isthemolecu-
droplets. Here, we consider only one soluble substance larweightofthesolute. Rv istheindividualgasconstant
as the CCN and M denotes the mass of solute dissolved forwatervapor,K isthecoefficientofthermalconductiv-
i
in the real-droplet i. We refer to these state variables, ityofair,Disthemoleculardiffusioncoefficient,Listhe
the velocity vi(t), the radius Ri(t), and the mass of so- latent heat of vaporization, and es(T) is the saturation
lute M (t), as the attributes of real-droplets denoted by vapor pressure.
i
a (t). Notethatwedonotregardthepositionx (t)asan
i i
attribute forconvenience. Assumingthatthe dropletve-
locity immediately reaches their terminal velocity, v (t) 4. stochastic coalescence of real-droplets
i
isdeterminedbytheirradiusR (t). Asaresult,thenum-
i
ber of independent attributes of our real-droplet is two, Two real-droplets may collide and coalesce into one
the radiusR (t)andthe massofsoluteM (t). Inthe fol- big real-droplet and this process is responsible for pre-
i i
lowing, we give the time evolution law of real-droplets. cipitation development. It is worth noting here that the
3
coalescenceprocessisofourparticularconcerninthispa- Here, D/Dt := ∂/∂t+U is the material derivative;
·∇
per because SDM is expected to overcome the difficulty ρ = ρ +ρ is the density of moist air, which is repre-
d v
in the numerical simulation of coalescence process. sented by the sum of dry air density ρ and vapor den-
d
Our primary postulate, which is commonly used in sity ρv; qv = ρv/ρ is the mixing ratio of vapor; U is the
many previous studies, is that the coalescence process wind velocity; T is the temperature; θ is the potential
can be described in a probabilistic way. Consider a re- temperature; Π = (P/P0)(Rd/cp) is the Exner function
gion with volume ∆V which is small enough compared with reference pressure P0; ρw is the density of liquid
to the space-scale of the atmospheric fluid. Then, we water; Sv is the source term of water vapor associated
assume that the real-droplets inside this region are well- with condensation/evaporationprocess;g is the gravita-
mixedandcoalescewitheachotherwithsomeprobability tionalconstant;λandκ arethe transportcoefficient;Rd
during a sufficiently short time interval (t,t+∆t), i.e., is the gas constant for dry air; cp is the specific heat of
there exists such C(R ,R ) that satisfies dryairatconstantpressure;Listhelatentheatofvapor-
j k
ization. Equations (5)-(9) represent equation of motion,
∆t equation of state, thermodynamic equation, mass conti-
P =C(R ,R )v v
jk j k | j − k|∆V nuity, and water continuity, respectively. Note here that
∆t it is assumed that qv 1 to derive these equations.
=K(R ,R ) ≪
j k ∆V (3) Therearethreecoupling termsfromthe microphysics:
=probability that real-droplet j and k ρw represents the momentum coupling; LSv/cpΠ is
the release of latent heat through the condensa-
inside the small region ∆V will coalesce
tion/evaporationprocess; S is the source term of vapor
v
in the short time interval (t,t+∆t).
through the condensation/evaporation process. These
sourcetermsareevaluatedbythemicrophysicsvariables,
Here,C(R ,R )canberegardedastheeffectivecollision
j k
cross-section of the real-droplets j and k and may be
Nr
evaluated as ρ (x,t):= m (t)δ3(x x (t)), (10)
w i i
−
C(R ,R )=E(R ,R )π(R +R )2, (4) Xi=1
j k j k j k 1 ∂ρ
S (x,t):= − w(x,t), (11)
v ρ(x,t) ∂t
where E(R ,R ) is the collection efficiency, which takes
j k
into account the effect that a smaller droplet would be
sweptaside bythe streamflow aroundthe largerdroplet where mi is the mass of real-droplet i.
or bounce on the surface (Davis 1972; Hall 1980; Jonas The principle model has defined completely now and
1972; Rogers and Yau 1989) . K(R ,R ) is defined by wecanseethatthedescriptionofthecloudmicrophysics
j k
K(R ,R ):=C(R ,R )v v anditiscalledthecoa- is very precise. In the next section, we review the tra-
j k j k j k
| − |
lescence kernel, which will later appear in the stochastic ditional methods from the point of view how they are
coalescence equation. related to the principle model.
Thecloudmicrophysicalprocessoftheprinciplemodel
has been defined. We will turn to determine the cloud
dynamics process of the principle model.
III. TRADITIONAL METHODS
B. cloud dynamics of the principle model In this section, we review the traditionalmethods and
see how they are derived from the principle model and
what are the problems involved.
We may apply non-hydrostatic model to describe the
clouddynamicsprocessoftheprinciplemodel. Thereare
several versions of non-hydrostatic model, but the basic
equations used in this paper is as the following.
A. The exact Monte Carlo method
DU
ρ = P (ρ+ρ )g+λ 2U, (5)
Dt −∇ − w ∇ The exact Monte Carlo method was developed by
P =ρR T, (6) Gillespie (1975) and improved by Seeßelberg et al.
d
(1996). This method simulates the principle modelvery
Dθ L
= Sv+κ 2θ, (7) directly. Their procedure repeatedly draws a random
Dt −c Π ∇
p waiting time for which the next one pair of real-droplets
Dρ
= ρ U, (8) will coalesce. Consequently, it is very rigorous, but not
Dt − ∇· so useful to simulate a cloud formation process in a suf-
Dqv =S +κ 2q . (9) ficiently large region because of the very high cost in
v v
Dt ∇ computation.
4
B. bulk parameterization method (12) reduces to
∂n(X,t) 1 X
Bulkparameterizationmethodismostcommonlyused = dX′n(X′)n(X′′)K(X′,X′′)
∂t 2
todealwithcloudmicrophysics. Thismethodempirically Z0 (13)
∞
evaluatesthesourcetermsfrommicrophysics,suchasρw n(X) dX′n(X′)K(X,X′).
and S , as a function of the field variables in the cloud −
v Z0
dynamics model (Kessler 1969; Ziegler 1985; Murakami
This is the simplest case of spectral (bin) method and
1990; Ferrier 1994; Meyers et al. 1997). Then, the gov-
(13) is the equation which is sometimes called the
erningequationsbecomeaclosedforminclouddynamics
stochasticcoalescenceequation(SCE),butinthispaper,
variablesandweneednotcopedirectlywithcloudmicro-
SCE is used as a general term refers to the governing
physics. Consequently,thenumberofdegreesoffreedom
equation of spectral (bin) method.
isreducedsomuchthatthecomputationaldemandisex-
These equations are derived from the probability
tremely relaxed. However, the accuracy seems to be not
law (3) under the assumption that the number density
always enough. This is very likely because bulk param-
has no correlation (Gillespie 1972). This assumption
eterization is justified under the assumption that there
seems to be correct (Seeßelberg et al. 1996), thus the
exists a closed-form equation in cloud dynamics, though
solutionofSCEcanbeveryaccurate. However,thecom-
this is not guaranteed.
putationalcosttosolveSCEisveryexpensive,especially
if we consider many types of microphysics processes and
the number of attributes d becomes large. We will come
C. spectral (bin) method back to this remark in Sec. VI, but in a word, this is
owing to the fact that the governingequation contains a
Another approach to cloud microphysics is spectral d-multiple integral.
(bin)method(Berry1967;Soong1974;Bott1998,2000).
Spectral (bin) method is characterized by solving the
IV. SUPER-DROPLET METHOD
time evolution equation of the number density distribu-
tion of real-droplets.
LetX =(4π/3)R3 bethevolumeofthereal-dropleti, In this section, we develop a simple SDM for warm
i i
and n(X,M,x,t) be the number density distribution of rain, which can be regarded as a coarse-grained model
real-droplets, i.e., n(X,M,x,t)∆X∆M∆V denotes the of the principle model, and discuss the theoretical foun-
dationandcharacteristicsofSDM.Numericalsimulation
number of real-droplets with their volume in the range
of (X,X+∆X), mass of solute in the range of (M,M+ scheme is also proposed for each process. Especially, a
∆M), inside the space region ∆V around x, at time t. novel Monte Carlo scheme for the stochastic coalescence
Then,thetimeevolutionequationforn(X,M,x,t)which process is developed and validated.
correspondsto the principle model canbe derivedas the
following.
A. cloud microphysics for Super-Droplet Method
∂n ∂
+ x vt(X,x)n + fg(X,M,x)n 1. definition of super-droplets
∂t ∇ ·{ } ∂X { }
1 M X
= dM′ dX′n(X′,M′)n(X′′,M′′)K(X′,X′′) First of all, let us define what a super-droplet is.
4
Z0 Z0 Eachsuper-dropletrepresentsamultiple number ofreal-
∞ ∞
n(X,M) dM′ dX′n(X′,M′)K(X,X′), droplets with the same attributes and position, and the
− multiplicity is denoted by the positive integer ξ (t)
Z0 Z0 (12) 1,2,... , which value can be different in each isuper∈-
{ }
droplet and time-dependent due to the definition of co-
where v (X,x) := U(x) zˆv (R) is the terminal ve- alescence introduced later. That is to say, each super-
t ∞
−
locity of a droplet, f (X,M,x) := 4πR (r.h.s of (2)) droplet has their own position x (t) and own attributes
g i
×
satisfies dX/dt=f , K is the coalescence kernel defined a (t), which characterize the ξ (t) number of same real-
g i i
in (3), X′′ := X X′, and M′′ := M M′. The l.h.s droplets represented by the super-droplet i. Here, the
− −
of (12) corresponds to the advection of real-droplets in attribute a (t) consists of the radius R (t) and the so-
i i
real-spaceandattribute-space,andthe r.h.scorresponds lute mass M (t) for this simplest case. Note here that
i
to the expected dynamics of the stochastic coalescence no two real-droplets have exactly the same position and
process. attributes. In this sense, super-droplet is a kind of
On the ideal case when there is only the stochastic coarse-grained view of real-droplets both in real-space
coalescence process and neither sedimentation nor con- and attribute-space. Let Ns(t) be the number of super-
densation/evaporation process exists, and further, only dropletsfloatingintheatmosphereattimet. Then,these
the droplet volume X is regarded as the attribute, then super-droplets represent N (t)= Ns(t)ξ (t) number of
r i=1 i
P
5
real-droplets in total. In the following, we give the time
evolution law of the super-droplets.
2. motions of super-droplets
Exceptthecoalescenceprocess,eachsuper-dropletbe-
haves just like a real-droplet. Thus, a super-droplet
moves according to the motion equation (5) and v (t)
i
is evaluated by the terminal velocity.
It is easy to simulate this process numerically. Define
the time step ∆t and evaluate the terminal velocity
m
v (t)andthe positionx (t+∆t )=x (t)+∆t v (t) at
i i m i m i
each time step.
FIG. 1: Schematic view of the coalescence of super-droplets.
Two super-droplets with multiplicity 2 and 3 undergo coa-
lescence (upper left). This represents the coalescence of 2
real-dropletpairs(lowerleftandright). Asaresultthesuper-
3. condensation/evaporation of super-droplets
dropletwithmultiplicity2becomeslargerandthemultiplicity
of the othersuper-droplet decreases 3→1 (upperright).
Condensation/evaporationprocessofsuper-dropletsis
also governedby the growth equation (2).
We mayadoptimplicitEulerdiscretizationschemefor generality, and
the numerical simulation of (2), yielding
ξ′ =ξ ξ , ξ′ =ξ , (14)
j j − k k k
Ri2,n+2∆1−tgRi2,n = (S−1)−FkRi+,an+F1d+ R3i,bn+1, RMj′j′==RMj,j, RMk′ =j′ =(R(Mj3+j +Rk3M)1k/)3,, ((1165))
x′ =x , x′ =x , (17)
j j k k
whereR R (t=n∆t ). Itmaybenotpossibletode-
i,n i g
≃
termine analytically the next step value R from the where the dashed valuables represent the updated
i,n+1
previous value R . However, using Newton-Raphson value after the coalesce.
i,n
scheme, R can be evaluated numerically.
i,n+1 2. if ξ =ξ ,
j k
ξ′ =[ξ /2], ξ′ =ξ [ξ /2], (18)
j j k j − j
4. stochastic coalescence of super-droplets R′ =R′ = R3+R3 1/3, (19)
j k j k
M′ =M′ =(M +M ), (20)
Stochastic coalescence process for the super-droplets j k (cid:0) j (cid:1)k
must be formulated carefully so that the super-droplets x′j =xj, x′k =xk, (21)
well-approximate the behavior of real-droplets. The for-
whereGauss’symbol[x]isthegreatestintegerthat
mulation procedure is in two steps: 1) define how a pair
is less than or equal to x.
of super-droplets coalesce; 2) determine the probability
that the super-droplet coalescence occurs. Note here that in the case ξ = ξ , no real-droplets will
j k
Consider the case when a pair of super-droplets (j,k) beleftafterthecoalescencebecauseallξ numberofreal-
j
will coalesce. These super-droplets represent ξj number droplets contribute to the coalescence. Thus, we divide
of real-droplets (aj,xj) and ξk number of real-droplets the resulting ξj number of real-droplets into two super-
(ak,xk). Let us define that exactly min(ξj,ξk) pairs droplets.
of real-droplets will contribute to the coalescence of the Thisdefinitionofsuper-dropletcoalescencepossesshas
super-droplet pair (j,k). Figure. 1 is a schematic view a favorable property that the number of super-droplets
of the coalescence of super-droplets. In this example, is unchanged in most cases though the number of real-
ξj =3andξk =2. We canseethatmin(ξj,ξk)=2pairs dropletsalwaysdecreases. The numberofsuper-droplets
ofreal-dropletsundergocoalescence,whichresultsinthe isdecreasedthroughthecoalescenceonlywhenξ =ξ =
j k
decreaseofmultiplicityξj :3 1andtheincreaseofthe 1,i.e.,bothsuper-dropletsarereal-droplets. Thisresults
→
size of the super-droplet k. in ξ′ =0 and ξ′ =1, and we remove the super-droplet j
j k
Belowisthecompletedefinitionofhowapairofsuper- outofthesystem. Becausethe numberofsuper-droplets
droplets (j,k) change their state after the coalesce: correspondstotheaccuracyofSDM,thenumberconser-
vation of super-droplets suggests the flexible response of
1. if ξ = ξ , we can choose ξ > ξ without losing SDMtothedrasticchangeofthenumberofreal-droplets.
j k j k
6
6
Having defined how a pair of super-droplets coalesces, compared to the space-scale of the atmospheric fluid so
weturntodeterminetheprobabilitythatthecoalescence that the coalescence occurs according to the probability
occurs. If we require the consistency of the expectation (22) in each cell. Let us make a list of super-droplets in
value, the probability is determined as the following. a certain cell at time t: I := i ,i ,...i , where n is
{ 1 2 ns} s
the total number of super-droplets in the cell. Then,
Pj(ks) =max(ξj,ξk)Pjk by definition, each pair of super-droplets (j,k) I2,
∈
=probability that super-droplet j and k j = k, coalesces with the probability P(s) within the
(22) 6 jk
inside the small region ∆V will coalesce short time interval (t,t+∆t ). Note that P(s) 1 be-
c jk ≪
in the short time interval (t,t+∆t ). cause ∆t is small. We have to examine all the possible
c c
pairs (j,k) I2, j = k, to know the random realization
Indeed, we can confirm that the expectation value is ∈ 6
at the time t+∆t . However, this yields computational
exactly the same to that of the principle model (real- c
cost O(n2) and this is not efficient. Instead, we propose
droplet view). The super-droplet j represents ξ num- s
j anovelMonteCarloschemewhichleadsustoO(n )cost
ber of real-droplets with attribute a and super-droplet s
j in computation.
k represents ξ number of real-droplets with attribute
k Let L := (j ,k ),(j ,k ),...,(j ,k ) be a
ak. From the real-droplet view, there are ξjξk number { 1 1 2 2 [ns/2] [ns/2] }
listofrandomlygenerated,non-overlapping,[n /2]num-
of real-droplet pairs which have the possibility to coa- s
ber of pairs. Non-overlapping means that no super-
lesce with the probability P , i.e., the number of real-
jk
droplet belongs to more than one pair. This list can
droplet pairs which will coalescence follows the binomial
be made by generating a random permutation of I and
distribution with ξ ξ trials and success probability P .
j k jk
makepairsfromthefront,whichcostO(n )incomputa-
Thus, from the real-droplet view, the expectation value s
tion (see the appendix). By examining only these [n /2]
of the number of coalesced pairs is ξ ξ P . From the s
j k jk
pairs instead of the whole possible n (n 1)/2 pairs,
super-dropletview,acoalescenceofsuper-dropletsj and s s −
the cost will be O(n ). In compensation for this simpli-
k represents the coalescence of min(ξ ,ξ ) pairs of real- s
j k
fication, all the coalescence probability is divided by the
droplets with attribute a and a . Thus, the coales-
j k
decreasingratioofpairnumber,andthecorrected,scaled
cence of super-droplets j and k is expected to represent
up probability for the αth pair is obtained as
min(ξ ,ξ )P(s) = min(ξ ,ξ )max(ξ ,ξ )P = ξ ξ P
j k jk j k j k jk j k jk
numberofcoalescenceofreal-dropletpairs,whichisequal n (n 1) n
to the value of the real-droplet view. It is worth notic- pα :=Pj(αs)kα s s2− 2s . (23)
inghere thatthe varianceofsuper-dropletview becomes (cid:30)h i
larger than that of real-droplet view. From real-droplet Thismaybejustifiedifthefollowingrelationholdsgood,
view, the variance of the number of coalesced pairs is
ξtijoξnkPljikm(i1t)−. PFjrko)m≃suξpjeξrk-Pdjrkop=let: vVirew(P, otihsesonvardiiasntrciebuis- [ns/2]min(ξ ,ξ )p 1 min(ξ ,ξ )P(s)
jα kα α ≃ 2 j k jk
{min(ξj,ξk)}2Pj(ks) −{ξjξkPjk}2 ≃ min(ξj,ξk)Vr, which αX=1 Xj,k,jX6=k
is min(ξ ,ξ ) times larger than that of the real-droplet 1
j k = ξ ξ P ,
view. 2 j k jk
j,k,j6=k
We have defined the cloud microphysics of super- XX
droplets. If all the multiplicity ξ = 1 then our model is
i which representthe consistency of the expectation num-
equivalent to the principle model. Thus, we can change
ber of coalesced real-droplet pairs in this cell.
the accuracy of our model by changing the initial dis-
Based on the above idea, below is the complete pro-
tribution of multiplicity ξ (t = 0) , or equivalently, the
i cedure of our Monte Carlo scheme to obtain one of the
{ }
initial number of super-droplets N (t = 0). If the con-
s stochastic realizations in a certain cell at time t+∆t
c
vergenceto the principle model is rapid, we cansaythat
from the sate of super-droplets at time t.
the use of SDM is beneficial.
The definition of the coalescence process of super- 1. Makethesuper-dropletlistattimetinthiscertain
dropletsisrathertheoreticalandwehavetodonumerical cell: I = i ,i ,...i .
{ 1 2 ns}
simulation on this model somehow. In the next section,
we develop a Monte Carlo scheme for the stochastic co- 2. Make the [ns/2] number of candidate-pairs
alescence, which can simulate the coalescence process of L = {(j1,k1),(j2,k2),...,(j[ns/2],k[ns/2])} from
super-droplets efficiently. the random permutation of I.
3. For eachpair ofsuper-droplets(j ,k ) L,gener-
α α
∈
ate a uniform random number φ (0,1). Then,
5. Monte Carlo Scheme for super-droplet coalescence α ∈
evaluate
Let us introduce a grid which covers the real-space.
[p ]+1 if φ <p [p ]
We can choose any kind of grids, which may be not uni- γα := α α α− α
form, but the volume of each cell must be small enough ([pα] if φα ≥pα−[pα]
7
4. Ifγα =0,theαthpair(jα,kα)isupdatedtot+∆tc 6. validation of our Monte Carlo scheme
without changing their state.
In this section, we validate our Monte Carlo scheme
5. If γ =0, choose ξ ξ without losing general-
α 6 jα ≥ kα for the stochastic coalescenceprocess by confirming that
ity and evaluate γ˜ :=min(γ ,[ξ /ξ ]).
α α jα kα the numerical results agree well with the solution of the
(a) if ξjα −γ˜αξkα >0, SCE (13).
TheSCE(13)isthetimeevolutionequationofn(X,t)
ξ′ =ξ γ˜ ξ , ξ′ =ξ ,
jα jα − α kα kα kα in the ideal system in which only the stochastic coales-
R′ =R , R′ =(γ˜ R3 +R3 )1/3, cenceprocessexists. Correspondingly,thesuper-droplets
jα jα kα α jα kα
M′ =M , M′ =(γ˜ M +M ), in this section undergo only the coalescence process ac-
jα jα jα α jα kα cording to our Monte Carlo scheme. The super-droplets
x′jα =xjα, x′kα =xkα, have velocities, but do not move out of the coalescence
cell,nordotheyundergocondensation/evaporation. The
(b) If ξ γ˜ ξ = 0, i.e., γ˜ = ξ /ξ =
[ξ /jαξ −] αγkα, α jα kα radius Ri is the only attribute.
jα kα ≤ α To compare the simulation result, we have to recon-
ξ′ =[ξ /2], ξ′ =ξ [ξ /2], struct the number density distribution n(X,t) from the
jα kα kα kα − kα
data of super-droplets (ξ ,R ) . There are severalways
Rj′α =Rk′α = γ˜αRj3α +Rk3α 1/3, to do this, but we adop{t kierniel}density estimate method
M′ =M′ =(γ˜ M +M ), with Gaussian kernel (Terrell and Scott 1992). Some
jα kα (cid:0) α jα (cid:1)kα brief discussions on the application of kernel density es-
x′ =x , x′ =x .
jα jα kα kα timate method to SDM is also developed in Sec. IV B.
Ifξ′ =0,thesuper-dropletj isremovedout It is convenient to plot the results by the mass density
of tjhαe system. α function over lnR, which is defined by g(lnR)dlnR :=
ρ Xn(X)dX. The corresponding estimator function
liq
Our Monte Carlo scheme numerically solves the g˜(lnR) becomes,
stochastic coalescence process of super-droplets defined
in Sec. IV A 4. We used a sort of random sampling of 1 Ns
candidate-pairstoachievethecomputationalcostO(ns). g˜(lnR):= ξimiWσ(lnR lnRi),
∆V −
Further extension are implemented in our Monte Carlo
i=1
X
schemeto dealwith multiple coalescence. The integerγα W (Y):= 1 exp Y2/2σ2 .
representshowmanytimesthepair(j ,k )willcoalesce, σ
α α √2πσ −
andγ˜ istherestrictedvalueofγ byitsmaximumvalue
α α (cid:0) (cid:1)
[ξjα/ξkα]. Primarily, pα must be smaller than 1 because Here, ∆V is the volume of the coalescence cell and
pα represents a probability, thus γα should be either 0 σ = σ N−1/5 with some constant σ . Theoretically,
0 s 0
or 1. However, if we take ∆t too much larger, p may
c α the most efficient choice of σ is estimated as σ =
0 0
become larger than 1. Though the situation γα > 1 is (2√π∆V2 g′′2dlnR/M N )−1/5,whereM istheto-
tot r tot
not consistent with the fundamental premise of our the-
talmassofliquidwater. Becauseg istheestimatedfunc-
ory, we formulate our Monte Carlo scheme to cope with R
tion itself, we have to estimate also σ from the data of
0
thissituationsothatSDMrobustlywell-approximatethe
super-dropletstofindthe maximumlikelihoodestimator
principle model irrespective of the choice of ∆t .
c function g˜, but in this paper we just choose σ empiri-
0
Alternatively be said that the most rigorous selection
cally.
criteriaofthesimulationtimestepforcoalescenceprocess
We have to give an explicit form of the coalescence
∆t is that p 1 holds good for all α. Let us estimate
c α ≪ kernel K(Rj,Rk) in (3) to determine the dynamics. For
∆t veryroughly. Thetypicalnumberdensityandradius
c Golovin’s kernel K(R ,R ) = b(X +X ), we have the
ofsmallclouddropletsmaybe 109 m−3 andR=10µm. j k j k
analyticsolutionofSCE(13)(Golovin1963)thoughthis
The typical terminal velocity of a rain droplet may be
choiceofcoalescencekernelis notrealistic. The resultof
1 m s−1. Then, we may roughly evaluate p as follows.
α our comparison simulation for Golovin’s kernel is shown
(s) in Fig. 2a. We set b = 1.5 103 s−1. We need only
pα ≃nsPjαkα one big coalescence cell, the ×volume of which is chosen
= nsmax(ξjα,ξkα)Eπ(R +R )2 v v ∆t to be ∆V = 106 m3. The time step is fixed as ∆tc =
∆V jα kα | jα − kα| c 1.0 s. The initial number density of real-droplets is set
109m−3 1 π(2 10−5m)2 1 m s−1 ∆tc to be n0 = 223 m−3. The initial size distribution of the
≃ · · × · ·
real-droplets follows an exponential distribution of the
=0.4π∆t <1.
c
dropletvolumeX whichisdeterminedbytheprobability
i
Thus, it is estimated that ∆t < 1/0.4π s 0.8 s. It is density p(X ) = (1/X )exp( X /X ), X = (4π/3)R3,
worth noticing here, that ourcestimation is≃independent R = 30.531i 10−6 m0. This−settiing0resul0ts in the tota0l
0
ofthenumberofsuper-dropletsn andthevolumeofthe amountofliq×uid water1.0 g m−3. The initial number of
s
coalescence cell ∆V. super-dropletN ischangedasthesimulationparameter.
s
8
The initial multiplicity ξ is determined by the equality
i
ξ =n ∆V/N . σ is fixed as σ =0.62.
i 0 s 0 0
Figures 2b and 2c show the results of the case of so-
calledhydrodynamickernel,whichisamuchrealisticco-
alescence kernel. Hydrodynamic kernel is defined by (3)
and (4). The same coalescence efficiency E(R ,R ) is
j k
adopted as described in Seeßelberg et al. (1996) and
Bott (1998). For small droplets the dataset of Davis
(1972) and Jonas (1972) is used, and for large droplets
the dataset of Hall (1980) is used. Because the analytic
solutionisnotavailableforhydrodynamickernel,wegen-
eratedthereferencesolutionnumericallyusingExponen-
tialFluxMethod(EFM)developedbyBott(1998,2000).
AlogarithmicallyequidistantradiusgridisusedinEFM.
We adopted 1000bins rangingfrom 0.62 µm to 6.34 cm.
The time stepwassetto 0.1s. Theseparametersettings
provide a sufficiently accurate numerical solution of the
SCE(13). InFig.2b,simulationsettingsaresameascase
(a) except the coalescence kernel. In Fig. 2c we changed
the initial size distribution by reducing R to one-third,
0
i.e.,R =10.177 10−6m,andincreasedtheinitialnum-
0
ber density of re×al-droplets as n = 33 223 m−3. Con-
0
sequently, the total amount of liquid wa·ter 1.0 g m−3 is
unchanged from case (b). ∆t = 0.1 s and σ = 1.5 for
c 0
case (c).
In case (a) and (b), we cansee that the result of SDM
with N = 213 agrees fairly well with the solution of
s
SCE (13), and much improved when N = 217. How-
s
ever, case (c) is not so good as the previous two cases.
Even N = 217 is not good and we need N = 221
s s
for good agreement. This situation typically reveals the
weak point of SDM. The initial size of the droplets are
comparatively small and the coalescence hardly occurs.
We can see that the thick solid line at t = 1200 s is not
so changed from the initial distribution. However, once
a few big droplets are created, then the coalescence pro-
cessisabruptlyaccelerated. Accordingly,foranaccurate
prediction, we have to resolve the right tail of the distri-
bution at the time t= 1200 s. However, existence prob-
ability of super-droplets is very low for such region and
samplingerroroccursifthenumberofsuper-dropletsN
s
is not sufficient. As a result, we need much more super-
dropletsforcase(c)comparedtothatof(a)and(b). For
the practical applications to cloud formation, this fact
wouldnotimposesevererestrictionsonSDMbecausethe
condensation/evaporationprocess dominates the growth
of small droplets in such cases like (c). FIG. 2: Time evolution of the mass density distribution
Anyway, we have confirmed that SDM reproduces the g(lnR,t), which is numerically obtained by the Monte Carlo
solution of the SCE (13) if the number of super-droplets scheme of SDM. (a): The case of Golovin’s kernel. The solid
linerepresentstheanalyticsolutionofSCE(13). (b)and(c):
N issufficientlylarge. Theseresultssupportthevalidity
s
The case of hydrodynamic kernel. The solid line represents
ofourMonteCarloschemeforthestochasticcoalescence.
the approximate solution of SCE (13), which is numerically
obtained byEFM. Ineach case, we can see that theresult of
SDM agrees fairly well with solution of SCE (13) when the
B. cloud dynamics for Super-Droplet Method numberof super-dropletsNs issufficientlylarge. ∆tc =1.0 s
for (a) and (b), and ∆tc =0.1 s for (c).
Forthe clouddynamicsofSDM,wecanalsoapplythe
non-hydrostatic model (5)-(9). Note here that we have
toevaluatethe sourcetermsfrommicrophysics(10)-(11)
9
not by real-droplets but by super-droplets: bubble is inserted at the center of the horizontal axis
according to the equation below.
Ns
ρ (x,t)= ξ m (t)δ3(x x (t)), (24)
w i i i
− θ(x,t=0)=θ (x,t=0)+1.0 K
i=1 b
X
S (x,t)= −1 ∂ρw(x,t). (25) (z 500.0 m) 2 (x 6.4 km) 2
v ρ(x,t) ∂t exp − − ,
× "−(cid:26) 400.0 m (cid:27) −(cid:26) 1200.0 m (cid:27) #
For the numericalsimulationof (5)-(9),we may adapt
4th-order Runge-Kutta scheme for the time derivative, where θ (x,t = 0) is the base state of potential temper-
b
second-order central difference scheme for the spatial ature.
derivatives, with artificial viscosity to stabilize the nu-
We simulatethe clouddynamicsaccordingtothe non-
mericalinstability. The source terms (24) and (25) must
hydrostatic model introduced in Sec. II B with the nu-
be evaluatedon the numericalgridpoint of fluid simula-
merical schemes proposed in Sec. IV B. The grid width
tion x∗. Generally, introducing an appropriately chosen
j is ∆x = 3.125 m and ∆z = 4.0 m with time step
coarsening function w(x), it can be evaluated as
∆t = 0.0125 s. We set the artificial viscosity ν = κ =
1.56 m2 s−1.
ρ (x∗,t)= w(x∗ x)ρ (x,t)d3x
w j j − w Initially,NaClCCNsareuniformlydistributedinspace
Z with number density 1.0 107 m−3. The solute mass
Ns ×
= ξ m (t)w(x∗ x (t)). Mi(t = 0) has an exponential distribution given by the
i i j − i probability density
i=1
X
Herewsatisfies w(x)d3x=1,andwouldbeapiecewise 1 M
linear,localizedfunctionwithlength-scalecomparableto p(Mi)= exp i , M0 =1.0 10−16 g.
M −M ×
the grid width.RIf a super-droplet is simply shared by 0 (cid:18) 0(cid:19)
the 8 adjacent grid points at the corner of the cell, then
Because it is unsaturated at the beginning, the growth
w(x) = 1/8∆V inside the cube with volume 8∆V, and
equation (2) has a stable and steady solution, which is
w(x)=0 outside the cube.
adopted as the initial value of the radius R (t=0). The
i
SDM has been defined completely and the numerical
coalescenceefficiency usedinSec.IV A6 is alsoadopted
simulation schemes have been also proposed. This is the
to this simulation.
main resultof the presentpaper. In the next section, we
For the numerical simulation of the microphysics pro-
show a preliminary result of SDM.
cesses, we also use the numerical scheme proposed in
Sec. IV. Initially, the super-droplets are uniformly dis-
tributed in space with number density 1.92 m−3 with
V. DEMONSTRATION
common multiplicity, i.e., ξ (t = 0) = INT(1.0
i
107/1.92) = 5208333. The grid for our Monte Car×lo
Averypreliminaryresultofashallowmaritimecumu-
scheme of stochastic coalescence is same to that of non-
lus formationinitiated by a warm bubble is presented in
hydrostatic model. The time steps are ∆t = ∆t =
this section to demonstrate the practicality of SDM. m g
1.0 s and ∆t =10.0 s.
Thesimulationdomainis2-dimensional(x-z),12.8km c
Cloud started to form at about 8 min and began to
in horizontaldirection and 5.12 km in vertical direction.
rain at about 20 min. After raining for half an hour,
Initially, a humid but not saturatedatmosphere is strat-
the cloud remained for a while but finally disappeared
ified, which is almost absolutely unstable under 2 km in
at about 90 min. The amount of precipitation is about
altitude. There is a temperature inversion layer above
1mmin total. Figure3 is asnapshotatthe time 1608s.
2 km. The equations below determines the initial base
We plot all the super-droplets with the color and the
sounding of the atmosphere.
alphatransparencywhicharedeterminedbytheirradius
U(x,t=0)=0, Ri and multiplicity ξi. The color map is indicated in
the figure and the alpha value is determined by α =
T(x,z =0,t=0)=305.0 K, i
7.8 10−5(R /10−3 m)2 ξ ,whichisproportionaltoR2ξ
P(x,z =0,t=0)=101325 Pa, × i i i i
Wecanseethatthereisaturbulentlikestructuresinside
∂T 9.5 K km−1 if z <2.0 km, the cloud. Note that there are super-droplets also in the
∂z =(+−3.0 K km−1 if z 2.0 km, black region,but it is not depicted in this figure because
≥ the radius is very small, which represent CCNs.
2
(z+200.0 m) Inthissection,weshowedaverypreliminaryresultofa
q (x,t=0)=0.022exp .
v "−(cid:26) 2000.0 m (cid:27) # shallowmaritimecumulusformationinitiatedbyawarm
bubble, and demonstrated that SDM actually works as
The upper and lower boundary is fixed to their initial a cloud resolving model. Application to more realistic
value and the horizontal boundary is periodic. A warm situations are our future works.
10
VI. DISCUSSIONS Breakup of droplets can be also taken into account.
This can be done if a governing law is given for the
A. extension of Super-Droplet Method breakup process. However, note here that breakup re-
sults in the increase of super-droplets, and this is not
preferablebecausethemorethesuper-dropletisthemore
In Sec. IV, SDM was developed only for warm rain
the computational cost will be. Breakup process may
withonlyonesolublesubstanceasthe CCN.However,it
result in the explosion of the number of super-droplets
is very conceivable that we can extend SDM to incorpo-
N (t) . In other words,SDM is formulatedsuitable
rate various properties of clouds, such as, several types s →∞
for the simulation of such phenomena in which the coa-
of ice crystals, several sorts of soluble/insoluble CCNs,
lescenceprocessis predominantandthe breakupprocess
their chemical reactions,electrification, and the breakup
has a secondary importance. This is exactly the case for
of droplets. Here, we give a rough overview of how to
cloud formation process.
achieve these extensions.
SDM can be regarded as a sort of Direct Simula-
For this purpose, let us extend the type of attributes
tion Monte Carlo (DSMC) method, which was initially
a. To treat several kinds of soluble/insoluble chemicals
proposed to simulate the Boltzmann equation for pre-
for the CCNs, new attributes are necessary to represent
dicting rarefied gas flows (Bird 1994). Particularly,
the mass of each chemicals. To treat several types of
SDM has similarity to the Extended version of the No-
ice crystals, one attribute may represent the state of the
Time-Counter (ENTC) method, which was developed
super-droplets. Additional attributes may represent the
by Schmidt and Rutland (2000) for the simulation of
electric charge and the temperature of the droplet. If
spray flows. Though the basic idea of using a compu-
we do not give a terminal velocity, then the velocity of a
tational particle with varying multiplicity is common in
droplet is also an independent attribute.
both SDM and ENTC, there is a difference in the proce-
Let us denote the attributes of our new super-droplet
dureoftheMonteCarloscheme. SDMisusingrandomly
as a(t) = (a(1),a(2),...,a(d)), where d is the number of
generated, non-overlapping, [n /2] number of candidate
independent attributes. Then, we have to give a time s
pairs, and allows multiple coalescence for each pair. On
evolution law for each attributes, which determines the
the other hand, ENTC method chooses the candidate
individual dynamics of super-droplets:
pairsrandomlyfromthesetofallthepossiblepairs. The
da number of candidate-pairs to be sampled varies propor-
i =f(a ,A(x )), i=1,2,...,N (t), (26)
dt i i s tionally to the maximum coalescence probability in the
cell. Both SDM and ENTC method result in the com-
where A(x ) represents the field variable which charac-
i putational cost O(N ). Probably, ENTC method is also
terize the state of ambient atmosphere. Equation (26) s
applicabletocloudsimulationsandSDMisapplicableto
may include the growth equation (2) and the time evo-
spray simulations, but further verifications are still nec-
lution of solute mass dM /dt = 0. For example, the
i essary to compare their capability.
chemicalreactionsandthe electrificationprocessmaybe
also incorporated in (26).
The coalescence process must be specified for our new
B. computational cost and accuracy
super-droplets. Atfirst, we haveto determine whatkind
of real-dropleta′ will be created after the coalescence of
Ourestimationsuggestthatthe computationalcostof
real-droplets a and a . This can be either determinis-
j k SDM becomes much lower than spectral (bin) method
ticorstochastic. Then,thecoalescenceofsuper-droplets
when the number of attributes d becomes large. We
will be determined similarly to (14)-(21). It is not the
wouldliketoshowthebasisofthisassertionverybriefly.
casein(17)and(21),butwemayalsochangetheposition
Let us estimate the computational cost and accuracy
of super-droplets after the coalescence if necessary. Sim-
ofspectral(bin) method. We have extendedSDM in the
ilar to (3), the coalescence probability for real-droplets
previous section, and the corresponding,general form of
must be given in the form
SCE can be derived as follows.
∆t
P =C(a ,a ) c v v ∂n(a,x,t) ∂n
jk j k ∆V | j − k| + x vn + a fn = , (27)
∂t ∇ ·{ } ∇ ·{ } ∂t
=probability that real-droplet j and k (cid:18) (cid:19)c
inside the small region ∆V will coalesce wheren(a,x,t)isthenumberdensitydistributionandf
is defined in (26). The term(∂n/∂t) correspondsto the
in the short time interval (t,t+∆t ), c
c
stochasticcoalescenceprocessandthiswillbed-multiple
which in general is a function of the attributes a and integralingeneralbecausewehavetosurveyalloverthe
j
a . Then the coalescence probability for super-droplets d-dimensional attribute-space to evaluate the variation
k
will be determined by (22). of droplet number density at a. Obviously, this term is
Thecoalescenceprocessisashorttimeandshortrange the most expensive for computation, so we neglect the
interaction between the droplets, but we may also con- advection terms in the l.h.s of (27) and focus our atten-
sider other types of interactions if any. tion on the simplified equation ∂n/∂t=(∂n/∂t) , which
c