Go to the documentation of this file.
38 namespace regionModels
40 namespace surfaceFilmModels
82 deltaMin_(coeffDict_.get<scalar>(
"deltaMin")),
83 L_(coeffDict_.get<scalar>(
"L")),
84 TbFactor_(coeffDict_.getOrDefault<scalar>(
"TbFactor", 1.1)),
85 YInfZero_(coeffDict_.getOrDefault<
Switch>(
"YInfZero",
false))
91 template<
class YInfType>
106 const label vapId =
thermo.carrierId(filmThermo.
name());
120 max(scalar(0), availableMass - deltaMin_*
rho*magSf)
124 const scalar Wvap =
thermo.carrier().W(vapId);
127 const scalar Wliq = filmThermo.
W();
133 if (
delta[celli] > deltaMin_)
136 const scalar pc = pInf[celli];
139 const scalar Tb = filmThermo.
Tb(pc);
142 const scalar Tloc =
min(TbFactor_*Tb,
max(200.0,
T[celli]));
145 const scalar pSat = filmThermo.
pv(pc, Tloc);
148 const scalar hVap = filmThermo.
hl(pc, Tloc);
154 const scalar
Cp = filmThermo.
Cp(pc, Tloc);
155 const scalar Tcorr =
max(0.0,
T[celli] - Tb);
156 const scalar qCorr = limMass[celli]*
Cp*(Tcorr);
162 const scalar rhoInfc = rhoInf[celli];
165 const scalar muInfc = muInf[celli];
168 const scalar
Re = rhoInfc*
mag(dU[celli])*L_/muInfc;
171 const scalar Ys = Wliq*pSat/(Wliq*pSat + Wvap*(pc - pSat));
174 const scalar Dab = filmThermo.
D(pc, Tloc);
177 const scalar Sc = muInfc/(rhoInfc*(Dab + ROOTVSMALL));
180 const scalar Sh = this->Sh(
Re, Sc);
183 const scalar hm = Sh*Dab/(L_ + ROOTVSMALL);
186 dm = dt*magSf[celli]*rhoInfc*hm*(Ys - YInf[celli])/(1.0 - Ys);
189 dMass[celli] +=
min(limMass[celli],
max(dm, 0));
193 dEnergy[celli] += dm*hs[celli];
210 correctModel(dt, availableMass, dMass, dEnergy,
zeroField());
218 correctModel(dt, availableMass, dMass, dEnergy, YInf);
virtual scalar Cp(const scalar p, const scalar T) const =0
Return specific heat capacity [J/kg/K].
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
virtual const volScalarField & hs() const
Return the film sensible enthalpy [J/kg].
virtual scalar hl(const scalar p, const scalar T) const =0
Return latent heat [J/kg].
const PtrList< volScalarField > & YPrimary() const
Specie mass fractions [0-1].
virtual const volVectorField & Us() const
Return the film surface velocity [m/s].
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package,...
virtual scalar W() const =0
Return molecular weight [kg/kmol].
virtual const word & name() const =0
Return the specie name.
const SLGThermo & thermo() const
Return const reference to the SLGThermo object.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
const volVectorField & UPrimary() const
Velocity [m/s].
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
standardPhaseChange(const standardPhaseChange &)=delete
No copy construct.
virtual scalar Tb(const scalar p) const =0
Return boiling temperature [K].
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
const volScalarField & rhoPrimary() const
Density [kg/m3].
Thermodynamic form of single-cell layer surface film model.
Base class for surface film models.
const filmThermoModel & filmThermo() const
Film thermo.
scalar Sh(const scalar Re, const scalar Sc) const
Return Sherwood number as a function of Reynolds and Schmidt numbers.
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void correctModel(const scalar dt, scalarField &availableMass, scalarField &dMass, scalarField &dEnergy, YInfType YInf)
const volScalarField & pPrimary() const
Pressure [Pa].
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.
virtual const volScalarField & T() const
Return the film mean temperature [K].
const volScalarField & muPrimary() const
Viscosity [Pa.s].
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
virtual scalar pv(const scalar p, const scalar T) const =0
Return vapour pressure [Pa].
virtual scalar D(const scalar p, const scalar T) const =0
Return diffusivity [m2/s].
defineTypeNameAndDebug(kinematicSingleLayer, 0)
Base class for film thermo models.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const volScalarField & delta() const
Return const access to the film thickness [m].
label carrierId(const word &cmptName, bool allowNotFound=false) const
Index of carrier component.
scalarField Re(const UList< complex > &cf)
Extract real component.
const volScalarField & Cp
virtual const volScalarField & rho() const
Return the film density [kg/m3].
dimensionedScalar cbrt(const dimensionedScalar &ds)
Base class for surface film phase change models.