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