35 template<
class CloudType>
40 "(absorb bounce splashBai)" 47 template<
class CloudType>
51 forAll(interactionTypeNames_, i)
53 if (interactionTypeNames_[i] == it)
60 <<
"Unknown interaction type " << it
61 <<
". Valid interaction types include: " << interactionTypeNames_
68 template<
class CloudType>
74 if (it >= interactionTypeNames_.size())
80 return interactionTypeNames_[it];
84 template<
class CloudType>
91 scalar magTangent = 0.0;
93 while (magTangent < SMALL)
96 tangent = vTest - (vTest & v)*v;
97 magTangent =
mag(tangent);
100 return tangent/magTangent;
104 template<
class CloudType>
113 const scalar phiSi =
twoPi*rndGen_.sample01<scalar>();
116 const scalar thetaSi =
pi/180.0*(rndGen_.sample01<scalar>()*(50 - 5) + 5);
120 const scalar dcorr =
cos(thetaSi);
125 return dirVec/
mag(dirVec);
129 template<
class CloudType>
142 Info<<
"Parcel " << p.origId() <<
" absorbInteraction" <<
endl;
149 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
155 const vector Un = nf*(Urel & nf);
158 const vector Ut = Urel - Un;
170 this->nParcelsTransferred()++;
172 keepParticle =
false;
176 template<
class CloudType>
187 Info<<
"Parcel " << p.origId() <<
" bounceInteraction" <<
endl;
194 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
200 p.U() -= 2.0*nf*(Urel & nf);
206 template<
class CloudType>
218 Info<<
"Parcel " << p.origId() <<
" drySplashInteraction" <<
endl;
224 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
228 const scalar pc = thermo_.thermo().p()[p.cell()];
231 const scalar m = p.mass()*p.nParticle();
232 const scalar
rho = p.rho();
233 const scalar d = p.d();
234 const scalar sigma = liq.
sigma(pc, p.T());
235 const scalar
mu = liq.
mu(pc, p.T());
237 const vector Un = nf*(Urel & nf);
240 const scalar La = rho*sigma*d/
sqr(mu);
246 const scalar Wec = Adry_*
pow(La, -0.183);
250 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
255 const scalar mRatio = 0.2 + 0.6*rndGen_.sample01<scalar>();
257 (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
262 template<
class CloudType>
274 Info<<
"Parcel " << p.origId() <<
" wetSplashInteraction" <<
endl;
280 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
284 const scalar pc = thermo_.thermo().p()[p.cell()];
287 const scalar m = p.mass()*p.nParticle();
288 const scalar
rho = p.rho();
289 const scalar d = p.d();
291 const scalar sigma = liq.
sigma(pc, p.T());
292 const scalar
mu = liq.
mu(pc, p.T());
294 const vector Un = nf*(Urel & nf);
295 const vector Ut = Urel - Un;
298 const scalar La = rho*sigma*d/
sqr(mu);
304 const scalar Wec = Awet_*
pow(La, -0.183);
308 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
310 else if ((We >= 1) && (We < 20))
313 const scalar theta =
pi/2 -
acos(U/
mag(U) & nf);
316 const scalar
epsilon = 0.993 - theta*(1.76 - theta*(1.56 - theta*0.49));
319 U = -epsilon*(Un) + 5/7*(Ut);
324 else if ((We >= 20) && (We < Wec))
326 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
332 const scalar mRatio = 0.2 + 0.9*rndGen_.sample01<scalar>();
334 (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
339 template<
class CloudType>
355 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
359 const vector tanVec1 = tangentVector(nf);
360 const vector tanVec2 = nf^tanVec1;
363 const scalar np = p.nParticle();
364 const scalar m = p.mass()*np;
365 const scalar d = p.d();
367 const vector Un = nf*(Urel & nf);
368 const vector Ut = Urel - Un;
369 const vector& posC = mesh.
C()[p.cell()];
373 const scalar mSplash = m*mRatio;
376 const scalar Ns = 5.0*(We/Wec - 1.0);
379 const scalar dBarSplash = 1/
cbrt(6.0)*
cbrt(mRatio/Ns)*d + ROOTVSMALL;
382 const scalar dMax = 0.9*
cbrt(mRatio)*d;
383 const scalar dMin = 0.1*dMax;
384 const scalar
K =
exp(-dMin/dBarSplash) -
exp(-dMax/dBarSplash);
387 scalar ESigmaSec = 0;
394 const scalar
y = rndGen_.sample01<scalar>();
395 dNew[i] = -dBarSplash*
log(
exp(-dMin/dBarSplash) - y*K);
396 npNew[i] = mRatio*np*
pow3(d)/
pow3(dNew[i])/parcelsPerSplash_;
397 ESigmaSec += npNew[i]*sigma*p.areaS(dNew[i]);
401 const scalar EKIn = 0.5*m*
magSqr(Urel);
404 const scalar ESigmaIn = np*sigma*p.areaS(d);
407 const scalar Ed =
max(0.8*EKIn, np*Wec/12*
pi*sigma*
sqr(d));
410 const scalar EKs = EKIn + ESigmaIn - ESigmaSec - Ed;
415 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
420 const scalar logD =
log(d);
421 const scalar coeff2 =
log(dNew[0]) - logD + ROOTVSMALL;
425 coeff1 +=
sqr(
log(dNew[i]) - logD);
429 const scalar magUns0 =
430 sqrt(2.0*parcelsPerSplash_*EKs/mSplash/(1.0 + coeff1/
sqr(coeff2)));
435 const vector dirVec = splashDirection(tanVec1, tanVec2, -nf);
440 pPtr->origId() = pPtr->getNewParticleID();
444 if (splashParcelType_ >= 0)
446 pPtr->typeId() = splashParcelType_;
450 pPtr->position() += 0.5*rndGen_.sample01<scalar>()*(posC - posCf);
452 pPtr->nParticle() = npNew[i];
456 pPtr->U() = dirVec*(
mag(Cf_*Ut) + magUns0*(
log(dNew[i]) - logD)/coeff2);
462 this->owner().addParticle(pPtr);
469 const scalar mDash = m - mSplash;
470 absorbInteraction(filmModel, p, pp, facei, mDash, keepParticle);
476 template<
class CloudType>
484 rndGen_(owner.rndGen()),
487 owner.db().objectRegistry::template lookupObject<SLGThermo>(
"SLGThermo")
493 interactionTypeEnum(this->coeffDict().
lookup(
"interactionType"))
496 splashParcelType_(0),
497 parcelsPerSplash_(0),
503 Info<<
" Applying " << interactionTypeStr(interactionType_)
504 <<
" interaction model" <<
endl;
506 if (interactionType_ == itSplashBai)
508 this->coeffDict().lookup(
"deltaWet") >> deltaWet_;
510 this->coeffDict().lookupOrDefault(
"splashParcelType", -1);
512 this->coeffDict().lookupOrDefault(
"parcelsPerSplash", 2);
513 this->coeffDict().lookup(
"Adry") >> Adry_;
514 this->coeffDict().lookup(
"Awet") >> Awet_;
515 this->coeffDict().lookup(
"Cf") >> Cf_;
520 template<
class CloudType>
544 template<
class CloudType>
551 template<
class CloudType>
563 this->owner().db().time().objectRegistry::template
564 lookupObject<regionModels::surfaceFilmModels::surfaceFilmModel>
566 "surfaceFilmProperties" 576 switch (interactionType_)
580 bounceInteraction(p, pp, facei, keepParticle);
586 const scalar m = p.nParticle()*p.mass();
587 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
593 bool dry = this->deltaFilmPatch_[
patchi][facei] < deltaWet_;
597 drySplashInteraction(filmModel, p, pp, facei, keepParticle);
601 wetSplashInteraction(filmModel, p, pp, facei, keepParticle);
609 <<
"Unknown interaction type enumeration" 623 template<
class CloudType>
626 const label filmPatchi,
627 const label primaryPatchi,
639 filmModel.
toPrimary(filmPatchi, TFilmPatch_);
642 filmModel.
toPrimary(filmPatchi, CpFilmPatch_);
646 template<
class CloudType>
650 const label filmFacei
656 p.T() = TFilmPatch_[filmFacei];
657 p.Cp() = CpFilmPatch_[filmFacei];
661 template<
class CloudType>
666 label nSplash0 = this->
template getModelProperty<label>(
"nParcelsSplashed");
670 os <<
" New film splash parcels = " << nSplashTotal <<
endl;
672 if (this->writeTime())
674 this->setModelProperty(
"nParcelsSplashed", nSplashTotal);
675 nParcelsSplashed_ = 0;
scalar Awet_
Wet surface roughness coefficient.
virtual void info(Ostream &os)
Write surface film info to stream.
Thermo parcel surface film model.
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
scalar Cf_
Skin friction typically in the range 0.6 < Cf < 0.8.
void wetSplashInteraction(regionModels::surfaceFilmModels::surfaceFilmModel &filmModel, parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with wetted surface.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void bounceInteraction(parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle) const
Bounce parcel (flip parcel normal velocity)
bool isRegionPatch(const label primaryPatchi) const
Return true if patchi on the primary region is a coupled.
label whichFace(const label l) const
Return label of face in patch from global face label.
void splashInteraction(regionModels::surfaceFilmModels::surfaceFilmModel &filmModel, const parcelType &p, const polyPatch &pp, const label facei, const scalar mRatio, const scalar We, const scalar Wec, const scalar sigma, bool &keepParticle)
Bai parcel splash interaction model.
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
ThermoSurfaceFilm(const dictionary &dict, CloudType &owner)
Construct from components.
scalar Adry_
Dry surface roughness coefficient.
virtual ~ThermoSurfaceFilm()
Destructor.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual scalar mu(scalar p, scalar T) const
Liquid viscosity [Pa s].
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
const volVectorField & C() const
Return cell centres as volVectorField.
scalarList TFilmPatch_
Film temperature / patch face.
CGAL::Exact_predicates_exact_constructions_kernel K
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmModel &filmModel)
Cache the film fields in preparation for injection.
Macros for easy insertion into run-time selection tables.
virtual void info(Ostream &os)
Write surface film info to stream.
vector splashDirection(const vector &tanVec1, const vector &tanVec2, const vector &nf) const
Return splashed parcel direction.
interactionType interactionType_
Interaction type enumeration.
stressControl lookup("compactNormalStress") >> compactNormalStress
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static wordList interactionTypeNames_
Word descriptions of interaction type names.
A class for handling words, derived from string.
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
dimensionedScalar cbrt(const dimensionedScalar &ds)
void absorbInteraction(regionModels::surfaceFilmModels::surfaceFilmModel &filmModel, const parcelType &p, const polyPatch &pp, const label facei, const scalar mass, bool &keepParticle)
Absorb parcel into film.
virtual scalar sigma(scalar p, scalar T) const
Surface tension [N/m].
The thermophysical properties of a liquidProperties.
virtual const volScalarField & Ts() const =0
Return the film surface temperature [K].
vector tangentVector(const vector &v) const
Return a vector tangential to input vector, v.
errorManip< error > abort(error &err)
const SLGThermo & thermo_
Reference to the cloud thermo package.
const scalar twoPi(2 *pi)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
dimensionedScalar sin(const dimensionedScalar &ds)
virtual const volScalarField & Cp() const =0
Return the film specific heat capacity [J/kg/K].
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource)=0
External hook to add sources to the film.
interactionType interactionTypeEnum(const word &it) const
Templated wall surface film model class.
cachedRandom & rndGen_
Reference to the cloud random number generator.
A normal distribution model.
word interactionTypeStr(const interactionType &it) const
const dimensionedScalar mu
Atomic mass unit.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
const Field< PointType > & faceNormals() const
Return face normals for patch.
void toPrimary(const label regionPatchi, List< Type > ®ionField) const
Convert a local region field to the primary region.
void drySplashInteraction(regionModels::surfaceFilmModels::surfaceFilmModel &filmModel, const parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with dry surface.
CloudType::parcelType parcelType
Convenience typedef to the cloud's parcel type.
Mesh data needed to do the Finite Volume discretisation.
label parcelsPerSplash_
Number of new parcels resulting from splash event.
Base class for surface film models.
dimensioned< scalar > mag(const dimensioned< Type > &)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
label nParcelsSplashed_
Counter for number of new splash parcels.
A patch is a list of labels that address the faces in the global face list.
virtual bool transferParcel(parcelType &p, const polyPatch &pp, bool &keepParticle)
Transfer parcel from cloud to surface film.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
scalar deltaWet_
Film thickness beyond which patch is assumed to be wet.
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
Templated base class for dsmc cloud.
label index() const
Return the index of this patch in the boundaryMesh.
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmModel &filmModel)
Cache the film fields in preparation for injection.
label splashParcelType_
Splash parcel type label - id assigned to identify parcel for.
scalarList CpFilmPatch_
Film specific heat capacity / patch face.