meshSearch.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 "meshSearch.H"
27 #include "polyMesh.H"
28 #include "indexedOctree.H"
29 #include "DynamicList.H"
30 #include "demandDrivenData.H"
31 #include "treeDataCell.H"
32 #include "treeDataFace.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(meshSearch, 0);
39 
40  scalar meshSearch::tol_ = 1e-3;
41 }
42 
43 
44 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
45 
46 bool Foam::meshSearch::findNearer
47 (
48  const point& sample,
49  const pointField& points,
50  label& nearestI,
51  scalar& nearestDistSqr
52 )
53 {
54  bool nearer = false;
55 
56  forAll(points, pointi)
57  {
58  scalar distSqr = magSqr(points[pointi] - sample);
59 
60  if (distSqr < nearestDistSqr)
61  {
62  nearestDistSqr = distSqr;
63  nearestI = pointi;
64  nearer = true;
65  }
66  }
67 
68  return nearer;
69 }
70 
71 
72 bool Foam::meshSearch::findNearer
73 (
74  const point& sample,
75  const pointField& points,
76  const labelList& indices,
77  label& nearestI,
78  scalar& nearestDistSqr
79 )
80 {
81  bool nearer = false;
82 
83  forAll(indices, i)
84  {
85  label pointi = indices[i];
86 
87  scalar distSqr = magSqr(points[pointi] - sample);
88 
89  if (distSqr < nearestDistSqr)
90  {
91  nearestDistSqr = distSqr;
92  nearestI = pointi;
93  nearer = true;
94  }
95  }
96 
97  return nearer;
98 }
99 
100 
101 // tree based searching
102 Foam::label Foam::meshSearch::findNearestCellTree(const point& location) const
103 {
104  const indexedOctree<treeDataCell>& tree = cellTree();
105 
106  pointIndexHit info = tree.findNearest
107  (
108  location,
109  magSqr(tree.bb().max()-tree.bb().min())
110  );
111 
112  if (!info.hit())
113  {
114  info = tree.findNearest(location, Foam::sqr(great));
115  }
116  return info.index();
117 }
118 
119 
120 Foam::label Foam::meshSearch::findNearestCellLinear(const point& location) const
121 {
122  const vectorField& centres = mesh_.cellCentres();
123 
124  label nearestIndex = 0;
125  scalar minProximity = magSqr(centres[nearestIndex] - location);
126 
127  findNearer
128  (
129  location,
130  centres,
131  nearestIndex,
132  minProximity
133  );
134 
135  return nearestIndex;
136 }
137 
138 
139 Foam::label Foam::meshSearch::findNearestCellWalk
140 (
141  const point& location,
142  const label seedCelli
143 ) const
144 {
145  if (seedCelli < 0)
146  {
148  << "illegal seedCell:" << seedCelli << exit(FatalError);
149  }
150 
151  // Walk in direction of face that decreases distance
152 
153  label curCelli = seedCelli;
154  scalar distanceSqr = magSqr(mesh_.cellCentres()[curCelli] - location);
155 
156  bool closer;
157 
158  do
159  {
160  // Try neighbours of curCelli
161  closer = findNearer
162  (
163  location,
164  mesh_.cellCentres(),
165  mesh_.cellCells()[curCelli],
166  curCelli,
167  distanceSqr
168  );
169  } while (closer);
170 
171  return curCelli;
172 }
173 
174 
175 Foam::label Foam::meshSearch::findNearestFaceTree(const point& location) const
176 {
177  // Search nearest cell centre.
178  const indexedOctree<treeDataCell>& tree = cellTree();
179 
180  // Search with decent span
181  pointIndexHit info = tree.findNearest
182  (
183  location,
184  magSqr(tree.bb().max()-tree.bb().min())
185  );
186 
187  if (!info.hit())
188  {
189  // Search with desperate span
190  info = tree.findNearest(location, Foam::sqr(great));
191  }
192 
193 
194  // Now check any of the faces of the nearest cell
195  const vectorField& centres = mesh_.faceCentres();
196  const cell& ownFaces = mesh_.cells()[info.index()];
197 
198  label nearestFacei = ownFaces[0];
199  scalar minProximity = magSqr(centres[nearestFacei] - location);
200 
201  findNearer
202  (
203  location,
204  centres,
205  ownFaces,
206  nearestFacei,
207  minProximity
208  );
209 
210  return nearestFacei;
211 }
212 
213 
214 Foam::label Foam::meshSearch::findNearestFaceLinear(const point& location) const
215 {
216  const vectorField& centres = mesh_.faceCentres();
217 
218  label nearestFacei = 0;
219  scalar minProximity = magSqr(centres[nearestFacei] - location);
220 
221  findNearer
222  (
223  location,
224  centres,
225  nearestFacei,
226  minProximity
227  );
228 
229  return nearestFacei;
230 }
231 
232 
233 Foam::label Foam::meshSearch::findNearestFaceWalk
234 (
235  const point& location,
236  const label seedFacei
237 ) const
238 {
239  if (seedFacei < 0)
240  {
242  << "illegal seedFace:" << seedFacei << exit(FatalError);
243  }
244 
245  const vectorField& centres = mesh_.faceCentres();
246 
247 
248  // Walk in direction of face that decreases distance
249 
250  label curFacei = seedFacei;
251  scalar distanceSqr = magSqr(centres[curFacei] - location);
252 
253  while (true)
254  {
255  label betterFacei = curFacei;
256 
257  findNearer
258  (
259  location,
260  centres,
261  mesh_.cells()[mesh_.faceOwner()[curFacei]],
262  betterFacei,
263  distanceSqr
264  );
265 
266  if (mesh_.isInternalFace(curFacei))
267  {
268  findNearer
269  (
270  location,
271  centres,
272  mesh_.cells()[mesh_.faceNeighbour()[curFacei]],
273  betterFacei,
274  distanceSqr
275  );
276  }
277 
278  if (betterFacei == curFacei)
279  {
280  break;
281  }
282 
283  curFacei = betterFacei;
284  }
285 
286  return curFacei;
287 }
288 
289 
290 Foam::label Foam::meshSearch::findCellLinear(const point& location) const
291 {
292  bool cellFound = false;
293  label n = 0;
294 
295  label celli = -1;
296 
297  while ((!cellFound) && (n < mesh_.nCells()))
298  {
299  if (mesh_.pointInCell(location, n, cellDecompMode_))
300  {
301  cellFound = true;
302  celli = n;
303  }
304  else
305  {
306  n++;
307  }
308  }
309  if (cellFound)
310  {
311  return celli;
312  }
313  else
314  {
315  return -1;
316  }
317 }
318 
319 
320 Foam::label Foam::meshSearch::findCellWalk
321 (
322  const point& location,
323  const label seedCelli
324 ) const
325 {
326  if (seedCelli < 0)
327  {
329  << "illegal seedCell:" << seedCelli << exit(FatalError);
330  }
331 
332  if (mesh_.pointInCell(location, seedCelli, cellDecompMode_))
333  {
334  return seedCelli;
335  }
336 
337  // Walk in direction of face that decreases distance
338  label curCelli = seedCelli;
339  scalar nearestDistSqr = magSqr(mesh_.cellCentres()[curCelli] - location);
340 
341  while(true)
342  {
343  // Try neighbours of curCelli
344 
345  const cell& cFaces = mesh_.cells()[curCelli];
346 
347  label nearestCelli = -1;
348 
349  forAll(cFaces, i)
350  {
351  label facei = cFaces[i];
352 
353  if (mesh_.isInternalFace(facei))
354  {
355  label celli = mesh_.faceOwner()[facei];
356  if (celli == curCelli)
357  {
358  celli = mesh_.faceNeighbour()[facei];
359  }
360 
361  // Check if this is the correct cell
362  if (mesh_.pointInCell(location, celli, cellDecompMode_))
363  {
364  return celli;
365  }
366 
367  // Also calculate the nearest cell
368  scalar distSqr = magSqr(mesh_.cellCentres()[celli] - location);
369 
370  if (distSqr < nearestDistSqr)
371  {
372  nearestDistSqr = distSqr;
373  nearestCelli = celli;
374  }
375  }
376  }
377 
378  if (nearestCelli == -1)
379  {
380  return -1;
381  }
382 
383  // Continue with the nearest cell
384  curCelli = nearestCelli;
385  }
386 
387  return -1;
388 }
389 
390 
391 Foam::label Foam::meshSearch::findNearestBoundaryFaceWalk
392 (
393  const point& location,
394  const label seedFacei
395 ) const
396 {
397  if (seedFacei < 0)
398  {
400  << "illegal seedFace:" << seedFacei << exit(FatalError);
401  }
402 
403  // Start off from seedFacei
404 
405  label curFacei = seedFacei;
406 
407  const face& f = mesh_.faces()[curFacei];
408 
409  scalar minDist = f.nearestPoint
410  (
411  location,
412  mesh_.points()
413  ).distance();
414 
415  bool closer;
416 
417  do
418  {
419  closer = false;
420 
421  // Search through all neighbouring boundary faces by going
422  // across edges
423 
424  label lastFacei = curFacei;
425 
426  const labelList& myEdges = mesh_.faceEdges()[curFacei];
427 
428  forAll(myEdges, myEdgeI)
429  {
430  const labelList& neighbours = mesh_.edgeFaces()[myEdges[myEdgeI]];
431 
432  // Check any face which uses edge, is boundary face and
433  // is not curFacei itself.
434 
435  forAll(neighbours, nI)
436  {
437  label facei = neighbours[nI];
438 
439  if
440  (
441  (facei >= mesh_.nInternalFaces())
442  && (facei != lastFacei)
443  )
444  {
445  const face& f = mesh_.faces()[facei];
446 
447  pointHit curHit = f.nearestPoint
448  (
449  location,
450  mesh_.points()
451  );
452 
453  // If the face is closer, reset current face and distance
454  if (curHit.distance() < minDist)
455  {
456  minDist = curHit.distance();
457  curFacei = facei;
458  closer = true; // a closer neighbour has been found
459  }
460  }
461  }
462  }
463  } while (closer);
464 
465  return curFacei;
466 }
467 
468 
469 Foam::vector Foam::meshSearch::offset
470 (
471  const point& bPoint,
472  const label bFacei,
473  const vector& dir
474 ) const
475 {
476  // Get the neighbouring cell
477  label ownerCelli = mesh_.faceOwner()[bFacei];
478 
479  const point& c = mesh_.cellCentres()[ownerCelli];
480 
481  // Typical dimension: distance from point on face to cell centre
482  scalar typDim = mag(c - bPoint);
483 
484  return tol_*typDim*dir;
485 }
486 
487 
488 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
489 
490 Foam::meshSearch::meshSearch
491 (
492  const polyMesh& mesh,
493  const polyMesh::cellDecomposition cellDecompMode
494 )
495 :
496  mesh_(mesh),
497  cellDecompMode_(cellDecompMode)
498 {
499  if
500  (
501  cellDecompMode_ == polyMesh::FACE_DIAG_TRIS
502  || cellDecompMode_ == polyMesh::CELL_TETS)
503  {
504  // Force construction of face diagonals
505  (void)mesh.tetBasePtIs();
506  }
507 }
508 
509 
510 Foam::meshSearch::meshSearch
511 (
512  const polyMesh& mesh,
513  const treeBoundBox& bb,
514  const polyMesh::cellDecomposition cellDecompMode
515 )
516 :
517  mesh_(mesh),
518  cellDecompMode_(cellDecompMode)
519 {
520  overallBbPtr_.reset(new treeBoundBox(bb));
521 
522  if
523  (
524  cellDecompMode_ == polyMesh::FACE_DIAG_TRIS
525  || cellDecompMode_ == polyMesh::CELL_TETS
526  )
527  {
528  // Force construction of face diagonals
529  (void)mesh.tetBasePtIs();
530  }
531 }
532 
533 
534 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
535 
537 {
538  clearOut();
539 }
540 
541 
542 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
543 
546 {
547  if (!boundaryTreePtr_.valid())
548  {
549  //
550  // Construct tree
551  //
552 
553  if (!overallBbPtr_.valid())
554  {
555  overallBbPtr_.reset
556  (
557  new treeBoundBox(mesh_.points())
558  );
559 
560  treeBoundBox& overallBb = overallBbPtr_();
561 
562  // Extend slightly and make 3D
563  overallBb = overallBb.extend(1e-4);
564  }
565 
566  // all boundary faces (not just walls)
567  labelList bndFaces(mesh_.nFaces()-mesh_.nInternalFaces());
568  forAll(bndFaces, i)
569  {
570  bndFaces[i] = mesh_.nInternalFaces() + i;
571  }
572 
573  boundaryTreePtr_.reset
574  (
576  (
577  treeDataFace // all information needed to search faces
578  (
579  false, // do not cache bb
580  mesh_,
581  bndFaces // boundary faces only
582  ),
583  overallBbPtr_(), // overall search domain
584  8, // maxLevel
585  10, // leafsize
586  3.0 // duplicity
587  )
588  );
589  }
590 
591  return boundaryTreePtr_();
592 }
593 
594 
597 {
598  if (!cellTreePtr_.valid())
599  {
600  //
601  // Construct tree
602  //
603 
604  if (!overallBbPtr_.valid())
605  {
606  overallBbPtr_.reset
607  (
608  new treeBoundBox(mesh_.points())
609  );
610 
611  treeBoundBox& overallBb = overallBbPtr_();
612 
613  // Extend slightly and make 3D
614  overallBb = overallBb.extend(1e-4);
615  }
616 
617  cellTreePtr_.reset
618  (
620  (
622  (
623  false, // not cache bb
624  mesh_,
625  cellDecompMode_ // cell decomposition mode for inside tests
626  ),
627  overallBbPtr_(),
628  8, // maxLevel
629  10, // leafsize
630  6.0 // duplicity
631  )
632  );
633  }
634 
635  return cellTreePtr_();
636 }
637 
638 
640 (
641  const point& location,
642  const label seedCelli,
643  const bool useTreeSearch
644 ) const
645 {
646  if (seedCelli == -1)
647  {
648  if (useTreeSearch)
649  {
650  return findNearestCellTree(location);
651  }
652  else
653  {
654  return findNearestCellLinear(location);
655  }
656  }
657  else
658  {
659  return findNearestCellWalk(location, seedCelli);
660  }
661 }
662 
663 
665 (
666  const point& location,
667  const label seedFacei,
668  const bool useTreeSearch
669 ) const
670 {
671  if (seedFacei == -1)
672  {
673  if (useTreeSearch)
674  {
675  return findNearestFaceTree(location);
676  }
677  else
678  {
679  return findNearestFaceLinear(location);
680  }
681  }
682  else
683  {
684  return findNearestFaceWalk(location, seedFacei);
685  }
686 }
687 
688 
690 (
691  const point& location,
692  const label seedCelli,
693  const bool useTreeSearch
694 ) const
695 {
696  // Find the nearest cell centre to this location
697  if (seedCelli == -1)
698  {
699  if (useTreeSearch)
700  {
701  return cellTree().findInside(location);
702  }
703  else
704  {
705  return findCellLinear(location);
706  }
707  }
708  else
709  {
710  return findCellWalk(location, seedCelli);
711  }
712 }
713 
714 
716 (
717  const point& location,
718  const label seedFacei,
719  const bool useTreeSearch
720 ) const
721 {
722  if (seedFacei == -1)
723  {
724  if (useTreeSearch)
725  {
726  const indexedOctree<treeDataFace>& tree = boundaryTree();
727 
728  pointIndexHit info = boundaryTree().findNearest
729  (
730  location,
731  magSqr(tree.bb().max()-tree.bb().min())
732  );
733 
734  if (!info.hit())
735  {
736  info = boundaryTree().findNearest
737  (
738  location,
739  Foam::sqr(great)
740  );
741  }
742 
743  return tree.shapes().faceLabels()[info.index()];
744  }
745  else
746  {
747  scalar minDist = great;
748 
749  label minFacei = -1;
750 
751  for
752  (
753  label facei = mesh_.nInternalFaces();
754  facei < mesh_.nFaces();
755  facei++
756  )
757  {
758  const face& f = mesh_.faces()[facei];
759 
760  pointHit curHit =
761  f.nearestPoint
762  (
763  location,
764  mesh_.points()
765  );
766 
767  if (curHit.distance() < minDist)
768  {
769  minDist = curHit.distance();
770  minFacei = facei;
771  }
772  }
773  return minFacei;
774  }
775  }
776  else
777  {
778  return findNearestBoundaryFaceWalk(location, seedFacei);
779  }
780 }
781 
782 
784 (
785  const point& pStart,
786  const point& pEnd
787 ) const
788 {
789  pointIndexHit curHit = boundaryTree().findLine(pStart, pEnd);
790 
791  if (curHit.hit())
792  {
793  // Change index into octreeData into face label
794  curHit.setIndex(boundaryTree().shapes().faceLabels()[curHit.index()]);
795  }
796  return curHit;
797 }
798 
799 
801 (
802  const point& pStart,
803  const point& pEnd
804 ) const
805 {
807 
808  vector edgeVec = pEnd - pStart;
809  edgeVec /= mag(edgeVec);
810 
811  point pt = pStart;
812 
813  pointIndexHit bHit;
814  do
815  {
816  bHit = intersection(pt, pEnd);
817 
818  if (bHit.hit())
819  {
820  hits.append(bHit);
821 
822  const vector& area = mesh_.faceAreas()[bHit.index()];
823 
824  scalar typDim = Foam::sqrt(mag(area));
825 
826  if ((mag(bHit.hitPoint() - pEnd)/typDim) < small)
827  {
828  break;
829  }
830 
831  // Restart from hitPoint shifted a little bit in the direction
832  // of the destination
833 
834  pt =
835  bHit.hitPoint()
836  + offset(bHit.hitPoint(), bHit.index(), edgeVec);
837  }
838 
839  } while (bHit.hit());
840 
841 
842  hits.shrink();
843 
844  return hits;
845 }
846 
847 
849 {
850  return (boundaryTree().getVolumeType(p) == volumeType::INSIDE);
851 }
852 
853 
855 {
856  boundaryTreePtr_.clear();
857  cellTreePtr_.clear();
858  overallBbPtr_.clear();
859 }
860 
861 
863 {
864  clearOut();
865 }
866 
867 
868 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const indexedOctree< treeDataCell > & cellTree() const
Get (demand driven) reference to octree holding all cells.
Definition: meshSearch.C:596
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
label findNearestCell(const point &location, const label seedCelli=-1, const bool useTreeSearch=true) const
Find nearest cell in terms of cell centre.
Definition: meshSearch.C:640
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:841
Encapsulation of data needed to search for faces.
Definition: treeDataFace.H:58
static scalar tol_
Tolerance on linear dimensions.
Definition: meshSearch.H:163
dimensionedSymmTensor sqr(const dimensionedVector &dv)
label findNearestFace(const point &location, const label seedFacei=-1, const bool useTreeSearch=true) const
Definition: meshSearch.C:665
void setIndex(const label index)
dimensionedScalar sqrt(const dimensionedScalar &ds)
scalar distance(const vector &p1, const vector &p2)
Definition: curveTools.C:12
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
cellDecomposition
Enumeration defining the decomposition of the cell for.
Definition: polyMesh.H:98
scalar minDist(const List< pointIndexHit > &hitList)
const Point & hitPoint() const
Return hit point.
Foam::intersection.
Definition: intersection.H:49
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
void correct()
Correct for mesh geom/topo changes.
Definition: meshSearch.C:862
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
label findCell(const point &location, const label seedCelli=-1, const bool useTreeSearch=true) const
Find cell containing location.
Definition: meshSearch.C:690
const indexedOctree< treeDataFace > & boundaryTree() const
Get (demand driven) reference to octree holding all.
Definition: meshSearch.C:545
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:292
bool hit() const
Is there a hit.
const Type & shapes() const
Reference to shape.
List< label > labelList
A List of labels.
Definition: labelList.H:56
const treeBoundBox & bb() const
Top bounding box.
Templated 3D Vector derived from VectorSpace adding construction from 3 components, element access using x(), y() and z() member functions and the inner-product (dot-product) and cross product operators.
Definition: Vector.H:57
dimensioned< scalar > magSqr(const dimensioned< Type > &)
DynamicList< T, SizeInc, SizeMult, SizeDiv > & shrink()
Shrink the allocated space to the number of elements used.
Definition: DynamicListI.H:240
defineTypeNameAndDebug(combustionModel, 0)
Encapsulation of data needed to search in/for cells. Used to find the cell containing a point (e...
Definition: treeDataCell.H:54
void clearOut()
Delete all storage.
Definition: meshSearch.C:854
label findNearestBoundaryFace(const point &location, const label seedFacei=-1, const bool useTreeSearch=true) const
Find nearest boundary face.
Definition: meshSearch.C:716
Template functions to aid in the implementation of demand driven data.
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
vector point
Point is a vector.
Definition: point.H:41
pointHit nearestPoint(const point &p, const pointField &) const
Return nearest point to face.
~meshSearch()
Destructor.
Definition: meshSearch.C:536
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
Field< vector > vectorField
Specialisation of Field<T> for vector.
bool isInside(const point &) const
Determine inside/outside status.
Definition: meshSearch.C:848
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
PointHit< point > pointHit
Definition: pointHit.H:41
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
volScalarField & p
List< pointIndexHit > intersections(const point &pStart, const point &pEnd) const
Find all intersections of boundary within segment pStart .. pEnd.
Definition: meshSearch.C:801
pointIndexHit intersection(const point &pStart, const point &pEnd) const
Find first intersection of boundary in segment [pStart, pEnd].
Definition: meshSearch.C:784
label index() const
Return index.
const labelList & faceLabels() const
Definition: treeDataFace.H:179
treeBoundBox extend(const scalar s) const
Return asymetrically extended bounding box, with guaranteed.
Namespace for OpenFOAM.