Go to the documentation of this file.
49 solidificationMeltingSource,
61 { thermoMode::mdThermo,
"thermo" },
62 { thermoMode::mdLookup,
"lookup" },
69 Foam::fv::solidificationMeltingSource::Cp()
const
75 const basicThermo&
thermo =
83 if (CpName_ ==
"CpRef")
85 const scalar CpRef =
coeffs_.
get<scalar>(
"CpRef");
104 extrapolatedCalculatedFvPatchScalarField::typeName
126 void Foam::fv::solidificationMeltingSource::update(
const volScalarField&
Cp)
128 if (curTimeIndex_ == mesh_.time().timeIndex())
135 Info<<
type() <<
": " << name_ <<
" - updating phase indicator" <<
endl;
145 label celli = cells_[i];
147 scalar Tc =
T[celli];
148 scalar Cpc =
Cp[celli];
149 scalar alpha1New = alpha1_[celli] + relax_*Cpc*(Tc - Tmelt_)/L_;
151 alpha1_[celli] =
max(0,
min(alpha1New, 1));
152 deltaT_[i] = Tc - Tmelt_;
155 alpha1_.correctBoundaryConditions();
157 curTimeIndex_ = mesh_.time().timeIndex();
163 Foam::fv::solidificationMeltingSource::solidificationMeltingSource
165 const word& sourceName,
166 const word& modelType,
172 Tmelt_(coeffs_.get<scalar>(
"Tmelt")),
173 L_(coeffs_.get<scalar>(
"L")),
174 relax_(coeffs_.lookupOrDefault(
"relax", 0.9)),
175 mode_(thermoModeTypeNames_.get(
"thermoMode", coeffs_)),
176 rhoRef_(coeffs_.get<scalar>(
"rhoRef")),
177 TName_(coeffs_.lookupOrDefault<
word>(
"T",
"T")),
178 CpName_(coeffs_.lookupOrDefault<
word>(
"Cp",
"Cp")),
179 UName_(coeffs_.lookupOrDefault<
word>(
"U",
"U")),
180 phiName_(coeffs_.lookupOrDefault<
word>(
"phi",
"phi")),
181 Cu_(coeffs_.lookupOrDefault<scalar>(
"Cu", 100000)),
182 q_(coeffs_.lookupOrDefault(
"q", 0.001)),
183 beta_(coeffs_.get<scalar>(
"beta")),
189 mesh.time().timeName(),
196 zeroGradientFvPatchScalarField::typeName
199 deltaT_(cells_.size(), 0)
201 fieldNames_.setSize(2);
202 fieldNames_[0] = UName_;
211 fieldNames_[1] =
thermo.he().name();
216 fieldNames_[1] = TName_;
222 <<
"Unhandled thermo mode: " << thermoModeTypeNames_[mode_]
227 applied_.setSize(fieldNames_.size(),
false);
277 label celli = cells_[i];
279 scalar Vc = V[celli];
280 scalar alpha1c = alpha1_[celli];
282 scalar S = -Cu_*
sqr(1.0 - alpha1c)/(
pow3(alpha1c) + q_);
283 vector Sb = rhoRef_*
g*beta_*deltaT_[i];
307 coeffs_.readEntry(
"Tmelt", Tmelt_);
308 coeffs_.readEntry(
"L", L_);
310 coeffs_.readIfPresent(
"relax", relax_);
312 thermoModeTypeNames_.readEntry(
"thermoMode", coeffs_, mode_);
314 coeffs_.readEntry(
"rhoRef", rhoRef_);
315 coeffs_.readIfPresent(
"T", TName_);
316 coeffs_.readIfPresent(
"U", UName_);
318 coeffs_.readIfPresent(
"Cu", Cu_);
319 coeffs_.readIfPresent(
"q", q_);
321 coeffs_.readIfPresent(
"beta", beta_);
int debug
Static debugging option.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
A class for handling words, derived from Foam::string.
Cell-set options abtract base class. Provides a base set of controls, e.g.:
A class for managing temporary objects.
static constexpr const zero Zero
Global zero.
const dimensionSet dimEnergy
static word timeName(const scalar t, const int precision=precision_)
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Abstract base-class for fluid and solid thermodynamic properties.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const word name_
Source name.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const GeometricField< Type, fvPatchField, volMesh > & psi() const
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const fvMesh & mesh_
Reference to the mesh database.
#define forAll(list, i)
Loop across all elements in list.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (uses stdout - output is on the master only)
dictionary coeffs_
Dictionary containing source coefficients.
virtual bool read(const dictionary &dict)
Read source dictionary.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Add explicit contribution to enthalpy equation.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const Type & lookupObject(const word &name, const bool recursive=false) const
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
word dictName() const
The local dictionary name (final part of scoped name)
static const Enum< thermoMode > thermoModeTypeNames_
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.
Mesh data needed to do the Finite Volume discretisation.
const uniformDimensionedVectorField & g
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
errorManip< error > abort(error &err)
defineTypeNameAndDebug(option, 0)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual bool read(const dictionary &dict)
Read source dictionary.
static void apply(bitSet &selection, const Detail::parcelSelection::actionType action, const Predicate &accept, const UList< Type > &list, const AccessOp &aop)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const volScalarField & Cp
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
const Time & time() const
Return the top-level database.
static const gravity & New(const Time &runTime)
Construct on Time.