51 mixedFvPatchScalarField(
p, iF),
58 TnbrName_(
"undefined-Tnbr")
60 this->refValue() = 0.0;
61 this->refGrad() = 0.0;
62 this->valueFraction() = 1.0;
75 mixedFvPatchScalarField(ptf,
p, iF, mapper),
83 TnbrName_(ptf.TnbrName_),
84 thicknessLayers_(ptf.thicknessLayers_),
85 thicknessLayer_(ptf.thicknessLayer_.clone(
p.
patch())),
86 kappaLayers_(ptf.kappaLayers_),
87 kappaLayer_(ptf.kappaLayer_.clone(
p.
patch()))
99 mixedFvPatchScalarField(
p, iF),
109 if (!isA<mappedPatchBase>(this->
patch().
patch()))
112 <<
"' not type '" << mappedPatchBase::typeName <<
"'" 113 <<
"\n for patch " <<
p.name()
114 <<
" of field " << internalField().name()
115 <<
" in file " << internalField().objectPath()
120 <<
"This BC has been superseded by " 121 <<
"compressible::turbulentTemperatureRadCoupledMixed " 122 <<
"which has more functionalities and it can handle " 123 <<
"the assemble coupled option for energy. " 127 if (
dict.readIfPresent(
"thicknessLayers", thicknessLayers_))
129 dict.readEntry(
"kappaLayers", kappaLayers_);
148 if (
dict.found(
"refValue"))
160 valueFraction() = 1.0;
182 const turbulentTemperatureCoupledBaffleMixedFvPatchScalarField& wtcsf,
183 const DimensionedField<scalar, volMesh>& iF
186 mixedFvPatchScalarField(wtcsf, iF),
187 temperatureCoupledBase(
patch(), wtcsf),
188 mappedPatchFieldBase<scalar>
190 mappedPatchFieldBase<scalar>::mapper(
patch(), iF),
194 TnbrName_(wtcsf.TnbrName_),
195 thicknessLayers_(wtcsf.thicknessLayers_),
196 thicknessLayer_(wtcsf.thicknessLayer_.clone(
patch().
patch())),
197 kappaLayers_(wtcsf.kappaLayers_),
198 kappaLayer_(wtcsf.kappaLayer_.clone(
patch().
patch()))
208 mixedFvPatchScalarField(wtcsf),
216 TnbrName_(wtcsf.TnbrName_),
217 thicknessLayers_(wtcsf.thicknessLayers_),
218 thicknessLayer_(wtcsf.thicknessLayer_.clone(
patch().
patch())),
219 kappaLayers_(wtcsf.kappaLayers_),
231 mixedFvPatchScalarField::autoMap(
mapper);
237 kappaLayer_().autoMap(
mapper);
248 mixedFvPatchScalarField::rmap(ptf, addr);
260 thicknessLayer_().rmap(tiptf.thicknessLayer_(), addr);
261 kappaLayer_().rmap(tiptf.kappaLayer_(), addr);
276 if (thicknessLayer_ || thicknessLayers_.
size())
285 const scalar t = db().time().timeOutputValue();
287 thicknessLayer_().value(t)
288 /kappaLayer_().value(t);
290 if (thicknessLayers_.
size())
292 forAll(thicknessLayers_, iLayer)
294 KDelta += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
321 const mappedPatchBase& mpp =
325 this->internalField()
330 const tmp<scalarField> myKDelta = kappaTp*
patch().deltaCoeffs();
338 const auto& nbrMesh = refCast<const fvMesh>(mpp.sampleMesh());
339 const label nbrPatchID = mpp.samplePolyPatch().index();
340 const auto& nbrPatch = nbrMesh.boundary()[nbrPatchID];
356 nbrIntFld = nbrField.patchInternalField();
357 nbrKDelta = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
363 nbrIntFld = patchInternalField();
364 nbrKDelta = myKDelta.ref();
385 this->refValue() = nbrIntFld;
386 this->refGrad() = 0.0;
387 this->valueFraction() = nbrKDelta/(nbrKDelta + myKDelta());
389 mixedFvPatchScalarField::updateCoeffs();
395 Info<<
patch().boundaryMesh().mesh().name() <<
':' 396 <<
patch().name() <<
':' 397 << this->internalField().name() <<
" <- " 398 << mpp.sampleRegion() <<
':' 399 << mpp.samplePatch() <<
':' 400 << this->internalField().name() <<
" :" 401 <<
" heat transfer rate:" << Q
402 <<
" walltemperature " 403 <<
" min:" <<
gMin(*
this)
404 <<
" max:" <<
gMax(*
this)
422 <<
"This BC does not support energy coupling " 423 <<
"Use compressible::turbulentTemperatureRadCoupledMixed " 424 <<
"which has more functionalities and it can handle " 425 <<
"the assemble coupled option for energy. " 431 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::coeffs
439 <<
"This BC does not support energy coupling " 440 <<
"Use compressible::turbulentTemperatureRadCoupledMixed " 441 <<
"which has more functionalities and it can handle " 442 <<
"the assemble coupled option for energy. " 467 os.writeEntry(
"Tnbr", TnbrName_);
470 thicknessLayer_().writeData(
os);
471 kappaLayer_().writeData(
os);
473 if (thicknessLayers_.
size())
488 turbulentTemperatureCoupledBaffleMixedFvPatchScalarField
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Mixed boundary condition for temperature, to be used for heat-transfer on back-to-back baffles...
void size(const label n)
Older name for setAddressableSize.
turbulentTemperatureCoupledBaffleMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
static const mappedPatchBase & mapper(const fvPatch &p, const DimensionedField< Type, volMesh > &iF)
Check that patch is of correct type.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Type gMin(const FieldField< Field, Type > &f)
virtual void write(Ostream &os) const
Write.
Type & refCast(U &obj)
A dynamic_cast (for references) that generates FatalError on failed casts, uses the virtual type() me...
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)=0
Reverse map the given fvPatchField onto this fvPatchField.
static int & msgType() noexcept
Message tag of standard messages.
virtual void autoMap(const fvPatchFieldMapper &)=0
Map (and resize as needed) from self given a mapping object.
void write(Ostream &os) const
Write.
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
makePatchTypeField(fvPatchScalarField, alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Type gSum(const FieldField< Field, Type > &f)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
fvPatchField< scalar > fvPatchScalarField
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A FieldMapper for finite-volume patch fields.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
errorManip< error > abort(error &err)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
int debug
Static debugging option.
Type gMax(const FieldField< Field, Type > &f)
OBJstream os(runTime.globalPath()/outputName)
Common functions used in temperature coupled boundaries.
#define WarningInFunction
Report a warning using Foam::Warning.
Type gAverage(const FieldField< Field, Type > &f)
virtual void operator=(const UList< scalar > &)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const std::string patch
OpenFOAM patch number as a std::string.
static autoPtr< PatchFunction1< Type > > NewIfPresent(const polyPatch &pp, const word &entryName, const dictionary &dict, const bool faceValues=true)
An optional selector.
messageStream Info
Information stream (stdout output on master, null elsewhere)
void writeEntry(Ostream &os) const
Write the UList with its compound type.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field. Override temperatureCoupledBase::kappa to in...
A class for managing temporary objects.
void distribute(const word &fieldName, Field< T > &newValues) const
Wrapper for mapDistribute::distribute that knows about dabase mapping.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void manipulateMatrix(fvMatrix< scalar > &m, const label iMatrix, const direction cmpt)
Manipulate matrix.
virtual void write(Ostream &os) const
Write.