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 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 // check mesh for two fields
37 #define checkField(df1, df2, op) \
38 if (&(df1).mesh() != &(df2).mesh()) \
39 { \
40  FatalErrorInFunction \
41  << "different mesh for fields " \
42  << (df1).name() << " and " << (df2).name() \
43  << " during operatrion " << op \
44  << abort(FatalError); \
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 template<class Type, class GeoMesh>
52 (
53  const IOobject& io,
54  const Mesh& mesh,
55  const dimensionSet& dims,
56  const Field<Type>& field
57 )
58 :
59  regIOobject(io),
60  Field<Type>(field),
61  mesh_(mesh),
62  dimensions_(dims)
63 {
64  if (field.size() && field.size() != GeoMesh::size(mesh))
65  {
67  << "size of field = " << field.size()
68  << " is not the same as the size of mesh = "
69  << GeoMesh::size(mesh)
70  << abort(FatalError);
71  }
72 }
73 
74 
75 template<class Type, class GeoMesh>
77 (
78  const IOobject& io,
79  const Mesh& mesh,
80  const dimensionSet& dims,
81  const bool checkIOFlags
82 )
83 :
84  regIOobject(io),
85  Field<Type>(GeoMesh::size(mesh)),
86  mesh_(mesh),
87  dimensions_(dims)
88 {
89  if (checkIOFlags)
90  {
91  readIfPresent();
92  }
93 }
94 
95 
96 template<class Type, class GeoMesh>
98 (
99  const IOobject& io,
100  const Mesh& mesh,
101  const dimensioned<Type>& dt,
102  const bool checkIOFlags
103 )
104 :
105  regIOobject(io),
106  Field<Type>(GeoMesh::size(mesh), dt.value()),
107  mesh_(mesh),
108  dimensions_(dt.dimensions())
109 {
110  if (checkIOFlags)
111  {
112  readIfPresent();
113  }
114 }
115 
116 
117 template<class Type, class GeoMesh>
119 (
121 )
122 :
123  regIOobject(df),
124  Field<Type>(df),
125  mesh_(df.mesh_),
126  dimensions_(df.dimensions_)
127 {}
128 
129 
130 template<class Type, class GeoMesh>
132 (
134  bool reuse
135 )
136 :
137  regIOobject(df, reuse),
138  Field<Type>(df, reuse),
139  mesh_(df.mesh_),
140  dimensions_(df.dimensions_)
141 {}
142 
143 
144 template<class Type, class GeoMesh>
146 (
148 )
149 :
150  regIOobject(move(df), true),
151  Field<Type>(move(df)),
152  mesh_(df.mesh_),
153  dimensions_(move(df.dimensions_))
154 {}
155 
156 
157 template<class Type, class GeoMesh>
159 (
161 )
162 :
163  regIOobject(tdf(), tdf.isTmp()),
165  (
166  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
167  tdf.isTmp()
168  ),
169  mesh_(tdf().mesh_),
170  dimensions_(tdf().dimensions_)
171 {
172  tdf.clear();
173 }
174 
175 
176 template<class Type, class GeoMesh>
178 (
179  const IOobject& io,
181 )
182 :
183  regIOobject(io),
184  Field<Type>(df),
185  mesh_(df.mesh_),
186  dimensions_(df.dimensions_)
187 {}
188 
189 
190 template<class Type, class GeoMesh>
192 (
193  const IOobject& io,
195  bool reuse
196 )
197 :
198  regIOobject(io, df),
199  Field<Type>(df, reuse),
200  mesh_(df.mesh_),
201  dimensions_(df.dimensions_)
202 {}
203 
204 
205 template<class Type, class GeoMesh>
207 (
208  const word& newName,
210 )
211 :
212  regIOobject(newName, df, newName != df.name()),
213  Field<Type>(df),
214  mesh_(df.mesh_),
215  dimensions_(df.dimensions_)
216 {}
217 
218 
219 template<class Type, class GeoMesh>
221 (
222  const word& newName,
224  bool reuse
225 )
226 :
227  regIOobject(newName, df, true),
228  Field<Type>(df, reuse),
229  mesh_(df.mesh_),
230  dimensions_(df.dimensions_)
231 {}
232 
233 
234 template<class Type, class GeoMesh>
236 (
237  const word& newName,
239 )
240 :
241  regIOobject(newName, tdf(), true),
243  (
244  const_cast<DimensionedField<Type, GeoMesh>&>(tdf()),
245  tdf.isTmp()
246  ),
247  mesh_(tdf().mesh_),
248  dimensions_(tdf().dimensions_)
249 {
250  tdf.clear();
251 }
252 
253 
254 template<class Type, class GeoMesh>
257 {
259  (
261  );
262 }
263 
264 
265 template<class Type, class GeoMesh>
268 (
269  const word& name,
270  const Mesh& mesh,
271  const dimensionSet& ds
272 )
273 {
275  (
277  (
278  IOobject
279  (
280  name,
281  mesh.time().timeName(),
282  mesh,
285  false
286  ),
287  mesh,
288  ds,
289  false
290  )
291  );
292 }
293 
294 
295 template<class Type, class GeoMesh>
298 (
299  const word& name,
300  const Mesh& mesh,
301  const dimensioned<Type>& dt
302 )
303 {
305  (
307  (
308  IOobject
309  (
310  name,
311  mesh.time().timeName(),
312  mesh,
315  false
316  ),
317  mesh,
318  dt,
319  false
320  )
321  );
322 }
323 
324 
325 template<class Type, class GeoMesh>
328 (
329  const word& newName,
331 )
332 {
334  (
336  (
337  IOobject
338  (
339  newName,
340  df.instance(),
341  df.local(),
342  df.db(),
345  false
346  ),
347  df
348  )
349  );
350 }
351 
352 
353 template<class Type, class GeoMesh>
356 (
357  const word& newName,
359 )
360 {
362  (
364  (
365  IOobject
366  (
367  newName,
368  tdf().instance(),
369  tdf().local(),
370  tdf().db(),
373  false
374  ),
375  tdf
376  )
377  );
378 }
379 
380 
381 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
382 
383 template<class Type, class GeoMesh>
385 {}
386 
387 
388 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
389 
390 template<class Type, class GeoMesh>
391 tmp
392 <
394  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>
395 >
397 (
398  const direction d
399 ) const
400 {
402  (
404  (
405  name() + ".component(" + ::Foam::name(d) + ')',
406  mesh_,
407  dimensions_
408  )
409  );
410 
411  Foam::component(result(), *this, d);
412 
413  return result;
414 }
415 
416 
417 template<class Type, class GeoMesh>
419 (
420  const direction d,
421  const DimensionedField
422  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>& df
423 )
424 {
425  Field<Type>::replace(d, df);
426 }
427 
428 
429 template<class Type, class GeoMesh>
431 (
432  const direction d,
433  const tmp
434  <
436  <typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh>
437  >& tdf
438 )
439 {
440  replace(d, tdf());
441  tdf.clear();
442 }
443 
444 
445 template<class Type, class GeoMesh>
448 {
450  (
452  (
453  name() + ".T()",
454  mesh_,
455  dimensions_
456  )
457  );
458 
459  Foam::T(result(), *this);
460 
461  return result;
462 }
463 
464 
465 template<class Type, class GeoMesh>
467 {
468  dimensioned<Type> Average
469  (
470  this->name() + ".average()",
471  this->dimensions(),
472  gAverage(field())
473  );
474 
475  return Average;
476 }
477 
478 
479 template<class Type, class GeoMesh>
481 (
482  const DimensionedField<scalar, GeoMesh>& weightField
483 ) const
484 {
485  return
486  (
488  (
489  this->name() + ".weightedAverage(weights)",
490  this->dimensions(),
491  gSum(weightField*field())/gSum(weightField)
492  )
493  );
494 }
495 
496 
497 template<class Type, class GeoMesh>
499 (
500  const tmp<DimensionedField<scalar, GeoMesh>>& tweightField
501 ) const
502 {
503  dimensioned<Type> wa = weightedAverage(tweightField());
504  tweightField.clear();
505  return wa;
506 }
507 
508 
509 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
510 
511 template<class Type, class GeoMesh>
512 void DimensionedField<Type, GeoMesh>::operator=
513 (
515 )
516 {
517  // Check for assignment to self
518  if (this == &df)
519  {
521  << "attempted assignment to self"
522  << abort(FatalError);
523  }
524 
525  checkField(*this, df, "=");
526 
527  dimensions_ = df.dimensions();
529 }
530 
531 
532 template<class Type, class GeoMesh>
533 void DimensionedField<Type, GeoMesh>::operator=
534 (
536 )
537 {
538  // Check for assignment to self
539  if (this == &df)
540  {
542  << "attempted assignment to self"
543  << abort(FatalError);
544  }
545 
546  checkField(*this, df, "=");
547 
548  dimensions_ = move(df.dimensions());
549  Field<Type>::operator=(move(df));
550 }
551 
552 
553 template<class Type, class GeoMesh>
554 void DimensionedField<Type, GeoMesh>::operator=
555 (
557 )
558 {
559  const DimensionedField<Type, GeoMesh>& df = tdf();
560 
561  // Check for assignment to self
562  if (this == &df)
563  {
565  << "attempted assignment to self"
566  << abort(FatalError);
567  }
568 
569  checkField(*this, df, "=");
570 
571  dimensions_ = df.dimensions();
572  this->transfer(const_cast<DimensionedField<Type, GeoMesh>&>(df));
573  tdf.clear();
574 }
575 
576 
577 template<class Type, class GeoMesh>
578 void DimensionedField<Type, GeoMesh>::operator=
579 (
580  const dimensioned<Type>& dt
581 )
582 {
583  dimensions_ = dt.dimensions();
585 }
586 
587 
588 template<class Type, class GeoMesh>
590 {
592 }
593 
594 
595 #define COMPUTED_ASSIGNMENT(TYPE, op) \
596  \
597 template<class Type, class GeoMesh> \
598 void DimensionedField<Type, GeoMesh>::operator op \
599 ( \
600  const DimensionedField<TYPE, GeoMesh>& df \
601 ) \
602 { \
603  checkField(*this, df, #op); \
604  \
605  dimensions_ op df.dimensions(); \
606  Field<Type>::operator op(df); \
607 } \
608  \
609 template<class Type, class GeoMesh> \
610 void DimensionedField<Type, GeoMesh>::operator op \
611 ( \
612  const tmp<DimensionedField<TYPE, GeoMesh>>& tdf \
613 ) \
614 { \
615  operator op(tdf()); \
616  tdf.clear(); \
617 } \
618  \
619 template<class Type, class GeoMesh> \
620 void DimensionedField<Type, GeoMesh>::operator op \
621 ( \
622  const dimensioned<TYPE>& dt \
623 ) \
624 { \
625  dimensions_ op dt.dimensions(); \
626  Field<Type>::operator op(dt.value()); \
627 }
628 
633 
634 #undef COMPUTED_ASSIGNMENT
635 
636 
637 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
638 
639 #undef checkField
640 
641 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
642 
643 } // End namespace Foam
644 
645 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
646 
647 #include "DimensionedFieldIO.C"
649 
650 // ************************************************************************* //
const word & name() const
Return name.
Definition: IOobject.H:295
#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:163
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:124
A class for handling words, derived from string.
Definition: word.H:59
const fileName & local() const
Definition: IOobject.H:388
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:378
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
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:354
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.