42Foam::humidityTemperatureCoupledMixedFvPatchScalarField::massModeTypeNames_
44 { massTransferMode::mtConstantMass,
"constantMass" },
45 { massTransferMode::mtCondensation,
"condensation" },
46 { massTransferMode::mtEvaporation,
"evaporation" },
48 massTransferMode::mtCondensationAndEvaporation,
49 "condensationAndEvaporation"
74Foam::humidityTemperatureCoupledMixedFvPatchScalarField::htcCondensation
81 if (Tsat > 295.15 && Tsat < 373.15)
83 return 51104 + 2044*(Tsat - 273.15);
93Foam::humidityTemperatureCoupledMixedFvPatchScalarField::thicknessField
95 const word& fieldName,
127Foam::humidityTemperatureCoupledMixedFvPatchScalarField::
128humidityTemperatureCoupledMixedFvPatchScalarField
134 mixedFvPatchScalarField(
p, iF),
143 mode_(mtConstantMass),
147 muName_(
"thermo:mu"),
153 liquidDict_(nullptr),
154 mass_(patch().size(),
Zero),
156 myKDelta_(patch().size(),
Zero),
157 dmHfg_(patch().size(),
Zero),
158 mpCpTp_(patch().size(),
Zero),
162 cp_(patch().size(),
Zero),
163 thickness_(patch().size(),
Zero),
164 rho_(patch().size(),
Zero),
168 this->refValue() = 0.0;
169 this->refGrad() = 0.0;
170 this->valueFraction() = 1.0;
174Foam::humidityTemperatureCoupledMixedFvPatchScalarField::
175humidityTemperatureCoupledMixedFvPatchScalarField
183 mixedFvPatchScalarField(psf,
p, iF, mapper),
188 rhoName_(psf.rhoName_),
189 muName_(psf.muName_),
190 TnbrName_(psf.TnbrName_),
191 qrNbrName_(psf.qrNbrName_),
192 qrName_(psf.qrName_),
193 specieName_(psf.specieName_),
194 liquid_(psf.liquid_.clone()),
195 liquidDict_(psf.liquidDict_),
196 mass_(psf.mass_, mapper),
198 myKDelta_(psf.myKDelta_, mapper),
199 dmHfg_(psf.dmHfg_, mapper),
200 mpCpTp_(psf.mpCpTp_, mapper),
204 cp_(psf.cp_, mapper),
205 thickness_(psf.thickness_, mapper),
206 rho_(psf.rho_, mapper),
207 thicknessLayers_(psf.thicknessLayers_),
208 kappaLayers_(psf.kappaLayers_)
212Foam::humidityTemperatureCoupledMixedFvPatchScalarField::
213humidityTemperatureCoupledMixedFvPatchScalarField
220 mixedFvPatchScalarField(
p, iF),
222 mode_(mtCondensationAndEvaporation),
223 pName_(
dict.getOrDefault<
word>(
"p",
"p")),
224 UName_(
dict.getOrDefault<
word>(
"U",
"U")),
225 rhoName_(
dict.getOrDefault<
word>(
"rho",
"rho")),
226 muName_(
dict.getOrDefault<
word>(
"mu",
"thermo:mu")),
227 TnbrName_(
dict.getOrDefault<
word>(
"Tnbr",
"T")),
228 qrNbrName_(
dict.getOrDefault<
word>(
"qrNbr",
"none")),
229 qrName_(
dict.getOrDefault<
word>(
"qr",
"none")),
230 specieName_(
dict.getOrDefault<
word>(
"specie",
"none")),
233 mass_(patch().size(),
Zero),
235 myKDelta_(patch().size(),
Zero),
236 dmHfg_(patch().size(),
Zero),
237 mpCpTp_(patch().size(),
Zero),
241 cp_(patch().size(),
Zero),
242 thickness_(patch().size(),
Zero),
243 rho_(patch().size(),
Zero),
247 if (!isA<mappedPatchBase>(this->patch().patch()))
250 <<
"\n patch type '" <<
p.type()
252 <<
"\n for patch " <<
p.
name()
253 <<
" of field " << internalField().name()
254 <<
" in file " << internalField().objectPath()
304 thickness_[i]*liquid_->rho(pf, Tp[i])*magSf[i];
314 <<
"Did not find mode " << mode_
315 <<
" on patch " << patch().name()
317 <<
"Please set 'mode' to one of "
338 valueFraction() = 1.0;
343Foam::humidityTemperatureCoupledMixedFvPatchScalarField::
344humidityTemperatureCoupledMixedFvPatchScalarField
350 mixedFvPatchScalarField(psf, iF),
355 rhoName_(psf.rhoName_),
356 muName_(psf.muName_),
357 TnbrName_(psf.TnbrName_),
358 qrNbrName_(psf.qrNbrName_),
359 qrName_(psf.qrName_),
360 specieName_(psf.specieName_),
361 liquid_(psf.liquid_.clone()),
362 liquidDict_(psf.liquidDict_),
365 myKDelta_(psf.myKDelta_),
367 mpCpTp_(psf.mpCpTp_),
372 thickness_(psf.thickness_),
374 thicknessLayers_(psf.thicknessLayers_),
375 kappaLayers_(psf.kappaLayers_)
386 mixedFvPatchScalarField::autoMap(m);
392 myKDelta_.autoMap(m);
396 thickness_.autoMap(m);
408 mixedFvPatchScalarField::rmap(ptf, addr);
411 refCast<const humidityTemperatureCoupledMixedFvPatchScalarField>
419 mass_.rmap(tiptf.mass_, addr);
420 myKDelta_.rmap(tiptf.myKDelta_, addr);
421 dmHfg_.rmap(tiptf.dmHfg_, addr);
422 mpCpTp_.rmap(tiptf.mpCpTp_, addr);
423 cp_.rmap(tiptf.cp_, addr);
424 thickness_.rmap(tiptf.thickness_, addr);
425 rho_.rmap(tiptf.rho_, addr);
439 refCast<const mappedPatchBase>(patch().patch());
447 refCast<const fvMesh>(nbrMesh).boundary()[nbrPatchI];
460 scalarField nbrIntFld(nbrField.patchInternalField());
486 myKDelta_ =
K*patch().deltaCoeffs();
488 if (thicknessLayers_.size() > 0)
490 myKDelta_ = 1.0/myKDelta_;
491 forAll(thicknessLayers_, iLayer)
493 myKDelta_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
495 myKDelta_ = 1.0/myKDelta_;
508 if (mode_ != mtConstantMass)
549 const scalar Tf = Tp[faceI];
550 const scalar Tint = Tin[faceI];
552 const scalar pf = pp[faceI];
554 const scalar muf = mup[faceI];
555 const scalar
rhof = rhop[faceI];
556 const scalar nuf = muf/
rhof;
557 const scalar pSat = liquid_->pv(pf, Tint);
558 const scalar Mv = liquid_->W();
559 const scalar TSat = liquid_->pvInvert(pSat);
560 const scalar
Re =
mag(
Uf)*L_/nuf;
562 cp[faceI] = liquid_->Cp(pf, Tf);
563 hfg[faceI] = liquid_->hl(pf, Tf);
566 const scalar invMwmean =
567 Yi[faceI]/Mv + (1.0 - Yi[faceI])/Mcomp_;
568 const scalar Xv = Yi[faceI]/invMwmean/Mv;
569 RH[faceI] =
min(Xv*pf/pSat, 1.0);
574 if (RH[faceI] > RHmin)
578 scalar TintDeg = Tint - 273;
580 b*(
log(RH[faceI]) + (c*TintDeg)/(
b + TintDeg))
581 /(c -
log(RH[faceI]) - ((c*TintDeg)/(
b + TintDeg)))
590 mode_ == mtCondensation
591 || mode_ == mtCondensationAndEvaporation
595 htc[faceI] = htcCondensation(TSat,
Re)*nbrK[faceI]/L_;
598 1.0/((1.0/myKDelta_[faceI]) + (1.0/htc[faceI]));
601 const scalar q = (Tint - Tf)*htcTotal*magSf[faceI];
604 dm[faceI] = q/hfg[faceI]/magSf[faceI];
606 mass_[faceI] += q/hfg[faceI]*dt;
609 const scalar Dab = liquid_->D(pf, Tf);
611 -
min(dm[faceI]/Dab/
rhof, Yi[faceI]*myDelta[faceI]);
616 && mass_[faceI] > 0.0
618 mode_ == mtEvaporation
619 || mode_ == mtCondensationAndEvaporation
623 const scalar Dab = liquid_->D(pf, Tf);
625 const scalar Sc = nuf/Dab;
626 const scalar Sh = this->Sh(
Re, Sc);
628 const scalar Ys = Mv*pSat/(Mv*pSat + Mcomp_*(pf - pSat));
631 const scalar hm = Dab*Sh/L_;
633 const scalar Yinf =
max(Yi[faceI], 0.0);
636 dm[faceI] = -
rhof*hm*
max((Ys - Yinf), 0.0)/(1.0 - Ys);
639 Yvp[faceI] = -dm[faceI]/Dab/
rhof;
642 mass_[faceI] += dm[faceI]*magSf[faceI]*dt;
644 htc[faceI] = htcCondensation(TSat,
Re)*nbrK[faceI]/L_;
646 else if (Tf > Tdew[faceI] && Tf < Tvap_ && mass_[faceI] > 0.0)
648 htc[faceI] = htcCondensation(TSat,
Re)*nbrK[faceI]/L_;
650 else if (mass_[faceI] == 0.0)
655 liquidRho[faceI] = liquid_->rho(pf, Tf);
658 mass_ =
max(mass_, scalar(0));
663 const word fieldName(specieName_ +
"Thickness");
669 refCast<const fvMesh>(
mesh)
670 ).boundaryFieldRef()[patch().index()];
673 pDelta = mass_/liquidRho/magSf;
676 myKDelta_ = 1.0/((1.0/myKDelta_) + (1.0/htc));
678 mpCpTp_ = mass_*
cp/dt/magSf;
690 refCast<const fvMesh>(
mesh)
691 ).boundaryFieldRef()[patch().index()];
698 refCast<const fvMesh>(
mesh)
699 ).boundaryFieldRef()[patch().index()];
706 mpCpTp_ = thickness_*rho_*cp_/dt;
715 mpCpTpNbr = nbrField.
mpCpTp();
718 dmHfgNbr = nbrField.
dmHfg();
724 if (qrName_ !=
"none")
730 if (qrNbrName_ !=
"none")
745 refValue() = (KDeltaNbr*nbrIntFld + mpCpdt*TpOld + dmHfg)/
alpha;
747 mixedFvPatchScalarField::updateCoeffs();
751 scalar Qdm =
gSum(dm*magSf);
752 scalar QMass =
gSum(mass_);
753 scalar Qt =
gSum(myKDelta_*(Tp - Tin)*magSf);
754 scalar QtSolid =
gSum(KDeltaNbr*(Tp - nbrIntFld)*magSf);
757 << patch().name() <<
':'
758 << internalField().name() <<
" <- "
759 << nbrMesh.
name() <<
':'
760 << nbrPatch.
name() <<
':'
761 << internalField().name() <<
" :" <<
nl
762 <<
" Total mass flux [Kg/s] : " << Qdm <<
nl
763 <<
" Total mass on the wall [Kg] : " << QMass <<
nl
764 <<
" Total heat (>0 leaving the wall to the fluid) [W] : "
766 <<
" Total heat (>0 leaving the wall to the solid) [W] : "
768 <<
" wall temperature "
769 <<
" min:" <<
gMin(*
this)
770 <<
" max:" <<
gMax(*
this)
782 mixedFvPatchScalarField::write(
os);
804 if (mode_ == mtConstantMass)
813 liquidDict_.write(
os);
816 if (thicknessLayers_.size())
818 thicknessLayers_.writeEntry(
"thicknessLayers",
os);
819 kappaLayers_.writeEntry(
"kappaLayers",
os);
CGAL::Exact_predicates_exact_constructions_kernel K
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
bool readIfPresent(const word &key, const dictionary &dict, EnumType &val) const
Find an entry if present, and assign to T val.
List< word > sortedToc() const
The sorted list of enum names.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const noexcept
Return the object name.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
scalar Sh() const
Sherwood number.
The TAB Method for Numerical Calculation of Spray Droplet Breakup.
scalar deltaTValue() const noexcept
Return time step value.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
void size(const label n)
Older name for setAddressableSize.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
This boundary condition supplies a fixed gradient condition, such that the patch values are calculate...
virtual Field< Type > & gradient()
Return gradient at boundary.
virtual bool write()
Write the output fields.
const Time & time() const
Return the top-level database.
const word & name() const
Return reference to name.
A FieldMapper for finite-volume patch fields.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
virtual void operator=(const UList< scalar > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual const word & name() const
Return name.
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
const scalarField & deltaCoeffs() const
const scalarField mpCpTp() const
Return mpCpTp.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
massTransferMode
Modes of mass transfer.
@ mtCondensationAndEvaporation
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
const scalarField dmHfg() const
Return dmHfg.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
const polyMesh & sampleMesh() const
Get the region mesh.
const polyPatch & samplePolyPatch() const
Get the patch on the region.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
Type * getObjectPtr(const word &name, const bool recursive=false) const
label index() const noexcept
The index of this patch in the boundaryMesh.
Mesh consisting of general polyhedral cells.
Common functions used in temperature coupled boundaries.
virtual void autoMap(const fvPatchFieldMapper &)=0
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)=0
Reverse map the given fvPatchField onto this fvPatchField.
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
A class for handling words, derived from Foam::string.
autoPtr< surfaceVectorField > Uf
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
surfaceScalarField rhof(fvc::interpolate(rho, "div(phi,rho)"))
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Type gSum(const FieldField< Field, Type > &f)
To & refCast(From &r)
Reference type cast template function.
GeometricField< vector, fvPatchField, volMesh > volVectorField
const dimensionSet dimless
Dimensionless.
scalarField Re(const UList< complex > &cf)
Extract real component.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
dimensionedScalar cbrt(const dimensionedScalar &ds)
Type gMin(const FieldField< Field, Type > &f)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Type gMax(const FieldField< Field, Type > &f)
constexpr char nl
The newline '\n' character (0x0a)
const volScalarField & cp
#define forAll(list, i)
Loop across all elements in list.
static const char *const typeName
The type name used in ensight case files.