Go to the documentation of this file.
42 tabulatedNTUHeatTransfer,
54 { geometryModeType::gmCalculated,
"calculated" },
55 { geometryModeType::gmUser,
"user" },
62 Foam::fv::tabulatedNTUHeatTransfer::ntuTable()
64 if (!ntuTable_.valid())
66 ntuTable_.reset(
new interpolation2DTable<scalar>(
coeffs_));
81 <<
" on mesh " <<
mesh.name()
90 void Foam::fv::tabulatedNTUHeatTransfer::initialiseGeometry()
94 geometryMode_ = geometryModelNames_.get(
"geometryMode", coeffs_);
96 Info<<
"Region " << mesh_.name() <<
" " <<
type() <<
" " << name_
97 <<
" " << geometryModelNames_[geometryMode_] <<
" geometry:" <<
nl;
99 switch (geometryMode_)
103 const fvMesh& nbrMesh =
104 mesh_.time().lookupObject<fvMesh>(nbrRegionName());
106 word inletPatchName(coeffs_.get<word>(
"inletPatch"));
107 word inletPatchNbrName(coeffs_.get<word>(
"inletPatchNbr"));
109 Info<<
" Inlet patch : " << inletPatchName <<
nl
110 <<
" Inlet patch neighbour : " << inletPatchNbrName
113 label patchI = mesh_.boundary().findPatchID(inletPatchName);
115 nbrMesh.boundary().findPatchID(inletPatchNbrName);
117 scalar
alpha(coeffs_.get<scalar>(
"inletBlockageRatio"));
119 if (alpha < 0 || alpha > 1)
122 <<
"Inlet patch blockage ratio must be between 0 and 1"
123 <<
". Current value: " <<
alpha
127 scalar alphaNbr(coeffs_.get<scalar>(
"inletBlockageRatioNbr"));
129 if (alphaNbr < 0 || alphaNbr > 1)
132 <<
"Inlet patch neighbour blockage ratio must be "
133 <<
"between 0 and 1. Current value: " << alphaNbr
138 <<
" Inlet blockage ratio neighbour : " << alphaNbr
143 *
gSum(mesh_.magSf().boundaryField()[patchI]);
146 (scalar(1) - alphaNbr)
147 *
gSum(nbrMesh.magSf().boundaryField()[patchINbr]);
149 scalar
beta(coeffs_.get<scalar>(
"coreBlockageRatio"));
151 if (beta < 0 || beta > 1)
154 <<
"Core volume blockage ratio must be between 0 and 1"
155 <<
". Current value: " <<
beta
159 Info<<
" Core volume blockage ratio : " <<
beta <<
nl;
161 Vcore_ = (scalar(1) -
beta)*meshInterp().V();
167 coeffs_.readEntry(
"Ain", Ain_);
168 coeffs_.readEntry(
"AinNbr", AinNbr_);
170 if (!coeffs_.readIfPresent(
"Vcore", Vcore_))
172 Vcore_ = meshInterp().V();
180 <<
"Unhandled enumeration " << geometryMode_
185 Info<<
" Inlet area local : " << Ain_ <<
nl
186 <<
" Inlet area neighbour : " << AinNbr_ <<
nl
187 <<
" Core volume : " << Vcore_ <<
nl
198 const word& modelType,
204 UName_(coeffs_.getOrDefault<
word>(
"U",
"U")),
205 UNbrName_(coeffs_.getOrDefault<
word>(
"UNbr",
"U")),
206 rhoName_(coeffs_.getOrDefault<
word>(
"rho",
"rho")),
207 rhoNbrName_(coeffs_.getOrDefault<
word>(
"rhoNbr",
"rho")),
209 geometryMode_(gmCalculated),
226 initialiseGeometry();
248 const scalarField mDotNbr(UMagNbrMapped*rhoNbrMapped*AinNbr_);
256 scalar Cpc =
Cp[cellI];
257 scalar CpcNbr = CpNbr[cellI];
258 scalar mDotc = mDot[cellI];
259 scalar mDotcNbr = mDotNbr[cellI];
260 scalar Cmin =
min(Cpc*mDotc, CpcNbr*mDotcNbr);
261 scalar ntu = ntuTable(mDotc, mDotcNbr);
263 htcc[cellI] = Cmin*ntu/Vcore_;
272 coeffs_.readIfPresent(
"U", UName_);
273 coeffs_.readIfPresent(
"UNbr", UNbrName_);
274 coeffs_.readIfPresent(
"rho", rhoName_);
275 coeffs_.readIfPresent(
"rhoNbr", rhoNbrName_);
279 initialiseGeometry();
static const Enum< geometryModeType > geometryModelNames_
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
virtual bool read(const dictionary &dict)
Read dictionary.
A class for handling words, derived from Foam::string.
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure [J/kg/K].
Base class for inter region heat exchange. The derived classes must provide the heat transfer coeffis...
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
psiReactionThermo & thermo
virtual void calculateHtc()
Calculate the heat transfer coefficient.
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)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
2D table interpolation. The data must be in ascending order in both dimensions x and y.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
dictionary coeffs_
Dictionary containing source coefficients.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
const Type & lookupObject(const word &name, const bool recursive=false) const
word dictName() const
The local dictionary name (final part of scoped name)
virtual ~tabulatedNTUHeatTransfer()
Destructor.
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.
Mesh data needed to do the Finite Volume discretisation.
errorManip< error > abort(error &err)
virtual bool read(const dictionary &dict)
Read source dictionary.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
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.
tabulatedNTUHeatTransfer(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from dictionary.
const volScalarField & Cp
const Time & time() const
Return the top-level database.
defineTypeNameAndDebug(atmAmbientTurbSource, 0)
addToRunTimeSelectionTable(option, atmAmbientTurbSource, dictionary)