treeDataCell.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-2022 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 "treeDataCell.H"
27 #include "indexedOctree.H"
28 #include "polyMesh.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(treeDataCell, 0);
35 }
36 
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
39 
40 void Foam::treeDataCell::update()
41 {
42  if (cacheBb_)
43  {
44  bbs_.setSize(cellLabels_.size());
45 
46  forAll(cellLabels_, i)
47  {
48  bbs_[i] =
49  treeBoundBox
50  (
51  mesh_.cells()[cellLabels_[i]].bb
52  (
53  mesh_.points(),
54  mesh_.faces()
55  )
56  );
57  }
58  }
59 }
60 
61 
62 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63 
65 (
66  const bool cacheBb,
67  const polyMesh& mesh,
68  const labelUList& cellLabels,
69  const polyMesh::cellDecomposition decompMode
70 )
71 :
72  mesh_(mesh),
73  cellLabels_(cellLabels),
74  cacheBb_(cacheBb),
75  decompMode_(decompMode)
76 {
77  update();
78 }
79 
80 
82 (
83  const bool cacheBb,
84  const polyMesh& mesh,
85  labelList&& cellLabels,
86  const polyMesh::cellDecomposition decompMode
87 )
88 :
89  mesh_(mesh),
90  cellLabels_(move(cellLabels)),
91  cacheBb_(cacheBb),
92  decompMode_(decompMode)
93 {
94  update();
95 }
96 
97 
99 (
100  const bool cacheBb,
101  const polyMesh& mesh,
102  const polyMesh::cellDecomposition decompMode
103 )
104 :
105  mesh_(mesh),
106  cellLabels_(identity(mesh_.nCells())),
107  cacheBb_(cacheBb),
108  decompMode_(decompMode)
109 {
110  update();
111 }
112 
113 
115 (
116  const indexedOctree<treeDataCell>& tree
117 )
118 :
119  tree_(tree)
120 {}
121 
122 
124 (
125  const indexedOctree<treeDataCell>& tree
126 )
127 :
128  tree_(tree)
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
133 
135 {
136  pointField cc(cellLabels_.size());
137 
138  forAll(cellLabels_, i)
139  {
140  cc[i] = mesh_.cellCentres()[cellLabels_[i]];
141  }
142 
143  return cc;
144 }
145 
146 
148 (
149  const label index,
150  const treeBoundBox& cubeBb
151 ) const
152 {
153  if (cacheBb_)
154  {
155  return cubeBb.overlaps(bbs_[index]);
156  }
157  else
158  {
159  return cubeBb.overlaps
160  (
161  mesh_.cells()[cellLabels_[index]].bb
162  (
163  mesh_.points(),
164  mesh_.faces()
165  )
166  );
167  }
168 }
169 
170 
172 (
173  const label index,
174  const point& sample
175 ) const
176 {
177  return mesh_.pointInCell(sample, cellLabels_[index], decompMode_);
178 }
179 
180 
181 void Foam::treeDataCell::findNearestOp::operator()
182 (
183  const labelUList& indices,
184  const point& sample,
185 
186  scalar& nearestDistSqr,
187  label& minIndex,
188  point& nearestPoint
189 ) const
190 {
191  const treeDataCell& shape = tree_.shapes();
192 
193  forAll(indices, i)
194  {
195  label index = indices[i];
196  label celli = shape.cellLabels()[index];
197  scalar distSqr = magSqr(sample - shape.mesh().cellCentres()[celli]);
198 
199  if (distSqr < nearestDistSqr)
200  {
201  nearestDistSqr = distSqr;
202  minIndex = index;
203  nearestPoint = shape.mesh().cellCentres()[celli];
204  }
205  }
206 }
207 
208 
209 void Foam::treeDataCell::findNearestOp::operator()
210 (
211  const labelUList& indices,
212  const linePointRef& ln,
213 
214  treeBoundBox& tightest,
215  label& minIndex,
216  point& linePoint,
217  point& nearestPoint
218 ) const
219 {
221 }
222 
223 
224 bool Foam::treeDataCell::findIntersectOp::operator()
225 (
226  const label index,
227  const point& start,
228  const point& end,
229  point& intersectionPoint
230 ) const
231 {
232  const treeDataCell& shape = tree_.shapes();
233 
234  // Do quick rejection test
235  if (shape.cacheBb_)
236  {
237  const treeBoundBox& cellBb = shape.bbs_[index];
238 
239  if ((cellBb.posBits(start) & cellBb.posBits(end)) != 0)
240  {
241  // Start and end in same block outside of cellBb.
242  return false;
243  }
244  }
245  else
246  {
247  const treeBoundBox cellBb
248  (
249  shape.mesh_.cells()[shape.cellLabels_[index]].bb
250  (
251  shape.mesh_.points(),
252  shape.mesh_.faces()
253  )
254  );
255 
256  if ((cellBb.posBits(start) & cellBb.posBits(end)) != 0)
257  {
258  // Start and end in same block outside of cellBb.
259  return false;
260  }
261  }
262 
263 
264  // Do intersection with all faces of cell
265  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
266 
267  // Disable picking up intersections behind us.
268  scalar oldTol = intersection::setPlanarTol(0.0);
269 
270  const cell& cFaces = shape.mesh_.cells()[shape.cellLabels_[index]];
271 
272  const vector dir(end - start);
273  scalar minDistSqr = magSqr(dir);
274  bool hasMin = false;
275 
276  forAll(cFaces, i)
277  {
278  const face& f = shape.mesh_.faces()[cFaces[i]];
279 
280  pointHit inter = f.ray
281  (
282  start,
283  dir,
284  shape.mesh_.points(),
286  );
287 
288  if (inter.hit() && sqr(inter.distance()) <= minDistSqr)
289  {
290  // Note: no extra test on whether intersection is in front of us
291  // since using half_ray AND zero tolerance. (note that tolerance
292  // is used to look behind us)
293  minDistSqr = sqr(inter.distance());
294  intersectionPoint = inter.hitPoint();
295  hasMin = true;
296  }
297  }
298 
299  // Restore picking tolerance
301 
302  return hasMin;
303 }
304 
305 
306 // ************************************************************************* //
bool overlaps(const boundBox &) const
Overlaps/touches boundingBox?
Definition: boundBoxI.H:120
bool overlaps(const label index, const treeBoundBox &sampleBb) const
Does (bb of) shape at index overlap bb.
Definition: treeDataCell.C:148
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A line primitive.
Definition: line.H:56
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
dimensionedSymmTensor sqr(const dimensionedVector &dv)
pointHit ray(const point &p, const vector &n, const pointField &, const intersection::algorithm alg=intersection::algorithm::fullRay, const intersection::direction dir=intersection::direction::vector) const
Return potential intersection with face with a ray starting.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
const cellList & cells() const
direction posBits(const point &) const
Position of point relative to bounding box.
Definition: treeBoundBox.C:465
cellDecomposition
Enumeration defining the decomposition of the cell for.
Definition: polyMesh.H:100
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1211
static scalar setPlanarTol(const scalar t)
Set the planar tolerance, returning the previous value.
Definition: intersection.H:89
const polyMesh & mesh() const
Definition: treeDataCell.H:174
findIntersectOp(const indexedOctree< treeDataCell > &tree)
Definition: treeDataCell.C:124
const Point & hitPoint() const
Return hit point.
Definition: PointHit.H:126
bool hit() const
Is there a hit.
Definition: PointHit.H:120
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1224
const vectorField & cellCentres() const
treeDataCell(const bool cacheBb, const polyMesh &, const labelUList &, const polyMesh::cellDecomposition decompMode)
Construct from mesh and subset of cells.
Definition: treeDataCell.C:65
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.
Definition: POSIX.C:908
dimensioned< scalar > magSqr(const dimensioned< Type > &)
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
Encapsulation of data needed to search in/for cells. Used to find the cell containing a point (e...
Definition: treeDataCell.H:54
pointField shapePoints() const
Get representative point cloud for all shapes inside.
Definition: treeDataCell.C:134
findNearestOp(const indexedOctree< treeDataCell > &tree)
Definition: treeDataCell.C:115
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:57
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
const labelList & cellLabels() const
Definition: treeDataCell.H:169
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:353
bool contains(const label index, const point &sample) const
Does shape at index contain sample.
Definition: treeDataCell.C:172
Namespace for OpenFOAM.