ForkMaster/netgen
1
0
Fork 0
mirror of https://github.com/NGSolve/netgen.git synced 2024-12-26 05:50:32 +05:00
Code Issues Packages Projects Releases Wiki Activity

delaunay optimizations

Browse Source
This commit is contained in:
Joachim Schöberl 2019-09-26 18:28:31 +02:00
parent 46116ce723
commit fa3b84038f
2 changed files with 52 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

55
libsrc/meshing/delaunay.cpp
View File

@ -173,7 +173,7 @@ namespace netgen
void Insert (const Point<dim> & pmin, const Point<dim> & pmax, T pi)
{
static Timer timer("BTree::Insert"); RegionTimer rt(timer);
// static Timer timer("BTree::Insert"); RegionTimer rt(timer);
int dir = 0;
Point<2*dim> p;
for (auto i : IntRange(dim))
@ -247,7 +247,7 @@ namespace netgen
void DeleteElement (T pi)
{
static Timer timer("BTree::DeleteElement"); RegionTimer rt(timer);
// static Timer timer("BTree::DeleteElement"); RegionTimer rt(timer);
Leaf *leaf = leaf_index[pi];
auto & n_elements = leaf->n_elements;
auto & index = leaf->index;
@ -487,11 +487,18 @@ namespace netgen
NgArray<int> & freelist, SphereList & list,
IndexSet & insphere, IndexSet & closesphere)
{
// static Timer t("Meshing3::AddDelaunayPoint"); RegionTimer reg(t);
static Timer t("Meshing3::AddDelaunayPoint"); RegionTimer reg(t);
static Timer tsearch("addpoint, search");
static Timer tinsert("addpoint, insert");
static Timer t0("addpoint, 0");
static Timer t1("addpoint, 1");
static Timer t2("addpoint, 2");
static Timer t3("addpoint, 3");
/*
find any sphere, such that newp is contained in
*/
DelaunayTet el;
// DelaunayTet el;
int cfelind = -1;
const Point<3> * pp[4];
@ -500,6 +507,7 @@ namespace netgen
/*
// stop search if intersecting point is close enough to center
tettree.GetFirstIntersecting (newp, newp, treesearch, [&](const auto pi)
{
@ -523,7 +531,28 @@ namespace netgen
break;
}
}
*/
tsearch.Start();
double minquot{1e20};
tettree.GetFirstIntersecting
(newp, newp, treesearch, [&](const auto pi)
{
double quot = Dist2 (centers.Get(pi), newp) / radi2.Get(pi);
if (quot < 0.917632)
{
cfelind = pi;
return true;
}
if((cfelind == -1 || quot < 0.99*minquot) && quot < 1)
{
minquot = quot;
cfelind = pi;
}
return false;
} );
tsearch.Stop();
if (cfelind == -1)
{
@ -549,6 +578,7 @@ namespace netgen
int changed = 1;
int nstarti = 1, starti;
t0.Start();
while (changed)
{
changed = 0;
@ -637,8 +667,11 @@ namespace netgen
}
} // while (changed)
t0.Stop();
t1.Start();
// NgArray<Element> newels;
NgArray<DelaunayTet> newels;
static NgArray<DelaunayTet> newels;
newels.SetSize(0);
Element2d face(TRIG);
@ -684,6 +717,8 @@ namespace netgen
}
}
t1.Stop();
t2.Start();
meshnb.ResetFaceHT (10*insphere.GetArray().Size()+1);
for (auto celind : insphere.GetArray())
@ -707,7 +742,8 @@ namespace netgen
fabs (Dist2 (centers.Get (ind), newp) - radi2.Get(ind)) < 1e-8 )
hasclose = true;
}
t2.Stop();
t3.Start();
for (int j = 1; j <= newels.Size(); j++)
{
const auto & newel = newels.Get(j);
@ -781,8 +817,11 @@ namespace netgen
tpmax.SetToMax (*pp[k]);
}
tpmax = tpmax + 0.01 * (tpmax - tpmin);
// tinsert.Start();
tettree.Insert (tpmin, tpmax, nelind);
// tinsert.Stop();
}
t3.Stop();
}
@ -893,7 +932,9 @@ namespace netgen
// "random" reordering of points (speeds a factor 3 - 5 !!!)
NgArray<PointIndex, PointIndex::BASE, PointIndex> mixed(np);
int prims[] = { 11, 13, 17, 19, 23, 29, 31, 37 };
// int prims[] = { 11, 13, 17, 19, 23, 29, 31, 37 };
// int prims[] = { 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97 };
int prims[] = { 211, 223, 227, 229, 233, 239, 241, 251, 257, 263 };
int prim;
{

2
libsrc/meshing/meshing2.cpp
View File

@ -506,7 +506,7 @@ namespace netgen
cout << "set debugflag" << endl;
}
if (debugparam.haltlargequalclass && qualclass > 50)
if (debugparam.haltlargequalclass && qualclass == 50)
debugflag = 1;
// problem recognition !