From 2abe97d55c5731aed80b74d0c9b5a7a0f320052d Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Thu, 22 Jan 2009 09:01:16 +0000
Subject: [PATCH] * added isInsideTolerance in MElement, too * point not found
should NOT be an error in the octree
---
Common/OctreeInternals.cpp | 3 ++-
Geo/MElement.cpp | 1 +
Geo/MElement.h | 32 +++++++++++++++++++++++---------
Post/OctreePost.cpp | 9 ++++++---
4 files changed, 32 insertions(+), 13 deletions(-)
diff --git a/Common/OctreeInternals.cpp b/Common/OctreeInternals.cpp
index d76f27b86a..4e45e79ff5 100644
--- a/Common/OctreeInternals.cpp
+++ b/Common/OctreeInternals.cpp
@@ -295,7 +295,8 @@ void *searchElement(octantBucket *_buckets_head, double *_pt, globalInfo *_globa
ptrBucket = findElementBucket(_buckets_head, _pt);
if (ptrBucket == NULL) {
- Msg::Error("The point is not in the domain");
+ // this is not an error
+ Msg::Debug("The point is not in the domain");
return NULL;
}
diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 1fa18ff896..57243bdafb 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -28,6 +28,7 @@
extern Context_T CTX;
int MElement::_globalNum = 0;
+double MElement::_isInsideTolerance = 1.e-6;
double MElementLessThanLexicographic::tolerance = 1.e-6;
void MElement::_getEdgeRep(MVertex *v0, MVertex *v1,
diff --git a/Geo/MElement.h b/Geo/MElement.h
index ab3475377b..3bfd5002f5 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -33,6 +33,9 @@ class MElement
// a visibility flag
char _visible;
protected:
+ // the tolerance used to determine if a point is inside an element,
+ // in parametric coordinates
+ static double _isInsideTolerance;
void _getEdgeRep(MVertex *v0, MVertex *v1,
double *x, double *y, double *z, SVector3 *n,
int faceIndex=-1);
@@ -55,6 +58,10 @@ class MElement
// reset the global node number
static void resetGlobalNumber(){ _globalNum = 0; }
+ // set/get the tolerance for isInside() test
+ static void setTolerance (const double tol){ _isInsideTolerance = tol; }
+ static double getTolerance () { return _isInsideTolerance; }
+
// return the tag of the element
virtual int getNum(){ return _num; }
@@ -187,7 +194,7 @@ class MElement
// test if a point, given in parametric coordinates, belongs to the
// element
- virtual bool isInside(double u, double v, double w, double tol=1.e-8) = 0;
+ virtual bool isInside(double u, double v, double w) = 0;
// interpolate the given nodal data (resp. its gradient, curl and
// divergence) at point (u,v,w) in parametric coordinates
@@ -296,7 +303,7 @@ class MPoint : public MElement {
{
s[0][0] = s[0][1] = s[0][2] = 0.;
}
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
return true;
}
@@ -357,8 +364,9 @@ class MLine : public MElement {
MVertex *tmp = _v[0]; _v[0] = _v[1]; _v[1] = tmp;
}
virtual const gmshFunctionSpace* getFunctionSpace(int o=-1) const;
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(u < -(1. + tol) || u > (1. + tol))
return false;
return true;
@@ -569,8 +577,9 @@ class MTriangle : public MElement {
MVertex *tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
}
virtual const gmshFunctionSpace* getFunctionSpace(int o=-1) const;
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(u < (-tol) || v < (-tol) || u > ((1. + tol) - v))
return false;
return true;
@@ -862,8 +871,9 @@ class MQuadrangle : public MElement {
s[2][0] = 0.25 * (1. + v); s[2][1] = 0.25 * (1. + u); s[2][2] = 0.;
s[3][0] = -0.25 * (1. + v); s[3][1] = 0.25 * (1. - u); s[3][2] = 0.;
}
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(u < -(1. + tol) || v < -(1. + tol) || u > (1. + tol) || v > (1. + tol))
return false;
return true;
@@ -1190,8 +1200,9 @@ class MTetrahedron : public MElement {
virtual double etaShapeMeasure();
void xyz2uvw(double xyz[3], double uvw[3]);
virtual const gmshFunctionSpace* getFunctionSpace(int o=-1) const;
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(u < (-tol) || v < (-tol) || w < (-tol) || u > ((1. + tol) - v - w))
return false;
return true;
@@ -1657,8 +1668,9 @@ class MHexahedron : public MElement {
s[7][1] = 0.125 * (1. - u) * (1. + w);
s[7][2] = 0.125 * (1. - u) * (1. + v);
}
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(u < -(1. + tol) || v < -(1. + tol) || w < -(1. + tol) ||
u > (1. + tol) || v > (1. + tol) || w > (1. + tol))
return false;
@@ -2131,8 +2143,9 @@ class MPrism : public MElement {
s[5][1] = 0.5 * (1. + w) ;
s[5][2] = 0.5 * v ;
}
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(w > (1. + tol) || w < -(1. + tol) || u < (1. + tol)
|| v < (1. + tol) || u > ((1. + tol) - v))
return false;
@@ -2577,8 +2590,9 @@ class MPyramid : public MElement {
s[4][1] = 0.;
s[4][2] = 1.;
}
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ virtual bool isInside(double u, double v, double w)
{
+ double tol = _isInsideTolerance;
if(u < (w - (1. + tol)) || u > ((1. + tol) - w) || v < (w - (1. + tol)) ||
v > ((1. + tol) - w) || w < (-tol) || w > (1. + tol))
return false;
diff --git a/Post/OctreePost.cpp b/Post/OctreePost.cpp
index 087ca75426..a7731aa4af 100644
--- a/Post/OctreePost.cpp
+++ b/Post/OctreePost.cpp
@@ -417,11 +417,14 @@ bool OctreePost::searchScalar(double x, double y, double z, double *values,
int step, double *size)
{
bool a = _searchScalar(x, y, z, values, step, size);
- if (!a){
- double oldeps = element::getTolerance();
+ if(!a){
+ double oldeps1 = element::getTolerance();
+ double oldeps2 = MElement::getTolerance();
element::setTolerance(10.);
+ MElement::setTolerance(10.);
a = _searchScalar(x, y, z, values, step, size);
- element::setTolerance(oldeps);
+ element::setTolerance(oldeps1);
+ MElement::setTolerance(oldeps2);
}
if (!a) Msg::Debug("No element found containing point (%g,%g,%g)", x, y, z);
return a;
--
GitLab