wallHeatTransferCoeff.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2017-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "wallHeatTransferCoeff.H"
28 #include "wallPolyPatch.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace functionObjects
36 {
37  defineTypeNameAndDebug(wallHeatTransferCoeff, 0);
39  (
40  functionObject,
41  wallHeatTransferCoeff,
42  dictionary
43  );
44 }
45 }
46 
47 
48 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
49 
51 (
52  const label i
53 )
54 {
55  // Add headers to output data
56  writeHeader(file(), "Wall heat transfer coefficient");
57  writeCommented(file(), "Time");
58  writeTabbed(file(), "patch");
59  writeTabbed(file(), "min");
60  writeTabbed(file(), "max");
61  writeTabbed(file(), "integral");
62  file() << endl;
63 }
64 
65 
68 (
69  const volScalarField& nu,
70  const volScalarField& nut
71 )
72 {
73  tmp<volScalarField> twallHeatTransferCoeff
74  (
76  (
77  type(),
78  mesh_,
80  (
82  0
83  )
84  )
85  );
86 
87  volScalarField::Boundary& wallHeatTransferCoeffBf =
88  twallHeatTransferCoeff.ref().boundaryFieldRef();
89 
90  const volScalarField::Boundary& nuBf = nu.boundaryField();
91  const volScalarField::Boundary& nutBf = nut.boundaryField();
92 
93  forAll(wallHeatTransferCoeffBf, patchi)
94  {
95  if (!wallHeatTransferCoeffBf[patchi].coupled())
96  {
97  wallHeatTransferCoeffBf[patchi] =
98  rho_*Cp_*(nuBf[patchi]/Prl_ + nutBf[patchi]/Prt_);
99  }
100  }
101 
102  return twallHeatTransferCoeff;
103 }
104 
105 
106 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
107 
109 (
110  const word& name,
111  const Time& runTime,
112  const dictionary& dict
113 )
114 :
115  fvMeshFunctionObject(name, runTime, dict),
116  logFiles(obr_, name),
117  writeLocalObjects(obr_, log),
118  patchSet_()
119 {
120  read(dict);
121  resetName(typeName);
122  resetLocalObjectName(typeName);
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
127 
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
133 
135 {
138 
139  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
140 
141  patchSet_ =
142  mesh_.boundaryMesh().patchSet
143  (
144  wordReList(dict.lookupOrDefault("patches", wordReList()))
145  );
146 
147  Info<< type() << " " << name() << ":" << nl;
148 
149  if (patchSet_.empty())
150  {
151  forAll(pbm, patchi)
152  {
153  if (isA<wallPolyPatch>(pbm[patchi]))
154  {
155  patchSet_.insert(patchi);
156  }
157  }
158 
159  Info<< " processing all wall patches" << nl << endl;
160  }
161  else
162  {
163  Info<< " processing wall patches: " << nl;
164  labelHashSet filteredPatchSet;
165  forAllConstIter(labelHashSet, patchSet_, iter)
166  {
167  label patchi = iter.key();
168  if (isA<wallPolyPatch>(pbm[patchi]))
169  {
170  filteredPatchSet.insert(patchi);
171  Info<< " " << pbm[patchi].name() << endl;
172  }
173  else
174  {
176  << "Requested wall heat-transferCoeff on non-wall boundary "
177  << "type patch: " << pbm[patchi].name() << endl;
178  }
179  }
180 
181  Info<< endl;
182 
183  patchSet_ = filteredPatchSet;
184  }
185 
186  dict.lookup("rho") >> rho_;
187  dict.lookup("Cp") >> Cp_;
188  dict.lookup("Prl") >> Prl_;
189  dict.lookup("Prt") >> Prt_;
190 
191  return true;
192 }
193 
194 
196 {
197  word name(type());
198 
199  if
200  (
201  foundObject<incompressible::momentumTransportModel>
202  (
203  momentumTransportModel::typeName
204  )
205  )
206  {
207  const incompressible::momentumTransportModel& turbModel =
208  lookupObject<incompressible::momentumTransportModel>
209  (
210  momentumTransportModel::typeName
211  );
212 
213  return store
214  (
215  name,
216  calcHeatTransferCoeff(turbModel.nu(), turbModel.nut())
217  );
218  }
219  else
220  {
222  << "Unable to find incompressible turbulence model in the "
223  << "database" << exit(FatalError);
224 
225  return false;
226  }
227 }
228 
229 
231 {
232  Log << type() << " " << name() << " write:" << nl;
233 
235 
236  logFiles::write();
237 
240 
241  const fvPatchList& patches = mesh_.boundary();
242 
243  const surfaceScalarField::Boundary& magSf =
244  mesh_.magSf().boundaryField();
245 
246  forAllConstIter(labelHashSet, patchSet_, iter)
247  {
248  label patchi = iter.key();
249  const fvPatch& pp = patches[patchi];
250 
251  const scalarField& hfp = wallHeatTransferCoeff.boundaryField()[patchi];
252 
253  const scalar minHtcp = gMin(hfp);
254  const scalar maxHtcp = gMax(hfp);
255  const scalar integralHtcp = gSum(magSf[patchi]*hfp);
256 
257  if (Pstream::master())
258  {
259  file()
260  << mesh_.time().value()
261  << tab << pp.name()
262  << tab << minHtcp
263  << tab << maxHtcp
264  << tab << integralHtcp
265  << endl;
266  }
267 
268  Log << " min/max/integ(" << pp.name() << ") = "
269  << minHtcp << ", " << maxHtcp << ", " << integralHtcp << endl;
270  }
271 
272  Log << endl;
273 
274  return true;
275 }
276 
277 
278 // ************************************************************************* //
virtual bool write()
Write function.
Calculates the natural logarithm of the specified scalar field.
Definition: log.H:103
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual bool write()
Write function.
Definition: logFiles.C:180
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
Definition: IOobject.H:303
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const char tab
Definition: Ostream.H:259
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual bool write()
Write the wall heat transfer coefficient.
Type gMin(const FieldField< Field, Type > &f)
wallHeatTransferCoeff(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
virtual bool read(const dictionary &)
Read the wallHeatTransferCoeff data.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
Calculates and write the estimated incompressible flow heat transfer coefficient at wall patches as t...
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
tmp< volScalarField > calcHeatTransferCoeff(const volScalarField &nu, const volScalarField &nut)
Calculate the heat transfer coefficient.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
patches[0]
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:111
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
const ObjectType & lookupObject(const word &fieldName) const
Lookup object from the objectRegistry.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual void writeFileHeader(const label i)
File header information.
Macros for easy insertion into run-time selection tables.
virtual bool read(const dictionary &)
Read optional controls.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
Type gSum(const FieldField< Field, Type > &f)
void writeHeader(std::ostream &, const bool isBinary, const std::string &title)
Write header.
A class for handling words, derived from string.
Definition: word.H:59
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Foam::polyBoundaryMesh.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
static const char nl
Definition: Ostream.H:260
Type gMax(const FieldField< Field, Type > &f)
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
defineTypeNameAndDebug(Qdot, 0)
dimensionedScalar pow3(const dimensionedScalar &ds)
Templated abstract base class for single-phase incompressible turbulence models.
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
FunctionObject base class for managing a list of objects on behalf of the inheriting function object...
virtual bool read(const dictionary &)
Read the list of objects to be written.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual bool execute()
Calculate the wall heat transfer coefficient.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
#define Log
Report write to Foam::Info if the local log switch is true.
messageStream Info
virtual const word & name() const
Return name.
Definition: fvPatch.H:143
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812
functionObject base class for creating, maintaining and writing log files e.g. integrated of averaged...
Definition: logFiles.H:57