src/conversion/meshReader/createPolyBoundary.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     boundary faces
  26     - use pointCells when searching for connectivity
  27     - initialize the cell connectivity with '-1'
  28     - find both cell faces corresponding to the baffles and mark them
  29       to prevent a connection
  30     - standard connectivity checks
  31
  32     - added baffle and monitoring support
  33
  34 \*---------------------------------------------------------------------------*/
  35
  36 #include "meshReader.H"
  37 #include "Time.H"
  38 #include "polyPatch.H"
  39 #include "emptyPolyPatch.H"
  40 #include "preservePatchTypes.H"
  41
  42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  43
  44 void Foam::meshReader::addPolyBoundaryFace
  45 (
  46     const label cellId,
  47     const label cellFaceId,
  48     const label nCreatedFaces
  49 )
  50 {
  51 #ifdef DEBUG_BOUNDARY
  52     Info<< nCreatedFaces
  53         << " add bnd for cell " << cellId
  54         << " face " << cellFaceId
  55         << " (original cell " << origCellId_[cellId] << ")"
  56         << endl;
  57 #endif
  58
  59     // standard case: volume cells
  60     const face& thisFace = cellFaces_[cellId][cellFaceId];
  61
  62     // Debugging
  63     if (cellPolys_[cellId][cellFaceId] > nInternalFaces_)
  64     {
  65         Info<< "meshReader::createPolyBoundary(): "
  66             << "Problem with face: " << thisFace << endl
  67             << "Probably multiple definitions "
  68             << "of a single boundary face." << endl
  69             << endl;
  70     }
  71     else if (cellPolys_[cellId][cellFaceId] >= 0)
  72     {
  73         Info<< "meshReader::createPolyBoundary(): "
  74             << "Problem with face: " << thisFace << endl
  75             << "Probably trying to define a boundary face "
  76             << "on a previously matched internal face." << endl
  77             << "Internal face: "
  78             << meshFaces_[cellPolys_[cellId][cellFaceId]]
  79             << endl;
  80     }
  81
  82     meshFaces_[nCreatedFaces] = thisFace;
  83     cellPolys_[cellId][cellFaceId] = nCreatedFaces;
  84 }
  85
  86
  87 void Foam::meshReader::addPolyBoundaryFace
  88 (
  89     const cellFaceIdentifier& identifier,
  90     const label nCreatedFaces
  91 )
  92 {
  93     addPolyBoundaryFace(identifier.cell, identifier.face, nCreatedFaces);
  94 }
  95
  96
  97 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  98
  99 void Foam::meshReader::createPolyBoundary()
 100 {
 101     label nBoundaryFaces = 0;
 102     label nMissingFaces = 0;
 103     label nInterfaces = 0;
 104
 105     const faceListList& cFaces = cellFaces();
 106
 107     // determine number of non-patched faces:
 108     forAll(cellPolys_, cellI)
 109     {
 110         cell& curCell = cellPolys_[cellI];
 111
 112         forAll(curCell, fI)
 113         {
 114             if (curCell[fI] < 0)
 115             {
 116                 nMissingFaces++;
 117             }
 118         }
 119     }
 120
 121     forAll(boundaryIds_, patchI)
 122     {
 123         nBoundaryFaces += boundaryIds_[patchI].size();
 124     }
 125
 126     Info<< nl
 127         << "There are " << nMissingFaces
 128         << " faces to be patched and " << nBoundaryFaces
 129         << " specified - collect missed boundaries to final patch" << endl;
 130
 131     patchStarts_.setSize(boundaryIds_.size());
 132     patchSizes_.setSize(boundaryIds_.size());
 133
 134     label nCreatedFaces = nInternalFaces_;
 135     label baffleOffset  = cFaces.size();
 136     interfaces_.setSize(baffleIds_.size());
 137     nBoundaryFaces = 0;
 138
 139     forAll(boundaryIds_, patchI)
 140     {
 141         const List<cellFaceIdentifier>& idList = boundaryIds_[patchI];
 142
 143         patchStarts_[patchI] = nCreatedFaces;
 144
 145         // write each baffle side separately
 146         if (patchPhysicalTypes_[patchI] == "baffle")
 147         {
 148             label count = 0;
 149
 150             for (label side = 0; side < 2; ++side)
 151             {
 152                 label position = nInterfaces;
 153
 154                 forAll(idList, bndI)
 155                 {
 156                     label baffleI = idList[bndI].cell - baffleOffset;
 157
 158                     if
 159                     (
 160                         baffleI >= 0
 161                      && baffleI < baffleFaces_.size()
 162                      && baffleIds_[baffleI].size()
 163                     )
 164                     {
 165                         addPolyBoundaryFace
 166                         (
 167                             baffleIds_[baffleI][side],
 168                             nCreatedFaces
 169                         );
 170
 171                         // remove applied boundaries (2nd pass)
 172                         if (side == 1)
 173                         {
 174                             baffleIds_[baffleI].clear();
 175                         }
 176
 177                         interfaces_[position][side] = nCreatedFaces;
 178
 179                         nBoundaryFaces++;
 180                         nCreatedFaces++;
 181                         position++;
 182                         count++;
 183                     }
 184                 }
 185             }
 186
 187             nInterfaces += (count - (count % 2)) / 2;
 188         }
 189         else if (patchPhysicalTypes_[patchI] == "monitoring")
 190         {
 191             // translate the "monitoring" pseudo-boundaries to face sets
 192             List<label> monitoring(idList.size());
 193
 194             label monitorI = 0;
 195             forAll(idList, bndI)
 196             {
 197                 label cellId = idList[bndI].cell;
 198                 label faceId = idList[bndI].face;
 199
 200                 // standard case: volume cells
 201                 if (cellId < baffleOffset)
 202                 {
 203                     label faceNr = cellPolys_[cellId][faceId];
 204                     if (faceNr >= 0)
 205                     {
 206                         monitoring[monitorI++] = faceNr;
 207                     }
 208                 }
 209             }
 210
 211             monitoringSets_.insert(patchNames_[patchI], monitoring);
 212         }
 213         else
 214         {
 215             forAll(idList, bndI)
 216             {
 217                 // standard case: volume cells
 218                 if (idList[bndI].cell < baffleOffset)
 219                 {
 220                     addPolyBoundaryFace
 221                     (
 222                         idList[bndI],
 223                         nCreatedFaces
 224                     );
 225
 226                     nBoundaryFaces++;
 227                     nCreatedFaces++;
 228                 }
 229             }
 230         }
 231
 232         patchSizes_[patchI] = nCreatedFaces - patchStarts_[patchI];
 233     }
 234
 235     // add in missing faces
 236     Info<< "Missing faces added to patch after face "
 237         << nCreatedFaces << ":" <<endl;
 238     nMissingFaces = 0;
 239
 240     // look for baffles first - keep them together at the start of the patch
 241     for (label side = 0; side < 2; ++side)
 242     {
 243         label position = nInterfaces;
 244
 245         forAll(baffleIds_, baffleI)
 246         {
 247             if (baffleIds_[baffleI].size())
 248             {
 249                 // add each side for each baffle
 250                 addPolyBoundaryFace
 251                 (
 252                     baffleIds_[baffleI][side],
 253                     nCreatedFaces
 254                 );
 255
 256                 interfaces_[position][side] = nCreatedFaces;
 257
 258                 // remove applied boundaries (2nd pass)
 259                 if (side == 1)
 260                 {
 261                     baffleIds_[baffleI].clear();
 262                 }
 263
 264                 nMissingFaces++;
 265                 nCreatedFaces++;
 266                 position++;
 267             }
 268         }
 269     }
 270
 271     nInterfaces += (nMissingFaces - (nMissingFaces % 2)) / 2;
 272
 273     // scan for any other missing faces
 274     forAll(cellPolys_, cellI)
 275     {
 276         const labelList& curFaces = cellPolys_[cellI];
 277
 278         forAll(curFaces, cellFaceI)
 279         {
 280             if (curFaces[cellFaceI] < 0)
 281             {
 282                 // just report the first few
 283                 if (nMissingFaces < 4)
 284                 {
 285                     const face& thisFace = cFaces[cellI][cellFaceI];
 286
 287                     Info<< "  cell " << cellI << " face " << cellFaceI
 288                         << " (original cell " << origCellId_[cellI] << ")"
 289                         << " face: " << thisFace
 290                         << endl;
 291                 }
 292                 else if (nMissingFaces == 5)
 293                 {
 294                     Info<< "  ..." << nl << endl;
 295                 }
 296
 297                 addPolyBoundaryFace(cellI, cellFaceI, nCreatedFaces);
 298                 nMissingFaces++;
 299                 nCreatedFaces++;
 300             }
 301         }
 302     }
 303
 304     Info<< "Added " << nMissingFaces << " unmatched faces" << endl;
 305
 306     if (nMissingFaces > 0)
 307     {
 308         patchSizes_.last() = nMissingFaces;
 309     }
 310     else
 311     {
 312         patchStarts_.setSize(patchStarts_.size() - 1);
 313     }
 314
 315     // reset the size of the face list
 316     meshFaces_.setSize(nCreatedFaces);
 317
 318     // check the mesh for face mismatch
 319     // (faces addressed once or more than twice)
 320     labelList markupFaces(meshFaces_.size(), 0);
 321
 322     forAll(cellPolys_, cellI)
 323     {
 324         const labelList& curFaces = cellPolys_[cellI];
 325
 326         forAll(curFaces, faceI)
 327         {
 328             markupFaces[curFaces[faceI]]++;
 329         }
 330     }
 331
 332     for (label i = nInternalFaces_; i < markupFaces.size(); i++)
 333     {
 334         markupFaces[i]++;
 335     }
 336
 337     label nProblemFaces = 0;
 338
 339     forAll(markupFaces, faceI)
 340     {
 341         if (markupFaces[faceI] != 2)
 342         {
 343             const face& problemFace = meshFaces_[faceI];
 344
 345             Info<< "meshReader::createPolyBoundary() : "
 346                 << "problem with face " << faceI << ": addressed "
 347                 << markupFaces[faceI] << " times (should be 2!). Face: "
 348                 << problemFace << endl;
 349
 350             nProblemFaces++;
 351         }
 352     }
 353
 354     if (nProblemFaces > 0)
 355     {
 356         Info<< "Number of incorrectly matched faces: "
 357             << nProblemFaces << endl;
 358     }
 359
 360     // adjust for missing members
 361     if (nInterfaces < interfaces_.size())
 362     {
 363         interfaces_.setSize(nInterfaces);
 364     }
 365
 366     Info<< "Number of boundary faces: " << nBoundaryFaces << nl
 367         << "Total number of faces: " << nCreatedFaces << nl
 368         << "Number of interfaces: " << nInterfaces << endl;
 369 }
 370
 371
 372 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 373
 374 Foam::List<Foam::polyPatch*>
 375 Foam::meshReader::polyBoundaryPatches(const polyMesh& mesh)
 376 {
 377     label nUsed = 0, nEmpty = 0;
 378     label nPatches = patchStarts_.size();
 379
 380     // avoid empty patches - move to the end of the lists and truncate
 381     labelList oldToNew = identity(nPatches);
 382     forAll(patchSizes_, patchI)
 383     {
 384         if (patchSizes_[patchI] > 0)
 385         {
 386             oldToNew[patchI] = nUsed++;
 387         }
 388         else
 389         {
 390             nEmpty++;
 391             oldToNew[patchI] = nPatches - nEmpty;
 392         }
 393     }
 394
 395     nPatches = nUsed;
 396
 397     if (nEmpty)
 398     {
 399         Info<< "Removing " << nEmpty << " empty patches" << endl;
 400
 401         inplaceReorder(oldToNew, patchTypes_);
 402         inplaceReorder(oldToNew, patchNames_);
 403         inplaceReorder(oldToNew, patchStarts_);
 404         inplaceReorder(oldToNew, patchSizes_);
 405     }
 406
 407     patchTypes_.setSize(nPatches);
 408     patchNames_.setSize(nPatches);
 409     patchStarts_.setSize(nPatches);
 410     patchSizes_.setSize(nPatches);
 411
 412
 413     List<polyPatch*> p(nPatches);
 414
 415     // All patch dictionaries
 416     PtrList<dictionary> patchDicts(patchNames_.size());
 417     // Default boundary patch types
 418     word defaultFacesType(emptyPolyPatch::typeName);
 419
 420     // we could consider dropping this entirely
 421     preservePatchTypes
 422     (
 423         mesh,
 424         mesh.instance(),
 425         mesh.meshDir(),
 426         patchNames_,
 427         patchDicts,
 428         "defaultFaces",
 429         defaultFacesType
 430     );
 431     forAll(patchDicts, patchI)
 432     {
 433         if (!patchDicts.set(patchI))
 434         {
 435             patchDicts.set(patchI, new dictionary());
 436         }
 437         dictionary& patchDict = patchDicts[patchI];
 438
 439         // add but not overwrite type
 440         patchDict.add("type", patchTypes_[patchI], false);
 441         if (patchPhysicalTypes_.size() && patchPhysicalTypes_[patchI].size())
 442         {
 443             patchDict.add("startFace", patchPhysicalTypes_[patchI], false);
 444         }
 445
 446         // overwrite sizes and start
 447         patchDict.add("nFaces", patchSizes_[patchI], true);
 448         patchDict.add("startFace", patchStarts_[patchI], true);
 449     }
 450
 451
 452     forAll(patchStarts_, patchI)
 453     {
 454         p[patchI] = polyPatch::New
 455         (
 456             patchNames_[patchI],
 457             patchDicts[patchI],
 458             patchI,
 459             mesh.boundaryMesh()
 460         ).ptr();
 461     }
 462
 463     return p;
 464 }
 465
 466
 467 // ************************************************************************* //