30 template<
class CloudType>
45 if (!this->defaultCoeffs(
true))
47 this->coeffDict().lookup(
"B0") >> b0_;
48 this->coeffDict().lookup(
"B1") >> b1_;
49 this->coeffDict().lookup(
"Ctau") >> cTau_;
50 this->coeffDict().lookup(
"CRT") >> cRT_;
51 this->coeffDict().lookup(
"msLimit") >> msLimit_;
52 this->coeffDict().lookup(
"WeberLimit") >> weberLimit_;
57 template<
class CloudType>
65 msLimit_(bum.msLimit_),
66 weberLimit_(bum.weberLimit_)
72 template<
class CloudType>
79 template<
class CloudType>
105 bool addParcel =
false;
107 const scalar averageParcelMass = this->
owner().averageParcelMass();
111 scalar d03 =
pow3(d0);
114 scalar mass = nParticle*d3*rhopi6;
115 scalar mass0 = nParticle*d03*rhopi6;
117 scalar weGas = 0.5*rhoc*
sqr(Urmag)*d/
sigma;
118 scalar weLiquid = 0.5*rho*
sqr(Urmag)*d/
sigma;
121 scalar reLiquid = rho*Urmag*r/
mu;
122 scalar ohnesorge =
sqrt(weLiquid)/(reLiquid + vSmall);
123 scalar taylor = ohnesorge*
sqrt(weGas);
125 vector acceleration = Urel/tMom;
127 scalar gt = (g + acceleration) & trajectory;
131 (0.34 + 0.38*
pow(weGas, 1.5))
132 /((1.0 + ohnesorge)*(1.0 + 1.4*
pow(taylor, 0.6)))
139 *(1.0 + 0.45*
sqrt(ohnesorge))
140 *(1.0 + 0.4*
pow(taylor, 0.7))
141 /
pow(1.0 + 0.865*
pow(weGas, 1.67), 0.6);
144 scalar tauKH = 3.726*b1_*r/(omegaKH*lambdaKH);
147 scalar dc = 2.0*b0_*lambdaKH;
150 scalar helpVariable =
mag(gt*(rho - rhoc));
151 scalar omegaRT =
sqrt 153 2.0*
pow(helpVariable, 1.5)
154 /(3.0*
sqrt(3.0*sigma)*(rhoc + rho))
158 scalar KRT =
sqrt(helpVariable/(3.0*sigma + vSmall));
165 if ((tc > 0) || (lambdaRT < d) )
171 scalar tauRT = cTau_/(omegaRT + vSmall);
174 if ((tc > tauRT) && (lambdaRT < d))
178 scalar nDrops = d/lambdaRT;
185 if (weGas > weberLimit_)
187 scalar fraction = dt/tauKH;
190 d = (fraction*dc + d)/(1.0 + fraction);
193 scalar ms0 = mass0*(1.0 -
pow3(d/d0));
196 if (ms/averageParcelMass > msLimit_)
207 scalar de3 = d*d*(dc - d);
209 pow3(be3/(3.0*ae3)) - be3*ce3/(6.0*ae3*ae3) + de3/(2.0*ae3);
210 scalar pe3 = (3.0*ae3*ce3 - be3*be3)/(9.0*ae3*ae3);
211 scalar D3 = qe3*qe3 + pe3*pe3*pe3;
213 if (D3 < 0) br3 =
false;
218 scalar ue3 =
cbrt(-qe3 + D3);
219 scalar ve3 =
cbrt(-qe3 - D3);
220 scalar dParenDrops = ue3 + ve3 - be3/3.;
221 scalar mc = nParticle*(
pow3(d) -
pow3(dParenDrops));
222 scalar nChildDrops = mc/
pow3(dc);
224 if (nChildDrops >= nParticle)
230 massChild = mc*rhopi6;
239 else if (KHindex < 0.5)
246 scalar diameterLargerDrop =
cbrt(1.5*d*d*lengthScale);
247 d = diameterLargerDrop;
253 scalar massDrop =
pow3(d)*rhopi6;
254 nParticle = mass/massDrop;
secondary breakup model which uses the Kelvin-Helmholtz instability theory to predict the 'stripped' ...
A list of keyword definitions, which are a keyword followed by any number of values (e...
virtual bool update(const scalar dt, const vector &g, scalar &d, scalar &tc, scalar &ms, scalar &nParticle, scalar &KHindex, scalar &y, scalar &yDot, const scalar d0, const scalar rho, const scalar mu, const scalar sigma, const vector &U, const scalar rhoc, const scalar muc, const vector &Urel, const scalar Urmag, const scalar tMom, scalar &dChild, scalar &massChild)
Update the parcel diameter.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionedScalar & sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
dimensionedScalar sqrt(const dimensionedScalar &ds)
const CloudType & owner() const
Return const access to the owner cloud.
dimensionedScalar cbrt(const dimensionedScalar &ds)
const scalar twoPi(2 *pi)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
ReitzKHRT(const dictionary &, CloudType &)
Construct from dictionary.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionedScalar & mu
Atomic mass unit.
dimensioned< scalar > mag(const dimensioned< Type > &)
Templated break-up model class.
Templated base class for dsmc cloud.
virtual ~ReitzKHRT()
Destructor.