40 template<
class Thermo,
class OtherThermo>
44 const fvMesh&
mesh = this->mesh_;
53 cutCellIso cutCell(
mesh, ap);
55 forAll(interfaceArea_, celli)
57 label status = cutCell.calcSubCell(celli, isoAlpha_);
58 interfaceArea_[celli] = 0;
61 interfaceArea_[celli] =
62 mag(cutCell.faceArea())/
mesh.V()[celli];
70 template<
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] " 147 template<
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_;
223 template<
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));
251 template<
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));
274 return tmp<volScalarField> ();
Different types of constants.
dimensionedScalar sign(const dimensionedScalar &ds)
const Type & value() const noexcept
Return const reference to value.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
const vector L(dict.get< vector >("L"))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static const volPointInterpolation & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject. Registered with typeName.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Base class for interface composition models, templated on the two thermodynamic models either side of...
dimensionedScalar sqrt(const dimensionedScalar &ds)
word member() const
Return member (name without the extension)
const dimensionedScalar R
Universal gas constant: default SI units: [J/mol/K].
const dimensionSet dimless
Dimensionless.
#define forAll(list, i)
Loop across all elements in list.
const OtherThermo & toThermo_
Other Thermo (to)
virtual const multiphaseInter::phaseModel & to() const
To phase.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar pos(const dimensionedScalar &ds)
const pureMixture< ThermoType >::thermoType & getLocalThermo(const word &speciesName, const pureMixture< ThermoType > &globalThermo) const
Get a reference to the local thermo for a pure mixture.
const word transferSpecie() const
Return the transferring species name.
const dimensionedScalar e
Elementary charge.
virtual const multiphaseInter::phaseModel & from() const
From phase.
Considering the Hertz Knudsen formula, which gives the evaporation-condensation flux based on the kin...
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
constexpr scalar pi(M_PI)
virtual tmp< volScalarField > KSu(label modelVariable, const volScalarField &field)
Explicit mass transfer coefficient.
errorManip< error > abort(error &err)
virtual tmp< volScalarField > Kexp(const volScalarField &field)
Explicit total mass transfer coefficient.
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
tmp< volScalarField > rho() const
Return the phase density.
const dimensionSet dimDensity
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual tmp< volScalarField > KSp(label modelVariable, const volScalarField &field)
Implicit mass transfer coefficient.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
kineticGasEvaporation(const dictionary &dict, const phasePair &pair)
Construct from components.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Defines the attributes of an object for which implicit objectRegistry management is supported...
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
const dimensionSet dimArea(sqr(dimLength))
static constexpr const zero Zero
Global zero (0)