Go to the documentation of this file.
35 #include "surfaceInterpolate.H"
47 effectivenessHeatExchangerSource,
56 void Foam::fv::effectivenessHeatExchangerSource::initialise()
63 <<
type() <<
" " << this->
name() <<
": "
78 label facei = fZone[i];
80 label facePatchId = -1;
90 if (isA<coupledPolyPatch>(pp))
92 if (refCast<const coupledPolyPatch>(pp).owner())
94 faceId = pp.whichFace(facei);
101 else if (!isA<emptyPolyPatch>(pp))
103 faceId = facei - pp.start();
114 if (fZone.flipMap()[i])
138 const word& modelType,
144 secondaryMassFlowRate_(0),
147 userPrimaryInletT_(
false),
148 targetQdotActive_(
false),
150 targetQdotCalcInterval_(5),
151 targetQdotRelax_(0.5),
156 faceZoneName_(
"unknown-faceZone"),
168 fieldNames_.setSize(1,
thermo.he().name());
170 applied_.setSize(1,
false);
197 scalar sumMagPhi = 0;
199 scalar primaryInletTfMean = 0;
202 label facei = faceId_[i];
203 if (facePatchId_[i] != -1)
205 label patchi = facePatchId_[i];
206 scalar phii =
phi.boundaryField()[patchi][facei]*faceSign_[i];
210 scalar Cpfi = Cpf.boundaryField()[patchi][facei];
211 scalar Tfi = Tf.boundaryField()[patchi][facei];
212 scalar magPhii =
mag(phii);
214 sumMagPhi += magPhii;
215 CpfMean += Cpfi*magPhii;
216 primaryInletTfMean += Tfi*magPhii;
220 scalar phii =
phi[facei]*faceSign_[i];
221 scalar magPhii =
mag(phii);
224 sumMagPhi += magPhii;
225 CpfMean += Cpf[facei]*magPhii;
226 primaryInletTfMean += Tf[facei]*magPhii;
232 CpfMean /= sumMagPhi + ROOTVSMALL;
234 scalar primaryInletT = primaryInletT_;
235 if (!userPrimaryInletT_)
238 primaryInletT = primaryInletTfMean/(sumMagPhi + ROOTVSMALL);
242 eTable_()(
mag(sumPhi), secondaryMassFlowRate_)
243 *CpfMean*
mag(sumPhi);
245 const scalar Qt =
alpha*(secondaryInletT_ - primaryInletT);
250 && (mesh_.time().timeIndex() % targetQdotCalcInterval_ == 0)
253 scalar dT = (targetQdot_ - Qt)/(
alpha + ROOTVSMALL);
254 secondaryInletT_ += targetQdotRelax_*dT;
265 deltaTCells[i] =
max(Tref - TCells[i], 0.0);
273 deltaTCells[i] =
max(TCells[i] - Tref, 0.0);
280 scalar sumWeight = 0;
283 label celli = cells_[i];
284 sumWeight += V[celli]*
mag(
U[celli])*deltaTCells[i];
288 if (this->V() > VSMALL &&
mag(Qt) > VSMALL)
294 label celli = cells_[i];
296 Qt*V[celli]*
mag(
U[celli])*deltaTCells[i]
297 /(sumWeight + ROOTVSMALL);
302 <<
indent <<
"Net mass flux [Kg/s] : " << sumPhi <<
nl
303 <<
indent <<
"Total heat exchange [W] : " << Qt <<
nl
304 <<
indent <<
"Secondary inlet T [K] : " << secondaryInletT_ <<
nl
305 <<
indent <<
"Tref [K] : " << Tref <<
nl
306 <<
indent <<
"Effectiveness : "
307 << eTable_()(
mag(sumPhi), secondaryMassFlowRate_) <<
decrIndent
316 UName_ = coeffs_.getOrDefault<
word>(
"U",
"U");
317 TName_ = coeffs_.getOrDefault<
word>(
"T",
"T");
318 phiName_ = coeffs_.getOrDefault<
word>(
"phi",
"phi");
319 coeffs_.readEntry(
"faceZone", faceZoneName_);
321 coeffs_.readEntry(
"secondaryMassFlowRate", secondaryMassFlowRate_);
322 coeffs_.readEntry(
"secondaryInletT", secondaryInletT_);
324 if (coeffs_.readIfPresent(
"primaryInletT", primaryInletT_))
326 userPrimaryInletT_ =
true;
328 <<
"employing user-specified primary flow inlet temperature: "
329 << primaryInletT_ <<
endl;
334 <<
"employing flux-weighted primary flow inlet temperature"
338 if (coeffs_.readIfPresent(
"targetQdot", targetQdot_))
340 targetQdotActive_ =
true;
341 Info<<
indent <<
"employing target heat rejection of "
342 << targetQdot_ <<
nl;
344 coeffs_.readIfPresent
346 "targetQdotCalcInterval",
347 targetQdotCalcInterval_
350 Info<<
indent <<
"updating secondary inlet temperature every "
351 << targetQdotCalcInterval_ <<
" iterations" <<
nl;
353 coeffs_.readIfPresent(
"targetQdotRelax", targetQdotRelax_);
356 << targetQdotRelax_ <<
endl;
A class for handling words, derived from Foam::string.
Intermediate abstract class for handling cell-set options for the derived fvOptions.
static constexpr const zero Zero
Global zero (0)
const word & name() const
Return const access to the source name.
virtual bool read(const dictionary &dict)
Read dictionary.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
bool read(const char *buf, int32_t &val)
Same as readInt32.
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Add explicit/implicit contribution to momentum equation.
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.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
const fvMesh & mesh_
Reference to the mesh database.
#define forAll(list, i)
Loop across all elements in list.
const faceZoneMesh & faceZones() const
Return face zone mesh.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
labelList faceId_
Local list of face IDs.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
virtual bool read(const dictionary &dict)
Read source dictionary.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
labelList faceSign_
List of +1/-1 representing face flip map (1 use as is, -1 negate)
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
word dictName() const
The local dictionary name (final part of scoped name)
word faceZoneName_
Name of the faceZone at the heat exchange inlet.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
const labelIOList & zoneID
Mesh data needed to do the Finite Volume discretisation.
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
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.
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
Ostream & indent(Ostream &os)
Indent stream.
wordList names() const
A list of the zone names.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
labelList facePatchId_
Local list of patch ID per face.
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)
Type gMin(const FieldField< Field, Type > &f)
void setSize(const label newSize)
Alias for resize(const label)
effectivenessHeatExchangerSource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
Type gMax(const FieldField< Field, Type > &f)