From 99f025f002a831a610ed452fabcd18cff5187807 Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Wed, 27 Aug 2008 20:27:32 +0000
Subject: [PATCH] add StopAtDistMax for Treshold field tmp fix build
---
Mesh/Field.cpp | 10 ++++++++--
Mesh/HighOrder.cpp | 28 +++++++++++++++-------------
2 files changed, 23 insertions(+), 15 deletions(-)
diff --git a/Mesh/Field.cpp b/Mesh/Field.cpp
index cbf92dde45..f149924c64 100644
--- a/Mesh/Field.cpp
+++ b/Mesh/Field.cpp
@@ -451,7 +451,7 @@ class ThresholdField : public Field
{
int iField;
double dmin, dmax, lcmin, lcmax;
- bool sigmoid;
+ bool sigmoid, stopAtDistMax;
public:
const char *get_name()
{
@@ -465,6 +465,7 @@ class ThresholdField : public Field
lcmin = 0.1;
lcmax = 1;
sigmoid = false;
+ stopAtDistMax = false;
options["IField"] = new FieldOptionInt(iField, "Field index");
options["DistMin"] = new FieldOptionDouble(dmin, "Distance from entity up to which "
"element size will be LcMin");
@@ -475,6 +476,8 @@ class ThresholdField : public Field
options["Sigmoid"] = new FieldOptionBool(sigmoid, "True to interpolate between LcMin "
"and LcMax using a sigmoid, false to "
"interpolate linearly");
+ options["StopAtDistMax"] = new FieldOptionBool(stopAtDistMax, "True to not impose element "
+ "size outside DistMax");
}
double operator() (double x, double y, double z)
{
@@ -482,7 +485,10 @@ class ThresholdField : public Field
double r = ((*field) (x, y, z) - dmin) / (dmax - dmin);
r = std::max(std::min(r, 1.), 0.);
double lc;
- if(sigmoid){
+ if(stopAtDistMax && r >= 1.){
+ lc = MAX_LC;
+ }
+ else if(sigmoid){
double s = exp(12. * r - 6.) / (1. + exp(12. * r - 6.));
lc = lcmin * (1. - s) + lcmax * s;
}
diff --git a/Mesh/HighOrder.cpp b/Mesh/HighOrder.cpp
index ecb99aa999..adc8fea423 100644
--- a/Mesh/HighOrder.cpp
+++ b/Mesh/HighOrder.cpp
@@ -861,21 +861,23 @@ void getRegionVertices(GRegion *gr,
const double t2 = points(k, 1);
const double t3 = points(k, 2);
- SPoint3 pos;
- incomplete->pnt(t1,t2,t3,pos);
+ // FIXME: KOEN - I had to comment this out (MElement does not have
+ // pnt() member) -- CG
- v = new MVertex(pos.x(),pos.y(),pos.z(),gr);
+ // SPoint3 pos;
+ // incomplete->pnt(t1,t2,t3,pos);
+ // v = new MVertex(pos.x(),pos.y(),pos.z(),gr);
-// double X(0),Y(0),Z(0);
-// for (int j=0; j<incomplete->getNumVertices(); j++){
-// double sf ; incomplete->getShapeFunction(j,t1,t2,t3,sf);
-// MVertex *vt = incomplete->getVertex(j);
-// X += sf * vt->x();
-// Y += sf * vt->y();
-// Z += sf * vt->z();
-// }
-//
- // v = new MVertex(X,Y,Z, gr);
+ double X(0),Y(0),Z(0);
+ for (int j=0; j<incomplete->getNumVertices(); j++){
+ double sf ; incomplete->getShapeFunction(j,t1,t2,t3,sf);
+ MVertex *vt = incomplete->getVertex(j);
+ X += sf * vt->x();
+ Y += sf * vt->y();
+ Z += sf * vt->z();
+ }
+ v = new MVertex(X,Y,Z, gr);
+
gr->mesh_vertices.push_back(v);
vr.push_back(v);
}
--
GitLab