GeometricField.H
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) 2011-2017 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 Class
25  Foam::GeometricField
26 
27 Description
28  Generic GeometricField class.
29 
30 SourceFiles
31  GeometricFieldI.H
32  GeometricField.C
33  Boundary.C
34  GeometricFieldFunctions.H
35  GeometricFieldFunctions.C
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #ifndef GeometricField_H
40 #define GeometricField_H
41 
42 #include "regIOobject.H"
43 #include "dimensionedTypes.H"
44 #include "DimensionedField.H"
45 #include "FieldField.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 
54 class dictionary;
55 
56 // Forward declaration of friend functions and operators
57 
58 template<class Type, template<class> class PatchField, class GeoMesh>
59 class GeometricField;
60 
61 template<class Type, template<class> class PatchField, class GeoMesh>
62 Ostream& operator<<
63 (
64  Ostream&,
65  const GeometricField<Type, PatchField, GeoMesh>&
66 );
67 
68 template<class Type, template<class> class PatchField, class GeoMesh>
69 Ostream& operator<<
70 (
71  Ostream&,
72  const tmp<GeometricField<Type, PatchField, GeoMesh>>&
73 );
74 
75 
76 /*---------------------------------------------------------------------------*\
77  Class GeometricField Declaration
78 \*---------------------------------------------------------------------------*/
79 
80 template<class Type, template<class> class PatchField, class GeoMesh>
81 class GeometricField
82 :
83  public DimensionedField<Type, GeoMesh>
84 {
85  // Private Member Functions
86 
87  //- Read from file if it is present
88  bool readIfPresent();
89 
90  //- Read old time field from file if it is present
91  bool readOldTimeIfPresent();
92 
93 
94 public:
95 
96  // Public typedefs
97 
98  //- Type of mesh on which this GeometricField is instantiated
99  typedef typename GeoMesh::Mesh Mesh;
100 
101  //- Type of boundary mesh on which this
102  // GeometricField::Boundary is instantiated
103  typedef typename GeoMesh::BoundaryMesh BoundaryMesh;
104 
105  //- Type of the internal field from which this GeometricField is derived
107 
108  //- Type of the patch field of which the
109  // GeometricField::Boundary is composed
110  typedef PatchField<Type> Patch;
111 
113  class Boundary
114  :
115  public FieldField<PatchField, Type>
116  {
117  // Private data
118 
119  //- Reference to BoundaryMesh for which this field is defined
120  const BoundaryMesh& bmesh_;
121 
122 
123  public:
124 
125  // Constructors
126 
127  //- Construct from a BoundaryMesh
128  Boundary(const BoundaryMesh&);
129 
130  //- Construct from a BoundaryMesh,
131  // reference to the internal field
132  // and a patch type
133  Boundary
134  (
135  const BoundaryMesh&,
136  const Internal&,
137  const word&
138  );
139 
140  //- Construct from a BoundaryMesh,
141  // reference to the internal field
142  // and a wordList of patch types and optional the actual patch
143  // types (to override constraint patches)
144  Boundary
145  (
146  const BoundaryMesh&,
147  const Internal&,
148  const wordList& wantedPatchTypes,
149  const wordList& actualPatchTypes = wordList()
150  );
151 
152  //- Construct from a BoundaryMesh,
153  // reference to the internal field
154  // and a PtrList<PatchField<Type>>
155  Boundary
156  (
157  const BoundaryMesh&,
158  const Internal&,
159  const PtrList<PatchField<Type>>&
160  );
161 
162  //- Construct as copy setting the reference to the internal field
163  Boundary
164  (
165  const Internal&,
166  const Boundary&
167  );
168 
169  //- Construct as copy
170  // Dangerous because Field may be set to a field which gets deleted
171  // Need new type of BoundaryField, one which is part of a geometric
172  // field for which snGrad etc. may be called and a free standing
173  // BoundaryField for which such operations are unavailable.
174  Boundary
175  (
176  const Boundary&
177  );
178 
179  //- Construct from dictionary
180  Boundary
181  (
182  const BoundaryMesh&,
183  const Internal&,
184  const dictionary&
185  );
186 
187 
188  // Member functions
189 
190  //- Read the boundary field
191  void readField
192  (
193  const Internal& field,
194  const dictionary& dict
195  );
196 
197  //- Update the boundary condition coefficients
198  void updateCoeffs();
199 
200  //- Evaluate boundary conditions
201  void evaluate();
202 
203  //- Return a list of the patch types
204  wordList types() const;
205 
206  //- Return BoundaryField of the cell values neighbouring
207  // the boundary
209 
210  //- Return a list of pointers for each patch field with only those
211  // pointing to interfaces being set
213 
214  //- Return a list of pointers for each patch field with only those
215  // pointing to interfaces being set
217 
218  //- Write boundary field as dictionary entry
219  void writeEntry(const word& keyword, Ostream& os) const;
220 
221 
222  // Member operators
223 
224  //- Assignment to BoundaryField<Type, PatchField, BoundaryMesh>
225  void operator=(const Boundary&);
226 
227  //- Assignment to FieldField<PatchField, Type>
229 
230  //- Assignment to Type
231  void operator=(const Type&);
232 
233 
234  //- Forced assignment to
235  // BoundaryField<Type, PatchField, BoundaryMesh>
236  void operator==(const Boundary&);
237 
238  //- Forced assignment to FieldField<PatchField, Type>
240 
241  //- Forced assignment to Type
242  void operator==(const Type&);
243  };
244 
245 
246 private:
247 
248  // Private data
249 
250  //- Current time index.
251  // Used to trigger the storing of the old-time value
252  mutable label timeIndex_;
253 
254  //- Pointer to old time field
255  mutable GeometricField<Type, PatchField, GeoMesh>* field0Ptr_;
256 
257  //- Pointer to previous iteration (used for under-relaxation)
258  mutable GeometricField<Type, PatchField, GeoMesh>* fieldPrevIterPtr_;
259 
260  //- Boundary Type field containing boundary field values
261  Boundary boundaryField_;
262 
263 
264  // Private Member Functions
265 
266  //- Read the field from the dictionary
267  void readFields(const dictionary&);
268 
269  //- Read the field - create the field dictionary on-the-fly
270  void readFields();
271 
272 
273 public:
274 
275  //- Runtime type information
276  TypeName("GeometricField");
277 
278 
279  // Public typedefs
281  typedef typename Field<Type>::cmptType cmptType;
282 
283  // Static Member Functions
284 
285  //- Return a null geometric field
286  inline static const GeometricField<Type, PatchField, GeoMesh>& null();
287 
288 
289  // Constructors
290 
291  //- Constructor given IOobject, mesh, dimensions and patch type.
292  // This allocates storage for the field but not values.
293  // Used only within this class to create TEMPORARY variables
295  (
296  const IOobject&,
297  const Mesh&,
298  const dimensionSet&,
299  const word& patchFieldType=PatchField<Type>::calculatedType()
300  );
301 
302  //- Constructor given IOobject, mesh, dimensions and patch types.
303  // This allocates storage for the field but not values.
304  // Used only within this class to create TEMPORARY variables
306  (
307  const IOobject&,
308  const Mesh&,
309  const dimensionSet&,
310  const wordList& wantedPatchTypes,
311  const wordList& actualPatchTypes = wordList()
312  );
313 
314  //- Constructor given IOobject, mesh, dimensioned<Type> and patch type.
316  (
317  const IOobject&,
318  const Mesh&,
319  const dimensioned<Type>&,
320  const word& patchFieldType=PatchField<Type>::calculatedType()
321  );
322 
323  //- Constructor given IOobject, mesh, dimensioned<Type> and patch types.
325  (
326  const IOobject&,
327  const Mesh&,
328  const dimensioned<Type>&,
329  const wordList& wantedPatchTypes,
330  const wordList& actualPatchTypes = wordList()
331  );
332 
333  //- Constructor from components
335  (
336  const IOobject&,
337  const Internal&,
338  const PtrList<PatchField<Type>>&
339  );
340 
341  //- Constructor from components
343  (
344  const IOobject&,
345  const Mesh&,
346  const dimensionSet&,
347  const Field<Type>&,
348  const PtrList<PatchField<Type>>&
349  );
350 
351  //- Construct and read given IOobject
353  (
354  const IOobject&,
355  const Mesh&
356  );
357 
358  //- Construct from dictionary
360  (
361  const IOobject&,
362  const Mesh&,
363  const dictionary&
364  );
365 
366  //- Construct as copy
368  (
370  );
371 
372  //- Construct as copy of tmp<GeometricField> deleting argument
373  #ifndef NoConstructFromTmp
375  (
377  );
378  #endif
379 
380  //- Construct as copy resetting IO parameters
382  (
383  const IOobject&,
385  );
386 
387  //- Construct as copy of tmp<GeometricField> resetting IO parameters
388  #ifndef NoConstructFromTmp
390  (
391  const IOobject&,
393  );
394  #endif
395 
396  //- Construct as copy resetting name
398  (
399  const word& newName,
401  );
402 
403  //- Construct as copy resetting name
404  #ifndef NoConstructFromTmp
406  (
407  const word& newName,
409  );
410  #endif
411 
412  //- Construct as copy resetting IO parameters and patch type
414  (
415  const IOobject&,
417  const word& patchFieldType
418  );
419 
420  //- Construct as copy resetting IO parameters and boundary types
422  (
423  const IOobject&,
425  const wordList& patchFieldTypes,
426  const wordList& actualPatchTypes = wordList()
427  );
428 
429  //- Construct as copy resetting IO parameters and boundary types
430  #ifndef NoConstructFromTmp
432  (
433  const IOobject&,
435  const wordList& patchFieldTypes,
436  const wordList& actualPatchTypes = wordList()
437  );
438  #endif
439 
440  //- Clone
442 
443 
444  //- Destructor
445  virtual ~GeometricField();
446 
447 
448  // Member Functions
449 
450  //- Return a reference to the dimensioned internal field
451  // Note: this increments the event counter and checks the
452  // old-time fields; avoid in loops.
453  Internal& ref();
454 
455  //- Return a const-reference to the dimensioned internal field
456  inline const Internal& internalField() const;
457 
458  //- Return a const-reference to the dimensioned internal field
459  // of a "vol" field. Useful in the formulation of source-terms
460  // for FV equations
461  inline const Internal& v() const;
462 
463  //- Return a reference to the internal field
464  // Note: this increments the event counter and checks the
465  // old-time fields; avoid in loops.
467 
468  //- Return a const-reference to the internal field
469  inline const typename Internal::FieldType& primitiveField() const;
470 
471  //- Return a reference to the boundary field
472  // Note: this increments the event counter and checks the
473  // old-time fields; avoid in loops.
475 
476  //- Return const-reference to the boundary field
477  inline const Boundary& boundaryField() const;
478 
479  //- Return the time index of the field
480  inline label timeIndex() const;
481 
482  //- Return the time index of the field
483  inline label& timeIndex();
484 
485  //- Store the old-time fields
486  void storeOldTimes() const;
487 
488  //- Store the old-time field
489  void storeOldTime() const;
490 
491  //- Return the number of old time fields stored
492  label nOldTimes() const;
493 
494  //- Return old time field
496 
497  //- Return non-const old time field
498  // (Not a good idea but it is used for sub-cycling)
500 
501  //- Store the field as the previous iteration value
502  void storePrevIter() const;
503 
504  //- Return previous iteration field
506 
507  //- Correct boundary field
509 
510  //- Does the field need a reference level for solution
511  bool needReference() const;
512 
513  //- Return a component of the field
515  (
516  const direction
517  ) const;
518 
519  //- WriteData member function required by regIOobject
520  bool writeData(Ostream&) const;
521 
522  //- Return transpose (only if it is a tensor field)
524 
525  //- Relax field (for steady-state solution).
526  // alpha = 1 : no relaxation
527  // alpha < 1 : relaxation
528  // alpha = 0 : do nothing
529  void relax(const scalar alpha);
530 
531  //- Relax field (for steady-state solution).
532  // alpha is read from controlDict
533  void relax();
534 
535  //- Select the final iteration parameters if `final' is true
536  // by returning the field name + "Final"
537  // otherwise the standard parameters by returning the field name
538  word select(bool final) const;
539 
540  //- Helper function to write the min and max to an Ostream
541  void writeMinMax(Ostream& os) const;
542 
543 
544  // Member function *this operators
545 
546  void negate();
547 
548  void replace
549  (
550  const direction,
552  );
553 
554  void replace
555  (
556  const direction,
557  const dimensioned<cmptType>&
558  );
559 
560  void max(const dimensioned<Type>&);
561  void min(const dimensioned<Type>&);
562 
563  void maxMin
564  (
565  const dimensioned<Type>& minDt,
566  const dimensioned<Type>& maxDt
567  );
568 
569  void max
570  (
572  const dimensioned<Type>&
573  );
574 
575  void min
576  (
578  const dimensioned<Type>&
579  );
580 
581  void scale
582  (
585  );
586 
587  void scale
588  (
590  const dimensioned<Type>&
591  );
592 
593 
594  // Member operators
595 
596  //- Return a const-reference to the dimensioned internal field
597  // Useful in the formulation of source-terms for FV equations
598  inline const Internal& operator()() const;
599 
602  void operator=(const dimensioned<Type>&);
603 
605  void operator==(const dimensioned<Type>&);
606 
609 
612 
615 
618 
619  void operator+=(const dimensioned<Type>&);
620  void operator-=(const dimensioned<Type>&);
621 
622  void operator*=(const dimensioned<scalar>&);
623  void operator/=(const dimensioned<scalar>&);
624 
625 
626  // Ostream operators
627 
628  friend Ostream& operator<< <Type, PatchField, GeoMesh>
629  (
630  Ostream&,
632  );
633 
634  friend Ostream& operator<< <Type, PatchField, GeoMesh>
635  (
636  Ostream&,
638  );
639 };
640 
641 
642 template<class Type, template<class> class PatchField, class GeoMesh>
643 Ostream& operator<<
644 (
645  Ostream&,
647  Boundary&
648 );
649 
650 
651 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
652 
653 } // End namespace Foam
654 
655 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
656 
657 #include "GeometricFieldI.H"
658 
659 #ifdef NoRepository
660  #include "GeometricField.C"
661 #endif
662 
663 #include "GeometricFieldFunctions.H"
664 
665 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
666 
667 #endif
668 
669 // ************************************************************************* //
void operator=(const Boundary &)
Assignment to BoundaryField<Type, PatchField, BoundaryMesh>
dictionary dict
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Boundary(const BoundaryMesh &)
Construct from a BoundaryMesh.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
GeoMesh::BoundaryMesh BoundaryMesh
Type of boundary mesh on which this.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
void writeEntry(const word &keyword, Ostream &os) const
Write boundary field as dictionary entry.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
bool writeData(Ostream &) const
WriteData member function required by regIOobject.
MESH::BoundaryMesh BoundaryMesh
Definition: GeoMesh.H:62
uint8_t direction
Definition: direction.H:45
pTraits< Type >::cmptType cmptType
Component type.
Definition: Field.H:86
void evaluate()
Evaluate boundary conditions.
GeoMesh::Mesh Mesh
Type of mesh on which this GeometricField is instantiated.
TypeName("GeometricField")
Runtime type information.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
Generic dimensioned Type class.
LduInterfaceFieldPtrsList< Type > interfaces() const
Return a list of pointers for each patch field with only those.
label nOldTimes() const
Return the number of old time fields stored.
void operator==(const Boundary &)
Forced assignment to.
void maxMin(const dimensioned< Type > &minDt, const dimensioned< Type > &maxDt)
DimensionedField< Type, GeoMesh > Internal
Type of the internal field from which this GeometricField is derived.
word select(bool final) const
Select the final iteration parameters if `final&#39; is true.
Generic field type.
Definition: FieldField.H:51
wordList types() const
Return a list of the patch types.
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
void readField(const Internal &field, const dictionary &dict)
Read the boundary field.
void updateCoeffs()
Update the boundary condition coefficients.
bool needReference() const
Does the field need a reference level for solution.
Dimension set for the base types.
Definition: dimensionSet.H:120
const Internal & v() const
Return a const-reference to the dimensioned internal field.
Definition: volFieldsI.H:31
void negate()
Negate this field.
Definition: FieldField.C:212
void storeOldTimes() const
Store the old-time fields.
Pre-declare SubField and related Field type.
Definition: Field.H:57
A class for handling words, derived from string.
Definition: word.H:59
void storeOldTime() const
Store the old-time field.
void writeMinMax(Ostream &os) const
Helper function to write the min and max to an Ostream.
void operator*=(const FieldField< PatchField, scalar > &)
Definition: FieldField.C:359
void min(const dimensioned< Type > &)
GeometricField(const IOobject &, const Mesh &, const dimensionSet &, const word &patchFieldType=PatchField< Type >::calculatedType())
Constructor given IOobject, mesh, dimensions and patch type.
void operator+=(const FieldField< PatchField, Type > &)
Definition: FieldField.C:357
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const GeometricField< Type, PatchField, GeoMesh > & prevIter() const
Return previous iteration field.
MESH Mesh
Definition: GeoMesh.H:61
void operator/=(const FieldField< PatchField, scalar > &)
Definition: FieldField.C:360
const Field< Type > & field() const
tmp< FieldField< PatchField, cmptType > > component(const direction) const
Return a component field of the field.
Definition: FieldField.C:224
virtual ~GeometricField()
Destructor.
void replace(const direction, const FieldField< PatchField, cmptType > &)
Replace a component field of the field.
Definition: FieldField.C:242
label timeIndex() const
Return the time index of the field.
Internal & ref()
Return a reference to the dimensioned internal field.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
List< word > wordList
A List of words.
Definition: fileName.H:54
pTraits< Type >::cmptType cmptType
Component type.
Definition: FieldField.H:82
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
PatchField< Type > Patch
Type of the patch field of which the.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
void max(const dimensioned< Type > &)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void correctBoundaryConditions()
Correct boundary field.
lduInterfaceFieldPtrsList scalarInterfaces() const
Return a list of pointers for each patch field with only those.
tmp< FieldField< PatchField, Type > > T() const
Return the field transpose (only defined for second rank tensors)
Definition: FieldField.C:269
void storePrevIter() const
Store the field as the previous iteration value.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
void operator-=(const FieldField< PatchField, Type > &)
Definition: FieldField.C:358
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const T * operator()(const label) const
Return element const pointer.
Definition: UPtrListI.H:135
void scale(const GeometricField< Type, PatchField, GeoMesh > &, const GeometricField< Type, PatchField, GeoMesh > &)
Boundary boundaryInternalField() const
Return BoundaryField of the cell values neighbouring.
List of coupled interface fields to be used in coupling.
void relax()
Relax field (for steady-state solution).
tmp< FieldField< PatchField, Type > > clone() const
Clone.
Definition: FieldField.C:182
Namespace for OpenFOAM.