FreeFOAM The Cross-Platform CFD Toolkit
faceOnlySet.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 "faceOnlySet.H"
27 #include <meshTools/meshSearch.H>
28 #include <OpenFOAM/DynamicList.H>
29 #include <OpenFOAM/polyMesh.H>
30 
31 #include <lagrangian/Cloud.H>
33 #include <OpenFOAM/IDLList.H>
34 
36 
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41  defineTypeNameAndDebug(faceOnlySet, 0);
42  addToRunTimeSelectionTable(sampledSet, faceOnlySet, word);
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 
49 // Sample singly connected segment. Returns false if end_ reached.
50 bool Foam::faceOnlySet::trackToBoundary
51 (
52  Particle<passiveParticle>& singleParticle,
53  DynamicList<point>& samplingPts,
54  DynamicList<label>& samplingCells,
55  DynamicList<label>& samplingFaces,
56  DynamicList<scalar>& samplingCurveDist
57 ) const
58 {
59  // distance vector between sampling points
60  const vector offset = end_ - start_;
61  const vector smallVec = tol*offset;
62  const scalar smallDist = mag(smallVec);
63 
64  // Alias
65  const point& trackPt = singleParticle.position();
66 
67  while(true)
68  {
69  point oldPoint = trackPt;
70 
71  singleParticle.trackToFace(end_);
72 
73  if (singleParticle.face() != -1 && mag(oldPoint - trackPt) > smallDist)
74  {
75  // Reached face. Sample.
76  samplingPts.append(trackPt);
77  samplingCells.append(singleParticle.cell());
78  samplingFaces.append(singleParticle.face());
79  samplingCurveDist.append(mag(trackPt - start_));
80  }
81 
82  if (mag(trackPt - end_) < smallDist)
83  {
84  // end reached
85  return false;
86  }
87  else if (singleParticle.onBoundary())
88  {
89  // Boundary reached.
90  return true;
91  }
92  }
93 }
94 
95 
96 void Foam::faceOnlySet::calcSamples
97 (
98  DynamicList<point>& samplingPts,
99  DynamicList<label>& samplingCells,
100  DynamicList<label>& samplingFaces,
101  DynamicList<label>& samplingSegments,
102  DynamicList<scalar>& samplingCurveDist
103 ) const
104 {
105  // distance vector between sampling points
106  if (mag(end_ - start_) < SMALL)
107  {
108  FatalErrorIn("faceOnlySet::calcSamples()")
109  << "Incorrect sample specification :"
110  << " start equals end point." << endl
111  << " start:" << start_
112  << " end:" << end_
113  << exit(FatalError);
114  }
115 
116  const vector offset = (end_ - start_);
117  const vector normOffset = offset/mag(offset);
118  const vector smallVec = tol*offset;
119  const scalar smallDist = mag(smallVec);
120 
121 
122  // Get all boundary intersections
123  List<pointIndexHit> bHits = searchEngine().intersections
124  (
125  start_ - smallVec,
126  end_ + smallVec
127  );
128 
129  point bPoint(GREAT, GREAT, GREAT);
130  label bFaceI = -1;
131 
132  if (bHits.size())
133  {
134  bPoint = bHits[0].hitPoint();
135  bFaceI = bHits[0].index();
136  }
137 
138  // Get first tracking point. Use bPoint, bFaceI if provided.
139 
140  point trackPt;
141  label trackCellI = -1;
142  label trackFaceI = -1;
143 
144  //Info<< "before getTrackingPoint : bPoint:" << bPoint
145  // << " bFaceI:" << bFaceI << endl;
146 
147  getTrackingPoint
148  (
149  offset,
150  start_,
151  bPoint,
152  bFaceI,
153 
154  trackPt,
155  trackCellI,
156  trackFaceI
157  );
158 
159  //Info<< "after getTrackingPoint : "
160  // << " trackPt:" << trackPt
161  // << " trackCellI:" << trackCellI
162  // << " trackFaceI:" << trackFaceI
163  // << endl;
164 
165  if (trackCellI == -1)
166  {
167  // Line start_ - end_ does not intersect domain at all.
168  // (or is along edge)
169  // Set points and cell/face labels to empty lists
170  //Info<< "calcSamples : Both start_ and end_ outside domain"
171  // << endl;
172 
173  return;
174  }
175 
176  if (trackFaceI == -1)
177  {
178  // No boundary face. Check for nearish internal face
179  trackFaceI = findNearFace(trackCellI, trackPt, smallDist);
180  }
181 
182  //Info<< "calcSamples : got first point to track from :"
183  // << " trackPt:" << trackPt
184  // << " trackCell:" << trackCellI
185  // << " trackFace:" << trackFaceI
186  // << endl;
187 
188  //
189  // Track until hit end of all boundary intersections
190  //
191 
192  // current segment number
193  label segmentI = 0;
194 
195  // starting index of current segment in samplePts
196  label startSegmentI = 0;
197 
198  // index in bHits; current boundary intersection
199  label bHitI = 1;
200 
201  while(true)
202  {
203  if (trackFaceI != -1)
204  {
205  //Info<< "trackPt:" << trackPt << " on face so use." << endl;
206  samplingPts.append(trackPt);
207  samplingCells.append(trackCellI);
208  samplingFaces.append(trackFaceI);
209  samplingCurveDist.append(mag(trackPt - start_));
210  }
211 
212  // Initialize tracking starting from trackPt
213  Cloud<passiveParticle> particles(mesh(), IDLList<passiveParticle>());
214 
215  passiveParticle singleParticle
216  (
217  particles,
218  trackPt,
219  trackCellI
220  );
221 
222  bool reachedBoundary = trackToBoundary
223  (
224  singleParticle,
225  samplingPts,
226  samplingCells,
227  samplingFaces,
228  samplingCurveDist
229  );
230 
231  // fill sampleSegments
232  for(label i = samplingPts.size() - 1; i >= startSegmentI; --i)
233  {
234  samplingSegments.append(segmentI);
235  }
236 
237 
238  if (!reachedBoundary)
239  {
240  //Info<< "calcSamples : Reached end of samples: "
241  // << " samplePt now:" << singleParticle.position()
242  // << endl;
243  break;
244  }
245 
246 
247  // Go past boundary intersection where tracking stopped
248  // Use coordinate comparison instead of face comparison for
249  // accuracy reasons
250 
251  bool foundValidB = false;
252 
253  while (bHitI < bHits.size())
254  {
255  scalar dist =
256  (bHits[bHitI].hitPoint() - singleParticle.position())
257  & normOffset;
258 
259  //Info<< "Finding next boundary : "
260  // << "bPoint:" << bHits[bHitI].hitPoint()
261  // << " tracking:" << singleParticle.position()
262  // << " dist:" << dist
263  // << endl;
264 
265  if (dist > smallDist)
266  {
267  // hitpoint is past tracking position
268  foundValidB = true;
269  break;
270  }
271  else
272  {
273  bHitI++;
274  }
275  }
276 
277  if (!foundValidB)
278  {
279  // No valid boundary intersection found beyond tracking position
280  break;
281  }
282 
283  // Update starting point for tracking
284  trackFaceI = bHits[bHitI].index();
285  trackPt = pushIn(bHits[bHitI].hitPoint(), trackFaceI);
286  trackCellI = getBoundaryCell(trackFaceI);
287 
288  segmentI++;
289 
290  startSegmentI = samplingPts.size();
291  }
292 }
293 
294 
295 void Foam::faceOnlySet::genSamples()
296 {
297  // Storage for sample points
298  DynamicList<point> samplingPts;
299  DynamicList<label> samplingCells;
300  DynamicList<label> samplingFaces;
301  DynamicList<label> samplingSegments;
302  DynamicList<scalar> samplingCurveDist;
303 
304  calcSamples
305  (
306  samplingPts,
307  samplingCells,
308  samplingFaces,
309  samplingSegments,
310  samplingCurveDist
311  );
312 
313  samplingPts.shrink();
314  samplingCells.shrink();
315  samplingFaces.shrink();
316  samplingSegments.shrink();
317  samplingCurveDist.shrink();
318 
319  // Copy into *this
320  setSamples
321  (
322  samplingPts,
323  samplingCells,
324  samplingFaces,
325  samplingSegments,
326  samplingCurveDist
327  );
328 }
329 
330 
331 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
332 
334 (
335  const word& name,
336  const polyMesh& mesh,
337  meshSearch& searchEngine,
338  const word& axis,
339  const point& start,
340  const point& end
341 )
342 :
343  sampledSet(name, mesh, searchEngine, axis),
344  start_(start),
345  end_(end)
346 {
347  genSamples();
348 
349  if (debug)
350  {
351  write(Info);
352  }
353 }
354 
355 
357 (
358  const word& name,
359  const polyMesh& mesh,
360  meshSearch& searchEngine,
361  const dictionary& dict
362 )
363 :
364  sampledSet(name, mesh, searchEngine, dict),
365  start_(dict.lookup("start")),
366  end_(dict.lookup("end"))
367 {
368  genSamples();
369 
370  if (debug)
371  {
372  write(Info);
373  }
374 }
375 
376 
377 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
378 
380 {}
381 
382 
383 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
384 
386 {
387  return start_;
388 }
389 
390 
391 // ************************ vim: set sw=4 sts=4 et: ************************ //