fvMeshSubset.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-2023 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::fvMeshSubset
26 
27 Description
28  Post-processing mesh subset tool. Given the original mesh and the
29  list of selected cells, it creates the mesh consisting only of the
30  desired cells, with the mapping list for points, faces, and cells.
31 
32  Puts all exposed internal faces into either
33  - a user supplied patch
34  - a newly created patch "oldInternalFaces"
35 
36  - setCellSubset is for small subsets. Uses Maps to minimise memory.
37  - setLargeCellSubset is for largish subsets (>10% of mesh).
38  Uses labelLists instead.
39 
40  - setLargeCellSubset does coupled patch subsetting as well. If it detects
41  a face on a coupled patch 'losing' its neighbour it will move the
42  face into the oldInternalFaces patch.
43 
44  - if a user supplied patch is used it is up to the destination
45  patchField to handle exposed internal faces (mapping from face -1).
46  If not provided the default is to assign the internalField. All the
47  basic patch field types (e.g. fixedValue) will give a warning and
48  preferably derived patch field types should be used that know how to
49  handle exposed faces (e.g. use uniformFixedValue instead of fixedValue)
50 
51 SourceFiles
52  fvMeshSubset.C
53 
54 \*---------------------------------------------------------------------------*/
55 
56 #ifndef fvMeshSubset_H
57 #define fvMeshSubset_H
58 
59 #include "fvMesh.H"
60 #include "HashSet.H"
61 
62 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
63 
64 namespace Foam
65 {
66 
67 class pointMesh;
68 class surfaceMesh;
69 
70 /*---------------------------------------------------------------------------*\
71  Class fvMeshSubset Declaration
72 \*---------------------------------------------------------------------------*/
73 
74 class fvMeshSubset
75 {
76  // Private Data
77 
78  //- Mesh to subset from
79  const fvMesh& baseMesh_;
80 
81  //- Subset mesh pointer
82  autoPtr<fvMesh> fvMeshSubsetPtr_;
83 
84  //- Point mapping array
85  labelList pointMap_;
86 
87  //- Face mapping array
88  labelList faceMap_;
89 
90  //- Cell mapping array
91  labelList cellMap_;
92 
93  //- Patch mapping array
94  labelList patchMap_;
95 
96  //- Optional face mapping array with flip encoded
97  mutable autoPtr<labelList> faceFlipMapPtr_;
98 
99 
100  // Private Member Functions
101 
102  //- Check if subset has been performed
103  bool checkCellSubset() const;
104 
105  //- Mark points in Map
106  static void markPoints(const labelList&, Map<label>&);
107 
108  //- Mark points (with 0) in labelList
109  static void markPoints(const labelList&, labelList&);
110 
111  //- Adapt nCellsUsingFace for coupled faces becoming 'uncoupled'.
112  void doCoupledPatches
113  (
114  const bool syncPar,
115  Map<label>& facesToSubset,
116  labelList& nCellsUsingFace
117  ) const;
118 
119  //- Subset of subset
120  static labelList subset
121  (
122  const label nElems,
123  const labelList& selectedElements,
124  const labelList& subsetMap
125  );
126 
127  //- Create zones for submesh
128  void subsetZones();
129 
130  //- Helper: extract cells-to-remove from cells-to-keep
131  labelList getCellsToRemove
132  (
133  const labelList& region,
134  const label currentRegion
135  ) const;
136 
137 
138 public:
139 
140  //- Runtime type information
141  TypeName("fvMeshSubset");
142 
143 
144  // Constructors
145 
146  //- Construct given a mesh to subset
147  explicit fvMeshSubset(const fvMesh&);
148 
149  //- Disallow default bitwise copy construction
150  fvMeshSubset(const fvMeshSubset&) = delete;
151 
152 
153  //- Destructor
154  virtual ~fvMeshSubset();
155 
156 
157  // Member Functions
158 
159  // Edit
160 
161  //- Set the subset. Create "oldInternalFaces" patch for exposed
162  // internal faces (patchID==-1) or use supplied patch.
163  // Does not handle coupled patches correctly if only one side
164  // gets deleted.
165  void setCellSubset
166  (
167  const labelHashSet& globalCellMap,
168  const label patchID = -1,
169  const bool syncPar = true
170  );
171 
172  //- Set the subset from all cells with region == currentRegion.
173  // Create "oldInternalFaces" patch for exposed
174  // internal faces (patchID==-1) or use supplied patch.
175  // Handles coupled patches by if necessary making coupled patch
176  // face part of patchID (so uncoupled)
177  void setLargeCellSubset
178  (
179  const labelList& region,
180  const label currentRegion,
181  const label patchID = -1,
182  const bool syncCouples = true
183  );
184 
185  //- setLargeCellSubset but with labelHashSet.
186  void setLargeCellSubset
187  (
188  const labelHashSet& globalCellMap,
189  const label patchID = -1,
190  const bool syncPar = true
191  );
192 
193 
194  // Two step subsetting
195 
196  //- Get labels of exposed faces.
197  // These are
198  // - internal faces that become boundary faces
199  // - coupled faces that become uncoupled (since one of the
200  // sides gets deleted)
202  (
203  const labelList& region,
204  const label currentRegion,
205  const bool syncCouples = true
206  ) const;
207 
208  //- For every exposed face (from above getExposedFaces)
209  // used supplied (existing!) patch
210  void setLargeCellSubset
211  (
212  const labelList& region,
213  const label currentRegion,
214  const labelList& exposedFaces,
215  const labelList& patchIDs,
216  const bool syncCouples = true
217  );
218 
219 
220  // Access
221 
222  //- Original mesh
223  const fvMesh& baseMesh() const
224  {
225  return baseMesh_;
226  }
227 
228  //- Have subMesh?
229  bool hasSubMesh() const;
230 
231  //- Return reference to subset mesh
232  const fvMesh& subMesh() const;
233 
234  fvMesh& subMesh();
235 
236  //- Return point map
237  const labelList& pointMap() const;
238 
239  //- Return face map
240  const labelList& faceMap() const;
241 
242  //- Return face map with sign to encode flipped faces
243  const labelList& faceFlipMap() const;
244 
245  //- Return cell map
246  const labelList& cellMap() const;
247 
248  //- Return patch map
249  const labelList& patchMap() const;
250 
251 
252  // Field mapping
253 
254  //- Map volume field
255  template<class Type>
258  (
259  const VolField<Type>&,
260  const fvMesh& sMesh,
261  const labelList& patchMap,
262  const labelList& cellMap,
263  const labelList& faceMap
264  );
265 
266  template<class Type>
269  (
270  const VolField<Type>&
271  ) const;
272 
273  //- Map surface field. Optionally negates value if flipping
274  // a face (from exposing an internal face)
275  template<class Type>
278  (
279  const SurfaceField<Type>&,
280  const fvMesh& sMesh,
281  const labelList& patchMap,
282  const labelList& cellMap,
283  const labelList& faceMap
284  );
285 
286  template<class Type>
289  (
290  const SurfaceField<Type>&
291  ) const;
292 
293  //- Map point field
294  template<class Type>
297  (
298  const PointField<Type>&,
299  const pointMesh& sMesh,
300  const labelList& patchMap,
301  const labelList& pointMap
302  );
303 
304  template<class Type>
307  (
308  const PointField<Type>&
309  ) const;
310 
311  //- Map dimensioned field
312  template<class Type>
315  (
317  const fvMesh& sMesh,
318  const labelList& cellMap
319  );
320 
321  template<class Type>
324 
325 
326  // Member Operators
327 
328  //- Disallow default bitwise assignment
329  void operator=(const fvMeshSubset&) = delete;
330 };
331 
332 
333 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
334 
335 } // End namespace Foam
336 
337 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
338 
339 #ifdef NoRepository
340  #include "fvMeshSubsetInterpolate.C"
341 #endif
342 
343 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
344 
345 #endif
346 
347 // ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Post-processing mesh subset tool. Given the original mesh and the list of selected cells,...
Definition: fvMeshSubset.H:74
const labelList & faceMap() const
Return face map.
const fvMesh & baseMesh() const
Original mesh.
Definition: fvMeshSubset.H:222
void setCellSubset(const labelHashSet &globalCellMap, const label patchID=-1, const bool syncPar=true)
Set the subset. Create "oldInternalFaces" patch for exposed.
Definition: fvMeshSubset.C:490
const labelList & cellMap() const
Return cell map.
void operator=(const fvMeshSubset &)=delete
Disallow default bitwise assignment.
fvMeshSubset(const fvMesh &)
Construct given a mesh to subset.
Definition: fvMeshSubset.C:469
bool hasSubMesh() const
Have subMesh?
TypeName("fvMeshSubset")
Runtime type information.
const labelList & faceFlipMap() const
Return face map with sign to encode flipped faces.
virtual ~fvMeshSubset()
Destructor.
Definition: fvMeshSubset.C:483
void setLargeCellSubset(const labelList &region, const label currentRegion, const label patchID=-1, const bool syncCouples=true)
Set the subset from all cells with region == currentRegion.
Definition: fvMeshSubset.C:893
const fvMesh & subMesh() const
Return reference to subset mesh.
const labelList & patchMap() const
Return patch map.
const labelList & pointMap() const
Return point map.
labelList getExposedFaces(const labelList &region, const label currentRegion, const bool syncCouples=true) const
Get labels of exposed faces.
static tmp< VolField< Type > > interpolate(const VolField< Type > &, const fvMesh &sMesh, const labelList &patchMap, const labelList &cellMap, const labelList &faceMap)
Map volume field.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:53
A class for managing temporary objects.
Definition: tmp.H:55
Namespace for OpenFOAM.
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