PatchToolsSortEdges.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-2021 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<class FaceList, class PointField>
34 (
36 )
37 {
38  typedef typename PrimitivePatch<FaceList, PointField>::FaceType FaceType;
39  typedef typename PrimitivePatch<FaceList, PointField>::PointType PointType;
40 
41  const edgeList& edges = p.edges();
42  const labelListList& edgeFaces = p.edgeFaces();
43  const List<FaceType>& localFaces = p.localFaces();
44  const Field<PointType>& localPoints = p.localPoints();
45 
46  // create the lists for the various results. (resized on completion)
47  labelListList sortedEdgeFaces(edgeFaces.size());
48 
49  forAll(edgeFaces, edgeI)
50  {
51  const labelList& faceNbs = edgeFaces[edgeI];
52 
53  if (faceNbs.size() > 2)
54  {
55  // Get point on edge and normalised direction of edge (= e2 base
56  // of our coordinate system)
57  const edge& e = edges[edgeI];
58 
59  const point& edgePt = localPoints[e.start()];
60 
61  vector e2 = e.vec(localPoints);
62  e2 /= mag(e2) + vSmall;
63 
64  // Get the vertex on 0th face that forms a vector with the first
65  // edge point that has the largest angle with the edge
66  const FaceType& f0 = localFaces[faceNbs[0]];
67 
68  scalar maxAngle = great;
69  vector maxAngleEdgeDir(vector::max);
70 
71  forAll(f0, fpI)
72  {
73  if (f0[fpI] != e.start())
74  {
75  vector faceEdgeDir = localPoints[f0[fpI]] - edgePt;
76  faceEdgeDir /= mag(faceEdgeDir) + vSmall;
77 
78  const scalar angle = e2 & faceEdgeDir;
79 
80  if (mag(angle) < maxAngle)
81  {
82  maxAngle = angle;
83  maxAngleEdgeDir = faceEdgeDir;
84  }
85  }
86  }
87 
88  // Get vector normal both to e2 and to edge from opposite vertex
89  // to edge (will be x-axis of our coordinate system)
90  vector e0 = e2 ^ maxAngleEdgeDir;
91  e0 /= mag(e0) + vSmall;
92 
93  // Get y-axis of coordinate system
94  vector e1 = e2 ^ e0;
95 
96  SortableList<scalar> faceAngles(faceNbs.size());
97 
98  // e0 is reference so angle is 0
99  faceAngles[0] = 0;
100 
101  for (label nbI = 1; nbI < faceNbs.size(); nbI++)
102  {
103  // Get the vertex on face that forms a vector with the first
104  // edge point that has the largest angle with the edge
105  const FaceType& f = localFaces[faceNbs[nbI]];
106 
107  maxAngle = great;
108  maxAngleEdgeDir = vector::max;
109 
110  forAll(f, fpI)
111  {
112  if (f[fpI] != e.start())
113  {
114  vector faceEdgeDir = localPoints[f[fpI]] - edgePt;
115  faceEdgeDir /= mag(faceEdgeDir) + vSmall;
116 
117  const scalar angle = e2 & faceEdgeDir;
118 
119  if (mag(angle) < maxAngle)
120  {
121  maxAngle = angle;
122  maxAngleEdgeDir = faceEdgeDir;
123  }
124  }
125  }
126 
127  vector vec = e2 ^ maxAngleEdgeDir;
128  vec /= mag(vec) + vSmall;
129 
130  faceAngles[nbI] = pseudoAngle
131  (
132  e0,
133  e1,
134  vec
135  );
136  }
137 
138  faceAngles.sort();
139 
140  sortedEdgeFaces[edgeI] = UIndirectList<label>
141  (
142  faceNbs,
143  faceAngles.indices()
144  );
145  }
146  else
147  {
148  // No need to sort. Just copy.
149  sortedEdgeFaces[edgeI] = faceNbs;
150  }
151  }
152 
153  return sortedEdgeFaces;
154 }
155 
156 
157 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void sort()
(stable) sort the list (if changed after construction time)
Definition: SortableList.C:112
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:306
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
A list that is sorted upon construction or when explicitly requested with the sort() method...
Definition: List.H:80
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
const Field< PointType > & localPoints() const
Return pointField of points in patch.
A list of faces which address into the list of points.
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:96
vector vec(const pointField &) const
Return the vector (end - start)
Definition: edgeI.H:175
labelList f(nPoints)
std::remove_reference< FaceList >::type::value_type FaceType
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
A List with indirect addressing.
Definition: fvMatrix.H:106
dimensioned< scalar > mag(const dimensioned< Type > &)
static labelListList sortedEdgeFaces(const PrimitivePatch< FaceList, PointField > &)
Return edge-face addressing sorted by angle around the edge.
std::remove_reference< PointField >::type::value_type PointType
const dimensionedScalar e
Elementary charge.
Definition: doubleScalar.H:105
label start() const
Return start vertex label.
Definition: edgeI.H:81