take care of tolerance in searchtree

This commit is contained in:
Joachim Schöberl 2019-02-05 09:02:31 +01:00
parent 1303e92379
commit e42f81b5d5
3 changed files with 31 additions and 17 deletions

View File

@ -2400,13 +2400,13 @@ namespace netgen
Array<T> & pis) const Array<T> & pis) const
{ {
Point<2*dim> tpmin, tpmax; Point<2*dim> tpmin, tpmax;
double tol = Tolerance();
for (size_t i = 0; i < dim; i++) for (size_t i = 0; i < dim; i++)
{ {
tpmin(i) = boxpmin(i); tpmin(i) = boxpmin(i);
tpmax(i) = pmax(i); tpmax(i) = pmax(i)+tol;
tpmin(i+dim) = pmin(i); tpmin(i+dim) = pmin(i)-tol;
tpmax(i+dim) = boxpmax(i); tpmax(i+dim) = boxpmax(i);
} }

View File

@ -570,9 +570,9 @@ public:
{ 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,
Array<T> & pis) const; Array<T> & pis) const;
double Tolerance() const { return 1e-7 * Dist(boxpmax, boxpmin); } // single precision
// const T_ADTree<2*dim> & Tree() const { return *tree; }; const auto & Tree() const { return *tree; };
// T_ADTree<2*dim> & Tree() { return *tree; }; auto & Tree() { return *tree; };
}; };
} }

View File

@ -1680,6 +1680,7 @@ namespace netgen
void Mesh :: CalcSurfacesOfNode () void Mesh :: CalcSurfacesOfNode ()
{ {
static Timer t("Mesh::CalcSurfacesOfNode"); RegionTimer reg (t);
// surfacesonnode.SetSize (GetNP()); // surfacesonnode.SetSize (GetNP());
TABLE<int,PointIndex::BASE> surfacesonnode(GetNP()); TABLE<int,PointIndex::BASE> surfacesonnode(GetNP());
@ -1864,8 +1865,7 @@ namespace netgen
void Mesh :: FindOpenElements (int dom) void Mesh :: FindOpenElements (int dom)
{ {
static int timer = NgProfiler::CreateTimer ("Mesh::FindOpenElements"); static Timer t("Mesh::FindOpenElements"); RegionTimer reg (t);
NgProfiler::RegionTimer reg (timer);
int np = GetNP(); int np = GetNP();
int ne = GetNE(); int ne = GetNE();
@ -2717,6 +2717,8 @@ namespace netgen
void Mesh :: CalcLocalH (double grading) void Mesh :: CalcLocalH (double grading)
{ {
static Timer t("Mesh::CalcLocalH"); RegionTimer reg(t);
if (!lochfunc) if (!lochfunc)
{ {
Point3d pmin, pmax; Point3d pmin, pmax;
@ -3242,6 +3244,8 @@ namespace netgen
void Mesh :: Compress () void Mesh :: Compress ()
{ {
static Timer t("Mesh::Compress"); RegionTimer reg(t);
Array<PointIndex,PointIndex::BASE,PointIndex> op2np(GetNP()); Array<PointIndex,PointIndex::BASE,PointIndex> op2np(GetNP());
Array<MeshPoint> hpoints; Array<MeshPoint> hpoints;
BitArrayChar<PointIndex::BASE> pused(GetNP()); BitArrayChar<PointIndex::BASE> pused(GetNP());
@ -3381,7 +3385,7 @@ namespace netgen
for (int i = 0; i < lockedpoints.Size(); i++) for (int i = 0; i < lockedpoints.Size(); i++)
lockedpoints[i] = op2np[lockedpoints[i]]; lockedpoints[i] = op2np[lockedpoints[i]];
/*
for (int i = 0; i < facedecoding.Size(); i++) for (int i = 0; i < facedecoding.Size(); i++)
facedecoding[i].firstelement = -1; facedecoding[i].firstelement = -1;
for (int i = surfelements.Size()-1; i >= 0; i--) for (int i = surfelements.Size()-1; i >= 0; i--)
@ -3390,7 +3394,8 @@ namespace netgen
surfelements[i].next = facedecoding[ind-1].firstelement; surfelements[i].next = facedecoding[ind-1].firstelement;
facedecoding[ind-1].firstelement = i; facedecoding[ind-1].firstelement = i;
} }
*/
RebuildSurfaceElementLists ();
CalcSurfacesOfNode(); CalcSurfacesOfNode();
@ -3508,6 +3513,8 @@ namespace netgen
int Mesh :: CheckOverlappingBoundary () int Mesh :: CheckOverlappingBoundary ()
{ {
static Timer t("Mesh::CheckOverlappingBoundary"); RegionTimer reg(t);
int i, j, k; int i, j, k;
Point3d pmin, pmax; Point3d pmin, pmax;
@ -4824,9 +4831,10 @@ namespace netgen
bool build_searchtree, bool build_searchtree,
const bool allowindex) const const bool allowindex) const
{ {
const double pointtol = 1e-12; // const double pointtol = 1e-12;
netgen::Point<3> pmin = p - Vec<3> (pointtol, pointtol, pointtol); // netgen::Point<3> pmin = p - Vec<3> (pointtol, pointtol, pointtol);
netgen::Point<3> pmax = p + Vec<3> (pointtol, pointtol, pointtol); // netgen::Point<3> pmax = p + Vec<3> (pointtol, pointtol, pointtol);
if (dimension == 2) if (dimension == 2)
{ {
int ne; int ne;
@ -4839,8 +4847,11 @@ namespace netgen
if (elementsearchtree || build_searchtree) if (elementsearchtree || build_searchtree)
{ {
// update if necessary: // update if necessary:
const_cast<Mesh&>(*this).BuildElementSearchTree (); const_cast<Mesh&>(*this).BuildElementSearchTree ();
elementsearchtree->GetIntersecting (pmin, pmax, locels); // double tol = elementsearchtree->Tolerance();
// netgen::Point<3> pmin = p - Vec<3> (tol, tol, tol);
// netgen::Point<3> pmax = p + Vec<3> (tol, tol, tol);
elementsearchtree->GetIntersecting (p, p, locels);
ne = locels.Size(); ne = locels.Size();
} }
else else
@ -4882,8 +4893,11 @@ namespace netgen
if (elementsearchtree || build_searchtree) if (elementsearchtree || build_searchtree)
{ {
// update if necessary: // update if necessary:
const_cast<Mesh&>(*this).BuildElementSearchTree (); const_cast<Mesh&>(*this).BuildElementSearchTree ();
elementsearchtree->GetIntersecting (pmin, pmax, locels); // double tol = elementsearchtree->Tolerance();
// netgen::Point<3> pmin = p - Vec<3> (tol, tol, tol);
// netgen::Point<3> pmax = p + Vec<3> (tol, tol, tol);
elementsearchtree->GetIntersecting (p, p, locels);
ne = locels.Size(); ne = locels.Size();
} }
else else