32 template<
class CloudType>
47 if (!this->defaultCoeffs(
true))
49 this->coeffDict().readEntry(
"B0", b0_);
50 this->coeffDict().readEntry(
"B1", b1_);
51 this->coeffDict().readEntry(
"Ctau", cTau_);
52 this->coeffDict().readEntry(
"CRT", cRT_);
53 this->coeffDict().readEntry(
"msLimit", msLimit_);
54 this->coeffDict().readEntry(
"WeberLimit", weberLimit_);
59 template<
class CloudType>
67 msLimit_(bum.msLimit_),
68 weberLimit_(bum.weberLimit_)
74 template<
class CloudType>
81 template<
class CloudType>
107 bool addParcel =
false;
109 const scalar averageParcelMass = this->owner().averageParcelMass();
113 scalar d03 =
pow3(d0);
116 scalar mass = nParticle*d3*rhopi6;
117 scalar mass0 = nParticle*d03*rhopi6;
119 scalar weGas = 0.5*rhoc*
sqr(Urmag)*d/
sigma;
123 scalar reLiquid =
rho*Urmag*r/
mu;
124 scalar ohnesorge =
sqrt(weLiquid)/(reLiquid + VSMALL);
125 scalar taylor = ohnesorge*
sqrt(weGas);
129 scalar gt = (
g + acceleration) & trajectory;
133 (0.34 + 0.38*
pow(weGas, 1.5))
134 /((1.0 + ohnesorge)*(1.0 + 1.4*
pow(taylor, 0.6)))
141 *(1.0 + 0.45*
sqrt(ohnesorge))
142 *(1.0 + 0.4*
pow(taylor, 0.7))
143 /
pow(1.0 + 0.865*
pow(weGas, 1.67), 0.6);
146 scalar tauKH = 3.726*b1_*r/(omegaKH*lambdaKH);
149 scalar dc = 2.0*b0_*lambdaKH;
152 scalar helpVariable =
mag(gt*(
rho - rhoc));
153 scalar omegaRT =
sqrt
155 2.0*
pow(helpVariable, 1.5)
160 scalar KRT =
sqrt(helpVariable/(3.0*
sigma + VSMALL));
167 if ((tc > 0) || (lambdaRT < d) )
173 scalar tauRT = cTau_/(omegaRT + VSMALL);
176 if ((tc > tauRT) && (lambdaRT < d))
180 scalar nDrops = d/lambdaRT;
187 if (weGas > weberLimit_)
189 scalar fraction = dt/tauKH;
192 d = (fraction*dc + d)/(1.0 + fraction);
195 scalar ms0 = mass0*(1.0 -
pow3(d/d0));
198 if (ms/averageParcelMass > msLimit_)
209 scalar de3 = d*d*(dc - d);
211 pow3(be3/(3.0*ae3)) - be3*ce3/(6.0*ae3*ae3) + de3/(2.0*ae3);
212 scalar pe3 = (3.0*ae3*ce3 - be3*be3)/(9.0*ae3*ae3);
213 scalar D3 = qe3*qe3 + pe3*pe3*pe3;
215 if (D3 < 0) br3 =
false;
220 scalar ue3 =
cbrt(-qe3 + D3);
221 scalar ve3 =
cbrt(-qe3 - D3);
222 scalar dParenDrops = ue3 + ve3 - be3/3.;
223 scalar mc = nParticle*(
pow3(d) -
pow3(dParenDrops));
224 scalar nChildDrops = mc/
pow3(dc);
226 if (nChildDrops >= nParticle)
232 massChild = mc*rhopi6;
241 else if (KHindex < 0.5)
248 scalar diameterLargerDrop =
cbrt(1.5*d*d*lengthScale);
249 d = diameterLargerDrop;
255 scalar massDrop =
pow3(d)*rhopi6;
256 nParticle = mass/massDrop;