ENH: cvMesh: Add quick rejection tests to octree and buildParallelInterface

applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C

@@ -724,6 +724,8 @@ void Foam::backgroundMeshDecomposition::buildPatchAndTree()
 
     globalBackgroundBounds_ = treeBoundBox(bbMin, bbMax);
 
+    octreeNearestDistances_ = bFTreePtr_().calcNearestDistance();
+
     if (cvMesh_.cvMeshControls().objOutput())
     {
         OFstream fStr
@@ -795,6 +797,7 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
     ),
     boundaryFacesPtr_(),
     bFTreePtr_(),
+    octreeNearestDistances_(),
     allBackgroundMeshBounds_(Pstream::nProcs()),
     globalBackgroundBounds_(),
     decomposeDict_
@@ -1150,6 +1153,8 @@ bool Foam::backgroundMeshDecomposition::positionOnThisProcessor
     const point& pt
 ) const
 {
+//    return bFTreePtr_().findAnyOverlap(pt, 0.0);
+
     return
         bFTreePtr_().getVolumeType(pt)
      == indexedOctree<treeDataBPatch>::INSIDE;
@@ -1176,6 +1181,7 @@ bool Foam::backgroundMeshDecomposition::overlapsThisProcessor
     const treeBoundBox& box
 ) const
 {
+//    return !procBounds().contains(box);
     return !bFTreePtr_().findBox(box).empty();
 }
 
@@ -1183,9 +1189,11 @@ bool Foam::backgroundMeshDecomposition::overlapsThisProcessor
 bool Foam::backgroundMeshDecomposition::overlapsThisProcessor
 (
     const point& centre,
-    scalar radiusSqr
+    const scalar radiusSqr
 ) const
 {
+    //return bFTreePtr_().findAnyOverlap(centre, radiusSqr);
+
     return bFTreePtr_().findNearest(centre, radiusSqr).hit();
 }
 
@@ -1645,6 +1653,7 @@ Foam::labelListList Foam::backgroundMeshDecomposition::overlapsProcessors
 
         // If the sphere finds a nearest element of the patch, then it overlaps
         sphereOverlapsCandidate[sI] = bFTreePtr_().findNearest(c, rSqr).hit();
+        //sphereOverlapsCandidate[sI] = bFTreePtr_().findAnyOverlap(c, rSqr);
     }
 
     map().reverseDistribute

applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H

@@ -111,6 +111,8 @@ class backgroundMeshDecomposition
         //- Search tree for the boundaryFaces_ patch
         autoPtr<indexedOctree<treeDataBPatch> > bFTreePtr_;
 
+        List<scalar> octreeNearestDistances_;
+
         //- The bounds of all background meshes on all processors
         treeBoundBoxList allBackgroundMeshBounds_;
 
@@ -225,16 +227,16 @@ public:
         //- Are the given positions inside the domain of this decomposition
         boolList positionOnThisProcessor(const List<point>& pts) const;
 
-        //- Does the given box overlap the faces of the bounday of this
+        //- Does the given box overlap the faces of the boundary of this
         //  processor
         bool overlapsThisProcessor(const treeBoundBox& box) const;
 
-        //- Does the given sphere overlap the faces of the bounday of this
+        //- Does the given sphere overlap the faces of the boundary of this
         //  processor
         bool overlapsThisProcessor
         (
             const point& centre,
-            scalar radiusSqr
+            const scalar radiusSqr
         ) const;
 
         //- Find nearest intersection of line between start and end, (exposing
@@ -289,6 +291,12 @@ public:
             //- Return access to the underlying mesh
             inline const fvMesh& mesh() const;
 
+            //- Return access to the underlying tree
+            inline const indexedOctree<treeDataBPatch>& tree() const;
+
+            //- Return access to the nearest distance of the octree nodes
+            inline const List<scalar>& octreeNearestDistances() const;
+
             //- Return the boundBox of this processor
             inline const treeBoundBox& procBounds() const;
 };

applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H

@@ -30,6 +30,18 @@ const Foam::fvMesh& Foam::backgroundMeshDecomposition::mesh() const
     return mesh_;
 }
 
+const Foam::indexedOctree<Foam::treeDataBPatch>&
+Foam::backgroundMeshDecomposition::tree() const
+{
+    return bFTreePtr_();
+}
+
+const Foam::List<Foam::scalar>&
+Foam::backgroundMeshDecomposition::octreeNearestDistances() const
+{
+    return octreeNearestDistances_;
+}
+
 
 const Foam::treeBoundBox& Foam::backgroundMeshDecomposition::procBounds() const
 {

src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C

@@ -379,6 +379,124 @@ Foam::label Foam::indexedOctree<Type>::compactContents
 }
 
 
