38 template<
class CloudType>
47 const scalar phiSi =
twoPi*rndGen_.sample01<scalar>();
50 const scalar thetaSi =
pi/180.0*(rndGen_.sample01<scalar>()*(50 - 5) + 5);
54 const scalar dcorr =
cos(thetaSi);
59 return dirVec/
mag(dirVec);
63 template<
class CloudType>
76 Info<<
"Parcel " <<
p.origId() <<
" absorbInteraction" <<
endl;
90 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
96 const vector Un = nf*(Urel & nf);
99 const vector Ut = Urel - Un;
104 const scalar pc = carrierThermo.
p()[
p.cell()];
112 mass*liq.
Hs(pc,
p.T())
115 this->nParcelsTransferred()++;
117 keepParticle =
false;
121 template<
class CloudType>
132 Info<<
"Parcel " <<
p.origId() <<
" bounceInteraction" <<
endl;
139 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
142 const vector Urel =
p.U() - Up;
145 p.U() -= 2.0*nf*(Urel & nf);
151 template<
class CloudType>
163 Info<<
"Parcel " <<
p.origId() <<
" drySplashInteraction" <<
endl;
176 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
180 const scalar pc = carrierThermo.
p()[
p.cell()];
183 const scalar m =
p.mass()*
p.nParticle();
184 const scalar
rho =
p.rho();
185 const scalar d =
p.d();
187 const scalar
mu = liq.
mu(pc,
p.T());
188 const vector Urel =
p.U() - Up;
189 const vector Un = nf*(Urel & nf);
198 const scalar Wec = Adry_*
pow(La, -0.183);
202 absorbInteraction(filmCloudTransfer,
p, pp, facei, m, keepParticle);
207 const scalar mRatio = 0.2 + 0.6*rndGen_.sample01<scalar>();
224 template<
class CloudType>
236 Info<<
"Parcel " <<
p.origId() <<
" wetSplashInteraction" <<
endl;
249 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
253 const scalar pc = carrierThermo.
p()[
p.cell()];
256 const scalar m =
p.mass()*
p.nParticle();
257 const scalar
rho =
p.rho();
258 const scalar d =
p.d();
261 const scalar
mu = liq.
mu(pc,
p.T());
262 const vector Urel =
p.U() - Up;
263 const vector Un = nf*(Urel & nf);
264 const vector Ut = Urel - Un;
273 const scalar Wec = Awet_*
pow(La, -0.183);
277 absorbInteraction(filmCloudTransfer,
p, pp, facei, m, keepParticle);
279 else if ((We >= 2) && (We < 20))
285 const scalar
epsilon = 0.993 - theta*(1.76 - theta*(1.56 - theta*0.49));
293 else if ((We >= 20) && (We < Wec))
295 absorbInteraction(filmCloudTransfer,
p, pp, facei, m, keepParticle);
301 const scalar mRatio = 0.2 + 0.9*rndGen_.sample01<scalar>();
318 template<
class CloudType>
333 const fvMesh& mesh = this->owner().mesh();
334 const vector& Up = this->owner().U().boundaryField()[pp.
index()][facei];
339 const vector tanVec2 = nf^tanVec1;
342 const scalar np =
p.nParticle();
343 const scalar m =
p.mass()*np;
344 const scalar d =
p.d();
345 const vector Urel =
p.U() - Up;
346 const vector Un = nf*(Urel & nf);
347 const vector Ut = Urel - Un;
348 const vector& posC = mesh.
C()[
p.cell()];
352 const scalar mSplash = m*mRatio;
355 const scalar Ns = 5.0*(We/Wec - 1.0);
358 const scalar dBarSplash = 1/
cbrt(6.0)*
cbrt(mRatio/Ns)*d + rootVSmall;
361 const scalar dMax = 0.9*
cbrt(mRatio)*d;
362 const scalar dMin = 0.1*dMax;
363 const scalar
K =
exp(-dMin/dBarSplash) -
exp(-dMax/dBarSplash);
366 scalar ESigmaSec = 0;
373 const scalar
y = rndGen_.sample01<scalar>();
374 dNew[i] = -dBarSplash*
log(
exp(-dMin/dBarSplash) -
y*
K);
375 npNew[i] = mRatio*np*
pow3(d)/
pow3(dNew[i])/parcelsPerSplash_;
376 ESigmaSec += npNew[i]*
sigma*
p.areaS(dNew[i]);
380 const scalar EKIn = 0.5*m*
magSqr(Un);
383 const scalar ESigmaIn = np*
sigma*
p.areaS(d);
389 const scalar EKs = EKIn + ESigmaIn - ESigmaSec - Ed;
394 absorbInteraction(filmCloudTransfer,
p, pp, facei, m, keepParticle);
399 const scalar logD =
log(d);
400 const scalar coeff2 =
log(dNew[0]) - logD + rootVSmall;
404 coeff1 +=
sqr(
log(dNew[i]) - logD);
408 const scalar magUns0 =
409 sqrt(2.0*parcelsPerSplash_*EKs/mSplash/(1.0 + coeff1/
sqr(coeff2)));
414 const vector dirVec = splashDirection(tanVec1, tanVec2, -nf);
419 pPtr->origId() = pPtr->getNewParticleID();
423 if (splashParcelType_ >= 0)
425 pPtr->typeId() = splashParcelType_;
429 pPtr->track(mesh, 0.5*rndGen_.sample01<scalar>()*(posC - posCf), 0);
431 pPtr->nParticle() = npNew[i];
435 pPtr->U() = dirVec*(
mag(Cf_*Ut) + magUns0*(
log(dNew[i]) - logD)/coeff2);
441 this->owner().addParticle(pPtr);
448 const scalar mDash = m - mSplash;
449 absorbInteraction(filmCloudTransfer,
p, pp, facei, mDash, keepParticle);
453 template<
class CloudType>
457 if (filmTransfers_.empty())
474 refCast<const fvMesh>(mpb.
nbrMesh())
482 if (isType<fv::filmCloudTransfer>(
fvModels[i]))
485 &refCast<fv::filmCloudTransfer>(
fvModels[i]);
491 Info<<
"Found filmCloudTransfer fvModel "
495 filmTransfers_.resize(nFilms + 1);
496 filmPatches_.resize(nFilms + 1);
498 filmTransfers_.set(nFilms, filmCloudPtr);
499 filmPatches_[nFilms] =
patchi;
507 return filmTransfers_;
511 template<
class CloudType>
521 template<
class CloudType>
528 this->deltaFilmPatch_ = filmCloudTransfer.
deltaToCloud();
534 UFilmPatch_ = filmCloudTransfer.
UToCloud();
535 rhoFilmPatch_ = filmCloudTransfer.
rhoToCloud();
536 TFilmPatch_ = filmCloudTransfer.
TToCloud();
537 CpFilmPatch_ = filmCloudTransfer.
CpToCloud();
542 this->massParcelPatch_.setSize(0);
547 template<
class CloudType>
551 const label filmFacei
557 p.d() = this->diameterParcelPatch_[filmFacei];
558 p.U() = UFilmPatch_[filmFacei];
559 p.rho() = rhoFilmPatch_[filmFacei];
561 p.nParticle() = this->massParcelPatch_[filmFacei]/
p.rho()/vol;
563 if (this->ejectedParcelType_ >= 0)
565 p.typeId() = this->ejectedParcelType_;
569 p.T() = TFilmPatch_[filmFacei];
570 p.Cp() = CpFilmPatch_[filmFacei];
576 template<
class CloudType>
591 interactionTypeNames_.
read(this->coeffDict().lookup(
"interactionType"))
594 splashParcelType_(0),
595 parcelsPerSplash_(0),
602 <<
" interaction model" <<
endl;
618 template<
class CloudType>
625 rndGen_(sfm.rndGen_),
626 UFilmPatch_(sfm.UFilmPatch_),
627 rhoFilmPatch_(sfm.rhoFilmPatch_),
628 TFilmPatch_(sfm.TFilmPatch_),
629 CpFilmPatch_(sfm.CpFilmPatch_),
630 interactionType_(sfm.interactionType_),
631 deltaWet_(sfm.deltaWet_),
632 splashParcelType_(sfm.splashParcelType_),
633 parcelsPerSplash_(sfm.parcelsPerSplash_),
637 nParcelsSplashed_(sfm.nParcelsSplashed_)
643 template<
class CloudType>
650 template<
class CloudType>
660 forAll(this->filmTransferPtrs(), filmi)
662 if (
patchi == filmPatches_[filmi])
665 filmTransferPtrs()[filmi];
669 switch (interactionType_)
671 case interactionType::bounce:
673 bounceInteraction(
p, pp, facei, keepParticle);
676 case interactionType::absorb:
684 p.nParticle()*
p.mass(),
689 case interactionType::splashBai:
691 if (this->deltaFilmPatch_[facei] < deltaWet_)
718 <<
"Unknown interaction type enumeration"
733 template<
class CloudType>
738 label nSplash0 = this->
template getModelProperty<label>(
"nParcelsSplashed");
742 os <<
" New film splash parcels = " << nSplashTotal <<
endl;
744 if (this->writeTime())
746 this->setModelProperty(
"nParcelsSplashed", nSplashTotal);
747 nParcelsSplashed_ = 0;
#define forAll(list, i)
Loop across all elements in list.
Macros for easy insertion into run-time selection tables.
static const NamedEnum< interactionType, 3 > interactionTypeNames_
Interaction type names.
Thermo parcel<->film transfer model.
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
label splashParcelType_
Splash parcel type label - id assigned to identify parcel for.
virtual void cacheFilmFields(const label filmi)
Cache the film fields in preparation for injection.
virtual const labelList & filmPatches() const
Return pointers to the films.
scalar Cf_
Skin friction typically in the range 0.6 < Cf < 0.8.
void drySplashInteraction(fv::filmCloudTransfer &, const parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with dry surface.
scalar Awet_
Wet surface roughness coefficient.
virtual void info(Ostream &os)
Write film info to stream.
void splashInteraction(fv::filmCloudTransfer &, 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.
vector splashDirection(const vector &tanVec1, const vector &tanVec2, const vector &nf) const
Return splashed parcel direction.
void bounceInteraction(parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle) const
Bounce parcel (flip parcel normal velocity)
void absorbInteraction(fv::filmCloudTransfer &, const parcelType &p, const polyPatch &pp, const label facei, const scalar mass, bool &keepParticle)
Absorb parcel into film.
scalar Adry_
Dry surface roughness coefficient.
scalar deltaWet_
Film thickness beyond which patch is assumed to be wet.
interactionType interactionType_
Interaction type enumeration.
virtual ~CloudFilmTransfer()
Destructor.
CloudFilmTransfer(const dictionary &dict, CloudType &owner)
Construct from components.
void wetSplashInteraction(fv::filmCloudTransfer &, parcelType &p, const polyPatch &pp, const label facei, bool &keepParticle)
Parcel interaction with wetted surface.
label parcelsPerSplash_
Number of new parcels resulting from splash event.
virtual bool transferParcel(parcelType &p, const polyPatch &pp, bool &keepParticle)
Transfer parcel from cloud to film.
Templated base class for dsmc cloud.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const
Return name.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const Field< PointType > & faceNormals() const
Return face normals for patch.
Templated wall surface film model class.
virtual void info(Ostream &os)
Write surface film info to stream.
CloudType::parcelType parcelType
Convenience typedef to the cloud's parcel type.
Templated base class for thermodynamic cloud.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
A list of keyword definitions, which are a keyword followed by any number of values (e....
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Base-class for fluid thermodynamic properties.
virtual volScalarField & p()=0
Pressure [Pa].
Mesh data needed to do the Finite Volume discretisation.
const volVectorField & C() const
Return cell centres.
const surfaceVectorField & Cf() const
Return face centres.
Film<->cloud transfer model.
tmp< Field< scalar > > deltaToCloud() const
Transfer the film delta field to the cloud.
void resetFromCloudFields()
Reset the fields accumulated cloud transfer fields.
tmp< Field< scalar > > TToCloud() const
Transfer the film temperature field to the cloud.
tmp< Field< scalar > > CpToCloud() const
Transfer the film heat capacity field to the cloud.
tmp< Field< scalar > > ejectedDiameterToCloud() const
Transfer the ejected droplet diameter to the cloud.
tmp< Field< vector > > UToCloud() const
Transfer the film velocity field to the cloud.
tmp< Field< scalar > > ejectedMassToCloud() const
Transfer the ejected mass to the cloud.
void parcelFromCloud(const label facei, const scalar mass, const vector &momentum, const scalar pressure, const scalar energy)
Transfer parcel properties from cloud to the film.
bool ejecting() const
Return true if the film is ejecting to the cloud.
tmp< Field< scalar > > rhoToCloud() const
Transfer the film density field to the cloud.
The thermophysical properties of a liquid.
virtual scalar mu(scalar p, scalar T) const =0
Liquid viscosity [Pa s].
virtual scalar sigma(scalar p, scalar T) const =0
Surface tension [N/m].
virtual scalar Hs(scalar p, scalar T) const =0
Liquid sensible enthalpy [J/kg].
Engine which provides mapping between two patches.
const polyMesh & nbrMesh() const
Get the mesh for the region to map from.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
label index() const
Return the index of this patch in the boundaryMesh.
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
A patch is a list of labels that address the faces in the global face list.
label whichFace(const label l) const
Return label of face in patch from global face label.
const dictionary & coeffDict() const
Return const access to the coefficients dictionary.
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const fvPatchList & patches
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const scalar twoPi(2 *pi)
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
const dimensionedScalar mu
Atomic mass unit.
dimensionedScalar exp(const dimensionedScalar &ds)
bool read(const char *, int32_t &)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pow3(const dimensionedScalar &ds)
errorManip< error > abort(error &err)
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Vector< Cmpt > perpendicular(const Vector< Cmpt > &v)
dimensionedScalar sqrt(const dimensionedScalar &ds)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Form normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
dimensionedScalar cbrt(const dimensionedScalar &ds)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
fluidMulticomponentThermo & thermo