diff --git a/applications/utilities/surface/surfaceSubset/surfaceSubset.C b/applications/utilities/surface/surfaceSubset/surfaceSubset.C
index 9df166e6e28ee194e01e067a32a994a3995d5999..195835523860a8ff9218ed7aaaef74a8a41a6970 100644
--- a/applications/utilities/surface/surfaceSubset/surfaceSubset.C
+++ b/applications/utilities/surface/surfaceSubset/surfaceSubset.C
@@ -31,6 +31,7 @@ Description
\*---------------------------------------------------------------------------*/
#include "triSurface.H"
+#include "triSurfaceSearch.H"
#include "argList.H"
#include "OFstream.H"
#include "IFstream.H"
@@ -40,6 +41,7 @@ Description
#include "indexedOctree.H"
#include "treeDataTriSurface.H"
#include "Random.H"
+#include "volumeType.H"
using namespace Foam;
@@ -242,26 +244,16 @@ int main(int argc, char *argv[])
// Read surface to select on
triSurface selectSurf(surfName);
- // bb of surface
- treeBoundBox bb(selectSurf.localPoints());
-
- // Random number generator
- Random rndGen(354543);
-
- // search engine
- indexedOctree<treeDataTriSurface> selectTree
+ triSurfaceSearch searchSelectSurf
(
- treeDataTriSurface
- (
- selectSurf,
- indexedOctree<treeDataTriSurface>::perturbTol()
- ),
- bb.extend(rndGen, 1e-4), // slightly randomize bb
- 8, // maxLevel
- 10, // leafsize
- 3.0 // duplicity
+ selectSurf,
+ indexedOctree<treeDataTriSurface>::perturbTol(),
+ 8
);
+ const indexedOctree<treeDataTriSurface>& selectTree =
+ searchSelectSurf.tree();
+
// Check if face (centre) is in outside or inside.
forAll(facesToSubset, faceI)
{
@@ -269,14 +261,13 @@ int main(int argc, char *argv[])
{
const point fc(surf1[faceI].centre(surf1.points()));
- indexedOctree<treeDataTriSurface>::volumeType t =
- selectTree.getVolumeType(fc);
+ volumeType t = selectTree.getVolumeType(fc);
if
(
outside
- ? (t == indexedOctree<treeDataTriSurface>::OUTSIDE)
- : (t == indexedOctree<treeDataTriSurface>::INSIDE)
+ ? (t == volumeType::OUTSIDE)
+ : (t == volumeType::INSIDE)
)
{
facesToSubset[faceI] = true;
diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index c7d000e8511396116aa6424da8cd63308809490b..19a8138fe5d1fc5f988befd6d5f388e639052c93 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -627,6 +627,7 @@ $(interpolationWeights)/splineInterpolationWeights/splineInterpolationWeights.C
algorithms/indexedOctree/indexedOctreeName.C
algorithms/indexedOctree/treeDataCell.C
+algorithms/indexedOctree/volumeType.C
algorithms/dynamicIndexedOctree/dynamicIndexedOctreeName.C
diff --git a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.C b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.C
index b18b2892cb8205db1ea80f0b522aa74aa8e3351e..92b2932510469a73455cdd89aedc3e43a8efb1fb 100644
--- a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.C
+++ b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.C
@@ -334,15 +334,14 @@ void Foam::dynamicIndexedOctree<Type>::recursiveSubDivision
// Recurses to determine status of lowest level boxes. Level above is
// combination of octants below.
template<class Type>
-typename Foam::dynamicIndexedOctree<Type>::volumeType
-Foam::dynamicIndexedOctree<Type>::calcVolumeType
+Foam::volumeType Foam::dynamicIndexedOctree<Type>::calcVolumeType
(
const label nodeI
) const
{
const node& nod = nodes_[nodeI];
- volumeType myType = UNKNOWN;
+ volumeType myType = volumeType::UNKNOWN;
for (direction octant = 0; octant < nod.subNodes_.size(); octant++)
{
@@ -359,7 +358,7 @@ Foam::dynamicIndexedOctree<Type>::calcVolumeType
{
// Contents. Depending on position in box might be on either
// side.
- subType = MIXED;
+ subType = volumeType::MIXED;
}
else
{
@@ -378,13 +377,13 @@ Foam::dynamicIndexedOctree<Type>::calcVolumeType
// Combine sub node types into type for treeNode. Result is 'mixed' if
// types differ among subnodes.
- if (myType == UNKNOWN)
+ if (myType == volumeType::UNKNOWN)
{
myType = subType;
}
else if (subType != myType)
{
- myType = MIXED;
+ myType = volumeType::MIXED;
}
}
return myType;
@@ -392,8 +391,7 @@ Foam::dynamicIndexedOctree<Type>::calcVolumeType
template<class Type>
-typename Foam::dynamicIndexedOctree<Type>::volumeType
-Foam::dynamicIndexedOctree<Type>::getVolumeType
+Foam::volumeType Foam::dynamicIndexedOctree<Type>::getVolumeType
(
const label nodeI,
const point& sample
@@ -403,22 +401,22 @@ Foam::dynamicIndexedOctree<Type>::getVolumeType
direction octant = nod.bb_.subOctant(sample);
- volumeType octantType = volumeType(nodeTypes_.get((nodeI<<3)+octant));
+ volumeType octantType = volumeType::type(nodeTypes_.get((nodeI<<3)+octant));
- if (octantType == INSIDE)
+ if (octantType == volumeType::INSIDE)
{
return octantType;
}
- else if (octantType == OUTSIDE)
+ else if (octantType == volumeType::OUTSIDE)
{
return octantType;
}
- else if (octantType == UNKNOWN)
+ else if (octantType == volumeType::UNKNOWN)
{
// Can happen for e.g. non-manifold surfaces.
return octantType;
}
- else if (octantType == MIXED)
+ else if (octantType == volumeType::MIXED)
{
labelBits index = nod.subNodes_[octant];
@@ -447,7 +445,7 @@ Foam::dynamicIndexedOctree<Type>::getVolumeType
<< "Empty subnode has invalid volume type MIXED."
<< abort(FatalError);
- return UNKNOWN;
+ return volumeType::UNKNOWN;
}
}
else
@@ -462,14 +460,13 @@ Foam::dynamicIndexedOctree<Type>::getVolumeType
<< "Node has invalid volume type " << octantType
<< abort(FatalError);
- return UNKNOWN;
+ return volumeType::UNKNOWN;
}
}
template<class Type>
-typename Foam::dynamicIndexedOctree<Type>::volumeType
-Foam::dynamicIndexedOctree<Type>::getSide
+Foam::volumeType Foam::dynamicIndexedOctree<Type>::getSide
(
const vector& outsideNormal,
const vector& vec
@@ -477,11 +474,11 @@ Foam::dynamicIndexedOctree<Type>::getSide
{
if ((outsideNormal&vec) >= 0)
{
- return OUTSIDE;
+ return volumeType::OUTSIDE;
}
else
{
- return INSIDE;
+ return volumeType::INSIDE;
}
}
@@ -2582,15 +2579,14 @@ const Foam::labelList& Foam::dynamicIndexedOctree<Type>::findIndices
// Determine type (inside/outside/mixed) per node.
template<class Type>
-typename Foam::dynamicIndexedOctree<Type>::volumeType
-Foam::dynamicIndexedOctree<Type>::getVolumeType
+Foam::volumeType Foam::dynamicIndexedOctree<Type>::getVolumeType
(
const point& sample
) const
{
if (nodes_.empty())
{
- return UNKNOWN;
+ return volumeType::UNKNOWN;
}
if (nodeTypes_.size() != 8*nodes_.size())
@@ -2598,7 +2594,7 @@ Foam::dynamicIndexedOctree<Type>::getVolumeType
// Calculate type for every octant of node.
nodeTypes_.setSize(8*nodes_.size());
- nodeTypes_ = UNKNOWN;
+ nodeTypes_ = volumeType::UNKNOWN;
calcVolumeType(0);
@@ -2611,21 +2607,21 @@ Foam::dynamicIndexedOctree<Type>::getVolumeType
forAll(nodeTypes_, i)
{
- volumeType type = volumeType(nodeTypes_.get(i));
+ volumeType type = volumeType::type(nodeTypes_.get(i));
- if (type == UNKNOWN)
+ if (type == volumeType::UNKNOWN)
{
nUNKNOWN++;
}
- else if (type == MIXED)
+ else if (type == volumeType::MIXED)
{
nMIXED++;
}
- else if (type == INSIDE)
+ else if (type == volumeType::INSIDE)
{
nINSIDE++;
}
- else if (type == OUTSIDE)
+ else if (type == volumeType::OUTSIDE)
{
nOUTSIDE++;
}
diff --git a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.H b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.H
index 0bbf7a0e2ea76808020ba80aa74af13abb807eec..e0a298855f9019869207fa84cde1047c9828fd44 100644
--- a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.H
+++ b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicIndexedOctree.H
@@ -43,6 +43,7 @@ SourceFiles
#include "HashSet.H"
#include "labelBits.H"
#include "PackedList.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -83,16 +84,6 @@ public:
// Data types
- //- volume types
- enum volumeType
- {
- UNKNOWN = 0,
- MIXED = 1,
- INSIDE = 2,
- OUTSIDE = 3
- };
-
-
//- Tree node. Has up pointer and down pointers.
class node
{
diff --git a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.C b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.C
index ebc93a63c7b5820efe1353835a1fb6fd567a5a1a..821e4ff5007bb35a1f2d07d45a98ef2a2ca0f01c 100644
--- a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.C
+++ b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -57,13 +57,13 @@ Foam::dynamicTreeDataPoint::shapePoints() const
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
-Foam::label Foam::dynamicTreeDataPoint::getVolumeType
+Foam::volumeType Foam::dynamicTreeDataPoint::getVolumeType
(
const dynamicIndexedOctree<dynamicTreeDataPoint>& oc,
const point& sample
) const
{
- return dynamicIndexedOctree<dynamicTreeDataPoint>::UNKNOWN;
+ return volumeType::UNKNOWN;
}
diff --git a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.H b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.H
index 2cbdc0005798ba88ff03bc96cee55aebbca51589..529d923b742cadf046b19348a8d3f1b25d23a31a 100644
--- a/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.H
+++ b/src/OpenFOAM/algorithms/dynamicIndexedOctree/dynamicTreeDataPoint.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -43,6 +43,7 @@ SourceFiles
#include "pointField.H"
#include "treeBoundBox.H"
#include "linePointRef.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -92,7 +93,7 @@ public:
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
- label getVolumeType
+ volumeType getVolumeType
(
const dynamicIndexedOctree<dynamicTreeDataPoint>&,
const point&
diff --git a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C
index 83ac739941284b74848944396492962e4d5f5964..945f9d989ce6a9f870c2a6ba7732a28decd2b62f 100644
--- a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C
+++ b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C
@@ -356,15 +356,14 @@ Foam::label Foam::indexedOctree<Type>::compactContents
// Recurses to determine status of lowest level boxes. Level above is
// combination of octants below.
template<class Type>
-typename Foam::indexedOctree<Type>::volumeType
-Foam::indexedOctree<Type>::calcVolumeType
+Foam::volumeType Foam::indexedOctree<Type>::calcVolumeType
(
const label nodeI
) const
{
const node& nod = nodes_[nodeI];
- volumeType myType = UNKNOWN;
+ volumeType myType = volumeType::UNKNOWN;
for (direction octant = 0; octant < nod.subNodes_.size(); octant++)
{
@@ -381,7 +380,7 @@ Foam::indexedOctree<Type>::calcVolumeType
{
// Contents. Depending on position in box might be on either
// side.
- subType = MIXED;
+ subType = volumeType::MIXED;
}
else
{
@@ -389,10 +388,7 @@ Foam::indexedOctree<Type>::calcVolumeType
// of its bounding box.
const treeBoundBox subBb = nod.bb_.subBbox(octant);
- subType = volumeType
- (
- shapes_.getVolumeType(*this, subBb.midpoint())
- );
+ subType = shapes_.getVolumeType(*this, subBb.midpoint());
}
// Store octant type
@@ -400,13 +396,13 @@ Foam::indexedOctree<Type>::calcVolumeType
// Combine sub node types into type for treeNode. Result is 'mixed' if
// types differ among subnodes.
- if (myType == UNKNOWN)
+ if (myType == volumeType::UNKNOWN)
{
myType = subType;
}
else if (subType != myType)
{
- myType = MIXED;
+ myType = volumeType::MIXED;
}
}
return myType;
@@ -414,8 +410,7 @@ Foam::indexedOctree<Type>::calcVolumeType
template<class Type>
-typename Foam::indexedOctree<Type>::volumeType
-Foam::indexedOctree<Type>::getVolumeType
+Foam::volumeType Foam::indexedOctree<Type>::getVolumeType
(
const label nodeI,
const point& sample
@@ -425,22 +420,22 @@ Foam::indexedOctree<Type>::getVolumeType
direction octant = nod.bb_.subOctant(sample);
- volumeType octantType = volumeType(nodeTypes_.get((nodeI<<3)+octant));
+ volumeType octantType = volumeType::type(nodeTypes_.get((nodeI<<3)+octant));
- if (octantType == INSIDE)
+ if (octantType == volumeType::INSIDE)
{
return octantType;
}
- else if (octantType == OUTSIDE)
+ else if (octantType == volumeType::OUTSIDE)
{
return octantType;
}
- else if (octantType == UNKNOWN)
+ else if (octantType == volumeType::UNKNOWN)
{
// Can happen for e.g. non-manifold surfaces.
return octantType;
}
- else if (octantType == MIXED)
+ else if (octantType == volumeType::MIXED)
{
labelBits index = nod.subNodes_[octant];
@@ -469,7 +464,7 @@ Foam::indexedOctree<Type>::getVolumeType
<< "Empty subnode has invalid volume type MIXED."
<< abort(FatalError);
- return UNKNOWN;
+ return volumeType::UNKNOWN;
}
}
else
@@ -484,14 +479,13 @@ Foam::indexedOctree<Type>::getVolumeType
<< "Node has invalid volume type " << octantType
<< abort(FatalError);
- return UNKNOWN;
+ return volumeType::UNKNOWN;
}
}
template<class Type>
-typename Foam::indexedOctree<Type>::volumeType
-Foam::indexedOctree<Type>::getSide
+Foam::volumeType Foam::indexedOctree<Type>::getSide
(
const vector& outsideNormal,
const vector& vec
@@ -499,11 +493,11 @@ Foam::indexedOctree<Type>::getSide
{
if ((outsideNormal&vec) >= 0)
{
- return OUTSIDE;
+ return volumeType::OUTSIDE;
}
else
{
- return INSIDE;
+ return volumeType::INSIDE;
}
}
@@ -576,6 +570,7 @@ Foam::indexedOctree<Type>::getSide
// Find nearest point starting from nodeI
template<class Type>
+template<class FindNearestOp>
void Foam::indexedOctree<Type>::findNearest
(
const label nodeI,
@@ -583,7 +578,9 @@ void Foam::indexedOctree<Type>::findNearest
scalar& nearestDistSqr,
label& nearestShapeI,
- point& nearestPoint
+ point& nearestPoint,
+
+ const FindNearestOp& fnOp
) const
{
const node& nod = nodes_[nodeI];
@@ -614,7 +611,9 @@ void Foam::indexedOctree<Type>::findNearest
nearestDistSqr,
nearestShapeI,
- nearestPoint
+ nearestPoint,
+
+ fnOp
);
}
}
@@ -631,7 +630,7 @@ void Foam::indexedOctree<Type>::findNearest
)
)
{
- shapes_.findNearest
+ fnOp
(
contents_[getContent(index)],
sample,
@@ -648,6 +647,7 @@ void Foam::indexedOctree<Type>::findNearest
// Find nearest point to line.
template<class Type>
+template<class FindNearestOp>
void Foam::indexedOctree<Type>::findNearest
(
const label nodeI,
@@ -656,7 +656,9 @@ void Foam::indexedOctree<Type>::findNearest
treeBoundBox& tightest,
label& nearestShapeI,
point& linePoint,
- point& nearestPoint
+ point& nearestPoint,
+
+ const FindNearestOp& fnOp
) const
{
const node& nod = nodes_[nodeI];
@@ -687,7 +689,9 @@ void Foam::indexedOctree<Type>::findNearest
tightest,
nearestShapeI,
linePoint,
- nearestPoint
+ nearestPoint,
+
+ fnOp
);
}
}
@@ -697,7 +701,7 @@ void Foam::indexedOctree<Type>::findNearest
if (subBb.overlaps(tightest))
{
- shapes_.findNearest
+ fnOp
(
contents_[getContent(index)],
ln,
@@ -1620,6 +1624,7 @@ Foam::word Foam::indexedOctree<Type>::faceString
// hitInfo.point = coordinate of intersection of ray with bounding box
// hitBits = posbits of point on bounding box
template<class Type>
+template<class FindIntersectOp>
void Foam::indexedOctree<Type>::traverseNode
(
const bool findAny,
@@ -1632,7 +1637,9 @@ void Foam::indexedOctree<Type>::traverseNode
const direction octant,
pointIndexHit& hitInfo,
- direction& hitBits
+ direction& hitBits,
+
+ const FindIntersectOp& fiOp
) const
{
if (debug)
@@ -1667,7 +1674,7 @@ void Foam::indexedOctree<Type>::traverseNode
label shapeI = indices[elemI];
point pt;
- bool hit = shapes_.intersects(shapeI, start, end, pt);
+ bool hit = fiOp(shapeI, start, end, pt);
// Note that intersection of shape might actually be
// in a neighbouring box. For findAny this is not important.
@@ -1695,13 +1702,7 @@ void Foam::indexedOctree<Type>::traverseNode
label shapeI = indices[elemI];
point pt;
- bool hit = shapes_.intersects
- (
- shapeI,
- start,
- nearestPoint,
- pt
- );
+ bool hit = fiOp(shapeI, start, nearestPoint, pt);
// Note that intersection of shape might actually be
// in a neighbouring box. Since we need to maintain strict
@@ -1774,7 +1775,9 @@ void Foam::indexedOctree<Type>::traverseNode
octant,
hitInfo,
- hitBits
+ hitBits,
+
+ fiOp
);
}
}
@@ -1782,6 +1785,7 @@ void Foam::indexedOctree<Type>::traverseNode
// Find first intersection
template<class Type>
+template<class FindIntersectOp>
Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
(
const bool findAny,
@@ -1789,6 +1793,7 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
const point& treeEnd,
const label startNodeI,
const direction startOctant,
+ const FindIntersectOp& fiOp,
const bool verbose
) const
{
@@ -1864,7 +1869,9 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
octant,
hitInfo,
- hitFaceID
+ hitFaceID,
+
+ fiOp
);
// Did we hit a triangle?
@@ -1948,7 +1955,8 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
treeEnd,
startNodeI,
startOctant,
- true //verbose
+ fiOp,
+ true //verbose,
);
}
if (debug)
@@ -2007,11 +2015,13 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
// Find first intersection
template<class Type>
+template<class FindIntersectOp>
Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
(
const bool findAny,
const point& start,
- const point& end
+ const point& end,
+ const FindIntersectOp& fiOp
) const
{
pointIndexHit hitInfo;
@@ -2069,7 +2079,8 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
trackStart,
trackEnd,
parentNodeI,
- octant
+ octant,
+ fiOp
);
}
@@ -2656,6 +2667,25 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
const point& sample,
const scalar startDistSqr
) const
+{
+ return findNearest
+ (
+ sample,
+ startDistSqr,
+ typename Type::findNearestOp(*this)
+ );
+}
+
+
+template <class Type>
+template <class FindNearestOp>
+Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
+(
+ const point& sample,
+ const scalar startDistSqr,
+
+ const FindNearestOp& fnOp
+) const
{
scalar nearestDistSqr = startDistSqr;
label nearestShapeI = -1;
@@ -2670,7 +2700,9 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
nearestDistSqr,
nearestShapeI,
- nearestPoint
+ nearestPoint,
+
+ fnOp
);
}
@@ -2685,6 +2717,27 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
treeBoundBox& tightest,
point& linePoint
) const
+{
+ return findNearest
+ (
+ ln,
+ tightest,
+ linePoint,
+ typename Type::findNearestOp(*this)
+ );
+}
+
+
+template <class Type>
+template <class FindNearestOp>
+Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
+(
+ const linePointRef& ln,
+ treeBoundBox& tightest,
+ point& linePoint,
+
+ const FindNearestOp& fnOp
+) const
{
label nearestShapeI = -1;
point nearestPoint = vector::zero;
@@ -2699,7 +2752,9 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
tightest,
nearestShapeI,
linePoint,
- nearestPoint
+ nearestPoint,
+
+ fnOp
);
}
@@ -2715,7 +2770,13 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
const point& end
) const
{
- return findLine(false, start, end);
+ return findLine
+ (
+ false,
+ start,
+ end,
+ typename Type::findIntersectOp(*this)
+ );
}
@@ -2727,7 +2788,41 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLineAny
const point& end
) const
{
- return findLine(true, start, end);
+ return findLine
+ (
+ true,
+ start,
+ end,
+ typename Type::findIntersectOp(*this)
+ );
+}
+
+
+// Find nearest intersection
+template <class Type>
+template <class FindIntersectOp>
+Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
+(
+ const point& start,
+ const point& end,
+ const FindIntersectOp& fiOp
+) const
+{
+ return findLine(false, start, end, fiOp);
+}
+
+
+// Find nearest intersection
+template <class Type>
+template <class FindIntersectOp>
+Foam::pointIndexHit Foam::indexedOctree<Type>::findLineAny
+(
+ const point& start,
+ const point& end,
+ const FindIntersectOp& fiOp
+) const
+{
+ return findLine(true, start, end, fiOp);
}
@@ -2871,23 +2966,29 @@ const Foam::labelList& Foam::indexedOctree<Type>::findIndices
// Determine type (inside/outside/mixed) per node.
template<class Type>
-typename Foam::indexedOctree<Type>::volumeType
-Foam::indexedOctree<Type>::getVolumeType
+Foam::volumeType Foam::indexedOctree<Type>::getVolumeType
(
const point& sample
) const
{
if (nodes_.empty())
{
- return UNKNOWN;
+ return volumeType::UNKNOWN;
}
+// // If the sample is not within the octree, then have to query shapes
+// // directly
+// if (!nodes_[0].bb_.contains(sample))
+// {
+// return volumeType(shapes_.getVolumeType(*this, sample));
+// }
+
if (nodeTypes_.size() != 8*nodes_.size())
{
// Calculate type for every octant of node.
nodeTypes_.setSize(8*nodes_.size());
- nodeTypes_ = UNKNOWN;
+ nodeTypes_ = volumeType::UNKNOWN;
calcVolumeType(0);
@@ -2900,21 +3001,21 @@ Foam::indexedOctree<Type>::getVolumeType
forAll(nodeTypes_, i)
{
- volumeType type = volumeType(nodeTypes_.get(i));
+ volumeType type = volumeType::type(nodeTypes_.get(i));
- if (type == UNKNOWN)
+ if (type == volumeType::UNKNOWN)
{
nUNKNOWN++;
}
- else if (type == MIXED)
+ else if (type == volumeType::MIXED)
{
nMIXED++;
}
- else if (type == INSIDE)
+ else if (type == volumeType::INSIDE)
{
nINSIDE++;
}
- else if (type == OUTSIDE)
+ else if (type == volumeType::OUTSIDE)
{
nOUTSIDE++;
}
diff --git a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H
index 6a5c56ac7bb2d7f7900677668d274beeffba177e..229993e406ca3d56ac06b6a2a610983901c260eb 100644
--- a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H
+++ b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H
@@ -42,6 +42,7 @@ SourceFiles
#include "HashSet.H"
#include "labelBits.H"
#include "PackedList.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -74,16 +75,6 @@ public:
// Data types
- //- volume types
- enum volumeType
- {
- UNKNOWN = 0,
- MIXED = 1,
- INSIDE = 2,
- OUTSIDE = 3
- };
-
-
//- Tree node. Has up pointer and down pointers.
class node
{
@@ -212,6 +203,7 @@ private:
// Query
//- Find nearest point to line.
+ template <class FindNearestOp>
void findNearest
(
const label nodeI,
@@ -220,7 +212,9 @@ private:
treeBoundBox& tightest,
label& nearestShapeI,
point& linePoint,
- point& nearestPoint
+ point& nearestPoint,
+
+ const FindNearestOp& fnOp
) const;
//- Return bbox of octant
@@ -294,6 +288,7 @@ private:
// intersection point.
// findAny=true : return any intersection
// findAny=false: return nearest (to start) intersection
+ template <class FindIntersectOp>
void traverseNode
(
const bool findAny,
@@ -306,10 +301,13 @@ private:
const direction octantI,
pointIndexHit& hitInfo,
- direction& faceID
+ direction& faceID,
+
+ const FindIntersectOp& fiOp
) const;
//- Find any or nearest intersection
+ template <class FindIntersectOp>
pointIndexHit findLine
(
const bool findAny,
@@ -317,6 +315,7 @@ private:
const point& treeEnd,
const label startNodeI,
const direction startOctantI,
+ const FindIntersectOp& fiOp,
const bool verbose = false
) const;
@@ -328,11 +327,13 @@ private:
// ) const;
//- Find any or nearest intersection of line between start and end.
+ template <class FindIntersectOp>
pointIndexHit findLine
(
const bool findAny,
const point& start,
- const point& end
+ const point& end,
+ const FindIntersectOp& fiOp
) const;
//- Find all elements intersecting box.
@@ -528,15 +529,24 @@ public:
// Queries
+ pointIndexHit findNearest
+ (
+ const point& sample,
+ const scalar nearestDistSqr
+ ) const;
+
//- Calculate nearest point on nearest shape.
// Returns
// - bool : any point found nearer than nearestDistSqr
// - label: index in shapes
// - point: actual nearest point found
+ template <class FindNearestOp>
pointIndexHit findNearest
(
const point& sample,
- const scalar nearestDistSqr
+ const scalar nearestDistSqr,
+
+ const FindNearestOp& fnOp
) const;
// bool findAnyOverlap
@@ -553,6 +563,7 @@ public:
// ) const;
//- Low level: calculate nearest starting from subnode.
+ template <class FindNearestOp>
void findNearest
(
const label nodeI,
@@ -560,7 +571,9 @@ public:
scalar& nearestDistSqr,
label& nearestShapeI,
- point& nearestPoint
+ point& nearestPoint,
+
+ const FindNearestOp& fnOp
) const;
//- Find nearest to line.
@@ -577,6 +590,16 @@ public:
point& linePoint
) const;
+ template <class FindNearestOp>
+ pointIndexHit findNearest
+ (
+ const linePointRef& ln,
+ treeBoundBox& tightest,
+ point& linePoint,
+
+ const FindNearestOp& fnOp
+ ) const;
+
//- Find nearest intersection of line between start and end.
pointIndexHit findLine
(
@@ -591,6 +614,24 @@ public:
const point& end
) const;
+ //- Find nearest intersection of line between start and end.
+ template <class FindIntersectOp>
+ pointIndexHit findLine
+ (
+ const point& start,
+ const point& end,
+ const FindIntersectOp& fiOp
+ ) const;
+
+ //- Find any intersection of line between start and end.
+ template <class FindIntersectOp>
+ pointIndexHit findLineAny
+ (
+ const point& start,
+ const point& end,
+ const FindIntersectOp& fiOp
+ ) const;
+
//- Find (in no particular order) indices of all shapes inside or
// overlapping bounding box (i.e. all shapes not outside box)
labelList findBox(const treeBoundBox& bb) const;
diff --git a/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.C b/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.C
index e4d78fc83d6d5e22bbad3598446f546dc28f391c..2ab17d0de855e57bbb6def51fefbc446da6fa3ce 100644
--- a/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.C
+++ b/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -133,6 +133,24 @@ Foam::treeDataCell::treeDataCell
}
+Foam::treeDataCell::findNearestOp::findNearestOp
+(
+ const indexedOctree<treeDataCell>& tree
+)
+:
+ tree_(tree)
+{}
+
+
+Foam::treeDataCell::findIntersectOp::findIntersectOp
+(
+ const indexedOctree<treeDataCell>& tree
+)
+:
+ tree_(tree)
+{}
+
+
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::pointField Foam::treeDataCell::shapePoints() const
@@ -175,7 +193,7 @@ bool Foam::treeDataCell::contains
}
-void Foam::treeDataCell::findNearest
+void Foam::treeDataCell::findNearestOp::operator()
(
const labelUList& indices,
const point& sample,
@@ -185,23 +203,51 @@ void Foam::treeDataCell::findNearest
point& nearestPoint
) const
{
+ const treeDataCell& shape = tree_.shapes();
+
forAll(indices, i)
{
label index = indices[i];
- label cellI = cellLabels_[index];
- scalar distSqr = magSqr(sample - mesh_.cellCentres()[cellI]);
+ label cellI = shape.cellLabels()[index];
+ scalar distSqr = magSqr(sample - shape.mesh().cellCentres()[cellI]);
if (distSqr < nearestDistSqr)
{
nearestDistSqr = distSqr;
minIndex = index;
- nearestPoint = mesh_.cellCentres()[cellI];
+ nearestPoint = shape.mesh().cellCentres()[cellI];
}
}
}
-bool Foam::treeDataCell::intersects
+void Foam::treeDataCell::findNearestOp::operator()
+(
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+) const
+{
+ notImplemented
+ (
+ "treeDataCell::findNearestOp::operator()"
+ "("
+ " const labelUList&,"
+ " const linePointRef&,"
+ " treeBoundBox&,"
+ " label&,"
+ " point&,"
+ " point&"
+ ") const"
+ );
+}
+
+
+bool Foam::treeDataCell::findIntersectOp::operator()
(
const label index,
const point& start,
@@ -209,10 +255,12 @@ bool Foam::treeDataCell::intersects
point& intersectionPoint
) const
{
+ const treeDataCell& shape = tree_.shapes();
+
// Do quick rejection test
- if (cacheBb_)
+ if (shape.cacheBb_)
{
- const treeBoundBox& cellBb = bbs_[index];
+ const treeBoundBox& cellBb = shape.bbs_[index];
if ((cellBb.posBits(start) & cellBb.posBits(end)) != 0)
{
@@ -222,7 +270,7 @@ bool Foam::treeDataCell::intersects
}
else
{
- const treeBoundBox cellBb = calcCellBb(cellLabels_[index]);
+ const treeBoundBox cellBb = shape.calcCellBb(shape.cellLabels_[index]);
if ((cellBb.posBits(start) & cellBb.posBits(end)) != 0)
{
@@ -238,7 +286,7 @@ bool Foam::treeDataCell::intersects
// Disable picking up intersections behind us.
scalar oldTol = intersection::setPlanarTol(0.0);
- const cell& cFaces = mesh_.cells()[cellLabels_[index]];
+ const cell& cFaces = shape.mesh_.cells()[shape.cellLabels_[index]];
const vector dir(end - start);
scalar minDistSqr = magSqr(dir);
@@ -246,13 +294,13 @@ bool Foam::treeDataCell::intersects
forAll(cFaces, i)
{
- const face& f = mesh_.faces()[cFaces[i]];
+ const face& f = shape.mesh_.faces()[cFaces[i]];
pointHit inter = f.ray
(
start,
dir,
- mesh_.points(),
+ shape.mesh_.points(),
intersection::HALF_RAY
);
diff --git a/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.H b/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.H
index 476266faf9a741f1c93eb6962876331fb67fa760..c302f6cbecedd3e1400f29ce40f531145f805dc1 100644
--- a/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.H
+++ b/src/OpenFOAM/algorithms/indexedOctree/treeDataCell.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -38,6 +38,7 @@ SourceFiles
#include "polyMesh.H"
#include "treeBoundBoxList.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -80,6 +81,56 @@ class treeDataCell
public:
+
+ class findNearestOp
+ {
+ const indexedOctree<treeDataCell>& tree_;
+
+ public:
+
+ findNearestOp(const indexedOctree<treeDataCell>& tree);
+
+ void operator()
+ (
+ const labelUList& indices,
+ const point& sample,
+
+ scalar& nearestDistSqr,
+ label& minIndex,
+ point& nearestPoint
+ ) const;
+
+ void operator()
+ (
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+ ) const;
+ };
+
+
+ class findIntersectOp
+ {
+ const indexedOctree<treeDataCell>& tree_;
+
+ public:
+
+ findIntersectOp(const indexedOctree<treeDataCell>& tree);
+
+ bool operator()
+ (
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ ) const;
+ };
+
+
// Declare name of the class and its debug switch
ClassName("treeDataCell");
@@ -146,7 +197,7 @@ public:
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
- label getVolumeType
+ volumeType getVolumeType
(
const indexedOctree<treeDataCell>&,
const point&
@@ -157,7 +208,7 @@ public:
"treeDataCell::getVolumeType"
"(const indexedOctree<treeDataCell>&, const point&)"
);
- return -1;
+ return volumeType::UNKNOWN;
}
//- Does (bb of) shape at index overlap bb
@@ -173,49 +224,6 @@ public:
const label index,
const point& sample
) const;
-
- //- Calculates nearest (to sample) point in shape.
- // Returns actual point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& nearestIndex,
- point& nearestPoint
- ) const;
-
- //- Calculates nearest (to line) point in shape.
- // Returns point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
- ) const
- {
- notImplemented
- (
- "treeDataCell::findNearest"
- "(const labelUList&, const linePointRef&, ..)"
- );
- }
-
- //- Calculate intersection of shape with ray. Sets result
- // accordingly
- bool intersects
- (
- const label index,
- const point& start,
- const point& end,
- point& result
- ) const;
-
};
diff --git a/src/OpenFOAM/algorithms/indexedOctree/volumeType.C b/src/OpenFOAM/algorithms/indexedOctree/volumeType.C
new file mode 100644
index 0000000000000000000000000000000000000000..95868f6c0c3f58e105b82ece76ecda59924daf50
--- /dev/null
+++ b/src/OpenFOAM/algorithms/indexedOctree/volumeType.C
@@ -0,0 +1,79 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "volumeType.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ template<>
+ const char* Foam::NamedEnum
+ <
+ Foam::volumeType,
+ 4
+ >::names[] =
+ {
+ "unknown",
+ "mixed",
+ "inside",
+ "outside"
+ };
+}
+
+const Foam::NamedEnum<Foam::volumeType, 4> Foam::volumeType::names;
+
+
+// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
+
+Foam::Istream& Foam::operator>>(Istream& is, volumeType& vt)
+{
+ // Read beginning of volumeType
+ is.readBegin("volumeType");
+
+ int type;
+ is >> type;
+
+ vt.t_ = static_cast<volumeType::type>(type);
+
+ // Read end of volumeType
+ is.readEnd("volumeType");
+
+ // Check state of Istream
+ is.check("operator>>(Istream&, volumeType&)");
+
+ return is;
+}
+
+
+Foam::Ostream& Foam::operator<<(Ostream& os, const volumeType& vt)
+{
+ os << static_cast<int>(vt.t_);
+
+ return os;
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/algorithms/indexedOctree/volumeType.H b/src/OpenFOAM/algorithms/indexedOctree/volumeType.H
new file mode 100644
index 0000000000000000000000000000000000000000..7bdfc043e6024e5b3a2417dbd4d8960b4927515a
--- /dev/null
+++ b/src/OpenFOAM/algorithms/indexedOctree/volumeType.H
@@ -0,0 +1,127 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::volumeType
+
+Description
+
+SourceFiles
+ volumeType.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef volumeType_H
+#define volumeType_H
+
+#include "NamedEnum.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward declaration of friend functions and operators
+
+class volumeType;
+Istream& operator>>(Istream& is, volumeType&);
+Ostream& operator<<(Ostream& os, const volumeType& C);
+
+
+/*---------------------------------------------------------------------------*\
+ Class volumeType Declaration
+\*---------------------------------------------------------------------------*/
+
+class volumeType
+{
+public:
+
+ //- Volume types
+ enum type
+ {
+ UNKNOWN = 0,
+ MIXED = 1,
+ INSIDE = 2,
+ OUTSIDE = 3
+ };
+
+
+private:
+
+ // Private data
+
+ //- Volume type
+ type t_;
+
+
+public:
+
+ // Static data
+
+ static const NamedEnum<volumeType, 4> names;
+
+
+ // Constructors
+
+ //- Construct null
+ volumeType()
+ :
+ t_(UNKNOWN)
+ {}
+
+ //- Construct from components
+ volumeType(type t)
+ :
+ t_(t)
+ {}
+
+
+ // Member Functions
+
+ operator type() const
+ {
+ return t_;
+ }
+
+
+ // IOstream operators
+
+ friend Istream& operator>>(Istream& is, volumeType& vt);
+ friend Ostream& operator<<(Ostream& os, const volumeType& vt);
+};
+
+
+//- Data associated with volumeType type are contiguous
+template<>
+inline bool contiguous<volumeType>() {return true;}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/containers/Lists/ListOps/ListOps.H b/src/OpenFOAM/containers/Lists/ListOps/ListOps.H
index f878a6f969deba67665afa3aab561d7a2426448f..c322c7856faae1accc1c36db8790131b1eb3e5c5 100644
--- a/src/OpenFOAM/containers/Lists/ListOps/ListOps.H
+++ b/src/OpenFOAM/containers/Lists/ListOps/ListOps.H
@@ -92,14 +92,23 @@ void inplaceMapKey(const labelUList& oldToNew, Container&);
template<class T>
void sortedOrder(const UList<T>&, labelList& order);
+template<class T, class Cmp>
+void sortedOrder(const UList<T>&, labelList& order, const Cmp& cmp);
+
//- Generate (sorted) indices corresponding to duplicate list values
template<class T>
void duplicateOrder(const UList<T>&, labelList& order);
+template<class T, class Cmp>
+void duplicateOrder(const UList<T>&, labelList& order, const Cmp& cmp);
+
//- Generate (sorted) indices corresponding to unique list values
template<class T>
void uniqueOrder(const UList<T>&, labelList& order);
+template<class T, class Cmp>
+void uniqueOrder(const UList<T>&, labelList& order, const Cmp& cmp);
+
//- Extract elements of List when select is a certain value.
// eg, to extract all selected elements:
// subset<bool, labelList>(selectedElems, true, lst);
diff --git a/src/OpenFOAM/containers/Lists/ListOps/ListOpsTemplates.C b/src/OpenFOAM/containers/Lists/ListOps/ListOpsTemplates.C
index ab6f9bf511e90a28abee4852d5e1d0ee4a40a53e..29ea65055121aedf96dc1ac262601134c97e6953 100644
--- a/src/OpenFOAM/containers/Lists/ListOps/ListOpsTemplates.C
+++ b/src/OpenFOAM/containers/Lists/ListOps/ListOpsTemplates.C
@@ -180,6 +180,18 @@ void Foam::sortedOrder
const UList<T>& lst,
labelList& order
)
+{
+ sortedOrder(lst, order, typename UList<T>::less(lst));
+}
+
+
+template<class T, class Cmp>
+void Foam::sortedOrder
+(
+ const UList<T>& lst,
+ labelList& order,
+ const Cmp& cmp
+)
{
// list lengths must be identical
if (order.size() != lst.size())
@@ -193,7 +205,7 @@ void Foam::sortedOrder
{
order[elemI] = elemI;
}
- Foam::stableSort(order, typename UList<T>::less(lst));
+ Foam::stableSort(order, cmp);
}
@@ -203,6 +215,18 @@ void Foam::duplicateOrder
const UList<T>& lst,
labelList& order
)
+{
+ duplicateOrder(lst, order, typename UList<T>::less(lst));
+}
+
+
+template<class T, class Cmp>
+void Foam::duplicateOrder
+(
+ const UList<T>& lst,
+ labelList& order,
+ const Cmp& cmp
+)
{
if (lst.size() < 2)
{
@@ -210,7 +234,7 @@ void Foam::duplicateOrder
return;
}
- sortedOrder(lst, order);
+ sortedOrder(lst, order, cmp);
label n = 0;
for (label i = 0; i < order.size() - 1; ++i)
@@ -231,7 +255,19 @@ void Foam::uniqueOrder
labelList& order
)
{
- sortedOrder(lst, order);
+ uniqueOrder(lst, order, typename UList<T>::less(lst));
+}
+
+
+template<class T, class Cmp>
+void Foam::uniqueOrder
+(
+ const UList<T>& lst,
+ labelList& order,
+ const Cmp& cmp
+)
+{
+ sortedOrder(lst, order, cmp);
if (order.size() > 1)
{
diff --git a/src/OpenFOAM/containers/Lists/SortableList/SortableList.C b/src/OpenFOAM/containers/Lists/SortableList/SortableList.C
index e7613ab1ab2c6fff2bdbb724e2f1c1018d631b83..d610a56a04cdf408325dcb8a8b1b805d53945563 100644
--- a/src/OpenFOAM/containers/Lists/SortableList/SortableList.C
+++ b/src/OpenFOAM/containers/Lists/SortableList/SortableList.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -23,26 +23,7 @@ License
\*---------------------------------------------------------------------------*/
-// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
-
-template<class T>
-void Foam::SortableList<T>::sortIndices(List<label>& order) const
-{
- // list lengths must be identical
- if (order.size() != this->size())
- {
- // avoid copying any elements, they are overwritten anyhow
- order.clear();
- order.setSize(this->size());
- }
-
- forAll(order, elemI)
- {
- order[elemI] = elemI;
- }
- Foam::stableSort(order, typename UList<T>::less(*this));
-}
-
+#include "ListOps.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
@@ -113,7 +94,7 @@ Foam::List<T>& Foam::SortableList<T>::shrink()
template<class T>
void Foam::SortableList<T>::sort()
{
- sortIndices(indices_);
+ sortedOrder(*this, indices_);
List<T> lst(this->size());
forAll(indices_, i)
@@ -128,13 +109,12 @@ void Foam::SortableList<T>::sort()
template<class T>
void Foam::SortableList<T>::reverseSort()
{
- sortIndices(indices_);
+ sortedOrder(*this, indices_, typename UList<T>::greater(*this));
List<T> lst(this->size());
- label endI = indices_.size();
forAll(indices_, i)
{
- lst[--endI] = this->operator[](indices_[i]);
+ lst[i] = this->operator[](indices_[i]);
}
List<T>::transfer(lst);
diff --git a/src/OpenFOAM/containers/Lists/SortableList/SortableList.H b/src/OpenFOAM/containers/Lists/SortableList/SortableList.H
index 1e62f39109e7952b9594441ba156db2bec227259..2d7b31a9e89ff315670d6608d23882528294d985 100644
--- a/src/OpenFOAM/containers/Lists/SortableList/SortableList.H
+++ b/src/OpenFOAM/containers/Lists/SortableList/SortableList.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -59,8 +59,6 @@ class SortableList
//- Original indices
labelList indices_;
- //- Resize, fill and sort the parameter according to the list values
- void sortIndices(List<label>&) const;
public:
diff --git a/src/OpenFOAM/containers/Lists/UList/UList.H b/src/OpenFOAM/containers/Lists/UList/UList.H
index 92557d3d352936354a00c01b8ca8cca78676cfaf..f61d6df7493ee2a7309aa2925b394002f236abd7 100644
--- a/src/OpenFOAM/containers/Lists/UList/UList.H
+++ b/src/OpenFOAM/containers/Lists/UList/UList.H
@@ -113,6 +113,24 @@ public:
}
};
+ //- Greater function class that can be used for sorting
+ class greater
+ {
+ const UList<T>& values_;
+
+ public:
+
+ greater(const UList<T>& values)
+ :
+ values_(values)
+ {}
+
+ bool operator()(const label a, const label b)
+ {
+ return values_[a] > values_[b];
+ }
+ };
+
// Constructors
diff --git a/src/OpenFOAM/meshes/meshShapes/triFace/triFace.H b/src/OpenFOAM/meshes/meshShapes/triFace/triFace.H
index 24ad5d3d1ccf920dd401a01b81edbe5354635d06..11f7d966bf7d1379fac7264751bf137131da1554 100644
--- a/src/OpenFOAM/meshes/meshShapes/triFace/triFace.H
+++ b/src/OpenFOAM/meshes/meshShapes/triFace/triFace.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -169,6 +169,16 @@ public:
const scalar tol = 0.0
) const;
+ inline pointHit intersection
+ (
+ const point& p,
+ const vector& q,
+ const point& ctr,
+ const pointField& points,
+ const intersection::algorithm alg,
+ const scalar tol = 0.0
+ ) const;
+
//- Return nearest point to face
inline pointHit nearestPoint
(
diff --git a/src/OpenFOAM/meshes/meshShapes/triFace/triFaceI.H b/src/OpenFOAM/meshes/meshShapes/triFace/triFaceI.H
index 6326c5cc6649246d7f67566ce51ed9ee1c36f3ed..b956780df4138ec0f0c4c8b7c52154bc292cdfd1 100644
--- a/src/OpenFOAM/meshes/meshShapes/triFace/triFaceI.H
+++ b/src/OpenFOAM/meshes/meshShapes/triFace/triFaceI.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -273,6 +273,20 @@ inline Foam::pointHit Foam::triFace::intersection
}
+inline Foam::pointHit Foam::triFace::intersection
+(
+ const point& p,
+ const vector& q,
+ const point& ctr,
+ const pointField& points,
+ const intersection::algorithm alg,
+ const scalar tol
+) const
+{
+ return intersection(p, q, points, alg, tol);
+}
+
+
inline Foam::pointHit Foam::triFace::nearestPoint
(
const point& p,
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchTools.H b/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchTools.H
index 36e792326dfaaf11e594396b1bae817e9953ec5f..e0800452a177d7cd26cdc1db1f9272e9e75a7d31 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchTools.H
+++ b/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchTools.H
@@ -54,6 +54,7 @@ namespace Foam
class polyMesh;
class PackedBoolList;
+class boundBox;
/*---------------------------------------------------------------------------*\
Class PatchTools Declaration
@@ -140,6 +141,21 @@ public:
labelList& faceMap
);
+ //-
+ template
+ <
+ class Face,
+ template<class> class FaceList,
+ class PointField,
+ class PointType
+ >
+ static void calcBounds
+ (
+ const PrimitivePatch<Face, FaceList, PointField, PointType>& p,
+ boundBox& bb,
+ label& nPoints
+ );
+
//- Return edge-face addressing sorted by angle around the edge.
// Orientation is anticlockwise looking from edge.vec(localPoints())
template
@@ -271,7 +287,6 @@ public:
List<Face>& mergedFaces,
labelList& pointMergeMap
);
-
};
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchToolsSearch.C b/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchToolsSearch.C
index 4fc6ca13c712f58c08906fcbc9272e80eba9932b..8ecd790b1ee6f1dd9cd5483b05d7f9c7e0922451 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchToolsSearch.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PatchTools/PatchToolsSearch.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -27,6 +27,8 @@ Description
\*---------------------------------------------------------------------------*/
#include "PatchTools.H"
+#include "PackedBoolList.H"
+#include "boundBox.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -222,4 +224,45 @@ Foam::PatchTools::subsetMap
}
+template
+<
+ class Face,
+ template<class> class FaceList,
+ class PointField,
+ class PointType
+>
+void Foam::PatchTools::calcBounds
+(
+ const PrimitivePatch<Face, FaceList, PointField, PointType>& p,
+ boundBox& bb,
+ label& nPoints
+)
+{
+ // Unfortunately nPoints constructs meshPoints() so do compact version
+ // ourselves
+ const PointField& points = p.points();
+
+ PackedBoolList pointIsUsed(points.size());
+
+ nPoints = 0;
+ bb = boundBox::invertedBox;
+
+ forAll(p, faceI)
+ {
+ const Face& f = p[faceI];
+
+ forAll(f, fp)
+ {
+ label pointI = f[fp];
+ if (pointIsUsed.set(pointI, 1u))
+ {
+ bb.min() = ::Foam::min(bb.min(), points[pointI]);
+ bb.max() = ::Foam::max(bb.max(), points[pointI]);
+ nPoints++;
+ }
+ }
+ }
+}
+
+
// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveShapes/objectHit/PointHit.H b/src/OpenFOAM/meshes/primitiveShapes/objectHit/PointHit.H
index 716bf88f734f9685622bf347b72b9c55ac293728..3e2e4dd6986d8cbfd8ef8df5edb6d0145fac3cdf 100644
--- a/src/OpenFOAM/meshes/primitiveShapes/objectHit/PointHit.H
+++ b/src/OpenFOAM/meshes/primitiveShapes/objectHit/PointHit.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -81,6 +81,15 @@ public:
// Constructors
+ //- Construct null
+ PointHit()
+ :
+ hit_(false),
+ hitPoint_(vector::zero),
+ distance_(GREAT),
+ eligibleMiss_(false)
+ {}
+
//- Construct from components
PointHit
(
diff --git a/src/OpenFOAM/meshes/primitiveShapes/plane/plane.H b/src/OpenFOAM/meshes/primitiveShapes/plane/plane.H
index 3d4e25211c10508a884ca36ee733b2946691c2f6..3a47e6bf2ea0c9e0b7fdacfcd39885163cedaede 100644
--- a/src/OpenFOAM/meshes/primitiveShapes/plane/plane.H
+++ b/src/OpenFOAM/meshes/primitiveShapes/plane/plane.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -26,7 +26,7 @@ Class
Description
Geometric class that creates a 2D plane and can return the intersection
- point between a line and the plane.
+ point between a line and the plane.
SourceFiles
plane.C
@@ -62,6 +62,14 @@ class plane
{
public:
+ //- Side of the plane
+ enum side
+ {
+ NORMAL,
+ FLIP
+ };
+
+
//- A direction and a reference point
class ray
{
@@ -182,6 +190,10 @@ public:
//- Return the cutting point between this plane and two other planes
point planePlaneIntersect(const plane&, const plane&) const;
+ //- Return the side of the plane that the point is on.
+ // If the point is on the plane, then returns NORMAL.
+ side sideOfPlane(const point& p) const;
+
//- Write to dictionary
void writeDict(Ostream&) const;
diff --git a/src/dynamicMesh/boundaryMesh/boundaryMesh.C b/src/dynamicMesh/boundaryMesh/boundaryMesh.C
index 1abf219dfb58dfd1306c8835ce0333c16952fb45..ad204d7c157fd9dad6c692a2b10447017a0fd924 100644
--- a/src/dynamicMesh/boundaryMesh/boundaryMesh.C
+++ b/src/dynamicMesh/boundaryMesh/boundaryMesh.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -924,6 +924,11 @@ Foam::labelList Foam::boundaryMesh::getNearest
bbMax.y() += 2*tol;
bbMax.z() += 2*tol;
+ const scalar planarTol =
+ indexedOctree<treeDataPrimitivePatch<uindirectPrimitivePatch> >::
+ perturbTol();
+
+
// Create the octrees
indexedOctree
<
@@ -933,7 +938,8 @@ Foam::labelList Foam::boundaryMesh::getNearest
treeDataPrimitivePatch<uindirectPrimitivePatch>
(
false, // cacheBb
- leftPatch
+ leftPatch,
+ planarTol
),
overallBb,
10, // maxLevel
@@ -948,7 +954,8 @@ Foam::labelList Foam::boundaryMesh::getNearest
treeDataPrimitivePatch<uindirectPrimitivePatch>
(
false, // cacheBb
- rightPatch
+ rightPatch,
+ planarTol
),
overallBb,
10, // maxLevel
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
index 1a699590953e2e94da8d701e8f989521169ebf47..f66c3acfa6e047eba3f95f563e74824a9a7cb3cd 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
@@ -241,7 +241,7 @@ surfaceInterpolationScheme<Type>::interpolate
//- Return the face-interpolate of the given cell field
-// with the given weigting factors
+// with the given weighting factors
template<class Type>
tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
surfaceInterpolationScheme<Type>::interpolate
diff --git a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
index a45a06a6b2c735d34770d51937a9d9f536ac75e3..f53968532c36a0b54c38a0a215bf92517b904793 100644
--- a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
+++ b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
@@ -31,6 +31,7 @@ License
#include "labelPair.H"
#include "searchableSurfacesQueries.H"
#include "UPtrList.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -1302,7 +1303,7 @@ void Foam::refinementSurfaces::findInside
if (allGeometry_[surfaces_[surfI]].hasVolumeType())
{
- List<searchableSurface::volumeType> volType;
+ List<volumeType> volType;
allGeometry_[surfaces_[surfI]].getVolumeType(pt, volType);
forAll(volType, pointI)
@@ -1312,11 +1313,11 @@ void Foam::refinementSurfaces::findInside
if
(
(
- volType[pointI] == triSurfaceMesh::INSIDE
+ volType[pointI] == volumeType::INSIDE
&& zoneInside_[surfI] == INSIDE
)
|| (
- volType[pointI] == triSurfaceMesh::OUTSIDE
+ volType[pointI] == volumeType::OUTSIDE
&& zoneInside_[surfI] == OUTSIDE
)
)
diff --git a/src/mesh/autoMesh/autoHexMesh/shellSurfaces/shellSurfaces.C b/src/mesh/autoMesh/autoHexMesh/shellSurfaces/shellSurfaces.C
index bc3289f18e13055f6189ffa357300e970e6ad499..a7a44f3fbfb388de08469b2a80f6616d4b31c685 100644
--- a/src/mesh/autoMesh/autoHexMesh/shellSurfaces/shellSurfaces.C
+++ b/src/mesh/autoMesh/autoHexMesh/shellSurfaces/shellSurfaces.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -31,6 +31,7 @@ License
#include "searchableSurfaces.H"
#include "orientedSurface.H"
#include "pointIndexHit.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -329,7 +330,7 @@ void Foam::shellSurfaces::findHigherLevel
candidateMap.setSize(candidateI);
// Do the expensive nearest test only for the candidate points.
- List<searchableSurface::volumeType> volType;
+ List<volumeType> volType;
allGeometry_[shells_[shellI]].getVolumeType(candidates, volType);
forAll(volType, i)
@@ -340,11 +341,11 @@ void Foam::shellSurfaces::findHigherLevel
(
(
modes_[shellI] == INSIDE
- && volType[i] == searchableSurface::INSIDE
+ && volType[i] == volumeType::INSIDE
)
|| (
modes_[shellI] == OUTSIDE
- && volType[i] == searchableSurface::OUTSIDE
+ && volType[i] == volumeType::OUTSIDE
)
)
{
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
index 920c78228682b85f6bcfbf1143e77eb5a85667b2..9c145081342bd8b36d0e3cf4a5a9d9b26ce89d03 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
@@ -210,7 +210,12 @@ void Foam::AMIMethod<SourcePatch, TargetPatch>::resetTree()
(
new indexedOctree<treeType>
(
- treeType(false, tgtPatch_),
+ treeType
+ (
+ false,
+ tgtPatch_,
+ indexedOctree<treeType>::perturbTol()
+ ),
bb, // overall search domain
8, // maxLevel
10, // leaf size
diff --git a/src/meshTools/Make/files b/src/meshTools/Make/files
index f11e6fc4712a237cb2bf553d83850e0d9b79d10e..91bbbc11fd2ea7df9c1a0fec04b4dc37ded5b369 100644
--- a/src/meshTools/Make/files
+++ b/src/meshTools/Make/files
@@ -160,6 +160,7 @@ $(intersectedSurface)/intersectedSurface.C
$(intersectedSurface)/edgeSurface.C
triSurface/triSurfaceSearch/triSurfaceSearch.C
+triSurface/triSurfaceSearch/triSurfaceRegionSearch.C
triSurface/triangleFuncs/triangleFuncs.C
triSurface/surfaceFeatures/surfaceFeatures.C
triSurface/triSurfaceTools/triSurfaceTools.C
diff --git a/src/meshTools/cellClassification/cellClassification.C b/src/meshTools/cellClassification/cellClassification.C
index 09e84ac93c0d1d032a240649ce14e8165d60da65..800d007f6fb66e7dc8c7400966a46724e3086503 100644
--- a/src/meshTools/cellClassification/cellClassification.C
+++ b/src/meshTools/cellClassification/cellClassification.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -34,6 +34,7 @@ License
#include "ListOps.H"
#include "meshTools.H"
#include "cpuTime.H"
+#include "triSurface.H"
#include "globalMeshData.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
diff --git a/src/meshTools/cellClassification/cellClassification.H b/src/meshTools/cellClassification/cellClassification.H
index aeb7a8102c59ff0e6bf3581237b0503319dd17b8..b040819fb2330ec5a0a73adde4c20a7d1e7ab5a0 100644
--- a/src/meshTools/cellClassification/cellClassification.H
+++ b/src/meshTools/cellClassification/cellClassification.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -107,9 +107,7 @@ namespace Foam
class triSurfaceSearch;
class meshSearch;
class polyMesh;
-class polyMesh;
class primitiveMesh;
-class triSurface;
/*---------------------------------------------------------------------------*\
Class cellClassification Declaration
diff --git a/src/meshTools/indexedOctree/treeDataEdge.C b/src/meshTools/indexedOctree/treeDataEdge.C
index 567fae04786363a0e404fbed217d38bc3cb30559..4d0b9b4d813624bb635424c125b1f61e7e3f0d07 100644
--- a/src/meshTools/indexedOctree/treeDataEdge.C
+++ b/src/meshTools/indexedOctree/treeDataEdge.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -96,6 +96,24 @@ Foam::treeDataEdge::treeDataEdge
}
+Foam::treeDataEdge::findNearestOp::findNearestOp
+(
+ const indexedOctree<treeDataEdge>& tree
+)
+:
+ tree_(tree)
+{}
+
+
+Foam::treeDataEdge::findIntersectOp::findIntersectOp
+(
+ const indexedOctree<treeDataEdge>& tree
+)
+:
+ tree_(tree)
+{}
+
+
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::pointField Foam::treeDataEdge::shapePoints() const
@@ -114,13 +132,13 @@ Foam::pointField Foam::treeDataEdge::shapePoints() const
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
-Foam::label Foam::treeDataEdge::getVolumeType
+Foam::volumeType Foam::treeDataEdge::getVolumeType
(
const indexedOctree<treeDataEdge>& oc,
const point& sample
) const
{
- return indexedOctree<treeDataEdge>::UNKNOWN;
+ return volumeType::UNKNOWN;
}
@@ -165,9 +183,7 @@ bool Foam::treeDataEdge::overlaps
}
-// Calculate nearest point to sample. Updates (if any) nearestDistSqr, minIndex,
-// nearestPoint.
-void Foam::treeDataEdge::findNearest
+void Foam::treeDataEdge::findNearestOp::operator()
(
const labelUList& indices,
const point& sample,
@@ -177,13 +193,15 @@ void Foam::treeDataEdge::findNearest
point& nearestPoint
) const
{
+ const treeDataEdge& shape = tree_.shapes();
+
forAll(indices, i)
{
const label index = indices[i];
- const edge& e = edges_[edgeLabels_[index]];
+ const edge& e = shape.edges()[shape.edgeLabels()[index]];
- pointHit nearHit = e.line(points_).nearestDist(sample);
+ pointHit nearHit = e.line(shape.points()).nearestDist(sample);
scalar distSqr = sqr(nearHit.distance());
@@ -197,9 +215,7 @@ void Foam::treeDataEdge::findNearest
}
-//- Calculates nearest (to line) point in shape.
-// Returns point and distance (squared)
-void Foam::treeDataEdge::findNearest
+void Foam::treeDataEdge::findNearestOp::operator()
(
const labelUList& indices,
const linePointRef& ln,
@@ -210,6 +226,8 @@ void Foam::treeDataEdge::findNearest
point& nearestPoint
) const
{
+ const treeDataEdge& shape = tree_.shapes();
+
// Best so far
scalar nearestDistSqr = magSqr(linePoint - nearestPoint);
@@ -217,13 +235,13 @@ void Foam::treeDataEdge::findNearest
{
const label index = indices[i];
- const edge& e = edges_[edgeLabels_[index]];
+ const edge& e = shape.edges()[shape.edgeLabels()[index]];
// Note: could do bb test ? Worthwhile?
// Nearest point on line
point ePoint, lnPt;
- scalar dist = e.line(points_).nearestDist(ln, ePoint, lnPt);
+ scalar dist = e.line(shape.points()).nearestDist(ln, ePoint, lnPt);
scalar distSqr = sqr(dist);
if (distSqr < nearestDistSqr)
@@ -252,4 +270,21 @@ void Foam::treeDataEdge::findNearest
}
+bool Foam::treeDataEdge::findIntersectOp::operator()
+(
+ const label index,
+ const point& start,
+ const point& end,
+ point& result
+) const
+{
+ notImplemented
+ (
+ "treeDataEdge::intersects(const label, const point&,"
+ "const point&, point&)"
+ );
+ return false;
+}
+
+
// ************************************************************************* //
diff --git a/src/meshTools/indexedOctree/treeDataEdge.H b/src/meshTools/indexedOctree/treeDataEdge.H
index 9a8deadb73c2883e9869bb80442ad2dd289e8b9a..12fd0d6aa10ba98c5db55977549b4775d58ca4d7 100644
--- a/src/meshTools/indexedOctree/treeDataEdge.H
+++ b/src/meshTools/indexedOctree/treeDataEdge.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -37,6 +37,7 @@ SourceFiles
#include "treeBoundBoxList.H"
#include "linePointRef.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -86,6 +87,58 @@ class treeDataEdge
public:
+
+ class findNearestOp
+ {
+ const indexedOctree<treeDataEdge>& tree_;
+
+ public:
+
+ findNearestOp(const indexedOctree<treeDataEdge>& tree);
+
+ void operator()
+ (
+ const labelUList& indices,
+ const point& sample,
+
+ scalar& nearestDistSqr,
+ label& minIndex,
+ point& nearestPoint
+ ) const;
+
+ void operator()
+ (
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+ ) const;
+ };
+
+
+ class findIntersectOp
+ {
+ const indexedOctree<treeDataEdge>& tree_;
+
+ public:
+
+ findIntersectOp(const indexedOctree<treeDataEdge>& tree);
+
+ //- Calculate intersection of triangle with ray. Sets result
+ // accordingly
+ bool operator()
+ (
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ ) const;
+ };
+
+
// Declare name of the class and its debug switch
ClassName("treeDataEdge");
@@ -116,6 +169,16 @@ public:
// Access
+ const edgeList& edges() const
+ {
+ return edges_;
+ }
+
+ const pointField& points() const
+ {
+ return points_;
+ }
+
const labelList& edgeLabels() const
{
return edgeLabels_;
@@ -135,7 +198,7 @@ public:
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
- label getVolumeType
+ volumeType getVolumeType
(
const indexedOctree<treeDataEdge>&,
const point&
@@ -155,50 +218,6 @@ public:
const point& centre,
const scalar radiusSqr
) const;
-
- //- Calculates nearest (to sample) point in shape.
- // Returns actual point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& nearestIndex,
- point& nearestPoint
- ) const;
-
- //- Calculates nearest (to line) point in shape.
- // Returns point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
- ) const;
-
- //- Calculate intersection of shape with ray. Sets result
- // accordingly
- bool intersects
- (
- const label index,
- const point& start,
- const point& end,
- point& result
- ) const
- {
- notImplemented
- (
- "treeDataEdge::intersects(const label, const point&,"
- "const point&, point&)"
- );
- return false;
- }
-
};
diff --git a/src/meshTools/indexedOctree/treeDataFace.C b/src/meshTools/indexedOctree/treeDataFace.C
index 9699416c13412d5e68cf62c0333500f5f9a16440..578a0de2ae056be937f97f226e847f10285ccbab 100644
--- a/src/meshTools/indexedOctree/treeDataFace.C
+++ b/src/meshTools/indexedOctree/treeDataFace.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -145,6 +145,24 @@ Foam::treeDataFace::treeDataFace
}
+Foam::treeDataFace::findNearestOp::findNearestOp
+(
+ const indexedOctree<treeDataFace>& tree
+)
+:
+ tree_(tree)
+{}
+
+
+Foam::treeDataFace::findIntersectOp::findIntersectOp
+(
+ const indexedOctree<treeDataFace>& tree
+)
+:
+ tree_(tree)
+{}
+
+
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::pointField Foam::treeDataFace::shapePoints() const
@@ -162,7 +180,7 @@ Foam::pointField Foam::treeDataFace::shapePoints() const
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
-Foam::label Foam::treeDataFace::getVolumeType
+Foam::volumeType Foam::treeDataFace::getVolumeType
(
const indexedOctree<treeDataFace>& oc,
const point& sample
@@ -415,7 +433,7 @@ Foam::label Foam::treeDataFace::getVolumeType
// - tolerances are wrong. (if e.g. face has zero area)
// - or (more likely) surface is not closed.
- return indexedOctree<treeDataFace>::UNKNOWN;
+ return volumeType::UNKNOWN;
}
@@ -477,9 +495,7 @@ bool Foam::treeDataFace::overlaps
}
-// Calculate nearest point to sample. Updates (if any) nearestDistSqr, minIndex,
-// nearestPoint.
-void Foam::treeDataFace::findNearest
+void Foam::treeDataFace::findNearestOp::operator()
(
const labelUList& indices,
const point& sample,
@@ -489,13 +505,15 @@ void Foam::treeDataFace::findNearest
point& nearestPoint
) const
{
+ const treeDataFace& shape = tree_.shapes();
+
forAll(indices, i)
{
const label index = indices[i];
- const face& f = mesh_.faces()[faceLabels_[index]];
+ const face& f = shape.mesh().faces()[shape.faceLabels()[index]];
- pointHit nearHit = f.nearestPoint(sample, mesh_.points());
+ pointHit nearHit = f.nearestPoint(sample, shape.mesh().points());
scalar distSqr = sqr(nearHit.distance());
if (distSqr < nearestDistSqr)
@@ -508,7 +526,33 @@ void Foam::treeDataFace::findNearest
}
-bool Foam::treeDataFace::intersects
+void Foam::treeDataFace::findNearestOp::operator()
+(
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+) const
+{
+ notImplemented
+ (
+ "treeDataFace::findNearestOp::operator()"
+ "("
+ " const labelUList&,"
+ " const linePointRef&,"
+ " treeBoundBox&,"
+ " label&,"
+ " point&,"
+ " point&"
+ ") const"
+ );
+}
+
+
+bool Foam::treeDataFace::findIntersectOp::operator()
(
const label index,
const point& start,
@@ -516,10 +560,12 @@ bool Foam::treeDataFace::intersects
point& intersectionPoint
) const
{
+ const treeDataFace& shape = tree_.shapes();
+
// Do quick rejection test
- if (cacheBb_)
+ if (shape.cacheBb_)
{
- const treeBoundBox& faceBb = bbs_[index];
+ const treeBoundBox& faceBb = shape.bbs_[index];
if ((faceBb.posBits(start) & faceBb.posBits(end)) != 0)
{
@@ -528,16 +574,16 @@ bool Foam::treeDataFace::intersects
}
}
- const label faceI = faceLabels_[index];
+ const label faceI = shape.faceLabels_[index];
const vector dir(end - start);
- pointHit inter = mesh_.faces()[faceI].intersection
+ pointHit inter = shape.mesh_.faces()[faceI].intersection
(
start,
dir,
- mesh_.faceCentres()[faceI],
- mesh_.points(),
+ shape.mesh_.faceCentres()[faceI],
+ shape.mesh_.points(),
intersection::HALF_RAY
);
diff --git a/src/meshTools/indexedOctree/treeDataFace.H b/src/meshTools/indexedOctree/treeDataFace.H
index 2b8ab258318cadc6d6e7a80168f2bc0938f47893..3867ba0206289f721d6bb037b70c251a515371fe 100644
--- a/src/meshTools/indexedOctree/treeDataFace.H
+++ b/src/meshTools/indexedOctree/treeDataFace.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -39,6 +39,8 @@ SourceFiles
#include "indexedOctree.H"
#include "treeBoundBoxList.H"
#include "PackedBoolList.H"
+#include "primitiveMesh.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -46,7 +48,7 @@ namespace Foam
{
// Forward declaration of classes
-class primitiveMesh;
+//class primitiveMesh;
//template<class Type> class indexedOctree;
class polyPatch;
@@ -90,6 +92,58 @@ class treeDataFace
public:
+
+ class findNearestOp
+ {
+ const indexedOctree<treeDataFace>& tree_;
+
+ public:
+
+ findNearestOp(const indexedOctree<treeDataFace>& tree);
+
+ void operator()
+ (
+ const labelUList& indices,
+ const point& sample,
+
+ scalar& nearestDistSqr,
+ label& minIndex,
+ point& nearestPoint
+ ) const;
+
+ void operator()
+ (
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+ ) const;
+ };
+
+
+ class findIntersectOp
+ {
+ const indexedOctree<treeDataFace>& tree_;
+
+ public:
+
+ findIntersectOp(const indexedOctree<treeDataFace>& tree);
+
+ //- Calculate intersection of triangle with ray. Sets result
+ // accordingly
+ bool operator()
+ (
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ ) const;
+ };
+
+
// Declare name of the class and its debug switch
ClassName("treeDataFace");
@@ -147,7 +201,7 @@ public:
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
- label getVolumeType
+ volumeType getVolumeType
(
const indexedOctree<treeDataFace>&,
const point&
@@ -159,49 +213,6 @@ public:
const label index,
const treeBoundBox& sampleBb
) const;
-
- //- Calculates nearest (to sample) point in shape.
- // Returns actual point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& nearestIndex,
- point& nearestPoint
- ) const;
-
- //- Calculates nearest (to line) point in shape.
- // Returns point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
- ) const
- {
- notImplemented
- (
- "treeDataFace::findNearest"
- "(const labelUList&, const linePointRef&, ..)"
- );
- }
-
- //- Calculate intersection of shape with ray. Sets result
- // accordingly
- bool intersects
- (
- const label index,
- const point& start,
- const point& end,
- point& result
- ) const;
-
};
diff --git a/src/meshTools/indexedOctree/treeDataPoint.C b/src/meshTools/indexedOctree/treeDataPoint.C
index 1d7f3d133ec2710e5ad64f1adb12723456db8cd0..fd7f4404b418cbab465ca6eb95b21505bdb8ddf5 100644
--- a/src/meshTools/indexedOctree/treeDataPoint.C
+++ b/src/meshTools/indexedOctree/treeDataPoint.C
@@ -57,6 +57,24 @@ Foam::treeDataPoint::treeDataPoint
{}
+Foam::treeDataPoint::findNearestOp::findNearestOp
+(
+ const indexedOctree<treeDataPoint>& tree
+)
+:
+ tree_(tree)
+{}
+
+
+Foam::treeDataPoint::findIntersectOp::findIntersectOp
+(
+ const indexedOctree<treeDataPoint>& tree
+)
+:
+ tree_(tree)
+{}
+
+
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::pointField Foam::treeDataPoint::shapePoints() const
@@ -74,13 +92,13 @@ Foam::pointField Foam::treeDataPoint::shapePoints() const
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
-Foam::label Foam::treeDataPoint::getVolumeType
+Foam::volumeType Foam::treeDataPoint::getVolumeType
(
const indexedOctree<treeDataPoint>& oc,
const point& sample
) const
{
- return indexedOctree<treeDataPoint>::UNKNOWN;
+ return volumeType::UNKNOWN;
}
@@ -115,9 +133,7 @@ bool Foam::treeDataPoint::overlaps
}
-// Calculate nearest point to sample. Updates (if any) nearestDistSqr, minIndex,
-// nearestPoint.
-void Foam::treeDataPoint::findNearest
+void Foam::treeDataPoint::findNearestOp::operator()
(
const labelUList& indices,
const point& sample,
@@ -127,12 +143,19 @@ void Foam::treeDataPoint::findNearest
point& nearestPoint
) const
{
+ const treeDataPoint& shape = tree_.shapes();
+
forAll(indices, i)
{
const label index = indices[i];
- label pointI = (useSubset_ ? pointLabels_[index] : index);
+ label pointI =
+ (
+ shape.useSubset()
+ ? shape.pointLabels()[index]
+ : index
+ );
- const point& pt = points_[pointI];
+ const point& pt = shape.points()[pointI];
scalar distSqr = magSqr(pt - sample);
@@ -146,9 +169,7 @@ void Foam::treeDataPoint::findNearest
}
-//- Calculates nearest (to line) point in shape.
-// Returns point and distance (squared)
-void Foam::treeDataPoint::findNearest
+void Foam::treeDataPoint::findNearestOp::operator()
(
const labelUList& indices,
const linePointRef& ln,
@@ -159,6 +180,8 @@ void Foam::treeDataPoint::findNearest
point& nearestPoint
) const
{
+ const treeDataPoint& shape = tree_.shapes();
+
// Best so far
scalar nearestDistSqr = GREAT;
if (minIndex >= 0)
@@ -169,9 +192,14 @@ void Foam::treeDataPoint::findNearest
forAll(indices, i)
{
const label index = indices[i];
- label pointI = (useSubset_ ? pointLabels_[index] : index);
+ label pointI =
+ (
+ shape.useSubset()
+ ? shape.pointLabels()[index]
+ : index
+ );
- const point& shapePt = points_[pointI];
+ const point& shapePt = shape.points()[pointI];
if (tightest.contains(shapePt))
{
@@ -206,4 +234,21 @@ void Foam::treeDataPoint::findNearest
}
+bool Foam::treeDataPoint::findIntersectOp::operator()
+(
+ const label index,
+ const point& start,
+ const point& end,
+ point& result
+) const
+{
+ notImplemented
+ (
+ "treeDataPoint::intersects(const label, const point&,"
+ "const point&, point&)"
+ );
+ return false;
+}
+
+
// ************************************************************************* //
diff --git a/src/meshTools/indexedOctree/treeDataPoint.H b/src/meshTools/indexedOctree/treeDataPoint.H
index f867a3274acf02f04c7786c623cf501b142ec152..020f60a11d14452cb21c0c678f7cca11486e1f51 100644
--- a/src/meshTools/indexedOctree/treeDataPoint.H
+++ b/src/meshTools/indexedOctree/treeDataPoint.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -43,6 +43,7 @@ SourceFiles
#include "pointField.H"
#include "treeBoundBox.H"
#include "linePointRef.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -69,6 +70,58 @@ class treeDataPoint
public:
+
+ class findNearestOp
+ {
+ const indexedOctree<treeDataPoint>& tree_;
+
+ public:
+
+ findNearestOp(const indexedOctree<treeDataPoint>& tree);
+
+ void operator()
+ (
+ const labelUList& indices,
+ const point& sample,
+
+ scalar& nearestDistSqr,
+ label& minIndex,
+ point& nearestPoint
+ ) const;
+
+ void operator()
+ (
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+ ) const;
+ };
+
+
+ class findIntersectOp
+ {
+ const indexedOctree<treeDataPoint>& tree_;
+
+ public:
+
+ findIntersectOp(const indexedOctree<treeDataPoint>& tree);
+
+ //- Calculate intersection of triangle with ray. Sets result
+ // accordingly
+ bool operator()
+ (
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ ) const;
+ };
+
+
// Declare name of the class and its debug switch
ClassName("treeDataPoint");
@@ -106,6 +159,11 @@ public:
return points_;
}
+ bool useSubset() const
+ {
+ return useSubset_;
+ }
+
//- Get representative point cloud for all shapes inside
// (one point per shape)
pointField shapePoints() const;
@@ -115,7 +173,7 @@ public:
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
- label getVolumeType
+ volumeType getVolumeType
(
const indexedOctree<treeDataPoint>&,
const point&
@@ -135,50 +193,6 @@ public:
const point& centre,
const scalar radiusSqr
) const;
-
- //- Calculates nearest (to sample) point in shape.
- // Returns actual point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& nearestIndex,
- point& nearestPoint
- ) const;
-
- //- Calculates nearest (to line) point in shape.
- // Returns point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
- ) const;
-
- //- Calculate intersection of shape with ray. Sets result
- // accordingly
- bool intersects
- (
- const label index,
- const point& start,
- const point& end,
- point& result
- ) const
- {
- notImplemented
- (
- "treeDataPoint::intersects(const label, const point&,"
- "const point&, point&)"
- );
- return false;
- }
-
};
diff --git a/src/meshTools/indexedOctree/treeDataPrimitivePatch.C b/src/meshTools/indexedOctree/treeDataPrimitivePatch.C
index 4cd29682323df7d68cb1259cdb9836ac7325749f..47d6f57d69e1a5d05a2b721289b94c2bd6c977ed 100644
--- a/src/meshTools/indexedOctree/treeDataPrimitivePatch.C
+++ b/src/meshTools/indexedOctree/treeDataPrimitivePatch.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -26,12 +26,8 @@ License
#include "treeDataPrimitivePatch.H"
#include "indexedOctree.H"
#include "triangleFuncs.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-template<class PatchType>
-Foam::scalar Foam::treeDataPrimitivePatch<PatchType>::tolSqr = sqr(1e-6);
-
+#include "triSurfaceTools.H"
+#include "triFace.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
@@ -70,6 +66,75 @@ void Foam::treeDataPrimitivePatch<PatchType>::update()
}
+template<class PatchType>
+bool Foam::treeDataPrimitivePatch<PatchType>::findIntersection
+(
+ const indexedOctree<treeDataPrimitivePatch<PatchType> >& tree,
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+)
+{
+ const treeDataPrimitivePatch<PatchType>& shape = tree.shapes();
+ const PatchType& patch = shape.patch();
+
+ const pointField& points = patch.points();
+ const typename PatchType::FaceType& f = patch[index];
+
+ // Do quick rejection test
+ if (shape.cacheBb_)
+ {
+ const treeBoundBox& faceBb = shape.bbs_[index];
+
+ if ((faceBb.posBits(start) & faceBb.posBits(end)) != 0)
+ {
+ // start and end in same block outside of faceBb.
+ return false;
+ }
+ }
+
+ const vector dir(end - start);
+ pointHit inter;
+
+ if (f.size() == 3)
+ {
+ inter = triPointRef
+ (
+ points[f[0]],
+ points[f[1]],
+ points[f[2]]
+ ).intersection(start, dir, intersection::HALF_RAY, shape.planarTol_);
+ }
+ else
+ {
+ const pointField& faceCentres = patch.faceCentres();
+
+ inter = f.intersection
+ (
+ start,
+ dir,
+ faceCentres[index],
+ points,
+ intersection::HALF_RAY
+ );
+ }
+
+ if (inter.hit() && inter.distance() <= 1)
+ {
+ // Note: no extra test on whether intersection is in front of us
+ // since using half_ray
+ intersectionPoint = inter.hitPoint();
+
+ return true;
+ }
+ else
+ {
+ return false;
+ }
+}
+
+
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from components
@@ -77,16 +142,50 @@ template<class PatchType>
Foam::treeDataPrimitivePatch<PatchType>::treeDataPrimitivePatch
(
const bool cacheBb,
- const PatchType& patch
+ const PatchType& patch,
+ const scalar planarTol
)
:
patch_(patch),
- cacheBb_(cacheBb)
+ cacheBb_(cacheBb),
+ planarTol_(planarTol)
{
update();
}
+template<class PatchType>
+Foam::treeDataPrimitivePatch<PatchType>::findNearestOp::findNearestOp
+(
+ const indexedOctree<treeDataPrimitivePatch<PatchType> >& tree
+)
+:
+ tree_(tree)
+{}
+
+
+template<class PatchType>
+Foam::treeDataPrimitivePatch<PatchType>::findIntersectOp::findIntersectOp
+(
+ const indexedOctree<treeDataPrimitivePatch<PatchType> >& tree
+)
+:
+ tree_(tree)
+{}
+
+
+template<class PatchType>
+Foam::treeDataPrimitivePatch<PatchType>::findAllIntersectOp::findAllIntersectOp
+(
+ const indexedOctree<treeDataPrimitivePatch<PatchType> >& tree,
+ DynamicList<label>& shapeMask
+)
+:
+ tree_(tree),
+ shapeMask_(shapeMask)
+{}
+
+
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class PatchType>
@@ -106,7 +205,7 @@ Foam::pointField Foam::treeDataPrimitivePatch<PatchType>::shapePoints() const
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
template<class PatchType>
-Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
+Foam::volumeType Foam::treeDataPrimitivePatch<PatchType>::getVolumeType
(
const indexedOctree<treeDataPrimitivePatch<PatchType> >& oc,
const point& sample
@@ -136,7 +235,6 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
<< abort(FatalError);
}
-
// Get actual intersection point on face
label faceI = info.index();
@@ -147,7 +245,7 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
}
const pointField& points = patch_.localPoints();
- const face& f = patch_.localFaces()[faceI];
+ const typename PatchType::FaceType& f = patch_.localFaces()[faceI];
// Retest to classify where on face info is. Note: could be improved. We
// already have point.
@@ -188,9 +286,10 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
const scalar typDimSqr = mag(area) + VSMALL;
+
forAll(f, fp)
{
- if ((magSqr(points[f[fp]] - curPt)/typDimSqr) < tolSqr)
+ if ((magSqr(points[f[fp]] - curPt)/typDimSqr) < planarTol_)
{
// Face intersection point equals face vertex fp
@@ -207,7 +306,7 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
const point fc(f.centre(points));
- if ((magSqr(fc - curPt)/typDimSqr) < tolSqr)
+ if ((magSqr(fc - curPt)/typDimSqr) < planarTol_)
{
// Face intersection point equals face centre. Normal at face centre
// is already average of face normals
@@ -240,7 +339,7 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
pointHit edgeHit = e.line(points).nearestDist(sample);
- if ((magSqr(edgeHit.rawPoint() - curPt)/typDimSqr) < tolSqr)
+ if ((magSqr(edgeHit.rawPoint() - curPt)/typDimSqr) < planarTol_)
{
// Face intersection point lies on edge e
@@ -285,7 +384,7 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
fc
).nearestDist(sample);
- if ((magSqr(edgeHit.rawPoint() - curPt)/typDimSqr) < tolSqr)
+ if ((magSqr(edgeHit.rawPoint() - curPt)/typDimSqr) < planarTol_)
{
// Face intersection point lies on edge between two face triangles
@@ -336,7 +435,7 @@ Foam::label Foam::treeDataPrimitivePatch<PatchType>:: getVolumeType
// - tolerances are wrong. (if e.g. face has zero area)
// - or (more likely) surface is not closed.
- return indexedOctree<treeDataPrimitivePatch>::UNKNOWN;
+ return volumeType::UNKNOWN;
}
@@ -368,7 +467,7 @@ bool Foam::treeDataPrimitivePatch<PatchType>::overlaps
// 2. Check if one or more face points inside
const pointField& points = patch_.points();
- const face& f = patch_[index];
+ const typename PatchType::FaceType& f = patch_[index];
if (cubeBb.containsAny(points, f))
{
@@ -379,21 +478,35 @@ bool Foam::treeDataPrimitivePatch<PatchType>::overlaps
// go through cube. Use triangle-bounding box intersection.
const point fc = f.centre(points);
- forAll(f, fp)
+ if (f.size() == 3)
{
- bool triIntersects = triangleFuncs::intersectBb
+ return triangleFuncs::intersectBb
(
- points[f[fp]],
- points[f[f.fcIndex(fp)]],
- fc,
+ points[f[0]],
+ points[f[1]],
+ points[f[2]],
cubeBb
);
-
- if (triIntersects)
+ }
+ else
+ {
+ forAll(f, fp)
{
- return true;
+ bool triIntersects = triangleFuncs::intersectBb
+ (
+ points[f[fp]],
+ points[f[f.fcIndex(fp)]],
+ fc,
+ cubeBb
+ );
+
+ if (triIntersects)
+ {
+ return true;
+ }
}
}
+
return false;
}
@@ -439,10 +552,8 @@ bool Foam::treeDataPrimitivePatch<PatchType>::overlaps
}
-// Calculate nearest point to sample. Updates (if any) nearestDistSqr, minIndex,
-// nearestPoint.
template<class PatchType>
-void Foam::treeDataPrimitivePatch<PatchType>::findNearest
+void Foam::treeDataPrimitivePatch<PatchType>::findNearestOp::operator()
(
const labelUList& indices,
const point& sample,
@@ -452,13 +563,15 @@ void Foam::treeDataPrimitivePatch<PatchType>::findNearest
point& nearestPoint
) const
{
- const pointField& points = patch_.points();
+ const treeDataPrimitivePatch<PatchType>& shape = tree_.shapes();
+ const PatchType& patch = shape.patch();
+
+ const pointField& points = patch.points();
forAll(indices, i)
{
const label index = indices[i];
-
- const face& f = patch_[index];
+ const typename PatchType::FaceType& f = patch[index];
pointHit nearHit = f.nearestPoint(sample, points);
scalar distSqr = sqr(nearHit.distance());
@@ -474,7 +587,34 @@ void Foam::treeDataPrimitivePatch<PatchType>::findNearest
template<class PatchType>
-bool Foam::treeDataPrimitivePatch<PatchType>::intersects
+void Foam::treeDataPrimitivePatch<PatchType>::findNearestOp::operator()
+(
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+) const
+{
+ notImplemented
+ (
+ "treeDataPrimitivePatch<PatchType>::findNearestOp::operator()"
+ "("
+ " const labelUList&,"
+ " const linePointRef&,"
+ " treeBoundBox&,"
+ " label&,"
+ " point&,"
+ " point&"
+ ") const"
+ );
+}
+
+
+template<class PatchType>
+bool Foam::treeDataPrimitivePatch<PatchType>::findIntersectOp::operator()
(
const label index,
const point& start,
@@ -482,43 +622,25 @@ bool Foam::treeDataPrimitivePatch<PatchType>::intersects
point& intersectionPoint
) const
{
- // Do quick rejection test
- if (cacheBb_)
- {
- const treeBoundBox& faceBb = bbs_[index];
-
- if ((faceBb.posBits(start) & faceBb.posBits(end)) != 0)
- {
- // start and end in same block outside of faceBb.
- return false;
- }
- }
-
- const pointField& points = patch_.points();
- const face& f = patch_[index];
- const point fc = f.centre(points);
- const vector dir(end - start);
+ return findIntersection(tree_, index, start, end, intersectionPoint);
+}
- pointHit inter = patch_[index].intersection
- (
- start,
- dir,
- fc,
- points,
- intersection::HALF_RAY
- );
- if (inter.hit() && inter.distance() <= 1)
- {
- // Note: no extra test on whether intersection is in front of us
- // since using half_ray
- intersectionPoint = inter.hitPoint();
- return true;
- }
- else
+template<class PatchType>
+bool Foam::treeDataPrimitivePatch<PatchType>::findAllIntersectOp::operator()
+(
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+) const
+{
+ if (!shapeMask_.empty() && findIndex(shapeMask_, index) != -1)
{
return false;
}
+
+ return findIntersection(tree_, index, start, end, intersectionPoint);
}
diff --git a/src/meshTools/indexedOctree/treeDataPrimitivePatch.H b/src/meshTools/indexedOctree/treeDataPrimitivePatch.H
index bbb3c749562c49d99ce168a8e76cbc0fc2988a01..9d43d2ba85154dcec67dbab742a7d7a7bc89b19c 100644
--- a/src/meshTools/indexedOctree/treeDataPrimitivePatch.H
+++ b/src/meshTools/indexedOctree/treeDataPrimitivePatch.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -35,9 +35,8 @@ SourceFiles
#ifndef treeDataPrimitivePatch_H
#define treeDataPrimitivePatch_H
-#include "PrimitivePatch.H"
-//#include "indexedOctree.H"
#include "treeBoundBoxList.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -64,11 +63,6 @@ class treeDataPrimitivePatch
:
public treeDataPrimitivePatchName
{
- // Static data
-
- //- tolerance on linear dimensions
- static scalar tolSqr;
-
// Private data
//- Underlying geometry
@@ -77,6 +71,9 @@ class treeDataPrimitivePatch
//- Whether to precalculate and store face bounding box
const bool cacheBb_;
+ //- Tolerance to use for intersection tests
+ const scalar planarTol_;
+
//- face bounding boxes (valid only if cacheBb_)
treeBoundBoxList bbs_;
@@ -89,15 +86,107 @@ class treeDataPrimitivePatch
//- Initialise all member data
void update();
+ //- Find intersection of line with shapes
+ static bool findIntersection
+ (
+ const indexedOctree<treeDataPrimitivePatch<PatchType> >& tree,
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ );
+
+
public:
+
+ class findNearestOp
+ {
+ const indexedOctree<treeDataPrimitivePatch>& tree_;
+
+ public:
+
+ findNearestOp(const indexedOctree<treeDataPrimitivePatch>& tree);
+
+ void operator()
+ (
+ const labelUList& indices,
+ const point& sample,
+
+ scalar& nearestDistSqr,
+ label& minIndex,
+ point& nearestPoint
+ ) const;
+
+ //- Calculates nearest (to line) point in shape.
+ // Returns point and distance (squared)
+ void operator()
+ (
+ const labelUList& indices,
+ const linePointRef& ln,
+
+ treeBoundBox& tightest,
+ label& minIndex,
+ point& linePoint,
+ point& nearestPoint
+ ) const;
+ };
+
+
+ class findIntersectOp
+ {
+ const indexedOctree<treeDataPrimitivePatch>& tree_;
+
+ public:
+
+ findIntersectOp(const indexedOctree<treeDataPrimitivePatch>& tree);
+
+ //- Calculate intersection of triangle with ray. Sets result
+ // accordingly
+ bool operator()
+ (
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ ) const;
+ };
+
+
+ class findAllIntersectOp
+ {
+ const indexedOctree<treeDataPrimitivePatch>& tree_;
+
+ DynamicList<label>& shapeMask_;
+
+ public:
+
+ findAllIntersectOp
+ (
+ const indexedOctree<treeDataPrimitivePatch>& tree,
+ DynamicList<label>& shapeMask
+ );
+
+ //- Calculate intersection of triangle with ray. Sets result
+ // accordingly
+ bool operator()
+ (
+ const label index,
+ const point& start,
+ const point& end,
+ point& intersectionPoint
+ ) const;
+ };
+
+
// Constructors
//- Construct from patch.
treeDataPrimitivePatch
(
const bool cacheBb,
- const PatchType&
+ const PatchType&,
+ const scalar planarTol
);
@@ -125,7 +214,7 @@ public:
//- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
// Only makes sense for closed surfaces.
- label getVolumeType
+ volumeType getVolumeType
(
const indexedOctree<treeDataPrimitivePatch<PatchType> >&,
const point&
@@ -145,49 +234,6 @@ public:
const point& centre,
const scalar radiusSqr
) const;
-
- //- Calculates nearest (to sample) point in shape.
- // Returns actual point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& nearestIndex,
- point& nearestPoint
- ) const;
-
- //- Calculates nearest (to line) point in shape.
- // Returns point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
- ) const
- {
- notImplemented
- (
- "treeDataPrimitivePatch::findNearest"
- "(const labelUList&, const linePointRef&, ..)"
- );
- }
-
- //- Calculate intersection of shape with ray. Sets result
- // accordingly
- bool intersects
- (
- const label index,
- const point& start,
- const point& end,
- point& result
- ) const;
-
};
diff --git a/src/meshTools/indexedOctree/treeDataTriSurface.C b/src/meshTools/indexedOctree/treeDataTriSurface.C
index 0aa567f68b0ce38a80b6bf4576dc25f25accaa6d..ef00434eb28b28116232681c4d37d1eb152e43e5 100644
--- a/src/meshTools/indexedOctree/treeDataTriSurface.C
+++ b/src/meshTools/indexedOctree/treeDataTriSurface.C
@@ -25,491 +25,58 @@ License
#include "treeDataTriSurface.H"
#include "triSurfaceTools.H"
-#include "triangleFuncs.H"
-#include "indexedOctree.H"
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-namespace Foam
-{
-defineTypeNameAndDebug(treeDataTriSurface, 0);
-}
-
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-// // Fast distance to triangle calculation. From
-// // "Distance Between Point and Triangle in 3D"
-// // David Eberly, Magic Software Inc. Aug. 2003.
-// // Works on function Q giving distance to point and tries to minimize this.
-// Foam::scalar Foam::treeDataTriSurface::nearestCoords
-// (
-// const point& base,
-// const point& E0,
-// const point& E1,
-// const scalar a,
-// const scalar b,
-// const scalar c,
-// const point& P,
-// scalar& s,
-// scalar& t
-// )
-// {
-// // distance vector
-// const vector D(base - P);
-
-// // Precalculate distance factors.
-// const scalar d = E0 & D;
-// const scalar e = E1 & D;
-
-// // Do classification
-// const scalar det = a*c - b*b;
-
-// s = b*e - c*d;
-// t = b*d - a*e;
-
-// if (s+t < det)
-// {
-// if (s < 0)
-// {
-// if (t < 0)
-// {
-// //region 4
-// if (e > 0)
-// {
-// //min on edge t = 0
-// t = 0;
-// s = (d >= 0 ? 0 : (-d >= a ? 1 : -d/a));
-// }
-// else
-// {
-// //min on edge s=0
-// s = 0;
-// t = (e >= 0 ? 0 : (-e >= c ? 1 : -e/c));
-// }
-// }
-// else
-// {
-// //region 3. Min on edge s = 0
-// s = 0;
-// t = (e >= 0 ? 0 : (-e >= c ? 1 : -e/c));
-// }
-// }
-// else if (t < 0)
-// {
-// //region 5
-// t = 0;
-// s = (d >= 0 ? 0 : (-d >= a ? 1 : -d/a));
-// }
-// else
-// {
-// //region 0
-// const scalar invDet = 1/det;
-// s *= invDet;
-// t *= invDet;
-// }
-// }
-// else
-// {
-// if (s < 0)
-// {
-// //region 2
-// const scalar tmp0 = b + d;
-// const scalar tmp1 = c + e;
-// if (tmp1 > tmp0)
-// {
-// //min on edge s+t=1
-// const scalar numer = tmp1 - tmp0;
-// const scalar denom = a-2*b+c;
-// s = (numer >= denom ? 1 : numer/denom);
-// t = 1 - s;
-// }
-// else
-// {
-// //min on edge s=0
-// s = 0;
-// t = (tmp1 <= 0 ? 1 : (e >= 0 ? 0 : - e/c));
-// }
-// }
-// else if (t < 0)
-// {
-// //region 6
-// const scalar tmp0 = b + d;
-// const scalar tmp1 = c + e;
-// if (tmp1 > tmp0)
-// {
-// //min on edge s+t=1
-// const scalar numer = tmp1 - tmp0;
-// const scalar denom = a-2*b+c;
-// s = (numer >= denom ? 1 : numer/denom);
-// t = 1 - s;
-// }
-// else
-// {
-// //min on edge t=0
-// t = 0;
-// s = (tmp1 <= 0 ? 1 : (d >= 0 ? 0 : - d/a));
-// }
-// }
-// else
-// {
-// //region 1
-// const scalar numer = c+e-(b+d);
-// if (numer <= 0)
-// {
-// s = 0;
-// }
-// else
-// {
-// const scalar denom = a-2*b+c;
-// s = (numer >= denom ? 1 : numer/denom);
-// }
-// }
-// t = 1 - s;
-// }
-
-
-// // Calculate distance.
-// // Note: abs should not be needed but truncation error causes problems
-// // with points very close to one of the triangle vertices.
-// // (seen up to -9e-15). Alternatively add some small value.
-
-// const scalar f = D & D;
-// return Foam::mag(a*s*s + 2*b*s*t + c*t*t + 2*d*s + 2*e*t + f);
-// }
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-// Construct from components
-Foam::treeDataTriSurface::treeDataTriSurface
+template<>
+Foam::volumeType Foam::treeDataPrimitivePatch<Foam::triSurface>::getVolumeType
(
- const triSurface& surface,
- const scalar planarTol
-)
-:
- surface_(surface),
- planarTol_(planarTol)
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-Foam::pointField Foam::treeDataTriSurface::shapePoints() const
-{
- const pointField& points = surface_.points();
-
- pointField centres(surface_.size());
-
- forAll(surface_, triI)
- {
- centres[triI] = surface_[triI].centre(points);
- }
- return centres;
-}
-
-
-//- Get type of sample (inside/outside/mixed) w.r.t. surface.
-Foam::label Foam::treeDataTriSurface::getVolumeType
-(
- const indexedOctree<treeDataTriSurface>& tree,
+ const indexedOctree<treeDataPrimitivePatch<triSurface> >& oc,
const point& sample
) const
{
- // Find nearest point
- const treeBoundBox& treeBb = tree.bb();
+ // Find nearest face to sample
+ pointIndexHit info = oc.findNearest(sample, sqr(GREAT));
- pointIndexHit pHit = tree.findNearest
- (
- sample,
- max
- (
- Foam::sqr(GREAT),
- Foam::magSqr(treeBb.span())
- )
- );
-
- if (!pHit.hit())
+ if (info.index() == -1)
{
- FatalErrorIn("treeDataTriSurface::getVolumeType(..)")
- << "treeBb:" << treeBb
- << " sample:" << sample
- << " pHit:" << pHit
+ FatalErrorIn
+ (
+ "treeDataPrimitivePatch::getSampleType"
+ "(indexedOctree<treeDataPrimitivePatch>&, const point&)"
+ ) << "Could not find " << sample << " in octree."
<< abort(FatalError);
}
+ // Get actual intersection point on face
+ label faceI = info.index();
+
triSurfaceTools::sideType t = triSurfaceTools::surfaceSide
(
- surface_,
+ patch_,
sample,
- pHit.index()
+ faceI
);
if (t == triSurfaceTools::UNKNOWN)
{
- return indexedOctree<treeDataTriSurface>::UNKNOWN;
+ return volumeType::UNKNOWN;
}
else if (t == triSurfaceTools::INSIDE)
{
- return indexedOctree<treeDataTriSurface>::INSIDE;
+ return volumeType::INSIDE;
}
else if (t == triSurfaceTools::OUTSIDE)
{
- return indexedOctree<treeDataTriSurface>::OUTSIDE;
+ return volumeType::OUTSIDE;
}
else
{
- FatalErrorIn("treeDataTriSurface::getVolumeType(..)")
+ FatalErrorIn("treeDataPrimitivePatch<PatchType>::getVolumeType(..)")
<< "problem" << abort(FatalError);
- return indexedOctree<treeDataTriSurface>::UNKNOWN;
- }
-}
-
-
-// Check if any point on triangle is inside cubeBb.
-bool Foam::treeDataTriSurface::overlaps
-(
- const label index,
- const treeBoundBox& cubeBb
-) const
-{
- const pointField& points = surface_.points();
- const labelledTri& f = surface_[index];
-
- treeBoundBox triBb(points, surface_[index]);
-
- //- For testing: robust one
- //return cubeBb.overlaps(triBb);
-
- //- Exact test of triangle intersecting bb
-
- // Quick rejection. If whole bounding box of tri is outside cubeBb then
- // there will be no intersection.
- if (!cubeBb.overlaps(triBb))
- {
- return false;
- }
-
- // Check if one or more triangle point inside
- if (cubeBb.containsAny(points, f))
- {
- return true;
- }
-
- // Triangle points
- const point& p0 = points[f[0]];
- const point& p1 = points[f[1]];
- const point& p2 = points[f[2]];
-
-
- // Now we have the difficult case: all points are outside but connecting
- // edges might go through cube. Use fast intersection of bounding box.
-
- //return triangleFuncs::intersectBbExact(p0, p1, p2, cubeBb);
- return triangleFuncs::intersectBb(p0, p1, p2, cubeBb);
-}
-
-
-// Calculate nearest point to sample. Updates (if any) nearestDistSqr, minIndex,
-// nearestPoint.
-void Foam::treeDataTriSurface::findNearest
-(
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& minIndex,
- point& nearestPoint
-) const
-{
- const pointField& points = surface_.points();
-
- forAll(indices, i)
- {
- label index = indices[i];
- const triSurface::FaceType& f = surface_[index];
-
- ////- Possible optimization: do quick rejection of triangle if bounding
- //// sphere does not intersect triangle bounding box. From simplistic
- //// test was not found to speed up things.
- //
- //// Triangle bounding box.
- //point triBbMin = min(p0, min(p1, p2));
- //point triBbMax = max(p0, max(p1, p2));
- //
- //if
- //(
- // indexedOctree<treeDataTriSurface>::intersects
- // (
- // triBbMin,
- // triBbMax,
- // nearestDistSqr,
- // sample
- // )
- //)
- {
- // // Get spanning vectors of triangle
- // vector base(p1);
- // vector E0(p0 - p1);
- // vector E1(p2 - p1);
-
- // scalar a(E0& E0);
- // scalar b(E0& E1);
- // scalar c(E1& E1);
-
- // // Get nearest point in s,t coordinates (s is along E0, t
- // // is along E1)
- // scalar s;
- // scalar t;
-
- // scalar distSqr = nearestCoords
- // (
- // base,
- // E0,
- // E1,
- // a,
- // b,
- // c,
- // sample,
-
- // s,
- // t
- // );
-
- pointHit pHit = f.nearestPoint(sample, points);
-
- scalar distSqr = sqr(pHit.distance());
-
- if (distSqr < nearestDistSqr)
- {
- nearestDistSqr = distSqr;
- minIndex = index;
- nearestPoint = pHit.rawPoint();
- }
- }
+ return volumeType::UNKNOWN;
}
}
-// Calculate nearest point to line. Updates (if any) nearestDistSqr, minIndex,
-// nearestPoint.
-void Foam::treeDataTriSurface::findNearest
-(
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
-) const
-{
- // Best so far
- scalar nearestDistSqr = VGREAT;
- if (minIndex >= 0)
- {
- nearestDistSqr = magSqr(linePoint - nearestPoint);
- }
-
- const pointField& points = surface_.points();
-
- forAll(indices, i)
- {
- label index = indices[i];
- const triSurface::FaceType& f = surface_[index];
-
- triPointRef tri(f.tri(points));
-
- pointHit lineInfo(point::max);
- pointHit pHit = tri.nearestPoint(ln, lineInfo);
-
- scalar distSqr = sqr(pHit.distance());
-
- if (distSqr < nearestDistSqr)
- {
- nearestDistSqr = distSqr;
- minIndex = index;
- linePoint = lineInfo.rawPoint();
- nearestPoint = pHit.rawPoint();
-
- {
- point& minPt = tightest.min();
- minPt = min(ln.start(), ln.end());
- minPt.x() -= pHit.distance();
- minPt.y() -= pHit.distance();
- minPt.z() -= pHit.distance();
- }
- {
- point& maxPt = tightest.max();
- maxPt = max(ln.start(), ln.end());
- maxPt.x() += pHit.distance();
- maxPt.y() += pHit.distance();
- maxPt.z() += pHit.distance();
- }
- }
- }
-}
-
-
-bool Foam::treeDataTriSurface::intersects
-(
- const label index,
- const point& start,
- const point& end,
- point& intersectionPoint
-) const
-{
- const pointField& points = surface_.points();
-
- const triSurface::FaceType& f = surface_[index];
-
- // Do quick rejection test
- treeBoundBox triBb(points, f);
-
- const direction startBits(triBb.posBits(start));
- const direction endBits(triBb.posBits(end));
-
- if ((startBits & endBits) != 0)
- {
- // start and end in same block outside of triBb.
- return false;
- }
-
- const vector dir(end - start);
-
- // Use relative tolerance (from octree) to determine intersection.
- pointHit inter = f.intersection
- (
- start,
- dir,
- points,
- intersection::HALF_RAY,
- planarTol_
- );
-
- if (inter.hit() && inter.distance() <= 1)
- {
- // Note: no extra test on whether intersection is in front of us
- // since using half_ray.
- intersectionPoint = inter.hitPoint();
-
- return true;
- }
- else
- {
- return false;
- }
-
-
- //- Using exact intersections
- //pointHit info = f.tri(points).intersectionExact(start, end);
- //
- //if (info.hit())
- //{
- // intersectionPoint = info.hitPoint();
- //}
- //return info.hit();
-}
-
-
// ************************************************************************* //
diff --git a/src/meshTools/indexedOctree/treeDataTriSurface.H b/src/meshTools/indexedOctree/treeDataTriSurface.H
index 9184294c0be0e211f5f5cd799727848ebaa959c9..fbc3d00b129311f4d390daea58204ec605fe219e 100644
--- a/src/meshTools/indexedOctree/treeDataTriSurface.H
+++ b/src/meshTools/indexedOctree/treeDataTriSurface.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -25,7 +25,7 @@ Class
Foam::treeDataTriSurface
Description
- Encapsulates data for (indexedOc)tree searches on triSurface.
+ Encapsulates data for (indexedOc)tree searches on a triSurface.
SourceFiles
treeDataTriSurface.C
@@ -35,139 +35,27 @@ SourceFiles
#ifndef treeDataTriSurface_H
#define treeDataTriSurface_H
+#include "treeDataPrimitivePatch.H"
#include "triSurface.H"
+#include "point.H"
+#include "indexedOctree.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
+ typedef treeDataPrimitivePatch<triSurface> treeDataTriSurface;
-// Forward declaration of classes
-class treeBoundBox;
-class treeDataTriSurface;
-template<class Type> class indexedOctree;
+ //- Template specialisation of getVolumeType for treeDataTriSurface
+ template<>
+ volumeType treeDataPrimitivePatch<triSurface>::getVolumeType
+ (
+ const indexedOctree<treeDataPrimitivePatch<triSurface> >& oc,
+ const point& sample
+ ) const;
+}
-/*---------------------------------------------------------------------------*\
- Class treeDataTriSurface Declaration
-\*---------------------------------------------------------------------------*/
-
-class treeDataTriSurface
-{
- // Private data
-
- //- Reference to triSurface
- const triSurface& surface_;
-
- //- Tolerance to use for intersection tests
- const scalar planarTol_;
-
- // Private Member Functions
-
- // //- fast triangle nearest point calculation. Returns point in E0, E1
- // // coordinate system: base + s*E0 + t*E1
- // static scalar nearestCoords
- // (
- // const point& base,
- // const point& E0,
- // const point& E1,
- // const scalar a,
- // const scalar b,
- // const scalar c,
- // const point& P,
- // scalar& s,
- // scalar& t
- // );
-
-public:
-
- // Declare name of the class and its debug switch
- ClassName("treeDataTriSurface");
-
-
- // Constructors
-
- //- Construct from triSurface and tolerance for intersection
- // tests. Holds reference.
- treeDataTriSurface(const triSurface&, const scalar planarTol);
-
-
- // Member Functions
-
- // Access
-
- inline const triSurface& surface() const
- {
- return surface_;
- }
-
- inline label size() const
- {
- return surface_.size();
- }
-
- //- Get representative point cloud for all shapes inside
- // (one point per shape)
- pointField shapePoints() const;
-
-
- // Search
-
- //- Get type (inside,outside,mixed,unknown) of point w.r.t. surface.
- // Only makes sense for closed surfaces.
- label getVolumeType
- (
- const indexedOctree<treeDataTriSurface>&,
- const point&
- ) const;
-
- //- Does (bb of) shape at index overlap bb
- bool overlaps
- (
- const label index,
- const treeBoundBox& sampleBb
- ) const;
-
- //- Calculates nearest (to sample) point in shape.
- // Returns actual point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const point& sample,
-
- scalar& nearestDistSqr,
- label& nearestIndex,
- point& nearestPoint
- ) const;
-
- //- Calculates nearest (to line) point in shape.
- // Returns point and distance (squared)
- void findNearest
- (
- const labelUList& indices,
- const linePointRef& ln,
-
- treeBoundBox& tightest,
- label& minIndex,
- point& linePoint,
- point& nearestPoint
- ) const;
-
- //- Calculate intersection of triangle with ray. Sets result
- // accordingly
- bool intersects
- (
- const label index,
- const point& start,
- const point& end,
- point& result
- ) const;
-
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/meshTools/meshSearch/meshSearch.C b/src/meshTools/meshSearch/meshSearch.C
index f05d31710435f424066c378c20d0ef5f0c605744..a28b88b625cbe627ac6b47f118d01452b4608110 100644
--- a/src/meshTools/meshSearch/meshSearch.C
+++ b/src/meshTools/meshSearch/meshSearch.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -947,9 +947,7 @@ Foam::List<Foam::pointIndexHit> Foam::meshSearch::intersections
bool Foam::meshSearch::isInside(const point& p) const
{
- return
- boundaryTree().getVolumeType(p)
- == indexedOctree<treeDataFace>::INSIDE;
+ return (boundaryTree().getVolumeType(p) == volumeType::INSIDE);
}
diff --git a/src/meshTools/searchableSurface/searchableBox.C b/src/meshTools/searchableSurface/searchableBox.C
index 3503308e422347620127daecb09212fa12e32175..7d639b99f8697205b13c4d6159a395b6630c65a9 100644
--- a/src/meshTools/searchableSurface/searchableBox.C
+++ b/src/meshTools/searchableSurface/searchableBox.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -583,7 +583,7 @@ void Foam::searchableBox::getVolumeType
) const
{
volType.setSize(points.size());
- volType = INSIDE;
+ volType = volumeType::INSIDE;
forAll(points, pointI)
{
@@ -593,7 +593,7 @@ void Foam::searchableBox::getVolumeType
{
if (pt[dir] < min()[dir] || pt[dir] > max()[dir])
{
- volType[pointI] = OUTSIDE;
+ volType[pointI] = volumeType::OUTSIDE;
break;
}
}
diff --git a/src/meshTools/searchableSurface/searchableCylinder.C b/src/meshTools/searchableSurface/searchableCylinder.C
index 3231dd95cd36104c06616ad0fa831a457c1f3baf..70d674d521b580acae0845f51890be762d07ce73 100644
--- a/src/meshTools/searchableSurface/searchableCylinder.C
+++ b/src/meshTools/searchableSurface/searchableCylinder.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -673,7 +673,7 @@ void Foam::searchableCylinder::getVolumeType
) const
{
volType.setSize(points.size());
- volType = INSIDE;
+ volType = volumeType::INSIDE;
forAll(points, pointI)
{
@@ -687,12 +687,12 @@ void Foam::searchableCylinder::getVolumeType
if (parallel < 0)
{
// left of point1 endcap
- volType[pointI] = OUTSIDE;
+ volType[pointI] = volumeType::OUTSIDE;
}
else if (parallel > magDir_)
{
// right of point2 endcap
- volType[pointI] = OUTSIDE;
+ volType[pointI] = volumeType::OUTSIDE;
}
else
{
@@ -701,11 +701,11 @@ void Foam::searchableCylinder::getVolumeType
if (mag(v) > radius_)
{
- volType[pointI] = OUTSIDE;
+ volType[pointI] = volumeType::OUTSIDE;
}
else
{
- volType[pointI] = INSIDE;
+ volType[pointI] = volumeType::INSIDE;
}
}
}
diff --git a/src/meshTools/searchableSurface/searchableSphere.C b/src/meshTools/searchableSurface/searchableSphere.C
index 5f3e98af90d52edb319a805e95b9bbd9cf4491e2..3d17b52d9fd2c6a392462b8a1c2ad73dfa46c18e 100644
--- a/src/meshTools/searchableSurface/searchableSphere.C
+++ b/src/meshTools/searchableSurface/searchableSphere.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -334,7 +334,7 @@ void Foam::searchableSphere::getVolumeType
) const
{
volType.setSize(points.size());
- volType = INSIDE;
+ volType = volumeType::INSIDE;
forAll(points, pointI)
{
@@ -342,11 +342,11 @@ void Foam::searchableSphere::getVolumeType
if (magSqr(pt - centre_) <= sqr(radius_))
{
- volType[pointI] = INSIDE;
+ volType[pointI] = volumeType::INSIDE;
}
else
{
- volType[pointI] = OUTSIDE;
+ volType[pointI] = volumeType::OUTSIDE;
}
}
}
diff --git a/src/meshTools/searchableSurface/searchableSurface.C b/src/meshTools/searchableSurface/searchableSurface.C
index decd4c7259ab41d187d438985f5517c2a40527e6..102b8cc8ef9d1b9eac3fe2336812ef39229944be 100644
--- a/src/meshTools/searchableSurface/searchableSurface.C
+++ b/src/meshTools/searchableSurface/searchableSurface.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -30,25 +30,11 @@ License
namespace Foam
{
- defineTypeNameAndDebug(searchableSurface, 0);
- defineRunTimeSelectionTable(searchableSurface, dict);
-
- template<>
- const char* Foam::NamedEnum
- <
- Foam::searchableSurface::volumeType,
- 4
- >::names[] =
- {
- "unknown",
- "mixed",
- "inside",
- "outside"
- };
-}
-const Foam::NamedEnum<Foam::searchableSurface::volumeType, 4>
- Foam::searchableSurface::volumeTypeNames;
+defineTypeNameAndDebug(searchableSurface, 0);
+defineRunTimeSelectionTable(searchableSurface, dict);
+
+}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/meshTools/searchableSurface/searchableSurface.H b/src/meshTools/searchableSurface/searchableSurface.H
index baeb970844dfe8e7757afa4c0f3679866a23d17a..90224ae9ce54c517b3cdf589467b3533607d930e 100644
--- a/src/meshTools/searchableSurface/searchableSurface.H
+++ b/src/meshTools/searchableSurface/searchableSurface.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -48,6 +48,7 @@ SourceFiles
#include "pointIndexHit.H"
#include "linePointRef.H"
#include "objectRegistry.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -67,25 +68,6 @@ class searchableSurface
:
public regIOobject
{
-public:
-
- // Data types
-
- //- Volume types
- enum volumeType
- {
- UNKNOWN = 0,
- MIXED = 1, // not used. only here to maintain consistency with
- // indexedOctree volumeType.
- INSIDE = 2,
- OUTSIDE = 3
- };
-
- static const NamedEnum<volumeType, 4> volumeTypeNames;
-
-
-private:
-
// Private data
const word name_;
@@ -278,6 +260,19 @@ public:
List<pointIndexHit>&
) const = 0;
+ //- Find the nearest locations for the supplied points to a
+ // particular region in the searchable surface.
+ virtual void findNearest
+ (
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ List<pointIndexHit>& info
+ ) const
+ {
+ findNearest(samples, nearestDistSqr, info);
+ }
+
//- Find first intersection on segment from start to end.
// Note: searchableSurfacesQueries expects no
// intersection to be found if start==end. Is problem?
diff --git a/src/meshTools/searchableSurface/searchableSurfaces.C b/src/meshTools/searchableSurface/searchableSurfaces.C
index 0612e20d80744494d533b0896240b8b79bed12b6..08569aba78c155437d42bbf086470fd11dedc419 100644
--- a/src/meshTools/searchableSurface/searchableSurfaces.C
+++ b/src/meshTools/searchableSurface/searchableSurfaces.C
@@ -293,6 +293,18 @@ Foam::label Foam::searchableSurfaces::findSurfaceID
}
+Foam::label Foam::searchableSurfaces::findSurfaceRegionID
+(
+ const word& surfaceName,
+ const word& regionName
+) const
+{
+ label surfaceIndex = findSurfaceID(surfaceName);
+
+ return findIndex(regionNames()[surfaceIndex], regionName);
+}
+
+
// Find any intersection
void Foam::searchableSurfaces::findAnyIntersection
(
@@ -382,6 +394,28 @@ void Foam::searchableSurfaces::findNearest
}
+// Find nearest. Return -1 or nearest point
+void Foam::searchableSurfaces::findNearest
+(
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ labelList& nearestSurfaces,
+ List<pointIndexHit>& nearestInfo
+) const
+{
+ searchableSurfacesQueries::findNearest
+ (
+ *this,
+ allSurfaces_,
+ samples,
+ nearestDistSqr,
+ regionIndices,
+ nearestSurfaces,
+ nearestInfo
+ );
+}
+
//- Calculate bounding box
Foam::boundBox Foam::searchableSurfaces::bounds() const
{
diff --git a/src/meshTools/searchableSurface/searchableSurfaces.H b/src/meshTools/searchableSurface/searchableSurfaces.H
index 107e343ef3d262f077a105dfce0441f757b435a5..22ba450e9d2db6aae1d08710d7849a650319fd11 100644
--- a/src/meshTools/searchableSurface/searchableSurfaces.H
+++ b/src/meshTools/searchableSurface/searchableSurfaces.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -140,6 +140,11 @@ public:
//- Find index of surface. Return -1 if not found.
label findSurfaceID(const word& name) const;
+ label findSurfaceRegionID
+ (
+ const word& surfaceName,
+ const word& regionName
+ ) const;
// Multiple point queries.
@@ -185,6 +190,15 @@ public:
List<pointIndexHit>&
) const;
+ void findNearest
+ (
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ labelList& nearestSurfaces,
+ List<pointIndexHit>& nearestInfo
+ ) const;
+
//- Calculate bounding box
boundBox bounds() const;
diff --git a/src/meshTools/searchableSurface/searchableSurfacesQueries.C b/src/meshTools/searchableSurface/searchableSurfacesQueries.C
index b2144cf24aa38c06393e9bc4ec561d4aa2c67e5f..2073de3d767d44972d07174948fafb552e286405 100644
--- a/src/meshTools/searchableSurface/searchableSurfacesQueries.C
+++ b/src/meshTools/searchableSurface/searchableSurfacesQueries.C
@@ -324,9 +324,8 @@ bool Foam::searchableSurfacesQueries::morphTet
void Foam::searchableSurfacesQueries::mergeHits
(
const point& start,
- const scalar mergeDist,
- const label testI, // index of surface
+ const label testI, // index of surface
const List<pointIndexHit>& surfHits, // hits on surface
labelList& allSurfaces,
@@ -338,7 +337,7 @@ void Foam::searchableSurfacesQueries::mergeHits
scalarList surfDistSqr(surfHits.size());
forAll(surfHits, i)
{
- surfDistSqr[i] = magSqr(surfHits[i].hitPoint()-start);
+ surfDistSqr[i] = magSqr(surfHits[i].hitPoint() - start);
}
forAll(surfDistSqr, i)
@@ -346,11 +345,7 @@ void Foam::searchableSurfacesQueries::mergeHits
label index = findLower(allDistSqr, surfDistSqr[i]);
// Check if equal to lower.
- if
- (
- index >= 0
- && (mag(allDistSqr[index]-surfDistSqr[i]) < mergeDist)
- )
+ if (index >= 0)
{
// Same. Do not count.
//Pout<< "point:" << surfHits[i].hitPoint()
@@ -361,12 +356,9 @@ void Foam::searchableSurfacesQueries::mergeHits
else
{
// Check if equal to higher
- label next = index+1;
- if
- (
- next < allDistSqr.size()
- && (mag(allDistSqr[next]-surfDistSqr[i]) < mergeDist)
- )
+ label next = index + 1;
+
+ if (next < allDistSqr.size())
{
//Pout<< "point:" << surfHits[i].hitPoint()
// << " considered same as:" << allInfo[next].hitPoint()
@@ -510,21 +502,6 @@ void Foam::searchableSurfacesQueries::findAllIntersections
if (surfacesToTest.size() > 1)
{
- // Small vector to increment start vector by to find next intersection
- // along line. Constant factor added to make sure that start==end still
- // ends iteration in findAllIntersections. Also SMALL is just slightly
- // too small.
- const vectorField smallVec
- (
- 1E2*SMALL*(end-start)
- + vector(ROOTVSMALL,ROOTVSMALL,ROOTVSMALL)
- );
-
- // Tolerance used to check whether points are equal. Note: used to
- // compare distance^2. Note that we use the maximum possible tolerance
- // (reached at intersections close to the end point)
- const scalarField mergeDist(2*mag(smallVec)*mag(end-start));
-
// Test the other surfaces and merge (according to distance from start).
for (label testI = 1; testI < surfacesToTest.size(); testI++)
{
@@ -541,7 +518,6 @@ void Foam::searchableSurfacesQueries::findAllIntersections
mergeHits
(
start[pointI], // Current segment
- mergeDist[pointI],
testI, // Surface and its hits
surfHits[pointI],
@@ -691,13 +667,75 @@ void Foam::searchableSurfacesQueries::findNearest
}
+// Find nearest. Return -1 or nearest point
+void Foam::searchableSurfacesQueries::findNearest
+(
+ const PtrList<searchableSurface>& allSurfaces,
+ const labelList& surfacesToTest,
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ labelList& nearestSurfaces,
+ List<pointIndexHit>& nearestInfo
+)
+{
+ if (regionIndices.empty())
+ {
+ findNearest
+ (
+ allSurfaces,
+ surfacesToTest,
+ samples,
+ nearestDistSqr,
+ nearestSurfaces,
+ nearestInfo
+ );
+ }
+
+ // Initialise
+ nearestSurfaces.setSize(samples.size());
+ nearestSurfaces = -1;
+ nearestInfo.setSize(samples.size());
+
+ // Work arrays
+ scalarField minDistSqr(nearestDistSqr);
+ List<pointIndexHit> hitInfo(samples.size());
+
+ forAll(surfacesToTest, testI)
+ {
+ allSurfaces[surfacesToTest[testI]].findNearest
+ (
+ samples,
+ minDistSqr,
+ regionIndices,
+ hitInfo
+ );
+
+ // Update minDistSqr and arguments
+ forAll(hitInfo, pointI)
+ {
+ if (hitInfo[pointI].hit())
+ {
+ minDistSqr[pointI] = magSqr
+ (
+ hitInfo[pointI].hitPoint()
+ - samples[pointI]
+ );
+ nearestInfo[pointI] = hitInfo[pointI];
+ nearestSurfaces[pointI] = testI;
+ }
+ }
+ }
+}
+
+
void Foam::searchableSurfacesQueries::signedDistance
(
const PtrList<searchableSurface>& allSurfaces,
const labelList& surfacesToTest,
const pointField& samples,
const scalarField& nearestDistSqr,
- const searchableSurface::volumeType illegalHandling,
+ const volumeType illegalHandling,
labelList& nearestSurfaces,
scalarField& distance
)
@@ -737,7 +775,7 @@ void Foam::searchableSurfacesQueries::signedDistance
}
// Calculate sideness of these surface points
- List<searchableSurface::volumeType> volType;
+ List<volumeType> volType;
allSurfaces[surfacesToTest[testI]].getVolumeType(surfPoints, volType);
// Push back to original
@@ -746,13 +784,13 @@ void Foam::searchableSurfacesQueries::signedDistance
label pointI = surfIndices[i];
scalar dist = mag(samples[pointI] - nearestInfo[pointI].hitPoint());
- searchableSurface::volumeType vT = volType[i];
+ volumeType vT = volType[i];
- if (vT == searchableSurface::OUTSIDE)
+ if (vT == volumeType::OUTSIDE)
{
distance[pointI] = dist;
}
- else if (vT == searchableSurface::INSIDE)
+ else if (vT == volumeType::INSIDE)
{
distance[i] = -dist;
}
@@ -760,12 +798,12 @@ void Foam::searchableSurfacesQueries::signedDistance
{
switch (illegalHandling)
{
- case searchableSurface::OUTSIDE:
+ case volumeType::OUTSIDE:
{
distance[pointI] = dist;
break;
}
- case searchableSurface::INSIDE:
+ case volumeType::INSIDE:
{
distance[pointI] = -dist;
break;
@@ -779,7 +817,7 @@ void Foam::searchableSurfacesQueries::signedDistance
<< " surface:"
<< allSurfaces[surfacesToTest[testI]].name()
<< " volType:"
- << searchableSurface::volumeTypeNames[vT]
+ << volumeType::names[vT]
<< exit(FatalError);
break;
}
diff --git a/src/meshTools/searchableSurface/searchableSurfacesQueries.H b/src/meshTools/searchableSurface/searchableSurfacesQueries.H
index 9f2ebfbadb6c575aaf0c5fb9a431acf6079df455..529d4b94fabb44b8d08695573fd5f0f60d204226 100644
--- a/src/meshTools/searchableSurface/searchableSurfacesQueries.H
+++ b/src/meshTools/searchableSurface/searchableSurfacesQueries.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -116,7 +116,6 @@ class searchableSurfacesQueries
static void mergeHits
(
const point& start,
- const scalar mergeDist,
const label surfI,
const List<pointIndexHit>& surfHits,
@@ -184,6 +183,18 @@ public:
List<pointIndexHit>&
);
+ //- Find nearest points to a specific region of the surface
+ static void findNearest
+ (
+ const PtrList<searchableSurface>& allSurfaces,
+ const labelList& surfacesToTest,
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ labelList& nearestSurfaces,
+ List<pointIndexHit>& nearestInfo
+ );
+
//- Find signed distance to nearest surface. Outside is positive.
// illegalHandling: how to handle non-inside or outside
// OUTSIDE : treat as outside
@@ -195,7 +206,7 @@ public:
const labelList& surfacesToTest,
const pointField& samples,
const scalarField& nearestDistSqr,
- const searchableSurface::volumeType illegalHandling,
+ const volumeType illegalHandling,
labelList& nearestSurfaces,
scalarField& distance
);
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C
index 4e42f5d3adce4ae3cbdf8486c76430f910a2574d..22229817932e31a1704fbc6608975d6e8b65ec48 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -29,7 +29,7 @@ License
#include "EdgeMap.H"
#include "triSurfaceFields.H"
#include "Time.H"
-#include "PackedBoolList.H"
+#include "PatchTools.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -218,99 +218,6 @@ bool Foam::triSurfaceMesh::isSurfaceClosed() const
}
-// Gets all intersections after initial one. Adds smallVec and starts tracking
-// from there.
-void Foam::triSurfaceMesh::getNextIntersections
-(
- const indexedOctree<treeDataTriSurface>& octree,
- const point& start,
- const point& end,
- const vector& smallVec,
- DynamicList<pointIndexHit, 1, 1>& hits
-)
-{
- const vector dirVec(end-start);
- const scalar magSqrDirVec(magSqr(dirVec));
-
- // Initial perturbation amount
- vector perturbVec(smallVec);
-
- while (true)
- {
- // Start tracking from last hit.
- point pt = hits.last().hitPoint() + perturbVec;
-
- if (((pt-start)&dirVec) > magSqrDirVec)
- {
- return;
- }
-
- // See if any intersection between pt and end
- pointIndexHit inter = octree.findLine(pt, end);
-
- if (!inter.hit())
- {
- return;
- }
-
- // Check if already found this intersection
- bool duplicateHit = false;
- forAllReverse(hits, i)
- {
- if (hits[i].index() == inter.index())
- {
- duplicateHit = true;
- break;
- }
- }
-
-
- if (duplicateHit)
- {
- // Hit same triangle again. Increase perturbVec and try again.
- perturbVec *= 2;
- }
- else
- {
- // Proper hit
- hits.append(inter);
- // Restore perturbVec
- perturbVec = smallVec;
- }
- }
-}
-
-
-void Foam::triSurfaceMesh::calcBounds(boundBox& bb, label& nPoints) const
-{
- // Unfortunately nPoints constructs meshPoints() so do compact version
- // ourselves
-
- const triSurface& s = static_cast<const triSurface&>(*this);
-
- PackedBoolList pointIsUsed(points()().size());
-
- nPoints = 0;
- bb = boundBox::invertedBox;
-
- forAll(s, faceI)
- {
- const triSurface::FaceType& f = s[faceI];
-
- forAll(f, fp)
- {
- label pointI = f[fp];
- if (pointIsUsed.set(pointI, 1u))
- {
- bb.min() = ::Foam::min(bb.min(), points()()[pointI]);
- bb.max() = ::Foam::max(bb.max(), points()()[pointI]);
- nPoints++;
- }
- }
- }
-}
-
-
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
@@ -330,9 +237,8 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
)
),
triSurface(s),
- tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ triSurfaceRegionSearch(s),
minQuality_(-1),
- maxTreeDepth_(10),
surfaceClosed_(-1)
{
bounds() = boundBox(points());
@@ -377,9 +283,8 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
searchableSurface::objectPath()
)
),
- tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ triSurfaceRegionSearch(static_cast<const triSurface&>(*this)),
minQuality_(-1),
- maxTreeDepth_(10),
surfaceClosed_(-1)
{
bounds() = boundBox(points());
@@ -427,9 +332,8 @@ Foam::triSurfaceMesh::triSurfaceMesh
searchableSurface::objectPath()
)
),
- tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ triSurfaceRegionSearch(static_cast<const triSurface&>(*this), dict),
minQuality_(-1),
- maxTreeDepth_(10),
surfaceClosed_(-1)
{
scalar scaleFactor = 0;
@@ -445,13 +349,6 @@ Foam::triSurfaceMesh::triSurfaceMesh
bounds() = boundBox(points());
- // Have optional non-standard search tolerance for gappy surfaces.
- if (dict.readIfPresent("tolerance", tolerance_) && tolerance_ > 0)
- {
- Info<< searchableSurface::name() << " : using intersection tolerance "
- << tolerance_ << endl;
- }
-
// Have optional minimum quality for normal calculation
if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
{
@@ -459,13 +356,6 @@ Foam::triSurfaceMesh::triSurfaceMesh
<< " : ignoring triangles with quality < "
<< minQuality_ << " for normals calculation." << endl;
}
-
- // Have optional non-standard tree-depth to limit storage.
- if (dict.readIfPresent("maxTreeDepth", maxTreeDepth_) && maxTreeDepth_ > 0)
- {
- Info<< searchableSurface::name() << " : using maximum tree depth "
- << maxTreeDepth_ << endl;
- }
}
@@ -479,7 +369,7 @@ Foam::triSurfaceMesh::~triSurfaceMesh()
void Foam::triSurfaceMesh::clearOut()
{
- tree_.clear();
+ triSurfaceRegionSearch::clearOut();
edgeTree_.clear();
triSurface::clearOut();
}
@@ -521,64 +411,12 @@ bool Foam::triSurfaceMesh::overlaps(const boundBox& bb) const
void Foam::triSurfaceMesh::movePoints(const pointField& newPoints)
{
- tree_.clear();
+ triSurfaceRegionSearch::clearOut();
edgeTree_.clear();
triSurface::movePoints(newPoints);
}
-const Foam::indexedOctree<Foam::treeDataTriSurface>&
-Foam::triSurfaceMesh::tree() const
-{
- if (tree_.empty())
- {
- // Calculate bb without constructing local point numbering.
- treeBoundBox bb;
- label nPoints;
- calcBounds(bb, nPoints);
-
- if (nPoints != points()().size())
- {
- WarningIn("triSurfaceMesh::tree() const")
- << "Surface " << searchableSurface::name()
- << " does not have compact point numbering."
- << " Of " << points()().size() << " only " << nPoints
- << " are used. This might give problems in some routines."
- << endl;
- }
-
-
- // Random number generator. Bit dodgy since not exactly random ;-)
- Random rndGen(65431);
-
- // Slightly extended bb. Slightly off-centred just so on symmetric
- // geometry there are less face/edge aligned items.
- bb = bb.extend(rndGen, 1e-4);
- bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
- bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-
- scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
-
- tree_.reset
- (
- new indexedOctree<treeDataTriSurface>
- (
- treeDataTriSurface(*this, tolerance_),
- bb,
- maxTreeDepth_, // maxLevel
- 10, // leafsize
- 3.0 // duplicity
- )
- );
-
- indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
- }
-
- return tree_();
-}
-
-
const Foam::indexedOctree<Foam::treeDataEdge>&
Foam::triSurfaceMesh::edgeTree() const
{
@@ -597,7 +435,12 @@ Foam::triSurfaceMesh::edgeTree() const
treeBoundBox bb;
label nPoints;
- calcBounds(bb, nPoints);
+ PatchTools::calcBounds
+ (
+ static_cast<const triSurface&>(*this),
+ bb,
+ nPoints
+ );
// Random number generator. Bit dodgy since not exactly random ;-)
Random rndGen(65431);
@@ -609,8 +452,8 @@ Foam::triSurfaceMesh::edgeTree() const
bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
- scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataEdge>::perturbTol() = tolerance_;
+ scalar oldTol = indexedOctree<treeDataEdge>::perturbTol();
+ indexedOctree<treeDataEdge>::perturbTol() = tolerance();
edgeTree_.reset
(
@@ -624,7 +467,7 @@ Foam::triSurfaceMesh::edgeTree() const
bEdges // selected edges
),
bb, // bb
- maxTreeDepth_, // maxLevel
+ maxTreeDepth(), // maxLevel
10, // leafsize
3.0 // duplicity
)
@@ -636,12 +479,6 @@ Foam::triSurfaceMesh::edgeTree() const
}
-Foam::scalar Foam::triSurfaceMesh::tolerance() const
-{
- return tolerance_;
-}
-
-
const Foam::wordList& Foam::triSurfaceMesh::regions() const
{
if (regions_.empty())
@@ -682,23 +519,25 @@ void Foam::triSurfaceMesh::findNearest
List<pointIndexHit>& info
) const
{
- const indexedOctree<treeDataTriSurface>& octree = tree();
-
- info.setSize(samples.size());
-
- scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
+ triSurfaceSearch::findNearest(samples, nearestDistSqr, info);
+}
- forAll(samples, i)
- {
- static_cast<pointIndexHit&>(info[i]) = octree.findNearest
- (
- samples[i],
- nearestDistSqr[i]
- );
- }
- indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
+void Foam::triSurfaceMesh::findNearest
+(
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ List<pointIndexHit>& info
+) const
+{
+ triSurfaceRegionSearch::findNearest
+ (
+ samples,
+ nearestDistSqr,
+ regionIndices,
+ info
+ );
}
@@ -709,23 +548,7 @@ void Foam::triSurfaceMesh::findLine
List<pointIndexHit>& info
) const
{
- const indexedOctree<treeDataTriSurface>& octree = tree();
-
- info.setSize(start.size());
-
- scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
-
- forAll(start, i)
- {
- static_cast<pointIndexHit&>(info[i]) = octree.findLine
- (
- start[i],
- end[i]
- );
- }
-
- indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
+ triSurfaceSearch::findLine(start, end, info);
}
@@ -736,23 +559,7 @@ void Foam::triSurfaceMesh::findLineAny
List<pointIndexHit>& info
) const
{
- const indexedOctree<treeDataTriSurface>& octree = tree();
-
- info.setSize(start.size());
-
- scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
-
- forAll(start, i)
- {
- static_cast<pointIndexHit&>(info[i]) = octree.findLineAny
- (
- start[i],
- end[i]
- );
- }
-
- indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
+ triSurfaceSearch::findLineAny(start, end, info);
}
@@ -763,57 +570,7 @@ void Foam::triSurfaceMesh::findLineAll
List<List<pointIndexHit> >& info
) const
{
- const indexedOctree<treeDataTriSurface>& octree = tree();
-
- info.setSize(start.size());
-
- scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
-
- // Work array
- DynamicList<pointIndexHit, 1, 1> hits;
-
- // Tolerances:
- // To find all intersections we add a small vector to the last intersection
- // This is chosen such that
- // - it is significant (SMALL is smallest representative relative tolerance;
- // we need something bigger since we're doing calculations)
- // - if the start-end vector is zero we still progress
- const vectorField dirVec(end-start);
- const vectorField smallVec
- (
- indexedOctree<treeDataTriSurface>::perturbTol()*dirVec
- + vector(ROOTVSMALL,ROOTVSMALL,ROOTVSMALL)
- );
-
- forAll(start, pointI)
- {
- // See if any intersection between pt and end
- pointIndexHit inter = octree.findLine(start[pointI], end[pointI]);
-
- if (inter.hit())
- {
- hits.clear();
- hits.append(inter);
-
- getNextIntersections
- (
- octree,
- start[pointI],
- end[pointI],
- smallVec[pointI],
- hits
- );
-
- info[pointI].transfer(hits);
- }
- else
- {
- info[pointI].clear();
- }
- }
-
- indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
+ triSurfaceSearch::findLineAll(start, end, info);
}
@@ -975,7 +732,7 @@ void Foam::triSurfaceMesh::getVolumeType
volType.setSize(points.size());
scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
- indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
+ indexedOctree<treeDataTriSurface>::perturbTol() = tolerance();
forAll(points, pointI)
{
@@ -983,7 +740,7 @@ void Foam::triSurfaceMesh::getVolumeType
// - use cached volume type per each tree node
// - cheat conversion since same values
- volType[pointI] = static_cast<volumeType>(tree().getVolumeType(pt));
+ volType[pointI] = tree().getVolumeType(pt);
}
indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H
index f27df49d7c85d89d876a63287277a344c8284f06..fa7e9e511c769a0d37b7eb6d9b9a885ecaa367cd 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -46,8 +46,11 @@ SourceFiles
#include "objectRegistry.H"
#include "indexedOctree.H"
#include "treeDataTriSurface.H"
+#include "treeDataPrimitivePatch.H"
#include "treeDataEdge.H"
#include "EdgeMap.H"
+#include "triSurface.H"
+#include "triSurfaceRegionSearch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -62,25 +65,17 @@ class triSurfaceMesh
:
public searchableSurface,
public objectRegistry, // so we can store fields
- public triSurface
+ public triSurface,
+ public triSurfaceRegionSearch
{
private:
// Private member data
- //- Optional tolerance to use in searches
- scalar tolerance_;
-
//- Optional min triangle quality. Triangles below this get
// ignored for normal calculation
scalar minQuality_;
- //- Optional max tree depth of octree
- label maxTreeDepth_;
-
- //- Search tree (triangles)
- mutable autoPtr<indexedOctree<treeDataTriSurface> > tree_;
-
//- Search tree for boundary edges.
mutable autoPtr<indexedOctree<treeDataEdge> > edgeTree_;
@@ -136,11 +131,6 @@ private:
void operator=(const triSurfaceMesh&);
-protected:
-
- //- Calculate (number of)used points and their bounding box
- void calcBounds(boundBox& bb, label& nPoints) const;
-
public:
//- Runtime type information
@@ -176,13 +166,9 @@ public:
//- Move points
virtual void movePoints(const pointField&);
- //- Demand driven contruction of octree
- const indexedOctree<treeDataTriSurface>& tree() const;
-
//- Demand driven contruction of octree for boundary edges
const indexedOctree<treeDataEdge>& edgeTree() const;
- scalar tolerance() const;
// searchableSurface implementation
@@ -214,6 +200,14 @@ public:
List<pointIndexHit>&
) const;
+ virtual void findNearest
+ (
+ const pointField& sample,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ List<pointIndexHit>&
+ ) const;
+
virtual void findLine
(
const pointField& start,
diff --git a/src/meshTools/sets/pointSources/surfaceToPoint/surfaceToPoint.C b/src/meshTools/sets/pointSources/surfaceToPoint/surfaceToPoint.C
index 4aaf9eee17664ef0ea78003a83f29ab9666ca9fa..34381e4a027f78c847581301a5fb7ceceedc1d41 100644
--- a/src/meshTools/sets/pointSources/surfaceToPoint/surfaceToPoint.C
+++ b/src/meshTools/sets/pointSources/surfaceToPoint/surfaceToPoint.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -26,6 +26,7 @@ License
#include "surfaceToPoint.H"
#include "polyMesh.H"
#include "triSurfaceSearch.H"
+#include "triSurface.H"
#include "cpuTime.H"
#include "addToRunTimeSelectionTable.H"
diff --git a/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C b/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C
index 9af93cf0b7be7b5b6962ded550c1f707903dcc06..b54ad7143d7adce09b0f47a6e392d347b3d054e4 100644
--- a/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C
+++ b/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -694,6 +694,7 @@ void Foam::surfaceIntersection::doCutEdges
}
while (doTrack);
}
+
intersection::setPlanarTol(oldTol);
}
diff --git a/src/meshTools/triSurface/orientedSurface/orientedSurface.C b/src/meshTools/triSurface/orientedSurface/orientedSurface.C
index 95acfa125c60ad32dbb8b5bc4ffc43f0d51939b3..0466343eb5307c5690dbdd76fe75ec6c90f372e7 100644
--- a/src/meshTools/triSurface/orientedSurface/orientedSurface.C
+++ b/src/meshTools/triSurface/orientedSurface/orientedSurface.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -258,7 +258,8 @@ void Foam::orientedSurface::findZoneSide
zoneFaceI = -1;
isOutside = false;
- List<pointIndexHit> hits;
+ pointField start(1, outsidePoint);
+ List<List<pointIndexHit> > hits(1, List<pointIndexHit>());
forAll(faceZone, faceI)
{
@@ -273,7 +274,7 @@ void Foam::orientedSurface::findZoneSide
// Check if normal different enough to decide upon
if (magD > SMALL && (mag(n & d/magD) > 1e-6))
{
- point end = fc + d;
+ pointField end(1, fc + d);
//Info<< "Zone " << zoneI << " : Shooting ray"
// << " from " << outsidePoint
@@ -281,12 +282,13 @@ void Foam::orientedSurface::findZoneSide
// << " to pierce triangle " << faceI
// << " with centre " << fc << endl;
- surfSearches.findLineAll(outsidePoint, end, hits);
+
+ surfSearches.findLineAll(start, end, hits);
label zoneIndex = -1;
- forAll(hits, i)
+ forAll(hits[0], i)
{
- if (hits[i].index() == faceI)
+ if (hits[0][i].index() == faceI)
{
zoneIndex = i;
break;
diff --git a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C
new file mode 100644
index 0000000000000000000000000000000000000000..738cb8dcdd183a5a4e9ecdc0f603bae68b73d855
--- /dev/null
+++ b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C
@@ -0,0 +1,246 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "triSurfaceRegionSearch.H"
+#include "indexedOctree.H"
+#include "triSurface.H"
+#include "PatchTools.H"
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::triSurfaceRegionSearch::triSurfaceRegionSearch(const triSurface& surface)
+:
+ triSurfaceSearch(surface),
+ indirectRegionPatches_(),
+ treeByRegion_()
+{}
+
+
+Foam::triSurfaceRegionSearch::triSurfaceRegionSearch
+(
+ const triSurface& surface,
+ const dictionary& dict
+)
+:
+ triSurfaceSearch(surface, dict),
+ indirectRegionPatches_(),
+ treeByRegion_()
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::triSurfaceRegionSearch::~triSurfaceRegionSearch()
+{
+ clearOut();
+}
+
+
+void Foam::triSurfaceRegionSearch::clearOut()
+{
+ triSurfaceSearch::clearOut();
+ treeByRegion_.clear();
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+const Foam::PtrList<Foam::triSurfaceRegionSearch::treeType>&
+Foam::triSurfaceRegionSearch::treeByRegion() const
+{
+ if (treeByRegion_.empty())
+ {
+ Map<label> regionSizes;
+ forAll(surface(), fI)
+ {
+ const label regionI = surface()[fI].region();
+
+ regionSizes(regionI)++;
+ }
+
+ label nRegions = regionSizes.size();
+
+ indirectRegionPatches_.setSize(nRegions);
+ treeByRegion_.setSize(nRegions);
+
+ labelListList regionsAddressing(nRegions);
+
+ forAll(regionsAddressing, regionI)
+ {
+ regionsAddressing[regionI] = labelList(regionSizes[regionI], -1);
+ }
+
+ labelList nFacesInRegions(nRegions, 0);
+
+ forAll(surface(), fI)
+ {
+ const label regionI = surface()[fI].region();
+
+ regionsAddressing[regionI][nFacesInRegions[regionI]++] = fI;
+ }
+
+ forAll(regionsAddressing, regionI)
+ {
+ scalar oldTol = treeType::perturbTol();
+ treeType::perturbTol() = tolerance();
+
+ indirectRegionPatches_.set
+ (
+ regionI,
+ new indirectTriSurface
+ (
+ IndirectList<labelledTri>
+ (
+ surface(),
+ regionsAddressing[regionI]
+ ),
+ surface().points()
+ )
+ );
+
+ // Calculate bb without constructing local point numbering.
+ treeBoundBox bb;
+ label nPoints;
+ PatchTools::calcBounds
+ (
+ indirectRegionPatches_[regionI],
+ bb,
+ nPoints
+ );
+
+// if (nPoints != surface().points().size())
+// {
+// WarningIn("triSurfaceRegionSearch::treeByRegion() const")
+// << "Surface does not have compact point numbering. "
+// << "Of " << surface().points().size()
+// << " only " << nPoints
+// << " are used. This might give problems in some routines."
+// << endl;
+// }
+
+ // Random number generator. Bit dodgy since not exactly random ;-)
+ Random rndGen(65431);
+
+ // Slightly extended bb. Slightly off-centred just so on symmetric
+ // geometry there are fewer face/edge aligned items.
+ bb = bb.extend(rndGen, 1e-4);
+ bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+ bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+
+ treeByRegion_.set
+ (
+ regionI,
+ new treeType
+ (
+ treeDataIndirectTriSurface
+ (
+ true,
+ indirectRegionPatches_[regionI],
+ tolerance()
+ ),
+ bb,
+ maxTreeDepth(), // maxLevel
+ 10, // leafsize
+ 3.0 // duplicity
+ )
+ );
+
+ treeType::perturbTol() = oldTol;
+ }
+ }
+
+ return treeByRegion_;
+}
+
+
+void Foam::triSurfaceRegionSearch::findNearest
+(
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ List<pointIndexHit>& info
+) const
+{
+ if (regionIndices.empty())
+ {
+ triSurfaceSearch::findNearest(samples, nearestDistSqr, info);
+ }
+ else
+ {
+ scalar oldTol = treeType::perturbTol();
+ treeType::perturbTol() = tolerance();
+
+ const PtrList<treeType>& octrees = treeByRegion();
+
+ info.setSize(samples.size());
+
+ forAll(octrees, treeI)
+ {
+ if (findIndex(regionIndices, treeI) == -1)
+ {
+ continue;
+ }
+
+ const treeType& octree = octrees[treeI];
+
+ forAll(samples, i)
+ {
+// if (!octree.bb().contains(samples[i]))
+// {
+// continue;
+// }
+
+ pointIndexHit currentRegionHit = octree.findNearest
+ (
+ samples[i],
+ nearestDistSqr[i],
+ treeDataIndirectTriSurface::findNearestOp(octree)
+ );
+
+ if
+ (
+ currentRegionHit.hit()
+ &&
+ (
+ !info[i].hit()
+ ||
+ (
+ magSqr(currentRegionHit.hitPoint() - samples[i])
+ < magSqr(info[i].hitPoint() - samples[i])
+ )
+ )
+ )
+ {
+ info[i] = currentRegionHit;
+ }
+ }
+ }
+
+ treeType::perturbTol() = oldTol;
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.H b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.H
new file mode 100644
index 0000000000000000000000000000000000000000..c8c421b900200e1c72d3ae62e23432e69b4254d8
--- /dev/null
+++ b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.H
@@ -0,0 +1,144 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::triSurfaceRegionSearch
+
+Description
+ Helper class to search on triSurface. Creates an octree for each region of
+ the surface and only searches on the specified regions.
+
+SourceFiles
+ triSurfaceRegionSearch.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef triSurfaceRegionSearch_H
+#define triSurfaceRegionSearch_H
+
+#include "pointField.H"
+#include "pointIndexHit.H"
+#include "triSurfaceSearch.H"
+#include "labelledTri.H"
+#include "IndirectList.H"
+#include "PtrList.H"
+#include "indexedOctree.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class triSurfaceRegionSearch Declaration
+\*---------------------------------------------------------------------------*/
+
+class triSurfaceRegionSearch
+:
+ public triSurfaceSearch
+{
+ // Private typedefs
+
+ typedef PrimitivePatch
+ <
+ labelledTri,
+ IndirectList,
+ const pointField&
+ > indirectTriSurface;
+
+ typedef treeDataPrimitivePatch<indirectTriSurface>
+ treeDataIndirectTriSurface;
+
+ typedef indexedOctree<treeDataIndirectTriSurface> treeType;
+
+
+ // Private data
+
+ //- Surface is split into patches by region
+ mutable PtrList<indirectTriSurface> indirectRegionPatches_;
+
+ //- Search tree for each region
+ mutable PtrList<treeType> treeByRegion_;
+
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ triSurfaceRegionSearch(const triSurfaceRegionSearch&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const triSurfaceRegionSearch&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from surface. Holds reference to surface!
+ explicit triSurfaceRegionSearch(const triSurface&);
+
+ //- Construct from surface and dictionary. Holds reference to surface!
+ triSurfaceRegionSearch(const triSurface&, const dictionary& dict);
+
+
+ //- Destructor
+ ~triSurfaceRegionSearch();
+
+ //- Clear storage
+ void clearOut();
+
+
+ // Member Functions
+
+ // Access
+
+ //- Demand driven construction of octree for each region.
+ // @todo Currently creates a tree for each region; could optimise
+ // by only constructing trees when they are in regionIndices
+ const PtrList<treeType>& treeByRegion() const;
+
+
+ // Query
+
+ //- Find the nearest point on the surface out of the regions
+ // supplied in the list regionIndices. Ignores regions that are
+ // not specified
+ void findNearest
+ (
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ const labelList& regionIndices,
+ List<pointIndexHit>& info
+ ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C
index 28f6393b2c28b84d209dc10b1c3d0ce7ccb7b468..aad1a824eb8985a9b0c7cddf94a2f1f74db51a9c 100644
--- a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C
+++ b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -24,62 +24,203 @@ License
\*---------------------------------------------------------------------------*/
#include "triSurfaceSearch.H"
-#include "indexedOctree.H"
-#include "boolList.H"
-#include "treeDataTriSurface.H"
#include "triSurface.H"
-#include "line.H"
-#include "cpuTime.H"
+#include "PatchTools.H"
+#include "volumeType.H"
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-const Foam::point Foam::triSurfaceSearch::greatPoint(GREAT, GREAT, GREAT);
+bool Foam::triSurfaceSearch::checkUniqueHit
+(
+ const pointIndexHit& currHit,
+ const DynamicList<pointIndexHit, 1, 1>& hits,
+ const vector& lineVec
+) const
+{
+ // Classify the hit
+ label nearType = -1;
+ label nearLabel = -1;
+
+ const labelledTri& f = surface()[currHit.index()];
+
+ f.nearestPointClassify
+ (
+ currHit.hitPoint(),
+ surface().points(),
+ nearType,
+ nearLabel
+ );
+
+ if (nearType == 1)
+ {
+ // near point
+ }
+ else if (nearType == 2)
+ {
+ // near edge
+ // check if the other face of the edge is already hit
+
+ const labelList& fEdges = surface().faceEdges()[currHit.index()];
+
+ const label edgeI = fEdges[nearLabel];
+
+ const labelList& edgeFaces = surface().edgeFaces()[edgeI];
+
+ forAll(edgeFaces, fI)
+ {
+ const label edgeFaceI = edgeFaces[fI];
+
+ if (edgeFaceI != currHit.index())
+ {
+ forAll(hits, hI)
+ {
+ const pointIndexHit& hit = hits[hI];
+
+ if (hit.index() == edgeFaceI)
+ {
+ // Check normals
+ const vector currHitNormal =
+ surface().faceNormals()[currHit.index()];
+
+ const vector existingHitNormal =
+ surface().faceNormals()[edgeFaceI];
+
+ const label signCurrHit =
+ pos(currHitNormal & lineVec);
+
+ const label signExistingHit =
+ pos(existingHitNormal & lineVec);
+
+ if (signCurrHit == signExistingHit)
+ {
+ return false;
+ }
+ }
+ }
+ }
+ }
+ }
+
+ return true;
+}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-// Construct from surface. Holds reference!
Foam::triSurfaceSearch::triSurfaceSearch(const triSurface& surface)
:
surface_(surface),
+ tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ maxTreeDepth_(10),
+ treePtr_(NULL)
+{}
+
+
+Foam::triSurfaceSearch::triSurfaceSearch
+(
+ const triSurface& surface,
+ const dictionary& dict
+)
+:
+ surface_(surface),
+ tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ maxTreeDepth_(10),
treePtr_(NULL)
{
- // Random number generator. Bit dodgy since not exactly random ;-)
- Random rndGen(65431);
+ // Have optional non-standard search tolerance for gappy surfaces.
+ if (dict.readIfPresent("tolerance", tolerance_) && tolerance_ > 0)
+ {
+ Info<< " using intersection tolerance " << tolerance_ << endl;
+ }
- // Slightly extended bb. Slightly off-centred just so on symmetric
- // geometry there are less face/edge aligned items.
- treeBoundBox treeBb
- (
- treeBoundBox(surface_.points(), surface_.meshPoints()).extend
- (
- rndGen,
- 1e-4
- )
- );
- treeBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
- treeBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+ // Have optional non-standard tree-depth to limit storage.
+ if (dict.readIfPresent("maxTreeDepth", maxTreeDepth_) && maxTreeDepth_ > 0)
+ {
+ Info<< " using maximum tree depth " << maxTreeDepth_ << endl;
+ }
+}
- treePtr_.reset
- (
- new indexedOctree<treeDataTriSurface>
- (
- treeDataTriSurface
- (
- surface_,
- indexedOctree<treeDataTriSurface>::perturbTol()
- ),
- treeBb,
- 8, // maxLevel
- 10, // leafsize
- 3.0 // duplicity
- )
- );
+
+Foam::triSurfaceSearch::triSurfaceSearch
+(
+ const triSurface& surface,
+ const scalar tolerance,
+ const label maxTreeDepth
+)
+:
+ surface_(surface),
+ tolerance_(tolerance),
+ maxTreeDepth_(maxTreeDepth),
+ treePtr_(NULL)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::triSurfaceSearch::~triSurfaceSearch()
+{
+ clearOut();
+}
+
+
+void Foam::triSurfaceSearch::clearOut()
+{
+ treePtr_.clear();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+const Foam::indexedOctree<Foam::treeDataTriSurface>&
+Foam::triSurfaceSearch::tree() const
+{
+ if (treePtr_.empty())
+ {
+ // Calculate bb without constructing local point numbering.
+ treeBoundBox bb;
+ label nPoints;
+ PatchTools::calcBounds(surface(), bb, nPoints);
+
+ if (nPoints != surface().points().size())
+ {
+ WarningIn("triSurfaceSearch::tree() const")
+ << "Surface does not have compact point numbering."
+ << " Of " << surface().points().size() << " only " << nPoints
+ << " are used. This might give problems in some routines."
+ << endl;
+ }
+
+ // Random number generator. Bit dodgy since not exactly random ;-)
+ Random rndGen(65431);
+
+ // Slightly extended bb. Slightly off-centred just so on symmetric
+ // geometry there are less face/edge aligned items.
+ bb = bb.extend(rndGen, 1e-4);
+ bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+ bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+
+ scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
+ indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
+
+ treePtr_.reset
+ (
+ new indexedOctree<treeDataTriSurface>
+ (
+ treeDataTriSurface(true, surface_, tolerance_),
+ bb,
+ maxTreeDepth_, // maxLevel
+ 10, // leafsize
+ 3.0 // duplicity
+ )
+ );
+
+ indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
+ }
+
+ return treePtr_();
+}
+
+
// Determine inside/outside for samples
Foam::boolList Foam::triSurfaceSearch::calcInside
(
@@ -96,11 +237,7 @@ Foam::boolList Foam::triSurfaceSearch::calcInside
{
inside[sampleI] = false;
}
- else if
- (
- tree().getVolumeType(sample)
- == indexedOctree<treeDataTriSurface>::INSIDE
- )
+ else if (tree().getVolumeType(sample) == volumeType::INSIDE)
{
inside[sampleI] = true;
}
@@ -113,159 +250,167 @@ Foam::boolList Foam::triSurfaceSearch::calcInside
}
-Foam::labelList Foam::triSurfaceSearch::calcNearestTri
+void Foam::triSurfaceSearch::findNearest
(
const pointField& samples,
- const vector& span
+ const scalarField& nearestDistSqr,
+ List<pointIndexHit>& info
) const
{
- labelList nearest(samples.size());
+ scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
+ indexedOctree<treeDataTriSurface>::perturbTol() = tolerance();
- const scalar nearestDistSqr = 0.25*magSqr(span);
+ const indexedOctree<treeDataTriSurface>& octree = tree();
- pointIndexHit hitInfo;
+ info.setSize(samples.size());
- forAll(samples, sampleI)
+ forAll(samples, i)
{
- hitInfo = tree().findNearest(samples[sampleI], nearestDistSqr);
-
- if (hitInfo.hit())
- {
- nearest[sampleI] = hitInfo.index();
- }
- else
- {
- nearest[sampleI] = -1;
- }
+ static_cast<pointIndexHit&>(info[i]) = octree.findNearest
+ (
+ samples[i],
+ nearestDistSqr[i],
+ treeDataTriSurface::findNearestOp(octree)
+ );
}
- return nearest;
+ indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
}
-// Nearest point on surface
-Foam::tmp<Foam::pointField> Foam::triSurfaceSearch::calcNearest
+Foam::pointIndexHit Foam::triSurfaceSearch::nearest
(
- const pointField& samples,
+ const point& pt,
const vector& span
-) const
+)
+const
{
const scalar nearestDistSqr = 0.25*magSqr(span);
- tmp<pointField> tnearest(new pointField(samples.size()));
- pointField& nearest = tnearest();
+ return tree().findNearest(pt, nearestDistSqr);
+}
- pointIndexHit hitInfo;
- forAll(samples, sampleI)
- {
- hitInfo = tree().findNearest(samples[sampleI], nearestDistSqr);
+void Foam::triSurfaceSearch::findLine
+(
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>& info
+) const
+{
+ const indexedOctree<treeDataTriSurface>& octree = tree();
- if (hitInfo.hit())
- {
- nearest[sampleI] = hitInfo.hitPoint();
- }
- else
- {
- nearest[sampleI] = greatPoint;
- }
+ info.setSize(start.size());
+
+ scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
+ indexedOctree<treeDataTriSurface>::perturbTol() = tolerance();
+
+ forAll(start, i)
+ {
+ static_cast<pointIndexHit&>(info[i]) = octree.findLine
+ (
+ start[i],
+ end[i]
+ );
}
- return tnearest;
+ indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
}
-Foam::pointIndexHit Foam::triSurfaceSearch::nearest
+void Foam::triSurfaceSearch::findLineAny
(
- const point& pt,
- const vector& span
-)
-const
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>& info
+) const
{
- const scalar nearestDistSqr = 0.25*magSqr(span);
+ const indexedOctree<treeDataTriSurface>& octree = tree();
- return tree().findNearest(pt, nearestDistSqr);
+ info.setSize(start.size());
+
+ scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
+ indexedOctree<treeDataTriSurface>::perturbTol() = tolerance();
+
+ forAll(start, i)
+ {
+ static_cast<pointIndexHit&>(info[i]) = octree.findLineAny
+ (
+ start[i],
+ end[i]
+ );
+ }
+
+ indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
}
void Foam::triSurfaceSearch::findLineAll
(
- const point& start,
- const point& end,
- List<pointIndexHit>& hits
-)
-const
+ const pointField& start,
+ const pointField& end,
+ List<List<pointIndexHit> >& info
+) const
{
- // See if any intersection between pt and end
- pointIndexHit inter = tree().findLine(start, end);
+ const indexedOctree<treeDataTriSurface>& octree = tree();
- if (inter.hit())
- {
- hits.setSize(1);
- hits[0] = inter;
+ info.setSize(start.size());
- const vector dirVec(end-start);
- const scalar magSqrDirVec(magSqr(dirVec));
- const vector smallVec
- (
- indexedOctree<treeDataTriSurface>::perturbTol()*dirVec
- + vector(ROOTVSMALL,ROOTVSMALL,ROOTVSMALL)
- );
+ scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
+ indexedOctree<treeDataTriSurface>::perturbTol() = tolerance();
+ // Work array
+ DynamicList<pointIndexHit, 1, 1> hits;
- // Initial perturbation amount
- vector perturbVec(smallVec);
+ DynamicList<label> shapeMask;
- while (true)
- {
- // Start tracking from last hit.
- point pt = hits.last().hitPoint() + perturbVec;
+ treeDataTriSurface::findAllIntersectOp allIntersectOp(octree, shapeMask);
- if (((pt-start)&dirVec) > magSqrDirVec)
- {
- return;
- }
+ forAll(start, pointI)
+ {
+ hits.clear();
+ shapeMask.clear();
+ while (true)
+ {
// See if any intersection between pt and end
- pointIndexHit inter = tree().findLine(pt, end);
-
- if (!inter.hit())
- {
- return;
- }
+ pointIndexHit inter = octree.findLine
+ (
+ start[pointI],
+ end[pointI],
+ allIntersectOp
+ );
- // Check if already found this intersection
- bool duplicateHit = false;
- forAllReverse(hits, i)
+ if (inter.hit())
{
- if (hits[i].index() == inter.index())
+ vector lineVec = end[pointI] - start[pointI];
+ lineVec /= mag(lineVec) + VSMALL;
+
+ if
+ (
+ checkUniqueHit
+ (
+ inter,
+ hits,
+ lineVec
+ )
+ )
{
- duplicateHit = true;
- break;
+ hits.append(inter);
}
- }
-
- if (duplicateHit)
- {
- // Hit same triangle again. Increase perturbVec and try again.
- perturbVec *= 2;
+ shapeMask.append(inter.index());
}
else
{
- // Proper hit
- label sz = hits.size();
- hits.setSize(sz+1);
- hits[sz] = inter;
- // Restore perturbVec
- perturbVec = smallVec;
+ break;
}
}
+
+ info[pointI].transfer(hits);
}
- else
- {
- hits.clear();
- }
+
+ indexedOctree<treeDataTriSurface>::perturbTol() = oldTol;
}
diff --git a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.H b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.H
index bbf7d1c6563c18af8002991c7a5aca34261beccc..bdca1bca83592031ea3d8cde238312ca337df8b8 100644
--- a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.H
+++ b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -48,11 +48,9 @@ namespace Foam
// Forward declaration of classes
class triSurface;
-class treeDataTriSurface;
-template<class Type> class indexedOctree;
/*---------------------------------------------------------------------------*\
- Class triSurfaceSearch Declaration
+ Class triSurfaceSearch Declaration
\*---------------------------------------------------------------------------*/
class triSurfaceSearch
@@ -62,12 +60,26 @@ class triSurfaceSearch
//- Reference to surface to work on
const triSurface& surface_;
+ //- Optional tolerance to use in searches
+ scalar tolerance_;
+
+ //- Optional max tree depth of octree
+ label maxTreeDepth_;
+
//- Octree for searches
- autoPtr<indexedOctree<treeDataTriSurface> > treePtr_;
+ mutable autoPtr<indexedOctree<treeDataTriSurface> > treePtr_;
// Private Member Functions
+ //- Check whether the current hit on the surface which lies on lineVec
+ // is unique.
+ bool checkUniqueHit
+ (
+ const pointIndexHit& currHit,
+ const DynamicList<pointIndexHit, 1, 1>& hits,
+ const vector& lineVec
+ ) const;
//- Disallow default bitwise copy construct
triSurfaceSearch(const triSurfaceSearch&);
@@ -75,51 +87,64 @@ class triSurfaceSearch
//- Disallow default bitwise assignment
void operator=(const triSurfaceSearch&);
+
public:
- // Static data members
+ // Constructors
- //- Point far away; used for illegal finds
- static const point greatPoint;
+ //- Construct from surface. Holds reference to surface!
+ explicit triSurfaceSearch(const triSurface&);
+ //- Construct from surface and dictionary.
+ triSurfaceSearch(const triSurface&, const dictionary& dict);
- // Constructors
+ //- Construct from components
+ triSurfaceSearch
+ (
+ const triSurface& surface,
+ const scalar tolerance,
+ const label maxTreeDepth
+ );
- //- Construct from surface. Holds reference to surface!
- triSurfaceSearch(const triSurface&);
+ //- Destructor
+ ~triSurfaceSearch();
+
+ //- Clear storage
+ void clearOut();
// Member Functions
- const indexedOctree<treeDataTriSurface>& tree() const
- {
- return treePtr_();
- }
+ //- Demand driven construction of the octree
+ const indexedOctree<treeDataTriSurface>& tree() const;
+ //- Return reference to the surface.
const triSurface& surface() const
{
return surface_;
}
+ //- Return tolerance to use in searches
+ scalar tolerance() const
+ {
+ return tolerance_;
+ }
+
+ //- Return max tree depth of octree
+ label maxTreeDepth() const
+ {
+ return maxTreeDepth_;
+ }
+
//- Calculate for each searchPoint inside/outside status.
boolList calcInside(const pointField& searchPoints) const;
- //- Calculate index of nearest triangle (or -1) for each sample.
- // Looks only in box of size 2*span around sample.
- labelList calcNearestTri
- (
- const pointField& samples,
- const vector& span
- ) const;
-
- //- Calculate nearest points (to searchPoints) on surface.
- // Looks only in box of size 2*span around sample. Returns greatPoint
- // if not found.
- tmp<pointField> calcNearest
+ void findNearest
(
const pointField& samples,
- const vector& span
+ const scalarField& nearestDistSqr,
+ List<pointIndexHit>& info
) const;
//- Calculate nearest point on surface for single searchPoint. Returns
@@ -129,14 +154,27 @@ public:
// - index() : surface triangle label
pointIndexHit nearest(const point&, const vector& span) const;
+ void findLine
+ (
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>& info
+ ) const;
+
+ void findLineAny
+ (
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>& info
+ ) const;
+
//- Calculate all intersections from start to end
void findLineAll
(
- const point& start,
- const point& end,
- List<pointIndexHit>&
+ const pointField& start,
+ const pointField& end,
+ List<List<pointIndexHit> >& info
) const;
-
};
diff --git a/src/meshTools/twoDPointCorrector/twoDPointCorrector.H b/src/meshTools/twoDPointCorrector/twoDPointCorrector.H
index d9e184347e6858ff7de85c82b3768e794f76e7b8..ed5e8dddf27a0cbf8a699e180f8b7df917ff7925 100644
--- a/src/meshTools/twoDPointCorrector/twoDPointCorrector.H
+++ b/src/meshTools/twoDPointCorrector/twoDPointCorrector.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -31,7 +31,7 @@ Description
and thus introduce a third dimension into a 2-D problem.
The operation is performed by looping through all edges approximately
- normal to the plane and enforcing their orthoginality onto the plane by
+ normal to the plane and enforcing their orthogonality onto the plane by
adjusting points on their ends.
SourceFiles
diff --git a/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C b/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C
index f0b13b4afb450df10063ec14c173076e5cefb02a..2dab63202beda8e0088bde24cce5127885dc0aaa 100644
--- a/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C
+++ b/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -37,6 +37,7 @@ License
#include "geomDecomp.H"
#include "vectorList.H"
#include "PackedBoolList.H"
+#include "PatchTools.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -2356,7 +2357,7 @@ void Foam::distributedTriSurfaceMesh::writeStats(Ostream& os) const
{
boundBox bb;
label nPoints;
- calcBounds(bb, nPoints);
+ PatchTools::calcBounds(static_cast<const triSurface&>(*this), bb, nPoints);
reduce(bb.min(), minOp<point>());
reduce(bb.max(), maxOp<point>());
diff --git a/src/sampling/sampledSurface/distanceSurface/distanceSurface.C b/src/sampling/sampledSurface/distanceSurface/distanceSurface.C
index 7d9c4eace454d7cd0bab0c351267530c03ab8cf5..2dafe3e98ad268312e11ca50d68c84a2669b83ec 100644
--- a/src/sampling/sampledSurface/distanceSurface/distanceSurface.C
+++ b/src/sampling/sampledSurface/distanceSurface/distanceSurface.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -29,6 +29,7 @@ License
#include "volPointInterpolation.H"
#include "addToRunTimeSelectionTable.H"
#include "fvMesh.H"
+#include "volumeType.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -94,19 +95,19 @@ void Foam::distanceSurface::createGeometry()
if (signed_)
{
- List<searchableSurface::volumeType> volType;
+ List<volumeType> volType;
surfPtr_().getVolumeType(cc, volType);
forAll(volType, i)
{
- searchableSurface::volumeType vT = volType[i];
+ volumeType vT = volType[i];
- if (vT == searchableSurface::OUTSIDE)
+ if (vT == volumeType::OUTSIDE)
{
fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
}
- else if (vT == searchableSurface::INSIDE)
+ else if (vT == volumeType::INSIDE)
{
fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint());
}
@@ -146,19 +147,19 @@ void Foam::distanceSurface::createGeometry()
if (signed_)
{
- List<searchableSurface::volumeType> volType;
+ List<volumeType> volType;
surfPtr_().getVolumeType(cc, volType);
forAll(volType, i)
{
- searchableSurface::volumeType vT = volType[i];
+ volumeType vT = volType[i];
- if (vT == searchableSurface::OUTSIDE)
+ if (vT == volumeType::OUTSIDE)
{
fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
}
- else if (vT == searchableSurface::INSIDE)
+ else if (vT == volumeType::INSIDE)
{
fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint());
}
@@ -203,20 +204,20 @@ void Foam::distanceSurface::createGeometry()
if (signed_)
{
- List<searchableSurface::volumeType> volType;
+ List<volumeType> volType;
surfPtr_().getVolumeType(pts, volType);
forAll(volType, i)
{
- searchableSurface::volumeType vT = volType[i];
+ volumeType vT = volType[i];
- if (vT == searchableSurface::OUTSIDE)
+ if (vT == volumeType::OUTSIDE)
{
pointDistance_[i] =
Foam::mag(pts[i] - nearest[i].hitPoint());
}
- else if (vT == searchableSurface::INSIDE)
+ else if (vT == volumeType::INSIDE)
{
pointDistance_[i] =
-Foam::mag(pts[i] - nearest[i].hitPoint());