+template <class Type>
+bool Foam::indexedOctree<Type>::quickCircumsphereRejection
+(
+    const label nodeI,
+    const point& cc,
+    const scalar crSqr,
+    const List<scalar>& nearestDistances
+) const
+{
+    const node& nod = nodes_[nodeI];
+
+    volumeType nodeType = volumeType(nodeTypes_.get(nodeI<<3));
+
+    //scalar boxDist = nearestDistances[nodeI] + 0.5*nod.bb_.mag();
+    scalar boxDist = crSqr + magSqr(cc - nod.bb_.midpoint());
+
+    if
+    (
+         nodeType == INSIDE
+      //&& (crSqr < sqr(boxDist))
+      && (boxDist < sqr(nearestDistances[nodeI]))
+    )
+    {
+        return true;
+    }
+    else
+    {
+        direction octant = nod.bb_.subOctant(cc);
+
+        labelBits index = nod.subNodes_[octant];
+
+        if (isNode(index))
+        {
+            return quickCircumsphereRejection
+            (
+                getNode(index),
+                cc,
+                crSqr,
+                nearestDistances
+            );
+        }
+        else
+        {
+            return false;
+        }
+    }
+}
+
+
+template <class Type>
+bool Foam::indexedOctree<Type>::quickCircumsphereRejection
+(
+    const point& cc,
+    const scalar crSqr,
+    const List<scalar>& nearestDistances
+) const
+{
+    if (nodes_.size())
+    {
+        return quickCircumsphereRejection
+        (
+            0,
+            cc,
+            crSqr,
+            nearestDistances
+        );
+    }
+
+    return false;
+}
+
+
+template <class Type>
+Foam::scalar
+Foam::indexedOctree<Type>::calcNearestDistance
+(
+    const label nodeI
+) const
+{
+    const node& nod = nodes_[nodeI];
+
+    const point& nodeCentre = nod.bb_.midpoint();
+
+    scalar nearestDistance = 0.0;
+
+    pointIndexHit pHit = findNearest(nodeCentre, sqr(GREAT));
+
+    if (pHit.hit())
+    {
+        nearestDistance = mag(pHit.hitPoint() - nodeCentre);
+    }
+    else
+    {
+        WarningIn("Foam::indexedOctree<Type>::calcNearestDistance(const label)")
+            << "Cannot calculate distance of nearest point on surface from "
+            << "the midpoint of the octree node. Returning distance of zero."
+            << endl;
+    }
+
+    return nearestDistance;
+}
+
+
+template <class Type>
+Foam::List<Foam::scalar>
+Foam::indexedOctree<Type>::calcNearestDistance() const
+{
+    List<scalar> nearestDistances(nodes_.size());
+
+    forAll(nearestDistances, nodeI)
+    {
+        nearestDistances[nodeI] = calcNearestDistance(nodeI);
+    }
+
+    return nearestDistances;
+}
+
+
 // Pre-calculates wherever possible the volume status per node/subnode.
 // Recurses to determine status of lowest level boxes. Level above is
 // combination of octants below.
@@ -540,6 +658,66 @@ Foam::indexedOctree<Type>::getSide
 // ~~~~~~~~~~~~~~
 //
 
+
+//template <class Type>
+//bool Foam::indexedOctree<Type>::findAnyOverlap
+//(
+//    const label nodeI,
+//    const point& sample,
+//    const scalar nearestDistSqr
+//) const
+//{
+//    const node& nod = nodes_[nodeI];
+//
+//    // Determine order to walk through octants
+//    FixedList<direction, 8> octantOrder;
+//    nod.bb_.searchOrder(sample, octantOrder);
+//
+//    // Go into all suboctants (one containing sample first) and update nearest.
+//    for (direction i = 0; i < 8; i++)
+//    {
+//        direction octant = octantOrder[i];
+//
+//        labelBits index = nod.subNodes_[octant];
+//
+//        if (isNode(index))
+//        {
+//            label subNodeI = getNode(index);
+//
+//            const treeBoundBox& subBb = nodes_[subNodeI].bb_;
+//
+//            if (overlaps(subBb.min(), subBb.max(), nearestDistSqr, sample))
+//            {
+//                return findAnyOverlap
+//                (
+//                    subNodeI,
+//                    sample,
+//                    nearestDistSqr
+//                );
+//            }
+//        }
+//        else if (isContent(index))
+//        {
+//            if
+//            (
+//                overlaps
+//                (
+//                    nod.bb_,
+//                    octant,
+//                    nearestDistSqr,
+//                    sample
+//                )
+//            )
+//            {
+//                return true;
+//            }
+//        }
+//    }
+//
+//    return false;
+//}
+
+
 // Find nearest point starting from nodeI
 template <class Type>
 void Foam::indexedOctree<Type>::findNearest
