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-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 "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),
61  Field<Type>(field),
62  OldTimeField<DimensionedField>(this->time().timeIndex()),
63  mesh_(mesh),
64  dimensions_(dims)
65 {
66  if (field.size() && field.size() != GeoMesh::size(mesh))
67  {
69  << "size of field = " << field.size()
70  << " is not the same as the size of mesh = "
71  << GeoMesh::size(mesh)
72  << abort(FatalError);
73  }
74 }
75 
76 
77 template<class Type, class GeoMesh>
79 (
80  const IOobject& io,
81  const Mesh& mesh,
82  const dimensionSet& dims,
83  const bool checkIOFlags
84 )
85 :
86  regIOobject(io),
87  Field<Type>(GeoMesh::size(mesh)),
88  OldTimeField<DimensionedField>(this->time().timeIndex()),
89  mesh_(mesh),
90  dimensions_(dims)
91 {
92  if (checkIOFlags)
93  {
94  readIfPresent();
95  }
96 }
97 
98 
99 template<class Type, class GeoMesh>
101 (
102  const IOobject& io,
103  const Mesh& mesh,
104  const dimensioned<Type>& dt,
105  const bool checkIOFlags
106 )
107 :
108  regIOobject(io),
109  Field<Type>(GeoMesh::size(mesh), dt.value()),
110  OldTimeField<DimensionedField>(this->time().timeIndex()),
111  mesh_(mesh),
112  dimensions_(dt.dimensions())
113 {
114  if (checkIOFlags)
115  {
116  readIfPresent();
117  }
118 }
119 
120 
121 template<class Type, class GeoMesh>
123 (
125 )
126 :
127  regIOobject(df),
128  Field<Type>(df),
130  mesh_(df.mesh_),
131  dimensions_(df.dimensions_)
132 {}
133 
134 
135 template<class Type, class GeoMesh>
137 (
139  bool reuse
140 )
141 :
142  regIOobject(df, reuse && df.registered()),
143  Field<Type>(df, reuse),
144  OldTimeField<DimensionedField>(this->time().timeIndex()),
145  mesh_(df.mesh_),
146  dimensions_(df.dimensions_)
147 {}
148 
149 
150 template<class Type, class GeoMesh>
152 (
154 )
155 :
156  regIOobject(move(df)),
157  Field<Type>(move(df)),
158  OldTimeField<DimensionedField>(move(df)),
159  mesh_(df.mesh_),
160  dimensions_(move(df.dimensions_))
161 {}
162 
163 
164 template<class Type, class GeoMesh>
166 (
168 )
169 :
170  regIOobject(tdf(), tdf.isTmp()),
171  Field<Type>
172  (
173  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
174  tdf.isTmp()
175  ),
176  OldTimeField<DimensionedField>(this->time().timeIndex()),
177  mesh_(tdf().mesh_),
178  dimensions_(tdf().dimensions_)
179 {
180  tdf.clear();
181 }
182 
183 
184 template<class Type, class GeoMesh>
186 (
187  const IOobject& io,
189  const bool checkIOFlags
190 )
191 :
193  Field<Type>(df),
195  mesh_(df.mesh_),
196  dimensions_(df.dimensions_)
197 {
198  if (!checkIOFlags || !readIfPresent())
199  {
200  copyOldTimes(io, df);
201  }
202 }
203 
204 
205 template<class Type, class GeoMesh>
207 (
208  const IOobject& io,
210  bool reuse,
211  const bool checkIOFlags
212 )
213 :
214  regIOobject(io, df),
215  Field<Type>(df, reuse),
216  OldTimeField<DimensionedField>(this->time().timeIndex()),
217  mesh_(df.mesh_),
218  dimensions_(df.dimensions_)
219 {
220  if (checkIOFlags)
221  {
222  readIfPresent();
223  }
224 }
225 
226 
227 template<class Type, class GeoMesh>
229 (
230  const IOobject& io,
232  const bool checkIOFlags
233 )
234 :
235  regIOobject(io),
236  Field<Type>
237  (
238  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
239  tdf.isTmp()
240  ),
241  OldTimeField<DimensionedField>(this->time().timeIndex()),
242  mesh_(tdf().mesh_),
243  dimensions_(tdf().dimensions_)
244 {
245  tdf.clear();
246 
247  if (checkIOFlags)
248  {
249  readIfPresent();
250  }
251 }
252 
253 
254 template<class Type, class GeoMesh>
256 (
257  const word& newName,
259 )
260 :
261  regIOobject(newName, df, newName != df.name()),
262  Field<Type>(df),
263  OldTimeField<DimensionedField>(this->time().timeIndex()),
264  mesh_(df.mesh_),
265  dimensions_(df.dimensions_)
266 {
267  copyOldTimes(newName, df);
268 }
269 
270 
271 template<class Type, class GeoMesh>
273 (
274  const word& newName,
276  bool reuse
277 )
278 :
279  regIOobject(newName, df, true),
280  Field<Type>(df, reuse),
281  OldTimeField<DimensionedField>(this->time().timeIndex()),
282  mesh_(df.mesh_),
283  dimensions_(df.dimensions_)
284 {}
285 
286 
287 template<class Type, class GeoMesh>
289 (
290  const word& newName,
292 )
293 :
294  regIOobject(newName, tdf(), true),
295  Field<Type>
296  (
297  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
298  tdf.isTmp()
299  ),
300  OldTimeField<DimensionedField>(this->time().timeIndex()),
301  mesh_(tdf().mesh_),
302  dimensions_(tdf().dimensions_)
303 {
304  tdf.clear();
305 }
306 
307 
308 template<class Type, class GeoMesh>
311 {
313  (
315  );
316 }
317 
318 
319 template<class Type, class GeoMesh>
322 (
323  const word& name,
324  const Mesh& mesh,
325  const dimensionSet& ds,
326  const Field<Type>& field
327 )
328 {
329  const bool cacheTmp = mesh.thisDb().cacheTemporaryObject(name);
330 
332  (
334  (
335  IOobject
336  (
337  name,
338  mesh.thisDb().time().name(),
339  mesh.thisDb(),
340  IOobject::NO_READ,
341  IOobject::NO_WRITE,
342  cacheTmp
343  ),
344  mesh,
345  ds,
346  field
347  ),
348  cacheTmp
349  );
350 }
351 
352 
353 template<class Type, class GeoMesh>
356 (
357  const word& name,
358  const Mesh& mesh,
359  const dimensionSet& ds
360 )
361 {
362  const bool cacheTmp = mesh.thisDb().cacheTemporaryObject(name);
363 
365  (
367  (
368  IOobject
369  (
370  name,
371  mesh.thisDb().time().name(),
372  mesh.thisDb(),
373  IOobject::NO_READ,
374  IOobject::NO_WRITE,
375  cacheTmp
376  ),
377  mesh,
378  ds,
379  false
380  ),
381  cacheTmp
382  );
383 }
384 
385 
386 template<class Type, class GeoMesh>
389 (
390  const word& name,
391  const Mesh& mesh,
393 )
394 {
395  const bool cacheTmp = mesh.thisDb().cacheTemporaryObject(name);
396 
398  (
400  (
401  IOobject
402  (
403  name,
404  mesh.thisDb().time().name(),
405  mesh.thisDb(),
406  IOobject::NO_READ,
407  IOobject::NO_WRITE,
408  cacheTmp
409  ),
410  mesh,
411  dt,
412  false
413  ),
414  cacheTmp
415  );
416 }
417 
418 
419 template<class Type, class GeoMesh>
422 (
423  const word& newName,
425 )
426 {
427  const bool cacheTmp = df.db().cacheTemporaryObject(newName);
428 
430  (
432  (
433  IOobject
434  (
435  newName,
436  df.instance(),
437  df.local(),
438  df.db(),
439  IOobject::NO_READ,
440  IOobject::NO_WRITE,
441  cacheTmp
442  ),
443  df
444  ),
445  cacheTmp
446  );
447 }
448 
449 
450 template<class Type, class GeoMesh>
453 (
454  const word& newName,
456 )
457 {
458  const bool cacheTmp = tdf().db().cacheTemporaryObject(newName);
459 
461  (
463  (
464  IOobject
465  (
466  newName,
467  tdf().instance(),
468  tdf().local(),
469  tdf().db(),
470  IOobject::NO_READ,
471  IOobject::NO_WRITE,
472  cacheTmp
473  ),
474  tdf
475  ),
476  cacheTmp
477  );
478 }
479 
480 
481 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
482 
483 template<class Type, class GeoMesh>
485 {
486  db().cacheTemporaryObject(*this);
487 }
488 
489 
490 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
491 
492 template<class Type, class GeoMesh>
494 {
495  this->setUpToDate();
496  storeOldTimes();
497  return *this;
498 }
499 
500 
501 template<class Type, class GeoMesh>
502 tmp
503 <
506 >
508 (
509  const direction d
510 ) const
511 {
513  (
515  (
516  name() + ".component(" + ::Foam::name(d) + ')',
517  mesh_,
518  dimensions_
519  )
520  );
521 
522  Foam::component(result.ref(), *this, d);
523 
524  return result;
525 }
526 
527 
528 template<class Type, class GeoMesh>
529 void DimensionedField<Type, GeoMesh>::replace
530 (
531  const direction d,
532  const DimensionedField
533  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>& df
534 )
535 {
536  Field<Type>::replace(d, df);
537 }
538 
539 
540 template<class Type, class GeoMesh>
541 void DimensionedField<Type, GeoMesh>::replace
542 (
543  const direction d,
544  const tmp
545  <
546  DimensionedField
547  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>
548  >& tdf
549 )
550 {
551  replace(d, tdf());
552  tdf.clear();
553 }
554 
555 
556 template<class Type, class GeoMesh>
557 tmp<DimensionedField<Type, GeoMesh>>
559 {
561  (
563  (
564  name() + ".T()",
565  mesh_,
566  dimensions_
567  )
568  );
569 
570  Foam::T(result.ref(), *this);
571 
572  return result;
573 }
574 
575 
576 template<class Type, class GeoMesh>
578 {
579  dimensioned<Type> Average
580  (
581  this->name() + ".average()",
582  this->dimensions(),
583  gAverage(primitiveField())
584  );
585 
586  return Average;
587 }
588 
589 
590 template<class Type, class GeoMesh>
592 (
593  const DimensionedField<scalar, GeoMesh>& weightField
594 ) const
595 {
596  return
597  (
599  (
600  this->name() + ".weightedAverage(weights)",
601  this->dimensions(),
602  gSum(weightField*primitiveField())/gSum(weightField)
603  )
604  );
605 }
606 
607 
608 template<class Type, class GeoMesh>
610 (
611  const tmp<DimensionedField<scalar, GeoMesh>>& tweightField
612 ) const
613 {
614  dimensioned<Type> wa = weightedAverage(tweightField());
615  tweightField.clear();
616  return wa;
617 }
618 
619 
620 template<class Type, class GeoMesh>
622 (
624 )
625 {
626  // Check for assignment to self
627  if (this == &df)
628  {
630  << "attempted assignment to self"
631  << abort(FatalError);
632  }
633 
634  dimensions_ = df.dimensions();
636 }
637 
638 
639 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
640 
641 template<class Type, class GeoMesh>
643 (
645 )
646 {
647  // Check for assignment to self
648  if (this == &df)
649  {
651  << "attempted assignment to self"
652  << abort(FatalError);
653  }
654 
655  checkField(*this, df, "=");
656 
657  dimensions_ = df.dimensions();
659 }
660 
661 
662 template<class Type, class GeoMesh>
664 (
666 )
667 {
668  // Check for assignment to self
669  if (this == &df)
670  {
672  << "attempted assignment to self"
673  << abort(FatalError);
674  }
675 
676  checkField(*this, df, "=");
677 
678  dimensions_ = move(df.dimensions());
679  Field<Type>::operator=(move(df));
680 }
681 
682 
683 template<class Type, class GeoMesh>
685 (
687 )
688 {
689  const DimensionedField<Type, GeoMesh>& df = tdf();
690 
691  // Check for assignment to self
692  if (this == &df)
693  {
695  << "attempted assignment to self"
696  << abort(FatalError);
697  }
698 
699  checkField(*this, df, "=");
700 
701  dimensions_ = df.dimensions();
702 
703  if (tdf.isTmp())
704  {
705  this->transfer(tdf.ref());
706  }
707  else
708  {
710  }
711 
712  tdf.clear();
713 }
714 
715 
716 template<class Type, class GeoMesh>
718 (
719  const dimensioned<Type>& dt
720 )
721 {
722  dimensions_ = dt.dimensions();
724 }
725 
726 
727 template<class Type, class GeoMesh>
729 {
731 }
732 
733 
734 template<class Type, class GeoMesh>
736 (
738 )
739 {
740  const DimensionedField<Type, GeoMesh>& df = tdf();
741 
742  // Check for assignment to self
743  if (this == &df)
744  {
746  << "attempted assignment to self"
747  << abort(FatalError);
748  }
749 
750  checkField(*this, df, "==");
751 
752  dimensions_ = df.dimensions();
753 
754  if (tdf.isTmp())
755  {
756  this->transfer(tdf.ref());
757  }
758  else
759  {
761  }
762 
763  tdf.clear();
764 }
765 
766 
767 template<class Type, class GeoMesh>
769 (
770  const dimensioned<Type>& dt
771 )
772 {
773  dimensions_ = dt.dimensions();
775 }
776 
777 
778 template<class Type, class GeoMesh>
780 {
782 }
783 
784 
785 #define COMPUTED_ASSIGNMENT(TYPE, op) \
786  \
787 template<class Type, class GeoMesh> \
788 void DimensionedField<Type, GeoMesh>::operator op \
789 ( \
790  const DimensionedField<TYPE, GeoMesh>& df \
791 ) \
792 { \
793  checkField(*this, df, #op); \
794  \
795  dimensions_ op df.dimensions(); \
796  Field<Type>::operator op(df); \
797 } \
798  \
799 template<class Type, class GeoMesh> \
800 void DimensionedField<Type, GeoMesh>::operator op \
801 ( \
802  const tmp<DimensionedField<TYPE, GeoMesh>>& tdf \
803 ) \
804 { \
805  operator op(tdf()); \
806  tdf.clear(); \
807 } \
808  \
809 template<class Type, class GeoMesh> \
810 void DimensionedField<Type, GeoMesh>::operator op \
811 ( \
812  const dimensioned<TYPE>& dt \
813 ) \
814 { \
815  dimensions_ op dt.dimensions(); \
816  Field<Type>::operator op(dt.value()); \
817 }
818 
823 
824 #undef COMPUTED_ASSIGNMENT
825 
826 
827 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
828 
829 #undef checkField
830 
831 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
832 
833 } // End namespace Foam
834 
835 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
836 
837 #include "DimensionedFieldIO.C"
839 
840 // ************************************************************************* //
#define checkField(df1, df2, op)
#define COMPUTED_ASSIGNMENT(TYPE, op)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Field< Type > & primitiveFieldRef()
Return a reference to the internal field.
GeoMesh::Mesh Mesh
Type of mesh on which this DimensionedField is instantiated.
const dimensionSet & dimensions() const
Return dimensions.
const Mesh & mesh() const
Return mesh.
Field< Type >::cmptType cmptType
Component type of the elements of the field.
DimensionedField(const IOobject &, const Mesh &mesh, const dimensionSet &, const Field< Type > &)
Construct from components.
Pre-declare SubField and related Field type.
Definition: Field.H:83
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:47
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const fileName & local() const
Definition: IOobject.H:400
fileName & instance() const
Return the instance directory, constant, system, <time> etc.
Definition: IOobject.C:355
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
Class to add into field types to provide old-time storage and retrieval.
Definition: OldTimeField.H:93
Dimension set for the base types.
Definition: dimensionSet.H:125
Generic dimensioned Type class.
const dimensionSet & dimensions() const
Return const reference to dimensions.
const Type & value() const
Return const reference to value.
bool cacheTemporaryObject(const word &name) const
Return true if given name is in the cacheTemporaryObjects set.
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:55
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:50
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< VolField< Type > > average(const SurfaceField< Type > &ssf)
Area-weighted average a surfaceField creating a volField.
Definition: fvcAverage.C:46
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
Type gSum(const FieldField< Field, Type > &f)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
T clone(const T &t)
Definition: List.H:55
Type gAverage(const FieldField< Field, Type > &f)
error FatalError
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
uint8_t direction
Definition: direction.H:45
label timeIndex
Definition: getTimeIndex.H:4