FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

faceZone.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 Description
00025     A subset of mesh faces.
00026 
00027 \*---------------------------------------------------------------------------*/
00028 
00029 #include "faceZone.H"
00030 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00031 #include <OpenFOAM/faceZoneMesh.H>
00032 #include <OpenFOAM/polyMesh.H>
00033 #include <OpenFOAM/primitiveMesh.H>
00034 #include <OpenFOAM/demandDrivenData.H>
00035 #include <OpenFOAM/mapPolyMesh.H>
00036 #include <OpenFOAM/syncTools.H>
00037 
00038 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00039 
00040 namespace Foam
00041 {
00042     defineTypeNameAndDebug(faceZone, 0);
00043     defineRunTimeSelectionTable(faceZone, dictionary);
00044     addToRunTimeSelectionTable(faceZone, faceZone, dictionary);
00045 }
00046 
00047 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00048 
00049 void Foam::faceZone::calcFaceZonePatch() const
00050 {
00051     if (debug)
00052     {
00053         Info<< "void faceZone::calcFaceZonePatch() const : "
00054             << "Calculating primitive patch"
00055             << endl;
00056     }
00057 
00058     if (patchPtr_)
00059     {
00060         FatalErrorIn
00061         (
00062             "void faceZone::calcFaceZonePatch() const"
00063         )   << "primitive face zone patch already calculated"
00064             << abort(FatalError);
00065     }
00066 
00067     patchPtr_ =
00068         new primitiveFacePatch
00069         (
00070             faceList(size()),
00071             zoneMesh().mesh().points()
00072         );
00073 
00074     primitiveFacePatch& patch = *patchPtr_;
00075 
00076     const faceList& f = zoneMesh().mesh().faces();
00077 
00078     const labelList& addr = *this;
00079     const boolList& flip = flipMap();
00080 
00081     forAll (addr, faceI)
00082     {
00083         if (flip[faceI])
00084         {
00085             patch[faceI] = f[addr[faceI]].reverseFace();
00086         }
00087         else
00088         {
00089             patch[faceI] = f[addr[faceI]];
00090         }
00091     }
00092 
00093     if (debug)
00094     {
00095         Info<< "void faceZone::calcFaceZonePatch() const : "
00096             << "Finished calculating primitive patch"
00097             << endl;
00098     }
00099 }
00100 
00101 
00102 const Foam::Map<Foam::label>& Foam::faceZone::faceLookupMap() const
00103 {
00104     if (!faceLookupMapPtr_)
00105     {
00106         calcFaceLookupMap();
00107     }
00108 
00109     return *faceLookupMapPtr_;
00110 }
00111 
00112 
00113 void Foam::faceZone::calcFaceLookupMap() const
00114 {
00115     if (debug)
00116     {
00117         Info<< "void faceZone::calcFaceLookupMap() const : "
00118             << "Calculating face lookup map"
00119             << endl;
00120     }
00121 
00122     if (faceLookupMapPtr_)
00123     {
00124         FatalErrorIn
00125         (
00126             "void faceZone::calcFaceLookupMap() const"
00127         )   << "face lookup map already calculated"
00128             << abort(FatalError);
00129     }
00130 
00131     const labelList& addr = *this;
00132 
00133     faceLookupMapPtr_ = new Map<label>(2*addr.size());
00134     Map<label>& flm = *faceLookupMapPtr_;
00135 
00136     forAll (addr, faceI)
00137     {
00138         flm.insert(addr[faceI], faceI);
00139     }
00140 
00141     if (debug)
00142     {
00143         Info<< "void faceZone::calcFaceLookupMap() const : "
00144             << "Finished calculating face lookup map"
00145             << endl;
00146     }
00147 }
00148 
00149 
00150 void Foam::faceZone::calcCellLayers() const
00151 {
00152     if (debug)
00153     {
00154         Info<< "void Foam::faceZone::calcCellLayers() const : "
00155             << "calculating master cells"
00156             << endl;
00157     }
00158 
00159     // It is an error to attempt to recalculate edgeCells
00160     // if the pointer is already set
00161     if (masterCellsPtr_ || slaveCellsPtr_)
00162     {
00163         FatalErrorIn("void faceZone::calcCellLayers() const")
00164             << "cell layers already calculated"
00165             << abort(FatalError);
00166     }
00167     else
00168     {
00169         // Go through all the faces in the master zone.  Choose the
00170         // master or slave cell based on the face flip
00171 
00172         const labelList& own = zoneMesh().mesh().faceOwner();
00173         const labelList& nei = zoneMesh().mesh().faceNeighbour();
00174 
00175         const labelList& mf = *this;
00176 
00177         const boolList& faceFlip = flipMap();
00178 
00179         masterCellsPtr_ = new labelList(mf.size());
00180         labelList& mc = *masterCellsPtr_;
00181 
00182         slaveCellsPtr_ = new labelList(mf.size());
00183         labelList& sc = *slaveCellsPtr_;
00184 
00185         forAll (mf, faceI)
00186         {
00187             label ownCellI = own[mf[faceI]];
00188             label neiCellI =
00189             (
00190                 zoneMesh().mesh().isInternalFace(mf[faceI])
00191               ? nei[mf[faceI]]
00192               : -1
00193             );
00194 
00195             if (!faceFlip[faceI])
00196             {
00197                 // Face is oriented correctly, no flip needed
00198                 mc[faceI] = neiCellI;
00199                 sc[faceI] = ownCellI;
00200             }
00201             else
00202             {
00203                 mc[faceI] = ownCellI;
00204                 sc[faceI] = neiCellI;
00205             }
00206         }
00207         //Info << "masterCells: " << mc << endl;
00208         //Info << "slaveCells: " << sc << endl;
00209     }
00210 }
00211 
00212 
00213 void Foam::faceZone::checkAddressing() const
00214 {
00215     if (size() != flipMap_.size())
00216     {
00217         FatalErrorIn("void Foam::faceZone::checkAddressing() const")
00218             << "Different sizes of the addressing and flipMap arrays.  "
00219             << "Size of addressing: " << size()
00220             << " size of flip map: " << flipMap_.size()
00221             << abort(FatalError);
00222     }
00223 }
00224 
00225 
00226 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00227 
00228 // Construct from components
00229 Foam::faceZone::faceZone
00230 (
00231     const word& name,
00232     const labelList& addr,
00233     const boolList& fm,
00234     const label index,
00235     const faceZoneMesh& zm
00236 )
00237 :
00238     labelList(addr),
00239     name_(name),
00240     flipMap_(fm),
00241     index_(index),
00242     zoneMesh_(zm),
00243     patchPtr_(NULL),
00244     masterCellsPtr_(NULL),
00245     slaveCellsPtr_(NULL),
00246     mePtr_(NULL),
00247     faceLookupMapPtr_(NULL)
00248 {
00249     checkAddressing();
00250 }
00251 
00252 
00253 Foam::faceZone::faceZone
00254 (
00255     const word& name,
00256     const Xfer<labelList>& addr,
00257     const Xfer<boolList>& fm,
00258     const label index,
00259     const faceZoneMesh& zm
00260 )
00261 :
00262     labelList(addr),
00263     name_(name),
00264     flipMap_(fm),
00265     index_(index),
00266     zoneMesh_(zm),
00267     patchPtr_(NULL),
00268     masterCellsPtr_(NULL),
00269     slaveCellsPtr_(NULL),
00270     mePtr_(NULL),
00271     faceLookupMapPtr_(NULL)
00272 {
00273     checkAddressing();
00274 }
00275 
00276 
00277 // Construct from dictionary
00278 Foam::faceZone::faceZone
00279 (
00280     const word& name,
00281     const dictionary& dict,
00282     const label index,
00283     const faceZoneMesh& zm
00284 )
00285 :
00286     labelList(dict.lookup("faceLabels")),
00287     name_(name),
00288     flipMap_(dict.lookup("flipMap")),
00289     index_(index),
00290     zoneMesh_(zm),
00291     patchPtr_(NULL),
00292     masterCellsPtr_(NULL),
00293     slaveCellsPtr_(NULL),
00294     mePtr_(NULL),
00295     faceLookupMapPtr_(NULL)
00296 {
00297     checkAddressing();
00298 }
00299 
00300 
00301 // Construct given the original zone and resetting the
00302 // face list and zone mesh information
00303 Foam::faceZone::faceZone
00304 (
00305     const faceZone& fz,
00306     const labelList& addr,
00307     const boolList& fm,
00308     const label index,
00309     const faceZoneMesh& zm
00310 )
00311 :
00312     labelList(addr),
00313     name_(fz.name()),
00314     flipMap_(fm),
00315     index_(index),
00316     zoneMesh_(zm),
00317     patchPtr_(NULL),
00318     masterCellsPtr_(NULL),
00319     slaveCellsPtr_(NULL),
00320     mePtr_(NULL),
00321     faceLookupMapPtr_(NULL)
00322 {
00323     checkAddressing();
00324 }
00325 
00326 
00327 Foam::faceZone::faceZone
00328 (
00329     const faceZone& fz,
00330     const Xfer<labelList>& addr,
00331     const Xfer<boolList>& fm,
00332     const label index,
00333     const faceZoneMesh& zm
00334 )
00335 :
00336     labelList(addr),
00337     name_(fz.name()),
00338     flipMap_(fm),
00339     index_(index),
00340     zoneMesh_(zm),
00341     patchPtr_(NULL),
00342     masterCellsPtr_(NULL),
00343     slaveCellsPtr_(NULL),
00344     mePtr_(NULL),
00345     faceLookupMapPtr_(NULL)
00346 {
00347     checkAddressing();
00348 }
00349 
00350 
00351 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00352 
00353 Foam::faceZone::~faceZone()
00354 {
00355     clearAddressing();
00356 }
00357 
00358 
00359 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00360 
00361 Foam::label Foam::faceZone::whichFace(const label globalFaceID) const
00362 {
00363     const Map<label>& flm = faceLookupMap();
00364 
00365     Map<label>::const_iterator flmIter = flm.find(globalFaceID);
00366 
00367     if (flmIter == flm.end())
00368     {
00369         return -1;
00370     }
00371     else
00372     {
00373         return flmIter();
00374     }
00375 }
00376 
00377 
00378 const Foam::faceZoneMesh& Foam::faceZone::zoneMesh() const
00379 {
00380     return zoneMesh_;
00381 }
00382 
00383 
00384 const Foam::primitiveFacePatch& Foam::faceZone::operator()() const
00385 {
00386     if (!patchPtr_)
00387     {
00388         calcFaceZonePatch();
00389     }
00390 
00391     return *patchPtr_;
00392 }
00393 
00394 
00395 const Foam::labelList& Foam::faceZone::masterCells() const
00396 {
00397     if (!masterCellsPtr_)
00398     {
00399         calcCellLayers();
00400     }
00401 
00402     return *masterCellsPtr_;
00403 }
00404 
00405 
00406 const Foam::labelList& Foam::faceZone::slaveCells() const
00407 {
00408     if (!slaveCellsPtr_)
00409     {
00410         calcCellLayers();
00411     }
00412 
00413     return *slaveCellsPtr_;
00414 }
00415 
00416 
00417 const Foam::labelList& Foam::faceZone::meshEdges() const
00418 {
00419     if (!mePtr_)
00420     {
00421         //labelList faceCells(size());
00422         //
00423         //const labelList& own = zoneMesh().mesh().faceOwner();
00424         //
00425         //const labelList& faceLabels = *this;
00426         //
00427         //forAll (faceCells, faceI)
00428         //{
00429         //    faceCells[faceI] = own[faceLabels[faceI]];
00430         //}
00431         //
00432         //mePtr_ =
00433         //    new labelList
00434         //    (
00435         //        operator()().meshEdges
00436         //        (
00437         //            zoneMesh().mesh().edges(),
00438         //            zoneMesh().mesh().cellEdges(),
00439         //            faceCells
00440         //        )
00441         //    );
00442 
00443         mePtr_ =
00444             new labelList
00445             (
00446                 operator()().meshEdges
00447                 (
00448                     zoneMesh().mesh().edges(),
00449                     zoneMesh().mesh().pointEdges()
00450                 )
00451             );
00452     }
00453 
00454     return *mePtr_;
00455 }
00456 
00457 
00458 void Foam::faceZone::clearAddressing()
00459 {
00460     deleteDemandDrivenData(patchPtr_);
00461 
00462     deleteDemandDrivenData(masterCellsPtr_);
00463     deleteDemandDrivenData(slaveCellsPtr_);
00464 
00465     deleteDemandDrivenData(mePtr_);
00466     deleteDemandDrivenData(faceLookupMapPtr_);
00467 }
00468 
00469 
00470 void Foam::faceZone::resetAddressing
00471 (
00472     const labelList& addr,
00473     const boolList& flipMap
00474 )
00475 {
00476     clearAddressing();
00477     labelList::operator=(addr);
00478     flipMap_ = flipMap;
00479 }
00480 
00481 
00482 void Foam::faceZone::updateMesh(const mapPolyMesh& mpm)
00483 {
00484     clearAddressing();
00485 
00486     labelList newAddressing(size());
00487     boolList newFlipMap(flipMap_.size());
00488     label nFaces = 0;
00489 
00490     const labelList& faceMap = mpm.reverseFaceMap();
00491 
00492     forAll(*this, i)
00493     {
00494         label faceI = operator[](i);
00495 
00496         if (faceMap[faceI] >= 0)
00497         {
00498             newAddressing[nFaces] = faceMap[faceI];
00499             newFlipMap[nFaces] = flipMap_[i];       // Keep flip map.
00500             nFaces++;
00501         }
00502     }
00503 
00504     newAddressing.setSize(nFaces);
00505     newFlipMap.setSize(nFaces);
00506 
00507     transfer(newAddressing);
00508     flipMap_.transfer(newFlipMap);
00509 }
00510 
00511 
00512 bool Foam::faceZone::checkDefinition(const bool report) const
00513 {
00514     const labelList& addr = *this;
00515 
00516     bool boundaryError = false;
00517 
00518     forAll(addr, i)
00519     {
00520         if (addr[i] < 0 || addr[i] >= zoneMesh().mesh().faces().size())
00521         {
00522             boundaryError = true;
00523 
00524             if (report)
00525             {
00526                 SeriousErrorIn
00527                 (
00528                     "bool faceZone::checkDefinition("
00529                     "const bool report) const"
00530                 )   << "Zone " << name()
00531                     << " contains invalid face label " << addr[i] << nl
00532                     << "Valid face labels are 0.."
00533                     << zoneMesh().mesh().faces().size()-1 << endl;
00534             }
00535         }
00536     }
00537     return boundaryError;
00538 }
00539 
00540 
00541 bool Foam::faceZone::checkParallelSync(const bool report) const
00542 {
00543     const polyMesh& mesh = zoneMesh().mesh();
00544     const polyBoundaryMesh& bm = mesh.boundaryMesh();
00545 
00546     bool boundaryError = false;
00547 
00548 
00549     // Check that zone faces are synced
00550     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00551 
00552     {
00553         boolList neiZoneFace(mesh.nFaces()-mesh.nInternalFaces(), false);
00554         boolList neiZoneFlip(mesh.nFaces()-mesh.nInternalFaces(), false);
00555         forAll(*this, i)
00556         {
00557             label faceI = operator[](i);
00558 
00559             if (!mesh.isInternalFace(faceI))
00560             {
00561                 neiZoneFace[faceI-mesh.nInternalFaces()] = true;
00562                 neiZoneFlip[faceI-mesh.nInternalFaces()] = flipMap()[i];
00563             }
00564         }
00565         boolList myZoneFace(neiZoneFace);
00566         syncTools::swapBoundaryFaceList(mesh, neiZoneFace, false);
00567         boolList myZoneFlip(neiZoneFlip);
00568         syncTools::swapBoundaryFaceList(mesh, neiZoneFlip, false);
00569 
00570         forAll(*this, i)
00571         {
00572             label faceI = operator[](i);
00573 
00574             label patchI = bm.whichPatch(faceI);
00575 
00576             if (patchI != -1 && bm[patchI].coupled())
00577             {
00578                 label bFaceI = faceI-mesh.nInternalFaces();
00579 
00580                 // Check face in zone on both sides
00581                 if (myZoneFace[bFaceI] != neiZoneFace[bFaceI])
00582                 {
00583                     boundaryError = true;
00584 
00585                     if (report)
00586                     {
00587                         Pout<< " ***Problem with faceZone " << index()
00588                             << " named " << name()
00589                             << ". Face " << faceI
00590                             << " on coupled patch "
00591                             << bm[patchI].name()
00592                             << " is not consistent with its coupled neighbour."
00593                             << endl;
00594                     }
00595                 }
00596                 else if (myZoneFlip[bFaceI] == neiZoneFlip[bFaceI])
00597                 {
00598                     // Flip state should be opposite.
00599                     boundaryError = true;
00600 
00601                     if (report)
00602                     {
00603                         Pout<< " ***Problem with faceZone " << index()
00604                             << " named " << name()
00605                             << ". Face " << faceI
00606                             << " on coupled patch "
00607                             << bm[patchI].name()
00608                             << " does not have consistent flipMap"
00609                             << " across coupled faces."
00610                             << endl;
00611                     }
00612                 }
00613             }
00614         }
00615     }
00616 
00617     return returnReduce(boundaryError, orOp<bool>());
00618 }
00619 
00620 
00621 void Foam::faceZone::movePoints(const pointField& p)
00622 {
00623     if (patchPtr_)
00624     {
00625         patchPtr_->movePoints(p);
00626     }
00627 }
00628 
00629 void Foam::faceZone::write(Ostream& os) const
00630 {
00631     os  << nl << name()
00632         << nl << static_cast<const labelList&>(*this)
00633         << nl << flipMap();
00634 }
00635 
00636 
00637 void Foam::faceZone::writeDict(Ostream& os) const
00638 {
00639     os  << nl << name() << nl << token::BEGIN_BLOCK << nl
00640         << "    type " << type() << token::END_STATEMENT << nl;
00641 
00642     writeEntry("faceLabels", os);
00643     flipMap().writeEntry("flipMap", os);
00644 
00645     os  << token::END_BLOCK << endl;
00646 }
00647 
00648 
00649 // * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
00650 
00651 Foam::Ostream& Foam::operator<<(Ostream& os, const faceZone& p)
00652 {
00653     p.write(os);
00654     os.check("Ostream& operator<<(Ostream& f, const faceZone& p");
00655     return os;
00656 }
00657 
00658 
00659 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines