FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/SortableList.H>
28 #include <OpenFOAM/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 = labelListList(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 opposite vertex for 0th face
70  const Face& f = localFaces[faceNbs[0]];
71 
72  label fp0 = findIndex(f, e[0]);
73  label fp1 = f.fcIndex(fp0);
74  label vertI = (f[fp1] != e[1] ? f[fp1] : f.fcIndex(fp1));
75 
76  // Get vector normal both to e2 and to edge from opposite vertex
77  // to edge (will be x-axis of our coordinate system)
78  vector e0 = e2 ^ (localPoints[vertI] - edgePt);
79  e0 /= mag(e0) + VSMALL;
80 
81  // Get y-axis of coordinate system
82  vector e1 = e2 ^ e0;
83 
84  SortableList<scalar> faceAngles(faceNbs.size());
85 
86  // e0 is reference so angle is 0
87  faceAngles[0] = 0;
88 
89  for (label nbI = 1; nbI < faceNbs.size(); nbI++)
90  {
91  // Get opposite vertex
92  const Face& f = localFaces[faceNbs[nbI]];
93  label fp0 = findIndex(f, e[0]);
94  label fp1 = f.fcIndex(fp0);
95  label vertI = (f[fp1] != e[1] ? f[fp1] : f.fcIndex(fp1));
96 
97  vector vec = e2 ^ (localPoints[vertI] - edgePt);
98  vec /= mag(vec) + VSMALL;
99 
100  faceAngles[nbI] = pseudoAngle
101  (
102  e0,
103  e1,
104  vec
105  );
106  }
107 
108  faceAngles.sort();
109 
110  sortedEdgeFaces[edgeI] = UIndirectList<label>
111  (
112  faceNbs,
113  faceAngles.indices()
114  );
115  }
116  else
117  {
118  // No need to sort. Just copy.
119  sortedEdgeFaces[edgeI] = faceNbs;
120  }
121  }
122 
123  return sortedEdgeFaces;
124 }
125 
126 
127 // ************************ vim: set sw=4 sts=4 et: ************************ //