PatchToolsCheck.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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<class FaceList, class PointField>
32 (
34  const bool report,
35  labelHashSet* setPtr
36 )
37 {
38  typedef typename PrimitivePatch<FaceList, PointField>::FaceType FaceType;
39 
40  bool foundError = false;
41 
42  // Check edge normals, face normals, point normals.
43  forAll(p.faceEdges(), facei)
44  {
45  const labelList& edgeLabels = p.faceEdges()[facei];
46  bool valid = true;
47 
48  if (edgeLabels.size() < 3)
49  {
50  if (report)
51  {
52  Info<< "Face[" << facei << "] " << p[facei]
53  << " has fewer than 3 edges. Edges: " << edgeLabels
54  << endl;
55  }
56  valid = false;
57  }
58  else
59  {
60  forAll(edgeLabels, i)
61  {
62  if (edgeLabels[i] < 0 || edgeLabels[i] >= p.nEdges())
63  {
64  if (report)
65  {
66  Info<< "edge number " << edgeLabels[i]
67  << " on face " << facei
68  << " out-of-range\n"
69  << "This usually means the input surface has "
70  << "edges with more than 2 faces connected."
71  << endl;
72  }
73  valid = false;
74  }
75  }
76  }
77 
78  if (!valid)
79  {
80  foundError = true;
81  continue;
82  }
83 
84 
85  //
86  //- Compute normal from 3 points, use the first as the origin
87  // minor warpage should not be a problem
88  const FaceType& f = p[facei];
89  const point& p0 = p.points()[f[0]];
90  const point& p1 = p.points()[f[1]];
91  const point& p2 = p.points()[f.last()];
92 
93  const vector pointNormal((p1 - p0) ^ (p2 - p0));
94  if ((pointNormal & p.faceNormals()[facei]) < 0)
95  {
96  foundError = false;
97 
98  if (report)
99  {
100  Info
101  << "Normal calculated from points inconsistent"
102  << " with faceNormal" << nl
103  << "face: " << f << nl
104  << "points: " << p0 << ' ' << p1 << ' ' << p2 << nl
105  << "pointNormal:" << pointNormal << nl
106  << "faceNormal:" << p.faceNormals()[facei] << endl;
107  }
108  }
109  }
110 
111 
112  forAll(p.edges(), edgeI)
113  {
114  const edge& e = p.edges()[edgeI];
115  const labelList& neighbouringFaces = p.edgeFaces()[edgeI];
116 
117  if (neighbouringFaces.size() == 2)
118  {
119  // we use localFaces() since edges() are LOCAL
120  // these are both already available
121  const FaceType& faceA = p.localFaces()[neighbouringFaces[0]];
122  const FaceType& faceB = p.localFaces()[neighbouringFaces[1]];
123 
124  // If the faces are correctly oriented, the edges must go in
125  // different directions on connected faces.
126  if (faceA.edgeDirection(e) == faceB.edgeDirection(e))
127  {
128  if (report)
129  {
130  Info<< "face orientation incorrect." << nl
131  << "localEdge[" << edgeI << "] " << e
132  << " between faces:" << nl
133  << " face[" << neighbouringFaces[0] << "] "
134  << p[neighbouringFaces[0]]
135  << " localFace: " << faceA
136  << nl
137  << " face[" << neighbouringFaces[1] << "] "
138  << p[neighbouringFaces[1]]
139  << " localFace: " << faceB
140  << endl;
141  }
142  if (setPtr)
143  {
144  setPtr->insert(edgeI);
145  }
146 
147  foundError = true;
148  }
149  }
150  else if (neighbouringFaces.size() != 1)
151  {
152  if (report)
153  {
154  Info
155  << "Wrong number of edge neighbours." << nl
156  << "edge[" << edgeI << "] " << e
157  << " with points:" << p.localPoints()[e.start()]
158  << ' ' << p.localPoints()[e.end()]
159  << " has neighbouringFaces:" << neighbouringFaces << endl;
160  }
161 
162  if (setPtr)
163  {
164  setPtr->insert(edgeI);
165  }
166 
167  foundError = true;
168  }
169  }
170 
171  return foundError;
172 }
173 
174 
175 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:109
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static bool checkOrientation(const PrimitivePatch< FaceList, PointField > &, const bool report=false, labelHashSet *marked=0)
Check for orientation issues.
A list of faces which address into the list of points.
std::remove_reference< FaceList >::type::value_type FaceType
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:61
bool valid(const PtrList< ModelType > &l)
const doubleScalar e
Definition: doubleScalar.H:106
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
static const char nl
Definition: Ostream.H:266
labelList f(nPoints)
volScalarField & p