src/dynamicMesh/slidingInterface/decoupleSlidingInterface.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 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 \*---------------------------------------------------------------------------*/
  25
  26 #include "slidingInterface.H"
  27 #include "polyMesh.H"
  28 #include "primitiveMesh.H"
  29 #include "polyTopoChange.H"
  30 #include "polyTopoChanger.H"
  31 #include "polyModifyFace.H"
  32 #include "polyModifyPoint.H"
  33
  34 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  35
  36 void Foam::slidingInterface::decoupleInterface
  37 (
  38     polyTopoChange& ref
  39 ) const
  40 {
  41     if (debug)
  42     {
  43         Pout<< "void slidingInterface::decoupleInterface("
  44             << "polyTopoChange& ref) const : "
  45             << "Decoupling sliding interface " << name() << endl;
  46     }
  47
  48     if (!attached_)
  49     {
  50         if (debug)
  51         {
  52             Pout<< "void slidingInterface::decoupleInterface("
  53                 << "polyTopoChange& ref) const : "
  54                 << "Interface already decoupled." << endl;
  55         }
  56
  57         return;
  58     }
  59
  60     // Clear previous couple
  61     clearCouple(ref);
  62
  63     const polyMesh& mesh = topoChanger().mesh();
  64     const faceList& faces = mesh.faces();
  65     const cellList& cells = mesh.cells();
  66
  67     const labelList& own = mesh.faceOwner();
  68     const labelList& nei = mesh.faceNeighbour();
  69
  70     // Master side
  71
  72     const primitiveFacePatch& masterPatch =
  73         mesh.faceZones()[masterFaceZoneID_.index()]();
  74
  75     const labelList& masterPatchAddr =
  76         mesh.faceZones()[masterFaceZoneID_.index()];
  77
  78     const boolList& masterPatchFlip =
  79         mesh.faceZones()[masterFaceZoneID_.index()].flipMap();
  80
  81     const labelList& masterFc = masterFaceCells();
  82
  83     // Recover faces in master patch
  84
  85     forAll(masterPatchAddr, faceI)
  86     {
  87         // Make a copy of the face and turn it if necessary
  88         face newFace = faces[masterPatchAddr[faceI]];
  89
  90         if (masterPatchFlip[faceI])
  91         {
  92             newFace.flip();
  93         }
  94
  95         ref.setAction
  96         (
  97             polyModifyFace
  98             (
  99                 newFace,                         // new face
 100                 masterPatchAddr[faceI],          // master face index
 101                 masterFc[faceI],                 // owner
 102                 -1,                              // neighbour
 103                 false,                           // flux flip
 104                 masterPatchID_.index(),          // patch ID
 105                 false,                           // remove from zone
 106                 masterFaceZoneID_.index(),       // zone ID
 107                 false                            // zone flip.  Face corrected
 108             )
 109         );
 110
 111         // Pout<< "Modifying master patch face no "
 112         //     << masterPatchAddr[faceI]
 113         //     << " face: " << faces[masterPatchAddr[faceI]]
 114         //     << " old owner: " << own[masterPatchAddr[faceI]]
 115         //     << " new owner: " << masterFc[faceI]
 116         //     << endl;
 117     }
 118
 119     // Slave side
 120
 121     const primitiveFacePatch& slavePatch =
 122         mesh.faceZones()[slaveFaceZoneID_.index()]();
 123
 124     const labelList& slavePatchAddr =
 125         mesh.faceZones()[slaveFaceZoneID_.index()];
 126
 127     const boolList& slavePatchFlip =
 128         mesh.faceZones()[slaveFaceZoneID_.index()].flipMap();
 129
 130     const labelList& slaveFc = slaveFaceCells();
 131
 132     // Grab retired point mapping
 133     const Map<label>& rpm = retiredPointMap();
 134
 135     // Recover faces in slave patch
 136
 137     forAll(slavePatchAddr, faceI)
 138     {
 139         // Make a copy of face and turn it if necessary
 140         face newFace = faces[slavePatchAddr[faceI]];
 141
 142         if (slavePatchFlip[faceI])
 143         {
 144             newFace.flip();
 145         }
 146
 147         // Recover retired points on the slave side
 148         forAll(newFace, pointI)
 149         {
 150             Map<label>::const_iterator rpmIter = rpm.find(newFace[pointI]);
 151             if (rpmIter != rpm.end())
 152             {
 153                 // Master of retired point; grab its original
 154                 // Pout<< "Reinstating retired point: " << newFace[pointI]
 155                 //     << " with old: " << rpm.find(newFace[pointI])()
 156                 //     << endl;
 157
 158                 newFace[pointI] = rpmIter();
 159             }
 160         }
 161
 162         ref.setAction
 163         (
 164             polyModifyFace
 165             (
 166                 newFace,                         // new face
 167                 slavePatchAddr[faceI],           // master face index
 168                 slaveFc[faceI],                  // owner
 169                 -1,                              // neighbour
 170                 false,                           // flux flip
 171                 slavePatchID_.index(),           // patch ID
 172                 false,                           // remove from zone
 173                 slaveFaceZoneID_.index(),        // zone ID
 174                 false                            // zone flip.  Face corrected
 175             )
 176         );
 177     }
 178
 179     // Re-create the master stick-out faces
 180
 181     // Grab the list of faces in the layer
 182     const labelList& masterStickOuts = masterStickOutFaces();
 183
 184     forAll(masterStickOuts, faceI)
 185     {
 186         // Renumber the face and remove additional points
 187
 188         const label curFaceID = masterStickOuts[faceI];
 189
 190         const face& oldFace = faces[curFaceID];
 191
 192         DynamicList<label> newFaceLabels(oldFace.size());
 193
 194         bool changed = false;
 195
 196         forAll(oldFace, pointI)
 197         {
 198             // Check if the point is removed
 199             if (ref.pointRemoved(oldFace[pointI]))
 200             {
 201                 // Point removed; skip it
 202                 changed = true;
 203             }
 204             else
 205             {
 206                 newFaceLabels.append(oldFace[pointI]);
 207             }
 208         }
 209
 210         if (changed)
 211         {
 212             if (newFaceLabels.size() < 3)
 213             {
 214                 FatalErrorIn
 215                 (
 216                     "void slidingInterface::decoupleInterface("
 217                     "polyTopoChange& ref) const"
 218                 )   << "Face " << curFaceID << " reduced to less than "
 219                     << "3 points.  Topological/cutting error." << nl
 220                     << "Old face: " << oldFace << " new face: " << newFaceLabels
 221                     << abort(FatalError);
 222             }
 223
 224             // Get face zone and its flip
 225             label modifiedFaceZone = mesh.faceZones().whichZone(curFaceID);
 226             bool modifiedFaceZoneFlip = false;
 227
 228             if (modifiedFaceZone >= 0)
 229             {
 230                 modifiedFaceZoneFlip =
 231                     mesh.faceZones()[modifiedFaceZone].flipMap()
 232                     [
 233                         mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
 234                     ];
 235             }
 236
 237             face newFace;
 238             newFace.transfer(newFaceLabels);
 239
 240             // Pout<< "Modifying master stick-out face " << curFaceID
 241             //     << " old face: " << oldFace
 242             //     << " new face: " << newFace
 243             //     << endl;
 244
 245             // Modify the face
 246             ref.setAction
 247             (
 248                 polyModifyFace
 249                 (
 250                     newFace,                // modified face
 251                     curFaceID,              // label of face being modified
 252                     own[curFaceID],         // owner
 253                     nei[curFaceID],         // neighbour
 254                     false,                  // face flip
 255                     mesh.boundaryMesh().whichPatch(curFaceID), // patch for face
 256                     false,                  // remove from zone
 257                     modifiedFaceZone,       // zone for face
 258                     modifiedFaceZoneFlip    // face flip in zone
 259                 )
 260             );
 261         }
 262     }
 263
 264     // Re-create the slave stick-out faces
 265
 266     labelHashSet slaveLayerCellFaceMap
 267     (
 268         primitiveMesh::facesPerCell_*(masterPatch.size() + slavePatch.size())
 269     );
 270
 271     forAll(slaveFc, faceI)
 272     {
 273         const labelList& curFaces = cells[slaveFc[faceI]];
 274
 275         forAll(curFaces, faceI)
 276         {
 277             // Check if the face belongs to the slave face zone; and
 278             // if it has been removed; if not add it
 279             if
 280             (
 281                 mesh.faceZones().whichZone(curFaces[faceI])
 282              != slaveFaceZoneID_.index()
 283              && !ref.faceRemoved(curFaces[faceI])
 284
 285             )
 286             {
 287                 slaveLayerCellFaceMap.insert(curFaces[faceI]);
 288             }
 289         }
 290     }
 291
 292     // Grab the list of faces in the layer
 293     const labelList& slaveStickOuts = slaveStickOutFaces();
 294
 295     // Grab master point mapping
 296     const Map<label>& masterPm = masterPatch.meshPointMap();
 297
 298     forAll(slaveStickOuts, faceI)
 299     {
 300         // Renumber the face and remove additional points
 301
 302         const label curFaceID = slaveStickOuts[faceI];
 303
 304         const face& oldFace = faces[curFaceID];
 305
 306         DynamicList<label> newFaceLabels(oldFace.size());
 307
 308         bool changed = false;
 309
 310         forAll(oldFace, pointI)
 311         {
 312             // Check if the point is removed or retired
 313             if (rpm.found(oldFace[pointI]))
 314             {
 315                 // Master of retired point; grab its original
 316                 changed = true;
 317
 318                 // Pout<< "Reinstating retired point: " << oldFace[pointI]
 319                 //     << " with old: " << rpm.find(oldFace[pointI])()
 320                 //     << endl;
 321
 322                 newFaceLabels.append(rpm.find(oldFace[pointI])());
 323             }
 324             else if (ref.pointRemoved(oldFace[pointI]))
 325             {
 326                 // Point removed; skip it
 327                 changed = true;
 328             }
 329             else if (masterPm.found(oldFace[pointI]))
 330             {
 331                 // Point from master patch only; skip it
 332                 changed = true;
 333             }
 334             else
 335             {
 336                 newFaceLabels.append(oldFace[pointI]);
 337             }
 338         }
 339
 340         if (changed)
 341         {
 342             if (newFaceLabels.size() < 3)
 343             {
 344                 FatalErrorIn
 345                 (
 346                     "void slidingInterface::decoupleInterface("
 347                     "polyTopoChange& ref) const"
 348                 )   << "Face " << curFaceID << " reduced to less than "
 349                     << "3 points.  Topological/cutting error." << nl
 350                     << "Old face: " << oldFace << " new face: " << newFaceLabels
 351                     << abort(FatalError);
 352             }
 353
 354             // Get face zone and its flip
 355             label modifiedFaceZone = mesh.faceZones().whichZone(curFaceID);
 356             bool modifiedFaceZoneFlip = false;
 357
 358             if (modifiedFaceZone >= 0)
 359             {
 360                 modifiedFaceZoneFlip =
 361                     mesh.faceZones()[modifiedFaceZone].flipMap()
 362                     [
 363                         mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
 364                     ];
 365             }
 366
 367             face newFace;
 368             newFace.transfer(newFaceLabels);
 369
 370             // Pout<< "Modifying slave stick-out face " << curFaceID
 371             //     << " old face: " << oldFace
 372             //     << " new face: " << newFace
 373             //     << endl;
 374
 375             // Modify the face
 376             ref.setAction
 377             (
 378                 polyModifyFace
 379                 (
 380                     newFace,                // modified face
 381                     curFaceID,              // label of face being modified
 382                     own[curFaceID],         // owner
 383                     nei[curFaceID],         // neighbour
 384                     false,                  // face flip
 385                     mesh.boundaryMesh().whichPatch(curFaceID), // patch for face
 386                     false,                  // remove from zone
 387                     modifiedFaceZone,       // zone for face
 388                     modifiedFaceZoneFlip    // face flip in zone
 389                 )
 390             );
 391         }
 392     }
 393
 394     // Bring all slave patch points back to life
 395     const pointField& points = mesh.points();
 396
 397     const labelList& slaveMeshPoints =
 398         mesh.faceZones()[slaveFaceZoneID_.index()]().meshPoints();
 399
 400     forAll(slaveMeshPoints, pointI)
 401     {
 402         ref.setAction
 403         (
 404             polyModifyPoint
 405             (
 406                 slaveMeshPoints[pointI],             // point ID
 407                 points[slaveMeshPoints[pointI]],     // point
 408                 false,                               // remove from zone
 409                 mesh.pointZones().whichZone(slaveMeshPoints[pointI]), // zone
 410                 true                                // in a cell
 411             )
 412         );
 413     }
 414
 415     // Clear the retired point numbering
 416     retiredPointMapPtr_->clear();
 417
 418     // Finished decoupling
 419     attached_ = false;
 420
 421     if (debug)
 422     {
 423         Pout<< "void slidingInterface::coupleInterface("
 424             << "polyTopoChange& ref) const : "
 425             << "Finished decoupling sliding interface " << name() << endl;
 426     }
 427 }
 428
 429
 430 // ************************************************************************* //