diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C index edafdb50f9..a2d44f0ee1 100644 --- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C +++ b/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::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 diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H index 87ab88af35..d8e6388cfb 100644 --- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H +++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H @@ -111,6 +111,8 @@ class backgroundMeshDecomposition //- Search tree for the boundaryFaces_ patch autoPtr > bFTreePtr_; + List 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& 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& tree() const; + + //- Return access to the nearest distance of the octree nodes + inline const List& octreeNearestDistances() const; + //- Return the boundBox of this processor inline const treeBoundBox& procBounds() const; }; diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H index df659cfd9c..2897d3fc58 100644 --- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H +++ b/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::backgroundMeshDecomposition::tree() const +{ + return bFTreePtr_(); +} + +const Foam::List& +Foam::backgroundMeshDecomposition::octreeNearestDistances() const +{ + return octreeNearestDistances_; +} + const Foam::treeBoundBox& Foam::backgroundMeshDecomposition::procBounds() const { diff --git a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C index 17ef34fb6f..bfde44fcf9 100644 --- a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C +++ b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.C @@ -379,6 +379,124 @@ Foam::label Foam::indexedOctree::compactContents } +template +bool Foam::indexedOctree::quickCircumsphereRejection +( + const label nodeI, + const point& cc, + const scalar crSqr, + const List& 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 +bool Foam::indexedOctree::quickCircumsphereRejection +( + const point& cc, + const scalar crSqr, + const List& nearestDistances +) const +{ + if (nodes_.size()) + { + return quickCircumsphereRejection + ( + 0, + cc, + crSqr, + nearestDistances + ); + } + + return false; +} + + +template +Foam::scalar +Foam::indexedOctree::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::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 +Foam::List +Foam::indexedOctree::calcNearestDistance() const +{ + List 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::getSide // ~~~~~~~~~~~~~~ // + +//template +//bool Foam::indexedOctree::findAnyOverlap +//( +// const label nodeI, +// const point& sample, +// const scalar nearestDistSqr +//) const +//{ +// const node& nod = nodes_[nodeI]; +// +// // Determine order to walk through octants +// FixedList 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 void Foam::indexedOctree::findNearest @@ -1614,7 +1792,6 @@ void Foam::indexedOctree::traverseNode } } - const node& nod = nodes_[nodeI]; labelBits index = nod.subNodes_[octant]; @@ -1781,6 +1958,11 @@ Foam::pointIndexHit Foam::indexedOctree::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::findLine } +//template +//bool Foam::indexedOctree::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 Foam::pointIndexHit Foam::indexedOctree::findLine @@ -2559,6 +2773,27 @@ Foam::scalar& Foam::indexedOctree::perturbTol() } +//template +//bool Foam::indexedOctree::findAnyOverlap +//( +// const point& sample, +// const scalar startDistSqr +//) const +//{ +// if (nodes_.size()) +// { +// return findAnyOverlap +// ( +// 0, +// sample, +// startDistSqr +// ); +// } +// +// return false; +//} + + template Foam::pointIndexHit Foam::indexedOctree::findNearest ( diff --git a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H b/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H index aa493f34e9..40ee7a8a32 100644 --- a/src/OpenFOAM/algorithms/indexedOctree/indexedOctree.H +++ b/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& 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 calcNearestDistance() const; + + bool quickCircumsphereRejection + ( + const point& cc, + const scalar crSqr, + const List& nearestDistances + ) const; // Write