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

cellDistFuncs.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 \*---------------------------------------------------------------------------*/
00025 
00026 #include "cellDistFuncs.H"
00027 #include <OpenFOAM/polyMesh.H>
00028 #include <OpenFOAM/wallPolyPatch.H>
00029 #include <OpenFOAM/polyBoundaryMesh.H>
00030 
00031 namespace Foam
00032 {
00033 
00034 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00035 
00036 defineTypeNameAndDebug(cellDistFuncs, 0);
00037 
00038 }
00039 
00040 
00041 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00042 
00043 // Find val in first nElems elements of elems.
00044 Foam::label Foam::cellDistFuncs::findIndex
00045 (
00046     const label nElems,
00047     const labelList& elems,
00048     const label val
00049 )
00050 {
00051     for (label i = 0; i < nElems; i++)
00052     {
00053         if (elems[i] == val)
00054         {
00055             return i;
00056         }
00057     }
00058     return -1;
00059 }
00060 
00061 
00062 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00063 
00064 // Construct from mesh
00065 Foam::cellDistFuncs::cellDistFuncs(const polyMesh& mesh)
00066 :
00067     mesh_(mesh)
00068 {}
00069 
00070 
00071 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00072 
00073 // Get patch ids of named patches
00074 Foam::labelHashSet Foam::cellDistFuncs::getPatchIDs
00075 (
00076     const wordList& patchNames
00077 ) const
00078 {
00079     const polyBoundaryMesh& bMesh = mesh().boundaryMesh();
00080 
00081     // Construct Set of patchNames for easy checking if included
00082     HashSet<word> patchNamesHash(patchNames.size());
00083 
00084     forAll(patchNames, patchI)
00085     {
00086         patchNamesHash.insert(patchNames[patchI]);
00087     }
00088 
00089     // Loop over all patches and check if patch name in hashset
00090 
00091     labelHashSet patchIDs(bMesh.size());
00092 
00093     forAll(bMesh, patchI)
00094     {
00095         const polyPatch& patch = bMesh[patchI];
00096 
00097         if (patchNamesHash.found(patch.name()))
00098         {
00099             patchIDs.insert(patchI);
00100         }
00101     }
00102     return patchIDs;
00103 }
00104 
00105 
00106 // Return smallest true distance from p to any of wallFaces.
00107 // Note that even if normal hits face we still check other faces.
00108 // Note that wallFaces is untruncated and we explicitly pass in size.
00109 Foam::scalar Foam::cellDistFuncs::smallestDist
00110 (
00111     const point& p,
00112     const polyPatch& patch,
00113     const label nWallFaces,
00114     const labelList& wallFaces,
00115     label& minFaceI
00116 ) const
00117 {
00118     const pointField& points = patch.points();
00119 
00120     scalar minDist = GREAT;
00121     minFaceI = -1;
00122 
00123     for(label wallFaceI = 0; wallFaceI < nWallFaces; wallFaceI++)
00124     {
00125         label patchFaceI = wallFaces[wallFaceI];
00126 
00127         pointHit curHit = patch[patchFaceI].nearestPoint(p, points);
00128 
00129         if (curHit.distance() < minDist)
00130         {
00131             minDist = curHit.distance();
00132             minFaceI = patch.start() + patchFaceI;
00133         }
00134     }
00135 
00136     return minDist;
00137 }
00138 
00139 
00140 // Get point neighbours of faceI (including faceI). Returns number of faces.
00141 // Note: does not allocate storage but does use linear search to determine
00142 // uniqueness. For polygonal faces this might be quite inefficient.
00143 Foam::label Foam::cellDistFuncs::getPointNeighbours
00144 (
00145     const primitivePatch& patch,
00146     const label patchFaceI,
00147     labelList& neighbours
00148 ) const
00149 {
00150     label nNeighbours = 0;
00151 
00152     // Add myself
00153     neighbours[nNeighbours++] = patchFaceI;
00154 
00155     // Add all face neighbours
00156     const labelList& faceNeighbours = patch.faceFaces()[patchFaceI];
00157 
00158     forAll(faceNeighbours, faceNeighbourI)
00159     {
00160         neighbours[nNeighbours++] = faceNeighbours[faceNeighbourI];
00161     }
00162 
00163     // Remember part of neighbours that contains edge-connected faces.
00164     label nEdgeNbs = nNeighbours;
00165 
00166 
00167     // Add all point-only neighbours by linear searching in edge neighbours.
00168     // Assumes that point-only neighbours are not using multiple points on
00169     // face.
00170 
00171     const face& f = patch.localFaces()[patchFaceI];
00172 
00173     forAll(f, fp)
00174     {
00175         label pointI = f[fp];
00176 
00177         const labelList& pointNbs = patch.pointFaces()[pointI];
00178 
00179         forAll(pointNbs, nbI)
00180         {
00181             label faceI = pointNbs[nbI];
00182 
00183             // Check for faceI in edge-neighbours part of neighbours
00184             if (findIndex(nEdgeNbs, neighbours, faceI) == -1)
00185             {
00186                 neighbours[nNeighbours++] = faceI;
00187             }
00188         }
00189     }
00190 
00191 
00192     if (debug)
00193     {
00194         // Check for duplicates
00195 
00196         // Use hashSet to determine nbs.
00197         labelHashSet nbs(4*f.size());
00198 
00199         forAll(f, fp)
00200         {
00201             const labelList& pointNbs = patch.pointFaces()[f[fp]];
00202 
00203             forAll(pointNbs, i)
00204             {
00205                 nbs.insert(pointNbs[i]);
00206             }
00207         }
00208 
00209         // Subtract ours.
00210         for (label i = 0; i < nNeighbours; i++)
00211         {
00212             label nb = neighbours[i];
00213 
00214             if (!nbs.found(nb))
00215             {
00216                 SeriousErrorIn("Foam::cellDistFuncs::getPointNeighbours")
00217                     << "getPointNeighbours : patchFaceI:" << patchFaceI
00218                     << " verts:" << f << endl;
00219 
00220                 forAll(f, fp)
00221                 {
00222                     SeriousErrorIn("Foam::cellDistFuncs::getPointNeighbours")
00223                         << "point:" << f[fp] << " pointFaces:"
00224                         << patch.pointFaces()[f[fp]] << endl;
00225                 }
00226 
00227                 for (label i = 0; i < nNeighbours; i++)
00228                 {
00229                     SeriousErrorIn("Foam::cellDistFuncs::getPointNeighbours")
00230                         << "fast nbr:" << neighbours[i]
00231                         << endl;
00232                 }
00233 
00234                 FatalErrorIn("getPointNeighbours")
00235                     << "Problem: fast pointNeighbours routine included " << nb
00236                     << " which is not in proper neigbour list " << nbs.toc()
00237                     << abort(FatalError);
00238             }
00239             nbs.erase(nb);
00240         }
00241 
00242         if (nbs.size())
00243         {
00244             FatalErrorIn("getPointNeighbours")
00245                 << "Problem: fast pointNeighbours routine did not find "
00246                 << nbs.toc() << abort(FatalError);
00247         }
00248     }
00249 
00250     return nNeighbours;
00251 }
00252 
00253 
00254 // size of largest patch (out of supplied subset of patches)
00255 Foam::label Foam::cellDistFuncs::maxPatchSize(const labelHashSet& patchIDs)
00256  const
00257 {
00258     label maxSize = 0;
00259 
00260     forAll(mesh().boundaryMesh(), patchI)
00261     {
00262         if (patchIDs.found(patchI))
00263         {
00264             const polyPatch& patch = mesh().boundaryMesh()[patchI];
00265 
00266             maxSize = Foam::max(maxSize, patch.size());
00267         }
00268     }
00269     return maxSize;
00270 }
00271 
00272 
00273 // sum of patch sizes (out of supplied subset of patches)
00274 Foam::label Foam::cellDistFuncs::sumPatchSize(const labelHashSet& patchIDs)
00275  const
00276 {
00277     label sum = 0;
00278 
00279     forAll(mesh().boundaryMesh(), patchI)
00280     {
00281         if (patchIDs.found(patchI))
00282         {
00283             const polyPatch& patch = mesh().boundaryMesh()[patchI];
00284 
00285             sum += patch.size();
00286         }
00287     }
00288     return sum;
00289 }
00290 
00291 
00292 // Gets nearest wall for cells next to wall
00293 void Foam::cellDistFuncs::correctBoundaryFaceCells
00294 (
00295     const labelHashSet& patchIDs,
00296     scalarField& wallDistCorrected,
00297     Map<label>& nearestFace
00298 ) const
00299 {
00300     // Size neighbours array for maximum possible (= size of largest patch)
00301     label maxPointNeighbours = maxPatchSize(patchIDs);
00302 
00303     labelList neighbours(maxPointNeighbours);
00304 
00305 
00306     // Correct all cells with face on wall
00307     const vectorField& cellCentres = mesh().cellCentres();
00308     const labelList& faceOwner = mesh().faceOwner();
00309 
00310     forAll(mesh().boundaryMesh(), patchI)
00311     {
00312         if (patchIDs.found(patchI))
00313         {
00314             const polyPatch& patch = mesh().boundaryMesh()[patchI];
00315 
00316             // Check cells with face on wall
00317             forAll(patch, patchFaceI)
00318             {
00319                 label nNeighbours = getPointNeighbours
00320                 (
00321                     patch,
00322                     patchFaceI,
00323                     neighbours
00324                 );
00325 
00326                 label cellI = faceOwner[patch.start() + patchFaceI];
00327 
00328                 label minFaceI = -1;
00329 
00330                 wallDistCorrected[cellI] = smallestDist
00331                 (
00332                     cellCentres[cellI],
00333                     patch,
00334                     nNeighbours,
00335                     neighbours,
00336                     minFaceI
00337                 );
00338 
00339                 // Store wallCell and its nearest neighbour
00340                 nearestFace.insert(cellI, minFaceI);
00341             }
00342         }
00343     }
00344 
00345 }
00346 
00347 
00348 // Correct all cells connected to wall (via point) and not in nearestFace
00349 void Foam::cellDistFuncs::correctBoundaryPointCells
00350 (
00351     const labelHashSet& patchIDs,
00352     scalarField& wallDistCorrected,
00353     Map<label>& nearestFace
00354 ) const
00355 {
00356     // Correct all (non-visited) cells with point on wall
00357 
00358     const labelListList& pointCells = mesh().pointCells();
00359     const vectorField& cellCentres = mesh().cellCentres();
00360 
00361     forAll(mesh().boundaryMesh(), patchI)
00362     {
00363         if (patchIDs.found(patchI))
00364         {
00365             const polyPatch& patch = mesh().boundaryMesh()[patchI];
00366 
00367             const labelList& meshPoints = patch.meshPoints();
00368             const labelListList& pointFaces = patch.pointFaces();
00369 
00370             forAll(meshPoints, meshPointI)
00371             {
00372                 label vertI = meshPoints[meshPointI];
00373 
00374                 const labelList& neighbours = pointCells[vertI];
00375 
00376                 forAll(neighbours, neighbourI)
00377                 {
00378                     label cellI = neighbours[neighbourI];
00379 
00380                     if (!nearestFace.found(cellI))
00381                     {
00382                         const labelList& wallFaces = pointFaces[meshPointI];
00383 
00384                         label minFaceI = -1;
00385 
00386                         wallDistCorrected[cellI] = smallestDist
00387                         (
00388                             cellCentres[cellI],
00389                             patch,
00390                             wallFaces.size(),
00391                             wallFaces,
00392                             minFaceI
00393                         );
00394 
00395                         // Store wallCell and its nearest neighbour
00396                         nearestFace.insert(cellI, minFaceI);
00397                     }
00398                 }
00399             }
00400         }
00401     }
00402 }
00403 
00404 
00405 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines