38namespace heatTransferCoeffModels
64 return rho.boundaryField()[patchi];
68 <<
"Unable to set rho for patch " << patchi
78 if (CpName_ ==
"CpInf")
80 const label
n = mesh_.boundary()[patchi].size();
91 return thermo.Cp(pp, Tp, patchi);
95 <<
"Unable to set Cp for patch " << patchi
108 if (mesh_.foundObject<cmpTurbModel>(cmpTurbModel::propertiesName))
111 mesh_.lookupObject<cmpTurbModel>(cmpTurbModel::propertiesName);
113 return turb.devRhoReff()/
turb.rho();
115 else if (mesh_.foundObject<icoTurbModel>(icoTurbModel::propertiesName))
118 mesh_.lookupObject<icoTurbModel>(icoTurbModel::propertiesName);
120 return turb.devReff();
131 else if (mesh_.foundObject<
transportModel>(
"transportProperties"))
133 const auto& laminarT =
140 else if (mesh_.foundObject<
dictionary>(
"transportProperties"))
143 mesh_.lookupObject<
dictionary>(
"transportProperties");
153 <<
"No valid model for viscous stress calculation"
167 auto& Cf = tCf.ref();
177 for (
const label patchi : patchSet_)
187 Cf[patchi] = 2*tauByRhop/
magSqr(URef_);
224 dict.readEntry(
"UInf", URef_);
226 dict.readIfPresent(
"Cp", CpName_);
227 if (CpName_ ==
"CpInf")
229 dict.readEntry(
"CpInf", CpRef_);
232 dict.readIfPresent(
"rho", rhoName_);
233 if (rhoName_ ==
"rhoInf")
235 dict.readEntry(
"rhoInf", rhoRef_);
252 const scalar magU =
mag(URef_);
256 for (
const label patchi : patchSet_)
261 htcBf[patchi] = 0.5*rhop*Cpp*magU*CfBf[patchi];
#define R(A, B, C, D, E, F, K, M)
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
compressible::turbulenceModel & turb
A field of fields is a PtrList of fields with reference counting.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Templated abstract base class for single-phase incompressible turbulence models.
virtual bool read()
Re-read model coefficients if they have changed.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
label size() const noexcept
The number of elements in the list.
static const word dictName
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Fundamental fluid thermodynamic properties.
Mesh data needed to do the Finite Volume discretisation.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
const fvPatch & patch() const
Return patch.
tmp< vectorField > nf() const
Return face normals.
An abstract base class for heat transfer coeffcient models.
static autoPtr< heatTransferCoeffModel > New(const dictionary &dict, const fvMesh &mesh, const word &TName)
Return a reference to the selected heat transfer coefficient model.
const fvMesh & mesh_
Mesh reference.
Heat transfer coefficient calculation based on Reynolds Analogy, which is used to relate turbulent mo...
scalar rhoRef_
Reference density.
virtual tmp< volSymmTensorField > devReff() const
virtual bool read(const dictionary &dict)
Read from dictionary.
word rhoName_
Name of density field.
tmp< FieldField< Field, scalar > > Cf() const
virtual void htc(volScalarField &htc, const FieldField< Field, scalar > &q)
Set the heat transfer coefficient.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
const Type & lookupObject(const word &name, const bool recursive=false) const
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Base-class for all transport models used by the incompressible turbulence models.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const volScalarField & Cp
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
const dimensionSet dimViscosity
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
IOdictionary transportProperties(IOobject("transportProperties", runTime.constant(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
#define forAll(list, i)
Loop across all elements in list.