Go to the documentation of this file.
42 namespace regionModels
44 namespace surfaceFilmModels
79 waxSolventEvaporation::waxSolventEvaporation
94 coeffDict_.get<scalar>(
"Wwax")
100 typeName +
":Wsolvent",
104 coeffDict_.get<scalar>(
"Wsolvent")
110 typeName +
":Ysolvent0",
121 typeName +
":Ysolvent",
129 deltaMin_(coeffDict_.get<scalar>(
"deltaMin")),
130 L_(coeffDict_.get<scalar>(
"L")),
131 TbFactor_(coeffDict_.lookupOrDefault<scalar>(
"TbFactor", 1.1)),
132 YInfZero_(coeffDict_.lookupOrDefault(
"YInfZero",
false)),
148 template<
class YInfType>
180 max(scalar(0), availableMass - deltaMin_*
rho*magSf)
184 const scalar Wvap =
thermo.carrier().W(vapId);
186 const scalar Wwax = Wwax_.value();
187 const scalar Wsolvent = Wsolvent_.value();
193 typeName +
":evapRateCoeff",
208 typeName +
":evapRateInf",
219 bool filmPresent =
false;
223 if (
delta[celli] > deltaMin_)
227 const scalar Ysolvent = Ysolvent_[celli];
230 const scalar Xsolvent
232 Ysolvent*Wsolvent/((1 - Ysolvent)*Wwax + Ysolvent*Wsolvent)
236 const scalar rhoInfc = rhoInf[celli];
239 const scalar pc = pInf[celli];
242 const scalar Tb = filmThermo.
Tb(pc);
245 const scalar Tloc =
min(TbFactor_*Tb,
max(200.0,
T[celli]));
247 const scalar pPartialCoeff
249 filmThermo.
pv(pc, Tloc)*activityCoeff_->value(Xsolvent)
252 scalar XsCoeff = pPartialCoeff/pc;
255 scalar Xs = XsCoeff*Xsolvent;
260 <<
"Solvent vapour pressure > ambient pressure"
270 XsCoeff/(XsCoeff*Xsolvent*Wsolvent + (1 - Xs)*Wvap)
274 const scalar muInfc = muInf[celli];
277 const scalar
Re = rhoInfc*
mag(dU[celli])*L_/muInfc;
280 const scalar Dab = filmThermo.
D(pc, Tloc);
283 const scalar Sc = muInfc/(rhoInfc*(Dab + ROOTVSMALL));
286 const scalar Sh = this->Sh(
Re, Sc);
289 evapRateCoeff[celli] = rhoInfc*Sh*Dab/(L_ + ROOTVSMALL);
297 *evapRateCoeff[celli]*(YsCoeff*Ysolvent - YInf[celli]),
302 if (dm > limMass[celli])
304 evapRateCoeff[celli] *= limMass[celli]/dm;
307 evapRateInf[celli] = evapRateCoeff[celli]*YInf[celli];
308 evapRateCoeff[celli] *= YsCoeff;
338 deltaRho0Bydt*Ysolvent_()
344 + impingementRate*Ysolvent0_
364 dt*magSf*rhoInf*(evapRateCoeff*Ysolvent_ + evapRateInf)
389 correctModel(dt, availableMass, dMass, dEnergy,
zeroField());
397 correctModel(dt, availableMass, dMass, dEnergy, YInf);
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual const volScalarField & hs() const
Return the film sensible enthalpy [J/kg].
Wax solvent mixture evaporation model.
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,...
static constexpr const zero Zero
Global zero.
const dimensionSet dimVelocity
const dimensionSet dimDensity
scalar Sh(const scalar Re, const scalar Sc) const
Return Sherwood number as a function of Reynolds and Schmidt numbers.
virtual const word & name() const =0
Return the specie name.
static word timeName(const scalar t, const int precision=precision_)
const SLGThermo & thermo() const
Return const reference to the SLGThermo object.
Calculate the divergence of the given field.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
volScalarField & rhoSp()
Mass [kg/m2/s].
Ostream & endl(Ostream &os)
Add newline and flush stream.
const volVectorField & UPrimary() const
Velocity [m/s].
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
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.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Calculate the matrix for the divergence of the given field and flux.
#define forAll(list, i)
Loop across all elements in list.
const volScalarField & rhoPrimary() const
Density [kg/m3].
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Thermodynamic form of single-cell layer surface film model.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Base class for surface film models.
const filmThermoModel & filmThermo() const
Film thermo.
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
const fvMesh & regionMesh() const
Return the region mesh database.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
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 surfaceScalarField & phi() const
Return the film flux [kg.m/s].
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].
Calculate the matrix for implicit and explicit sources.
defineTypeNameAndDebug(kinematicSingleLayer, 0)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
void correctModel(const scalar dt, scalarField &availableMass, scalarField &dMass, scalarField &dEnergy, YInfType YInf)
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.
virtual ~waxSolventEvaporation()
Destructor.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const Time & time() const
Return the top-level database.
virtual const volScalarField & rho() const
Return the film density [kg/m3].
dimensionedScalar cbrt(const dimensionedScalar &ds)
virtual const volScalarField & deltaRho() const
Return the film thickness*density (helper field) [kg/m3].
const word & constant() const
Return constant name.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
#define WarningInFunction
Report a warning using Foam::Warning.
Calulate the matrix for the first temporal derivative.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Base class for surface film phase change models.