45 VoFSolidificationMeltingSource,
54void Foam::fv::VoFSolidificationMeltingSource::update()
64 <<
" - updating solid phase fraction" <<
endl;
69 const twoPhaseMixtureThermo&
thermo
83 const label celli =
cells_[i];
85 alphaSolid_[celli] =
min
87 relax_*alphaVoF[celli]*alphaSolidT_->value(TVoF[celli])
88 + (1 - relax_)*alphaSolid_[celli],
99Foam::word Foam::fv::VoFSolidificationMeltingSource::alphaSolidName()
const
101 const twoPhaseMixtureThermo&
thermo
103 mesh_.lookupObject<twoPhaseMixtureThermo>
119 const word& sourceName,
120 const word& modelType,
121 const dictionary&
dict,
125 fv::cellSetOption(sourceName, modelType,
dict,
mesh),
126 alphaSolidT_(Function1<scalar>::
New(
"alphaSolidT", coeffs_, &
mesh)),
128 relax_(coeffs_.getOrDefault(
"relax", 0.9)),
129 Cu_(coeffs_.getOrDefault<scalar>(
"Cu", 100000)),
130 q_(coeffs_.getOrDefault<scalar>(
"q", 0.001)),
138 IOobject::READ_IF_PRESENT,
143 zeroGradientFvPatchScalarField::typeName
147 fieldNames_.resize(2);
148 fieldNames_[0] =
"U";
149 fieldNames_[1] =
"T";
151 fv::option::resetApplied();
159 fvMatrix<scalar>& eqn,
163 apply(geometricOneField(), eqn);
170 fvMatrix<scalar>& eqn,
180 fvMatrix<vector>& eqn,
186 Info<<
type() <<
": applying source to " << eqn.psi().name() <<
endl;
196 const label celli = cells_[i];
197 const scalar Vc = V[celli];
198 const scalar alphaFluid = 1 - alphaSolid_[celli];
200 const scalar S = Cu_*
sqr(1 - alphaFluid)/(
pow3(alphaFluid) + q_);
210 fvMatrix<vector>& eqn,
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
void correctBoundaryConditions()
Correct boundary field.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
label timeIndex() const noexcept
Return current time index.
static const word dictName
const Time & time() const
Return the top-level database.
Solidification and melting model for VoF simulations.
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Add explicit contribution to enthalpy equation.
labelList cells_
Set of cells to apply source to.
const fvMesh & mesh_
Reference to the mesh database.
const word name_
Source name.
const Type & lookupObject(const word &name, const bool recursive=false) const
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimEnergy
dimensionedScalar pow3(const dimensionedScalar &ds)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
static void apply(bitSet &selection, const Detail::parcelSelection::actionType action, const Predicate &accept, const UList< Type > &list, const AccessOp &aop)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
static constexpr const zero Zero
Global zero (0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
#define forAll(list, i)
Loop across all elements in list.