MapGeometricFields.H
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-2018 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::MapInternalField
26 
27 Description
28  Generic internal field mapper. For "real" mapping, add template
29  specialisations for mapping of internal fields depending on mesh
30  type.
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #ifndef MapGeometricFields_H
35 #define MapGeometricFields_H
36 
37 #include "polyMesh.H"
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43 
44 template<class Type, class MeshMapper, class GeoMesh>
45 class MapInternalField
46 {
47 public:
48 
50  {}
51 
52  void operator()
53  (
54  Field<Type>& field,
55  const MeshMapper& mapper
56  ) const;
57 };
58 
59 
60 //- Generic Geometric field mapper.
61 // For "real" mapping, add template specialisations
62 // for mapping of internal fields depending on mesh type.
63 template
64 <
65  class Type,
66  template<class> class PatchField,
67  class MeshMapper,
68  class GeoMesh
69 >
71 (
72  const MeshMapper& mapper
73 )
74 {
76  (
77  mapper.thisDb().objectRegistry::template
78  lookupClass<GeometricField<Type, PatchField, GeoMesh>>()
79  );
80 
81  // It is necessary to enforce that all old-time fields are stored
82  // before the mapping is performed. Otherwise, if the
83  // old-time-level field is mapped before the field itself, sizes
84  // will not match.
85 
86  for
87  (
89  iterator fieldIter = fields.begin();
90  fieldIter != fields.end();
91  ++fieldIter
92  )
93  {
96  (*fieldIter());
97 
98  // Note: check can be removed once pointFields are actually stored on
99  // the pointMesh instead of now on the polyMesh!
100  if (&field.mesh() == &mapper.mesh())
101  {
102  field.storeOldTimes();
103  }
104  }
105 
106  for
107  (
109  iterator fieldIter = fields.begin();
110  fieldIter != fields.end();
111  ++fieldIter
112  )
113  {
116  (*fieldIter());
117 
118  if (&field.mesh() == &mapper.mesh())
119  {
120  if (polyMesh::debug)
121  {
122  Info<< "Mapping " << field.typeName << ' ' << field.name()
123  << endl;
124  }
125 
126  // Map the internal field
128  (
129  field.primitiveFieldRef(),
130  mapper
131  );
132 
133  // Map the patch fields
135  ::Boundary& bfield = field.boundaryFieldRef();
136  forAll(bfield, patchi)
137  {
138  // Cannot check sizes for patch fields because of
139  // empty fields in FV and because point fields get their size
140  // from the patch which has already been resized
141  //
142 
143  bfield[patchi].autoMap(mapper.boundaryMap()[patchi]);
144  }
145 
146  field.instance() = field.time().timeName();
147  }
148  else if (polyMesh::debug)
149  {
150  Info<< "Not mapping " << field.typeName << ' ' << field.name()
151  << " since originating mesh differs from that of mapper."
152  << endl;
153  }
154  }
155 }
156 
157 
158 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
159 
160 } // End namespace Foam
161 
162 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
163 
164 #endif
165 
166 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const word & name() const
Return name.
Definition: IOobject.H:297
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:110
static const char *const typeName
Definition: Field.H:94
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Generic GeometricField class.
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
void MapGeometricFields(const MeshMapper &mapper)
Generic Geometric field mapper.
void storeOldTimes() const
Store the old-time fields.
Pre-declare SubField and related Field type.
Definition: Field.H:57
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:419
An STL-conforming hash table.
Definition: HashTable.H:62
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const Mesh & mesh() const
Return mesh.
Generic internal field mapper. For "real" mapping, add template specialisations for mapping of intern...
const fileName & instance() const
Definition: IOobject.H:392
label patchi
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const Time & time() const
Return time.
Definition: IOobject.C:367
messageStream Info
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
Namespace for OpenFOAM.