src/foam/meshes/polyMesh/zones/faceZone/faceZone.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend 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     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25     A subset of mesh faces.
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "faceZone.H"
  30 #include "addToRunTimeSelectionTable.H"
  31 #include "faceZoneMesh.H"
  32 #include "polyMesh.H"
  33 #include "primitiveMesh.H"
  34 #include "demandDrivenData.H"
  35 #include "mapPolyMesh.H"
  36 #include "syncTools.H"
  37
  38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  39
  40 namespace Foam
  41 {
  42     defineTypeNameAndDebug(faceZone, 0);
  43     defineRunTimeSelectionTable(faceZone, dictionary);
  44     addToRunTimeSelectionTable(faceZone, faceZone, dictionary);
  45 }
  46
  47 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  48
  49 void Foam::faceZone::calcFaceZonePatch() const
  50 {
  51     if (debug)
  52     {
  53         Info<< "void faceZone::calcFaceZonePatch() const : "
  54             << "Calculating primitive patch"
  55             << endl;
  56     }
  57
  58     if (patchPtr_)
  59     {
  60         FatalErrorIn
  61         (
  62             "void faceZone::calcFaceZonePatch() const"
  63         )   << "primitive face zone patch already calculated"
  64             << abort(FatalError);
  65     }
  66
  67     patchPtr_ =
  68         new primitiveFacePatch
  69         (
  70             faceList(size()),
  71             zoneMesh().mesh().allPoints()
  72         );
  73
  74     primitiveFacePatch& patch = *patchPtr_;
  75
  76     const faceList& f = zoneMesh().mesh().allFaces();
  77
  78     const labelList& addr = *this;
  79     const boolList& flip = flipMap();
  80
  81     forAll (addr, faceI)
  82     {
  83         if (flip[faceI])
  84         {
  85             patch[faceI] = f[addr[faceI]].reverseFace();
  86         }
  87         else
  88         {
  89             patch[faceI] = f[addr[faceI]];
  90         }
  91     }
  92
  93     if (debug)
  94     {
  95         Info<< "void faceZone::calcFaceZonePatch() const : "
  96             << "Finished calculating primitive patch"
  97             << endl;
  98     }
  99 }
 100
 101
 102 const Foam::Map<Foam::label>& Foam::faceZone::faceLookupMap() const
 103 {
 104     if (!faceLookupMapPtr_)
 105     {
 106         calcFaceLookupMap();
 107     }
 108
 109     return *faceLookupMapPtr_;
 110 }
 111
 112
 113 void Foam::faceZone::calcFaceLookupMap() const
 114 {
 115     if (debug)
 116     {
 117         Info<< "void faceZone::calcFaceLookupMap() const : "
 118             << "Calculating face lookup map"
 119             << endl;
 120     }
 121
 122     if (faceLookupMapPtr_)
 123     {
 124         FatalErrorIn
 125         (
 126             "void faceZone::calcFaceLookupMap() const"
 127         )   << "face lookup map already calculated"
 128             << abort(FatalError);
 129     }
 130
 131     const labelList& addr = *this;
 132
 133     faceLookupMapPtr_ = new Map<label>(2*addr.size());
 134     Map<label>& flm = *faceLookupMapPtr_;
 135
 136     forAll (addr, faceI)
 137     {
 138         flm.insert(addr[faceI], faceI);
 139     }
 140
 141     if (debug)
 142     {
 143         Info<< "void faceZone::calcFaceLookupMap() const : "
 144             << "Finished calculating face lookup map"
 145             << endl;
 146     }
 147 }
 148
 149
 150 void Foam::faceZone::calcCellLayers() const
 151 {
 152     if (debug)
 153     {
 154         Info<< "void Foam::faceZone::calcCellLayers() const : "
 155             << "calculating master cells"
 156             << endl;
 157     }
 158
 159     // It is an error to attempt to recalculate edgeCells
 160     // if the pointer is already set
 161     if (masterCellsPtr_ || slaveCellsPtr_)
 162     {
 163         FatalErrorIn("void faceZone::calcCellLayers() const")
 164             << "cell layers already calculated"
 165             << abort(FatalError);
 166     }
 167     else
 168     {
 169         // Go through all the faces in the master zone.  Choose the
 170         // master or slave cell based on the face flip
 171         const polyMesh& mesh = zoneMesh().mesh();
 172
 173         const labelList& own = mesh.faceOwner();
 174         const labelList& nei = mesh.faceNeighbour();
 175
 176         const labelList& mf = *this;
 177
 178         const boolList& faceFlip = flipMap();
 179
 180         masterCellsPtr_ = new labelList(mf.size());
 181         labelList& mc = *masterCellsPtr_;
 182
 183         slaveCellsPtr_ = new labelList(mf.size());
 184         labelList& sc = *slaveCellsPtr_;
 185
 186         forAll (mf, faceI)
 187         {
 188             label curMc = -1;
 189             label curSc = -1;
 190
 191             if (!faceFlip[faceI])
 192             {
 193                 // Face is oriented correctly, no flip needed
 194                 // Master is the neighbour, with normal pointing into it
 195                 // Bug fix, HJ, 7/Mar/2011
 196                 if (mesh.isInternalFace(mf[faceI]))
 197                 {
 198                     curMc = nei[mf[faceI]];
 199                     curSc = own[mf[faceI]];
 200                 }
 201                 else if (mf[faceI] < mesh.nFaces())
 202                 {
 203                     curMc = -1;
 204                     curSc = own[mf[faceI]];
 205                 }
 206             }
 207             else
 208             {
 209                 // Face flip
 210                 // Master is the owner, with normal pointing into it
 211                 // Bug fix, HJ, 7/Mar/2011
 212                 if (mesh.isInternalFace(mf[faceI]))
 213                 {
 214                     curMc = own[mf[faceI]];
 215                     curSc = nei[mf[faceI]];
 216                 }
 217                 else if (mf[faceI] < mesh.nFaces())
 218                 {
 219                     curMc = own[mf[faceI]];
 220                     curSc = -1;
 221                 }
 222
 223             }
 224
 225             mc[faceI] = curMc;
 226             sc[faceI] = curSc;
 227         }
 228     }
 229 }
 230
 231
 232 void Foam::faceZone::checkAddressing() const
 233 {
 234     if (size() != flipMap_.size())
 235     {
 236         FatalErrorIn("void Foam::faceZone::checkAddressing() const")
 237             << "Different sizes of the addressing and flipMap arrays.  "
 238             << "Size of addressing: " << size()
 239             << " size of flip map: " << flipMap_.size()
 240             << abort(FatalError);
 241     }
 242 }
 243
 244
 245 void Foam::faceZone::clearAddressing()
 246 {
 247     deleteDemandDrivenData(patchPtr_);
 248
 249     deleteDemandDrivenData(masterCellsPtr_);
 250     deleteDemandDrivenData(slaveCellsPtr_);
 251
 252     deleteDemandDrivenData(mePtr_);
 253     deleteDemandDrivenData(faceLookupMapPtr_);
 254 }
 255
 256
 257 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 258
 259 // Construct from components
 260 Foam::faceZone::faceZone
 261 (
 262     const word& name,
 263     const labelList& addr,
 264     const boolList& fm,
 265     const label index,
 266     const faceZoneMesh& zm
 267 )
 268 :
 269     labelList(addr),
 270     name_(name),
 271     flipMap_(fm),
 272     index_(index),
 273     zoneMesh_(zm),
 274     patchPtr_(NULL),
 275     masterCellsPtr_(NULL),
 276     slaveCellsPtr_(NULL),
 277     mePtr_(NULL),
 278     faceLookupMapPtr_(NULL)
 279 {
 280     checkAddressing();
 281 }
 282
 283
 284 Foam::faceZone::faceZone
 285 (
 286     const word& name,
 287     const Xfer<labelList>& addr,
 288     const Xfer<boolList>& fm,
 289     const label index,
 290     const faceZoneMesh& zm
 291 )
 292 :
 293     labelList(addr),
 294     name_(name),
 295     flipMap_(fm),
 296     index_(index),
 297     zoneMesh_(zm),
 298     patchPtr_(NULL),
 299     masterCellsPtr_(NULL),
 300     slaveCellsPtr_(NULL),
 301     mePtr_(NULL),
 302     faceLookupMapPtr_(NULL)
 303 {
 304     checkAddressing();
 305 }
 306
 307
 308 // Construct from dictionary
 309 Foam::faceZone::faceZone
 310 (
 311     const word& name,
 312     const dictionary& dict,
 313     const label index,
 314     const faceZoneMesh& zm
 315 )
 316 :
 317     labelList(dict.lookup("faceLabels")),
 318     name_(name),
 319     flipMap_(dict.lookup("flipMap")),
 320     index_(index),
 321     zoneMesh_(zm),
 322     patchPtr_(NULL),
 323     masterCellsPtr_(NULL),
 324     slaveCellsPtr_(NULL),
 325     mePtr_(NULL),
 326     faceLookupMapPtr_(NULL)
 327 {
 328     checkAddressing();
 329 }
 330
 331
 332 // Construct given the original zone and resetting the
 333 // face list and zone mesh information
 334 Foam::faceZone::faceZone
 335 (
 336     const faceZone& fz,
 337     const labelList& addr,
 338     const boolList& fm,
 339     const label index,
 340     const faceZoneMesh& zm
 341 )
 342 :
 343     labelList(addr),
 344     name_(fz.name()),
 345     flipMap_(fm),
 346     index_(index),
 347     zoneMesh_(zm),
 348     patchPtr_(NULL),
 349     masterCellsPtr_(NULL),
 350     slaveCellsPtr_(NULL),
 351     mePtr_(NULL),
 352     faceLookupMapPtr_(NULL)
 353 {
 354     checkAddressing();
 355 }
 356
 357
 358 Foam::faceZone::faceZone
 359 (
 360     const faceZone& fz,
 361     const Xfer<labelList>& addr,
 362     const Xfer<boolList>& fm,
 363     const label index,
 364     const faceZoneMesh& zm
 365 )
 366 :
 367     labelList(addr),
 368     name_(fz.name()),
 369     flipMap_(fm),
 370     index_(index),
 371     zoneMesh_(zm),
 372     patchPtr_(NULL),
 373     masterCellsPtr_(NULL),
 374     slaveCellsPtr_(NULL),
 375     mePtr_(NULL),
 376     faceLookupMapPtr_(NULL)
 377 {
 378     checkAddressing();
 379 }
 380
 381
 382 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 383
 384 Foam::faceZone::~faceZone()
 385 {
 386     clearAddressing();
 387 }
 388
 389
 390 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 391
 392 Foam::label Foam::faceZone::whichFace(const label globalFaceID) const
 393 {
 394     const Map<label>& flm = faceLookupMap();
 395
 396     Map<label>::const_iterator flmIter = flm.find(globalFaceID);
 397
 398     if (flmIter == flm.end())
 399     {
 400         return -1;
 401     }
 402     else
 403     {
 404         return flmIter();
 405     }
 406 }
 407
 408
 409 const Foam::faceZoneMesh& Foam::faceZone::zoneMesh() const
 410 {
 411     return zoneMesh_;
 412 }
 413
 414
 415 const Foam::primitiveFacePatch& Foam::faceZone::operator()() const
 416 {
 417     if (!patchPtr_)
 418     {
 419         calcFaceZonePatch();
 420     }
 421
 422     return *patchPtr_;
 423 }
 424
 425
 426 const Foam::labelList& Foam::faceZone::masterCells() const
 427 {
 428     if (!masterCellsPtr_)
 429     {
 430         calcCellLayers();
 431     }
 432
 433     return *masterCellsPtr_;
 434 }
 435
 436
 437 const Foam::labelList& Foam::faceZone::slaveCells() const
 438 {
 439     if (!slaveCellsPtr_)
 440     {
 441         calcCellLayers();
 442     }
 443
 444     return *slaveCellsPtr_;
 445 }
 446
 447
 448 const Foam::labelList& Foam::faceZone::meshEdges() const
 449 {
 450     if (!mePtr_)
 451     {
 452         // Old form.  Under testing: merge, HJ, 17/Aug/2010
 453         //labelList faceCells(size());
 454         //
 455         //const labelList& own = zoneMesh().mesh().faceOwner();
 456         //
 457         //const labelList& faceLabels = *this;
 458         //
 459         //forAll (faceCells, faceI)
 460         //{
 461         //    faceCells[faceI] = own[faceLabels[faceI]];
 462         //}
 463         //
 464         //mePtr_ =
 465         //    new labelList
 466         //    (
 467         //        operator()().meshEdges
 468         //        (
 469         //            zoneMesh().mesh().edges(),
 470         //            zoneMesh().mesh().cellEdges(),
 471         //            faceCells
 472         //        )
 473         //    );
 474
 475         mePtr_ =
 476             new labelList
 477             (
 478                 operator()().meshEdges
 479                 (
 480                     zoneMesh().mesh().edges(),
 481                     zoneMesh().mesh().pointEdges()
 482                 )
 483             );
 484     }
 485
 486     return *mePtr_;
 487 }
 488
 489
 490 void Foam::faceZone::resetAddressing
 491 (
 492     const labelList& addr,
 493     const boolList& flipMap
 494 )
 495 {
 496     clearAddressing();
 497     labelList::operator=(addr);
 498     flipMap_ = flipMap;
 499 }
 500
 501
 502 void Foam::faceZone::updateMesh()
 503 {
 504     clearAddressing();
 505 }
 506
 507
 508 bool Foam::faceZone::checkDefinition(const bool report) const
 509 {
 510     const labelList& addr = *this;
 511
 512     bool boundaryError = false;
 513
 514     forAll(addr, i)
 515     {
 516         if (addr[i] < 0 || addr[i] >= zoneMesh().mesh().allFaces().size())
 517         {
 518             boundaryError = true;
 519
 520             if (report)
 521             {
 522                 SeriousErrorIn
 523                 (
 524                     "bool faceZone::checkDefinition("
 525                     "const bool report) const"
 526                 )   << "Zone " << name()
 527                     << " contains invalid face label " << addr[i] << nl
 528                     << "Valid face labels are 0.."
 529                     << zoneMesh().mesh().allFaces().size() - 1 << endl;
 530             }
 531         }
 532     }
 533     return boundaryError;
 534 }
 535
 536
 537 bool Foam::faceZone::checkParallelSync(const bool report) const
 538 {
 539     const polyMesh& mesh = zoneMesh().mesh();
 540     const polyBoundaryMesh& bm = mesh.boundaryMesh();
 541
 542     bool boundaryError = false;
 543
 544
 545     // Check that zone faces are synced
 546     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 547
 548     {
 549         boolList neiZoneFace(mesh.nFaces() - mesh.nInternalFaces(), false);
 550         boolList neiZoneFlip(mesh.nFaces() - mesh.nInternalFaces(), false);
 551
 552         forAll(*this, i)
 553         {
 554             label faceI = operator[](i);
 555
 556             // Only check live faces
 557             if (faceI < mesh.nFaces())
 558             {
 559                 if (!mesh.isInternalFace(faceI))
 560                 {
 561                     neiZoneFace[faceI - mesh.nInternalFaces()] = true;
 562                     neiZoneFlip[faceI - mesh.nInternalFaces()] = flipMap()[i];
 563                 }
 564             }
 565         }
 566
 567         boolList myZoneFace(neiZoneFace);
 568         syncTools::swapBoundaryFaceList(mesh, neiZoneFace, false);
 569         boolList myZoneFlip(neiZoneFlip);
 570         syncTools::swapBoundaryFaceList(mesh, neiZoneFlip, false);
 571
 572         forAll(*this, i)
 573         {
 574             label faceI = operator[](i);
 575
 576             // Only check live faces
 577             if (faceI < mesh.nFaces())
 578             {
 579                 label patchI = bm.whichPatch(faceI);
 580
 581                 if (patchI != -1 && isA<cyclicPolyPatch>(bm[patchI]))
 582                 {
 583                     label bFaceI = faceI-mesh.nInternalFaces();
 584
 585                     // Check face in zone on both sides
 586                     if (myZoneFace[bFaceI] != neiZoneFace[bFaceI])
 587                     {
 588                         boundaryError = true;
 589
 590                         if (report)
 591                         {
 592                             Pout<< " ***Problem with faceZone " << index()
 593                                 << " named " << name()
 594                                 << ". Face " << faceI
 595                                 << " on coupled patch "
 596                                 << bm[patchI].name()
 597                                 << " is not consistent with its "
 598                                 << "coupled neighbour."
 599                                 << endl;
 600                         }
 601                     }
 602
 603                     // Flip state should be opposite.
 604                     if (myZoneFlip[bFaceI] == neiZoneFlip[bFaceI])
 605                     {
 606                         boundaryError = true;
 607
 608                         if (report)
 609                         {
 610                             Pout<< " ***Problem with faceZone " << index()
 611                                 << " named " << name()
 612                                 << ". Face " << faceI
 613                                 << " on coupled patch "
 614                                 << bm[patchI].name()
 615                                 << " does not have consistent flipMap"
 616                                 << " across coupled faces."
 617                                 << endl;
 618                         }
 619                     }
 620                 }
 621             }
 622         }
 623     }
 624
 625     return returnReduce(boundaryError, orOp<bool>());
 626 }
 627
 628
 629 void Foam::faceZone::movePoints(const pointField& p)
 630 {
 631     if (patchPtr_)
 632     {
 633         patchPtr_->movePoints(p);
 634     }
 635 }
 636
 637 void Foam::faceZone::write(Ostream& os) const
 638 {
 639     os  << nl << name()
 640         << nl << static_cast<const labelList&>(*this)
 641         << nl << flipMap();
 642 }
 643
 644
 645 void Foam::faceZone::writeDict(Ostream& os) const
 646 {
 647     os  << nl << name() << nl << token::BEGIN_BLOCK << incrIndent << nl
 648         << indent << "type " << type() << token::END_STATEMENT << nl;
 649
 650     writeEntry("faceLabels", os);
 651     flipMap().writeEntry("flipMap", os);
 652
 653     os  << decrIndent << token::END_BLOCK << endl;
 654 }
 655
 656
 657 // * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
 658
 659 Foam::Ostream& Foam::operator<<(Ostream& os, const faceZone& p)
 660 {
 661     p.write(os);
 662     os.check("Ostream& operator<<(Ostream& f, const faceZone& p");
 663     return os;
 664 }
 665
 666
 667 // ************************************************************************* //