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>
93 scalar magTangent = 0.0;
95 while (magTangent < SMALL)
98 tangent = vTest - (vTest & v)*v;
99 magTangent =
mag(tangent);
102 return tangent/magTangent;
106 template<
class CloudType>
115 const scalar phiSi =
twoPi*rndGen_.sample01<scalar>();
118 const scalar thetaSi =
pi/180.0*(rndGen_.sample01<scalar>()*(50 - 5) + 5);
122 const scalar dcorr =
cos(thetaSi);
127 return dirVec/
mag(dirVec);
131 template<
class CloudType>
144 Info<<
"Parcel " << p.origId() <<
" absorbInteraction" <<
endl;
151 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
157 const vector Un = nf*(Urel & nf);
160 const vector Ut = Urel - Un;
172 this->nParcelsTransferred()++;
174 keepParticle =
false;
178 template<
class CloudType>
189 Info<<
"Parcel " << p.origId() <<
" bounceInteraction" <<
endl;
196 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
202 p.U() -= 2.0*nf*(Urel & nf);
208 template<
class CloudType>
220 Info<<
"Parcel " << p.origId() <<
" drySplashInteraction" <<
endl;
226 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
230 const scalar pc = thermo_.thermo().p()[p.cell()];
233 const scalar m = p.mass()*p.nParticle();
234 const scalar
rho = p.rho();
235 const scalar d = p.d();
236 const scalar sigma = liq.
sigma(pc, p.T());
237 const scalar
mu = liq.
mu(pc, p.T());
239 const vector Un = nf*(Urel & nf);
242 const scalar La = rho*sigma*d/
sqr(mu);
248 const scalar Wec = Adry_*
pow(La, -0.183);
252 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
257 const scalar mRatio = 0.2 + 0.6*rndGen_.sample01<scalar>();
259 (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
264 template<
class CloudType>
276 Info<<
"Parcel " << p.origId() <<
" wetSplashInteraction" <<
endl;
282 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
286 const scalar pc = thermo_.thermo().p()[p.cell()];
289 const scalar m = p.mass()*p.nParticle();
290 const scalar
rho = p.rho();
291 const scalar d = p.d();
293 const scalar sigma = liq.
sigma(pc, p.T());
294 const scalar
mu = liq.
mu(pc, p.T());
296 const vector Un = nf*(Urel & nf);
297 const vector Ut = Urel - Un;
300 const scalar La = rho*sigma*d/
sqr(mu);
306 const scalar Wec = Awet_*
pow(La, -0.183);
310 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
312 else if ((We >= 2) && (We < 20))
315 const scalar theta =
pi/2 -
acos(U/
mag(U) & nf);
318 const scalar
epsilon = 0.993 - theta*(1.76 - theta*(1.56 - theta*0.49));
321 U = -epsilon*(Un) + 5.0/7.0*(Ut);
326 else if ((We >= 20) && (We < Wec))
328 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
334 const scalar mRatio = 0.2 + 0.9*rndGen_.sample01<scalar>();
336 (filmModel, p, pp, facei, mRatio, We, Wec, sigma, keepParticle);
341 template<
class CloudType>
357 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
361 const vector tanVec1 = tangentVector(nf);
362 const vector tanVec2 = nf^tanVec1;
365 const scalar np = p.nParticle();
366 const scalar m = p.mass()*np;
367 const scalar d = p.d();
369 const vector Un = nf*(Urel & nf);
370 const vector Ut = Urel - Un;
371 const vector& posC = mesh.
C()[p.cell()];
375 const scalar mSplash = m*mRatio;
378 const scalar Ns = 5.0*(We/Wec - 1.0);
381 const scalar dBarSplash = 1/
cbrt(6.0)*
cbrt(mRatio/Ns)*d + ROOTVSMALL;
384 const scalar dMax = 0.9*
cbrt(mRatio)*d;
385 const scalar dMin = 0.1*dMax;
386 const scalar
K =
exp(-dMin/dBarSplash) -
exp(-dMax/dBarSplash);
389 scalar ESigmaSec = 0;
396 const scalar
y = rndGen_.sample01<scalar>();
397 dNew[i] = -dBarSplash*
log(
exp(-dMin/dBarSplash) - y*K);
398 npNew[i] = mRatio*np*
pow3(d)/
pow3(dNew[i])/parcelsPerSplash_;
399 ESigmaSec += npNew[i]*sigma*p.areaS(dNew[i]);
403 const scalar EKIn = 0.5*m*
magSqr(Urel);
406 const scalar ESigmaIn = np*sigma*p.areaS(d);
409 const scalar Ed =
max(0.8*EKIn, np*Wec/12*
pi*sigma*
sqr(d));
412 const scalar EKs = EKIn + ESigmaIn - ESigmaSec - Ed;
417 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
422 const scalar logD =
log(d);
423 const scalar coeff2 =
log(dNew[0]) - logD + ROOTVSMALL;
427 coeff1 +=
sqr(
log(dNew[i]) - logD);
431 const scalar magUns0 =
432 sqrt(2.0*parcelsPerSplash_*EKs/mSplash/(1.0 + coeff1/
sqr(coeff2)));
437 const vector dirVec = splashDirection(tanVec1, tanVec2, -nf);
442 pPtr->origId() = pPtr->getNewParticleID();
446 if (splashParcelType_ >= 0)
448 pPtr->typeId() = splashParcelType_;
452 pPtr->track(0.5*rndGen_.sample01<scalar>()*(posC - posCf), 0);
454 pPtr->nParticle() = npNew[i];
458 pPtr->U() = dirVec*(
mag(Cf_*Ut) + magUns0*(
log(dNew[i]) - logD)/coeff2);
464 this->owner().addParticle(pPtr);
471 const scalar mDash = m - mSplash;
472 absorbInteraction(filmModel, p, pp, facei, mDash, keepParticle);
478 template<
class CloudType>
486 rndGen_(owner.rndGen()),
489 owner.db().objectRegistry::template lookupObject<SLGThermo>(
"SLGThermo")
495 interactionTypeEnum(this->coeffDict().
lookup(
"interactionType"))
498 splashParcelType_(0),
499 parcelsPerSplash_(0),
505 Info<<
" Applying " << interactionTypeStr(interactionType_)
506 <<
" interaction model" <<
endl;
508 if (interactionType_ == itSplashBai)
510 this->coeffDict().lookup(
"deltaWet") >> deltaWet_;
512 this->coeffDict().lookupOrDefault(
"splashParcelType", -1);
514 this->coeffDict().lookupOrDefault(
"parcelsPerSplash", 2);
515 this->coeffDict().lookup(
"Adry") >> Adry_;
516 this->coeffDict().lookup(
"Awet") >> Awet_;
517 this->coeffDict().lookup(
"Cf") >> Cf_;
522 template<
class CloudType>
546 template<
class CloudType>
553 template<
class CloudType>
565 this->owner().db().time().objectRegistry::template
566 lookupObject<regionModels::surfaceFilmModels::surfaceFilmModel>
568 "surfaceFilmProperties" 578 switch (interactionType_)
582 bounceInteraction(p, pp, facei, keepParticle);
588 const scalar m = p.nParticle()*p.mass();
589 absorbInteraction(filmModel, p, pp, facei, m, keepParticle);
595 bool dry = this->deltaFilmPatch_[
patchi][facei] < deltaWet_;
599 drySplashInteraction(filmModel, p, pp, facei, keepParticle);
603 wetSplashInteraction(filmModel, p, pp, facei, keepParticle);
611 <<
"Unknown interaction type enumeration" 625 template<
class CloudType>
628 const label filmPatchi,
629 const label primaryPatchi,
641 filmModel.
toPrimary(filmPatchi, TFilmPatch_);
644 filmModel.
toPrimary(filmPatchi, CpFilmPatch_);
648 template<
class CloudType>
652 const label filmFacei
658 p.T() = TFilmPatch_[filmFacei];
659 p.Cp() = CpFilmPatch_[filmFacei];
663 template<
class CloudType>
668 label nSplash0 = this->
template getModelProperty<label>(
"nParcelsSplashed");
672 os <<
" New film splash parcels = " << nSplashTotal <<
endl;
674 if (this->writeTime())
676 this->setModelProperty(
"nParcelsSplashed", nSplashTotal);
677 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 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 > &)
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/m2/K4].
const Field< PointType > & faceNormals() const
Return face normals for patch.
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.
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.
stressControl lookup("compactNormalStress") >> compactNormalStress
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)
void absorbInteraction(regionModels::surfaceFilmModels::surfaceFilmModel &filmModel, const parcelType &p, const polyPatch &pp, const label facei, const scalar mass, bool &keepParticle)
Absorb parcel into film.
The thermophysical properties of a liquid.
virtual const volScalarField & Ts() const =0
Return the film surface temperature [K].
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 & 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.
Templated wall surface film model class.
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
cachedRandom & rndGen_
Reference to the cloud random number generator.
void bounceInteraction(parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle) const
Bounce parcel (flip parcel normal velocity)
A normal distribution model.
const dimensionedScalar mu
Atomic mass unit.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
vector tangentVector(const vector &v) const
Return a vector tangential to input vector, v.
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 index() const
Return the index of this patch in the boundaryMesh.
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.
const volVectorField & C() const
Return cell centres as volVectorField.
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.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
scalar deltaWet_
Film thickness beyond which patch is assumed to be wet.
Templated base class for dsmc cloud.
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.
label whichFace(const label l) const
Return label of face in patch from global face label.