volFieldValue.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) 2011-2024 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 "volFieldValue.H"
27 #include "fvMesh.H"
28 #include "volFields.H"
29 #include "writeFile.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace functionObjects
37 {
38 namespace fieldValues
39 {
43 }
44 }
45 }
46 
47 template<>
48 const char*
50 <
52  11
53 >::names[] =
54 {
55  "none",
56  "sum",
57  "sumMag",
58  "average",
59  "volAverage",
60  "volIntegrate",
61  "min",
62  "max",
63  "minMag",
64  "maxMag",
65  "CoV"
66 };
67 
68 const Foam::NamedEnum
69 <
71  11
73 
74 
75 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
76 
78 (
79  const dictionary& dict
80 )
81 {
82  dict.readIfPresent<Switch>("writeLocation", writeLocation_);
83 
84  if (dict.readIfPresent("weightFields", weightFieldNames_))
85  {
87  << " weight fields " << weightFieldNames_;
88  }
89  else if (dict.found("weightField"))
90  {
92  dict.lookup("weightField") >> weightFieldNames_[0];
93 
95  << " weight field " << weightFieldNames_[0];
96  }
97 
98  if (dict.readIfPresent("scaleFactor", scaleFactor_))
99  {
100  Info<< " scale factor = " << scaleFactor_ << nl;
101  }
102 
103  Info<< nl << endl;
104 }
105 
106 
107 template<class Type>
110 {
111  switch (operation_)
112  {
113  case operationType::minMag:
114  case operationType::maxMag:
115  file() << tab << "location" << tab << "cell";
116  if (Pstream::parRun())
117  {
118  file() << tab << "processor";
119  }
120  break;
121  default:
122  break;
123  }
124 }
125 
126 
127 template<>
128 void Foam::functionObjects::fieldValues::volFieldValue::
129 writeFileHeaderLocation<Foam::scalar>()
130 {
131  switch (operation_)
132  {
133  case operationType::min:
134  case operationType::max:
135  case operationType::minMag:
136  case operationType::maxMag:
137  file() << tab << "location" << tab << "cell";
138  if (Pstream::parRun())
139  {
140  file() << tab << "processor";
141  }
142  break;
143  default:
144  break;
145  }
146 }
147 
148 
150 (
151  const label i
152 )
153 {
155 
156  writeCommented(file(), "Time");
157 
158  forAll(fields_, fieldi)
159  {
160  file() << tab << operationTypeNames_[operation_] << "(";
161 
162  forAll(weightFieldNames_, i)
163  {
164  file() << weightFieldNames_[i] << ',';
165  }
166 
167  file() << fields_[fieldi] << ")";
168 
169  if (writeLocation_)
170  {
171  #define writeFileHeaderLocationFieldType(fieldType, none) \
172  if (validField<fieldType>(fields_[fieldi])) \
173  { \
174  writeFileHeaderLocation<fieldType>(); \
175  }
177  #undef writeHeaderLocationFieldType
178  }
179  }
180 
181  file() << endl;
182 }
183 
184 
186 (
187  const Field<scalar>& values,
188  const scalarField& weights,
189  const scalarField& V,
190  Result<scalar>& result
191 ) const
192 {
193  switch (operation_)
194  {
195  case operationType::min:
196  {
197  compareScalars(values, vGreat, result, lessOp<scalar>());
198  return true;
199  }
200  case operationType::minMag:
201  {
202  compareScalars(mag(values), vGreat, result, lessOp<scalar>());
203  return true;
204  }
205  case operationType::max:
206  {
207  compareScalars(values, -vGreat, result, greaterOp<scalar>());
208  return true;
209  }
210  case operationType::maxMag:
211  {
212  compareScalars(mag(values), -vGreat, result, greaterOp<scalar>());
213  return true;
214  }
215  default:
216  {
217  // Fall through to same-type operations
218  return processValuesTypeType(values, weights, V, result);
219  }
220  }
221 }
222 
223 
224 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
225 
227 (
228  const word& name,
229  const Time& runTime,
230  const dictionary& dict
231 )
232 :
233  fieldValue(name, runTime, dict, typeName),
234  fvCellSet(fieldValue::mesh_, dict),
235  writeLocation_(false),
236  operation_(operationTypeNames_.read(dict.lookup("operation"))),
237  scaleFactor_(1)
238 {
239  read(dict);
240 }
241 
242 
244 (
245  const word& name,
246  const objectRegistry& obr,
247  const dictionary& dict
248 )
249 :
250  fieldValue(name, obr, dict, typeName),
251  fvCellSet(fieldValue::mesh_, dict),
252  writeLocation_(false),
253  operation_(operationTypeNames_.read(dict.lookup("operation"))),
254  scaleFactor_(1)
255 {
256  read(dict);
257 }
258 
259 
260 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
261 
263 {}
264 
265 
266 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
267 
269 (
270  const dictionary& dict
271 )
272 {
274 
275  // No additional info to read
276  initialise(dict);
277 
278  return true;
279 }
280 
281 
283 {
284  // Look to see if any fields exist. Use the flag to suppress output later.
285  bool anyFields = false;
286  forAll(fields_, i)
287  {
288  #define validFieldType(fieldType, none) \
289  anyFields = anyFields || validField<fieldType>(fields_[i]);
291  #undef validFieldType
292  }
293  if (!anyFields && fields_.size() > 1) // (error for 1 will happen below)
294  {
295  cannotFindObjects(fields_);
296  }
297 
298  // Initialise the file, write the header, etc...
299  if (anyFields && operation_ != operationType::none)
300  {
302  }
303 
304  // Write the time
305  if (anyFields && operation_ != operationType::none && Pstream::master())
306  {
307  writeTime(file());
308  }
309 
310  // Construct the weight field and the volumes
311  scalarField weights(nCells(), 1);
312  forAll(weightFieldNames_, i)
313  {
314  weights *= getFieldValues<scalar>(weightFieldNames_[i]);
315  }
316  const scalarField V(filterField(fieldValue::mesh_.V()));
317 
318  // Process the fields
319  forAll(fields_, i)
320  {
321  const word& fieldName = fields_[i];
322  bool ok = false;
323 
324  #define writeValuesFieldType(fieldType, none) \
325  ok = ok || writeValues<fieldType>(fieldName, weights, V);
327  #undef writeValuesFieldType
328 
329  if (!ok)
330  {
331  cannotFindObject(fieldName);
332  }
333  }
334 
335  // Finalise
336  if (anyFields && operation_ != operationType::none && Pstream::master())
337  {
338  file() << endl;
339  }
340  Log << endl;
341 
342  return true;
343 }
344 
345 
346 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
void setSize(const label)
Reset size of List.
Definition: List.C:281
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:509
Abstract base-class for Time/database functionObjects.
const word & name() const
Return the name of this functionObject.
Base class for field value -based function objects.
Definition: fieldValue.H:66
virtual bool read(const dictionary &dict)
Read from dictionary.
Definition: fieldValue.C:88
const dictionary & dict() const
Return the reference to the construction dictionary.
Definition: fieldValueI.H:31
virtual bool write()
Write.
Definition: fieldValue.C:112
Provides a 'fvCellSet' specialisation of the fieldValue function object.
scalar scaleFactor_
Scale factor - optional.
static const NamedEnum< operationType, 11 > operationTypeNames_
Operation type names.
Switch writeLocation_
Write the location if available for this operation - optional.
void initialise(const dictionary &dict)
Initialise, e.g. cell addressing.
Definition: volFieldValue.C:78
void writeFileHeaderLocation()
Output file header location information for a given type.
operationType operation_
Operation to apply to values.
volFieldValue(const word &name, const Time &runTime, const dictionary &dict)
Construct from name, Time and dictionary.
virtual void writeFileHeader(const label i)
Output file header information.
virtual bool write()
Calculate and write.
bool processValues(const Field< Type > &values, const scalarField &weights, const scalarField &V, Result< ResultType > &result) const
Apply the operation to the values, and return true if successful.
virtual bool read(const dictionary &)
Read from dictionary.
const fvMesh & mesh_
Reference to the fvMesh.
General run-time selected cell set selection class for fvMesh.
Definition: fvCellSet.H:88
void writeFileHeader(const functionObjects::writeFile &wf, Ostream &file)
Output file header information.
Definition: fvCellSet.C:55
Registry of regIOobjects.
A class for handling words, derived from string.
Definition: word.H:62
#define Log
Report write to Foam::Info if the local log switch is true.
defineTypeNameAndDebug(fieldValueDelta, 0)
addToRunTimeSelectionTable(functionObject, fieldValueDelta, dictionary)
Namespace for OpenFOAM.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
static const char tab
Definition: Ostream.H:265
messageStream Info
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
FOR_ALL_FIELD_TYPES(makeFieldSourceTypedef)
static const char nl
Definition: Ostream.H:266
dictionary dict
#define writeFileHeaderLocationFieldType(fieldType, none)
#define writeValuesFieldType(fieldType, none)
#define validFieldType(fieldType, none)