wallHeatFlux.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) 2016-2017 OpenFOAM Foundation
9  Copyright (C) 2016-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "wallHeatFlux.H"
31 #include "solidThermo.H"
32 #include "surfaceInterpolate.H"
33 #include "fvcSnGrad.H"
34 #include "wallPolyPatch.H"
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
42 namespace functionObjects
43 {
44  defineTypeNameAndDebug(wallHeatFlux, 0);
45  addToRunTimeSelectionTable(functionObject, wallHeatFlux, dictionary);
46 }
47 }
48 
49 
50 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
51 
53 {
54  writeHeader(os, "Wall heat-flux");
55  writeCommented(os, "Time");
56  writeTabbed(os, "patch");
57  writeTabbed(os, "min");
58  writeTabbed(os, "max");
59  writeTabbed(os, "integral");
60  os << endl;
61 }
62 
63 
65 (
66  const volScalarField& alpha,
67  const volScalarField& he,
69 )
70 {
71  volScalarField::Boundary& wallHeatFluxBf = wallHeatFlux.boundaryFieldRef();
72 
73  const volScalarField::Boundary& heBf = he.boundaryField();
74 
75  const volScalarField::Boundary& alphaBf = alpha.boundaryField();
76 
77  for (const label patchi : patchSet_)
78  {
79  wallHeatFluxBf[patchi] = alphaBf[patchi]*heBf[patchi].snGrad();
80  }
81 
82 
83  const auto* qrPtr = cfindObject<volScalarField>(qrName_);
84 
85  if (qrPtr)
86  {
87  const volScalarField::Boundary& radHeatFluxBf = qrPtr->boundaryField();
88 
89  for (const label patchi : patchSet_)
90  {
91  wallHeatFluxBf[patchi] -= radHeatFluxBf[patchi];
92  }
93  }
94 }
95 
96 
97 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
98 
100 (
101  const word& name,
102  const Time& runTime,
103  const dictionary& dict
104 )
105 :
107  writeFile(obr_, name, typeName, dict),
108  patchSet_(),
109  qrName_("qr")
110 {
111  volScalarField* wallHeatFluxPtr
112  (
113  new volScalarField
114  (
115  IOobject
116  (
117  scopedName(typeName),
118  mesh_.time().timeName(),
119  mesh_,
122  ),
123  mesh_,
125  )
126  );
127 
128  mesh_.objectRegistry::store(wallHeatFluxPtr);
129 
130  read(dict);
131 
132  writeFileHeader(file());
133 }
134 
135 
136 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
137 
139 {
142 
143  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
144 
145  patchSet_ =
146  mesh_.boundaryMesh().patchSet
147  (
148  dict.getOrDefault<wordRes>("patches", wordRes())
149  );
150 
151  dict.readIfPresent("qr", qrName_);
152 
153  Info<< type() << " " << name() << ":" << nl;
154 
155  if (patchSet_.empty())
156  {
157  forAll(pbm, patchi)
158  {
159  if (isA<wallPolyPatch>(pbm[patchi]))
160  {
161  patchSet_.insert(patchi);
162  }
163  }
164 
165  Info<< " processing all wall patches" << nl << endl;
166  }
167  else
168  {
169  Info<< " processing wall patches: " << nl;
170  labelHashSet filteredPatchSet;
171  for (const label patchi : patchSet_)
172  {
173  if (isA<wallPolyPatch>(pbm[patchi]))
174  {
175  filteredPatchSet.insert(patchi);
176  Info<< " " << pbm[patchi].name() << endl;
177  }
178  else
179  {
181  << "Requested wall heat-flux on non-wall boundary "
182  << "type patch: " << pbm[patchi].name() << endl;
183  }
184  }
185 
186  Info<< endl;
187 
188  patchSet_ = filteredPatchSet;
189  }
190 
191  return true;
192 }
193 
194 
196 {
197  auto& wallHeatFlux = lookupObjectRef<volScalarField>(scopedName(typeName));
198 
199  if
200  (
201  foundObject<compressible::turbulenceModel>
202  (
204  )
205  )
206  {
207  const compressible::turbulenceModel& turbModel =
208  lookupObject<compressible::turbulenceModel>
209  (
211  );
212 
213  calcHeatFlux
214  (
215  turbModel.alphaEff()(),
216  turbModel.transport().he(),
218  );
219  }
220  else if (foundObject<fluidThermo>(fluidThermo::dictName))
221  {
222  const fluidThermo& thermo =
223  lookupObject<fluidThermo>(fluidThermo::dictName);
224 
225  calcHeatFlux
226  (
227  thermo.alpha(),
228  thermo.he(),
230  );
231  }
232  else if (foundObject<solidThermo>(solidThermo::dictName))
233  {
234  const solidThermo& thermo =
235  lookupObject<solidThermo>(solidThermo::dictName);
236 
237  calcHeatFlux(thermo.alpha(), thermo.he(), wallHeatFlux);
238  }
239  else
240  {
242  << "Unable to find compressible turbulence model in the "
243  << "database" << exit(FatalError);
244  }
245 
246  const fvPatchList& patches = mesh_.boundary();
247 
248  const surfaceScalarField::Boundary& magSf = mesh_.magSf().boundaryField();
249 
250  for (const label patchi : patchSet_)
251  {
252  const fvPatch& pp = patches[patchi];
253 
254  const scalarField& hfp = wallHeatFlux.boundaryField()[patchi];
255 
256  const scalar minHfp = gMin(hfp);
257  const scalar maxHfp = gMax(hfp);
258  const scalar integralHfp = gSum(magSf[patchi]*hfp);
259 
260  if (Pstream::master())
261  {
262  writeCurrentTime(file());
263 
264  file()
265  << token::TAB << pp.name()
266  << token::TAB << minHfp
267  << token::TAB << maxHfp
268  << token::TAB << integralHfp
269  << endl;
270  }
271 
272  Log << " min/max/integ(" << pp.name() << ") = "
273  << minHfp << ", " << maxHfp << ", " << integralHfp << endl;
274 
275  this->setResult("min(" + pp.name() + ")", minHfp);
276  this->setResult("max(" + pp.name() + ")", maxHfp);
277  this->setResult("int(" + pp.name() + ")", integralHfp);
278  }
279 
280 
281  return true;
282 }
283 
284 
286 {
287  const auto& wallHeatFlux =
288  lookupObject<volScalarField>(scopedName(typeName));
289 
290  Log << type() << " " << name() << " write:" << nl
291  << " writing field " << wallHeatFlux.name() << endl;
292 
294 
295  return true;
296 }
297 
298 
299 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::functionObjects::wallHeatFlux
Computes the wall-heat flux at selected wall patches.
Definition: wallHeatFlux.H:160
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Log
#define Log
Definition: PDRblock.C:35
wallHeatFlux.H
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::polyBoundaryMesh
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Definition: polyBoundaryMesh.H:63
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
Foam::functionObjects::wallHeatFlux::calcHeatFlux
void calcHeatFlux(const volScalarField &alpha, const volScalarField &he, volScalarField &wallHeatFlux)
Calculate the heat-flux.
Definition: wallHeatFlux.C:65
Foam::read
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
fvcSnGrad.H
Calculate the snGrad of the given volField.
Foam::turbulenceModel::propertiesName
static const word propertiesName
Default name of the turbulence properties dictionary.
Definition: turbulenceModel.H:100
wallPolyPatch.H
Foam::fluidThermo
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:52
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Foam::functionObjects::wallHeatFlux::writeFileHeader
virtual void writeFileHeader(Ostream &os) const
File header information.
Definition: wallHeatFlux.C:52
Foam::UPstream::master
static bool master(const label communicator=worldComm)
Am I the master process.
Definition: UPstream.H:457
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
Foam::fvPatch::name
virtual const word & name() const
Return name.
Definition: fvPatch.H:167
Foam::HashSet< label, Hash< label > >
Foam::functionObjects::fvMeshFunctionObject
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Definition: fvMeshFunctionObject.H:64
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
solidThermo.H
Foam::functionObjects::wallHeatFlux::execute
virtual bool execute()
Calculate the wall heat-flux.
Definition: wallHeatFlux.C:195
Foam::functionObjects::writeFile::writeHeader
virtual void writeHeader(Ostream &os, const string &str) const
Write a commented header to stream.
Definition: writeFile.C:298
Foam::solidThermo
Fundamental solid thermodynamic properties.
Definition: solidThermo.H:52
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::functionObjects::writeFile::read
virtual bool read(const dictionary &dict)
Read.
Definition: writeFile.C:213
Foam::Field< scalar >
Foam::pow3
dimensionedScalar pow3(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:89
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::functionObjects::regionFunctionObject::read
virtual bool read(const dictionary &dict)
Read optional controls.
Definition: regionFunctionObject.C:173
os
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::functionObjects::wallHeatFlux::read
virtual bool read(const dictionary &dict)
Read the wallHeatFlux data.
Definition: wallHeatFlux.C:138
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::dimMass
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::IOobject::name
const word & name() const noexcept
Return name.
Definition: IOobjectI.H:65
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::functionObject::name
const word & name() const noexcept
Return the name of this functionObject.
Definition: functionObject.C:143
Foam::functionObjects::writeFile::writeCommented
virtual void writeCommented(Ostream &os, const string &str) const
Write a commented string to stream.
Definition: writeFile.C:272
he
volScalarField & he
Definition: YEEqn.H:52
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::ThermalDiffusivity::alphaEff
virtual tmp< volScalarField > alphaEff() const
Return the effective turbulent thermal diffusivity for enthalpy.
Definition: ThermalDiffusivity.H:160
Foam::functionObjects::addToRunTimeSelectionTable
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
Foam::ThermalDiffusivity
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
Definition: phaseCompressibleTurbulenceModelFwd.H:47
Foam::token::TAB
Tab [isspace].
Definition: token.H:123
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
Foam::functionObjects::wallHeatFlux::wallHeatFlux
wallHeatFlux(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
Definition: wallHeatFlux.C:100
Foam::wordRes
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:51
Foam::polyBoundaryMesh::patchSet
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool useGroups=true) const
Return the set of patch IDs corresponding to the given names.
Definition: polyBoundaryMesh.C:864
Foam::HashSet::insert
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
Foam::functionObjects::defineTypeNameAndDebug
defineTypeNameAndDebug(ObukhovLength, 0)
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
Foam::functionObjects::wallHeatFlux::write
virtual bool write()
Write the wall heat-flux.
Definition: wallHeatFlux.C:285
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::functionObjects::writeFile
Base class for writing single files from the function objects.
Definition: writeFile.H:119
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::gMin
Type gMin(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:593
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::functionObjects::writeFile::writeTabbed
virtual void writeTabbed(Ostream &os, const string &str) const
Write a tabbed string to stream.
Definition: writeFile.C:288
Foam::gMax
Type gMax(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:592
turbulentFluidThermoModel.H
Foam::dictionary::dictName
word dictName() const
The local dictionary name (final part of scoped name)
Definition: dictionaryI.H:60