Go to the documentation of this file.
36 bool Foam::triSurfaceSearch::checkUniqueHit
39 const UList<pointIndexHit>& hits,
47 const labelledTri&
f =
surface()[currHit.index()];
49 f.nearestPointClassify
61 const label nearPointi =
f[nearLabel];
68 const label pointFacei = pointFaces[pI];
70 if (pointFacei != currHit.index())
76 if (hit.index() == pointFacei)
91 const label edgeI = fEdges[nearLabel];
97 const label edgeFacei = edgeFaces[fI];
99 if (edgeFacei != currHit.index())
105 if (hit.index() == edgeFacei)
108 const vector currHitNormal =
111 const vector existingHitNormal =
114 const label signCurrHit =
115 pos0(currHitNormal & lineVec);
117 const label signExistingHit =
118 pos0(existingHitNormal & lineVec);
120 if (signCurrHit == signExistingHit)
136 Foam::triSurfaceSearch::triSurfaceSearch(
const triSurface& surface)
145 Foam::triSurfaceSearch::triSurfaceSearch
159 Info<<
" using intersection tolerance " << tolerance_ <<
endl;
165 Info<<
" using maximum tree depth " << maxTreeDepth_ <<
endl;
170 Foam::triSurfaceSearch::triSurfaceSearch
173 const scalar tolerance,
174 const label maxTreeDepth
178 tolerance_(tolerance),
179 maxTreeDepth_(maxTreeDepth),
213 if (surface().size())
221 <<
"Surface does not have compact point numbering."
222 <<
" Of " << surface().
points().size()
224 <<
" are used. This might give problems in some routines."
272 if (!tree().bb().contains(
sample))
274 inside[sampleI] =
false;
278 inside[sampleI] =
true;
282 inside[sampleI] =
false;
307 info[i] = octree.findNearest
326 const scalar nearestDistSqr = 0.25*
magSqr(span);
328 return tree().findNearest(pt, nearestDistSqr);
348 info[i] = octree.findLine(start[i],
end[i]);
const Field< point_type > & points() const noexcept
Return reference to global points.
List< label > labelList
A List of labels.
treeBoundBox extend(Random &rndGen, const scalar s) const
Return slightly wider bounding box.
const labelListList & pointFaces() const
Return point-face addressing.
const labelListList & edgeFaces() const
Return edge-face addressing.
static Vector< Cmpt > uniform(const Cmpt &s)
Return a VectorSpace with all elements = s.
pointIndexHit nearest(const point &pt, const vector &span) const
Calculate nearest point on surface for single searchPoint. Returns.
static constexpr const zero Zero
Global zero (0)
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit >> &info) const
Calculate all intersections from start to end.
Standard boundBox with extra functionality for use in octree.
void findNearest(const pointField &samples, const scalarField &nearestDistSqr, List< pointIndexHit > &info) const
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pos0(const dimensionedScalar &ds)
const point & max() const
Maximum describing the bounding box.
const labelListList & faceEdges() const
Return face-edge addressing.
pointIndexHit findLineAny(const point &start, const point &end) const
Find any intersection of line between start and end.
void clear() noexcept
Clear the addressed list, i.e. set the size to zero.
const point & min() const
Minimum describing the bounding box.
#define forAll(list, i)
Loop across all elements in list.
static scalar & perturbTol()
Get the perturbation tolerance.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
bool hit() const noexcept
Is there a hit?
Triangulated surface description with patch information.
messageStream Info
Information stream (stdout output on master, null elsewhere)
void setSize(const label n)
Alias for resize()
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
const indexedOctree< treeDataTriSurface > & tree() const
Demand driven construction of the octree.
~triSurfaceSearch()
Destructor.
void transfer(List< T > &list)
Non-pointer based hierarchical recursive searching.
scalarField samples(nIntervals, Zero)
void clearOut()
Clear storage.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Vector< scalar > vector
A scalar version of the templated Vector.
boolList calcInside(const pointField &searchPoints) const
Calculate for each searchPoint inside/outside status.
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
treeDataPrimitivePatch< triSurface > treeDataTriSurface
void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &info) const
PointIndexHit< point > pointIndexHit
A PointIndexHit for 3D points.
Encapsulation of data needed to search on PrimitivePatches.
label index() const noexcept
Return the hit index.
const triSurface & surface() const
Return reference to the surface.
const dimensionedScalar e
Elementary charge.
A location inside the volume.
const Field< point_type > & faceNormals() const
Return face unit normals for patch.
void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &info) const
const Map< label > & meshPointMap() const
Mesh point map.
#define WarningInFunction
Report a warning using Foam::Warning.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Minimal example by using system/controlDict.functions: