PatchToolsSortEdges.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2013 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 #include "SortableList.H"
28 #include "transform.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 template
33 <
34  class Face,
35  template<class> class FaceList,
36  class PointField,
37  class PointType
38 >
39 
42 (
44 )
45 {
46  const edgeList& edges = p.edges();
47  const labelListList& edgeFaces = p.edgeFaces();
48  const List<Face>& localFaces = p.localFaces();
49  const Field<PointType>& localPoints = p.localPoints();
50 
51  // create the lists for the various results. (resized on completion)
52  labelListList sortedEdgeFaces(edgeFaces.size());
53 
54  forAll(edgeFaces, edgeI)
55  {
56  const labelList& faceNbs = edgeFaces[edgeI];
57 
58  if (faceNbs.size() > 2)
59  {
60  // Get point on edge and normalized direction of edge (= e2 base
61  // of our coordinate system)
62  const edge& e = edges[edgeI];
63 
64  const point& edgePt = localPoints[e.start()];
65 
66  vector e2 = e.vec(localPoints);
67  e2 /= mag(e2) + VSMALL;
68 
69  // Get the vertex on 0th face that forms a vector with the first
70  // edge point that has the largest angle with the edge
71  const Face& f0 = localFaces[faceNbs[0]];
72 
73  scalar maxAngle = GREAT;
74  vector maxAngleEdgeDir(vector::max);
75 
76  forAll(f0, fpI)
77  {
78  if (f0[fpI] != e.start())
79  {
80  vector faceEdgeDir = localPoints[f0[fpI]] - edgePt;
81  faceEdgeDir /= mag(faceEdgeDir) + VSMALL;
82 
83  const scalar angle = e2 & faceEdgeDir;
84 
85  if (mag(angle) < maxAngle)
86  {
87  maxAngle = angle;
88  maxAngleEdgeDir = faceEdgeDir;
89  }
90  }
91  }
92 
93  // Get vector normal both to e2 and to edge from opposite vertex
94  // to edge (will be x-axis of our coordinate system)
95  vector e0 = e2 ^ maxAngleEdgeDir;
96  e0 /= mag(e0) + VSMALL;
97 
98  // Get y-axis of coordinate system
99  vector e1 = e2 ^ e0;
100 
101  SortableList<scalar> faceAngles(faceNbs.size());
102 
103  // e0 is reference so angle is 0
104  faceAngles[0] = 0;
105 
106  for (label nbI = 1; nbI < faceNbs.size(); nbI++)
107  {
108  // Get the vertex on face that forms a vector with the first
109  // edge point that has the largest angle with the edge
110  const Face& f = localFaces[faceNbs[nbI]];
111 
112  maxAngle = GREAT;
113  maxAngleEdgeDir = vector::max;
114 
115  forAll(f, fpI)
116  {
117  if (f[fpI] != e.start())
118  {
119  vector faceEdgeDir = localPoints[f[fpI]] - edgePt;
120  faceEdgeDir /= mag(faceEdgeDir) + VSMALL;
121 
122  const scalar angle = e2 & faceEdgeDir;
123 
124  if (mag(angle) < maxAngle)
125  {
126  maxAngle = angle;
127  maxAngleEdgeDir = faceEdgeDir;
128  }
129  }
130  }
131 
132  vector vec = e2 ^ maxAngleEdgeDir;
133  vec /= mag(vec) + VSMALL;
134 
135  faceAngles[nbI] = pseudoAngle
136  (
137  e0,
138  e1,
139  vec
140  );
141  }
142 
143  faceAngles.sort();
144 
145  sortedEdgeFaces[edgeI] = UIndirectList<label>
146  (
147  faceNbs,
148  faceAngles.indices()
149  );
150  }
151  else
152  {
153  // No need to sort. Just copy.
154  sortedEdgeFaces[edgeI] = faceNbs;
155  }
156  }
157 
158  return sortedEdgeFaces;
159 }
160 
161 
162 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
static labelListList sortedEdgeFaces(const PrimitivePatch< Face, FaceList, PointField, PointType > &)
Return edge-face addressing sorted by angle around the edge.
void sort()
(stable) sort the list (if changed after construction time)
Definition: SortableList.C:104
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
scalar pseudoAngle(const vector &e0, const vector &e1, const vector &vec)
Estimate angle of vec in coordinate system (e0, e1, e0^e1).
Definition: transform.H:289
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
A list that is sorted upon construction or when explicitly requested with the sort() method...
Definition: List.H:73
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
const List< Face > & localFaces() const
Return patch faces addressing into local point list.
A list of faces which address into the list of points.
const dimensionedScalar e
Elementary charge.
Definition: doubleFloat.H:78
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
3D tensor transformation operations.
const labelListList & edgeFaces() const
Return edge-face addressing.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
const labelList & indices() const
Return the list of sorted indices. Updated every sort.
Definition: SortableList.H:93
const Field< PointType > & localPoints() const
Return pointField of points in patch.
vector vec(const pointField &) const
Return the vector (end - start)
Definition: edgeI.H:157
labelList f(nPoints)
A List with indirect addressing.
Definition: fvMatrix.H:106
dimensioned< scalar > mag(const dimensioned< Type > &)
label start() const
Return start vertex label.
Definition: edgeI.H:81