Loading...
Searching...
No Matches
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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2013 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "PatchTools.H"
31#include "transform.H"
32
33// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34
35template<class FaceList, class PointField>
38(
40)
41{
42 const edgeList& edges = p.edges();
43 const labelListList& edgeFaces = p.edgeFaces();
44 const auto& localFaces = p.localFaces();
45 const auto& localPoints = p.localPoints();
46
47 // create the lists for the various results. (resized on completion)
49
50 forAll(edgeFaces, edgeI)
51 {
52 const labelList& faceNbs = edgeFaces[edgeI];
53
54 if (faceNbs.size() > 2)
55 {
56 // Get point on edge and normalized direction of edge (= e2 base
57 // of our coordinate system)
58 const edge& e = edges[edgeI];
59
60 const point& edgePt = localPoints[e.start()];
61
62 const vector e2 = e.unitVec(localPoints);
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 auto& 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 const vector faceEdgeDir =
77 (
78 localPoints[f0[fpI]] - edgePt
79 );
80
81 const scalar angle = e2 & faceEdgeDir;
82
83 if (mag(angle) < maxAngle)
84 {
85 maxAngle = angle;
86 maxAngleEdgeDir = faceEdgeDir;
87 }
88 }
89 }
90
91 // Get vector normal both to e2 and to edge from opposite vertex
92 // to edge (will be x-axis of our coordinate system)
93 const vector e0 = normalised(e2 ^ maxAngleEdgeDir);
94
95 // Get y-axis of coordinate system
96 const vector e1 = e2 ^ e0;
97
98 SortableList<scalar> faceAngles(faceNbs.size());
99
100 // e0 is reference so angle is 0
101 faceAngles[0] = 0;
102
103 for (label nbI = 1; nbI < faceNbs.size(); nbI++)
104 {
105 // Get the vertex on face that forms a vector with the first
106 // edge point that has the largest angle with the edge
107 const auto& f = localFaces[faceNbs[nbI]];
108
109 maxAngle = GREAT;
110 maxAngleEdgeDir = vector::max;
111
112 forAll(f, fpI)
113 {
114 if (f[fpI] != e.start())
115 {
116 const vector faceEdgeDir =
118 (
119 localPoints[f[fpI]] - edgePt
120 );
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 const vector vec = normalised(e2 ^ maxAngleEdgeDir);
133
134 faceAngles[nbI] = pseudoAngle
135 (
136 e0,
137 e1,
138 vec
139 );
140 }
141
142 faceAngles.sort();
143
145 (
146 faceNbs,
147 faceAngles.indices()
148 );
149 }
150 else
151 {
152 // No need to sort. Just copy.
153 sortedEdgeFaces[edgeI] = faceNbs;
154 }
155 }
156
157 return sortedEdgeFaces;
158}
159
160
161// ************************************************************************* //
static labelListList sortedEdgeFaces(const PrimitivePatch< FaceList, PointField > &)
Return edge-face addressing sorted by angle around the edge.
A list of faces which address into the list of points.
A list that is sorted upon construction or when explicitly requested with the sort() method.
const labelList & indices() const noexcept
Return the list of sorted indices. Updated every sort.
void sort()
Forward (stable) sort the list (if changed after construction).
void size(const label n)
Older name for setAddressableSize.
Definition UList.H:118
static const Form max
An edge is a list of two vertex labels. This can correspond to a directed graph edge or an edge on a ...
Definition edge.H:62
volScalarField & p
List< edge > edgeList
List of edge.
Definition edgeList.H:32
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:469
List< labelList > labelListList
List of labelList.
Definition labelList.H:38
List< label > labelList
A List of labels.
Definition List.H:62
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
UIndirectList< label > labelUIndList
UIndirectList of labels.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
vector point
Point is a vector.
Definition point.H:37
Vector< scalar > vector
Definition vector.H:57
labelList f(nPoints)
volScalarField & e
#define forAll(list, i)
Loop across all elements in list.
Definition stdFoam.H:299
3D tensor transformation operations.