Use searchtree in STLBoundary

This commit is contained in:
Matthias Hochsteger 2019-10-01 11:23:22 +02:00
parent ef64a5e7eb
commit 7ced41e56f
5 changed files with 97 additions and 103 deletions

View File

@ -741,97 +741,98 @@ public:
const ADTree3 & Tree() const { return *tree; }; const ADTree3 & Tree() const { return *tree; };
}; };
template<int dim, typename T=INDEX> template<int dim, typename T=INDEX>
class BoxTree class BoxTree
{
public:
// Number of entries per leaf
static constexpr int N = 100;
struct Node;
struct Leaf
{ {
public: Point<2*dim> p[N];
// Number of entries per leaf T index[N];
static constexpr int N = 100; int n_elements;
struct Node; Leaf() : n_elements(0)
{ }
struct Leaf void Add( ClosedHashTable<T, Leaf*> &leaf_index, const Point<2*dim> &ap, T aindex )
{
Point<2*dim> p[N];
T index[N];
int n_elements;
Leaf() : n_elements(0) {}
void Add( ClosedHashTable<T, Leaf*> &leaf_index, const Point<2*dim> &ap, T aindex )
{ {
p[n_elements] = ap; p[n_elements] = ap;
index[n_elements] = aindex; index[n_elements] = aindex;
n_elements++; n_elements++;
leaf_index[aindex] = this; leaf_index[aindex] = this;
} }
}; };
struct Node struct Node
{
union
{ {
union Node *children[2];
{ Leaf *leaf;
Node *children[2];
Leaf *leaf;
};
double sep;
int level;
Node()
: children{nullptr,nullptr}
{ }
~Node()
{ }
Leaf *GetLeaf() const
{
return children[1] ? nullptr : leaf;
}
}; };
double sep;
int level;
private: Node()
Node root; : children{nullptr,nullptr}
// Array<Leaf*, INDEX> leaf_index;
ClosedHashTable<T, Leaf*> leaf_index;
Point<dim> global_min, global_max;
double tol;
size_t n_leaves;
size_t n_nodes;
BlockAllocator ball_nodes;
BlockAllocator ball_leaves;
public:
BoxTree (const Box<dim> & abox)
: BoxTree( abox.PMin(), abox.PMax() )
{ } { }
BoxTree (const Point<dim> & pmin, const Point<dim> & pmax) ~Node()
: global_min(pmin), global_max(pmax), n_leaves(1), n_nodes(1), ball_nodes(sizeof(Node)), ball_leaves(sizeof(Leaf)) { }
Leaf *GetLeaf() const
{
return children[1] ? nullptr : leaf;
}
};
private:
Node root;
ClosedHashTable<T, Leaf*> leaf_index;
Point<dim> global_min, global_max;
double tol;
size_t n_leaves;
size_t n_nodes;
BlockAllocator ball_nodes;
BlockAllocator ball_leaves;
public:
BoxTree (const Point<dim> & pmin, const Point<dim> & pmax)
: global_min(pmin), global_max(pmax), n_leaves(1), n_nodes(1), ball_nodes(sizeof(Node)), ball_leaves(sizeof(Leaf))
{ {
root.leaf = (Leaf*) ball_leaves.Alloc(); new (root.leaf) Leaf(); root.leaf = (Leaf*) ball_leaves.Alloc(); new (root.leaf) Leaf();
root.level = 0; root.level = 0;
tol = 1e-7 * Dist(pmax, pmin); tol = 1e-7 * Dist(pmax, pmin);
} }
size_t GetNLeaves() BoxTree (const Box<dim> & box)
: BoxTree(box.PMin(), box.PMax())
{ }
size_t GetNLeaves()
{ {
return n_leaves; return n_leaves;
} }
size_t GetNNodes() size_t GetNNodes()
{ {
return n_nodes; return n_nodes;
} }
template<typename TFunc> template<typename TFunc>
void GetFirstIntersecting (const Point<dim> & pmin, const Point<dim> & pmax, void GetFirstIntersecting (const Point<dim> & pmin, const Point<dim> & pmax,
TFunc func=[](auto pi){return false;}) const TFunc func=[](auto pi){return false;}) const
{ {
// static Timer timer("BTree::GetIntersecting"); RegionTimer rt(timer); // static Timer timer("BoxTree::GetIntersecting"); RegionTimer rt(timer);
// static Timer timer1("BTree::GetIntersecting-LinearSearch"); // static Timer timer1("BoxTree::GetIntersecting-LinearSearch");
ArrayMem<const Node*, 100> stack; ArrayMem<const Node*, 100> stack;
ArrayMem<int, 100> dir_stack; ArrayMem<int, 100> dir_stack;
@ -842,7 +843,7 @@ public:
{ {
tpmin(i) = global_min(i); tpmin(i) = global_min(i);
tpmax(i) = pmax(i)+tol; tpmax(i) = pmax(i)+tol;
tpmin(i+dim) = pmin(i)-tol; tpmin(i+dim) = pmin(i)-tol;
tpmax(i+dim) = global_max(i); tpmax(i+dim) = global_max(i);
} }
@ -862,7 +863,7 @@ public:
if(Leaf *leaf = node->GetLeaf()) if(Leaf *leaf = node->GetLeaf())
{ {
// RegionTimer rt1(timer1); // RegionTimer rt1(timer1);
for (auto i : IntRange(leaf->n_elements)) for (auto i : IntRange(leaf->n_elements))
{ {
bool intersect = true; bool intersect = true;
@ -875,7 +876,7 @@ public:
if (p[d] < tpmin[d]) if (p[d] < tpmin[d])
intersect = false; intersect = false;
if(intersect) if(intersect)
if(func(leaf->index[i])) return; if(func(leaf->index[i])) return;
} }
} }
else else
@ -896,21 +897,21 @@ public:
} }
} }
void GetIntersecting (const Point<dim> & pmin, const Point<dim> & pmax, void GetIntersecting (const Point<dim> & pmin, const Point<dim> & pmax,
NgArray<T> & pis) const NgArray<T> & pis) const
{ {
pis.SetSize(0); pis.SetSize(0);
GetFirstIntersecting(pmin, pmax, [&pis](auto pi) { pis.Append(pi); return false;}); GetFirstIntersecting(pmin, pmax, [&pis](auto pi) { pis.Append(pi); return false;});
} }
void Insert (const Box<dim> & box, T pi) void Insert (const Box<dim> & box, T pi)
{ {
Insert (box.PMin(), box.PMax(), pi); Insert (box.PMin(), box.PMax(), pi);
} }
void Insert (const Point<dim> & pmin, const Point<dim> & pmax, T pi) void Insert (const Point<dim> & pmin, const Point<dim> & pmax, T pi)
{ {
// static Timer timer("BTree::Insert"); RegionTimer rt(timer); // static Timer timer("BoxTree::Insert"); RegionTimer rt(timer);
int dir = 0; int dir = 0;
Point<2*dim> p; Point<2*dim> p;
for (auto i : IntRange(dim)) for (auto i : IntRange(dim))
@ -982,9 +983,9 @@ public:
} }
} }
void DeleteElement (T pi) void DeleteElement (T pi)
{ {
// static Timer timer("BTree::DeleteElement"); RegionTimer rt(timer); // static Timer timer("BoxTree::DeleteElement"); RegionTimer rt(timer);
Leaf *leaf = leaf_index[pi]; Leaf *leaf = leaf_index[pi];
leaf_index.Delete(pi); leaf_index.Delete(pi);
auto & n_elements = leaf->n_elements; auto & n_elements = leaf->n_elements;
@ -1005,7 +1006,7 @@ public:
} }
} }
} }
}; };
// template <int dim, typename T = INDEX> // template <int dim, typename T = INDEX>
// class BoxTree // class BoxTree
@ -1025,7 +1026,7 @@ public:
// } // }
// tree = new T_ADTree<2*dim,T> (tpmin, tpmax); // tree = new T_ADTree<2*dim,T> (tpmin, tpmax);
// } // }
// //
// BoxTree (const Point<dim> & apmin, const Point<dim> & apmax) // BoxTree (const Point<dim> & apmin, const Point<dim> & apmax)
// { // {
// boxpmin = apmin; // boxpmin = apmin;
@ -1038,36 +1039,36 @@ public:
// } // }
// tree = new T_ADTree<2*dim,T> (tpmin, tpmax); // tree = new T_ADTree<2*dim,T> (tpmin, tpmax);
// } // }
// //
// ~BoxTree () // ~BoxTree ()
// { // {
// delete tree; // delete tree;
// } // }
// //
// void Insert (const Point<dim> & bmin, const Point<dim> & bmax, T pi) // void Insert (const Point<dim> & bmin, const Point<dim> & bmax, T pi)
// { // {
// Point<2*dim> tp; // Point<2*dim> tp;
// //
// for (size_t i = 0; i < dim; i++) // for (size_t i = 0; i < dim; i++)
// { // {
// tp(i) = bmin(i); // tp(i) = bmin(i);
// tp(i+dim) = bmax(i); // tp(i+dim) = bmax(i);
// } // }
// //
// tree->Insert (tp, pi); // tree->Insert (tp, pi);
// } // }
// //
// void Insert (const Box<dim> & box, T pi) // void Insert (const Box<dim> & box, T pi)
// { // {
// Insert (box.PMin(), box.PMax(), pi); // Insert (box.PMin(), box.PMax(), pi);
// } // }
// //
// void DeleteElement (T pi) // void DeleteElement (T pi)
// { // {
// tree->DeleteElement(pi); // tree->DeleteElement(pi);
// } // }
// //
// void GetIntersecting (const Point<dim> & pmin, const Point<dim> & pmax, // void GetIntersecting (const Point<dim> & pmin, const Point<dim> & pmax,
// NgArray<T> & pis) const // NgArray<T> & pis) const
// { // {
// Point<2*dim> tpmin, tpmax; // Point<2*dim> tpmin, tpmax;
@ -1076,15 +1077,15 @@ public:
// { // {
// tpmin(i) = boxpmin(i); // tpmin(i) = boxpmin(i);
// tpmax(i) = pmax(i)+tol; // tpmax(i) = pmax(i)+tol;
// //
// tpmin(i+dim) = pmin(i)-tol; // tpmin(i+dim) = pmin(i)-tol;
// tpmax(i+dim) = boxpmax(i); // tpmax(i+dim) = boxpmax(i);
// } // }
// //
// tree->GetIntersecting (tpmin, tpmax, pis); // tree->GetIntersecting (tpmin, tpmax, pis);
// } // }
// //
// //
// double Tolerance() const { return 1e-7 * Dist(boxpmax, boxpmin); } // single precision // double Tolerance() const { return 1e-7 * Dist(boxpmax, boxpmin); } // single precision
// const auto & Tree() const { return *tree; }; // const auto & Tree() const { return *tree; };
// auto & Tree() { return *tree; }; // auto & Tree() { return *tree; };

