43template<
class sol
idType>
52 mixedFvPatchScalarField(
p, iF),
54 baffleActivated_(true),
59 qrPrevious_(
p.size()),
61 qrName_(
"undefined-qr")
65template<
class sol
idType>
76 mixedFvPatchScalarField(ptf,
p, iF, mapper),
78 baffleActivated_(ptf.baffleActivated_),
79 thickness_(ptf.thickness_, mapper),
81 solidDict_(ptf.solidDict_),
82 solidPtr_(ptf.solidPtr_),
83 qrPrevious_(ptf.qrPrevious_, mapper),
84 qrRelaxation_(ptf.qrRelaxation_),
89template<
class sol
idType>
99 mixedFvPatchScalarField(
p, iF),
101 baffleActivated_(
dict.getOrDefault(
"baffleActivated", true)),
106 qrPrevious_(
p.size(),
Zero),
109 dict.getOrDefaultCompat(
"qrRelaxation", {{
"relaxation", 1712}}, 1)
111 qrName_(
dict.getOrDefault<
word>(
"qr",
"none"))
115 if (
dict.found(
"thickness"))
120 if (
dict.found(
"qs"))
125 if (
dict.found(
"qrPrevious"))
130 if (
dict.found(
"refValue") && baffleActivated_)
142 valueFraction() = 0.0;
148template<
class sol
idType>
156 mixedFvPatchScalarField(ptf),
158 baffleActivated_(ptf.baffleActivated_),
159 thickness_(ptf.thickness_),
161 solidDict_(ptf.solidDict_),
162 solidPtr_(ptf.solidPtr_),
163 qrPrevious_(ptf.qrPrevious_),
164 qrRelaxation_(ptf.qrRelaxation_),
169template<
class sol
idType>
178 mixedFvPatchScalarField(ptf, iF),
180 baffleActivated_(ptf.baffleActivated_),
181 thickness_(ptf.thickness_),
183 solidDict_(ptf.solidDict_),
184 solidPtr_(ptf.solidPtr_),
185 qrPrevious_(ptf.qrPrevious_),
186 qrRelaxation_(ptf.qrRelaxation_),
193template<
class sol
idType>
196 const label patchi = patch().index();
198 const label nbrPatchi = samplePolyPatch().index();
200 return (patchi < nbrPatchi);
204template<
class sol
idType>
205const solidType& thermalBaffle1DFvPatchScalarField<solidType>::solid()
const
211 solidPtr_.reset(
new solidType(solidDict_));
217 const fvPatch& nbrPatch =
218 patch().boundaryMesh()[samplePolyPatch().index()];
220 const thermalBaffle1DFvPatchScalarField& nbrField =
221 refCast<const thermalBaffle1DFvPatchScalarField>
223 nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
226 return nbrField.solid();
231template<
class sol
idType>
232tmp<scalarField> thermalBaffle1DFvPatchScalarField<solidType>::
233baffleThickness()
const
237 if (thickness_.size() !=
patch().size())
240 <<
"Field thickness has not been specified"
241 " for patch " << this->
patch().name()
251 const fvPatch& nbrPatch =
252 patch().boundaryMesh()[samplePolyPatch().index()];
253 const thermalBaffle1DFvPatchScalarField& nbrField =
254 refCast<const thermalBaffle1DFvPatchScalarField>
256 nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
259 tmp<scalarField> tthickness
264 mapDist.distribute(thickness);
270template<
class sol
idType>
271tmp<scalarField> thermalBaffle1DFvPatchScalarField<solidType>::qs()
const
281 const fvPatch& nbrPatch =
282 patch().boundaryMesh()[samplePolyPatch().index()];
284 const thermalBaffle1DFvPatchScalarField& nbrField =
285 refCast<const thermalBaffle1DFvPatchScalarField>
287 nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
290 tmp<scalarField> tqs(
new scalarField(nbrField.qs()));
292 mapDist.distribute(qs);
298template<
class sol
idType>
306 mixedFvPatchScalarField::autoMap(m);
310 thickness_.autoMap(m);
316template<
class sol
idType>
323 mixedFvPatchScalarField::rmap(ptf, addr);
326 refCast<const thermalBaffle1DFvPatchScalarField>(ptf);
330 thickness_.rmap(tiptf.thickness_, addr);
331 qs_.rmap(tiptf.qs_, addr);
336template<
class sol
idType>
350 const label patchi = patch().index();
352 const label nbrPatchi = samplePolyPatch().index();
354 if (baffleActivated_)
356 const fvPatch& nbrPatch = patch().boundaryMesh()[nbrPatchi];
359 db().template lookupObject<compressible::turbulenceModel>
368 patch().template lookupPatchField<volScalarField, scalar>(TName_);
373 if (qrName_ !=
"none")
375 qr = patch().template lookupPatchField<volScalarField, scalar>
378 qr = qrRelaxation_*qr + (1.0 - qrRelaxation_)*qrPrevious_;
384 scalarField myKDelta(patch().deltaCoeffs()*kappaw);
388 turbModel.transport().T().boundaryField()[nbrPatchi];
395 kappas[i] = solid().kappa(0.0, (Tp[i] + nbrTp[i])/2.0);
404 refValue() = (KDeltaSolid*nbrTp + qs()/2.0)/
alpha;
408 scalar Q =
gAverage(kappaw*snGrad());
409 Info<< patch().boundaryMesh().mesh().name() <<
':'
410 << patch().name() <<
':'
411 << this->internalField().name() <<
" <- "
412 << nbrPatch.
name() <<
':'
413 << this->internalField().name() <<
" :"
415 <<
" walltemperature "
416 <<
" min:" <<
gMin(*
this)
417 <<
" max:" <<
gMax(*
this)
426 mixedFvPatchScalarField::updateCoeffs();
429template<
class sol
idType>
432 mixedFvPatchScalarField::write(
os);
437 baffleThickness()().writeEntry(
"thickness",
os);
438 qs()().writeEntry(
"qs",
os);
442 qrPrevious_.writeEntry(
"qrPrevious",
os);
443 os.writeEntry(
"qr", qrName_);
444 os.writeEntry(
"qrRelaxation", qrRelaxation_);
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual tmp< volScalarField > kappaEff() const
Return the effective turbulent thermal diffusivity for temperature.
void size(const label n)
Older name for setAddressableSize.
static int & msgType() noexcept
Message tag of standard messages.
This BC solves a steady 1D thermal baffle.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
virtual void operator=(const UList< scalar > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual const word & name() const
Return name.
Class containing processor-to-processor mapping information.
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
const mapDistribute & map() const
Return reference to the parallel distribution map.
static const word propertiesName
Default name of the phase properties dictionary.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
OBJstream os(runTime.globalPath()/outputName)
const std::string patch
OpenFOAM patch number as a std::string.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
Type gMin(const FieldField< Field, Type > &f)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Type gMax(const FieldField< Field, Type > &f)
#define forAll(list, i)
Loop across all elements in list.