@@ -1614,7 +1792,6 @@ void Foam::indexedOctree<Type>::traverseNode
         }
     }
 
-
     const node& nod = nodes_[nodeI];
 
     labelBits index = nod.subNodes_[octant];
@@ -1781,6 +1958,11 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
     label i = 0;
     for (; i < 100000; i++)
     {
+//        if (isLineInsideOrOutside(nodeI, treeStart, treeEnd))
+//        {
+//            return hitInfo;
+//        }
+
         // Ray-trace to end of current node. Updates point (either on triangle
         // in case of hit or on node bounding box in case of miss)
 
@@ -1935,6 +2117,38 @@ Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
 }
 
 
+//template <class Type>
+//bool Foam::indexedOctree<Type>::isLineInsideOrOutside
+//(
+//    const label nodeI,
+//    const point& start,
+//    const point& end
+//) const
+//{
+//    const node& nod = nodes_[nodeI];
+//
+//    direction startOctant = nod.bb_.subOctant(start);
+//    direction endOctant = nod.bb_.subOctant(end);
+//
+//    if (startOctant == endOctant)
+//    {
+//        volumeType startOctantType
+//            = volumeType(nodeTypes_.get((nodeI<<3) + startOctant));
+//
+//        if
+//        (
+//            startOctantType == INSIDE || startOctantType == OUTSIDE
+//        )
+//        {
+//            //Info<< nodeI << " | " << start << " " << end << endl;
+//            return true;
+//        }
+//    }
+//
+//    return false;
+//}
+
+
 // Find first intersection
 template <class Type>
 Foam::pointIndexHit Foam::indexedOctree<Type>::findLine
@@ -2559,6 +2773,27 @@ Foam::scalar& Foam::indexedOctree<Type>::perturbTol()
 }
 
 
+//template <class Type>
+//bool Foam::indexedOctree<Type>::findAnyOverlap
+//(
+//    const point& sample,
+//    const scalar startDistSqr
+//) const
+//{
+//    if (nodes_.size())
+//    {
+//        return findAnyOverlap
+//        (
+//            0,
+//            sample,
+//            startDistSqr
+//        );
+//    }
+//
+//    return false;
+//}
+
+
 template <class Type>
 Foam::pointIndexHit Foam::indexedOctree<Type>::findNearest
 (

src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H

@@ -201,6 +201,16 @@ private:
                 label& compactI
             );
 
+            scalar calcNearestDistance(const label nodeI) const;
+
+            bool quickCircumsphereRejection
+            (
+                const label nodeI,
+                const point& cc,
+                const scalar crSqr,
+                const List<scalar>& nearestDistances
+            ) const;
+
             //- Determine inside/outside per node (mixed if cannot be
             //  determined). Only valid for closed shapes.
             volumeType calcVolumeType(const label nodeI) const;
@@ -320,6 +330,13 @@ private:
                 const bool verbose = false
             ) const;
 
+//            bool isLineInsideOrOutside
+//            (
+//                const label nodeI,
+//                const point& start,
+//                const point& end
+//            ) const;
+
             //- Find any or nearest intersection of line between start and end.
             pointIndexHit findLine
             (
@@ -532,6 +549,19 @@ public:
                 const scalar nearestDistSqr
             ) const;
 
+//            bool findAnyOverlap
+//            (
+//                const point& sample,
+//                const scalar nearestDistSqr
+//            ) const;
+//
+//            bool findAnyOverlap
+//            (
+//                const label nodeI,
+//                const point& sample,
+//                const scalar nearestDistSqr
+//            ) const;
+
             //- Low level: calculate nearest starting from subnode.
             void findNearest
             (
@@ -628,6 +658,16 @@ public:
                 CompareOp& cop
             ) const;
 
+            //- Return a list containing the nearest distance of nodes to any
+            //  shapes
+            List<scalar> calcNearestDistance() const;
+
+            bool quickCircumsphereRejection
+            (
+                const point& cc,
+                const scalar crSqr,
+                const List<scalar>& nearestDistances
+            ) const;
 
         // Write