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 | Copyright (C) 2016-2017 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 "wallHeatFlux.H"
28 #include "solidThermo.H"
29 #include "surfaceInterpolate.H"
30 #include "fvcSnGrad.H"
31 #include "wallPolyPatch.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace functionObjects
39 {
40  defineTypeNameAndDebug(wallHeatFlux, 0);
41  addToRunTimeSelectionTable(functionObject, wallHeatFlux, dictionary);
42 }
43 }
44 
45 
46 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
47 
49 {
50  // Add headers to output data
51  writeHeader(file(), "Wall heat-flux");
52  writeCommented(file(), "Time");
53  writeTabbed(file(), "patch");
54  writeTabbed(file(), "min");
55  writeTabbed(file(), "max");
56  writeTabbed(file(), "integral");
57  file() << endl;
58 }
59 
60 
62 (
63  const volScalarField& alpha,
64  const volScalarField& he,
65  volScalarField& wallHeatFlux
66 )
67 {
68  surfaceScalarField heatFlux
69  (
70  fvc::interpolate(alpha)*fvc::snGrad(he)
71  );
72 
73  volScalarField::Boundary& wallHeatFluxBf =
74  wallHeatFlux.boundaryFieldRef();
75 
76  const surfaceScalarField::Boundary& heatFluxBf =
77  heatFlux.boundaryField();
78 
79  forAll(wallHeatFluxBf, patchi)
80  {
81  wallHeatFluxBf[patchi] = heatFluxBf[patchi];
82  }
83 
84  if (foundObject<volScalarField>("qr"))
85  {
86  const volScalarField& qr = lookupObject<volScalarField>("qr");
87 
88  const volScalarField::Boundary& radHeatFluxBf =
89  qr.boundaryField();
90 
91  forAll(wallHeatFluxBf, patchi)
92  {
93  wallHeatFluxBf[patchi] += radHeatFluxBf[patchi];
94  }
95  }
96 }
97 
98 
99 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
100 
101 Foam::functionObjects::wallHeatFlux::wallHeatFlux
102 (
103  const word& name,
104  const Time& runTime,
105  const dictionary& dict
106 )
107 :
108  fvMeshFunctionObject(name, runTime, dict),
109  logFiles(obr_, name),
110  writeLocalObjects(obr_, log),
111  patchSet_()
112 {
113  volScalarField* wallHeatFluxPtr
114  (
115  new volScalarField
116  (
117  IOobject
118  (
119  type(),
120  mesh_.time().timeName(),
121  mesh_,
124  ),
125  mesh_,
127  )
128  );
129 
130  mesh_.objectRegistry::store(wallHeatFluxPtr);
131 
132  read(dict);
133  resetName(typeName);
134  resetLocalObjectName(typeName);
135 }
136 
137 
138 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
139 
141 {}
142 
143 
144 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
145 
147 {
150 
151  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
152 
153  patchSet_ =
154  mesh_.boundaryMesh().patchSet
155  (
156  wordReList(dict.lookupOrDefault("patches", wordReList()))
157  );
158 
159  Info<< type() << " " << name() << ":" << nl;
160 
161  if (patchSet_.empty())
162  {
163  forAll(pbm, patchi)
164  {
165  if (isA<wallPolyPatch>(pbm[patchi]))
166  {
167  patchSet_.insert(patchi);
168  }
169  }
170 
171  Info<< " processing all wall patches" << nl << endl;
172  }
173  else
174  {
175  Info<< " processing wall patches: " << nl;
176  labelHashSet filteredPatchSet;
177  forAllConstIter(labelHashSet, patchSet_, iter)
178  {
179  label patchi = iter.key();
180  if (isA<wallPolyPatch>(pbm[patchi]))
181  {
182  filteredPatchSet.insert(patchi);
183  Info<< " " << pbm[patchi].name() << endl;
184  }
185  else
186  {
188  << "Requested wall heat-flux on non-wall boundary "
189  << "type patch: " << pbm[patchi].name() << endl;
190  }
191  }
192 
193  Info<< endl;
194 
195  patchSet_ = filteredPatchSet;
196  }
197 
198  return true;
199 }
200 
201 
203 {
204  volScalarField& wallHeatFlux = lookupObjectRef<volScalarField>(type());
205 
206  if
207  (
208  foundObject<compressible::turbulenceModel>
209  (
211  )
212  )
213  {
214  const compressible::turbulenceModel& turbModel =
215  lookupObject<compressible::turbulenceModel>
216  (
218  );
219 
220  calcHeatFlux
221  (
222  turbModel.alphaEff(),
223  turbModel.transport().he(),
224  wallHeatFlux
225  );
226  }
227  else if (foundObject<solidThermo>(solidThermo::dictName))
228  {
229  const solidThermo& thermo =
230  lookupObject<solidThermo>(solidThermo::dictName);
231 
232  calcHeatFlux(thermo.alpha(), thermo.he(), wallHeatFlux);
233  }
234  else
235  {
237  << "Unable to find compressible turbulence model in the "
238  << "database" << exit(FatalError);
239  }
240 
241  return true;
242 }
243 
244 
246 {
247  Log << type() << " " << name() << " write:" << nl;
248 
250 
251  logFiles::write();
252 
253  const volScalarField& wallHeatFlux =
254  obr_.lookupObject<volScalarField>(type());
255 
256  const fvPatchList& patches = mesh_.boundary();
257 
258  const surfaceScalarField::Boundary& magSf =
259  mesh_.magSf().boundaryField();
260 
261  forAllConstIter(labelHashSet, patchSet_, iter)
262  {
263  label patchi = iter.key();
264  const fvPatch& pp = patches[patchi];
265 
266  const scalarField& hfp = wallHeatFlux.boundaryField()[patchi];
267 
268  const scalar minHfp = gMin(hfp);
269  const scalar maxHfp = gMax(hfp);
270  const scalar integralHfp = gSum(magSf[patchi]*hfp);
271 
272  if (Pstream::master())
273  {
274  file()
275  << mesh_.time().value()
276  << token::TAB << pp.name()
277  << token::TAB << minHfp
278  << token::TAB << maxHfp
279  << token::TAB << integralHfp
280  << endl;
281  }
282 
283  Log << " min/max/integ(" << pp.name() << ") = "
284  << minHfp << ", " << maxHfp << ", " << integralHfp << endl;
285  }
286 
287  Log << endl;
288 
289  return true;
290 }
291 
292 
293 // ************************************************************************* //
virtual bool write()
Write function.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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:291
virtual bool execute()
Calculate the wall heat-flux.
Definition: wallHeatFlux.C:202
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual bool write()
Write the wall heat-flux.
Definition: wallHeatFlux.C:245
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual ~wallHeatFlux()
Destructor.
Definition: wallHeatFlux.C:140
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:412
virtual const volScalarField & alpha() const
Thermal diffusivity for enthalpy of mixture [kg/m/s].
Definition: basicThermo.C:501
Calculate the snGrad of the given volField.
patches[0]
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:116
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.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const word & name() const
Return name.
Definition: fvPatch.H:149
Macros for easy insertion into run-time selection tables.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
virtual tmp< volScalarField > alphaEff() const
Return the effective turbulent thermal diffusivity for enthalpy.
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)
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Foam::polyBoundaryMesh.
static const char nl
Definition: Ostream.H:262
Type gMax(const FieldField< Field, Type > &f)
virtual void writeFileHeader(const label i)
File header information.
Definition: wallHeatFlux.C:48
Fundamental solid thermodynamic properties.
Definition: solidThermo.H:54
void calcHeatFlux(const volScalarField &alpha, const volScalarField &he, volScalarField &wallHeatFlux)
Calculate the heat-flux.
Definition: wallHeatFlux.C:62
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensionedScalar pow3(const dimensionedScalar &ds)
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.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
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.
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 bool read(const dictionary &)
Read the wallHeatFlux data.
Definition: wallHeatFlux.C:146
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
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
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
addToRunTimeSelectionTable(functionObject, add, dictionary)
Namespace for OpenFOAM.
functionObject base class for creating, maintaining and writing log files e.g. integrated of averaged...
Definition: logFiles.H:57