fieldAverage.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-2018 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 "fieldAverage.H"
27 #include "volFields.H"
28 #include "fieldAverageItem.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace functionObjects
36 {
37  defineTypeNameAndDebug(fieldAverage, 0);
38  addToRunTimeSelectionTable(functionObject, fieldAverage, dictionary);
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
44 
46 {
47  forAll(faItems_, i)
48  {
49  if (faItems_[i].mean())
50  {
51  if (obr_.found(faItems_[i].meanFieldName()))
52  {
53  obr_.checkOut(*obr_[faItems_[i].meanFieldName()]);
54  }
55  }
56 
57  if (faItems_[i].prime2Mean())
58  {
59  if (obr_.found(faItems_[i].prime2MeanFieldName()))
60  {
61  obr_.checkOut(*obr_[faItems_[i].prime2MeanFieldName()]);
62  }
63  }
64  }
65 }
66 
67 
69 {
70  if (!totalIter_.size())
71  {
72  totalIter_.setSize(faItems_.size(), 1);
73  }
74 
75  if (!totalTime_.size())
76  {
77  totalTime_.setSize(faItems_.size(), obr_.time().deltaTValue());
78  }
79  else
80  {
81  // Check if totalTime_ has been set otherwise initialize
82  forAll(totalTime_, fieldi)
83  {
84  if (totalTime_[fieldi] < 0)
85  {
86  totalTime_[fieldi] = obr_.time().deltaTValue();
87  }
88  }
89  }
90 
91  resetFields();
92 
93  // Add mean fields to the field lists
94  forAll(faItems_, fieldi)
95  {
96  addMeanField<scalar>(fieldi);
97  addMeanField<vector>(fieldi);
98  addMeanField<sphericalTensor>(fieldi);
99  addMeanField<symmTensor>(fieldi);
100  addMeanField<tensor>(fieldi);
101  }
102 
103  // Add prime-squared mean fields to the field lists
104  forAll(faItems_, fieldi)
105  {
106  addPrime2MeanField<scalar, scalar>(fieldi);
107  addPrime2MeanField<vector, symmTensor>(fieldi);
108  }
109 
110  // ensure first averaging works unconditionally
111  prevTimeIndex_ = -1;
112 
113  initialised_ = true;
114 }
115 
116 
118 {
119  Log << " Restarting averaging at time " << obr_.time().timeName()
120  << nl << endl;
121 
122  totalIter_.clear();
123  totalTime_.clear();
124 
125  initialize();
126 }
127 
128 
130 {
131  if (!initialised_)
132  {
133  initialize();
134  }
135 
136  const label currentTimeIndex = obr_.time().timeIndex();
137  const scalar currentTime = obr_.time().value();
138 
139  if (prevTimeIndex_ == currentTimeIndex)
140  {
141  return;
142  }
143  else
144  {
145  prevTimeIndex_ = currentTimeIndex;
146  }
147 
148  if (periodicRestart_ && currentTime > restartPeriod_*periodIndex_)
149  {
150  restart();
151  periodIndex_++;
152  }
153 
154  Log << type() << " " << name() << " write:" << nl
155  << " Calculating averages" << nl;
156 
157  addMeanSqrToPrime2Mean<scalar, scalar>();
158  addMeanSqrToPrime2Mean<vector, symmTensor>();
159 
160  calculateMeanFields<scalar>();
161  calculateMeanFields<vector>();
162  calculateMeanFields<sphericalTensor>();
163  calculateMeanFields<symmTensor>();
164  calculateMeanFields<tensor>();
165 
166  calculatePrime2MeanFields<scalar, scalar>();
167  calculatePrime2MeanFields<vector, symmTensor>();
168 
169  forAll(faItems_, fieldi)
170  {
171  totalIter_[fieldi]++;
172  totalTime_[fieldi] += obr_.time().deltaTValue();
173  }
174 
175  Log << endl;
176 }
177 
178 
180 {
181  Log << " Writing average fields" << endl;
182 
183  writeFields<scalar>();
184  writeFields<vector>();
185  writeFields<sphericalTensor>();
186  writeFields<symmTensor>();
187  writeFields<tensor>();
188 
189  Log << endl;
190 }
191 
192 
194 {
196  (
197  IOobject
198  (
199  name() + "Properties",
200  obr_.time().timeName(),
201  "uniform",
202  obr_,
205  false
206  )
207  );
208 
209  forAll(faItems_, fieldi)
210  {
211  const word& fieldName = faItems_[fieldi].fieldName();
212  propsDict.add(fieldName, dictionary());
213  propsDict.subDict(fieldName).add("totalIter", totalIter_[fieldi]);
214  propsDict.subDict(fieldName).add("totalTime", totalTime_[fieldi]);
215  }
216 
217  propsDict.regIOobject::write();
218 
219  Log << endl;
220 }
221 
222 
224 {
225  if ((restartOnRestart_ || restartOnOutput_) && log)
226  {
227  Info<< " Starting averaging at time " << obr_.time().timeName()
228  << nl;
229  }
230  else
231  {
232  IOobject propsDictHeader
233  (
234  name() + "Properties",
235  obr_.time().timeName(obr_.time().startTime().value()),
236  "uniform",
237  obr_,
240  false
241  );
242 
243  if (!propsDictHeader.typeHeaderOk<IOdictionary>())
244  {
245  Log << " Starting averaging at time "
246  << obr_.time().timeName() << nl;
247 
248  return;
249  }
250 
251  IOdictionary propsDict(propsDictHeader);
252 
253  Log << " Restarting averaging for fields:" << nl;
254 
255  totalIter_.setSize(faItems_.size(), 1);
256 
257  // Initialize totalTime with negative values
258  // to indicate that it has not been set
259  totalTime_.setSize(faItems_.size(), -1);
260 
261  forAll(faItems_, fieldi)
262  {
263  const word& fieldName = faItems_[fieldi].fieldName();
264  if (propsDict.found(fieldName))
265  {
266  dictionary fieldDict(propsDict.subDict(fieldName));
267 
268  totalIter_[fieldi] = readLabel(fieldDict.lookup("totalIter"));
269  totalTime_[fieldi] = readScalar(fieldDict.lookup("totalTime"));
270 
271  Log << " " << fieldName
272  << " iters = " << totalIter_[fieldi]
273  << " time = " << totalTime_[fieldi] << nl;
274  }
275  }
276  }
277 }
278 
279 
280 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
281 
283 (
284  const word& name,
285  const Time& runTime,
286  const dictionary& dict
287 )
288 :
289  fvMeshFunctionObject(name, runTime, dict),
290  prevTimeIndex_(-1),
291  restartOnRestart_(false),
292  restartOnOutput_(false),
293  periodicRestart_(false),
294  restartPeriod_(great),
295  initialised_(false),
296  faItems_(),
297  totalIter_(),
298  totalTime_(),
299  periodIndex_(1)
300 {
301  read(dict);
302 }
303 
304 
305 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
306 
308 {}
309 
310 
311 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
312 
314 {
316 
317  initialised_ = false;
318 
319  Log << type() << " " << name() << ":" << nl;
320 
321  dict.readIfPresent("restartOnRestart", restartOnRestart_);
322  dict.readIfPresent("restartOnOutput", restartOnOutput_);
323  dict.readIfPresent("periodicRestart", periodicRestart_);
324  dict.lookup("fields") >> faItems_;
325 
326  if (periodicRestart_)
327  {
328  dict.lookup("restartPeriod") >> restartPeriod_;
329  }
330 
331  readAveragingProperties();
332 
333  Log << endl;
334 
335  return true;
336 }
337 
338 
340 {
341  calcAverages();
342 
343  return true;
344 }
345 
346 
348 {
349  writeAverages();
350  writeAveragingProperties();
351 
352  if (restartOnOutput_)
353  {
354  restart();
355  }
356 
357  return true;
358 }
359 
360 
361 // ************************************************************************* //
void writeAveragingProperties() const
Write averaging properties - steps and time.
Definition: fieldAverage.C:193
Calculates the natural logarithm of the specified scalar field.
Definition: log.H:57
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
virtual bool read(const dictionary &)
Read the field average data.
Definition: fieldAverage.C:313
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
bool typeHeaderOk(const bool checkType=true)
Read header (uses typeFilePath to find file) and check header.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
virtual bool write()
Write the field averages.
Definition: fieldAverage.C:347
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
virtual bool execute()
Calculate the field averages.
Definition: fieldAverage.C:339
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
IOdictionary propsDict(IOobject("particleTrackProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED))
Macros for easy insertion into run-time selection tables.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:821
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:699
virtual bool read(const dictionary &)
Read optional controls.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
fieldAverage(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
Definition: fieldAverage.C:283
void readAveragingProperties()
Read averaging properties - steps and time.
Definition: fieldAverage.C:223
A class for handling words, derived from string.
Definition: word.H:59
virtual void calcAverages()
Main calculation routine.
Definition: fieldAverage.C:129
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
label readLabel(Istream &is)
Definition: label.H:64
static const char nl
Definition: Ostream.H:265
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
defineTypeNameAndDebug(Qdot, 0)
void initialize()
Reset lists (clear existing values) and initialize averaging.
Definition: fieldAverage.C:68
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
void resetFields()
Checkout fields (causes deletion) from the database.
Definition: fieldAverage.C:45
void restart()
Restart averaging for restartOnOutput.
Definition: fieldAverage.C:117
#define Log
Report write to Foam::Info if the local log switch is true.
messageStream Info
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual ~fieldAverage()
Destructor.
Definition: fieldAverage.C:307
virtual void writeAverages() const
Write averages.
Definition: fieldAverage.C:179
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583