.. _program_listing_file_cellFinders_CrawlingCellFinder_CrawlingCellFinder.C:

Program Listing for File CrawlingCellFinder.C
=============================================

|exhale_lsh| :ref:`Return to documentation for file <file_cellFinders_CrawlingCellFinder_CrawlingCellFinder.C>` (``cellFinders/CrawlingCellFinder/CrawlingCellFinder.C``)

.. |exhale_lsh| unicode:: U+021B0 .. UPWARDS ARROW WITH TIP LEFTWARDS

.. code-block:: cpp

   /*---------------------------------------------------------------------------* \
   License
       This file is part of libWallModelledLES.
   
       libWallModelledLES is free software: you can redistribute it and/or modify
       it under the terms of the GNU General Public License as published by
       the Free Software Foundation, either version 3 of the License, or
       (at your option) any later version.
   
       libWallModelledLES is distributed in the hope that it will be useful, but
       WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
       or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
       for more details.
   
       You should have received a copy of the GNU General Public License
       along with libWallModelledLES. 
       If not, see <http://www.gnu.org/licenses/>.
   
   \*---------------------------------------------------------------------------*/
   
   #include "CrawlingCellFinder.H"
   #include "volFields.H"
   #include "surfaceFields.H"
   #include "codeRules.H"
   #include <cmath>
   
   
   // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
   
   #if !defined(DOXYGEN_SHOULD_SKIP_THIS)
   namespace Foam
   {
       defineTypeNameAndDebug(CrawlingCellFinder, 0);
       addToRunTimeSelectionTable(CellFinder, CrawlingCellFinder, Patch);
   }
   #endif
   
   
   // * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
   
   
   
   // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
   
   Foam::CrawlingCellFinder::CrawlingCellFinder
   (
       const fvPatch & p
   )
   :
       CellFinder(p)
   {
       if (debug)
       {
           Info << "CrawlingCellFinder: Constructing from patch" << nl;
       }
   
       createFields();
   }
   
   Foam::CrawlingCellFinder::CrawlingCellFinder
   (
       const word & cellFinderName,
       const fvPatch & p
   )
   :
      CellFinder(p)
   {
   }
   
   // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
   
   
   void Foam::CrawlingCellFinder::findCellIndices
   (
       labelList & indexList,
       const scalarField & h,
       const bool hIsIndex
   ) const
   {
       const labelList & owner = mesh_.faceOwner();
       const labelList & neighbour = mesh_.faceNeighbour();
       const vectorField & faceCentres = mesh().Cf().primitiveField();
       const List<cell> & cells = mesh().cells();
   
       const vectorField & patchFaceCentres = patch().Cf();
       const volVectorField & C = mesh_.C();
       const UList<label> & faceCells = patch().faceCells();
       const List<face> & faces = mesh().faces();
       const label patchStart = patch().start();
       const polyBoundaryMesh & boundaryMesh = mesh().boundaryMesh();
   
       forAll(patch(), patchFaceI)
       {
   
           label nLayers = 0;
   
           if (hIsIndex)
           {
               nLayers = std::round(h[patchFaceI]);
               if (nLayers < 1)
               {
                   Warning
                       << "CrawlingCellFinder: Could not round "
                       << h[patchFaceI] << " to a valid cell layer index. "
                       << "Will fall back to wall-adjacent cell for face "
                       <<  patchFaceI << " on patch " << patch().name() << nl;
                   indexList[patchFaceI] = faceCells[patchFaceI];
                   continue;
               }
               else if (nLayers == 1)
               {
                   indexList[patchFaceI] = faceCells[patchFaceI];
                   continue;
               }
           }
           else
           {
               if (h[patchFaceI] < 0)
               {
                   Warning
                       << "CrawlingCellFinder: " << h[patchFaceI]
                       << " is negative and thus not a valid distance. "
                       << "Will fall back to wall-adjacent cell for face "
                       <<  patchFaceI << " on patch " << patch().name() << nl;
                   indexList[patchFaceI] = faceCells[patchFaceI];
                   continue;
               }
               else if (h[patchFaceI] == 0)
               {
                   indexList[patchFaceI] = faceCells[patchFaceI];
                   continue;
               }
   
               nLayers = 100; // arbitrary big number
           }
   
   
           label startCellIndex = faceCells[patchFaceI];
           cell startCell = cells[startCellIndex];
           label startFaceLabel = patchStart + patchFaceI;
           label nextCellIndex = -1;
   
           for (label layer=0; layer < nLayers-1; layer++)
           {
               label opposingFace = 
                   startCell.opposingFaceLabel(startFaceLabel, faces);
   
               if (opposingFace == -1)
               {
                   Warning
                       << "CrawlingCellFinder: Could not find opposing face"
                       << "for cell " << startCellIndex << " with cell center "
                       << C[startCellIndex] << " and face " << startFaceLabel << nl
                       << "Will stop crawling and use the last valid cell " 
                       << "corresponding to index " << layer + 1
                       << nl;
                   
                   indexList[patchFaceI] = startCellIndex;
                   break;
               }
               else if (opposingFace >= mesh().nInternalFaces())
               {
                   label opposingPatchInd = boundaryMesh.whichPatch(opposingFace);
                   const fvPatch & opposingPatch =
                       mesh().boundary()[opposingPatchInd];
                   word opposingPatchName = opposingPatch.name();
   
                   // distance that will be used due to abortion in crawling
                   scalar distance =
                       mag(C[startCellIndex] - patchFaceCentres[patchFaceI]);
   
                   // If distance-based we might actually get the right h
                   // in the last cell, so need to check for that 
                   // before issuing the warning
                   if (hIsIndex || mag(distance - h[patchFaceI])/distance > SMALL)
                   {
                       Warning
                           << "CrawlingCellFinder: The opposing face for cell "
                           << startCellIndex << " with cell center "
                           << C[startCellIndex] << " and face " << startFaceLabel
                           << " belongs to patch " << opposingPatchName << nl 
                           << "Will stop crawling and use the last valid cell " 
                           << "corresponding to index " << layer + 1
                           << " and cell centre at distance " << distance
                           << ". Requested index/distance is " << h[patchFaceI]
                           << " which may or may not correspond to this cell." 
                           << nl;
                   }
   
                   break;
               }
   
               scalar distance =
                   mag(faceCentres[opposingFace] - patchFaceCentres[patchFaceI]);
   
               // if the opposing face is above h, stop and grab the cell
               if ((!hIsIndex) && (distance > h[patchFaceI]))
               {
                   break; 
               }
   
               if (owner[opposingFace] == startCellIndex)
               {
                   nextCellIndex = neighbour[opposingFace];
               }
               else
               {
                   nextCellIndex = owner[opposingFace];
               }
   
               startFaceLabel = opposingFace;
               startCellIndex = nextCellIndex;
               startCell = cells[startCellIndex];
           }
   
           indexList[patchFaceI] = startCellIndex;
       }
   }
   
   
   void Foam::CrawlingCellFinder::findCellIndices
   (
       labelListList & indexList,
       const scalarField & h,
       const bool hIsIndex,
       const bool excludeWallAdjacent
   ) const
   {
       const labelList & owner = mesh_.faceOwner();
       const labelList & neighbour = mesh_.faceNeighbour();
       const vectorField & faceCentres = mesh().Cf().primitiveField();
       const List<cell> & cells = mesh().cells();
   
       const vectorField & patchFaceCentres = patch().Cf();
       const volVectorField & C = mesh_.C();
       const UList<label> & faceCells = patch().faceCells();
       const List<face> & faces = mesh().faces();
       const label patchStart = patch().start();
       const polyBoundaryMesh & boundaryMesh = mesh().boundaryMesh();
   
       forAll(patch(), patchFaceI)
       {
   
           label nLayers = 0;
   
           if (hIsIndex)
           {
               nLayers = std::round(h[patchFaceI]);
               if (nLayers < 1)
               {
                   Warning
                       << "CrawlingCellFinder: Could not round "
                       << h[patchFaceI] << " to a valid cell layer index. "
                       << "Will fall back to wall-adjacent cell for face "
                       <<  patchFaceI << " on patch " << patch().name() << nl;
                   indexList[patchFaceI].setSize(1, faceCells[patchFaceI]);
                   continue;
               }
               else if (nLayers == 1)
               {
                   indexList[patchFaceI].setSize(1, faceCells[patchFaceI]);
                   continue;
               }
           }
           else
           {
               if (h[patchFaceI] < 0)
               {
                   Warning
                       << "CrawlingCellFinder: " << h[patchFaceI]
                       << " is negative and thus not a valid distance. "
                       << "Will fall back to wall-adjacent cell for face "
                       <<  patchFaceI << " on patch " << patch().name() << nl;
                   indexList[patchFaceI].setSize(1, faceCells[patchFaceI]);
                   continue;
               }
               else if (h[patchFaceI] == 0)
               {
                   indexList[patchFaceI].setSize(1, faceCells[patchFaceI]);
                   continue;
               }
   
               nLayers = 100; // arbitrary big number
           }
           indexList[patchFaceI].setSize(nLayers, -1);
   
           label startCellIndex = faceCells[patchFaceI];
           cell startCell = cells[startCellIndex];
           label startFaceLabel = patchStart + patchFaceI;
           label nextCellIndex = -1;
   
           label layerCounter = 0;
   
           for (label layer=0; layer < nLayers; layer++)
           {
               layerCounter++;
   
               // Grab the next cell
               indexList[patchFaceI][layer] = startCellIndex;
   
               label opposingFace = 
                   startCell.opposingFaceLabel(startFaceLabel, faces);
   
               if (opposingFace == -1)
               {
                   Warning
                       << "CrawlingCellFinder: Could not find opposing face"
                       << "for cell " << startCellIndex << " with cell center "
                       << C[startCellIndex] << " and face " << startFaceLabel << nl
                       << "Will stop crawling and use the last valid cell " 
                       << "corresponding to index " << layer + 1
                       << nl;
                   indexList[patchFaceI].setSize(layerCounter);
                   break;
               }
               // if we hit a patch
               else if (opposingFace >= mesh().nInternalFaces())
               {
                   label opposingPatchInd = boundaryMesh.whichPatch(opposingFace);
                   const fvPatch & opposingPatch =
                       mesh().boundary()[opposingPatchInd];
                   word opposingPatchName = opposingPatch.name();
   
                   // Recompute as the actual distance that will be used
                   // due to abortion in crawling
                   scalar distance =
                       mag(C[startCellIndex] - patchFaceCentres[patchFaceI]);
   
                   Warning
                       << "CrawlingCellFinder: The opposing face for cell "
                       << startCellIndex << " with cell center "
                       << C[startCellIndex] << " and face " << startFaceLabel
                       << " belongs to patch " << opposingPatchName << nl 
                       << "Will stop crawling and use the last valid cell " 
                       << "corresponding to index " << layer + 1
                       << " and distance " << distance
                       << nl;
                       
                   // Need this when hIsIndex
                   indexList[patchFaceI].setSize(layerCounter);
                   break;
               }
   
               scalar distance =
                   mag(faceCentres[opposingFace] - patchFaceCentres[patchFaceI]);
   
               // if the opposing face is above h, stop 
               if ((!hIsIndex) && (distance > h[patchFaceI]))
               {
                   indexList[patchFaceI].setSize(layerCounter);
                   break;
               }
   
               // find index of the next cell
               if (owner[opposingFace] == startCellIndex)
               {
                   nextCellIndex = neighbour[opposingFace];
               }
               else
               {
                   nextCellIndex = owner[opposingFace];
               }
   
               startFaceLabel = opposingFace;
               startCellIndex = nextCellIndex;
               startCell = cells[startCellIndex];
           }
   
           if ((indexList[patchFaceI].size() != 1) && excludeWallAdjacent)
           {
               for(int j=0; j<indexList[patchFaceI].size()-1; j++)
               {
                   indexList[patchFaceI][j] = indexList[patchFaceI][j+1];
               }
               indexList[patchFaceI].setSize(indexList[patchFaceI].size()-1);
   
           }
       }
   }
   
   // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //