surfaceIntersectionFuncs.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-2018 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 "surfaceIntersection.H"
27 #include "triSurfaceSearch.H"
28 #include "labelPairLookup.H"
29 #include "OFstream.H"
30 #include "HashSet.H"
31 #include "triSurface.H"
32 #include "pointIndexHit.H"
33 #include "meshTools.H"
34 
35 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
36 
37 void Foam::surfaceIntersection::writeOBJ(const point& pt, Ostream& os)
38 {
39  os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
40 }
41 
42 
43 void Foam::surfaceIntersection::writeOBJ
44 (
45  const List<point>& pts,
46  const List<edge>& edges,
47  Ostream& os
48 )
49 {
50  forAll(pts, i)
51  {
52  writeOBJ(pts[i], os);
53  }
54  forAll(edges, i)
55  {
56  const edge& e = edges[i];
57 
58  os << "l " << e.start()+1 << ' ' << e.end()+1 << endl;
59  }
60 }
61 
62 
63 // Get minimum length of all edges connected to point
64 Foam::scalar Foam::surfaceIntersection::minEdgeLen
65 (
66  const triSurface& surf,
67  const label pointi
68 )
69 {
70  const labelList& pEdges = surf.pointEdges()[pointi];
71 
72  scalar minLen = great;
73 
74  forAll(pEdges, pEdgeI)
75  {
76  const edge& e = surf.edges()[pEdges[pEdgeI]];
77 
78  minLen = min(minLen, e.mag(surf.localPoints()));
79  }
80 
81  return minLen;
82 }
83 
84 
85 // Get edge between fp and fp+1 on facei.
86 Foam::label Foam::surfaceIntersection::getEdge
87 (
88  const triSurface& surf,
89  const label facei,
90  const label fp
91 )
92 {
93  const edge faceEdge = surf.localFaces()[facei].faceEdge(fp);
94 
95  const labelList& eLabels = surf.faceEdges()[facei];
96 
97  forAll(eLabels, eI)
98  {
99  const label edgeI = eLabels[eI];
100 
101  if (surf.edges()[edgeI] == faceEdge)
102  {
103  return edgeI;
104  }
105  }
106 
108  << "Problem:: Cannot find edge with vertices " << faceEdge
109  << " in face " << facei
110  << abort(FatalError);
111 
112  return -1;
113 }
114 
115 
116 // Given a map remove all consecutive duplicate elements.
117 void Foam::surfaceIntersection::removeDuplicates
118 (
119  const labelList& map,
120  labelList& elems
121 )
122 {
123  bool hasDuplicate = false;
124 
125  label prevVertI = -1;
126 
127  forAll(elems, elemI)
128  {
129  label newVertI = map[elems[elemI]];
130 
131  if (newVertI == prevVertI)
132  {
133  hasDuplicate = true;
134 
135  break;
136  }
137  prevVertI = newVertI;
138  }
139 
140  if (hasDuplicate)
141  {
142  // Create copy
143  labelList oldElems(elems);
144 
145  label elemI = 0;
146 
147  // Insert first
148  elems[elemI++] = map[oldElems[0]];
149 
150  for (label vertI = 1; vertI < oldElems.size(); vertI++)
151  {
152  // Insert others only if they differ from one before
153  label newVertI = map[oldElems[vertI]];
154 
155  if (newVertI != elems.last())
156  {
157  elems[elemI++] = newVertI;
158  }
159  }
160  elems.setSize(elemI);
161  }
162 }
163 
164 
165 // Remap.
166 void Foam::surfaceIntersection::inlineRemap
167 (
168  const labelList& map,
169  labelList& elems
170 )
171 {
172  forAll(elems, elemI)
173  {
174  elems[elemI] = map[elems[elemI]];
175  }
176 }
177 
178 
179 // Remove all duplicate and degenerate elements. Return unique elements and
180 // map from old to new.
181 Foam::edgeList Foam::surfaceIntersection::filterEdges
182 (
183  const edgeList& edges,
184  labelList& map
185 )
186 {
187  HashSet<edge, Hash<edge>> uniqueEdges(10*edges.size());
188 
189  edgeList newEdges(edges.size());
190 
191  map.setSize(edges.size());
192  map = -1;
193 
194  label newEdgeI = 0;
195 
196  forAll(edges, edgeI)
197  {
198  const edge& e = edges[edgeI];
199 
200  if
201  (
202  (e.start() != e.end())
203  && (uniqueEdges.find(e) == uniqueEdges.end())
204  )
205  {
206  // Edge is -non degenerate and -not yet seen.
207  uniqueEdges.insert(e);
208 
209  map[edgeI] = newEdgeI;
210 
211  newEdges[newEdgeI++] = e;
212  }
213  }
214 
215  newEdges.setSize(newEdgeI);
216 
217  return newEdges;
218 }
219 
220 
221 // Remove all duplicate elements.
222 Foam::labelList Foam::surfaceIntersection::filterLabels
223 (
224  const labelList& elems,
225  labelList& map
226 )
227 {
228  labelHashSet uniqueElems(10*elems.size());
229 
230  labelList newElems(elems.size());
231 
232  map.setSize(elems.size());
233  map = -1;
234 
235  label newElemI = 0;
236 
237  forAll(elems, elemI)
238  {
239  label elem = elems[elemI];
240 
241  if (uniqueElems.find(elem) == uniqueElems.end())
242  {
243  // First time elem is seen
244  uniqueElems.insert(elem);
245 
246  map[elemI] = newElemI;
247 
248  newElems[newElemI++] = elem;
249  }
250  }
251 
252  newElems.setSize(newElemI);
253 
254  return newElems;
255 }
256 
257 
258 void Foam::surfaceIntersection::writeIntersectedEdges
259 (
260  const triSurface& surf,
261  const labelListList& edgeCutVerts,
262  Ostream& os
263 ) const
264 {
265  // Dump all points (surface followed by cutPoints)
266  const pointField& pts = surf.localPoints();
267 
268  forAll(pts, pointi)
269  {
270  writeOBJ(pts[pointi], os);
271  }
272  forAll(cutPoints(), cutPointi)
273  {
274  writeOBJ(cutPoints()[cutPointi], os);
275  }
276 
277  forAll(edgeCutVerts, edgeI)
278  {
279  const labelList& extraVerts = edgeCutVerts[edgeI];
280 
281  if (extraVerts.size())
282  {
283  const edge& e = surf.edges()[edgeI];
284 
285  // Start of original edge to first extra point
286  os << "l " << e.start()+1 << ' '
287  << extraVerts[0] + surf.nPoints() + 1 << endl;
288 
289  for (label i = 1; i < extraVerts.size(); i++)
290  {
291  os << "l " << extraVerts[i-1] + surf.nPoints() + 1 << ' '
292  << extraVerts[i] + surf.nPoints() + 1 << endl;
293  }
294 
295  os << "l " << extraVerts.last() + surf.nPoints() + 1
296  << ' ' << e.end()+1 << endl;
297  }
298  }
299 }
300 
301 
302 // Return 0 (p close to start), 1(close to end) or -1.
303 Foam::label Foam::surfaceIntersection::classify
304 (
305  const scalar startTol,
306  const scalar endTol,
307  const point& p,
308  const edge& e,
309  const pointField& points
310 )
311 {
312  if (mag(p - points[e.start()]) < startTol)
313  {
314  return 0;
315  }
316  else if (mag(p - points[e.end()]) < endTol)
317  {
318  return 1;
319  }
320  else
321  {
322  return -1;
323  }
324 }
325 
326 
327 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
const pointField & points
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
const doubleScalar e
Definition: doubleScalar.H:106
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
errorManip< error > abort(error &err)
Definition: errorManip.H:131
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
vector point
Point is a vector.
Definition: point.H:41
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
dimensioned< scalar > mag(const dimensioned< Type > &)
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:213
error FatalError
List< edge > edgeList
Definition: edgeList.H:38
volScalarField & p