DimensionedField.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-2019 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 "DimensionedField.H"
27 #include "dimensionedType.H"
28 #include "Time.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 // check mesh for two fields
38 #define checkField(df1, df2, op) \
39 if (&(df1).mesh() != &(df2).mesh()) \
40 { \
41  FatalErrorInFunction \
42  << "different mesh for fields " \
43  << (df1).name() << " and " << (df2).name() \
44  << " during operatrion " << op \
45  << abort(FatalError); \
46 }
47 
48 
49 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
50 
51 template<class Type, class GeoMesh>
53 (
54  const IOobject& io,
55  const Mesh& mesh,
56  const dimensionSet& dims,
57  const Field<Type>& field
58 )
59 :
60  regIOobject(io),
62  mesh_(mesh),
63  dimensions_(dims)
64 {
65  if (field.size() && field.size() != GeoMesh::size(mesh))
66  {
68  << "size of field = " << field.size()
69  << " is not the same as the size of mesh = "
70  << GeoMesh::size(mesh)
71  << abort(FatalError);
72  }
73 }
74 
75 
76 template<class Type, class GeoMesh>
78 (
79  const IOobject& io,
80  const Mesh& mesh,
81  const dimensionSet& dims,
82  const bool checkIOFlags
83 )
84 :
85  regIOobject(io),
86  Field<Type>(GeoMesh::size(mesh)),
87  mesh_(mesh),
88  dimensions_(dims)
89 {
90  if (checkIOFlags)
91  {
92  readIfPresent();
93  }
94 }
95 
96 
97 template<class Type, class GeoMesh>
99 (
100  const IOobject& io,
101  const Mesh& mesh,
102  const dimensioned<Type>& dt,
103  const bool checkIOFlags
104 )
105 :
106  regIOobject(io),
107  Field<Type>(GeoMesh::size(mesh), dt.value()),
108  mesh_(mesh),
109  dimensions_(dt.dimensions())
110 {
111  if (checkIOFlags)
112  {
113  readIfPresent();
114  }
115 }
116 
117 
118 template<class Type, class GeoMesh>
120 (
122 )
123 :
124  regIOobject(df),
125  Field<Type>(df),
126  mesh_(df.mesh_),
127  dimensions_(df.dimensions_)
128 {}
129 
130 
131 template<class Type, class GeoMesh>
133 (
135  bool reuse
136 )
137 :
138  regIOobject(df, reuse),
139  Field<Type>(df, reuse),
140  mesh_(df.mesh_),
141  dimensions_(df.dimensions_)
142 {}
143 
144 
145 template<class Type, class GeoMesh>
147 (
149 )
150 :
151  regIOobject(move(df), true),
152  Field<Type>(move(df)),
153  mesh_(df.mesh_),
154  dimensions_(move(df.dimensions_))
155 {}
156 
157 
158 template<class Type, class GeoMesh>
160 (
162 )
163 :
164  regIOobject(tdf(), tdf.isTmp()),
166  (
167  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
168  tdf.isTmp()
169  ),
170  mesh_(tdf().mesh_),
171  dimensions_(tdf().dimensions_)
172 {
173  tdf.clear();
174 }
175 
176 
177 template<class Type, class GeoMesh>
179 (
180  const IOobject& io,
182 )
183 :
184  regIOobject(io),
185  Field<Type>(df),
186  mesh_(df.mesh_),
187  dimensions_(df.dimensions_)
188 {}
189 
190 
191 template<class Type, class GeoMesh>
193 (
194  const IOobject& io,
196  bool reuse
197 )
198 :
199  regIOobject(io, df),
200  Field<Type>(df, reuse),
201  mesh_(df.mesh_),
202  dimensions_(df.dimensions_)
203 {}
204 
205 
206 template<class Type, class GeoMesh>
208 (
209  const word& newName,
211 )
212 :
213  regIOobject(newName, df, newName != df.name()),
214  Field<Type>(df),
215  mesh_(df.mesh_),
216  dimensions_(df.dimensions_)
217 {}
218 
219 
220 template<class Type, class GeoMesh>
222 (
223  const word& newName,
225  bool reuse
226 )
227 :
228  regIOobject(newName, df, true),
229  Field<Type>(df, reuse),
230  mesh_(df.mesh_),
231  dimensions_(df.dimensions_)
232 {}
233 
234 
235 template<class Type, class GeoMesh>
237 (
238  const word& newName,
240 )
241 :
242  regIOobject(newName, tdf(), true),
244  (
245  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
246  tdf.isTmp()
247  ),
248  mesh_(tdf().mesh_),
249  dimensions_(tdf().dimensions_)
250 {
251  tdf.clear();
252 }
253 
254 
255 template<class Type, class GeoMesh>
258 {
260  (
262  );
263 }
264 
265 
266 template<class Type, class GeoMesh>
269 (
270  const word& name,
271  const Mesh& mesh,
272  const dimensionSet& ds
273 )
274 {
276  (
278  (
279  IOobject
280  (
281  name,
282  mesh.time().timeName(),
283  mesh,
286  mesh.cacheTemporaryObject(name)
287  ),
288  mesh,
289  ds,
290  false
291  )
292  );
293 }
294 
295 
296 template<class Type, class GeoMesh>
299 (
300  const word& name,
301  const Mesh& mesh,
302  const dimensioned<Type>& dt
303 )
304 {
306  (
308  (
309  IOobject
310  (
311  name,
312  mesh.time().timeName(),
313  mesh,
316  mesh.cacheTemporaryObject(name)
317  ),
318  mesh,
319  dt,
320  false
321  )
322  );
323 }
324 
325 
326 template<class Type, class GeoMesh>
329 (
330  const word& newName,
332 )
333 {
335  (
337  (
338  IOobject
339  (
340  newName,
341  df.instance(),
342  df.local(),
343  df.db(),
346  df.db().cacheTemporaryObject(newName)
347  ),
348  df
349  )
350  );
351 }
352 
353 
354 template<class Type, class GeoMesh>
357 (
358  const word& newName,
360 )
361 {
363  (
365  (
366  IOobject
367  (
368  newName,
369  tdf().instance(),
370  tdf().local(),
371  tdf().db(),
374  tdf().db().cacheTemporaryObject(newName)
375  ),
376  tdf
377  )
378  );
379 }
380 
381 
382 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
383 
384 template<class Type, class GeoMesh>
386 {
387  db().cacheTemporaryObject(*this);
388 }
389 
390 
391 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
392 
393 template<class Type, class GeoMesh>
394 tmp
395 <
397  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>
398 >
400 (
401  const direction d
402 ) const
403 {
405  (
407  (
408  name() + ".component(" + ::Foam::name(d) + ')',
409  mesh_,
410  dimensions_
411  )
412  );
413 
414  Foam::component(result(), *this, d);
415 
416  return result;
417 }
418 
419 
420 template<class Type, class GeoMesh>
422 (
423  const direction d,
424  const DimensionedField
425  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>& df
426 )
427 {
428  Field<Type>::replace(d, df);
429 }
430 
431 
432 template<class Type, class GeoMesh>
434 (
435  const direction d,
436  const tmp
437  <
439  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>
440  >& tdf
441 )
442 {
443  replace(d, tdf());
444  tdf.clear();
445 }
446 
447 
448 template<class Type, class GeoMesh>
451 {
453  (
455  (
456  name() + ".T()",
457  mesh_,
458  dimensions_
459  )
460  );
461 
462  Foam::T(result(), *this);
463 
464  return result;
465 }
466 
467 
468 template<class Type, class GeoMesh>
470 {
471  dimensioned<Type> Average
472  (
473  this->name() + ".average()",
474  this->dimensions(),
475  gAverage(field())
476  );
477 
478  return Average;
479 }
480 
481 
482 template<class Type, class GeoMesh>
484 (
485  const DimensionedField<scalar, GeoMesh>& weightField
486 ) const
487 {
488  return
489  (
491  (
492  this->name() + ".weightedAverage(weights)",
493  this->dimensions(),
494  gSum(weightField*field())/gSum(weightField)
495  )
496  );
497 }
498 
499 
500 template<class Type, class GeoMesh>
502 (
503  const tmp<DimensionedField<scalar, GeoMesh>>& tweightField
504 ) const
505 {
506  dimensioned<Type> wa = weightedAverage(tweightField());
507  tweightField.clear();
508  return wa;
509 }
510 
511 
512 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
513 
514 template<class Type, class GeoMesh>
515 void DimensionedField<Type, GeoMesh>::operator=
516 (
518 )
519 {
520  // Check for assignment to self
521  if (this == &df)
522  {
524  << "attempted assignment to self"
525  << abort(FatalError);
526  }
527 
528  checkField(*this, df, "=");
529 
530  dimensions_ = df.dimensions();
532 }
533 
534 
535 template<class Type, class GeoMesh>
536 void DimensionedField<Type, GeoMesh>::operator=
537 (
539 )
540 {
541  // Check for assignment to self
542  if (this == &df)
543  {
545  << "attempted assignment to self"
546  << abort(FatalError);
547  }
548 
549  checkField(*this, df, "=");
550 
551  dimensions_ = move(df.dimensions());
552  Field<Type>::operator=(move(df));
553 }
554 
555 
556 template<class Type, class GeoMesh>
557 void DimensionedField<Type, GeoMesh>::operator=
558 (
560 )
561 {
562  const DimensionedField<Type, GeoMesh>& df = tdf();
563 
564  // Check for assignment to self
565  if (this == &df)
566  {
568  << "attempted assignment to self"
569  << abort(FatalError);
570  }
571 
572  checkField(*this, df, "=");
573 
574  dimensions_ = df.dimensions();
575  this->transfer(const_cast<DimensionedField<Type, GeoMesh>&>(df));
576  tdf.clear();
577 }
578 
579 
580 template<class Type, class GeoMesh>
581 void DimensionedField<Type, GeoMesh>::operator=
582 (
583  const dimensioned<Type>& dt
584 )
585 {
586  dimensions_ = dt.dimensions();
588 }
589 
590 
591 template<class Type, class GeoMesh>
593 {
595 }
596 
597 
598 #define COMPUTED_ASSIGNMENT(TYPE, op) \
599  \
600 template<class Type, class GeoMesh> \
601 void DimensionedField<Type, GeoMesh>::operator op \
602 ( \
603  const DimensionedField<TYPE, GeoMesh>& df \
604 ) \
605 { \
606  checkField(*this, df, #op); \
607  \
608  dimensions_ op df.dimensions(); \
609  Field<Type>::operator op(df); \
610 } \
611  \
612 template<class Type, class GeoMesh> \
613 void DimensionedField<Type, GeoMesh>::operator op \
614 ( \
615  const tmp<DimensionedField<TYPE, GeoMesh>>& tdf \
616 ) \
617 { \
618  operator op(tdf()); \
619  tdf.clear(); \
620 } \
621  \
622 template<class Type, class GeoMesh> \
623 void DimensionedField<Type, GeoMesh>::operator op \
624 ( \
625  const dimensioned<TYPE>& dt \
626 ) \
627 { \
628  dimensions_ op dt.dimensions(); \
629  Field<Type>::operator op(dt.value()); \
630 }
631 
636 
637 #undef COMPUTED_ASSIGNMENT
638 
639 
640 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
641 
642 #undef checkField
643 
644 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
645 
646 } // End namespace Foam
647 
648 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
649 
650 #include "DimensionedFieldIO.C"
652 
653 // ************************************************************************* //
bool cacheTemporaryObject(const word &name) const
Return true if given name is in the cacheTemporaryObjects set.
const word & name() const
Return name.
Definition: IOobject.H:303
#define checkField(df1, df2, op)
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &)
Return a temporary field constructed from name, mesh.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
tmp< DimensionedField< Type, GeoMesh > > clone() const
Clone.
uint8_t direction
Definition: direction.H:45
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
tmp< DimensionedField< Type, GeoMesh > > T() const
Return the field transpose (only defined for second rank tensors)
dimensioned< Type > weightedAverage(const DimensionedField< scalar, GeoMesh > &) const
Calculate and return weighted average.
Generic dimensioned Type class.
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
Definition: Field.C:472
const dimensionSet & dimensions() const
Return dimensions.
Dimension set for the base types.
Definition: dimensionSet.H:120
void operator=(const DimensionedField< Type, GeoMesh > &)
dynamicFvMesh & mesh
Type gSum(const FieldField< Field, Type > &f)
Pre-declare SubField and related Field type.
Definition: Field.H:56
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
A class for handling words, derived from string.
Definition: word.H:59
const fileName & local() const
Definition: IOobject.H:400
const Type & value() const
Return const reference to value.
static const zero Zero
Definition: zero.H:97
virtual ~DimensionedField()
Destructor.
Foam::pointMesh ::Mesh Mesh
Type of mesh on which this DimensionedField is instantiated.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define COMPUTED_ASSIGNMENT(TYPE, op)
void operator=(const Field< Type > &)
Definition: Field.C:531
DimensionedField(const IOobject &, const Mesh &mesh, const dimensionSet &, const Field< Type > &)
Construct from components.
void replace(const direction, const DimensionedField< cmptType, GeoMesh > &)
Replace a component field of the field.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
dimensioned< Type > average() const
Calculate and return arithmetic average.
const fileName & instance() const
Definition: IOobject.H:390
Type gAverage(const FieldField< Field, Type > &f)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return const reference to dimensions.
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:55
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:49
tmp< DimensionedField< cmptType, GeoMesh > > component(const direction) const
Return a component field of the field.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
rDeltaTY field()
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:322
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
Namespace for OpenFOAM.