29 #include "surfaceInterpolate.H" 36 namespace heatTransferCoeffModels
51 void Foam::heatTransferCoeffModels::faceZoneReferenceTemperature::
57 const word faceZoneName(
dict.get<
word>(
"referenceFaceZone"));
64 <<
"referenceFaceZone: " << faceZoneName
65 <<
" does not exist in referenceRegion: " << refRegionName_
71 label numFaces = fZone.
size();
76 <<
"referenceFaceZone: " << faceZoneName
77 <<
" contains no faces." 90 const label meshFacei = fZone[i];
94 label facePatchId = -1;
102 if (isA<emptyPolyPatch>(
pp))
107 const auto* cpp = isA<coupledPolyPatch>(
pp);
109 if (cpp && !cpp->owner())
119 faceId_[numFaces] =
faceId;
120 facePatchId_[numFaces] = facePatchId;
129 facePatchId_.
resize(numFaces);
133 Foam::scalar Foam::heatTransferCoeffModels::faceZoneReferenceTemperature::
134 faceZoneAverageTemperature()
137 mesh_.objectRegistry::db().lookupObject<fvMesh>(refRegionName_);
149 const label facei = faceId_[i];
150 if (facePatchId_[i] != -1)
152 const label patchi = facePatchId_[i];
153 const scalar sf = magSf.boundaryField()[patchi][facei];
155 Tmean += Tf.boundaryField()[patchi][facei]*sf;
160 const scalar sf = magSf[facei];
161 Tmean += Tf[facei]*sf;
165 reduce(Tmean, sumOp<scalar>());
166 reduce(sumMagSf, sumOp<scalar>());
174 void Foam::heatTransferCoeffModels::faceZoneReferenceTemperature::htc
177 const FieldField<Field, scalar>& q
190 const scalar Tref = faceZoneAverageTemperature();
193 for (
const label patchi : patchIDs_)
195 htcBf[patchi] = q[patchi]/(Tref - Tbf[patchi] + ROOTVSMALL);
212 refRegionName_(
polyMesh::defaultRegion),
232 dict.readIfPresent(
"referenceRegion", refRegionName_);
234 setFaceZoneFaces(
dict);
const Type & lookupObject(const word &name, const bool recursive=false) const
Lookup and return const reference to the object of the given Type. Fatal if not found or the wrong ty...
errorManipArg< error, int > exit(error &err, const int errNo=1)
void resize(const label len)
Adjust allocated size of list.
defineTypeNameAndDebug(faceZoneReferenceTemperature, 0)
A list of keyword definitions, which are a keyword followed by a number of values (eg...
faceZoneReferenceTemperature(const dictionary &dict, const fvMesh &mesh, const word &TName)
Construct from components.
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
GeometricBoundaryField< scalar, fvPatchField, volMesh > Boundary
Type of boundary fields.
virtual bool read(const dictionary &dict)
Read from dictionary.
Macros for easy insertion into run-time selection tables.
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
#define forAll(list, i)
Loop across all elements in list.
A base class for heat transfer coefficient models.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
label size() const noexcept
The number of entries in the list.
const fvMesh & mesh_
Const reference to the mesh.
label whichPatch(const label meshFacei) const
Return patch index for a given mesh face index. Uses binary search.
const fvMesh & mesh() const noexcept
Return const reference to the mesh.
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Mesh data needed to do the Finite Volume discretisation.
Heat transfer coefficient calculation that employs the area-average temperature of a specified face z...
addToRunTimeSelectionTable(heatTransferCoeffModel, faceZoneReferenceTemperature, dictionary)
Mesh consisting of general polyhedral cells.
virtual bool read(const dictionary &dict)
Read from dictionary.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
bool returnReduceOr(const bool value, const label comm=UPstream::worldComm)
Perform logical (or) MPI Allreduce on a copy. Uses UPstream::reduceOr.
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...