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];
86 template<
class CloudType>
95 const scalar phiSi =
twoPi*rndGen_.sample01<scalar>();
98 const scalar thetaSi =
pi/180.0*(rndGen_.sample01<scalar>()*(50 - 5) + 5);
102 const scalar dcorr =
cos(thetaSi);
103 const vector normal = alpha*(tanVec1*
cos(phiSi) + tanVec2*
sin(phiSi));
107 return dirVec/
mag(dirVec);
111 template<
class CloudType>
124 Info<<
"Parcel " << p.origId() <<
" absorbInteraction" <<
endl;
131 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
137 const vector Un = nf*(Urel & nf);
140 const vector Ut = Urel - Un;
152 this->nParcelsTransferred()++;
154 keepParticle =
false;
158 template<
class CloudType>
169 Info<<
"Parcel " << p.origId() <<
" bounceInteraction" <<
endl;
176 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
182 p.U() -= 2.0*nf*(Urel & nf);
188 template<
class CloudType>
200 Info<<
"Parcel " << p.origId() <<
" drySplashInteraction" <<
endl;
206 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
210 const scalar pc = thermo_.thermo().p()[p.cell()];
213 const scalar m = p.mass()*p.nParticle();
214 const scalar
rho = p.rho();
215 const scalar d = p.d();
216 const scalar sigma = liq.
sigma(pc, p.T());
217 const scalar
mu = liq.
mu(pc, p.T());
219 const vector Un = nf*(Urel & nf);
222 const scalar La = rho*sigma*d/
sqr(mu);
228 const scalar Wec = Adry_*
pow(La, -0.183);
232 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
237 const scalar mRatio = 0.2 + 0.6*rndGen_.sample01<scalar>();
239 (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
244 template<
class CloudType>
256 Info<<
"Parcel " << p.origId() <<
" wetSplashInteraction" <<
endl;
262 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
266 const scalar pc = thermo_.thermo().p()[p.cell()];
269 const scalar m = p.mass()*p.nParticle();
270 const scalar
rho = p.rho();
271 const scalar d = p.d();
273 const scalar sigma = liq.
sigma(pc, p.T());
274 const scalar
mu = liq.
mu(pc, p.T());
276 const vector Un = nf*(Urel & nf);
277 const vector Ut = Urel - Un;
280 const scalar La = rho*sigma*d/
sqr(mu);
286 const scalar Wec = Awet_*
pow(La, -0.183);
290 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
292 else if ((We >= 2) && (We < 20))
295 const scalar theta =
pi/2 -
acos(U/
mag(U) & nf);
298 const scalar
epsilon = 0.993 - theta*(1.76 - theta*(1.56 - theta*0.49));
301 U = -epsilon*(Un) + 5.0/7.0*(Ut);
306 else if ((We >= 20) && (We < Wec))
308 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
314 const scalar mRatio = 0.2 + 0.9*rndGen_.sample01<scalar>();
316 (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
321 template<
class CloudType>
337 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
342 const vector tanVec2 = nf^tanVec1;
345 const scalar np = p.nParticle();
346 const scalar m = p.mass()*np;
347 const scalar d = p.d();
349 const vector Un = nf*(Urel & nf);
350 const vector Ut = Urel - Un;
351 const vector& posC = mesh.
C()[p.cell()];
355 const scalar mSplash = m*mRatio;
358 const scalar Ns = 5.0*(We/Wec - 1.0);
361 const scalar dBarSplash = 1/
cbrt(6.0)*
cbrt(mRatio/Ns)*d + rootVSmall;
364 const scalar dMax = 0.9*
cbrt(mRatio)*d;
365 const scalar dMin = 0.1*dMax;
366 const scalar
K =
exp(-dMin/dBarSplash) -
exp(-dMax/dBarSplash);
369 scalar ESigmaSec = 0;
376 const scalar
y = rndGen_.sample01<scalar>();
377 dNew[i] = -dBarSplash*
log(
exp(-dMin/dBarSplash) - y*K);
378 npNew[i] = mRatio*np*
pow3(d)/
pow3(dNew[i])/parcelsPerSplash_;
379 ESigmaSec += npNew[i]*sigma*p.areaS(dNew[i]);
383 const scalar EKIn = 0.5*m*
magSqr(Un);
386 const scalar ESigmaIn = np*sigma*p.areaS(d);
389 const scalar Ed =
max(0.8*EKIn, np*Wec/12*
pi*sigma*
sqr(d));
392 const scalar EKs = EKIn + ESigmaIn - ESigmaSec - Ed;
397 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
402 const scalar logD =
log(d);
403 const scalar coeff2 =
log(dNew[0]) - logD + rootVSmall;
407 coeff1 +=
sqr(
log(dNew[i]) - logD);
411 const scalar magUns0 =
412 sqrt(2.0*parcelsPerSplash_*EKs/mSplash/(1.0 + coeff1/
sqr(coeff2)));
417 const vector dirVec = splashDirection(tanVec1, tanVec2, -nf);
422 pPtr->origId() = pPtr->getNewParticleID();
426 if (splashParcelType_ >= 0)
428 pPtr->typeId() = splashParcelType_;
432 pPtr->track(0.5*rndGen_.sample01<scalar>()*(posC - posCf), 0);
434 pPtr->nParticle() = npNew[i];
438 pPtr->U() = dirVec*(
mag(Cf_*Ut) + magUns0*(
log(dNew[i]) - logD)/coeff2);
444 this->owner().addParticle(pPtr);
451 const scalar mDash = m - mSplash;
452 absorbInteraction(filmModel, p, pp, facei, mDash, keepParticle);
458 template<
class CloudType>
466 rndGen_(owner.rndGen()),
469 owner.db().objectRegistry::template lookupObject<SLGThermo>(
"SLGThermo")
475 interactionTypeEnum(this->coeffDict().
lookup(
"interactionType"))
478 splashParcelType_(0),
479 parcelsPerSplash_(0),
485 Info<<
" Applying " << interactionTypeStr(interactionType_)
486 <<
" interaction model" <<
endl;
488 if (interactionType_ == itSplashBai)
490 this->coeffDict().lookup(
"deltaWet") >> deltaWet_;
492 this->coeffDict().lookupOrDefault(
"splashParcelType", -1);
494 this->coeffDict().lookupOrDefault(
"parcelsPerSplash", 2);
495 this->coeffDict().lookup(
"Adry") >> Adry_;
496 this->coeffDict().lookup(
"Awet") >> Awet_;
497 this->coeffDict().lookup(
"Cf") >> Cf_;
502 template<
class CloudType>
526 template<
class CloudType>
533 template<
class CloudType>
545 this->owner().db().time().objectRegistry::template
549 "surfaceFilmProperties" 559 switch (interactionType_)
563 bounceInteraction(p, pp, facei, keepParticle);
569 const scalar m = p.nParticle()*p.mass();
570 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
576 bool dry = this->deltaFilmPatch_[
patchi][facei] < deltaWet_;
580 drySplashInteraction(filmModel, p, pp, facei, keepParticle);
584 wetSplashInteraction(filmModel, p, pp, facei, keepParticle);
592 <<
"Unknown interaction type enumeration" 606 template<
class CloudType>
609 const label filmPatchi,
610 const label primaryPatchi,
622 filmModel.
toPrimary(filmPatchi, TFilmPatch_);
625 filmModel.
toPrimary(filmPatchi, CpFilmPatch_);
629 template<
class CloudType>
633 const label filmFacei
639 p.T() = TFilmPatch_[filmFacei];
640 p.Cp() = CpFilmPatch_[filmFacei];
644 template<
class CloudType>
649 label nSplash0 = this->
template getModelProperty<label>(
"nParcelsSplashed");
653 os <<
" New film splash parcels = " << nSplashTotal <<
endl;
655 if (this->writeTime())
657 this->setModelProperty(
"nParcelsSplashed", nSplashTotal);
658 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.
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.
scalar Cf_
Skin friction typically in the range 0.6 < Cf < 0.8.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedScalar log(const dimensionedScalar &ds)
Base class for surface film models.
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 > &)
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
bool isRegionPatch(const label primaryPatchi) const
Return true if patchi on the primary region is a coupled.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
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.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
scalarList TFilmPatch_
Film temperature / patch face.
CGAL::Exact_predicates_exact_constructions_kernel K
void drySplashInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, const parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with dry surface.
Macros for easy insertion into run-time selection tables.
virtual void info(Ostream &os)
Write surface film info to stream.
virtual scalar sigma(scalar p, scalar T) const =0
Surface tension [N/m].
void toPrimary(const label regionPatchi, List< Type > ®ionField) const
Convert a local region field to the primary region.
interactionType interactionType_
Interaction type enumeration.
void absorbInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, const parcelType &p, const polyPatch &pp, const label facei, const scalar mass, bool &keepParticle)
Absorb parcel into film.
stressControl lookup("compactNormalStress") >> compactNormalStress
Form normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
interactionType interactionTypeEnum(const word &it) const
static wordList interactionTypeNames_
Word descriptions of interaction type names.
A class for handling words, derived from string.
dimensionedScalar cbrt(const dimensionedScalar &ds)
The thermophysical properties of a liquid.
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmRegionModel &filmModel)
Cache the film fields in preparation for injection.
const Field< PointType > & faceNormals() const
Return face normals for patch.
errorManip< error > abort(error &err)
const SLGThermo & thermo_
Reference to the cloud thermo package.
const scalar twoPi(2 *pi)
vector splashDirection(const vector &tanVec1, const vector &tanVec2, const vector &nf) const
Return splashed parcel direction.
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 & Ts() const =0
Return the film surface temperature [K].
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
void bounceInteraction(parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle) const
Bounce parcel (flip parcel normal velocity)
Vector< Cmpt > perpendicular(const Vector< Cmpt > &v)
const dimensionedScalar mu
Atomic mass unit.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
CloudType::parcelType parcelType
Convenience typedef to the cloud's parcel type.
void splashInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, 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.
Mesh data needed to do the Finite Volume discretisation.
label index() const
Return the index of this patch in the boundaryMesh.
label parcelsPerSplash_
Number of new parcels resulting from splash event.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmRegionModel &)
Cache the film fields in preparation for injection.
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.
Random & rndGen_
Reference to the cloud random number generator.
virtual const volScalarField & Cp() const =0
Return the film specific heat capacity [J/kg/K].
const volVectorField & C() const
Return cell centres as volVectorField.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A patch is a list of labels that address the faces in the global face list.
word interactionTypeStr(const interactionType &it) const
virtual scalar mu(scalar p, scalar T) const =0
Liquid viscosity [Pa s].
virtual bool transferParcel(parcelType &p, const polyPatch &pp, bool &keepParticle)
Transfer parcel from cloud to surface film.
void wetSplashInteraction(regionModels::surfaceFilmModels::surfaceFilmRegionModel &, parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with wetted surface.
scalar deltaWet_
Film thickness beyond which patch is assumed to be wet.
Templated base class for dsmc cloud.
label splashParcelType_
Splash parcel type label - id assigned to identify parcel for.
scalarList CpFilmPatch_
Film specific heat capacity / patch face.
label whichFace(const label l) const
Return label of face in patch from global face label.