FreeFOAM The Cross-Platform CFD Toolkit
detachInterface.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 \*---------------------------------------------------------------------------*/
25 
26 #include "attachDetach.H"
27 #include <OpenFOAM/polyMesh.H>
28 #include <OpenFOAM/primitiveMesh.H>
34 
35 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
36 
37 void Foam::attachDetach::detachInterface
38 (
39  polyTopoChange& ref
40 ) const
41 {
42  // Algorithm:
43  // 1. Create new points for points of the master face zone
44  // 2. Modify all faces of the master zone, by putting them into the master
45  // patch (look for orientation) and their renumbered mirror images
46  // into the slave patch
47  // 3. Create a point renumbering list, giving a new point index for original
48  // points in the face patch
49  // 4. Grab all faces in cells on the master side and renumber them
50  // using the point renumbering list. Exclude the ones that belong to
51  // the master face zone
52  //
53  // Note on point creation:
54  // In order to take into account the issues related to partial
55  // blocking in an attach/detach mesh modifier, special treatment
56  // is required for the duplication of points on the edge of the
57  // face zone. Points which are shared only by internal edges need
58  // not to be duplicated, as this would propagate the discontinuity
59  // in the mesh beyond the face zone. Therefore, before creating
60  // the new points, check the external edge loop. For each edge
61  // check if the edge is internal (i.e. does not belong to a
62  // patch); if so, exclude both of its points from duplication.
63 
64  if (debug)
65  {
66  Pout<< "void attachDetach::detachInterface("
67  << "polyTopoChange& ref) const "
68  << " for object " << name() << " : "
69  << "Detaching interface" << endl;
70  }
71 
72  const polyMesh& mesh = topoChanger().mesh();
73  const faceZoneMesh& zoneMesh = mesh.faceZones();
74 
75  // Check that zone is in increasing order (needed since adding faces
76  // in same order - otherwise polyTopoChange face ordering will mess up
77  // correspondence)
78  if (debug)
79  {
80  const labelList& faceLabels = zoneMesh[faceZoneID_.index()];
81  if (faceLabels.size() > 0)
82  {
83  for (label i = 1; i < faceLabels.size(); i++)
84  {
85  if (faceLabels[i] <= faceLabels[i-1])
86  {
88  (
89  "attachDetach::detachInterface"
90  "(polyTopoChange&) const"
91  ) << "faceZone " << zoneMesh[faceZoneID_.index()].name()
92  << " does not have mesh face labels in"
93  << " increasing order." << endl
94  << "Face label " << faceLabels[i]
95  << " at position " << i
96  << " is smaller than the previous value "
97  << faceLabels[i-1]
98  << exit(FatalError);
99  }
100  }
101  }
102  }
103 
104 
105 
106  const primitiveFacePatch& masterFaceLayer = zoneMesh[faceZoneID_.index()]();
107  const pointField& points = mesh.points();
108  const labelListList& meshEdgeFaces = mesh.edgeFaces();
109 
110  const labelList& mp = masterFaceLayer.meshPoints();
111  const edgeList& zoneLocalEdges = masterFaceLayer.edges();
112 
113  const labelList& meshEdges = zoneMesh[faceZoneID_.index()].meshEdges();
114 
115  // Create the points
116 
117  labelList addedPoints(mp.size(), -1);
118 
119  // Go through boundary edges of the master patch. If all the faces from
120  // this patch are internal, mark the points in the addedPoints lookup
121  // with their original labels to stop duplication
122  label nIntEdges = masterFaceLayer.nInternalEdges();
123 
124  for (label curEdgeID = nIntEdges; curEdgeID < meshEdges.size(); curEdgeID++)
125  {
126  const labelList& curFaces = meshEdgeFaces[meshEdges[curEdgeID]];
127 
128  bool edgeIsInternal = true;
129 
130  forAll (curFaces, faceI)
131  {
132  if (!mesh.isInternalFace(curFaces[faceI]))
133  {
134  // The edge belongs to a boundary face
135  edgeIsInternal = false;
136  break;
137  }
138  }
139 
140  if (edgeIsInternal)
141  {
142  const edge& e = zoneLocalEdges[curEdgeID];
143 
144  // Reset the point creation
145  addedPoints[e.start()] = mp[e.start()];
146  addedPoints[e.end()] = mp[e.end()];
147  }
148  }
149 // Pout << "addedPoints before point creation: " << addedPoints << endl;
150 
151  // Create new points for face zone
152  forAll (addedPoints, pointI)
153  {
154  if (addedPoints[pointI] < 0)
155  {
156  addedPoints[pointI] =
157  ref.setAction
158  (
159  polyAddPoint
160  (
161  points[mp[pointI]], // point
162  mp[pointI], // master point
163  -1, // zone ID
164  true // supports a cell
165  )
166  );
167  //Pout<< "Adding point " << addedPoints[pointI]
168  // << " coord1:" << points[mp[pointI]]
169  // << " coord2:" << masterFaceLayer.localPoints()[pointI]
170  // << " for original point " << mp[pointI] << endl;
171  }
172  }
173 
174  // Modify faces in the master zone and duplicate for the slave zone
175 
176  const labelList& mf = zoneMesh[faceZoneID_.index()];
177  const boolList& mfFlip = zoneMesh[faceZoneID_.index()].flipMap();
178  const faceList& zoneFaces = masterFaceLayer.localFaces();
179 
180  const faceList& faces = mesh.faces();
181  const labelList& own = mesh.faceOwner();
182  const labelList& nei = mesh.faceNeighbour();
183 
184  forAll (mf, faceI)
185  {
186  const label curFaceID = mf[faceI];
187 
188  // Build the face for the slave patch by renumbering
189  const face oldFace = zoneFaces[faceI].reverseFace();
190 
191  face newFace(oldFace.size());
192 
193  forAll (oldFace, pointI)
194  {
195  newFace[pointI] = addedPoints[oldFace[pointI]];
196  }
197 
198  if (mfFlip[faceI])
199  {
200  // Face needs to be flipped for the master patch
201  ref.setAction
202  (
203  polyModifyFace
204  (
205  faces[curFaceID].reverseFace(), // modified face
206  curFaceID, // label of face being modified
207  nei[curFaceID], // owner
208  -1, // neighbour
209  true, // face flip
210  masterPatchID_.index(), // patch for face
211  false, // remove from zone
212  faceZoneID_.index(), // zone for face
213  !mfFlip[faceI] // face flip in zone
214  )
215  );
216 
217  // Add renumbered face into the slave patch
218  //label addedFaceI =
219  ref.setAction
220  (
221  polyAddFace
222  (
223  newFace, // face
224  own[curFaceID], // owner
225  -1, // neighbour
226  -1, // master point
227  -1, // master edge
228  curFaceID, // master face
229  false, // flip flux
230  slavePatchID_.index(), // patch to add the face to
231  -1, // zone for face
232  false // zone flip
233  )
234  );
235  //{
236  // pointField newPts(ref.points());
237  //Pout<< "Flip. Modifying face: " << ref.faces()[curFaceID]
238  // << " fc:" << ref.faces()[curFaceID].centre(newPts)
239  // << " next to cell: " << nei[curFaceID]
240  // << " and adding face: " << newFace
241  // << " fc:" << ref.faces()[addedFaceI].centre(newPts)
242  // << " next to cell: " << own[curFaceID] << endl;
243  //}
244  }
245  else
246  {
247  // No flip
248  ref.setAction
249  (
250  polyModifyFace
251  (
252  faces[curFaceID], // modified face
253  curFaceID, // label of face being modified
254  own[curFaceID], // owner
255  -1, // neighbour
256  false, // face flip
257  masterPatchID_.index(), // patch for face
258  false, // remove from zone
259  faceZoneID_.index(), // zone for face
260  mfFlip[faceI] // face flip in zone
261  )
262  );
263 
264  // Add renumbered face into the slave patch
265  //label addedFaceI =
266  ref.setAction
267  (
268  polyAddFace
269  (
270  newFace, // face
271  nei[curFaceID], // owner
272  -1, // neighbour
273  -1, // master point
274  -1, // master edge
275  curFaceID, // master face
276  true, // flip flux
277  slavePatchID_.index(), // patch to add the face to
278  -1, // zone for face
279  false // face flip in zone
280  )
281  );
282  //{
283  // pointField newPts(ref.points());
284  //Pout<< "No flip. Modifying face: " << ref.faces()[curFaceID]
285  // << " fc:" << ref.faces()[curFaceID].centre(newPts)
286  // << " next to cell: " << own[curFaceID]
287  // << " and adding face: " << newFace
288  // << " fc:" << ref.faces()[addedFaceI].centre(newPts)
289  // << " next to cell: " << nei[curFaceID] << endl;
290  //}
291  }
292  }
293 
294  // Modify the remaining faces of the master cells to reconnect to the new
295  // layer of faces.
296 
297  // Algorithm: Go through all the cells of the master zone and make
298  // a map of faces to avoid duplicates. Do not insert the faces in
299  // the master patch (as they have already been dealt with). Make
300  // a master layer point renumbering map, which for every point in
301  // the master layer gives its new label. Loop through all faces in
302  // the map and attempt to renumber them using the master layer
303  // point renumbering map. Once the face is renumbered, compare it
304  // with the original face; if they are the same, the face has not
305  // changed; if not, modify the face but keep all of its old
306  // attributes (apart from the vertex numbers).
307 
308  // Create the map of faces in the master cell layer
309  const labelList& mc =
310  mesh.faceZones()[faceZoneID_.index()].masterCells();
311 
312  labelHashSet masterCellFaceMap(6*mc.size());
313 
314  const cellList& cells = mesh.cells();
315 
316  forAll (mc, cellI)
317  {
318  const labelList& curFaces = cells[mc[cellI]];
319 
320  forAll (curFaces, faceI)
321  {
322  // Check if the face belongs to the master patch; if not add it
323  if (zoneMesh.whichZone(curFaces[faceI]) != faceZoneID_.index())
324  {
325  masterCellFaceMap.insert(curFaces[faceI]);
326  }
327  }
328  }
329 
330  // Extend the map to include first neighbours of the master cells to
331  // deal with multiple corners.
332  { // Protection and memory management
333  // Make a map of master cells for quick reject
334  labelHashSet mcMap(2*mc.size());
335 
336  forAll (mc, mcI)
337  {
338  mcMap.insert(mc[mcI]);
339  }
340 
341  // Go through all the faces in the masterCellFaceMap. If the
342  // cells around them are not already used, add all of their
343  // faces to the map
344  const labelList mcf = masterCellFaceMap.toc();
345 
346  forAll (mcf, mcfI)
347  {
348  // Do the owner side
349  const label ownCell = own[mcf[mcfI]];
350 
351  if (!mcMap.found(ownCell))
352  {
353  // Cell not found. Add its faces to the map
354  const cell& curFaces = cells[ownCell];
355 
356  forAll (curFaces, faceI)
357  {
358  masterCellFaceMap.insert(curFaces[faceI]);
359  }
360  }
361 
362  // Do the neighbour side if face is internal
363  if (mesh.isInternalFace(mcf[mcfI]))
364  {
365  const label neiCell = nei[mcf[mcfI]];
366 
367  if (!mcMap.found(neiCell))
368  {
369  // Cell not found. Add its faces to the map
370  const cell& curFaces = cells[neiCell];
371 
372  forAll (curFaces, faceI)
373  {
374  masterCellFaceMap.insert(curFaces[faceI]);
375  }
376  }
377  }
378  }
379  }
380 
381  // Create the master layer point map
382  Map<label> masterLayerPointMap(2*mp.size());
383 
384  forAll (mp, pointI)
385  {
386  masterLayerPointMap.insert
387  (
388  mp[pointI],
389  addedPoints[pointI]
390  );
391  }
392 
393  // Grab the list of faces of the master layer
394  const labelList masterCellFaces = masterCellFaceMap.toc();
395 
396  forAll (masterCellFaces, faceI)
397  {
398  // Attempt to renumber the face using the masterLayerPointMap.
399  // Missing point remain the same
400 
401  const label curFaceID = masterCellFaces[faceI];
402 
403  const face& oldFace = faces[curFaceID];
404 
405  face newFace(oldFace.size());
406 
407  bool changed = false;
408 
409  forAll (oldFace, pointI)
410  {
411  if (masterLayerPointMap.found(oldFace[pointI]))
412  {
413  changed = true;
414 
415  newFace[pointI] = masterLayerPointMap.find(oldFace[pointI])();
416  }
417  else
418  {
419  newFace[pointI] = oldFace[pointI];
420  }
421  }
422 
423  // If the face has changed, create a modification entry
424  if (changed)
425  {
426  if (mesh.isInternalFace(curFaceID))
427  {
428  ref.setAction
429  (
430  polyModifyFace
431  (
432  newFace, // face
433  curFaceID, // master face
434  own[curFaceID], // owner
435  nei[curFaceID], // neighbour
436  false, // flip flux
437  -1, // patch for face
438  false, // remove from zone
439  -1, // zone for face
440  false // face zone flip
441  )
442  );
443 // Pout << "modifying stick-out face. Internal Old face: " << oldFace << " new face: " << newFace << " own: " << own[curFaceID] << " nei: " << nei[curFaceID] << endl;
444  }
445  else
446  {
447  ref.setAction
448  (
449  polyModifyFace
450  (
451  newFace, // face
452  curFaceID, // master face
453  own[curFaceID], // owner
454  -1, // neighbour
455  false, // flip flux
456  mesh.boundaryMesh().whichPatch(curFaceID), // patch
457  false, // remove from zone
458  -1, // zone for face
459  false // face zone flip
460  )
461  );
462 // Pout << "modifying stick-out face. Boundary Old face: " << oldFace << " new face: " << newFace << " own: " << own[curFaceID] << " patch: " << mesh.boundaryMesh().whichPatch(curFaceID) << endl;
463  }
464  }
465  }
466 
467  if (debug)
468  {
469  Pout<< "void attachDetach::detachInterface("
470  << "polyTopoChange& ref) const "
471  << " for object " << name() << " : "
472  << "Finished detaching interface" << endl;
473  }
474 }
475 
476 
477 // ************************ vim: set sw=4 sts=4 et: ************************ //