FreeFOAM The Cross-Platform CFD Toolkit
fvMeshSubset.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) 1991-2010 OpenCFD Ltd.
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 minimize 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 the mapping becomes a problem.
45  Do the new faces get the value of the internal face they came from?
46  What if e.g. the user supplied patch is a fixedValue 0? So for now
47  they get the face of existing patch face 0.
48 
49 SourceFiles
50  fvMeshSubset.C
51 
52 \*---------------------------------------------------------------------------*/
53 
54 #ifndef fvMeshSubset_H
55 #define fvMeshSubset_H
56 
57 #include <finiteVolume/fvMesh.H>
58 #include <OpenFOAM/pointMesh.H>
62 #include <OpenFOAM/HashSet.H>
64 
65 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
66 
67 namespace Foam
68 {
69 
70 /*---------------------------------------------------------------------------*\
71  Class fvMeshSubset Declaration
72 \*---------------------------------------------------------------------------*/
73 
75 {
76 
77 public:
78 
79  //- Patch-field subset interpolation class
81  :
82  public fvPatchFieldMapper
83  {
84  const labelList& directAddressing_;
85 
86  public:
87 
88  // Constructors
89 
90  //- Construct given addressing
92  :
93  directAddressing_(directAddressing)
94  {}
95 
96  // Destructor
97 
99  {}
100 
101 
102  // Member Functions
103 
104  label size() const
105  {
106  return directAddressing_.size();
107  }
108 
109  bool direct() const
110  {
111  return true;
112  }
113 
115  {
116  return directAddressing_;
117  }
118  };
119 
120 
121  //- Patch-field subset interpolation class
123  :
124  public pointPatchFieldMapper
125  {
126  const labelList& directAddressing_;
127 
128  public:
129 
130  // Constructors
131 
132  //- Construct given addressing
134  :
135  directAddressing_(directAddressing)
136  {}
137 
138  // Destructor
139 
141  {}
142 
143 
144  // Member Functions
145 
146  label size() const
147  {
148  return directAddressing_.size();
149  }
150 
151  bool direct() const
152  {
153  return true;
154  }
155 
157  {
158  return directAddressing_;
159  }
160  };
161 
162 
163 private:
164 
165  // Private data
166 
167  //- Mesh to subset from
168  const fvMesh& baseMesh_;
169 
170  //- Subset mesh pointer
171  autoPtr<fvMesh> fvMeshSubsetPtr_;
172 
173  //- Point mapping array
174  labelList pointMap_;
175 
176  //- Face mapping array
177  labelList faceMap_;
178 
179  //- Cell mapping array
180  labelList cellMap_;
181 
182  //- Patch mapping array
183  labelList patchMap_;
184 
185 
186  // Private Member Functions
187 
188  //- Check if subset has been performed
189  bool checkCellSubset() const;
190 
191  //- Mark points in Map
192  static void markPoints(const labelList&, Map<label>&);
193 
194  //- Mark points (with 0) in labelList
195  static void markPoints(const labelList&, labelList&);
196 
197  //- Adapt nCellsUsingFace for coupled faces becoming 'uncoupled'.
198  void doCoupledPatches
199  (
200  const bool syncPar,
201  labelList& nCellsUsingFace
202  ) const;
203 
204  //- Subset of subset
205  static labelList subset
206  (
207  const label nElems,
208  const labelList& selectedElements,
209  const labelList& subsetMap
210  );
211 
212  //- Create zones for submesh
213  void subsetZones();
214 
215  //- Disallow default bitwise copy construct
216  fvMeshSubset(const fvMeshSubset&);
217 
218  //- Disallow default bitwise assignment
219  void operator=(const fvMeshSubset&);
220 
221 public:
222 
223  // Constructors
224 
225  //- Construct given a mesh to subset
226  explicit fvMeshSubset(const fvMesh&);
227 
228 
229  // Member Functions
230 
231  // Edit
232 
233  //- Set the subset. Create "oldInternalFaces" patch for exposed
234  // internal faces (patchID==-1) or use supplied patch.
235  // Does not handle coupled patches correctly if only one side
236  // gets deleted.
237  void setCellSubset
238  (
239  const labelHashSet& globalCellMap,
240  const label patchID = -1
241  );
242 
243  //- Set the subset from all cells with region == currentRegion.
244  // Create "oldInternalFaces" patch for exposed
245  // internal faces (patchID==-1) or use supplied patch.
246  // Handles coupled patches by if nessecary making coupled patch
247  // face part of patchID (so uncoupled)
248  void setLargeCellSubset
249  (
250  const labelList& region,
251  const label currentRegion,
252  const label patchID = -1,
253  const bool syncCouples = true
254  );
255 
256  //- setLargeCellSubset but with labelHashSet.
257  void setLargeCellSubset
258  (
259  const labelHashSet& globalCellMap,
260  const label patchID = -1,
261  const bool syncPar = true
262  );
263 
264 
265  // Access
266 
267  //- Original mesh
268  const fvMesh& baseMesh() const
269  {
270  return baseMesh_;
271  }
272 
273  //- Return reference to subset mesh
274  const fvMesh& subMesh() const;
275 
276  fvMesh& subMesh();
277 
278  //- Return point map
279  const labelList& pointMap() const;
280 
281  //- Return face map
282  const labelList& faceMap() const;
283 
284  //- Return cell map
285  const labelList& cellMap() const;
286 
287  //- Return patch map
288  const labelList& patchMap() const;
289 
290 
291  // Field mapping
292 
293  //- Map volume field
294  template<class Type>
297  (
299  const fvMesh& sMesh,
300  const labelList& patchMap,
301  const labelList& cellMap,
302  const labelList& faceMap
303  );
304 
305  template<class Type>
308  (
310  ) const;
311 
312  //- Map surface field
313  template<class Type>
316  (
318  const fvMesh& sMesh,
319  const labelList& patchMap,
320  const labelList& faceMap
321  );
322 
323  template<class Type>
326  (
328  ) const;
329 
330  //- Map point field
331  template<class Type>
334  (
336  const pointMesh& sMesh,
337  const labelList& patchMap,
338  const labelList& pointMap
339  );
340 
341  template<class Type>
344  (
346  ) const;
347 };
348 
349 
350 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
351 
352 } // End namespace Foam
353 
354 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
355 
356 #ifdef NoRepository
357 # include "fvMeshSubsetInterpolate.C"
358 #endif
359 
360 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
361 
362 #endif
363 
364 // ************************ vim: set sw=4 sts=4 et: ************************ //