32 template<
class CloudType>
40 weCorrCoeff_(this->coeffDict().getScalar(
"weCorrCoeff")),
41 weBuCrit_(this->coeffDict().getScalar(
"weBuCrit")),
42 weBuBag_(this->coeffDict().getScalar(
"weBuBag")),
43 weBuMM_(this->coeffDict().getScalar(
"weBuMM")),
44 ohnCoeffCrit_(this->coeffDict().getScalar(
"ohnCoeffCrit")),
45 ohnCoeffBag_(this->coeffDict().getScalar(
"ohnCoeffBag")),
46 ohnCoeffMM_(this->coeffDict().getScalar(
"ohnCoeffMM")),
47 ohnExpCrit_(this->coeffDict().getScalar(
"ohnExpCrit")),
48 ohnExpBag_(this->coeffDict().getScalar(
"ohnExpBag")),
49 ohnExpMM_(this->coeffDict().getScalar(
"ohnExpMM")),
50 cInit_(this->coeffDict().getScalar(
"Cinit")),
51 c1_(this->coeffDict().getScalar(
"C1")),
52 c2_(this->coeffDict().getScalar(
"C2")),
53 c3_(this->coeffDict().getScalar(
"C3")),
54 cExp1_(this->coeffDict().getScalar(
"Cexp1")),
55 cExp2_(this->coeffDict().getScalar(
"Cexp2")),
56 cExp3_(this->coeffDict().getScalar(
"Cexp3")),
57 weConst_(this->coeffDict().getScalar(
"Weconst")),
58 weCrit1_(this->coeffDict().getScalar(
"Wecrit1")),
59 weCrit2_(this->coeffDict().getScalar(
"Wecrit2")),
60 coeffD_(this->coeffDict().getScalar(
"CoeffD")),
61 onExpD_(this->coeffDict().getScalar(
"OnExpD")),
62 weExpD_(this->coeffDict().getScalar(
"WeExpD")),
63 mu_(this->coeffDict().getScalar(
"mu")),
64 sigma_(this->coeffDict().getScalar(
"sigma")),
65 d32Coeff_(this->coeffDict().getScalar(
"d32Coeff")),
66 cDmaxBM_(this->coeffDict().getScalar(
"cDmaxBM")),
67 cDmaxS_(this->coeffDict().getScalar(
"cDmaxS")),
68 corePerc_(this->coeffDict().getScalar(
"corePerc"))
72 template<
class CloudType>
76 weCorrCoeff_(bum.weCorrCoeff_),
77 weBuCrit_(bum.weBuCrit_),
78 weBuBag_(bum.weBuBag_),
80 ohnCoeffCrit_(bum.ohnCoeffCrit_),
81 ohnCoeffBag_(bum.ohnCoeffBag_),
82 ohnCoeffMM_(bum.ohnCoeffMM_),
83 ohnExpCrit_(bum.ohnExpCrit_),
84 ohnExpBag_(bum.ohnExpBag_),
85 ohnExpMM_(bum.ohnExpMM_),
93 weConst_(bum.weConst_),
94 weCrit1_(bum.weCrit1_),
95 weCrit2_(bum.weCrit2_),
101 d32Coeff_(bum.d32Coeff_),
102 cDmaxBM_(bum.cDmaxBM_),
103 cDmaxS_(bum.cDmaxS_),
104 corePerc_(bum.corePerc_)
110 template<
class CloudType>
117 template<
class CloudType>
145 bool addChild =
false;
147 scalar d03 =
pow3(d);
149 scalar mass0 = nParticle*rhopi6*d03;
152 scalar weGas = 0.5*rhoc*
sqr(Urmag)*d/
sigma;
156 scalar reLiquid = 0.5*Urmag*d/
mu;
157 scalar ohnesorge =
sqrt(weLiquid)/(reLiquid + VSMALL);
159 scalar weGasCorr = weGas/(1.0 + weCorrCoeff_*ohnesorge);
165 scalar rChar = Urmag/d*
sqrt(rhoc/
rho);
169 if (tc*rChar < SMALL)
175 scalar tChar = 1/rChar;
177 scalar tFirst = cInit_*tChar;
180 scalar tCharSecond = 0;
183 bool multimode =
false;
188 if (weGas > weConst_)
190 if (weGas < weCrit1_)
192 tCharSecond = c1_*
pow((weGas - weConst_), cExp1_);
194 else if (weGas >= weCrit1_ && weGas <= weCrit2_)
196 tCharSecond = c2_*
pow((weGas - weConst_), cExp2_);
200 tCharSecond = c3_*
pow((weGas - weConst_), cExp3_);
204 scalar weC = weBuCrit_*(1.0 + ohnCoeffCrit_*
pow(ohnesorge, ohnExpCrit_));
205 scalar weB = weBuBag_*(1.0 + ohnCoeffBag_*
pow(ohnesorge, ohnExpBag_));
206 scalar weMM = weBuMM_*(1.0 + ohnCoeffMM_*
pow(ohnesorge, ohnExpMM_));
208 if (weGas > weC && weGas < weB)
213 if (weGas >= weB && weGas <= weMM)
223 tSecond = tCharSecond*tChar;
225 scalar tBreakUP = tFirst + tSecond;
228 scalar d32 = coeffD_*d*
pow(ohnesorge, onExpD_)*
pow(weGasCorr, weExpD_);
230 if (bag || multimode)
232 scalar d05 = d32Coeff_ * d32;
240 x = cDmaxBM_*
rndGen.sample01<scalar>();
242 yGuess =
rndGen.sample01<scalar>();
247 *
exp(-0.5*
sqr((
x - mu_)/sigma_));
261 scalar dC = weConst_*
sigma/(rhoc*
sqr(Urmag));
262 scalar d32Red = 4.0*(d32*dC)/(5.0*dC - d32);
264 scalar d05 = d32Coeff_ * d32Red;
272 x = cDmaxS_*
rndGen.sample01<scalar>();
274 yGuess =
rndGen.sample01<scalar>();
279 *
exp(-0.5*
sqr((
x - mu_)/sigma_));
289 massChild = corePerc_*mass0;
298 nParticle = mass/(rhopi6*
pow3(d));
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 properties.
SHF(const dictionary &, CloudType &)
Construct from dictionary.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
Secondary Breakup Model to take account of the different breakup regimes, bag, molutimode, shear....
dimensionedScalar exp(const dimensionedScalar &ds)
constexpr scalar twoPi(2 *M_PI)
constexpr scalar pi(M_PI)
const uniformDimensionedVectorField & g
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
const dimensionedScalar mu
Atomic mass unit.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
Templated break-up model class.
Templated base class for dsmc cloud.
virtual ~SHF()
Destructor.