36 namespace compressible
41 template<
class sol
idType>
50 mixedFvPatchScalarField(
p, iF,
dict, false),
52 baffleActivated_(
dict.lookupOrDefault<bool>(
"baffleActivated", true)),
57 qrPrevious_(
p.size(), 0.0),
58 qrRelaxation_(
dict.lookupOrDefault<scalar>(
"qrRelaxation", 1)),
59 qrName_(
dict.lookupOrDefault<
word>(
"qr",
"none"))
72 if (
dict.found(
"thickness"))
81 else if (
dict.found(
"Qs"))
86 if (
dict.found(
"qrPrevious"))
91 if (
dict.found(
"refValue") && baffleActivated_)
103 valueFraction() = 0.0;
108 template<
class sol
idType>
118 mixedFvPatchScalarField(ptf,
p, iF, mapper),
120 baffleActivated_(ptf.baffleActivated_),
121 thickness_(mapper(ptf.thickness_)),
122 qs_(mapper(ptf.qs_)),
123 solidDict_(ptf.solidDict_),
124 solidPtr_(ptf.solidPtr_),
125 qrPrevious_(mapper(ptf.qrPrevious_)),
126 qrRelaxation_(ptf.qrRelaxation_),
131 template<
class sol
idType>
139 mixedFvPatchScalarField(ptf, iF),
141 baffleActivated_(ptf.baffleActivated_),
142 thickness_(ptf.thickness_),
144 solidDict_(ptf.solidDict_),
145 solidPtr_(ptf.solidPtr_),
146 qrPrevious_(ptf.qrPrevious_),
147 qrRelaxation_(ptf.qrRelaxation_),
154 template<
class sol
idType>
162 template<
class sol
idType>
163 const thermalBaffle1DFvPatchScalarField<solidType>&
164 thermalBaffle1DFvPatchScalarField<solidType>::nbrField()
const
170 refCast<const fvMesh>(nbrMesh).boundary()[nbrPatchi];
173 refCast<const thermalBaffle1DFvPatchScalarField>
175 nbrPatch.template lookupPatchField<volScalarField, scalar>
183 template<
class sol
idType>
184 const solidType& thermalBaffle1DFvPatchScalarField<solidType>::solid()
const
188 if (solidPtr_.empty())
190 solidPtr_.reset(
new solidType(
"solid", solidDict_));
197 return nbrField().solid();
202 template<
class sol
idType>
203 tmp<scalarField> thermalBaffle1DFvPatchScalarField<solidType>::
204 baffleThickness()
const
208 if (thickness_.size() != patch().size())
211 <<
" Field thickness has not been specified "
212 <<
" for patch " << this->patch().name()
221 return mpp.fromNeighbour(nbrField().baffleThickness());
226 template<
class sol
idType>
227 tmp<scalarField> thermalBaffle1DFvPatchScalarField<solidType>::qs()
const
236 return mpp.fromNeighbour(nbrField().qs());
241 template<
class sol
idType>
248 mixedFvPatchScalarField::map(ptf, mapper);
251 refCast<const thermalBaffle1DFvPatchScalarField>(ptf);
255 mapper(thickness_, tiptf.thickness_);
256 mapper(qs_, tiptf.qs_);
261 template<
class sol
idType>
267 mixedFvPatchScalarField::reset(ptf);
270 refCast<const thermalBaffle1DFvPatchScalarField>(ptf);
274 thickness_.reset(tiptf.thickness_);
275 qs_.reset(tiptf.qs_);
280 template<
class sol
idType>
295 if (baffleActivated_)
298 db().objectRegistry::template lookupObject
305 thermophysicalTransportModel::typeName,
306 internalField().
group()
312 patch().template lookupPatchField<volScalarField, scalar>(TName_);
318 if (qrName_ !=
"none")
320 qr = patch().template lookupPatchField<volScalarField, scalar>
323 qr = qrRelaxation_*qr + (1.0 - qrRelaxation_)*qrPrevious_;
327 const scalarField kappaDelta(kappap*patch().deltaCoeffs());
336 kappas[i] = solid().kappa(0.0, (Tp[i] + nbrTp[i])/2.0);
339 const scalarField kappaByDeltaSolid(kappas/baffleThickness());
345 refValue() = (kappaByDeltaSolid*nbrTp + qs()/2.0)/
alpha;
350 Info<< patch().boundaryMesh().mesh().name() <<
':'
351 << patch().name() <<
':'
352 << this->internalField().name() <<
" <- "
353 << nbrField().patch().name() <<
':'
354 << this->internalField().name() <<
" :"
356 <<
" walltemperature "
357 <<
" min:" <<
gMin(*
this)
358 <<
" max:" <<
gMax(*
this)
367 mixedFvPatchScalarField::updateCoeffs();
371 template<
class sol
idType>
378 writeEntry(os,
"thickness", baffleThickness()());
385 writeEntry(os,
"qrRelaxation", qrRelaxation_);
#define forAll(list, i)
Loop across all elements in list.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static word groupName(Name name, const word &group)
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
static int & msgType()
Message tag of standard messages.
This BC solves a steady 1D thermal baffle.
virtual void write(Ostream &) const
Write.
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchScalarField &, const fvPatchFieldMapper &)
Map the given fvPatchField onto this fvPatchField.
thermalBaffle1DFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Foam::fvPatchFieldMapper.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void operator=(const UList< Type > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Engine which provides mapping between two patches.
const polyPatch & nbrPolyPatch() const
Get the patch to map from.
static const mappedPatchBase & getMap(const polyPatch &patch)
Cast the given polyPatch to a mappedPatchBase. Handle errors.
const polyMesh & nbrMesh() const
Get the mesh for the region to map from.
static void validateMapForField(const PatchFieldType &field, const FieldType &iF, const dictionary &context, const label froms=from::any)
Validate that the map exists and is appropriate for the given.
tmp< Field< Type > > fromNeighbour(const Field< Type > &nbrFld) const
Map/interpolate the neighbour patch field to this patch.
label index() const
Return the index of this patch in the boundaryMesh.
Mesh consisting of general polyhedral cells.
Abstract base class for all fluid and solid thermophysical transport models.
virtual tmp< volScalarField > kappaEff() const =0
Effective thermal turbulent conductivity.
A class for handling words, derived from string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const char *const group
Group name for atomic constants.
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
errorManipArg< error, int > exit(error &err, const int errNo=1)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Type gAverage(const FieldField< Field, Type > &f)
Type gMin(const FieldField< Field, Type > &f)
Type gMax(const FieldField< Field, Type > &f)
static const label differentPatch
static const label sameRegion