41filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
44 HashTable<const filmModelType*> models
49 if (iter()->regionMesh().
name() == filmRegionName_)
55 DynamicList<word> modelNames;
58 modelNames.append(iter()->regionMesh().
name());
62 <<
"Unable to locate film region " << filmRegionName_
63 <<
". Available regions include: " << modelNames
66 return **models.begin();
72filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
75 HashTable<const pyrolysisModelType*> models =
80 if (iter()->regionMesh().
name() == pyrolysisRegionName_)
86 DynamicList<word> modelNames;
89 modelNames.append(iter()->regionMesh().
name());
93 <<
"Unable to locate pyrolysis region " << pyrolysisRegionName_
94 <<
". Available regions include: " << modelNames
97 return **models.begin();
103filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
104filmPyrolysisRadiativeCoupledMixedFvPatchScalarField
110 mixedFvPatchScalarField(
p, iF),
119 filmRegionName_(
"surfaceFilmProperties"),
120 pyrolysisRegionName_(
"pyrolysisProperties"),
121 TnbrName_(
"undefined-Tnbr"),
122 qrName_(
"undefined-qr"),
123 convectiveScaling_(1.0),
127 this->refValue() = 0.0;
128 this->refGrad() = 0.0;
129 this->valueFraction() = 1.0;
133filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
134filmPyrolysisRadiativeCoupledMixedFvPatchScalarField
142 mixedFvPatchScalarField(psf,
p, iF, mapper),
144 filmRegionName_(psf.filmRegionName_),
145 pyrolysisRegionName_(psf.pyrolysisRegionName_),
146 TnbrName_(psf.TnbrName_),
147 qrName_(psf.qrName_),
148 convectiveScaling_(psf.convectiveScaling_),
149 filmDeltaDry_(psf.filmDeltaDry_),
150 filmDeltaWet_(psf.filmDeltaWet_)
154filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
155filmPyrolysisRadiativeCoupledMixedFvPatchScalarField
162 mixedFvPatchScalarField(
p, iF),
166 dict.getOrDefault<
word>(
"filmRegion",
"surfaceFilmProperties")
170 dict.getOrDefault<
word>(
"pyrolysisRegion",
"pyrolysisProperties")
174 convectiveScaling_(
dict.getOrDefault<scalar>(
"convectiveScaling", 1)),
175 filmDeltaDry_(
dict.get<scalar>(
"filmDeltaDry")),
176 filmDeltaWet_(
dict.get<scalar>(
"filmDeltaWet"))
178 if (!isA<mappedPatchBase>(this->patch().patch()))
182 <<
"\n for patch " <<
p.
name()
183 <<
" of field " << internalField().name()
184 <<
" in file " << internalField().objectPath()
202 valueFraction() = 1.0;
207filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
208filmPyrolysisRadiativeCoupledMixedFvPatchScalarField
214 mixedFvPatchScalarField(psf, iF),
216 filmRegionName_(psf.filmRegionName_),
217 pyrolysisRegionName_(psf.pyrolysisRegionName_),
218 TnbrName_(psf.TnbrName_),
219 qrName_(psf.qrName_),
220 convectiveScaling_(psf.convectiveScaling_),
221 filmDeltaDry_(psf.filmDeltaDry_),
222 filmDeltaWet_(psf.filmDeltaWet_)
233 mixedFvPatchScalarField::autoMap(mapper);
244 mixedFvPatchScalarField::rmap(ptf, addr);
265 refCast<const mappedPatchBase>(patch().patch());
267 const label patchi = patch().index();
272 refCast<const fvMesh>(nbrMesh).boundary()[nbrPatchi];
287 scalarField nbrIntFld(nbrField.patchInternalField());
311 label coupledPatchi = -1;
312 if (pyrolysisRegionName_ ==
mesh.
name())
314 coupledPatchi = patchi;
315 if (qrName_ !=
"none")
323 coupledPatchi = nbrPatch.
index();
324 if (qrName_ !=
"none")
332 <<
type() <<
" condition is intended to be applied to either the "
333 <<
"primary or pyrolysis regions only"
375 (filmDelta - filmDeltaDry_)/(filmDeltaWet_ - filmDeltaDry_),
382 scalarField qConv(ratio*htcwfilm*(Tfilm - Tp)*convectiveScaling_);
388 valueFraction() =
alpha/(
alpha + (1.0 - ratio)*myKDelta);
390 refValue() = ratio*Tfilm + (1.0 - ratio)*(KDeltaNbr*nbrIntFld)/
alpha;
392 mixedFvPatchScalarField::updateCoeffs();
396 scalar Qc =
gSum(qConv*patch().magSf());
397 scalar qr =
gSum(qRad*patch().magSf());
398 scalar Qt =
gSum((qConv + qRad)*patch().magSf());
401 << patch().name() <<
':'
402 << this->internalField().name() <<
" <- "
403 << nbrMesh.
name() <<
':'
404 << nbrPatch.
name() <<
':'
405 << this->internalField().name() <<
" :" <<
nl
406 <<
" convective heat[W] : " << Qc <<
nl
407 <<
" radiative heat [W] : " << qr <<
nl
408 <<
" total heat [W] : " << Qt <<
nl
409 <<
" wall temperature "
410 <<
" min:" <<
gMin(*
this)
411 <<
" max:" <<
gMax(*
this)
423 mixedFvPatchScalarField::write(
os);
427 "surfaceFilmProperties",
433 "pyrolysisProperties",
CGAL::Exact_predicates_exact_constructions_kernel K
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const noexcept
Return the object name.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
void size(const label n)
Older name for setAddressableSize.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Mixed boundary condition for temperature, to be used in the flow and pyrolysis regions when a film re...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
tmp< scalarField > K() const
Get corresponding K field.
Foam::regionModels::pyrolysisModels::pyrolysisModel pyrolysisModelType
Foam::regionModels::surfaceFilmModels::thermoSingleLayer filmModelType
virtual bool write()
Write the output fields.
const word & name() const
Return reference to name.
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void operator=(const UList< scalar > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual const word & name() const
Return name.
label index() const
Return the index of this patch in the fvBoundaryMesh.
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
const scalarField & deltaCoeffs() const
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
const polyMesh & sampleMesh() const
Get the region mesh.
const polyPatch & samplePolyPatch() const
Get the patch on the region.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
HashTable< const Type * > lookupClass(const bool strict=false) const
Return all objects with a class satisfying isA<Type>
label index() const noexcept
The index of this patch in the boundaryMesh.
Mesh consisting of general polyhedral cells.
Lookup type of boundary radiation properties.
Base class for pyrolysis models.
const Time & time() const
Return the reference to the time database.
label nbrCoupledPatchID(const regionModel &nbrRegion, const label regionPatchi) const
Return the coupled patch ID paired with coupled patch.
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
void toPrimary(const label regionPatchi, List< Type > ®ionField) const
Convert a local region field to the primary region.
tmp< Foam::Field< Type > > mapRegionPatchField(const regionModel &nbrRegion, const label regionPatchi, const label nbrPatchi, const Field< Type > &nbrField, const bool flip=false) const
Map patch field from another region model to local patch.
Thermodynamic form of single-cell layer surface film model.
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
const heatTransferModel & htcw() const
Return const access to the (wall) heat transfer model.
Common functions used in temperature coupled boundaries.
virtual void autoMap(const fvPatchFieldMapper &)=0
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)=0
Reverse map the given fvPatchField onto this fvPatchField.
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
A class for handling words, derived from Foam::string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Type gSum(const FieldField< Field, Type > &f)
To & refCast(From &r)
Reference type cast template function.
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.
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.
errorManip< error > abort(error &err)
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
Type gMin(const FieldField< Field, Type > &f)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Type gMax(const FieldField< Field, Type > &f)
constexpr char nl
The newline '\n' character (0x0a)
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
static const char *const typeName
The type name used in ensight case files.