DimensionedField.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2015-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "DimensionedField.H"
30 #include "dimensionedType.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 // Check mesh for two fields
35 #define checkField(df1, df2, op) \
36 if (&(df1).mesh() != &(df2).mesh()) \
37 { \
38  FatalErrorInFunction \
39  << "different mesh for fields " \
40  << (df1).name() << " and " << (df2).name() \
41  << " during operation " << op \
42  << abort(FatalError); \
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 template<class Type, class GeoMesh>
50 {
51  const label fieldSize = this->size();
52  if (fieldSize)
53  {
54  const label meshSize = GeoMesh::size(this->mesh_);
55  if (fieldSize != meshSize)
56  {
58  << "size of field = " << fieldSize
59  << " is not the same as the size of mesh = "
60  << meshSize
61  << abort(FatalError);
62  }
63  }
64 }
65 
66 
67 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
68 
69 template<class Type, class GeoMesh>
71 (
72  const IOobject& io,
73  const Mesh& mesh,
74  const dimensionSet& dims,
75  const Field<Type>& field
76 )
77 :
78  regIOobject(io),
80  mesh_(mesh),
81  dimensions_(dims),
82  oriented_()
83 {
84  checkFieldSize();
85 }
86 
87 
88 template<class Type, class GeoMesh>
90 (
91  const IOobject& io,
92  const Mesh& mesh,
93  const dimensionSet& dims,
95 )
96 :
97  regIOobject(io),
98  Field<Type>(std::move(field)),
99  mesh_(mesh),
100  dimensions_(dims)
101 {
102  //Info<<"Move construct dimensioned for " << io.name() << nl;
103  checkFieldSize();
104 }
105 
106 
107 template<class Type, class GeoMesh>
109 (
110  const IOobject& io,
111  const Mesh& mesh,
112  const dimensionSet& dims,
113  List<Type>&& field
114 )
115 :
116  regIOobject(io),
117  Field<Type>(std::move(field)),
118  mesh_(mesh),
119  dimensions_(dims)
120 {
121  //Info<<"Move construct dimensioned for " << io.name() << nl;
122  checkFieldSize();
123 }
124 
125 
126 template<class Type, class GeoMesh>
128 (
129  const IOobject& io,
130  const Mesh& mesh,
131  const dimensionSet& dims,
132  const tmp<Field<Type>>& tfield
133 )
134 :
135  regIOobject(io),
136  Field<Type>(tfield.constCast(), tfield.movable()),
137  mesh_(mesh),
138  dimensions_(dims),
139  oriented_()
140 {
141  tfield.clear();
142  checkFieldSize();
143 }
144 
145 
146 template<class Type, class GeoMesh>
148 (
149  const IOobject& io,
150  const Mesh& mesh,
151  const dimensionSet& dims,
152  const bool checkIOFlags
153 )
154 :
155  regIOobject(io),
156  Field<Type>(GeoMesh::size(mesh)),
157  mesh_(mesh),
158  dimensions_(dims),
159  oriented_()
160 {
161  if (checkIOFlags)
162  {
163  readIfPresent();
164  }
165 }
166 
167 
168 template<class Type, class GeoMesh>
170 (
171  const IOobject& io,
172  const Mesh& mesh,
173  const dimensioned<Type>& dt,
174  const bool checkIOFlags
175 )
176 :
177  regIOobject(io),
178  Field<Type>(GeoMesh::size(mesh), dt.value()),
179  mesh_(mesh),
180  dimensions_(dt.dimensions()),
181  oriented_()
182 {
183  if (checkIOFlags)
184  {
185  readIfPresent();
186  }
187 }
188 
189 
190 template<class Type, class GeoMesh>
192 (
194 )
195 :
196  regIOobject(df),
197  Field<Type>(df),
198  mesh_(df.mesh_),
199  dimensions_(df.dimensions_),
200  oriented_(df.oriented_)
201 {}
202 
203 
204 template<class Type, class GeoMesh>
206 (
208 )
209 :
211 {}
212 
213 
214 template<class Type, class GeoMesh>
216 (
218  bool reuse
219 )
220 :
221  regIOobject(df, reuse),
222  Field<Type>(df, reuse),
223  mesh_(df.mesh_),
224  dimensions_(df.dimensions_),
225  oriented_(df.oriented_)
226 {}
227 
228 
229 template<class Type, class GeoMesh>
231 (
233 )
234 :
235  DimensionedField<Type, GeoMesh>(tdf.constCast(), tdf.movable())
236 {
237  tdf.clear();
238 }
239 
240 
241 template<class Type, class GeoMesh>
243 (
244  const IOobject& io,
246 )
247 :
248  regIOobject(io),
249  Field<Type>(df),
250  mesh_(df.mesh_),
251  dimensions_(df.dimensions_),
252  oriented_(df.oriented_)
253 {}
254 
255 
256 template<class Type, class GeoMesh>
258 (
259  const IOobject& io,
261 )
262 :
263  DimensionedField<Type, GeoMesh>(io, df, true)
264 {}
265 
266 
267 template<class Type, class GeoMesh>
269 (
270  const IOobject& io,
272  bool reuse
273 )
274 :
275  regIOobject(io, df),
276  Field<Type>(df, reuse),
277  mesh_(df.mesh_),
278  dimensions_(df.dimensions_),
279  oriented_(df.oriented_)
280 {}
281 
282 
283 template<class Type, class GeoMesh>
285 (
286  const IOobject& io,
288 )
289 :
290  DimensionedField<Type, GeoMesh>(io, tdf.constCast(), tdf.movable())
291 {
292  tdf.clear();
293 }
294 
295 
296 template<class Type, class GeoMesh>
298 (
299  const word& newName,
301 )
302 :
303  regIOobject(newName, df, newName != df.name()),
304  Field<Type>(df),
305  mesh_(df.mesh_),
306  dimensions_(df.dimensions_),
307  oriented_(df.oriented_)
308 {}
309 
310 
311 template<class Type, class GeoMesh>
313 (
314  const word& newName,
316 )
317 :
318  DimensionedField<Type, GeoMesh>(newName, df, true)
319 {}
320 
321 
322 template<class Type, class GeoMesh>
324 (
325  const word& newName,
327  bool reuse
328 )
329 :
330  regIOobject(newName, df, true),
331  Field<Type>(df, reuse),
332  mesh_(df.mesh_),
333  dimensions_(df.dimensions_),
334  oriented_(df.oriented_)
335 {}
336 
337 
338 template<class Type, class GeoMesh>
340 (
341  const word& newName,
343 )
344 :
345  DimensionedField<Type, GeoMesh>(newName, tdf.constCast(), tdf.movable())
346 {
347  tdf.clear();
348 }
349 
350 
351 template<class Type, class GeoMesh>
354 {
356 }
357 
358 
359 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
360 
361 template<class Type, class GeoMesh>
362 Foam::tmp
363 <
365  <
367  >
368 >
370 (
371  const direction d
372 ) const
373 {
375  (
376  IOobject
377  (
378  name() + ".component(" + ::Foam::name(d) + ')',
379  instance(),
380  db()
381  ),
382  mesh_,
383  dimensions_
384  );
385 
386  Foam::component(tresult.ref(), *this, d);
387 
388  return tresult;
389 }
390 
391 
392 template<class Type, class GeoMesh>
394 (
395  const direction d,
396  const DimensionedField
397  <
398  typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh
399  >& df
400 )
401 {
402  Field<Type>::replace(d, df);
403 }
404 
405 
406 template<class Type, class GeoMesh>
408 (
409  const direction d,
410  const tmp
411  <
412  DimensionedField
413  <
414  typename DimensionedField<Type, GeoMesh>::cmptType, GeoMesh
415  >
416  >& tdf
417 )
418 {
419  replace(d, tdf());
420  tdf.clear();
421 }
422 
423 
424 template<class Type, class GeoMesh>
427 {
429  (
430  IOobject
431  (
432  name() + ".T()",
433  instance(),
434  db()
435  ),
436  mesh_,
437  dimensions_
438  );
439 
440  Foam::T(tresult.ref(), *this);
441 
442  return tresult;
443 }
444 
445 
446 template<class Type, class GeoMesh>
448 {
449  return
451  (
452  this->name() + ".average()",
453  this->dimensions(),
454  gAverage(field())
455  );
456 }
457 
458 
459 template<class Type, class GeoMesh>
461 (
462  const DimensionedField<scalar, GeoMesh>& weightField
463 ) const
464 {
465  return
467  (
468  this->name() + ".weightedAverage(weights)",
469  this->dimensions(),
470  gSum(weightField*field())/gSum(weightField)
471  );
472 }
473 
474 
475 template<class Type, class GeoMesh>
477 (
478  const tmp<DimensionedField<scalar, GeoMesh>>& tweightField
479 ) const
480 {
481  dimensioned<Type> result = weightedAverage(tweightField());
482  tweightField.clear();
483  return result;
484 }
485 
486 
487 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
488 
489 template<class Type, class GeoMesh>
491 (
493 )
494 {
495  if (this == &df)
496  {
497  return; // Self-assignment is a no-op
498  }
499 
500  checkField(*this, df, "=");
501 
502  dimensions_ = df.dimensions();
503  oriented_ = df.oriented();
505 }
506 
507 
508 template<class Type, class GeoMesh>
510 (
512 )
513 {
514  auto& df = tdf.constCast();
515 
516  if (this == &df)
517  {
518  return; // Self-assignment is a no-op
519  }
520 
521  checkField(*this, df, "=");
522 
523  dimensions_ = df.dimensions();
524  oriented_ = df.oriented();
525  this->transfer(df);
526  tdf.clear();
527 }
528 
529 
530 template<class Type, class GeoMesh>
532 (
533  const dimensioned<Type>& dt
534 )
535 {
536  dimensions_ = dt.dimensions();
538 }
539 
540 
541 #define COMPUTED_ASSIGNMENT(TYPE, op) \
542  \
543 template<class Type, class GeoMesh> \
544 void Foam::DimensionedField<Type, GeoMesh>::operator op \
545 ( \
546  const DimensionedField<TYPE, GeoMesh>& df \
547 ) \
548 { \
549  checkField(*this, df, #op); \
550  \
551  dimensions_ op df.dimensions(); \
552  oriented_ op df.oriented(); \
553  Field<Type>::operator op(df); \
554 } \
555  \
556 template<class Type, class GeoMesh> \
557 void Foam::DimensionedField<Type, GeoMesh>::operator op \
558 ( \
559  const tmp<DimensionedField<TYPE, GeoMesh>>& tdf \
560 ) \
561 { \
562  operator op(tdf()); \
563  tdf.clear(); \
564 } \
565  \
566 template<class Type, class GeoMesh> \
567 void Foam::DimensionedField<Type, GeoMesh>::operator op \
568 ( \
569  const dimensioned<TYPE>& dt \
570 ) \
571 { \
572  dimensions_ op dt.dimensions(); \
573  Field<Type>::operator op(dt.value()); \
574 }
575 
580 
581 #undef COMPUTED_ASSIGNMENT
582 
583 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
584 
585 #undef checkField
586 
587 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
588 
589 #include "DimensionedFieldIO.C"
591 
592 // ************************************************************************* //
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::component
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
Definition: FieldFieldFunctions.C:44
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
checkField
#define checkField(df1, df2, op)
Definition: DimensionedField.C:35
Foam::gAverage
Type gAverage(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:604
DimensionedField.H
Foam::DimensionedField::DimensionedField
DimensionedField(const IOobject &io, const Mesh &mesh, const dimensionSet &dims, const Field< Type > &field)
Construct from components.
Definition: DimensionedField.C:71
Foam::DimensionedField::replace
void replace(const direction d, const DimensionedField< cmptType, GeoMesh > &df)
Replace a component field of the field.
Foam::DimensionedField::average
dimensioned< Type > average() const
Calculate and return arithmetic average.
Definition: DimensionedField.C:447
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
DimensionedFieldIO.C
Foam::dimensionSet
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Definition: dimensionSet.H:108
Foam::DimensionedField::clone
tmp< DimensionedField< Type, GeoMesh > > clone() const
Clone.
Definition: DimensionedField.C:353
Foam::DimensionedField< Type, Foam::pointMesh >::Mesh
Foam::pointMesh ::Mesh Mesh
Type of mesh on which this DimensionedField is instantiated.
Definition: DimensionedField.H:86
Foam::Field
Generic templated field type.
Definition: Field.H:63
Foam::DimensionedField::cmptType
Field< Type >::cmptType cmptType
Component type of the elements of the field.
Definition: DimensionedField.H:92
Foam::T
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Definition: FieldFieldFunctions.C:58
field
rDeltaTY field()
Foam::FatalError
error FatalError
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::dimensioned
Generic dimensioned Type class.
Definition: dimensionedScalarFwd.H:42
COMPUTED_ASSIGNMENT
#define COMPUTED_ASSIGNMENT(TYPE, op)
Definition: DimensionedField.C:541
Foam::GeoMesh
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:48
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::DimensionedField::oriented
const orientedType & oriented() const noexcept
Return oriented type.
Definition: DimensionedFieldI.H:64
Foam::DimensionedField::dimensions
const dimensionSet & dimensions() const
Return dimensions.
Definition: DimensionedFieldI.H:49
Foam::DimensionedField::T
tmp< DimensionedField< Type, GeoMesh > > T() const
Return the field transpose (only defined for second rank tensors)
Definition: DimensionedField.C:426
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::regIOobject
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:73
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::DimensionedField::weightedAverage
dimensioned< Type > weightedAverage(const DimensionedField< scalar, GeoMesh > &weightField) const
Calculate and return weighted average.
Definition: DimensionedField.C:461
Foam::List< Type >
Foam::direction
uint8_t direction
Definition: direction.H:52
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
DimensionedFieldFunctions.C
Foam::DimensionedField::component
tmp< DimensionedField< cmptType, GeoMesh > > component(const direction d) const
Return a component field of the field.
Foam::dimensioned::dimensions
const dimensionSet & dimensions() const
Return const reference to dimensions.
Definition: dimensionedType.C:420
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
dimensionedType.H