treeDataCell.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2019 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "treeDataCell.H"
30 #include "indexedOctree.H"
31 #include "polyMesh.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(treeDataCell, 0);
38 }
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 Foam::treeBoundBox Foam::treeDataCell::calcCellBb(const label celli) const
44 {
45  const cellList& cells = mesh_.cells();
46  const faceList& faces = mesh_.faces();
47  const pointField& points = mesh_.points();
48 
49  treeBoundBox cellBb
50  (
51  vector(GREAT, GREAT, GREAT),
52  vector(-GREAT, -GREAT, -GREAT)
53  );
54 
55  const cell& cFaces = cells[celli];
56 
57  forAll(cFaces, cFacei)
58  {
59  const face& f = faces[cFaces[cFacei]];
60 
61  forAll(f, fp)
62  {
63  const point& p = points[f[fp]];
64 
65  cellBb.min() = min(cellBb.min(), p);
66  cellBb.max() = max(cellBb.max(), p);
67  }
68  }
69  return cellBb;
70 }
71 
72 
73 void Foam::treeDataCell::update()
74 {
75  if (cacheBb_)
76  {
77  bbs_.setSize(cellLabels_.size());
78 
79  forAll(cellLabels_, i)
80  {
81  bbs_[i] = calcCellBb(cellLabels_[i]);
82  }
83  }
84 }
85 
86 
87 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
88 
90 (
91  const bool cacheBb,
92  const polyMesh& mesh,
93  const labelUList& cellLabels,
94  const polyMesh::cellDecomposition decompMode
95 )
96 :
97  mesh_(mesh),
98  cellLabels_(cellLabels),
99  cacheBb_(cacheBb),
100  decompMode_(decompMode)
101 {
102  update();
103 }
104 
105 
107 (
108  const bool cacheBb,
109  const polyMesh& mesh,
110  labelList&& cellLabels,
111  const polyMesh::cellDecomposition decompMode
112 )
113 :
114  mesh_(mesh),
115  cellLabels_(std::move(cellLabels)),
116  cacheBb_(cacheBb),
117  decompMode_(decompMode)
118 {
119  update();
120 }
121 
122 
124 (
125  const bool cacheBb,
126  const polyMesh& mesh,
127  const polyMesh::cellDecomposition decompMode
128 )
129 :
130  mesh_(mesh),
131  cellLabels_(identity(mesh_.nCells())),
132  cacheBb_(cacheBb),
133  decompMode_(decompMode)
134 {
135  update();
136 }
137 
138 
140 (
141  const indexedOctree<treeDataCell>& tree
142 )
143 :
144  tree_(tree)
145 {}
146 
147 
149 (
150  const indexedOctree<treeDataCell>& tree
151 )
152 :
153  tree_(tree)
154 {}
155 
156 
157 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
158 
160 {
161  pointField cc(cellLabels_.size());
162 
163  forAll(cellLabels_, i)
164  {
165  cc[i] = mesh_.cellCentres()[cellLabels_[i]];
166  }
167 
168  return cc;
169 }
170 
171 
173 (
174  const label index,
175  const treeBoundBox& cubeBb
176 ) const
177 {
178  if (cacheBb_)
179  {
180  return cubeBb.overlaps(bbs_[index]);
181  }
182 
183  return cubeBb.overlaps(calcCellBb(cellLabels_[index]));
184 }
185 
186 
188 (
189  const label index,
190  const point& sample
191 ) const
192 {
193  return mesh_.pointInCell(sample, cellLabels_[index], decompMode_);
194 }
195 
196 
197 void Foam::treeDataCell::findNearestOp::operator()
198 (
199  const labelUList& indices,
200  const point& sample,
201 
202  scalar& nearestDistSqr,
203  label& minIndex,
204  point& nearestPoint
205 ) const
206 {
207  const treeDataCell& shape = tree_.shapes();
208 
209  forAll(indices, i)
210  {
211  label index = indices[i];
212  label celli = shape.cellLabels()[index];
213  scalar distSqr = magSqr(sample - shape.mesh().cellCentres()[celli]);
214 
215  if (distSqr < nearestDistSqr)
216  {
217  nearestDistSqr = distSqr;
218  minIndex = index;
219  nearestPoint = shape.mesh().cellCentres()[celli];
220  }
221  }
222 }
223 
224 
225 void Foam::treeDataCell::findNearestOp::operator()
226 (
227  const labelUList& indices,
228  const linePointRef& ln,
229 
230  treeBoundBox& tightest,
231  label& minIndex,
232  point& linePoint,
233  point& nearestPoint
234 ) const
235 {
237 }
238 
239 
240 bool Foam::treeDataCell::findIntersectOp::operator()
241 (
242  const label index,
243  const point& start,
244  const point& end,
245  point& intersectionPoint
246 ) const
247 {
248  const treeDataCell& shape = tree_.shapes();
249 
250  // Do quick rejection test
251  if (shape.cacheBb_)
252  {
253  const treeBoundBox& cellBb = shape.bbs_[index];
254 
255  if ((cellBb.posBits(start) & cellBb.posBits(end)) != 0)
256  {
257  // start and end in same block outside of cellBb.
258  return false;
259  }
260  }
261  else
262  {
263  const treeBoundBox cellBb = shape.calcCellBb(shape.cellLabels_[index]);
264 
265  if ((cellBb.posBits(start) & cellBb.posBits(end)) != 0)
266  {
267  // start and end in same block outside of cellBb.
268  return false;
269  }
270  }
271 
272 
273  // Do intersection with all faces of cell
274  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
275 
276  // Disable picking up intersections behind us.
277  scalar oldTol = intersection::setPlanarTol(0.0);
278 
279  const cell& cFaces = shape.mesh_.cells()[shape.cellLabels_[index]];
280 
281  const vector dir(end - start);
282  scalar minDistSqr = magSqr(dir);
283  bool hasMin = false;
284 
285  forAll(cFaces, i)
286  {
287  const face& f = shape.mesh_.faces()[cFaces[i]];
288 
289  pointHit inter = f.ray
290  (
291  start,
292  dir,
293  shape.mesh_.points(),
295  );
296 
297  if (inter.hit() && sqr(inter.distance()) <= minDistSqr)
298  {
299  // Note: no extra test on whether intersection is in front of us
300  // since using half_ray AND zero tolerance. (note that tolerance
301  // is used to look behind us)
302  minDistSqr = sqr(inter.distance());
303  intersectionPoint = inter.hitPoint();
304  hasMin = true;
305  }
306  }
307 
308  // Restore picking tolerance
310 
311  return hasMin;
312 }
313 
314 
315 // ************************************************************************* //
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Foam::polyMesh::points
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1069
Foam::treeBoundBox::overlaps
bool overlaps(const boundBox &bb) const
Overlaps other bounding box?
Definition: boundBoxI.H:221
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::polyMesh::cellDecomposition
cellDecomposition
Enumeration defining the decomposition of the cell for.
Definition: polyMesh.H:100
Foam::treeDataCell::shapePoints
pointField shapePoints() const
Get representative point cloud for all shapes inside.
Definition: treeDataCell.C:159
Foam::intersection::HALF_RAY
Definition: intersection.H:75
update
mesh update()
Foam::PointHit
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
Definition: PointHit.H:53
Foam::intersection::setPlanarTol
static scalar setPlanarTol(const scalar t)
Set the planar tolerance, returning the previous value.
Definition: intersection.H:94
Foam::treeBoundBox
Standard boundBox with extra functionality for use in octree.
Definition: treeBoundBox.H:86
indexedOctree.H
Foam::primitiveMesh::cells
const cellList & cells() const
Definition: primitiveMeshCells.C:138
Foam::treeDataCell::findNearestOp::findNearestOp
findNearestOp(const indexedOctree< treeDataCell > &tree)
Definition: treeDataCell.C:140
Foam::PointHit::distance
scalar distance() const noexcept
Return distance to hit.
Definition: PointHit.H:139
Foam::treeDataCell::mesh
const polyMesh & mesh() const
Definition: treeDataCell.H:177
Foam::PointHit::hitPoint
const point_type & hitPoint() const
Return the hit point. Fatal if not hit.
Definition: PointHit.H:145
polyMesh.H
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::magSqr
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Foam::treeBoundBox::posBits
direction posBits(const point &pt) const
Position of point relative to bounding box.
Definition: treeBoundBox.C:393
NotImplemented
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:517
Foam::treeDataCell
Encapsulation of data needed to search in/for cells. Used to find the cell containing a point (e....
Definition: treeDataCell.H:56
Foam::Field< vector >
Foam::treeDataCell::cellLabels
const labelList & cellLabels() const
Definition: treeDataCell.H:172
Foam::cellList
List< cell > cellList
A List of cells.
Definition: cellListFwd.H:47
Foam::indexedOctree
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:50
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
stdFoam::end
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Definition: stdFoam.H:121
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
treeDataCell.H
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::treeDataCell::overlaps
bool overlaps(const label index, const treeBoundBox &sampleBb) const
Does (bb of) shape at index overlap bb.
Definition: treeDataCell.C:173
Foam::primitiveMesh::cellCentres
const vectorField & cellCentres() const
Definition: primitiveMeshCellCentresAndVols.C:84
Foam::polyMesh::faces
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1094
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
f
labelList f(nPoints)
Foam::treeDataCell::treeDataCell
treeDataCell(const bool cacheBb, const polyMesh &mesh, const labelUList &cellLabels, const polyMesh::cellDecomposition decompMode)
Construct from mesh, copying subset of cells.
Definition: treeDataCell.C:90
Foam::faceList
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
Foam::Vector< scalar >
Foam::List< label >
Foam::UList< label >
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::identity
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
Definition: labelList.C:38
Foam::treeDataCell::contains
bool contains(const label index, const point &sample) const
Does shape at index contain sample.
Definition: treeDataCell.C:188
Foam::line
A line primitive.
Definition: line.H:53
Foam::ln
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.
Definition: MSwindows.C:925
Foam::face
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:72
cells
const cellShapeList & cells
Definition: gmvOutputHeader.H:3
Foam::treeDataCell::findIntersectOp::findIntersectOp
findIntersectOp(const indexedOctree< treeDataCell > &tree)
Definition: treeDataCell.C:149
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::PointHit::hit
bool hit() const noexcept
Is there a hit.
Definition: PointHit.H:121
Foam::cell
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:54
sample
Minimal example by using system/controlDict.functions: