src/OpenFOAM/meshes/polyMesh/polyMeshFromShapeMesh.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 Description
  25     Create polyMesh from cell and patch shapes
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "polyMesh.H"
  30 #include "Time.H"
  31 #include "primitiveMesh.H"
  32 #include "DynamicList.H"
  33
  34 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  35
  36 Foam::labelListList Foam::polyMesh::cellShapePointCells
  37 (
  38     const cellShapeList& c
  39 ) const
  40 {
  41     List<DynamicList<label, primitiveMesh::cellsPerPoint_> >
  42         pc(points().size());
  43
  44     // For each cell
  45     forAll(c, i)
  46     {
  47         // For each vertex
  48         const labelList& labels = c[i];
  49
  50         forAll(labels, j)
  51         {
  52             // Set working point label
  53             label curPoint = labels[j];
  54             DynamicList<label, primitiveMesh::cellsPerPoint_>& curPointCells =
  55                 pc[curPoint];
  56
  57             // Enter the cell label in the point's cell list
  58             curPointCells.append(i);
  59         }
  60     }
  61
  62     labelListList pointCellAddr(pc.size());
  63
  64     forAll(pc, pointI)
  65     {
  66         pointCellAddr[pointI].transfer(pc[pointI]);
  67     }
  68
  69     return pointCellAddr;
  70 }
  71
  72
  73 Foam::labelList Foam::polyMesh::facePatchFaceCells
  74 (
  75     const faceList& patchFaces,
  76     const labelListList& pointCells,
  77     const faceListList& cellsFaceShapes,
  78     const label patchID
  79 ) const
  80 {
  81     register bool found;
  82
  83     labelList FaceCells(patchFaces.size());
  84
  85     forAll(patchFaces, fI)
  86     {
  87         found = false;
  88
  89         const face& curFace = patchFaces[fI];
  90         const labelList& facePoints = patchFaces[fI];
  91
  92         forAll(facePoints, pointI)
  93         {
  94             const labelList& facePointCells = pointCells[facePoints[pointI]];
  95
  96             forAll(facePointCells, cellI)
  97             {
  98                 faceList cellFaces = cellsFaceShapes[facePointCells[cellI]];
  99
 100                 forAll(cellFaces, cellFace)
 101                 {
 102                     if (cellFaces[cellFace] == curFace)
 103                     {
 104                         // Found the cell corresponding to this face
 105                         FaceCells[fI] = facePointCells[cellI];
 106
 107                         found = true;
 108                     }
 109                     if (found) break;
 110                 }
 111                 if (found) break;
 112             }
 113             if (found) break;
 114         }
 115
 116         if (!found)
 117         {
 118             FatalErrorIn
 119             (
 120                 "polyMesh::facePatchFaceCells(const faceList& patchFaces,"
 121                 "const labelListList& pointCells,"
 122                 "const faceListList& cellsFaceShapes,"
 123                 "const label patchID)"
 124             )   << "face " << fI << " in patch " << patchID
 125                 << " does not have neighbour cell"
 126                 << " face: " << patchFaces[fI]
 127                 << abort(FatalError);
 128         }
 129     }
 130
 131     return FaceCells;
 132 }
 133
 134
 135 //- Set faces_, calculate cells and patchStarts.
 136 void Foam::polyMesh::setTopology
 137 (
 138     const cellShapeList& cellsAsShapes,
 139     const faceListList& boundaryFaces,
 140     const wordList& boundaryPatchNames,
 141     labelList& patchSizes,
 142     labelList& patchStarts,
 143     label& defaultPatchStart,
 144     label& nFaces,
 145     cellList& cells
 146 )
 147 {
 148     // Calculate the faces of all cells
 149     // Initialise maximum possible numer of mesh faces to 0
 150     label maxFaces = 0;
 151
 152     // Set up a list of face shapes for each cell
 153     faceListList cellsFaceShapes(cellsAsShapes.size());
 154     cells.setSize(cellsAsShapes.size());
 155
 156     forAll(cellsFaceShapes, cellI)
 157     {
 158         cellsFaceShapes[cellI] = cellsAsShapes[cellI].faces();
 159
 160         cells[cellI].setSize(cellsFaceShapes[cellI].size());
 161
 162         // Initialise cells to -1 to flag undefined faces
 163         static_cast<labelList&>(cells[cellI]) = -1;
 164
 165         // Count maximum possible numer of mesh faces
 166         maxFaces += cellsFaceShapes[cellI].size();
 167     }
 168
 169     // Set size of faces array to maximum possible number of mesh faces
 170     faces_.setSize(maxFaces);
 171
 172     // Initialise number of faces to 0
 173     nFaces = 0;
 174
 175     // set reference to point-cell addressing
 176     labelListList PointCells = cellShapePointCells(cellsAsShapes);
 177
 178     bool found = false;
 179
 180     forAll(cells, cellI)
 181     {
 182         // Note:
 183         // Insertion cannot be done in one go as the faces need to be
 184         // added into the list in the increasing order of neighbour
 185         // cells.  Therefore, all neighbours will be detected first
 186         // and then added in the correct order.
 187
 188         const faceList& curFaces = cellsFaceShapes[cellI];
 189
 190         // Record the neighbour cell
 191         labelList neiCells(curFaces.size(), -1);
 192
 193         // Record the face of neighbour cell
 194         labelList faceOfNeiCell(curFaces.size(), -1);
 195
 196         label nNeighbours = 0;
 197
 198         // For all faces ...
 199         forAll(curFaces, faceI)
 200         {
 201             // Skip faces that have already been matched
 202             if (cells[cellI][faceI] >= 0) continue;
 203
 204             found = false;
 205
 206             const face& curFace = curFaces[faceI];
 207
 208             // Get the list of labels
 209             const labelList& curPoints = curFace;
 210
 211             // For all points
 212             forAll(curPoints, pointI)
 213             {
 214                 // dGget the list of cells sharing this point
 215                 const labelList& curNeighbours =
 216                     PointCells[curPoints[pointI]];
 217
 218                 // For all neighbours
 219                 forAll(curNeighbours, neiI)
 220                 {
 221                     label curNei = curNeighbours[neiI];
 222
 223                     // Reject neighbours with the lower label
 224                     if (curNei > cellI)
 225                     {
 226                         // Get the list of search faces
 227                         const faceList& searchFaces = cellsFaceShapes[curNei];
 228
 229                         forAll(searchFaces, neiFaceI)
 230                         {
 231                             if (searchFaces[neiFaceI] == curFace)
 232                             {
 233                                 // Match!!
 234                                 found = true;
 235
 236                                 // Record the neighbour cell and face
 237                                 neiCells[faceI] = curNei;
 238                                 faceOfNeiCell[faceI] = neiFaceI;
 239                                 nNeighbours++;
 240
 241                                 break;
 242                             }
 243                         }
 244                         if (found) break;
 245                     }
 246                     if (found) break;
 247                 }
 248                 if (found) break;
 249             } // End of current points
 250         }  // End of current faces
 251
 252         // Add the faces in the increasing order of neighbours
 253         for (label neiSearch = 0; neiSearch < nNeighbours; neiSearch++)
 254         {
 255             // Find the lowest neighbour which is still valid
 256             label nextNei = -1;
 257             label minNei = cells.size();
 258
 259             forAll(neiCells, ncI)
 260             {
 261                 if (neiCells[ncI] > -1 && neiCells[ncI] < minNei)
 262                 {
 263                     nextNei = ncI;
 264                     minNei = neiCells[ncI];
 265                 }
 266             }
 267
 268             if (nextNei > -1)
 269             {
 270                 // Add the face to the list of faces
 271                 faces_[nFaces] = curFaces[nextNei];
 272
 273                 // Set cell-face and cell-neighbour-face to current face label
 274                 cells[cellI][nextNei] = nFaces;
 275                 cells[neiCells[nextNei]][faceOfNeiCell[nextNei]] = nFaces;
 276
 277                 // Stop the neighbour from being used again
 278                 neiCells[nextNei] = -1;
 279
 280                 // Increment number of faces counter
 281                 nFaces++;
 282             }
 283             else
 284             {
 285                 FatalErrorIn
 286                 (
 287                     "polyMesh::setTopology\n"
 288                     "(\n"
 289                     "    const cellShapeList& cellsAsShapes,\n"
 290                     "    const faceListList& boundaryFaces,\n"
 291                     "    const wordList& boundaryPatchNames,\n"
 292                     "    labelList& patchSizes,\n"
 293                     "    labelList& patchStarts,\n"
 294                     "    label& defaultPatchStart,\n"
 295                     "    label& nFaces,\n"
 296                     "    cellList& cells\n"
 297                     ")"
 298                 )   << "Error in internal face insertion"
 299                     << abort(FatalError);
 300             }
 301         }
 302     }
 303
 304     // Do boundary faces
 305
 306     patchSizes.setSize(boundaryFaces.size(), -1);
 307     patchStarts.setSize(boundaryFaces.size(), -1);
 308
 309     forAll(boundaryFaces, patchI)
 310     {
 311         const faceList& patchFaces = boundaryFaces[patchI];
 312
 313         labelList curPatchFaceCells =
 314             facePatchFaceCells
 315             (
 316                 patchFaces,
 317                 PointCells,
 318                 cellsFaceShapes,
 319                 patchI
 320             );
 321
 322         // Grab the start label
 323         label curPatchStart = nFaces;
 324
 325         forAll(patchFaces, faceI)
 326         {
 327             const face& curFace = patchFaces[faceI];
 328
 329             const label cellInside = curPatchFaceCells[faceI];
 330
 331             faces_[nFaces] = curFace;
 332
 333             // get faces of the cell inside
 334             const faceList& facesOfCellInside = cellsFaceShapes[cellInside];
 335
 336             bool found = false;
 337
 338             forAll(facesOfCellInside, cellFaceI)
 339             {
 340                 if (facesOfCellInside[cellFaceI] == curFace)
 341                 {
 342                     if (cells[cellInside][cellFaceI] >= 0)
 343                     {
 344                         FatalErrorIn
 345                         (
 346                             "polyMesh::setTopology\n"
 347                             "(\n"
 348                             "    const cellShapeList& cellsAsShapes,\n"
 349                             "    const faceListList& boundaryFaces,\n"
 350                             "    const wordList& boundaryPatchNames,\n"
 351                             "    labelList& patchSizes,\n"
 352                             "    labelList& patchStarts,\n"
 353                             "    label& defaultPatchStart,\n"
 354                             "    label& nFaces,\n"
 355                             "    cellList& cells\n"
 356                             ")"
 357                         )   << "Trying to specify a boundary face " << curFace
 358                             << " on the face on cell " << cellInside
 359                             << " which is either an internal face or already "
 360                             << "belongs to some other patch.  This is face "
 361                             << faceI << " of patch "
 362                             << patchI << " named "
 363                             << boundaryPatchNames[patchI] << "."
 364                             << abort(FatalError);
 365                     }
 366
 367                     found = true;
 368
 369                     cells[cellInside][cellFaceI] = nFaces;
 370
 371                     break;
 372                 }
 373             }
 374
 375             if (!found)
 376             {
 377                 FatalErrorIn("polyMesh::polyMesh(... construct from shapes...)")
 378                     << "face " << faceI << " of patch " << patchI
 379                     << " does not seem to belong to cell " << cellInside
 380                     << " which, according to the addressing, "
 381                     << "should be next to it."
 382                     << abort(FatalError);
 383             }
 384
 385             // increment the counter of faces
 386             nFaces++;
 387         }
 388
 389         patchSizes[patchI] = nFaces - curPatchStart;
 390         patchStarts[patchI] = curPatchStart;
 391     }
 392
 393     // Grab "non-existing" faces and put them into a default patch
 394
 395     defaultPatchStart = nFaces;
 396
 397     forAll(cells, cellI)
 398     {
 399         labelList& curCellFaces = cells[cellI];
 400
 401         forAll(curCellFaces, faceI)
 402         {
 403             if (curCellFaces[faceI] == -1) // "non-existent" face
 404             {
 405                 curCellFaces[faceI] = nFaces;
 406                 faces_[nFaces] = cellsFaceShapes[cellI][faceI];
 407
 408                 nFaces++;
 409             }
 410         }
 411     }
 412
 413     // Reset the size of the face list
 414     faces_.setSize(nFaces);
 415
 416     return ;
 417 }
 418
 419
 420 Foam::polyMesh::polyMesh
 421 (
 422     const IOobject& io,
 423     const Xfer<pointField>& points,
 424     const cellShapeList& cellsAsShapes,
 425     const faceListList& boundaryFaces,
 426     const wordList& boundaryPatchNames,
 427     const wordList& boundaryPatchTypes,
 428     const word& defaultBoundaryPatchName,
 429     const word& defaultBoundaryPatchType,
 430     const wordList& boundaryPatchPhysicalTypes,
 431     const bool syncPar
 432 )
 433 :
 434     objectRegistry(io),
 435     primitiveMesh(),
 436     points_
 437     (
 438         IOobject
 439         (
 440             "points",
 441             instance(),
 442             meshSubDir,
 443             *this,
 444             IOobject::NO_READ,
 445             IOobject::AUTO_WRITE
 446         ),
 447         points
 448     ),
 449     faces_
 450     (
 451         IOobject
 452         (
 453             "faces",
 454             instance(),
 455             meshSubDir,
 456             *this,
 457             IOobject::NO_READ,
 458             IOobject::AUTO_WRITE
 459         ),
 460         0
 461     ),
 462     owner_
 463     (
 464         IOobject
 465         (
 466             "owner",
 467             instance(),
 468             meshSubDir,
 469             *this,
 470             IOobject::NO_READ,
 471             IOobject::AUTO_WRITE
 472         ),
 473         0
 474     ),
 475     neighbour_
 476     (
 477         IOobject
 478         (
 479             "neighbour",
 480             instance(),
 481             meshSubDir,
 482             *this,
 483             IOobject::NO_READ,
 484             IOobject::AUTO_WRITE
 485         ),
 486         0
 487     ),
 488     clearedPrimitives_(false),
 489     boundary_
 490     (
 491         IOobject
 492         (
 493             "boundary",
 494             instance(),
 495             meshSubDir,
 496             *this,
 497             IOobject::NO_READ,
 498             IOobject::AUTO_WRITE
 499         ),
 500         *this,
 501         boundaryFaces.size() + 1    // add room for a default patch
 502     ),
 503     bounds_(points_, syncPar),
 504     geometricD_(Vector<label>::zero),
 505     solutionD_(Vector<label>::zero),
 506     tetBasePtIsPtr_(NULL),
 507     pointZones_
 508     (
 509         IOobject
 510         (
 511             "pointZones",
 512             instance(),
 513             meshSubDir,
 514             *this,
 515             IOobject::NO_READ,
 516             IOobject::NO_WRITE
 517         ),
 518         *this,
 519         0
 520     ),
 521     faceZones_
 522     (
 523         IOobject
 524         (
 525             "faceZones",
 526             instance(),
 527             meshSubDir,
 528             *this,
 529             IOobject::NO_READ,
 530             IOobject::NO_WRITE
 531         ),
 532         *this,
 533         0
 534     ),
 535     cellZones_
 536     (
 537         IOobject
 538         (
 539             "cellZones",
 540             instance(),
 541             meshSubDir,
 542             *this,
 543             IOobject::NO_READ,
 544             IOobject::NO_WRITE
 545         ),
 546         *this,
 547         0
 548     ),
 549     globalMeshDataPtr_(NULL),
 550     moving_(false),
 551     curMotionTimeIndex_(time().timeIndex()),
 552     oldPointsPtr_(NULL)
 553 {
 554     if (debug)
 555     {
 556         Info<<"Constructing polyMesh from cell and boundary shapes." << endl;
 557     }
 558
 559     // Remove all of the old mesh files if they exist
 560     removeFiles(instance());
 561
 562     // Calculate faces and cells
 563     labelList patchSizes;
 564     labelList patchStarts;
 565     label defaultPatchStart;
 566     label nFaces;
 567     cellList cells;
 568     setTopology
 569     (
 570         cellsAsShapes,
 571         boundaryFaces,
 572         boundaryPatchNames,
 573         patchSizes,
 574         patchStarts,
 575         defaultPatchStart,
 576         nFaces,
 577         cells
 578     );
 579
 580     // Warning: Patches can only be added once the face list is
 581     // completed, as they hold a subList of the face list
 582     forAll(boundaryFaces, patchI)
 583     {
 584         // add the patch to the list
 585         boundary_.set
 586         (
 587             patchI,
 588             polyPatch::New
 589             (
 590                 boundaryPatchTypes[patchI],
 591                 boundaryPatchNames[patchI],
 592                 patchSizes[patchI],
 593                 patchStarts[patchI],
 594                 patchI,
 595                 boundary_
 596             )
 597         );
 598
 599         if
 600         (
 601             boundaryPatchPhysicalTypes.size()
 602          && boundaryPatchPhysicalTypes[patchI].size()
 603         )
 604         {
 605             boundary_[patchI].physicalType() =
 606                 boundaryPatchPhysicalTypes[patchI];
 607         }
 608     }
 609
 610     label nAllPatches = boundaryFaces.size();
 611
 612     if (nFaces > defaultPatchStart)
 613     {
 614         WarningIn("polyMesh::polyMesh(... construct from shapes...)")
 615             << "Found " << nFaces - defaultPatchStart
 616             << " undefined faces in mesh; adding to default patch." << endl;
 617
 618         // Check if there already exists a defaultFaces patch as last patch
 619         // and reuse it.
 620         label patchI = findIndex(boundaryPatchNames, defaultBoundaryPatchName);
 621
 622         if (patchI != -1)
 623         {
 624             if (patchI != boundaryFaces.size()-1 || boundary_[patchI].size())
 625             {
 626                 FatalErrorIn("polyMesh::polyMesh(... construct from shapes...)")
 627                     << "Default patch " << boundary_[patchI].name()
 628                     << " already has faces in it or is not"
 629                     << " last in list of patches." << exit(FatalError);
 630             }
 631
 632             WarningIn("polyMesh::polyMesh(... construct from shapes...)")
 633                 << "Reusing existing patch " << patchI
 634                 << " for undefined faces." << endl;
 635
 636             boundary_.set
 637             (
 638                 patchI,
 639                 polyPatch::New
 640                 (
 641                     boundaryPatchTypes[patchI],
 642                     boundaryPatchNames[patchI],
 643                     nFaces - defaultPatchStart,
 644                     defaultPatchStart,
 645                     patchI,
 646                     boundary_
 647                 )
 648             );
 649         }
 650         else
 651         {
 652             boundary_.set
 653             (
 654                 nAllPatches,
 655                 polyPatch::New
 656                 (
 657                     defaultBoundaryPatchType,
 658                     defaultBoundaryPatchName,
 659                     nFaces - defaultPatchStart,
 660                     defaultPatchStart,
 661                     boundary_.size() - 1,
 662                     boundary_
 663                 )
 664             );
 665
 666             nAllPatches++;
 667         }
 668     }
 669
 670     // Reset the size of the boundary
 671     boundary_.setSize(nAllPatches);
 672
 673     // Set the primitive mesh
 674     initMesh(cells);
 675
 676     if (syncPar)
 677     {
 678         // Calculate topology for the patches (processor-processor comms etc.)
 679         boundary_.updateMesh();
 680
 681         // Calculate the geometry for the patches (transformation tensors etc.)
 682         boundary_.calcGeometry();
 683     }
 684
 685     if (debug)
 686     {
 687         if (checkMesh())
 688         {
 689             Info<< "Mesh OK" << endl;
 690         }
 691     }
 692 }
 693
 694
 695 Foam::polyMesh::polyMesh
 696 (
 697     const IOobject& io,
 698     const Xfer<pointField>& points,
 699     const cellShapeList& cellsAsShapes,
 700     const faceListList& boundaryFaces,
 701     const wordList& boundaryPatchNames,
 702     const PtrList<dictionary>& boundaryDicts,
 703     const word& defaultBoundaryPatchName,
 704     const word& defaultBoundaryPatchType,
 705     const bool syncPar
 706 )
 707 :
 708     objectRegistry(io),
 709     primitiveMesh(),
 710     points_
 711     (
 712         IOobject
 713         (
 714             "points",
 715             instance(),
 716             meshSubDir,
 717             *this,
 718             IOobject::NO_READ,
 719             IOobject::AUTO_WRITE
 720         ),
 721         points
 722     ),
 723     faces_
 724     (
 725         IOobject
 726         (
 727             "faces",
 728             instance(),
 729             meshSubDir,
 730             *this,
 731             IOobject::NO_READ,
 732             IOobject::AUTO_WRITE
 733         ),
 734         0
 735     ),
 736     owner_
 737     (
 738         IOobject
 739         (
 740             "owner",
 741             instance(),
 742             meshSubDir,
 743             *this,
 744             IOobject::NO_READ,
 745             IOobject::AUTO_WRITE
 746         ),
 747         0
 748     ),
 749     neighbour_
 750     (
 751         IOobject
 752         (
 753             "neighbour",
 754             instance(),
 755             meshSubDir,
 756             *this,
 757             IOobject::NO_READ,
 758             IOobject::AUTO_WRITE
 759         ),
 760         0
 761     ),
 762     clearedPrimitives_(false),
 763     boundary_
 764     (
 765         IOobject
 766         (
 767             "boundary",
 768             instance(),
 769             meshSubDir,
 770             *this,
 771             IOobject::NO_READ,
 772             IOobject::AUTO_WRITE
 773         ),
 774         *this,
 775         boundaryFaces.size() + 1    // add room for a default patch
 776     ),
 777     bounds_(points_, syncPar),
 778     geometricD_(Vector<label>::zero),
 779     solutionD_(Vector<label>::zero),
 780     tetBasePtIsPtr_(NULL),
 781     pointZones_
 782     (
 783         IOobject
 784         (
 785             "pointZones",
 786             instance(),
 787             meshSubDir,
 788             *this,
 789             IOobject::NO_READ,
 790             IOobject::NO_WRITE
 791         ),
 792         *this,
 793         0
 794     ),
 795     faceZones_
 796     (
 797         IOobject
 798         (
 799             "faceZones",
 800             instance(),
 801             meshSubDir,
 802             *this,
 803             IOobject::NO_READ,
 804             IOobject::NO_WRITE
 805         ),
 806         *this,
 807         0
 808     ),
 809     cellZones_
 810     (
 811         IOobject
 812         (
 813             "cellZones",
 814             instance(),
 815             meshSubDir,
 816             *this,
 817             IOobject::NO_READ,
 818             IOobject::NO_WRITE
 819         ),
 820         *this,
 821         0
 822     ),
 823     globalMeshDataPtr_(NULL),
 824     moving_(false),
 825     curMotionTimeIndex_(time().timeIndex()),
 826     oldPointsPtr_(NULL)
 827 {
 828     if (debug)
 829     {
 830         Info<<"Constructing polyMesh from cell and boundary shapes." << endl;
 831     }
 832
 833     // Remove all of the old mesh files if they exist
 834     removeFiles(instance());
 835
 836     // Calculate faces and cells
 837     labelList patchSizes;
 838     labelList patchStarts;
 839     label defaultPatchStart;
 840     label nFaces;
 841     cellList cells;
 842     setTopology
 843     (
 844         cellsAsShapes,
 845         boundaryFaces,
 846         boundaryPatchNames,
 847         patchSizes,
 848         patchStarts,
 849         defaultPatchStart,
 850         nFaces,
 851         cells
 852     );
 853
 854     // Warning: Patches can only be added once the face list is
 855     // completed, as they hold a subList of the face list
 856     forAll(boundaryDicts, patchI)
 857     {
 858         dictionary patchDict(boundaryDicts[patchI]);
 859
 860         patchDict.set("nFaces", patchSizes[patchI]);
 861         patchDict.set("startFace", patchStarts[patchI]);
 862
 863         // add the patch to the list
 864         boundary_.set
 865         (
 866             patchI,
 867             polyPatch::New
 868             (
 869                 boundaryPatchNames[patchI],
 870                 patchDict,
 871                 patchI,
 872                 boundary_
 873             )
 874         );
 875     }
 876
 877     label nAllPatches = boundaryFaces.size();
 878
 879     if (nFaces > defaultPatchStart)
 880     {
 881         WarningIn("polyMesh::polyMesh(... construct from shapes...)")
 882             << "Found " << nFaces - defaultPatchStart
 883             << " undefined faces in mesh; adding to default patch." << endl;
 884
 885         // Check if there already exists a defaultFaces patch as last patch
 886         // and reuse it.
 887         label patchI = findIndex(boundaryPatchNames, defaultBoundaryPatchName);
 888
 889         if (patchI != -1)
 890         {
 891             if (patchI != boundaryFaces.size()-1 || boundary_[patchI].size())
 892             {
 893                 FatalErrorIn("polyMesh::polyMesh(... construct from shapes...)")
 894                     << "Default patch " << boundary_[patchI].name()
 895                     << " already has faces in it or is not"
 896                     << " last in list of patches." << exit(FatalError);
 897             }
 898
 899             WarningIn("polyMesh::polyMesh(... construct from shapes...)")
 900                 << "Reusing existing patch " << patchI
 901                 << " for undefined faces." << endl;
 902
 903             boundary_.set
 904             (
 905                 patchI,
 906                 polyPatch::New
 907                 (
 908                     boundary_[patchI].type(),
 909                     boundary_[patchI].name(),
 910                     nFaces - defaultPatchStart,
 911                     defaultPatchStart,
 912                     patchI,
 913                     boundary_
 914                 )
 915             );
 916         }
 917         else
 918         {
 919             boundary_.set
 920             (
 921                 nAllPatches,
 922                 polyPatch::New
 923                 (
 924                     defaultBoundaryPatchType,
 925                     defaultBoundaryPatchName,
 926                     nFaces - defaultPatchStart,
 927                     defaultPatchStart,
 928                     boundary_.size() - 1,
 929                     boundary_
 930                 )
 931             );
 932
 933             nAllPatches++;
 934         }
 935     }
 936
 937     // Reset the size of the boundary
 938     boundary_.setSize(nAllPatches);
 939
 940     // Set the primitive mesh
 941     initMesh(cells);
 942
 943     if (syncPar)
 944     {
 945         // Calculate topology for the patches (processor-processor comms etc.)
 946         boundary_.updateMesh();
 947
 948         // Calculate the geometry for the patches (transformation tensors etc.)
 949         boundary_.calcGeometry();
 950     }
 951
 952     if (debug)
 953     {
 954         if (checkMesh())
 955         {
 956             Info << "Mesh OK" << endl;
 957         }
 958     }
 959 }
 960
 961
 962 // ************************************************************************* //