From 3e7964f18a13e2c34acae450b044d0a6ae605ae1 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 3 Jan 2014 15:43:45 +0000
Subject: [PATCH] ENH: refinementFeatures: searching for nearest region edge
---
.../refinementFeatures/refinementFeatures.C | 116 +++++++++++++++++-
.../refinementFeatures/refinementFeatures.H | 23 +++-
2 files changed, 137 insertions(+), 2 deletions(-)
diff --git a/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.C b/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.C
index b2b974de02f..320f3b3d1c3 100644
--- a/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.C
+++ b/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -313,6 +313,58 @@ void Foam::refinementFeatures::buildTrees(const label featI)
}
+const Foam::PtrList<Foam::indexedOctree<Foam::treeDataEdge> >&
+Foam::refinementFeatures::regionEdgeTrees() const
+{
+ if (!regionEdgeTreesPtr_.valid())
+ {
+ regionEdgeTreesPtr_.reset
+ (
+ new PtrList<indexedOctree<treeDataEdge> >(size())
+ );
+ PtrList<indexedOctree<treeDataEdge> >& trees = regionEdgeTreesPtr_();
+
+ forAll(*this, featI)
+ {
+ const extendedEdgeMesh& eMesh = operator[](featI);
+ const pointField& points = eMesh.points();
+ const edgeList& edges = eMesh.edges();
+
+ // Calculate bb of all points
+ treeBoundBox bb(points);
+
+ // 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);
+
+ trees.set
+ (
+ featI,
+ new indexedOctree<treeDataEdge>
+ (
+ treeDataEdge
+ (
+ false, // do not cache bb
+ edges,
+ points,
+ eMesh.regionEdges()
+ ),
+ bb, // overall search domain
+ 8, // maxLevel
+ 10, // leafsize
+ 3.0 // duplicity
+ )
+ );
+ }
+ }
+ return regionEdgeTreesPtr_();
+}
+
// Find maximum level of a shell.
void Foam::refinementFeatures::findHigherLevel
@@ -543,6 +595,68 @@ void Foam::refinementFeatures::findNearestEdge
}
+void Foam::refinementFeatures::findNearestRegionEdge
+(
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ labelList& nearFeature,
+ List<pointIndexHit>& nearInfo,
+ vectorField& nearNormal
+) const
+{
+ nearFeature.setSize(samples.size());
+ nearFeature = -1;
+ nearInfo.setSize(samples.size());
+ nearInfo = pointIndexHit();
+ nearNormal.setSize(samples.size());
+ nearNormal = vector::zero;
+
+
+ const PtrList<indexedOctree<treeDataEdge> >& regionTrees =
+ regionEdgeTrees();
+
+ forAll(regionTrees, featI)
+ {
+ const indexedOctree<treeDataEdge>& regionTree = regionTrees[featI];
+
+ forAll(samples, sampleI)
+ {
+ const point& sample = samples[sampleI];
+
+ scalar distSqr;
+ if (nearInfo[sampleI].hit())
+ {
+ distSqr = magSqr(nearInfo[sampleI].hitPoint()-sample);
+ }
+ else
+ {
+ distSqr = nearestDistSqr[sampleI];
+ }
+
+ // Find anything closer than current best
+ pointIndexHit info = regionTree.findNearest(sample, distSqr);
+
+ if (info.hit())
+ {
+ const treeDataEdge& td = regionTree.shapes();
+
+ nearFeature[sampleI] = featI;
+ nearInfo[sampleI] = pointIndexHit
+ (
+ info.hit(),
+ info.hitPoint(),
+ regionTree.shapes().edgeLabels()[info.index()]
+ );
+
+ const edge& e = td.edges()[nearInfo[sampleI].index()];
+ nearNormal[sampleI] = e.vec(td.points());
+ nearNormal[sampleI] /= mag(nearNormal[sampleI])+VSMALL;
+ }
+ }
+ }
+}
+
+
//void Foam::refinementFeatures::findNearestPoint
//(
// const pointField& samples,
diff --git a/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.H b/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.H
index 7d27090b147..5c923508b8b 100644
--- a/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.H
+++ b/src/mesh/autoMesh/autoHexMesh/refinementFeatures/refinementFeatures.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -69,6 +69,10 @@ private:
//- Features points
PtrList<indexedOctree<treeDataPoint> > pointTrees_;
+ //- Region edge trees (demand driven)
+ mutable autoPtr<PtrList<indexedOctree<treeDataEdge> > >
+ regionEdgeTreesPtr_;
+
// Private Member Functions
@@ -100,6 +104,8 @@ protected:
}
+ const PtrList<indexedOctree<treeDataEdge> >& regionEdgeTrees() const;
+
public:
// Constructors
@@ -149,6 +155,21 @@ public:
vectorField& nearNormal
) const;
+ //- Find nearest point on nearest region edge. Sets
+ // - nearFeature: index of feature mesh
+ // - nearInfo : location on feature edge and edge index
+ // (note: not feature edge index but index into
+ // edges() directly)
+ // - nearNormal : local feature edge normal
+ void findNearestRegionEdge
+ (
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ labelList& nearFeature,
+ List<pointIndexHit>& nearInfo,
+ vectorField& nearNormal
+ ) const;
+
//- Find nearest feature point. Sets
// - nearFeature: index of feature mesh
// - nearInfo : location on feature point and point index.
--
GitLab