View File

@ -5,10 +5,6 @@
namespace netgen namespace netgen
{ {
// typedef BTree<3> TBoxTree;
typedef BoxTree<3> TBoxTree;
static const int deltetfaces[][3] = static const int deltetfaces[][3] =
{ { 1, 2, 3 }, { { 1, 2, 3 },
{ 2, 0, 3 }, { 2, 0, 3 },
@ -231,7 +227,7 @@ namespace netgen
void AddDelaunayPoint (PointIndex newpi, const Point3d & newp, void AddDelaunayPoint (PointIndex newpi, const Point3d & newp,
NgArray<DelaunayTet> & tempels, NgArray<DelaunayTet> & tempels,
Mesh & mesh, Mesh & mesh,
TBoxTree & tettree, BoxTree<3> & tettree,
MeshNB & meshnb, MeshNB & meshnb,
NgArray<Point<3> > & centers, NgArray<double> & radi2, NgArray<Point<3> > & centers, NgArray<double> & radi2,
NgArray<int> & connected, NgArray<int> & treesearch, NgArray<int> & connected, NgArray<int> & treesearch,
@ -639,7 +635,7 @@ namespace netgen
pmin2 = pmin2 + 0.1 * (pmin2 - pmax2); pmin2 = pmin2 + 0.1 * (pmin2 - pmax2);
pmax2 = pmax2 + 0.1 * (pmax2 - pmin2); pmax2 = pmax2 + 0.1 * (pmax2 - pmin2);
TBoxTree tettree(pmin2, pmax2); BoxTree<3> tettree(pmin2, pmax2);
tempels.Append (startel); tempels.Append (startel);

View File

@ -146,8 +146,6 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
chartbound.Clear(); chartbound.Clear();
chartbound.SetChart(&chart); chartbound.SetChart(&chart);
chartbound.BuildSearchTree(); // different !!!
if (!found) throw Exception("Make Atlas, no starttrig found"); if (!found) throw Exception("Make Atlas, no starttrig found");
//find surrounding trigs //find surrounding trigs
@ -324,7 +322,6 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
innerchartpts.SetSize(innerchartpoints.Size()); innerchartpts.SetSize(innerchartpoints.Size());
for (size_t i = 0; i < innerchartpoints.Size(); i++) for (size_t i = 0; i < innerchartpoints.Size(); i++)
innerchartpts[i] = GetPoint(innerchartpoints[i]); innerchartpts[i] = GetPoint(innerchartpoints[i]);
// chartbound.BuildSearchTree(); // different !!!
// NgProfiler::StopTimer (timer2); // NgProfiler::StopTimer (timer2);
// NgProfiler::StartTimer (timer3); // NgProfiler::StartTimer (timer3);

View File

@ -1095,6 +1095,9 @@ void STLBoundary ::AddTriangle(const STLTriangle & t)
segs[1] = INDEX_2(t[1], t[2]); segs[1] = INDEX_2(t[1], t[2]);
segs[2] = INDEX_2(t[2], t[0]); segs[2] = INDEX_2(t[2], t[0]);
if(!searchtree)
BuildSearchTree();
for (auto seg : segs) for (auto seg : segs)
{ {
STLBoundarySeg bseg(seg[0], seg[1], geometry->GetPoints(), chart); STLBoundarySeg bseg(seg[0], seg[1], geometry->GetPoints(), chart);
@ -1312,21 +1315,18 @@ bool STLBoundary :: TestSeg(const Point<3>& p1, const Point<3> & p2, const Vec<3
void STLBoundary :: BuildSearchTree() void STLBoundary :: BuildSearchTree()
{ {
delete searchtree;
Box<2> box2d(Box<2>::EMPTY_BOX); Box<2> box2d(Box<2>::EMPTY_BOX);
Box<3> box3d = geometry->GetBoundingBox(); Box<3> box3d = geometry->GetBoundingBox();
for (size_t i = 0; i < 8; i++) for (size_t i = 0; i < 8; i++)
box2d.Add ( chart->Project2d (box3d.GetPointNr(i))); box2d.Add ( chart->Project2d (box3d.GetPointNr(i)));
// comment to enable searchtree: searchtree = make_unique<BoxTree<2,INDEX_2>> (box2d);
// searchtree = new BoxTree<2,INDEX_2> (box2d); // searchtree = nullptr;
searchtree = nullptr;
} }
void STLBoundary :: DeleteSearchTree() void STLBoundary :: DeleteSearchTree()
{ {
delete searchtree;
searchtree = nullptr; searchtree = nullptr;
} }

View File

@ -191,10 +191,10 @@ private:
const STLChart * chart; const STLChart * chart;
// NgArray<STLBoundarySeg> boundary; // NgArray<STLBoundarySeg> boundary;
ClosedHashTable<INDEX_2, STLBoundarySeg> boundary_ht; ClosedHashTable<INDEX_2, STLBoundarySeg> boundary_ht;
BoxTree<2,INDEX_2> * searchtree = nullptr; unique_ptr<BoxTree<2,INDEX_2>> searchtree;
public: public:
STLBoundary(STLGeometry * ageometry); STLBoundary(STLGeometry * ageometry);
~STLBoundary() { delete searchtree; } ~STLBoundary() {}
void Clear() { /* boundary.SetSize(0); */ boundary_ht = ClosedHashTable<INDEX_2,STLBoundarySeg>(); } void Clear() { /* boundary.SetSize(0); */ boundary_ht = ClosedHashTable<INDEX_2,STLBoundarySeg>(); }
void SetChart (const STLChart * achart) { chart = achart; } void SetChart (const STLChart * achart) { chart = achart; }