src/mesh/cfMesh/utilities/meshes/polyMeshGenModifier/polyMeshGenModifierReorderBoundaryFaces.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 // #define DEBUGSearch
  36
  37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  38
  39 namespace Foam
  40 {
  41
  42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  43
  44 void polyMeshGenModifier::reorderBoundaryFaces()
  45 {
  46     Info << "Reordering boundary faces " << endl;
  47
  48     if( Pstream::parRun() )
  49         reorderProcBoundaryFaces();
  50
  51     faceListPMG& faces = mesh_.faces_;
  52     cellListPMG& cells = mesh_.cells_;
  53
  54     const labelList& neighbour = mesh_.neighbour();
  55     const label nInternalFaces = mesh_.nInternalFaces();
  56
  57     //- count internal and boundary faces
  58     const label numBFaces = faces.size() - nInternalFaces;
  59
  60     labelLongList newFaceLabel(faces.size(), -1);
  61
  62     //- find faces which should be repositioned
  63     label nReplaced(0);
  64     labelList internalToChange;
  65     labelList boundaryToChange;
  66
  67     # ifdef USE_OMP
  68     const label nThreads = 3 * omp_get_num_procs();
  69     # else
  70     const label nThreads(1);
  71     # endif
  72     labelList nInternalToChangeThread(nThreads);
  73     labelList nBoundaryToChangeThread(nThreads);
  74
  75     # ifdef USE_OMP
  76     # pragma omp parallel num_threads(nThreads)
  77     # endif
  78     {
  79         # ifdef USE_OMP
  80         const label threadI = omp_get_thread_num();
  81         # else
  82         const label threadI = 0;
  83         # endif
  84
  85         label& nItc = nInternalToChangeThread[threadI];
  86         label& nBtc = nBoundaryToChangeThread[threadI];
  87
  88         labelLongList internalToChangeLocal, boundaryToChangeLocal;
  89
  90         //- find the boundary faces within the range of internal faces
  91         # ifdef USE_OMP
  92         # pragma omp for schedule(static)
  93         # endif
  94         for(label faceI=0;faceI<nInternalFaces;++faceI)
  95         {
  96             if( neighbour[faceI] == -1 )
  97                 internalToChangeLocal.append(faceI);
  98         }
  99
 100         nItc = internalToChangeLocal.size();
 101
 102         //- find the internal faces within the range of boundary faces
 103         # ifdef USE_OMP
 104         # pragma omp for schedule(static)
 105         # endif
 106         for(label faceI=nInternalFaces;faceI<faces.size();++faceI)
 107         {
 108             if( neighbour[faceI] != -1 )
 109                 boundaryToChangeLocal.append(faceI);
 110         }
 111
 112         nBtc = boundaryToChangeLocal.size();
 113
 114         //- perform reduction such that all threads know how many faces
 115         //- need to be swapped
 116         # ifdef USE_OMP
 117         # pragma omp critical
 118         # endif
 119         nReplaced += nBtc;
 120
 121         # ifdef USE_OMP
 122         # pragma omp barrier
 123
 124         # pragma omp master
 125         # endif
 126         {
 127             internalToChange.setSize(nReplaced);
 128             boundaryToChange.setSize(nReplaced);
 129         }
 130
 131         # ifdef USE_OMP
 132         # pragma omp barrier
 133         # endif
 134
 135         label localStart(0);
 136         for(label i=0;i<threadI;++i)
 137             localStart += nInternalToChangeThread[i];
 138
 139         forAll(internalToChangeLocal, i)
 140             internalToChange[localStart++] = internalToChangeLocal[i];
 141
 142         localStart = 0;
 143         for(label i=0;i<threadI;++i)
 144             localStart += nBoundaryToChangeThread[i];
 145
 146         forAll(boundaryToChangeLocal, i)
 147             boundaryToChange[localStart++] = boundaryToChangeLocal[i];
 148
 149         # ifdef USE_OMP
 150         # pragma omp barrier
 151
 152         //- start moving positions of faces
 153         # pragma omp for schedule(static)
 154         # endif
 155         forAll(internalToChange, fI)
 156         {
 157             //- swap with the face at the location the face should be
 158             face f;
 159             f.transfer(faces[internalToChange[fI]]);
 160             faces[internalToChange[fI]].transfer(faces[boundaryToChange[fI]]);
 161             faces[boundaryToChange[fI]].transfer(f);
 162             newFaceLabel[internalToChange[fI]] = boundaryToChange[fI];
 163             newFaceLabel[boundaryToChange[fI]] = internalToChange[fI];
 164         }
 165
 166         # ifdef USE_OMP
 167         # pragma omp barrier
 168
 169         //- renumber cells
 170         # pragma omp for schedule(dynamic, 40)
 171         # endif
 172         forAll(cells, cellI)
 173         {
 174             cell& c = cells[cellI];
 175
 176             forAll(c, fI)
 177                 if( newFaceLabel[c[fI]] != -1 )
 178                     c[fI] = newFaceLabel[c[fI]];
 179         }
 180     }
 181
 182     //- re-create boundary data
 183     PtrList<boundaryPatch>& boundaries = mesh_.boundaries_;
 184     if( boundaries.size() != 1 )
 185     {
 186         boundaries.clear();
 187         boundaries.setSize(1);
 188         boundaries.set
 189         (
 190             0,
 191             new boundaryPatch
 192             (
 193                 "defaultFaces",
 194                 "patch",
 195                 numBFaces,
 196                 nInternalFaces
 197             )
 198         );
 199     }
 200     else
 201     {
 202         boundaries[0].patchStart() = nInternalFaces;
 203         boundaries[0].patchSize() = numBFaces;
 204     }
 205
 206     if( Pstream::parRun() )
 207     {
 208         //- processor boundary faces must be contained at the end
 209         label nProcFaces(0);
 210         forAll(mesh_.procBoundaries_, procPatchI)
 211             nProcFaces += mesh_.procBoundaries_[procPatchI].patchSize();
 212
 213         boundaries[0].patchSize() -= nProcFaces;
 214     }
 215
 216     //- update face subsets
 217     mesh_.updateFaceSubsets(newFaceLabel);
 218
 219     //- delete invalid data
 220     mesh_.clearOut();
 221     this->clearOut();
 222
 223     Info << "Finished reordering boundary faces" << endl;
 224 }
 225
 226 void polyMeshGenModifier::reorderProcBoundaryFaces()
 227 {
 228     PtrList<processorBoundaryPatch>& procBoundaries = mesh_.procBoundaries_;
 229     if( procBoundaries.size() == 0 )
 230     {
 231         Warning << "Processor " << Pstream::myProcNo() << " has no "
 232             << "processor boundaries!" << endl;
 233         return;
 234     }
 235
 236     //- check if there exist any internal or ordinary bnd faces
 237     //- which appear after processor bnd faces. Move those faces before
 238     //- the processor boundary
 239     const label origProcStart = procBoundaries[0].patchStart();
 240     label nProcFaces(0);
 241     forAll(procBoundaries, patchI)
 242         nProcFaces += procBoundaries[patchI].patchSize();
 243
 244     faceListPMG& faces = mesh_.faces_;
 245     cellListPMG& cells = mesh_.cells_;
 246
 247     const label shift = faces.size() - (origProcStart + nProcFaces);
 248     if( shift == 0 )
 249         return;
 250     if( shift < 0 )
 251         FatalErrorIn
 252         (
 253             "void polyMeshGenModifier::reorderProcBoundaryFaces()"
 254         ) << "Missing some faces!" << abort(FatalError);
 255
 256     labelLongList newFaceLabel(faces.size(), -1);
 257
 258     //- faces added after processor boundaries should be moved up front
 259     faceList facesAtEnd(shift);
 260     label counter(0);
 261     for(label faceI=(origProcStart + nProcFaces);faceI<faces.size();++faceI)
 262     {
 263         facesAtEnd[counter].transfer(faces[faceI]);
 264         newFaceLabel[faceI] = origProcStart + counter;
 265         ++counter;
 266     }
 267
 268     //- shift proc faces
 269     forAllReverse(procBoundaries, patchI)
 270     {
 271         const label start = procBoundaries[patchI].patchStart();
 272         const label end = start + procBoundaries[patchI].patchSize();
 273
 274         //- set patch start to the new value
 275         procBoundaries[patchI].patchStart() += shift;
 276
 277         for(label faceI=end-1;faceI>=start;--faceI)
 278         {
 279             faces[faceI+shift].transfer(faces[faceI]);
 280             newFaceLabel[faceI] = faceI + shift;
 281         }
 282     }
 283
 284     //- store faces taken from the end
 285     forAll(facesAtEnd, fI)
 286     {
 287         faces[origProcStart+fI].transfer(facesAtEnd[fI]);
 288     }
 289
 290     //- set correct patch size
 291     PtrList<boundaryPatch>& boundaries = mesh_.boundaries_;
 292     if( boundaries.size() == 1 )
 293     {
 294         boundaries[0].patchSize() =
 295             procBoundaries[0].patchStart() - boundaries[0].patchStart();
 296     }
 297     else
 298     {
 299         const label start = boundaries[0].patchStart();
 300
 301         boundaries.clear();
 302         boundaries.setSize(1);
 303         boundaries.set
 304         (
 305             0,
 306             new boundaryPatch
 307             (
 308                 "defaultFaces",
 309                 "patch",
 310                 procBoundaries[0].patchStart() - start,
 311                 start
 312             )
 313         );
 314     }
 315
 316     //- renumber cells
 317     # ifdef USE_OMP
 318     # pragma omp parallel for schedule(dynamic, 40)
 319     # endif
 320     forAll(cells, cellI)
 321     {
 322         cell& c = cells[cellI];
 323
 324         forAll(c, fI)
 325             if( newFaceLabel[c[fI]] != -1 )
 326                 c[fI] = newFaceLabel[c[fI]];
 327     }
 328
 329     //- update face subsets
 330     mesh_.updateFaceSubsets(newFaceLabel);
 331
 332     //- delete invalid data
 333     mesh_.clearOut();
 334     this->clearOut();
 335 }
 336
 337 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 338
 339 } // End namespace Foam
 340
 341 // ************************************************************************* //