Go to the documentation of this file.
42 Foam::humidityTemperatureCoupledMixedFvPatchScalarField::massModeTypeNames_
44 { massTransferMode::mtConstantMass,
"constantMass" },
45 { massTransferMode::mtCondensation,
"condensation" },
46 { massTransferMode::mtEvaporation,
"evaporation" },
48 massTransferMode::mtCondensationAndEvaporation,
49 "condensationAndEvaporation"
56 Foam::scalar Foam::humidityTemperatureCoupledMixedFvPatchScalarField::Sh
74 Foam::humidityTemperatureCoupledMixedFvPatchScalarField::htcCondensation
80 if (Tsat > 295 && Tsat < 373)
82 return 51104 + 2044*Tsat;
92 Foam::humidityTemperatureCoupledMixedFvPatchScalarField::thicknessField
94 const word& fieldName,
132 mixedFvPatchScalarField(
p, iF),
141 mode_(mtConstantMass),
145 muName_(
"thermo:mu"),
151 liquidDict_(
nullptr),
164 this->refValue() = 0.0;
165 this->refGrad() = 0.0;
166 this->valueFraction() = 1.0;
179 mixedFvPatchScalarField(psf,
p, iF, mapper),
184 rhoName_(psf.rhoName_),
185 muName_(psf.muName_),
186 TnbrName_(psf.TnbrName_),
187 qrNbrName_(psf.qrNbrName_),
188 qrName_(psf.qrName_),
189 specieName_(psf.specieName_),
190 liquid_(psf.liquid_.clone()),
191 liquidDict_(psf.liquidDict_),
192 mass_(psf.mass_, mapper),
194 myKDelta_(psf.myKDelta_, mapper),
195 dmHfg_(psf.dmHfg_, mapper),
196 mpCpTp_(psf.mpCpTp_, mapper),
200 cp_(psf.cp_, mapper),
201 thickness_(psf.thickness_, mapper),
202 rho_(psf.rho_, mapper)
214 mixedFvPatchScalarField(
p, iF),
216 mode_(mtCondensationAndEvaporation),
239 if (!isA<mappedPatchBase>(this->
patch().
patch()))
242 <<
"\n patch type '" <<
p.type()
243 <<
"' not type '" << mappedPatchBase::typeName <<
"'"
244 <<
"\n for patch " <<
p.name()
245 <<
" of field " << internalField().name()
246 <<
" in file " << internalField().objectPath()
252 if (massModeTypeNames_.readIfPresent(
"mode",
dict, mode_))
271 case mtCondensationAndEvaporation:
291 thickness_[i]*liquid_->rho(pf, Tp[i])*magSf[i];
301 <<
"Did not find mode " << mode_
302 <<
" on patch " <<
patch().name()
304 <<
"Please set 'mode' to one of "
305 << massModeTypeNames_.sortedToc()
325 valueFraction() = 1.0;
337 mixedFvPatchScalarField(psf, iF),
342 rhoName_(psf.rhoName_),
343 muName_(psf.muName_),
344 TnbrName_(psf.TnbrName_),
345 qrNbrName_(psf.qrNbrName_),
346 qrName_(psf.qrName_),
347 specieName_(psf.specieName_),
348 liquid_(psf.liquid_.clone()),
349 liquidDict_(psf.liquidDict_),
352 myKDelta_(psf.myKDelta_),
354 mpCpTp_(psf.mpCpTp_),
359 thickness_(psf.thickness_),
371 mixedFvPatchScalarField::autoMap(m);
376 myKDelta_.autoMap(m);
380 thickness_.autoMap(m);
392 mixedFvPatchScalarField::rmap(ptf, addr);
395 refCast<const humidityTemperatureCoupledMixedFvPatchScalarField>
402 mass_.rmap(tiptf.mass_, addr);
403 myKDelta_.rmap(tiptf.myKDelta_, addr);
404 dmHfg_.rmap(tiptf.dmHfg_, addr);
405 mpCpTp_.rmap(tiptf.mpCpTp_, addr);
406 cp_.rmap(tiptf.cp_, addr);
407 thickness_.rmap(tiptf.thickness_, addr);
408 rho_.rmap(tiptf.rho_, addr);
422 refCast<const mappedPatchBase>(
patch().
patch());
430 refCast<const fvMesh>(nbrMesh).boundary()[nbrPatchI];
443 scalarField nbrIntFld(nbrField.patchInternalField());
469 myKDelta_ =
K*
patch().deltaCoeffs();
481 if (mode_ != mtConstantMass)
520 const scalar Tf = Tp[faceI];
521 const scalar Tint = Tin[faceI];
523 const scalar pf = pp[faceI];
525 const scalar muf = mup[faceI];
526 const scalar
rhof = rhop[faceI];
527 const scalar nuf = muf/
rhof;
528 const scalar pSat = liquid_->pv(pf, Tint);
529 const scalar Mv = liquid_->W();
530 const scalar TSat = liquid_->pvInvert(pSat);
531 const scalar
Re =
mag(
Uf)*L_/nuf;
533 cp[faceI] = liquid_->Cp(pf, Tf);
534 hfg[faceI] = liquid_->hl(pf, Tf);
537 const scalar invMwmean =
538 Yi[faceI]/Mv + (1.0 - Yi[faceI])/Mcomp_;
539 const scalar Xv = Yi[faceI]/invMwmean/Mv;
540 const scalar RH =
min(Xv*pf/pSat, 1.0);
543 scalar Tdew = -GREAT;
549 scalar TintDeg = Tint - 273;
551 b*(
log(RH) + (
c*TintDeg)/(
b + TintDeg))
552 /(
c -
log(RH) - ((
c*TintDeg)/(
b + TintDeg))) + 273;
560 mode_ == mtCondensation
561 || mode_ == mtCondensationAndEvaporation
565 htc[faceI] = htcCondensation(TSat,
Re)*nbrK[faceI]/L_;
568 1.0/((1.0/myKDelta_[faceI]) + (1.0/htc[faceI]));
571 const scalar q = (Tint - Tf)*htcTotal*magSf[faceI];
574 dm[faceI] = q/hfg[faceI]/magSf[faceI];
576 mass_[faceI] += q/hfg[faceI]*dt;
579 const scalar Dab = liquid_->D(pf, Tf);
581 -
min(dm[faceI]/Dab/
rhof, Yi[faceI]*myDelta[faceI]);
586 && mass_[faceI] > 0.0
588 mode_ == mtEvaporation
589 || mode_ == mtCondensationAndEvaporation
593 const scalar Dab = liquid_->D(pf, Tf);
595 const scalar Sc = nuf/Dab;
596 const scalar Sh = this->Sh(
Re, Sc);
598 const scalar Ys = Mv*pSat/(Mv*pSat + Mcomp_*(pf - pSat));
601 const scalar hm = Dab*Sh/L_;
603 const scalar Yinf =
max(Yi[faceI], 0.0);
606 dm[faceI] = -
rhof*hm*
max((Ys - Yinf), 0.0)/(1.0 - Ys);
609 Yvp[faceI] = -dm[faceI]/Dab/
rhof;
612 mass_[faceI] += dm[faceI]*magSf[faceI]*dt;
614 htc[faceI] = htcCondensation(TSat,
Re)*nbrK[faceI]/L_;
616 else if (Tf > Tdew && Tf < Tvap_ && mass_[faceI] > 0.0)
618 htc[faceI] = htcCondensation(TSat,
Re)*nbrK[faceI]/L_;
620 else if (mass_[faceI] == 0.0)
625 liquidRho[faceI] = liquid_->rho(pf, Tf);
628 mass_ =
max(mass_, scalar(0));
633 const word fieldName(specieName_ +
"Thickness");
639 refCast<const fvMesh>(
mesh)
640 ).boundaryFieldRef()[
patch().index()];
643 pDelta = mass_/liquidRho/magSf;
646 myKDelta_ = 1.0/((1.0/myKDelta_) + (1.0/htc));
648 mpCpTp_ = mass_*
cp/dt/magSf;
656 mpCpTp_ = thickness_*rho_*cp_/dt;
665 mpCpTpNbr = nbrField.
mpCpTp();
668 dmHfgNbr = nbrField.
dmHfg();
674 if (qrName_ !=
"none")
680 if (qrNbrName_ !=
"none")
695 refValue() = (KDeltaNbr*nbrIntFld + mpCpdt*TpOld + dmHfg)/
alpha;
697 mixedFvPatchScalarField::updateCoeffs();
701 scalar Qdm =
gSum(dm);
702 scalar QMass =
gSum(mass_);
703 scalar Qt =
gSum(myKDelta_*(Tp - Tin)*magSf);
704 scalar QtSolid =
gSum(KDeltaNbr*(Tp - nbrIntFld)*magSf);
707 <<
patch().name() <<
':'
708 << internalField().name() <<
" <- "
709 << nbrMesh.
name() <<
':'
710 << nbrPatch.
name() <<
':'
711 << internalField().name() <<
" :" <<
nl
712 <<
" Total mass flux [Kg/s] : " << Qdm <<
nl
713 <<
" Total mass on the wall [Kg] : " << QMass <<
nl
714 <<
" Total heat (>0 leaving the wall to the fluid) [W] : "
716 <<
" Total heat (>0 leaving the wall to the solid) [W] : "
718 <<
" wall temperature "
719 <<
" min:" <<
gMin(*
this)
720 <<
" max:" <<
gMax(*
this)
743 os.
writeEntry(
"mode", massModeTypeNames_[mode_]);
754 if (mode_ == mtConstantMass)
763 liquidDict_.write(os);
int debug
Static debugging option.
Type * getObjectPtr(const word &name, const bool recursive=false) const
const scalarField mpCpTp() const
Return mpCpTp.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
const word & name() const
Return name.
A class for handling words, derived from Foam::string.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
static constexpr const zero Zero
Global zero (0)
static autoPtr< liquidProperties > New(const word &name)
Return a pointer to a new liquidProperties created from name.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Type gAverage(const FieldField< Field, Type > &f)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const polyPatch & samplePolyPatch() const
Get the patch on the region.
static word timeName(const scalar t, const int precision=precision_)
autoPtr< surfaceVectorField > Uf
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Common functions used in temperature coupled boundaries.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
virtual void write(Ostream &) const
Write.
virtual const word & name() const
Return name.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
virtual Field< Type > & gradient()
Return gradient at boundary.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
CGAL::Exact_predicates_exact_constructions_kernel K
const polyMesh & sampleMesh() const
Get the region mesh.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
const scalarField dmHfg() const
Return dmHfg.
messageStream Info
Information stream (uses stdout - output is on the master only)
humidityTemperatureCoupledMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dimensionedScalar log(const dimensionedScalar &ds)
Macros for easy insertion into run-time selection tables.
To & refCast(From &r)
Reference type cast template function.
errorManipArg< error, int > exit(error &err, const int errNo=1)
GeometricField< vector, fvPatchField, volMesh > volVectorField
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy the source to the destination (recursively if necessary).
const std::string patch
OpenFOAM patch number as a std::string.
scalar deltaTValue() const
Return time step value.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const scalarField & deltaCoeffs() const
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
scalarField Re(const UList< complex > &cf)
Extract real component.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
const dimensionedScalar c
Speed of light in a vacuum.
const Time & time() const
Return the top-level database.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Foam::fvPatchFieldMapper.
void write(Ostream &os) const
Write.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Type gMin(const FieldField< Field, Type > &f)
massTransferMode
Modes of mass transfer.
dimensionedScalar cbrt(const dimensionedScalar &ds)
This boundary condition supplies a fixed gradient condition, such that the patch values are calculate...
surfaceScalarField rhof(fvc::interpolate(rho, "div(phi,rho)"))
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
virtual void operator=(const UList< scalar > &)
label index() const
The index of this patch in the boundaryMesh.
Type gMax(const FieldField< Field, Type > &f)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const word & name() const
Return reference to name.