src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchToolsCheck.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 \*---------------------------------------------------------------------------*/
  25
  26 #include "PatchTools.H"
  27
  28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  29
  30 template
  31 <
  32     class Face,
  33     template<class> class FaceList,
  34     class PointField,
  35     class PointType
  36 >
  37
  38 bool
  39 Foam::PatchTools::checkOrientation
  40 (
  41     const PrimitivePatch<Face, FaceList, PointField, PointType>& p,
  42     const bool report,
  43     labelHashSet* setPtr
  44 )
  45 {
  46     bool foundError = false;
  47
  48     // Check edge normals, face normals, point normals.
  49     forAll(p.faceEdges(), faceI)
  50     {
  51         const labelList& edgeLabels = p.faceEdges()[faceI];
  52         bool valid = true;
  53
  54         if (edgeLabels.size() < 3)
  55         {
  56             if (report)
  57             {
  58                 Info<< "Face[" << faceI << "] " << p[faceI]
  59                     << " has fewer than 3 edges. Edges: " << edgeLabels
  60                     << endl;
  61             }
  62             valid = false;
  63         }
  64         else
  65         {
  66             forAll(edgeLabels, i)
  67             {
  68                 if (edgeLabels[i] < 0 || edgeLabels[i] >= p.nEdges())
  69                 {
  70                     if (report)
  71                     {
  72                         Info<< "edge number " << edgeLabels[i]
  73                             << " on face " << faceI
  74                             << " out-of-range\n"
  75                             << "This usually means the input surface has "
  76                             << "edges with more than 2 faces connected."
  77                             << endl;
  78                     }
  79                     valid = false;
  80                 }
  81             }
  82         }
  83
  84         if (!valid)
  85         {
  86             foundError = true;
  87             continue;
  88         }
  89
  90
  91         //
  92         //- Compute normal from 3 points, use the first as the origin
  93         // minor warpage should not be a problem
  94         const Face& f = p[faceI];
  95         const point& p0 = p.points()[f[0]];
  96         const point& p1 = p.points()[f[1]];
  97         const point& p2 = p.points()[f.last()];
  98
  99         const vector pointNormal((p1 - p0) ^ (p2 - p0));
 100         if ((pointNormal & p.faceNormals()[faceI]) < 0)
 101         {
 102             foundError = false;
 103
 104             if (report)
 105             {
 106                 Info
 107                     << "Normal calculated from points inconsistent"
 108                     << " with faceNormal" << nl
 109                     << "face: " << f << nl
 110                     << "points: " << p0 << ' ' << p1 << ' ' << p2 << nl
 111                     << "pointNormal:" << pointNormal << nl
 112                     << "faceNormal:" << p.faceNormals()[faceI] << endl;
 113             }
 114         }
 115     }
 116
 117
 118     forAll(p.edges(), edgeI)
 119     {
 120         const edge& e = p.edges()[edgeI];
 121         const labelList& neighbouringFaces = p.edgeFaces()[edgeI];
 122
 123         if (neighbouringFaces.size() == 2)
 124         {
 125             // we use localFaces() since edges() are LOCAL
 126             // these are both already available
 127             const Face& faceA = p.localFaces()[neighbouringFaces[0]];
 128             const Face& faceB = p.localFaces()[neighbouringFaces[1]];
 129
 130             // If the faces are correctly oriented, the edges must go in
 131             // different directions on connected faces.
 132             if (faceA.edgeDirection(e) == faceB.edgeDirection(e))
 133             {
 134                 if (report)
 135                 {
 136                     Info<< "face orientation incorrect." << nl
 137                         << "localEdge[" << edgeI << "] " << e
 138                         << " between faces:" << nl
 139                         << "  face[" << neighbouringFaces[0] << "] "
 140                         << p[neighbouringFaces[0]]
 141                         << "   localFace: " << faceA
 142                         << nl
 143                         << "  face[" << neighbouringFaces[1] << "] "
 144                         << p[neighbouringFaces[1]]
 145                         << "   localFace: " << faceB
 146                         << endl;
 147                 }
 148                 if (setPtr)
 149                 {
 150                     setPtr->insert(edgeI);
 151                 }
 152
 153                 foundError = true;
 154             }
 155         }
 156         else if (neighbouringFaces.size() != 1)
 157         {
 158             if (report)
 159             {
 160                 Info
 161                     << "Wrong number of edge neighbours." << nl
 162                     << "edge[" << edgeI << "] " << e
 163                     << " with points:" << p.localPoints()[e.start()]
 164                     << ' ' << p.localPoints()[e.end()]
 165                     << " has neighbouringFaces:" << neighbouringFaces << endl;
 166             }
 167
 168             if (setPtr)
 169             {
 170                 setPtr->insert(edgeI);
 171             }
 172
 173             foundError = true;
 174         }
 175     }
 176
 177     return foundError;
 178 }
 179
 180
 181 // ************************************************************************* //