Go to the documentation of this file.
51 mixedFvPatchScalarField(
p, iF),
65 TnbrName_(
"undefined-Tnbr")
67 this->refValue() = 0.0;
68 this->refGrad() = 0.0;
69 this->valueFraction() = 1.0;
82 mixedFvPatchScalarField(ptf,
p, iF, mapper),
90 TnbrName_(ptf.TnbrName_),
91 thicknessLayers_(ptf.thicknessLayers_),
92 thicknessLayer_(ptf.thicknessLayer_.clone(
p.patch())),
93 kappaLayers_(ptf.kappaLayers_),
94 kappaLayer_(ptf.kappaLayer_.clone(
p.patch()))
106 mixedFvPatchScalarField(
p, iF),
116 if (!isA<mappedPatchBase>(this->
patch().
patch()))
119 <<
"' not type '" << mappedPatchBase::typeName <<
"'"
120 <<
"\n for patch " <<
p.name()
121 <<
" of field " << internalField().name()
122 <<
" in file " << internalField().objectPath()
127 <<
"This BC has been superseded by "
128 <<
"compressible::turbulentTemperatureRadCoupledMixed "
129 <<
"which has more functionalities and it can handle "
130 <<
"the assemble coupled option for energy. "
134 if (
dict.readIfPresent(
"thicknessLayers", thicknessLayers_))
136 dict.readEntry(
"kappaLayers", kappaLayers_);
155 if (
dict.found(
"refValue"))
167 valueFraction() = 1.0;
193 mixedFvPatchScalarField(wtcsf, iF),
201 TnbrName_(wtcsf.TnbrName_),
202 thicknessLayers_(wtcsf.thicknessLayers_),
203 thicknessLayer_(wtcsf.thicknessLayer_.clone(
patch().
patch())),
204 kappaLayers_(wtcsf.kappaLayers_),
205 kappaLayer_(wtcsf.kappaLayer_.clone(
patch().
patch()))
215 mixedFvPatchScalarField(wtcsf),
223 TnbrName_(wtcsf.TnbrName_),
224 thicknessLayers_(wtcsf.thicknessLayers_),
225 thicknessLayer_(wtcsf.thicknessLayer_.clone(
patch().
patch())),
226 kappaLayers_(wtcsf.kappaLayers_),
227 kappaLayer_(wtcsf.kappaLayer_.clone(
patch().
patch()))
238 mixedFvPatchScalarField::autoMap(mapper);
243 thicknessLayer_().autoMap(mapper);
244 kappaLayer_().autoMap(mapper);
255 mixedFvPatchScalarField::rmap(ptf, addr);
267 thicknessLayer_().rmap(tiptf.thicknessLayer_(), addr);
268 kappaLayer_().rmap(tiptf.kappaLayer_(), addr);
283 if (thicknessLayer_ || thicknessLayers_.size())
292 const scalar t = db().time().timeOutputValue();
294 thicknessLayer_().value(t)
295 /kappaLayer_().value(t);
297 if (thicknessLayers_.size())
299 forAll(thicknessLayers_, iLayer)
301 KDelta += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
308 tk = KDelta/
patch().deltaCoeffs();
332 this->internalField()
345 const auto& nbrMesh = refCast<const fvMesh>(mpp.
sampleMesh());
347 const auto& nbrPatch = nbrMesh.boundary()[nbrPatchID];
363 nbrIntFld = nbrField.patchInternalField();
364 nbrKDelta = nbrField.
kappa(nbrField)*nbrPatch.deltaCoeffs();
370 nbrIntFld = patchInternalField();
371 nbrKDelta = myKDelta.
ref();
392 this->refValue() = nbrIntFld;
393 this->refGrad() = 0.0;
394 this->valueFraction() = nbrKDelta/(nbrKDelta + myKDelta());
396 mixedFvPatchScalarField::updateCoeffs();
402 Info<<
patch().boundaryMesh().mesh().name() <<
':'
403 <<
patch().name() <<
':'
404 << this->internalField().name() <<
" <- "
407 << this->internalField().name() <<
" :"
408 <<
" heat transfer rate:" << Q
409 <<
" walltemperature "
410 <<
" min:" <<
gMin(*
this)
411 <<
" max:" <<
gMax(*
this)
429 <<
"This BC does not support energy coupling "
430 <<
"Use compressible::turbulentTemperatureRadCoupledMixed "
431 <<
"which has more functionalities and it can handle "
432 <<
"the assemble coupled option for energy. "
438 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::coeffs
446 <<
"This BC does not support energy coupling "
447 <<
"Use compressible::turbulentTemperatureRadCoupledMixed "
448 <<
"which has more functionalities and it can handle "
449 <<
"the assemble coupled option for energy. "
474 os.writeEntry(
"Tnbr", TnbrName_);
477 thicknessLayer_().writeData(
os);
478 kappaLayer_().writeData(
os);
480 if (thicknessLayers_.size())
482 thicknessLayers_.writeEntry(
"thicknessLayers",
os);
483 kappaLayers_.writeEntry(
"kappaLayers",
os);
int debug
Static debugging option.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual void manipulateMatrix(fvMatrix< scalar > &m, const label iMatrix, const direction cmpt)
Manipulate matrix.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
Type gAverage(const FieldField< Field, Type > &f)
const polyPatch & samplePolyPatch() const
Get the patch on the region.
turbulentTemperatureCoupledBaffleMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Common functions used in temperature coupled boundaries.
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool sameWorld() const
Is sample world the local world?
Type gSum(const FieldField< Field, Type > &f)
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
virtual void autoMap(const fvPatchFieldMapper &)=0
Map (and resize as needed) from self given a mapping object.
makePatchTypeField(fvPatchScalarField, alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField)
#define forAll(list, i)
Loop across all elements in list.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const polyMesh & sampleMesh() const
Get the region mesh.
virtual tmp< scalarField > kappa(const scalarField &Tp) const
messageStream Info
Information stream (stdout output on master, null elsewhere)
void distribute(const word &fieldName, Field< T > &newValues) const
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Mixed boundary condition for temperature, to be used for heat-transfer on back-to-back baffles....
GeometricField< scalar, fvPatchField, volMesh > volScalarField
virtual void write(Ostream &os) const
Write.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
static const mappedPatchBase & mapper(const fvPatch &p, const DimensionedField< Type, volMesh > &iF)
Check that patch is of correct type.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)=0
Reverse map the given fvPatchField onto this fvPatchField.
label index() const noexcept
The index of this patch in the boundaryMesh.
errorManip< error > abort(error &err)
To & refCast(From &r)
Reference type cast template function.
const word & sampleRegion() const
Region to sample.
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void write(Ostream &os) const
Write.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static int & msgType() noexcept
Message tag of standard messages.
const std::string patch
OpenFOAM patch number as a std::string.
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Foam::fvPatchFieldMapper.
void write(Ostream &os) const
Write.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Type gMin(const FieldField< Field, Type > &f)
virtual void operator=(const UList< scalar > &)
#define WarningInFunction
Report a warning using Foam::Warning.
Type gMax(const FieldField< Field, Type > &f)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static autoPtr< PatchFunction1< Type > > NewIfPresent(const polyPatch &pp, const word &entryName, const dictionary &dict, const bool faceValues=true)
An optional selector.
const word & samplePatch() const
Patch (only if NEARESTPATCHFACE)