src/mesh/cfMesh/meshLibrary/utilities/meshes/polyMeshGenModifier/polyMeshGenModifierRemoveFaces.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | cfMesh: A library for mesh generation
   4    \\    /   O peration     |
   5     \\  /    A nd           | Author: Franjo Juretic (franjo.juretic@c-fields.com)
   6      \\/     M anipulation  | Copyright (C) Creative Fields, Ltd.
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of cfMesh.
  10
  11     cfMesh is free software; you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation; either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     cfMesh 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 cfMesh.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25
  26 \*---------------------------------------------------------------------------*/
  27
  28 #include "polyMeshGenModifier.H"
  29 #include "demandDrivenData.H"
  30
  31 # ifdef USE_OMP
  32 #include <omp.h>
  33 # endif
  34
  35 namespace Foam
  36 {
  37
  38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  39
  40 void polyMeshGenModifier::removeFaces(const boolList& removeFace)
  41 {
  42     //mesh_.clearOut();
  43
  44     Info << "Removing faces" << endl;
  45     faceListPMG& faces = mesh_.faces_;
  46     cellListPMG& cells = mesh_.cells_;
  47
  48     label nFaces(0);
  49     labelLongList newFaceLabel(faces.size(), -1);
  50
  51     //- copy internal faces
  52     const label nInternalFaces = mesh_.nInternalFaces();
  53     for(label faceI=0;faceI<nInternalFaces;++faceI)
  54         if( !removeFace[faceI] )
  55         {
  56             if( nFaces < faceI )
  57                 faces[nFaces].transfer(faces[faceI]);
  58
  59             newFaceLabel[faceI] = nFaces;
  60             ++nFaces;
  61         }
  62
  63     //- copy boundary faces
  64     PtrList<boundaryPatch>& boundaries = mesh_.boundaries_;
  65     labelList patchStart(boundaries.size());
  66     labelList nFacesInPatch(boundaries.size());
  67     label npI(0);
  68     forAll(boundaries, patchI)
  69     {
  70         const label oldStart = boundaries[patchI].patchStart();
  71         const label oldNumFacesInPatch = boundaries[patchI].patchSize();
  72
  73         patchStart[npI] = nFaces;
  74         nFacesInPatch[npI] = 0;
  75
  76         for(label faceI=0;faceI<oldNumFacesInPatch;++faceI)
  77             if( !removeFace[oldStart+faceI] )
  78             {
  79                 ++nFacesInPatch[npI];
  80                 if( nFaces < (oldStart+faceI) )
  81                     faces[nFaces].transfer(faces[oldStart+faceI]);
  82                 newFaceLabel[oldStart+faceI] = nFaces++;
  83             }
  84
  85         ++npI;
  86     }
  87
  88     forAll(boundaries, patchI)
  89     {
  90         boundaries[patchI].patchStart() = patchStart[patchI];
  91         boundaries[patchI].patchSize() = nFacesInPatch[patchI];
  92     }
  93
  94     if( Pstream::parRun() )
  95     {
  96         PtrList<processorBoundaryPatch>& procBoundaries =
  97             mesh_.procBoundaries_;
  98
  99         label nProcFaces(0);
 100         forAll(procBoundaries, patchI)
 101             nProcFaces += procBoundaries[patchI].patchSize();
 102
 103         //- copy processor faces into a separate list
 104         //- this preserves the order of faces after the modification is finished
 105         faceList procFaces(nProcFaces);
 106         nProcFaces = 0;
 107         forAll(procBoundaries, patchI)
 108         {
 109             const label start = procBoundaries[patchI].patchStart();
 110             const label end = start + procBoundaries[patchI].patchSize();
 111             for(label faceI=start;faceI<end;++faceI)
 112                 procFaces[nProcFaces++].transfer(faces[faceI]);
 113         }
 114
 115         //- create ordinary boundary faces from processor faces
 116         //- which are removed from one processor only
 117         //- these faces are stored in the latest ordinary patch
 118         forAll(procBoundaries, patchI)
 119         {
 120             //- send data to other proc
 121             const label start = procBoundaries[patchI].patchStart();
 122             boolList removeProcFace(procBoundaries[patchI].patchSize(), false);
 123             forAll(removeProcFace, faceI)
 124                 if( removeFace[start+faceI] )
 125                     removeProcFace[faceI] = true;
 126
 127             OPstream toOtherProc
 128             (
 129                 Pstream::blocking,
 130                 procBoundaries[patchI].neiProcNo(),
 131                 removeProcFace.byteSize()
 132             );
 133             toOtherProc << removeProcFace;
 134         }
 135
 136         --npI;
 137         boolList ordinaryBoundaryFace(procFaces.size(), false);
 138         nProcFaces = 0;
 139         forAll(procBoundaries, patchI)
 140         {
 141             const label start = procBoundaries[patchI].patchStart();
 142             boolList removeOtherProcFace;
 143             IPstream fromOtherProc
 144             (
 145                 Pstream::blocking,
 146                 procBoundaries[patchI].neiProcNo()
 147             );
 148             fromOtherProc >> removeOtherProcFace;
 149
 150             forAll(removeOtherProcFace, faceI)
 151             {
 152                 if( !removeFace[start+faceI] && removeOtherProcFace[faceI] )
 153                 {
 154                     //- this face becomes an ordinary boundary face
 155                     //- because the face on the other side of the processor
 156                     //- boundary has been deleted
 157                     ordinaryBoundaryFace[nProcFaces] = true;
 158                     ++boundaries[npI].patchSize();
 159                     if( nFaces < (start+faceI) )
 160                         faces[nFaces].transfer(procFaces[nProcFaces]);
 161                     newFaceLabel[start+faceI] = nFaces++;
 162                 }
 163
 164                 ++nProcFaces;
 165             }
 166         }
 167
 168         if( nProcFaces != procFaces.size() )
 169             FatalErrorIn
 170             (
 171                 "void polyMeshGenModifier::removeFaces()"
 172             ) << "Invalid number of processor faces!" << abort(FatalError);
 173
 174         //- remove faces from processor patches
 175         npI = 0;
 176         nFacesInPatch.setSize(procBoundaries.size());
 177         patchStart.setSize(procBoundaries.size());
 178         labelList neiProcNo(procBoundaries.size());
 179         nProcFaces = 0;
 180
 181         forAll(procBoundaries, patchI)
 182         {
 183             const label oldStart = procBoundaries[patchI].patchStart();
 184             const label oldNumFacesInPatch = procBoundaries[patchI].patchSize();
 185
 186             patchStart[npI] = nFaces;
 187             neiProcNo[npI] = procBoundaries[patchI].neiProcNo();
 188             nFacesInPatch[npI] = 0;
 189
 190             for(register label faceI=0;faceI<oldNumFacesInPatch;++faceI)
 191             {
 192                 if(
 193                     !removeFace[oldStart+faceI] &&
 194                     !ordinaryBoundaryFace[nProcFaces]
 195                 )
 196                 {
 197                     ++nFacesInPatch[npI];
 198                     faces[nFaces].transfer(procFaces[nProcFaces]);
 199                     newFaceLabel[oldStart+faceI] = nFaces++;
 200                 }
 201
 202                 ++nProcFaces;
 203             }
 204
 205             ++npI;
 206         }
 207
 208         //- all patches still exist
 209         forAll(procBoundaries, patchI)
 210         {
 211             procBoundaries[patchI].patchSize() = nFacesInPatch[patchI];
 212             procBoundaries[patchI].patchStart() = patchStart[patchI];
 213         }
 214     }
 215
 216     faces.setSize(nFaces);
 217
 218     //- update face subsets in the mesh
 219     mesh_.updateFaceSubsets(newFaceLabel);
 220
 221     //- change cells
 222     # ifdef USE_OMP
 223     # pragma omp parallel for if( cells.size() > 1000 ) schedule(dynamic, 40)
 224     # endif
 225     forAll(cells, cellI)
 226     {
 227         cell& c = cells[cellI];
 228
 229         DynList<label> newC;
 230
 231         forAll(c, fI)
 232             if( newFaceLabel[c[fI]] != -1 )
 233                 newC.append(newFaceLabel[c[fI]]);
 234
 235         c.setSize(newC.size());
 236
 237         forAll(c, fI)
 238             c[fI] = newC[fI];
 239     }
 240
 241     mesh_.clearOut();
 242
 243     Info << "Finished removing faces" << endl;
 244 }
 245
 246 void polyMeshGenModifier::removeDuplicateFaces()
 247 {
 248     VRWGraph& pointFaces = this->pointFaces();
 249
 250     faceListPMG& faces = mesh_.faces_;
 251
 252     labelLongList newFaceLabel(faces.size(), -1);
 253
 254     label nFaces(0);
 255     forAll(faces, faceI)
 256     {
 257         if( newFaceLabel[faceI] != -1 )
 258             continue;
 259
 260         const face& f = faces[faceI];
 261         DynList<label> duplicates;
 262         forAllRow(pointFaces, f[0], pfI)
 263         {
 264             const label faceJ = pointFaces(f[0], pfI);
 265
 266             if( faceI == faceJ )
 267                 continue;
 268
 269             if( faces[faceJ] == f )
 270                 duplicates.append(faceJ);
 271         }
 272
 273         newFaceLabel[faceI] = nFaces;
 274         forAll(duplicates, i)
 275             newFaceLabel[duplicates[i]] = nFaces;
 276         ++nFaces;
 277     }
 278
 279     label nInternalFaces(faces.size());
 280     if( mesh_.boundaries_.size() != 0 )
 281         nInternalFaces = mesh_.boundaries_[0].patchStart();
 282
 283     //- remove internal faces
 284     for(label faceI=0;faceI<nInternalFaces;++faceI)
 285     {
 286         if( newFaceLabel[faceI] >= faceI )
 287             continue;
 288
 289         faces[newFaceLabel[faceI]].transfer(faces[faceI]);
 290     }
 291
 292     //- update boundary faces (they cannot be duplicated)
 293     forAll(mesh_.boundaries_, patchI)
 294     {
 295         boundaryPatch& patch = mesh_.boundaries_[patchI];
 296
 297         const label start = patch.patchStart();
 298         const label end = start + patch.patchSize();
 299
 300         patch.patchStart() = newFaceLabel[start];
 301         for(label faceI=start;faceI<end;++faceI)
 302         {
 303             if( newFaceLabel[faceI] <  faceI )
 304                 faces[newFaceLabel[faceI]].transfer(faces[faceI]);
 305         }
 306     }
 307
 308     //- update processor faces (they cannot be duplicated)
 309     forAll(mesh_.procBoundaries_, patchI)
 310     {
 311         processorBoundaryPatch& patch = mesh_.procBoundaries_[patchI];
 312
 313         const label start = patch.patchStart();
 314         const label end = start + patch.patchSize();
 315
 316         patch.patchStart() = newFaceLabel[start];
 317         for(label faceI=start;faceI<end;++faceI)
 318         {
 319             if( newFaceLabel[faceI] < faceI )
 320                 faces[newFaceLabel[faceI]].transfer(faces[faceI]);
 321         }
 322     }
 323
 324     faces.setSize(nFaces);
 325     mesh_.updateFaceSubsets(newFaceLabel);
 326
 327     //- change cells
 328     cellListPMG& cells = mesh_.cells_;
 329     # ifdef USE_OMP
 330     # pragma omp parallel for if( cells.size() > 1000 ) schedule(dynamic, 40)
 331     # endif
 332     forAll(cells, cellI)
 333     {
 334         cell& c = cells[cellI];
 335
 336         DynList<label> newC;
 337
 338         forAll(c, fI)
 339             if( newFaceLabel[c[fI]] != -1 )
 340                 newC.append(newFaceLabel[c[fI]]);
 341
 342         c.setSize(newC.size());
 343
 344         forAll(c, fI)
 345             c[fI] = newC[fI];
 346     }
 347
 348     mesh_.clearOut();
 349 }
 350
 351 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 352
 353 } // End namespace Foam
 354
 355 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //