FreeFOAM The Cross-Platform CFD Toolkit
slidingInterfaceAttachedAddressing.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 "slidingInterface.H"
27 #include <OpenFOAM/polyMesh.H>
28 #include <OpenFOAM/mapPolyMesh.H>
30 
31 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
32 
33 void Foam::slidingInterface::calcAttachedAddressing() const
34 {
35  if (debug)
36  {
37  Pout<< "void Foam::slidingInterface::calcAttachedAddressing() const "
38  << " for object " << name() << " : "
39  << "Calculating zone face-cell addressing."
40  << endl;
41  }
42 
43  if (!attached_)
44  {
45  // Clear existing addressing
46  clearAttachedAddressing();
47 
48  const polyMesh& mesh = topoChanger().mesh();
49  const labelList& own = mesh.faceOwner();
50  const labelList& nei = mesh.faceNeighbour();
51  const faceZoneMesh& faceZones = mesh.faceZones();
52 
53  // Master side
54 
55  const primitiveFacePatch& masterPatch =
56  faceZones[masterFaceZoneID_.index()]();
57 
58  const labelList& masterPatchFaces =
59  faceZones[masterFaceZoneID_.index()];
60 
61  const boolList& masterFlip =
62  faceZones[masterFaceZoneID_.index()].flipMap();
63 
64  masterFaceCellsPtr_ = new labelList(masterPatchFaces.size());
65  labelList& mfc = *masterFaceCellsPtr_;
66 
67  forAll (masterPatchFaces, faceI)
68  {
69  if (masterFlip[faceI])
70  {
71  mfc[faceI] = nei[masterPatchFaces[faceI]];
72  }
73  else
74  {
75  mfc[faceI] = own[masterPatchFaces[faceI]];
76  }
77  }
78 
79  // Slave side
80 
81  const primitiveFacePatch& slavePatch =
82  faceZones[slaveFaceZoneID_.index()]();
83 
84  const labelList& slavePatchFaces =
85  faceZones[slaveFaceZoneID_.index()];
86 
87  const boolList& slaveFlip =
88  faceZones[slaveFaceZoneID_.index()].flipMap();
89 
90  slaveFaceCellsPtr_ = new labelList(slavePatchFaces.size());
91  labelList& sfc = *slaveFaceCellsPtr_;
92 
93  forAll (slavePatchFaces, faceI)
94  {
95  if (slaveFlip[faceI])
96  {
97  sfc[faceI] = nei[slavePatchFaces[faceI]];
98  }
99  else
100  {
101  sfc[faceI] = own[slavePatchFaces[faceI]];
102  }
103  }
104 
105  // Check that the addressing is valid
106  if (min(mfc) < 0 || min(sfc) < 0)
107  {
108  if (debug)
109  {
110  forAll (mfc, faceI)
111  {
112  if (mfc[faceI] < 0)
113  {
114  Pout << "No cell next to master patch face " << faceI
115  << ". Global face no: " << mfc[faceI]
116  << " own: " << own[masterPatchFaces[faceI]]
117  << " nei: " << nei[masterPatchFaces[faceI]]
118  << " flip: " << masterFlip[faceI] << endl;
119  }
120  }
121 
122  forAll (sfc, faceI)
123  {
124  if (sfc[faceI] < 0)
125  {
126  Pout << "No cell next to slave patch face " << faceI
127  << ". Global face no: " << sfc[faceI]
128  << " own: " << own[slavePatchFaces[faceI]]
129  << " nei: " << nei[slavePatchFaces[faceI]]
130  << " flip: " << slaveFlip[faceI] << endl;
131  }
132  }
133  }
134 
136  (
137  "void slidingInterface::calcAttachedAddressing()"
138  "const"
139  ) << "Error is zone face-cell addressing. Probable error in "
140  << "decoupled mesh or sliding interface definition."
141  << abort(FatalError);
142  }
143 
144  // Calculate stick-out faces
145  const labelListList& pointFaces = mesh.pointFaces();
146 
147  // Master side
148  labelHashSet masterStickOutFaceMap
149  (
150  primitiveMesh::facesPerCell_*(masterPatch.size())
151  );
152 
153  const labelList& masterMeshPoints = masterPatch.meshPoints();
154 
155  forAll (masterMeshPoints, pointI)
156  {
157  const labelList& curFaces = pointFaces[masterMeshPoints[pointI]];
158 
159  forAll (curFaces, faceI)
160  {
161  // Check if the face belongs to the master face zone;
162  // if not add it
163  if
164  (
165  faceZones.whichZone(curFaces[faceI])
166  != masterFaceZoneID_.index()
167  )
168  {
169  masterStickOutFaceMap.insert(curFaces[faceI]);
170  }
171  }
172  }
173 
174  masterStickOutFacesPtr_ = new labelList(masterStickOutFaceMap.toc());
175 
176  // Slave side
177  labelHashSet slaveStickOutFaceMap
178  (
179  primitiveMesh::facesPerCell_*(slavePatch.size())
180  );
181 
182  const labelList& slaveMeshPoints = slavePatch.meshPoints();
183 
184  forAll (slaveMeshPoints, pointI)
185  {
186  const labelList& curFaces = pointFaces[slaveMeshPoints[pointI]];
187 
188  forAll (curFaces, faceI)
189  {
190  // Check if the face belongs to the slave face zone;
191  // if not add it
192  if
193  (
194  faceZones.whichZone(curFaces[faceI])
195  != slaveFaceZoneID_.index()
196  )
197  {
198  slaveStickOutFaceMap.insert(curFaces[faceI]);
199  }
200  }
201  }
202 
203  slaveStickOutFacesPtr_ = new labelList(slaveStickOutFaceMap.toc());
204 
205 
206  // Retired point addressing does not exist at this stage.
207  // It will be filled when the interface is coupled.
208  retiredPointMapPtr_ =
209  new Map<label>
210  (
211  2*faceZones[slaveFaceZoneID_.index()]().nPoints()
212  );
213 
214  // Ditto for cut point edge map. This is a rough guess of its size
215  //
216  cutPointEdgePairMapPtr_ =
217  new Map<Pair<edge> >
218  (
219  faceZones[slaveFaceZoneID_.index()]().nEdges()
220  );
221  }
222  else
223  {
225  (
226  "void slidingInterface::calcAttachedAddressing() const"
227  ) << "The interface is attached. The zone face-cell addressing "
228  << "cannot be assembled for object " << name()
229  << abort(FatalError);
230  }
231 
232  if (debug)
233  {
234  Pout<< "void Foam::slidingInterface::calcAttachedAddressing() const "
235  << " for object " << name() << " : "
236  << "Finished calculating zone face-cell addressing."
237  << endl;
238  }
239 }
240 
241 
242 void Foam::slidingInterface::clearAttachedAddressing() const
243 {
244  deleteDemandDrivenData(masterFaceCellsPtr_);
245  deleteDemandDrivenData(slaveFaceCellsPtr_);
246 
247  deleteDemandDrivenData(masterStickOutFacesPtr_);
248  deleteDemandDrivenData(slaveStickOutFacesPtr_);
249 
250  deleteDemandDrivenData(retiredPointMapPtr_);
251  deleteDemandDrivenData(cutPointEdgePairMapPtr_);
252 }
253 
254 
255 void Foam::slidingInterface::renumberAttachedAddressing
256 (
257  const mapPolyMesh& m
258 ) const
259 {
260  // Get reference to reverse cell renumbering
261  // The renumbering map is needed the other way around, i.e. giving
262  // the new cell number for every old cell next to the interface.
263  const labelList& reverseCellMap = m.reverseCellMap();
264 
265  const labelList& mfc = masterFaceCells();
266  const labelList& sfc = slaveFaceCells();
267 
268  // Master side
269  labelList* newMfcPtr = new labelList(mfc.size(), -1);
270  labelList& newMfc = *newMfcPtr;
271 
272  const labelList& mfzRenumber =
273  m.faceZoneFaceMap()[masterFaceZoneID_.index()];
274 
275  forAll (mfc, faceI)
276  {
277  label newCellI = reverseCellMap[mfc[mfzRenumber[faceI]]];
278 
279  if (newCellI >= 0)
280  {
281  newMfc[faceI] = newCellI;
282  }
283  }
284 
285  // Slave side
286  labelList* newSfcPtr = new labelList(sfc.size(), -1);
287  labelList& newSfc = *newSfcPtr;
288 
289  const labelList& sfzRenumber =
290  m.faceZoneFaceMap()[slaveFaceZoneID_.index()];
291 
292  forAll (sfc, faceI)
293  {
294  label newCellI = reverseCellMap[sfc[sfzRenumber[faceI]]];
295 
296  if (newCellI >= 0)
297  {
298  newSfc[faceI] = newCellI;
299  }
300  }
301 
302  if (debug)
303  {
304  // Check if all the mapped cells are live
305  if (min(newMfc) < 0 || min(newSfc) < 0)
306  {
308  (
309  "void slidingInterface::renumberAttachedAddressing("
310  "const mapPolyMesh& m) const"
311  ) << "Error in cell renumbering for object " << name()
312  << ". Some of master cells next "
313  << "to the interface have been removed."
314  << abort(FatalError);
315  }
316  }
317 
318  // Renumber stick-out faces
319 
320  const labelList& reverseFaceMap = m.reverseFaceMap();
321 
322  // Master side
323  const labelList& msof = masterStickOutFaces();
324 
325  labelList* newMsofPtr = new labelList(msof.size(), -1);
326  labelList& newMsof = *newMsofPtr;
327 
328  forAll (msof, faceI)
329  {
330  label newFaceI = reverseFaceMap[msof[faceI]];
331 
332  if (newFaceI >= 0)
333  {
334  newMsof[faceI] = newFaceI;
335  }
336  }
337 // Pout << "newMsof: " << newMsof << endl;
338  // Slave side
339  const labelList& ssof = slaveStickOutFaces();
340 
341  labelList* newSsofPtr = new labelList(ssof.size(), -1);
342  labelList& newSsof = *newSsofPtr;
343 
344  forAll (ssof, faceI)
345  {
346  label newFaceI = reverseFaceMap[ssof[faceI]];
347 
348  if (newFaceI >= 0)
349  {
350  newSsof[faceI] = newFaceI;
351  }
352  }
353 // Pout << "newSsof: " << newSsof << endl;
354  if (debug)
355  {
356  // Check if all the mapped cells are live
357  if (min(newMsof) < 0 || min(newSsof) < 0)
358  {
360  (
361  "void slidingInterface::renumberAttachedAddressing("
362  "const mapPolyMesh& m) const"
363  ) << "Error in face renumbering for object " << name()
364  << ". Some of stick-out next "
365  << "to the interface have been removed."
366  << abort(FatalError);
367  }
368  }
369 
370  // Renumber the retired point map. Need to take a copy!
371  const Map<label> rpm = retiredPointMap();
372 
373  Map<label>* newRpmPtr = new Map<label>(rpm.size());
374  Map<label>& newRpm = *newRpmPtr;
375 
376  const labelList rpmToc = rpm.toc();
377 
378  // Get reference to point renumbering
379  const labelList& reversePointMap = m.reversePointMap();
380 
381  label key, value;
382 
383  forAll (rpmToc, rpmTocI)
384  {
385  key = reversePointMap[rpmToc[rpmTocI]];
386 
387  value = reversePointMap[rpm.find(rpmToc[rpmTocI])()];
388 
389  if (debug)
390  {
391  // Check if all the mapped cells are live
392  if (key < 0 || value < 0)
393  {
395  (
396  "void slidingInterface::renumberAttachedAddressing("
397  "const mapPolyMesh& m) const"
398  ) << "Error in retired point numbering for object "
399  << name() << ". Some of master "
400  << "points have been removed."
401  << abort(FatalError);
402  }
403  }
404 
405  newRpm.insert(key, value);
406  }
407 
408  // Renumber the cut point edge pair map. Need to take a copy!
409  const Map<Pair<edge> > cpepm = cutPointEdgePairMap();
410 
411  Map<Pair<edge> >* newCpepmPtr = new Map<Pair<edge> >(cpepm.size());
412  Map<Pair<edge> >& newCpepm = *newCpepmPtr;
413 
414  const labelList cpepmToc = cpepm.toc();
415 
416  forAll (cpepmToc, cpepmTocI)
417  {
418  key = reversePointMap[cpepmToc[cpepmTocI]];
419 
420  const Pair<edge>& oldPe = cpepm.find(cpepmToc[cpepmTocI])();
421 
422  // Re-do the edges in global addressing
423  const label ms = reversePointMap[oldPe.first().start()];
424  const label me = reversePointMap[oldPe.first().end()];
425 
426  const label ss = reversePointMap[oldPe.second().start()];
427  const label se = reversePointMap[oldPe.second().end()];
428 
429  if (debug)
430  {
431  // Check if all the mapped cells are live
432  if (key < 0 || ms < 0 || me < 0 || ss < 0 || se < 0)
433  {
435  (
436  "void slidingInterface::renumberAttachedAddressing("
437  "const mapPolyMesh& m) const"
438  ) << "Error in cut point edge pair map numbering for object "
439  << name() << ". Some of master points have been removed."
440  << abort(FatalError);
441  }
442  }
443 
444  newCpepm.insert(key, Pair<edge>(edge(ms, me), edge(ss, se)));
445  }
446 
447  if (!projectedSlavePointsPtr_)
448  {
450  (
451  "void slidingInterface::renumberAttachedAddressing("
452  "const mapPolyMesh& m) const"
453  ) << "Error in projected point numbering for object " << name()
454  << abort(FatalError);
455  }
456 
457  // Renumber the projected slave zone points
458  const pointField& projectedSlavePoints = *projectedSlavePointsPtr_;
459 
460  pointField* newProjectedSlavePointsPtr
461  (
462  new pointField(projectedSlavePoints.size())
463  );
464  pointField& newProjectedSlavePoints = *newProjectedSlavePointsPtr;
465 
466  const labelList& sfzPointRenumber =
467  m.faceZonePointMap()[slaveFaceZoneID_.index()];
468 
469  forAll (newProjectedSlavePoints, pointI)
470  {
471  if (sfzPointRenumber[pointI] > -1)
472  {
473  newProjectedSlavePoints[pointI] =
474  projectedSlavePoints[sfzPointRenumber[pointI]];
475  }
476  }
477 
478  // Re-set the lists
479  clearAttachedAddressing();
480 
481  deleteDemandDrivenData(projectedSlavePointsPtr_);
482 
483  masterFaceCellsPtr_ = newMfcPtr;
484  slaveFaceCellsPtr_ = newSfcPtr;
485 
486  masterStickOutFacesPtr_ = newMsofPtr;
487  slaveStickOutFacesPtr_ = newSsofPtr;
488 
489  retiredPointMapPtr_ = newRpmPtr;
490  cutPointEdgePairMapPtr_ = newCpepmPtr;
491  projectedSlavePointsPtr_ = newProjectedSlavePointsPtr;
492 }
493 
494 
495 const Foam::labelList& Foam::slidingInterface::masterFaceCells() const
496 {
497  if (!masterFaceCellsPtr_)
498  {
500  (
501  "const labelList& slidingInterface::masterFaceCells() const"
502  ) << "Master zone face-cell addressing not available for object "
503  << name()
504  << abort(FatalError);
505  }
506 
507  return *masterFaceCellsPtr_;
508 }
509 
510 
511 const Foam::labelList& Foam::slidingInterface::slaveFaceCells() const
512 {
513  if (!slaveFaceCellsPtr_)
514  {
516  (
517  "const labelList& slidingInterface::slaveFaceCells() const"
518  ) << "Slave zone face-cell addressing not available for object "
519  << name()
520  << abort(FatalError);
521  }
522 
523  return *slaveFaceCellsPtr_;
524 }
525 
526 
527 const Foam::labelList& Foam::slidingInterface::masterStickOutFaces() const
528 {
529  if (!masterStickOutFacesPtr_)
530  {
532  (
533  "const labelList& slidingInterface::masterStickOutFaces() const"
534  ) << "Master zone stick-out face addressing not available for object "
535  << name()
536  << abort(FatalError);
537  }
538 
539  return *masterStickOutFacesPtr_;
540 }
541 
542 
543 const Foam::labelList& Foam::slidingInterface::slaveStickOutFaces() const
544 {
545  if (!slaveStickOutFacesPtr_)
546  {
548  (
549  "const labelList& slidingInterface::slaveStickOutFaces() const"
550  ) << "Slave zone stick-out face addressing not available for object "
551  << name()
552  << abort(FatalError);
553  }
554 
555  return *slaveStickOutFacesPtr_;
556 }
557 
558 
559 const Foam::Map<Foam::label>& Foam::slidingInterface::retiredPointMap() const
560 {
561  if (!retiredPointMapPtr_)
562  {
564  (
565  "const Map<label>& slidingInterface::retiredPointMap() const"
566  ) << "Retired point map not available for object " << name()
567  << abort(FatalError);
568  }
569 
570  return *retiredPointMapPtr_;
571 }
572 
573 
575 Foam::slidingInterface::cutPointEdgePairMap() const
576 {
577  if (!cutPointEdgePairMapPtr_)
578  {
580  (
581  "const Map<Pair<edge> >& slidingInterface::"
582  "cutPointEdgePairMap() const"
583  ) << "Retired point map not available for object " << name()
584  << abort(FatalError);
585  }
586 
587  return *cutPointEdgePairMapPtr_;
588 }
589 
590 
591 // ************************ vim: set sw=4 sts=4 et: ************************ //