Go to the documentation of this file.
41 filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
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();
72 filmPyrolysisRadiativeCoupledMixedFvPatchScalarField::
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();
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;
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_)
162 mixedFvPatchScalarField(
p, iF),
175 filmDeltaDry_(
dict.
get<scalar>(
"filmDeltaDry")),
176 filmDeltaWet_(
dict.
get<scalar>(
"filmDeltaWet"))
178 if (!isA<mappedPatchBase>(this->
patch().
patch()))
181 <<
"' not type '" << mappedPatchBase::typeName <<
"'"
182 <<
"\n for patch " <<
p.name()
183 <<
" of field " << internalField().name()
184 <<
" in file " << internalField().objectPath()
202 valueFraction() = 1.0;
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_)
237 refCast<const mappedPatchBase>(
patch().
patch());
239 const label patchi =
patch().index();
244 refCast<const fvMesh>(nbrMesh).boundary()[nbrPatchi];
259 scalarField nbrIntFld(nbrField.patchInternalField());
283 label coupledPatchi = -1;
284 if (pyrolysisRegionName_ ==
mesh.name())
286 coupledPatchi = patchi;
287 if (qrName_ !=
"none")
295 coupledPatchi = nbrPatch.
index();
296 if (qrName_ !=
"none")
304 <<
type() <<
" condition is intended to be applied to either the "
305 <<
"primary or pyrolysis regions only"
347 (filmDelta - filmDeltaDry_)/(filmDeltaWet_ - filmDeltaDry_),
354 scalarField qConv(ratio*htcwfilm*(Tfilm - Tp)*convectiveScaling_);
360 valueFraction() =
alpha/(
alpha + (1.0 - ratio)*myKDelta);
362 refValue() = ratio*Tfilm + (1.0 - ratio)*(KDeltaNbr*nbrIntFld)/
alpha;
364 mixedFvPatchScalarField::updateCoeffs();
370 scalar Qt =
gSum((qConv + qRad)*
patch().magSf());
373 <<
patch().name() <<
':'
374 << this->internalField().name() <<
" <- "
375 << nbrMesh.
name() <<
':'
376 << nbrPatch.
name() <<
':'
377 << this->internalField().name() <<
" :" <<
nl
378 <<
" convective heat[W] : " << Qc <<
nl
379 <<
" radiative heat [W] : " << qr <<
nl
380 <<
" total heat [W] : " << Qt <<
nl
381 <<
" wall temperature "
382 <<
" min:" <<
gMin(*
this)
383 <<
" max:" <<
gMax(*
this)
399 "surfaceFilmProperties",
405 "pyrolysisProperties",
410 os.
writeEntry(
"convectiveScaling", convectiveScaling_);
Base class for pyrolysis models.
int debug
Static debugging option.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const word & name() const
Return name.
A class for handling words, derived from Foam::string.
static constexpr const zero Zero
Global zero (0)
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Type gAverage(const FieldField< Field, Type > &f)
Mixed boundary condition for temperature, to be used in the flow and pyrolysis regions when a film re...
virtual void write(Ostream &) const
Write.
const polyPatch & samplePolyPatch() const
Get the patch on the region.
Common functions used in temperature coupled boundaries.
const Time & time() const
Return the reference to the time database.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
virtual const word & name() const
Return name.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
Foam::regionModels::pyrolysisModels::pyrolysisModel pyrolysisModelType
CGAL::Exact_predicates_exact_constructions_kernel K
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.
HashTable< const Type * > lookupClass(const bool strict=false) const
Return all objects with a class satisfying isA<Type>
const polyMesh & sampleMesh() const
Get the region mesh.
tmp< scalarField > K() const
Get corresponding K field.
Thermodynamic form of single-cell layer surface film model.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
virtual tmp< volScalarField > h() const =0
Return the heat transfer coefficient [W/m2/K].
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual tmp< scalarField > alpha(const scalarField &Tp) const
Given patch temperature calculate corresponding alphaEff field.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
label nbrCoupledPatchID(const regionModel &nbrRegion, const label regionPatchi) const
Return the coupled patch ID paired with coupled patch.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
label index() const
Return the index of this patch in the fvBoundaryMesh.
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
To & refCast(From &r)
Reference type cast template function.
errorManipArg< error, int > exit(error &err, const int errNo=1)
filmPyrolysisRadiativeCoupledMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
forAllConstIters(mixture.phases(), phase)
const heatTransferModel & htcw() const
Return const access to the (wall) heat transfer model.
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.
const std::string patch
OpenFOAM patch number as a std::string.
const scalarField & deltaCoeffs() const
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Foam::regionModels::surfaceFilmModels::thermoSingleLayer filmModelType
Foam::fvPatchFieldMapper.
void write(Ostream &os) const
Write.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Type gMin(const FieldField< Field, Type > &f)
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
virtual void operator=(const UList< scalar > &)
label index() const
The index of this patch in the boundaryMesh.
Type gMax(const FieldField< Field, Type > &f)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const word & name() const
Return reference to name.