FreeFOAM The Cross-Platform CFD Toolkit
removeCellLayer.C
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 Description
25  Remove a layer of cells and prepare addressing data
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "layerAdditionRemoval.H"
30 #include <OpenFOAM/polyMesh.H>
31 #include <OpenFOAM/primitiveMesh.H>
33 #include <OpenFOAM/oppositeFace.H>
39 
40 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41 
42 bool Foam::layerAdditionRemoval::validCollapse() const
43 {
44  // Check for valid layer collapse
45  // - no boundary-to-boundary collapse
46 
47  if (debug)
48  {
49  Pout << "Checking layer collapse for object " << name() << endl;
50  }
51 
52  // Grab the face collapse mapping
53  const polyMesh& mesh = topoChanger().mesh();
54 
55  const labelList& ftc = facesPairing();
56  const labelList& mf = mesh.faceZones()[faceZoneID_.index()];
57 
58  label nBoundaryHits = 0;
59 
60  forAll (mf, faceI)
61  {
62  if
63  (
64  !mesh.isInternalFace(mf[faceI])
65  && !mesh.isInternalFace(ftc[faceI])
66  )
67  {
68  nBoundaryHits++;
69  }
70  }
71 
72 
73  if (debug)
74  {
75  Pout<< "Finished checking layer collapse for object "
76  << name() <<". Number of boundary-on-boundary hits: "
77  << nBoundaryHits << endl;
78  }
79 
80  if (returnReduce(nBoundaryHits, sumOp<label>()) > 0)
81  {
82  return false;
83  }
84  else
85  {
86  return true;
87  }
88 }
89 
90 void Foam::layerAdditionRemoval::removeCellLayer
91 (
92  polyTopoChange& ref
93 ) const
94 {
95  // Algorithm for layer removal. Second phase: topological change
96  // 2) Go through all the faces of the master cell layer and remove
97  // the ones that are not in the master face zone.
98  //
99  // 3) Grab all the faces coming out of points that are collapsed
100  // and select the ones that are not marked for removal. Go
101  // through the remaining faces and replace the point to remove by
102  // the equivalent point in the master face zone.
103  if (debug)
104  {
105  Pout << "Removing the cell layer for object " << name() << endl;
106  }
107 
108  const polyMesh& mesh = topoChanger().mesh();
109 
110  const labelList& ptc = pointsPairing();
111  const labelList& ftc = facesPairing();
112 
113  // Remove all the cells from the master layer
114  const labelList& mc =
115  topoChanger().mesh().faceZones()[faceZoneID_.index()].masterCells();
116 
117  forAll (mc, cellI)
118  {
119  ref.setAction(polyRemoveCell(mc[cellI]));
120  }
121 
122  // Remove all the faces from the master layer cells which are not in
123  // the master face layer
124  labelHashSet facesToRemoveMap(mc.size()*primitiveMesh::facesPerCell_);
125 
126  const cellList& cells = mesh.cells();
127 
128  forAll (mc, cellI)
129  {
130  const cell& curCell = cells[mc[cellI]];
131 
132  forAll (curCell, faceI)
133  {
134  // Check if the face is in the master zone. If not, remove it
135  if
136  (
137  mesh.faceZones().whichZone(curCell[faceI])
138  != faceZoneID_.index()
139  )
140  {
141  facesToRemoveMap.insert(curCell[faceI]);
142  }
143  }
144  }
145 
146  forAllConstIter(labelHashSet, facesToRemoveMap, iter)
147  {
148  ref.setAction(polyRemoveFace(iter.key()));
149  }
150 
151  // Remove all points that will be collapsed
152  forAll (ptc, pointI)
153  {
154  ref.setAction(polyRemovePoint(ptc[pointI]));
155  }
156 
157  // Grab all faces coming off points to be deleted. If the face
158  // has not been deleted, replace the removed point with the
159  // equivalent from the master layer.
160 
161  // Make a map of all point to be removed, giving the master point label
162  // for each of them
163 
164  Map<label> removedPointMap(2*ptc.size());
165 
166  const labelList& meshPoints =
167  mesh.faceZones()[faceZoneID_.index()]().meshPoints();
168 
169  forAll (ptc, pointI)
170  {
171  removedPointMap.insert(ptc[pointI], meshPoints[pointI]);
172  }
173 
174  const labelListList& pf = mesh.pointFaces();
175 
176  const faceList& faces = mesh.faces();
177  const labelList& own = mesh.faceOwner();
178  const labelList& nei = mesh.faceNeighbour();
179 
180  // Make a list of faces to be modified using the map to avoid duplicates
181  labelHashSet facesToModify(ptc.size()*primitiveMesh::facesPerPoint_);
182 
183  forAll (ptc, pointI)
184  {
185  const labelList& curFaces = pf[ptc[pointI]];
186 
187  forAll (curFaces, faceI)
188  {
189  if (!facesToRemoveMap.found(curFaces[faceI]))
190  {
191  facesToModify.insert(curFaces[faceI]);
192  }
193  }
194  }
195 
196  labelList ftm = facesToModify.toc();
197 
198 //Pout << "faces to modify: " << ftm << endl;
199 
200  forAll (ftm, faceI)
201  {
202  // For every face to modify, copy the face and re-map the vertices.
203  // It is known all the faces will be changed since they hang off
204  // re-mapped vertices
205  label curFaceID = ftm[faceI];
206 
207  face newFace(faces[curFaceID]);
208 
209  forAll (newFace, pointI)
210  {
211  Map<label>::iterator rpmIter =
212  removedPointMap.find(newFace[pointI]);
213 
214  if (rpmIter != removedPointMap.end())
215  {
216  // Point mapped. Replace it
217  newFace[pointI] = rpmIter();
218  }
219  }
220 
221 //Pout<< "face label: " << curFaceID
222 // << " old face: " << faces[curFaceID]
223 // << " new face: " << newFace << endl;
224 
225  // Get face zone and its flip
226  label modifiedFaceZone = mesh.faceZones().whichZone(curFaceID);
227  bool modifiedFaceZoneFlip = false;
228 
229  if (modifiedFaceZone >= 0)
230  {
231  modifiedFaceZoneFlip =
232  mesh.faceZones()[modifiedFaceZone].flipMap()
233  [
234  mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
235  ];
236  }
237 
238  label newNei;
239 
240  if (curFaceID < mesh.nInternalFaces())
241  {
242  newNei = nei[curFaceID];
243  }
244  else
245  {
246  newNei = -1;
247  }
248 
249  // Modify the face
250  ref.setAction
251  (
252  polyModifyFace
253  (
254  newFace, // modified face
255  curFaceID, // label of face being modified
256  own[curFaceID], // owner
257  newNei, // neighbour
258  false, // face flip
259  mesh.boundaryMesh().whichPatch(curFaceID),// patch for face
260  false, // remove from zone
261  modifiedFaceZone, // zone for face
262  modifiedFaceZoneFlip // face flip in zone
263  )
264  );
265  }
266 
267  // Modify the faces in the master layer to point past the removed cells
268 
269  const labelList& mf = mesh.faceZones()[faceZoneID_.index()];
270  const boolList& mfFlip = mesh.faceZones()[faceZoneID_.index()].flipMap();
271 
272  forAll (mf, faceI)
273  {
274  // Grab the owner and neighbour of the faces to be collapsed and get rid
275  // of the cell to be removed
276  label masterSideCell = own[mf[faceI]];
277 
278  if (masterSideCell == mc[faceI])
279  {
280  // Owner cell of the face is being removed.
281  // Grab the neighbour instead
282  masterSideCell = nei[mf[faceI]];
283  }
284 
285  label slaveSideCell = own[ftc[faceI]];
286 
287  if (slaveSideCell == mc[faceI])
288  {
289  // Owner cell of the face is being removed.
290  // Grab the neighbour instead
291  slaveSideCell = nei[ftc[faceI]];
292  }
293 
294  // Find out if the face needs to be flipped
295  label newOwner = -1;
296  label newNeighbour = -1;
297  bool flipFace = false;
298  label newPatchID = -1;
299  label newZoneID = -1;
300 
301  // Is any of the faces a boundary face? If so, grab the patch
302  // A boundary-to-boundary collapse is checked for in validCollapse()
303  // and cannot happen here.
304 
305  if (!mesh.isInternalFace(mf[faceI]))
306  {
307  // Master is the boundary face: it gets a new owner but no flip
308  newOwner = slaveSideCell;
309  newNeighbour = -1;
310  flipFace = false;
311  newPatchID = mesh.boundaryMesh().whichPatch(mf[faceI]);
312  newZoneID = mesh.faceZones().whichZone(mf[faceI]);
313  }
314  else if (!mesh.isInternalFace(ftc[faceI]))
315  {
316  // Slave is the boundary face: grab its patch
317  newOwner = slaveSideCell;
318  newNeighbour = -1;
319 
320  // Find out if the face flip is necessary
321  if (own[mf[faceI]] == slaveSideCell)
322  {
323  flipFace = false;
324  }
325  else
326  {
327  flipFace = true;
328  }
329 
330  newPatchID = mesh.boundaryMesh().whichPatch(ftc[faceI]);
331 
332  // The zone of the master face is preserved
333  newZoneID = mesh.faceZones().whichZone(mf[faceI]);
334  }
335  else
336  {
337  // Both faces are internal: flip is decided based on which of the
338  // new cells around it has a lower label
339  newOwner = min(masterSideCell, slaveSideCell);
340  newNeighbour = max(masterSideCell, slaveSideCell);
341 
342  if (newOwner == own[mf[faceI]] || newNeighbour == nei[mf[faceI]])
343  {
344  flipFace = false;
345  }
346  else
347  {
348  flipFace = true;
349  }
350 
351  newPatchID = -1;
352 
353  // The zone of the master face is preserved
354  newZoneID = mesh.faceZones().whichZone(mf[faceI]);
355  }
356 
357  // Modify the face and flip if necessary
358  face newFace = faces[mf[faceI]];
359  bool zoneFlip = mfFlip[faceI];
360 
361  if (flipFace)
362  {
363  newFace = newFace.reverseFace();
364  zoneFlip = !zoneFlip;
365  }
366 
367 // Pout<< "Modifying face " << mf[faceI]
368 // << " newFace: " << newFace << nl
369 // << " newOwner: " << newOwner
370 // << " newNeighbour: " << newNeighbour
371 // << " flipFace: " << flipFace
372 // << " newPatchID: " << newPatchID
373 // << " newZoneID: " << newZoneID << nl
374 // << " oldOwn: " << own[mf[faceI]]
375 // << " oldNei: " << nei[mf[faceI]] << endl;
376 
377  ref.setAction
378  (
379  polyModifyFace
380  (
381  newFace, // modified face
382  mf[faceI], // label of face being modified
383  newOwner, // owner
384  newNeighbour, // neighbour
385  flipFace, // flip
386  newPatchID, // patch for face
387  false, // remove from zone
388  newZoneID, // new zone
389  zoneFlip // face flip in zone
390  )
391  );
392  }
393 }
394 
395 
396 // ************************ vim: set sw=4 sts=4 et: ************************ //