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-2015 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  GeometricBoundaryField.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  typedef typename GeoMesh::Mesh Mesh;
99  typedef typename GeoMesh::BoundaryMesh BoundaryMesh;
103  typedef PatchField<Type> PatchFieldType;
104 
107  :
108  public FieldField<PatchField, Type>
109  {
110  // Private data
111 
112  //- Reference to BoundaryMesh for which this field is defined
113  const BoundaryMesh& bmesh_;
114 
115 
116  public:
117 
118  // Constructors
119 
120  //- Construct from a BoundaryMesh
121  GeometricBoundaryField(const BoundaryMesh&);
122 
123  //- Construct from a BoundaryMesh,
124  // reference to the internal field
125  // and a patch type
127  (
128  const BoundaryMesh&,
129  const DimensionedInternalField&,
130  const word&
131  );
132 
133  //- Construct from a BoundaryMesh,
134  // reference to the internal field
135  // and a wordList of patch types and optional the actual patch
136  // types (to override constraint patches)
138  (
139  const BoundaryMesh&,
140  const DimensionedInternalField&,
141  const wordList& wantedPatchTypes,
142  const wordList& actualPatchTypes = wordList()
143  );
144 
145  //- Construct from a BoundaryMesh,
146  // reference to the internal field
147  // and a PtrList<PatchField<Type> >
149  (
150  const BoundaryMesh&,
151  const DimensionedInternalField&,
152  const PtrList<PatchField<Type> >&
153  );
154 
155  //- Construct as copy setting the reference to the internal field
157  (
158  const DimensionedInternalField&,
160  );
161 
162  //- Construct as copy
163  // Dangerous because Field may be set to a field which gets deleted
164  // Need new type of BoundaryField, one which is part of a geometric
165  // field for which snGrad etc. may be called and a free standing
166  // BoundaryField for which such operations are unavailable.
168  (
170  );
171 
172  //- Construct from dictionary
174  (
175  const BoundaryMesh&,
176  const DimensionedInternalField&,
177  const dictionary&
178  );
179 
180 
181  // Member functions
182 
183  //- Read the boundary field
184  void readField
185  (
187  const dictionary& dict
188  );
189 
190  //- Update the boundary condition coefficients
191  void updateCoeffs();
192 
193  //- Evaluate boundary conditions
194  void evaluate();
195 
196  //- Return a list of the patch types
197  wordList types() const;
198 
199  //- Return BoundaryField of the cell values neighbouring
200  // the boundary
202 
203  //- Return a list of pointers for each patch field with only those
204  // pointing to interfaces being set
206 
207  //- Return a list of pointers for each patch field with only those
208  // pointing to interfaces being set
210 
211  //- Write boundary field as dictionary entry
212  void writeEntry(const word& keyword, Ostream& os) const;
213 
214 
215  // Member operators
216 
217  //- Assignment to BoundaryField<Type, PatchField, BoundaryMesh>
218  void operator=(const GeometricBoundaryField&);
219 
220  //- Assignment to FieldField<PatchField, Type>
222 
223  //- Assignment to Type
224  void operator=(const Type&);
225 
226 
227  //- Forced assignment to
228  // BoundaryField<Type, PatchField, BoundaryMesh>
229  void operator==(const GeometricBoundaryField&);
230 
231  //- Forced assignment to FieldField<PatchField, Type>
233 
234  //- Forced assignment to Type
235  void operator==(const Type&);
236  };
237 
238 
239 private:
240 
241  // Private data
242 
243  //- Current time index.
244  // Used to trigger the storing of the old-time value
245  mutable label timeIndex_;
246 
247  //- Pointer to old time field
248  mutable GeometricField<Type, PatchField, GeoMesh>* field0Ptr_;
249 
250  //- Pointer to previous iteration (used for under-relaxation)
251  mutable GeometricField<Type, PatchField, GeoMesh>* fieldPrevIterPtr_;
252 
253  //- Boundary Type field containing boundary field values
254  GeometricBoundaryField boundaryField_;
255 
256 
257  // Private Member Functions
258 
259  //- Read the field from the dictionary
260  void readFields(const dictionary&);
261 
262  //- Read the field - create the field dictionary on-the-fly
263  void readFields();
264 
265 
266 public:
267 
268  //- Runtime type information
269  TypeName("GeometricField");
270 
271 
272  // Public typedefs
274  typedef typename Field<Type>::cmptType cmptType;
275 
276  // Static Member Functions
277 
278  //- Return a null geometric field
279  inline static const GeometricField<Type, PatchField, GeoMesh>& null();
280 
281 
282  // Constructors
283 
284  //- Constructor given IOobject, mesh, dimensions and patch type.
285  // This allocates storage for the field but not values.
286  // Used only within this class to create TEMPORARY variables
288  (
289  const IOobject&,
290  const Mesh&,
291  const dimensionSet&,
292  const word& patchFieldType=PatchField<Type>::calculatedType()
293  );
294 
295  //- Constructor given IOobject, mesh, dimensions and patch types.
296  // This allocates storage for the field but not values.
297  // Used only within this class to create TEMPORARY variables
299  (
300  const IOobject&,
301  const Mesh&,
302  const dimensionSet&,
303  const wordList& wantedPatchTypes,
304  const wordList& actualPatchTypes = wordList()
305  );
306 
307  //- Constructor given IOobject, mesh, dimensioned<Type> and patch type.
309  (
310  const IOobject&,
311  const Mesh&,
312  const dimensioned<Type>&,
313  const word& patchFieldType=PatchField<Type>::calculatedType()
314  );
315 
316  //- Constructor given IOobject, mesh, dimensioned<Type> and patch types.
318  (
319  const IOobject&,
320  const Mesh&,
321  const dimensioned<Type>&,
322  const wordList& wantedPatchTypes,
323  const wordList& actualPatchTypes = wordList()
324  );
325 
326  //- Constructor from components
328  (
329  const IOobject&,
330  const Mesh&,
331  const dimensionSet&,
332  const Field<Type>&,
333  const PtrList<PatchField<Type> >&
334  );
335 
336  //- Construct and read given IOobject
338  (
339  const IOobject&,
340  const Mesh&
341  );
342 
343  //- Construct from dictionary
345  (
346  const IOobject&,
347  const Mesh&,
348  const dictionary&
349  );
350 
351  //- Construct as copy
353  (
355  );
356 
357  //- Construct as copy of tmp<GeometricField> deleting argument
358  #ifndef NoConstructFromTmp
360  (
362  );
363  #endif
364 
365  //- Construct as copy resetting IO parameters
367  (
368  const IOobject&,
370  );
371 
372  //- Construct as copy of tmp<GeometricField> resetting IO parameters
373  #ifndef NoConstructFromTmp
375  (
376  const IOobject&,
378  );
379  #endif
380 
381  //- Construct as copy resetting name
383  (
384  const word& newName,
386  );
387 
388  //- Construct as copy resetting name
389  #ifndef NoConstructFromTmp
391  (
392  const word& newName,
394  );
395  #endif
396 
397  //- Construct as copy resetting IO parameters and patch type
399  (
400  const IOobject&,
402  const word& patchFieldType
403  );
404 
405  //- Construct as copy resetting IO parameters and boundary types
407  (
408  const IOobject&,
410  const wordList& patchFieldTypes,
411  const wordList& actualPatchTypes = wordList()
412  );
413 
414 
415  //- Destructor
416  virtual ~GeometricField();
417 
418 
419  // Member Functions
420 
421  //- Return dimensioned internal field
422  DimensionedInternalField& dimensionedInternalField();
423 
424  //- Return dimensioned internal field
425  inline const DimensionedInternalField& dimensionedInternalField() const;
426 
427  //- Return internal field
428  InternalField& internalField();
429 
430  //- Return internal field
431  inline const InternalField& internalField() const;
432 
433  //- Return reference to GeometricBoundaryField
435 
436  //- Return reference to GeometricBoundaryField for const field
437  inline const GeometricBoundaryField& boundaryField() const;
438 
439  //- Return the time index of the field
440  inline label timeIndex() const;
441 
442  //- Return the time index of the field
443  inline label& timeIndex();
444 
445  //- Store the old-time fields
446  void storeOldTimes() const;
447 
448  //- Store the old-time field
449  void storeOldTime() const;
450 
451  //- Return the number of old time fields stored
452  label nOldTimes() const;
453 
454  //- Return old time field
456 
457  //- Return non-const old time field
458  // (Not a good idea but it is used for sub-cycling)
460 
461  //- Store the field as the previous iteration value
462  void storePrevIter() const;
463 
464  //- Return previous iteration field
466 
467  //- Correct boundary field
469 
470  //- Does the field need a reference level for solution
471  bool needReference() const;
472 
473  //- Return a component of the field
475  (
476  const direction
477  ) const;
478 
479  //- WriteData member function required by regIOobject
480  bool writeData(Ostream&) const;
481 
482  //- Return transpose (only if it is a tensor field)
484 
485  //- Relax field (for steady-state solution).
486  // alpha = 1 : no relaxation
487  // alpha < 1 : relaxation
488  // alpha = 0 : do nothing
489  void relax(const scalar alpha);
490 
491  //- Relax field (for steady-state solution).
492  // alpha is read from controlDict
493  void relax();
494 
495  //- Select the final iteration parameters if `final' is true
496  // by returning the field name + "Final"
497  // otherwise the standard parameters by returning the field name
498  word select(bool final) const;
499 
500  //- Helper function to write the min and max to an Ostream
501  void writeMinMax(Ostream& os) const;
502 
503 
504  // Member function *this operators
505 
506  void negate();
507 
508  void replace
509  (
510  const direction,
512  );
513 
514  void replace
515  (
516  const direction,
517  const dimensioned<cmptType>&
518  );
519 
520  void max(const dimensioned<Type>&);
521  void min(const dimensioned<Type>&);
522 
523  void max
524  (
526  const dimensioned<Type>&
527  );
528 
529  void min
530  (
532  const dimensioned<Type>&
533  );
534 
535  void scale
536  (
539  );
540 
541  void scale
542  (
544  const dimensioned<Type>&
545  );
546 
547 
548  // Member operators
549 
552  void operator=(const dimensioned<Type>&);
553 
555  void operator==(const dimensioned<Type>&);
556 
559 
562 
565 
568 
569  void operator+=(const dimensioned<Type>&);
570  void operator-=(const dimensioned<Type>&);
571 
572  void operator*=(const dimensioned<scalar>&);
573  void operator/=(const dimensioned<scalar>&);
574 
575 
576  // Ostream operators
577 
578  friend Ostream& operator<< <Type, PatchField, GeoMesh>
579  (
580  Ostream&,
582  );
583 
584  friend Ostream& operator<< <Type, PatchField, GeoMesh>
585  (
586  Ostream&,
588  );
589 };
590 
591 
592 template<class Type, template<class> class PatchField, class GeoMesh>
593 Ostream& operator<<
594 (
595  Ostream&,
597  GeometricBoundaryField&
598 );
599 
600 
601 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
602 
603 } // End namespace Foam
604 
605 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
606 
607 #include "GeometricFieldI.H"
608 
609 #ifdef NoRepository
610 # include "GeometricField.C"
611 #endif
612 
613 #include "GeometricFieldFunctions.H"
614 
615 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
616 
617 #endif
618 
619 // ************************************************************************* //
tmp< FieldField< PatchField, cmptType > > component(const direction) const
Return a component field of the field.
unsigned char direction
Definition: direction.H:43
void evaluate()
Evaluate boundary conditions.
virtual ~GeometricField()
Destructor.
bool needReference() const
Does the field need a reference level for solution.
void replace(const direction, const FieldField< PatchField, cmptType > &)
Replace a component field of the field.
Generic field type.
Definition: FieldField.H:51
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void writeEntry(const word &keyword, Ostream &os) const
Write boundary field as dictionary entry.
const Field< Type > & field() const
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
MESH::BoundaryMesh BoundaryMesh
Definition: GeoMesh.H:62
void operator-=(const FieldField< PatchField, Type > &)
PatchField< Type > PatchFieldType
A class for handling words, derived from string.
Definition: word.H:59
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
InternalField & internalField()
Return internal field.
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
tmp< FieldField< PatchField, Type > > T() const
Return the field transpose (only defined for second rank tensors)
LduInterfaceFieldPtrsList< Type > interfaces() const
Return a list of pointers for each patch field with only those.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
void writeMinMax(Ostream &os) const
Helper function to write the min and max to an Ostream.
void readField(const DimensionedField< Type, GeoMesh > &field, const dictionary &dict)
Read the boundary field.
GeometricBoundaryField boundaryInternalField() const
Return BoundaryField of the cell values neighbouring.
dictionary dict
Generic dimensioned Type class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
GeoMesh::BoundaryMesh BoundaryMesh
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, scalar > &)
void updateCoeffs()
Update the boundary condition coefficients.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
pTraits< Type >::cmptType cmptType
Component type.
Definition: Field.H:86
pTraits< Type >::cmptType cmptType
Component type.
Definition: FieldField.H:82
Dimension set for the base types.
Definition: dimensionSet.H:116
const GeometricField< Type, PatchField, GeoMesh > & prevIter() const
Return previous iteration field.
MESH Mesh
Definition: GeoMesh.H:61
label nOldTimes() const
Return the number of old time fields stored.
GeoMesh::Mesh Mesh
TypeName("GeometricField")
Runtime type information.
Pre-declare SubField and related Field type.
Definition: Field.H:57
void operator==(const GeometricBoundaryField &)
Forced assignment to.
wordList types() const
Return a list of the patch types.
void operator*=(const FieldField< PatchField, scalar > &)
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
void scale(const GeometricField< Type, PatchField, GeoMesh > &, const GeometricField< Type, PatchField, GeoMesh > &)
Generic GeometricField class.
void operator+=(const FieldField< PatchField, Type > &)
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
List< word > wordList
A List of words.
Definition: fileName.H:54
void correctBoundaryConditions()
Correct boundary field.
DimensionedField< Type, GeoMesh > DimensionedInternalField
void negate()
Negate this field.
lduInterfaceFieldPtrsList scalarInterfaces() const
Return a list of pointers for each patch field with only those.
void max(const dimensioned< Type > &)
label timeIndex() const
Return the time index of the field.
GeometricBoundaryField(const BoundaryMesh &)
Construct from a BoundaryMesh.
void operator=(const GeometricBoundaryField &)
Assignment to BoundaryField<Type, PatchField, BoundaryMesh>
void storeOldTimes() const
Store the old-time fields.
bool writeData(Ostream &) const
WriteData member function required by regIOobject.
void storePrevIter() const
Store the field as the previous iteration value.
void relax()
Relax field (for steady-state solution).
void min(const dimensioned< Type > &)
A class for managing temporary objects.
Definition: PtrList.H:118
List of coupled interface fields to be used in coupling.
Field< Type > InternalField
void storeOldTime() const
Store the old-time field.
word select(bool final) const
Select the final iteration parameters if `final&#39; is true.