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