heatTransferCoeffModel.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017-2020 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
29#include "fvMesh.H"
30#include "fluidThermo.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
40}
41
42
43// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
44
47{
50
51 auto tq = tmp<FieldField<Field, scalar>>::New(Tbf.size());
52 auto& q = tq.ref();
53
54 forAll(q, patchi)
55 {
56 q.set(patchi, new Field<scalar>(Tbf[patchi].size(), Zero));
57 }
58
59 typedef compressible::turbulenceModel cmpTurbModel;
60
61 if (mesh_.foundObject<cmpTurbModel>(cmpTurbModel::propertiesName))
62 {
63 const auto& turb =
64 mesh_.lookupObject<cmpTurbModel>(cmpTurbModel::propertiesName);
65
66 // Note: calling he(p,T) instead of he()
67 const volScalarField he(turb.transport().he(turb.transport().p(), T));
69
70 const volScalarField alphaEff(turb.alphaEff());
72
73 for (const label patchi : patchSet_)
74 {
75 q[patchi] = alphaEffbf[patchi]*hebf[patchi].snGrad();
76 }
77 }
79 {
80 const auto& thermo =
82
83 // Note: calling he(p,T) instead of he()
84 const volScalarField he(thermo.he(thermo.p(), T));
86
87 const volScalarField& alpha(thermo.alpha());
89
90 for (const label patchi : patchSet_)
91 {
92 q[patchi] = alphabf[patchi]*hebf[patchi].snGrad();
93 }
94 }
95 else
96 {
98 << "Unable to find a valid thermo model to evaluate q. " << nl
99 << "Database contents are: " << mesh_.objectRegistry::sortedToc()
100 << exit(FatalError);
101 }
102
103 // Add radiative heat flux contribution if present
104
105 const auto* qrPtr = mesh_.cfindObject<volScalarField>(qrName_);
106
107 if (qrPtr)
108 {
109 const volScalarField::Boundary& qrbf = qrPtr->boundaryField();
110
111 for (const label patchi : patchSet_)
112 {
113 q[patchi] += qrbf[patchi];
114 }
115 }
116
117 return tq;
118}
119
120
121// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
122
124(
125 const dictionary& dict,
126 const fvMesh& mesh,
127 const word& TName
128)
129:
130 mesh_(mesh),
131 patchSet_(),
132 TName_(TName),
133 qrName_("qr")
134{}
135
136
137// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
138
140{
141 patchSet_ = mesh_.boundaryMesh().patchSet(dict.get<wordRes>("patches"));
142
143 dict.readIfPresent("qr", qrName_);
144
145 return true;
146}
147
148
150(
151 volScalarField& result,
153)
154{
155 htc(result, q);
156
157 return true;
158}
159
160
161// ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:52
compressible::turbulenceModel & turb
A field of fields is a PtrList of fields with reference counting.
Definition: FieldField.H:80
Generic templated field type.
Definition: Field.H:82
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual bool read()
Re-read model coefficients if they have changed.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
static const word dictName
Definition: basicThermo.H:256
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:56
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
An abstract base class for heat transfer coeffcient models.
const word TName_
Temperature name.
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.
word qrName_
Name of radiative heat flux.
tmp< FieldField< Field, scalar > > q() const
Return q boundary fields.
labelHashSet patchSet_
Optional list of (wall) patches to process.
virtual bool calc(volScalarField &result, const FieldField< Field, scalar > &q)
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
const Type * cfindObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
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.
Definition: tmp.H:65
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:54
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
const volScalarField & T
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
volScalarField alphaEff("alphaEff", turbulence->nu()/Pr+alphat)
Namespace for OpenFOAM.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
volScalarField & alpha
#define defineRunTimeSelectionTable(baseType, argNames)
Define run-time selection table.
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333