Go to the documentation of this file.
52 mixedFvPatchScalarField(
p, iF),
61 TnbrName_(
"undefined-Tnbr"),
62 qrNbrName_(
"undefined-qrNbr"),
63 qrName_(
"undefined-qr"),
67 thermalInertia_(
false)
69 this->refValue() = 0.0;
70 this->refGrad() = 0.0;
71 this->valueFraction() = 1.0;
84 mixedFvPatchScalarField(psf,
p, iF, mapper),
86 TnbrName_(psf.TnbrName_),
87 qrNbrName_(psf.qrNbrName_),
89 thicknessLayers_(psf.thicknessLayers_),
90 kappaLayers_(psf.kappaLayers_),
91 contactRes_(psf.contactRes_),
92 thermalInertia_(psf.thermalInertia_)
104 mixedFvPatchScalarField(
p, iF),
106 TnbrName_(
dict.getOrDefault<
word>(
"Tnbr",
"T")),
107 qrNbrName_(
dict.getOrDefault<
word>(
"qrNbr",
"none")),
108 qrName_(
dict.getOrDefault<
word>(
"qr",
"none")),
112 thermalInertia_(
dict.getOrDefault<
Switch>(
"thermalInertia",
false))
114 if (!isA<mappedPatchBase>(this->
patch().
patch()))
117 <<
"' not type '" << mappedPatchBase::typeName <<
"'"
118 <<
"\n for patch " <<
p.name()
119 <<
" of field " << internalField().name()
120 <<
" in file " << internalField().objectPath()
124 if (
dict.readIfPresent(
"thicknessLayers", thicknessLayers_))
126 dict.readEntry(
"kappaLayers", kappaLayers_);
128 if (thicknessLayers_.size() > 0)
131 forAll(thicknessLayers_, iLayer)
133 contactRes_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
135 contactRes_ = 1.0/contactRes_;
141 if (
dict.found(
"refValue"))
153 valueFraction() = 1.0;
165 mixedFvPatchScalarField(psf, iF),
167 TnbrName_(psf.TnbrName_),
168 qrNbrName_(psf.qrNbrName_),
169 qrName_(psf.qrName_),
170 thicknessLayers_(psf.thicknessLayers_),
171 kappaLayers_(psf.kappaLayers_),
172 contactRes_(psf.contactRes_),
173 thermalInertia_(psf.thermalInertia_)
194 const label patchi =
patch().index();
196 refCast<const mappedPatchBase>(
patch().
patch());
200 refCast<const fvMesh>(nbrMesh).boundary()[samplePatchi];
216 if (contactRes_ == 0.0)
218 TcNbr.
ref() = nbrField.patchInternalField();
223 TcNbr.
ref() = nbrField;
224 KDeltaNbr.setSize(nbrField.size(), contactRes_);
232 if (qrName_ !=
"none")
238 if (qrNbrName_ !=
"none")
247 const scalar dt =
mesh.time().deltaTValue();
257 thermo->p().boundaryField()[samplePatchi];
259 thermo->T().boundaryField()[samplePatchi];
263 thermo->Cp(ppn, Tpn, samplePatchi)
264 *
thermo->rho(samplePatchi)
272 mCpDtNbr.setSize(Tp.size(),
Zero);
289 thermo->Cp(pp, Tp, patchi)
291 /
patch().deltaCoeffs()/dt
297 mCpDt.setSize(Tp.size(),
Zero);
304 this->internalField().name()
312 scalarField c(KDeltaNbr*TcNbr + (mCpDt + mCpDtNbr)*TpOld);
314 refGrad() = (qr + qrNbr)/
kappa(Tp);
318 valueFraction() = KDeltaNbr/(KDeltaNbr + KDelta);
320 refGrad() = (qr + qrNbr)/
kappa(Tp);
323 mixedFvPatchScalarField::updateCoeffs();
329 Info<<
patch().boundaryMesh().mesh().name() <<
':'
330 <<
patch().name() <<
':'
331 << this->internalField().name() <<
" <- "
332 << nbrMesh.
name() <<
':'
333 << nbrPatch.
name() <<
':'
334 << this->internalField().name() <<
" :"
335 <<
" heat transfer rate:" << Q
336 <<
" walltemperature "
337 <<
" min:" <<
gMin(Tp)
338 <<
" max:" <<
gMax(Tp)
358 os.
writeEntry(
"thermalInertia", thermalInertia_);
360 thicknessLayers_.
writeEntry(
"thicknessLayers", 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.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
const word & name() const
Return name.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Type gAverage(const FieldField< Field, Type > &f)
const polyPatch & samplePolyPatch() const
Get the patch on the region.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Common functions used in temperature coupled boundaries.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Abstract base-class for fluid and solid thermodynamic properties.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
virtual const word & name() const
Return name.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
Mesh consisting of general polyhedral cells.
makePatchTypeField(fvPatchScalarField, alphatFixedDmdtWallBoilingWallFunctionFvPatchScalarField)
#define forAll(list, i)
Loop across all elements in list.
const polyMesh & sampleMesh() const
Get the region mesh.
messageStream Info
Information stream (uses stdout - output is on the master only)
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< scalarField > alpha(const scalarField &Tp) const
Given patch temperature calculate corresponding alphaEff field.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
word dictName() const
The local dictionary name (final part of scoped name)
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
turbulentTemperatureRadCoupledMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual void write(Ostream &) const
Write.
Macros for easy insertion into run-time selection tables.
To & refCast(From &r)
Reference type cast template function.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static int & msgType()
Message tag of standard messages.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const Type * findObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
Mixed boundary condition for temperature and radiation heat transfer to be used for in multiregion ca...
const std::string patch
OpenFOAM patch number as a std::string.
const scalarField & deltaCoeffs() const
virtual tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K 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)
const dimensionedScalar c
Speed of light in a vacuum.
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 > &)
label index() const
The index of this patch in the boundaryMesh.
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...