FreeFOAM The Cross-Platform CFD Toolkit
fieldAverage.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) 1991-2010 OpenCFD Ltd.
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 <finiteVolume/volFields.H>
28 #include <OpenFOAM/Time.H>
29 
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
35 
38 
39 
40 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41 
43 {
44  forAll(names, fieldI)
45  {
46  if (names[fieldI].size())
47  {
48  obr_.checkOut(*obr_[names[fieldI]]);
49  }
50  }
51 
52  names.clear();
53  names.setSize(faItems_.size());
54 }
55 
56 
58 {
59  resetFields(meanScalarFields_);
60  resetFields(meanVectorFields_);
61  resetFields(meanSphericalTensorFields_);
62  resetFields(meanSymmTensorFields_);
63  resetFields(meanTensorFields_);
64 
65  resetFields(prime2MeanScalarFields_);
66  resetFields(prime2MeanSymmTensorFields_);
67 
68  totalIter_.clear();
69  totalIter_.setSize(faItems_.size(), 1);
70 
71  totalTime_.clear();
72  totalTime_.setSize(faItems_.size(), obr_.time().deltaT().value());
73 
74 
75  // Add mean fields to the field lists
76  forAll(faItems_, fieldI)
77  {
78  const word& fieldName = faItems_[fieldI].fieldName();
79  if (obr_.foundObject<volScalarField>(fieldName))
80  {
81  addMeanField<scalar>(fieldI, meanScalarFields_);
82  }
83  else if (obr_.foundObject<volVectorField>(fieldName))
84  {
85  addMeanField<vector>(fieldI, meanVectorFields_);
86  }
87  else if (obr_.foundObject<volSphericalTensorField>(fieldName))
88  {
89  addMeanField<sphericalTensor>(fieldI, meanSphericalTensorFields_);
90  }
91  else if (obr_.foundObject<volSymmTensorField>(fieldName))
92  {
93  addMeanField<symmTensor>(fieldI, meanSymmTensorFields_);
94  }
95  else if (obr_.foundObject<volTensorField>(fieldName))
96  {
97  addMeanField<tensor>(fieldI, meanTensorFields_);
98  }
99  else
100  {
101  FatalErrorIn("Foam::fieldAverage::initialize()")
102  << "Requested field " << faItems_[fieldI].fieldName()
103  << " does not exist in the database" << nl
104  << exit(FatalError);
105  }
106  }
107 
108  // Add prime-squared mean fields to the field lists
109  forAll(faItems_, fieldI)
110  {
111  if (faItems_[fieldI].prime2Mean())
112  {
113  const word& fieldName = faItems_[fieldI].fieldName();
114  if (!faItems_[fieldI].mean())
115  {
116  FatalErrorIn("Foam::fieldAverage::initialize()")
117  << "To calculate the prime-squared average, the "
118  << "mean average must also be selected for field "
119  << fieldName << nl << exit(FatalError);
120  }
121 
122  if (obr_.foundObject<volScalarField>(fieldName))
123  {
124  addPrime2MeanField<scalar, scalar>
125  (
126  fieldI,
127  meanScalarFields_,
128  prime2MeanScalarFields_
129  );
130  }
131  else if (obr_.foundObject<volVectorField>(fieldName))
132  {
133  addPrime2MeanField<vector, symmTensor>
134  (
135  fieldI,
136  meanVectorFields_,
137  prime2MeanSymmTensorFields_
138  );
139  }
140  else
141  {
142  FatalErrorIn("Foam::fieldAverage::initialize()")
143  << "prime2Mean average can only be applied to "
144  << "volScalarFields and volVectorFields"
145  << nl << " Field: " << fieldName << nl
146  << exit(FatalError);
147  }
148  }
149  }
150 }
151 
152 
154 {
155  const label currentTimeIndex =
156  static_cast<const fvMesh&>(obr_).time().timeIndex();
157 
158  if (prevTimeIndex_ == currentTimeIndex)
159  {
160  return;
161  }
162  else
163  {
164  prevTimeIndex_ = currentTimeIndex;
165  }
166 
167 
168  Info<< "Calculating averages" << nl << endl;
169  forAll(faItems_, fieldI)
170  {
171  totalIter_[fieldI]++;
172  totalTime_[fieldI] += obr_.time().deltaT().value();
173  }
174 
175  addMeanSqrToPrime2Mean<scalar, scalar>
176  (
177  meanScalarFields_,
178  prime2MeanScalarFields_
179  );
180  addMeanSqrToPrime2Mean<vector, symmTensor>
181  (
182  meanVectorFields_,
183  prime2MeanSymmTensorFields_
184  );
185 
186  calculateMeanFields<scalar>(meanScalarFields_);
187  calculateMeanFields<vector>(meanVectorFields_);
188  calculateMeanFields<sphericalTensor>(meanSphericalTensorFields_);
189  calculateMeanFields<symmTensor>(meanSymmTensorFields_);
190  calculateMeanFields<tensor>(meanTensorFields_);
191 
192  calculatePrime2MeanFields<scalar, scalar>
193  (
194  meanScalarFields_,
195  prime2MeanScalarFields_
196  );
197  calculatePrime2MeanFields<vector, symmTensor>
198  (
199  meanVectorFields_,
200  prime2MeanSymmTensorFields_
201  );
202 }
203 
204 
206 {
207  writeFieldList<scalar>(meanScalarFields_);
208  writeFieldList<vector>(meanVectorFields_);
209  writeFieldList<sphericalTensor>(meanSphericalTensorFields_);
210  writeFieldList<symmTensor>(meanSymmTensorFields_);
211  writeFieldList<tensor>(meanTensorFields_);
212 
213  writeFieldList<scalar>(prime2MeanScalarFields_);
214  writeFieldList<symmTensor>(prime2MeanSymmTensorFields_);
215 }
216 
217 
219 {
221  (
222  IOobject
223  (
224  "fieldAveragingProperties",
225  obr_.time().timeName(),
226  "uniform",
227  obr_,
230  false
231  )
232  );
233 
234  forAll(faItems_, fieldI)
235  {
236  const word& fieldName = faItems_[fieldI].fieldName();
237  propsDict.add(fieldName, dictionary());
238  propsDict.subDict(fieldName).add("totalIter", totalIter_[fieldI]);
239  propsDict.subDict(fieldName).add("totalTime", totalTime_[fieldI]);
240  }
241 
242  propsDict.regIOobject::write();
243 }
244 
245 
247 {
248  if (cleanRestart_)
249  {
250  Info<< "fieldAverage: starting averaging at time "
251  << obr_.time().timeName() << nl << endl;
252  }
253  else
254  {
255  IOobject propsDictHeader
256  (
257  "fieldAveragingProperties",
258  obr_.time().timeName(),
259  "uniform",
260  obr_,
263  false
264  );
265 
266  if (!propsDictHeader.headerOk())
267  {
268  Info<< "fieldAverage: starting averaging at time "
269  << obr_.time().timeName() << nl << endl;
270  return;
271  }
272 
273  IOdictionary propsDict(propsDictHeader);
274 
275  Info<< "fieldAverage: restarting averaging for fields:" << endl;
276  forAll(faItems_, fieldI)
277  {
278  const word& fieldName = faItems_[fieldI].fieldName();
279  if (propsDict.found(fieldName))
280  {
281  dictionary fieldDict(propsDict.subDict(fieldName));
282 
283  totalIter_[fieldI] = readLabel(fieldDict.lookup("totalIter"));
284  totalTime_[fieldI] = readScalar(fieldDict.lookup("totalTime"));
285  Info<< " " << fieldName
286  << " iters = " << totalIter_[fieldI]
287  << " time = " << totalTime_[fieldI] << endl;
288  }
289  }
290  Info<< endl;
291  }
292 }
293 
294 
295 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
296 
298 (
299  const word& name,
300  const objectRegistry& obr,
301  const dictionary& dict,
302  const bool loadFromFiles
303 )
304 :
305  name_(name),
306  obr_(obr),
307  active_(true),
308  prevTimeIndex_(-1),
309  cleanRestart_(false),
310  resetOnOutput_(false),
311  faItems_(),
312  meanScalarFields_(),
313  meanVectorFields_(),
314  meanSphericalTensorFields_(),
315  meanSymmTensorFields_(),
316  meanTensorFields_(),
317  prime2MeanScalarFields_(),
318  prime2MeanSymmTensorFields_(),
319  totalIter_(),
320  totalTime_()
321 {
322  // Only active if a fvMesh is available
323  if (isA<fvMesh>(obr_))
324  {
325  read(dict);
326  }
327  else
328  {
329  active_ = false;
330  WarningIn
331  (
332  "fieldAverage::fieldAverage\n"
333  "(\n"
334  "const word&,\n"
335  "const objectRegistry&,\n"
336  "const dictionary&,\n"
337  "const bool\n"
338  ")"
339  ) << "No fvMesh available, deactivating."
340  << nl << endl;
341  }
342 }
343 
344 
345 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
346 
348 {}
349 
350 
351 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
352 
354 {
355  if (active_)
356  {
357  dict.readIfPresent("cleanRestart", cleanRestart_);
358  dict.readIfPresent("resetOnOutput", resetOnOutput_);
359  dict.lookup("fields") >> faItems_;
360 
361  initialize();
362  readAveragingProperties();
363 
364  // ensure first averaging works unconditionally
365  prevTimeIndex_ = -1;
366  }
367 }
368 
369 
371 {
372  if (active_)
373  {
374  calcAverages();
375  }
376 }
377 
378 
380 {
381 }
382 
383 
385 {
386  if (active_)
387  {
388  calcAverages();
389  writeAverages();
390  writeAveragingProperties();
391 
392  if (resetOnOutput_)
393  {
394  Info<< "fieldAverage: restarting averaging at time "
395  << obr_.time().timeName() << nl << endl;
396 
397  initialize();
398 
399  // ensure first averaging works unconditionally
400  prevTimeIndex_ = -1;
401  }
402  }
403 }
404 
405 
407 {
408  // Do nothing
409 }
410 
411 
413 {
414  // Do nothing
415 }
416 
417 
418 // ************************ vim: set sw=4 sts=4 et: ************************ //