40template<
class Thermo,
class OtherThermo>
55 forAll(interfaceArea_, celli)
57 label status =
cutCell.calcSubCell(celli, isoAlpha_);
58 interfaceArea_[celli] = 0;
61 interfaceArea_[celli] =
70template<
class Thermo,
class OtherThermo>
121 isoAlpha_(
dict.getOrDefault<scalar>(
"isoAlpha", 0.5))
126 const typename OtherThermo::thermoType& toThermo =
134 Mv_.
value() = toThermo.W()*1
e-3;
136 if (Mv_.
value() == -1)
139 <<
" Please provide the molar weight (Mv) of vapour [g/mol] "
147template<
class Thermo,
class OtherThermo>
177 mesh.time().timeName(),
193 mesh.time().timeName(),
204 if (
sign(C_.value()) > 0)
206 rhov = this->pair().to().rho();
207 deltaT =
max(
T - Tactivate_,
T0);
211 rhov = this->pair().from().rho();
212 deltaT =
max(Tactivate_ -
T,
T0);
215 htc_ = 2*
mag(C_)/(2-
mag(C_))*(
L()*rhov/HerztKnudsConst);
217 mDotc_ = htc_*deltaT*interfaceArea_;
223template<
class Thermo,
class OtherThermo>
231 if (this->modelVariable_ == variable)
235 if (
sign(C_.value()) > 0)
237 return(coeff*
pos(refValue - Tactivate_));
241 return(coeff*
pos(Tactivate_ - refValue));
251template<
class Thermo,
class OtherThermo>
259 if (this->modelVariable_ == variable)
263 if (
sign(C_.value()) > 0)
265 return(-coeff*
pos(refValue - Tactivate_));
269 return(coeff*
pos(Tactivate_ - refValue));
Defines the attributes of an object for which implicit objectRegistry management is supported,...
word member() const
Return member (name without the extension)
Base class for interface composition models, templated on the two thermodynamic models either side of...
const OtherThermo & toThermo_
Other Thermo (to)
const pureMixture< ThermoType >::thermoType & getLocalThermo(const word &speciesName, const pureMixture< ThermoType > &globalThermo) const
Get a reference to the local thermo for a pure mixture.
Class for cutting a cell, celli, of an fvMesh, mesh_, at its intersection with an isosurface defined ...
Service routines for cutting a cell, celli, of an fvMesh, mesh_, at its intersection with a surface.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const Type & value() const
Return const reference to value.
Mesh data needed to do the Finite Volume discretisation.
Considering the Hertz Knudsen formula, which gives the evaporation-condensation flux based on the kin...
virtual tmp< volScalarField > KSu(label modelVariable, const volScalarField &field)
Explicit mass transfer coefficient.
virtual tmp< volScalarField > KSp(label modelVariable, const volScalarField &field)
Implicit mass transfer coefficient.
const word transferSpecie() const
Return the transferring species name.
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
constexpr scalar pi(M_PI)
const dimensionedScalar R
Universal gas constant: default SI units: [J/mol/K].
Different types of constants.
dimensionedScalar pos(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const dimensionSet dimless
Dimensionless.
dimensionedScalar sign(const dimensionedScalar &ds)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
const dimensionSet dimArea(sqr(dimLength))
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManip< error > abort(error &err)
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimDensity
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
#define forAll(list, i)
Loop across all elements in list.
const vector L(dict.get< vector >("L"))