From ad08dd6a579b094410329fa758d6de6738709e0a Mon Sep 17 00:00:00 2001
From: Jean-Francois Remacle <jean-francois.remacle@uclouvain.be>
Date: Wed, 18 Apr 2012 07:02:27 +0000
Subject: [PATCH] removed bugs in BL
---
Fltk/highOrderToolsWindow.cpp | 1 +
Mesh/Field.cpp | 10 ++++--
Mesh/meshGFaceBoundaryLayers.cpp | 6 ++--
contrib/HighOrderMeshOptimizer/OptHOM.cpp | 4 +--
contrib/HighOrderMeshOptimizer/OptHOM.h | 2 +-
contrib/HighOrderMeshOptimizer/OptHomRun.cpp | 36 +++++++++++---------
6 files changed, 35 insertions(+), 24 deletions(-)
diff --git a/Fltk/highOrderToolsWindow.cpp b/Fltk/highOrderToolsWindow.cpp
index cf40d0e66b..de8d0eac70 100644
--- a/Fltk/highOrderToolsWindow.cpp
+++ b/Fltk/highOrderToolsWindow.cpp
@@ -79,6 +79,7 @@ static void highordertools_runelas_cb(Fl_Widget *w, void *data)
p.BARRIER = threshold;
p.onlyVisible = onlyVisible;
p.dim = GModel::current()->getNumRegions() ? 3 : 2;
+ p.itMax = (int) o->value[3]->value();
HighOrderMeshOptimizer (GModel::current(),p);
}
diff --git a/Mesh/Field.cpp b/Mesh/Field.cpp
index ae3ead8182..addd52d47e 100644
--- a/Mesh/Field.cpp
+++ b/Mesh/Field.cpp
@@ -1898,7 +1898,13 @@ void BoundaryLayerField::operator() (AttractorField *cc, double dist,
std::pair<AttractorInfo,SPoint3> pp = cc->getAttractorInfo();
if (pp.first.dim ==0){
GVertex *v = GModel::current()->getVertexByTag(pp.first.ent);
- SVector3 t1 = SVector3(v->x() -x,v->y() -y,v->z() -z);
+ SVector3 t1;
+ if (dist < thickness){
+ t1 = SVector3(1,0,0);
+ }
+ else {
+ t1 = SVector3(v->x() -x,v->y() -y,v->z() -z);
+ }
metr = buildMetricTangentToCurve(t1,lc_n,lc_n);
return;
}
@@ -1955,7 +1961,7 @@ void BoundaryLayerField::operator() (double x, double y, double z,
}
for(std::list<int>::iterator it = edges_id.begin();
it != edges_id.end(); ++it) {
- _att_fields.push_back(new AttractorField(1,*it,300000));
+ _att_fields.push_back(new AttractorField(1,*it,3000));
}
for(std::list<int>::iterator it = faces_id.begin();
it != faces_id.end(); ++it) {
diff --git a/Mesh/meshGFaceBoundaryLayers.cpp b/Mesh/meshGFaceBoundaryLayers.cpp
index f74e209a3d..cd830de491 100644
--- a/Mesh/meshGFaceBoundaryLayers.cpp
+++ b/Mesh/meshGFaceBoundaryLayers.cpp
@@ -276,6 +276,7 @@ BoundaryLayerColumns* buidAdditionalPoints2D (GFace *gf)
}
_columns->addFan (*it,ee[0],ee[1],true);
+ // printf("fansize = %d\n",fanSize);
for (int i=-1; i<=fanSize; i++){
double t = (double)(i+1)/ (fanSize+1);
double alpha = t * AMAX + (1.-t)* AMIN;
@@ -330,14 +331,14 @@ BoundaryLayerColumns* buidAdditionalPoints2D (GFace *gf)
// if (_dirs.size() > 1) printf("l = %g metric = %g %g %g dim %d tag %d \n",l,metric[0],metric[1],metric[2],current->onWhat()->dim(),current->onWhat()->tag());
// printf("%g %g\n",l,LL);
if (l >= blf->hfar){
- // printf("stopping %g %g\n",l,LL);
+ // printf("stopping %g %g\n",l,LL);
break;
}
// printf("%g %g %g \n",current->x(),current->y(),blf->current_distance);
if (blf->current_closest != catt || blf -> current_distance < _current_distance){
SVector3 aaa (_close- blf->_closest_point);
if (aaa.norm() > 8*blf->hwall_n || blf -> current_distance < _current_distance){
- printf("reaching the skelton %d\n", (int) _column.size());
+ // printf("reaching the skelton %d\n", (int) _column.size());
delete _column[_column.size()-1];
_column.erase(--_column.end());
_metrics.erase(--_metrics.end());
@@ -366,6 +367,7 @@ BoundaryLayerColumns* buidAdditionalPoints2D (GFace *gf)
if (_column.size() > nbCol)break; // FIXME
p = pnew;
}
+ // if (_dirs.size() > 1)printf("adding column with %d nodes\n",_column.size());
_columns->addColumn(n,*it, _column, _metrics);
}
}
diff --git a/contrib/HighOrderMeshOptimizer/OptHOM.cpp b/contrib/HighOrderMeshOptimizer/OptHOM.cpp
index ce38c6b41d..4c7dd27eb2 100644
--- a/contrib/HighOrderMeshOptimizer/OptHOM.cpp
+++ b/contrib/HighOrderMeshOptimizer/OptHOM.cpp
@@ -200,7 +200,7 @@ void OptHOM::OptimPass(alglib::real_1d_array &x, const alglib::real_1d_array &in
-int OptHOM::optimize(double weightFixed, double weightFree, double barrier_, int pInt, int itMax)
+int OptHOM::optimize(double weightFixed, double weightFree, double barrier_, int pInt, int itMax, double &minJ, double &maxJ)
{
barrier = barrier_;
progressInterv = pInt;
@@ -221,7 +221,7 @@ int OptHOM::optimize(double weightFixed, double weightFree, double barrier_, int
mesh.getUvw(x.getcontent());
// Calculate initial performance
- double minJ, maxJ;
+ // double minJ, maxJ;
getDistances(initMaxD, initAvgD, minJ, maxJ);
const double jacBarStart = (minJ > 0.) ? 0.9*minJ : 1.1*minJ;
diff --git a/contrib/HighOrderMeshOptimizer/OptHOM.h b/contrib/HighOrderMeshOptimizer/OptHOM.h
index 575a5eb3d3..23c7efaab9 100644
--- a/contrib/HighOrderMeshOptimizer/OptHOM.h
+++ b/contrib/HighOrderMeshOptimizer/OptHOM.h
@@ -23,7 +23,7 @@ public:
OptHOM(GEntity *gf, std::set<MVertex*> & toFix, int method);
void getDistances(double &distMaxBND, double &distAvgBND, double &minJac, double &maxJac);
- int optimize(double lambda, double lambda2, double barrier, int pInt, int itMax); // optimize one list of elements
+ int optimize(double lambda, double lambda2, double barrier, int pInt, int itMax, double &minJ, double &maxJ); // optimize one list of elements
void updateMesh(const alglib::real_1d_array &x);
void evalObjGrad(const alglib::real_1d_array &x, double &Obj, alglib::real_1d_array &gradObj);
void printProgress(const alglib::real_1d_array &x, double Obj);
diff --git a/contrib/HighOrderMeshOptimizer/OptHomRun.cpp b/contrib/HighOrderMeshOptimizer/OptHomRun.cpp
index b0d03813eb..ef79978e23 100644
--- a/contrib/HighOrderMeshOptimizer/OptHomRun.cpp
+++ b/contrib/HighOrderMeshOptimizer/OptHomRun.cpp
@@ -55,8 +55,8 @@ std::set<MVertex*> filter2D(GFace *gf, int nbLayers, double _qual, bool onlythew
}
}
}
- printf("FILTER2D : %d bad triangles found among %6d\n", ts.size(), gf->triangles.size());
- printf("FILTER2D : %d bad quads found among %6d\n", qs.size(), gf->quadrangles.size());
+ // printf("FILTER2D : %d bad triangles found among %6d\n", ts.size(), gf->triangles.size());
+ // printf("FILTER2D : %d bad quads found among %6d\n", qs.size(), gf->quadrangles.size());
// add layers of elements around bad elements
for (int layer = 0; layer < nbLayers; layer++) {
@@ -113,7 +113,7 @@ std::set<MVertex*> filter2D(GFace *gf, int nbLayers, double _qual, bool onlythew
gf->triangles.insert(gf->triangles.begin(), ts.begin(), ts.end());
gf->quadrangles.insert(gf->quadrangles.begin(), qs.begin(), qs.end());
- printf("FILTER2D : AFTER FILTERING %d triangles %d quads remain %d nodes on the boundary\n", gf->triangles.size(), gf->quadrangles.size(), not_touched.size());
+ // printf("FILTER2D : AFTER FILTERING %d triangles %d quads remain %d nodes on the boundary\n", gf->triangles.size(), gf->quadrangles.size(), not_touched.size());
return not_touched;
}
@@ -132,7 +132,7 @@ std::set<MVertex*> filter3D(GRegion *gr, int nbLayers, double _qual)
if (ts.size() == 1) break;
}
- printf("FILTER3D : %d bad tets found among %6d\n", ts.size(), gr->tetrahedra.size());
+ // printf("FILTER3D : %d bad tets found among %6d\n", ts.size(), gr->tetrahedra.size());
// add layers of elements around bad elements
for (int layer = 0; layer < nbLayers; layer++) {
@@ -165,7 +165,7 @@ std::set<MVertex*> filter3D(GRegion *gr, int nbLayers, double _qual)
gr->tetrahedra.clear();
gr->tetrahedra.insert(gr->tetrahedra.begin(), ts.begin(), ts.end());
- printf("FILTER3D : AFTER FILTERING %d tets remain %d nodes on the boundary\n", gr->tetrahedra.size(), not_touched.size());
+ // printf("FILTER3D : AFTER FILTERING %d tets remain %d nodes on the boundary\n", gr->tetrahedra.size(), not_touched.size());
return not_touched;
}
@@ -187,6 +187,7 @@ void HighOrderMeshOptimizer (GModel *gm, OptHomParameters &p)
if (p.dim == 2) {
for (GModel::fiter itf = gm->firstFace(); itf != gm->lastFace(); ++itf) {
+ if (p.onlyVisible && !(*itf)->getVisibility())continue;
std::vector<MTriangle*> tris = (*itf)->triangles;
std::vector<MQuadrangle*> quas = (*itf)->quadrangles;
std::set<MVertex*> toFix = filter2D(*itf, p.nbLayers, p.BARRIER);
@@ -196,33 +197,34 @@ void HighOrderMeshOptimizer (GModel *gm, OptHomParameters &p)
double distMaxBND, distAvgBND, minJac, maxJac;
temp->getDistances(distMaxBND, distAvgBND, minJac, maxJac);
- std::ostringstream ossI;
- ossI << "initial_" << (*itf)->tag() << ".msh";
- temp->mesh.writeMSH(ossI.str().c_str());
- printf("--> Mesh Loaded for GFace %d : minJ %12.5E -- maxJ %12.5E\n", (*itf)->tag(), minJac, maxJac);
+ // std::ostringstream ossI;
+ // ossI << "initial_" << (*itf)->tag() << ".msh";
+ // temp->mesh.writeMSH(ossI.str().c_str());
+ // printf("--> Mesh Loaded for GFace %d : minJ %12.5E -- maxJ %12.5E\n", (*itf)->tag(), minJac, maxJac);
if (distMaxBND < 1.e-8 * SIZE && minJac > p.BARRIER) continue;
- int result = temp->optimize(p.weightFixed, p.weightFree, p.BARRIER, samples, p.itMax);
+ p.SUCCESS = temp->optimize(p.weightFixed, p.weightFree, p.BARRIER, samples, p.itMax, p.minJac, p.maxJac);
temp->mesh.updateGEntityPositions();
- std::ostringstream ossF;
- ossF << "final_" << (*itf)->tag() << ".msh";
- temp->mesh.writeMSH(ossF.str().c_str());
+ // std::ostringstream ossF;
+ // ossF << "final_" << (*itf)->tag() << ".msh";
+ // temp->mesh.writeMSH(ossF.str().c_str());
}
}
else if (p.dim == 3) {
for (GModel::riter itr = gm->firstRegion(); itr != gm->lastRegion(); ++itr) {
+ if (p.onlyVisible && !(*itr)->getVisibility())continue;
std::vector<MTetrahedron*> tets = (*itr)->tetrahedra;
std::set<MVertex*> toFix;
filter3D(*itr, p.nbLayers, p.BARRIER);
OptHOM *temp = new OptHOM(*itr, toFix, method);
double distMaxBND, distAvgBND, minJac, maxJac;
temp->getDistances(distMaxBND, distAvgBND, minJac, maxJac);
- temp->mesh.writeMSH("initial.msh");
- printf("--> Mesh Loaded for GRegion %d : minJ %12.5E -- maxJ %12.5E\n", (*itr)->tag(), minJac, maxJac);
+ // temp->mesh.writeMSH("initial.msh");
+ // printf("--> Mesh Loaded for GRegion %d : minJ %12.5E -- maxJ %12.5E\n", (*itr)->tag(), minJac, maxJac);
if (distMaxBND < 1.e-8 * SIZE && minJac > p.BARRIER) continue;
- int result = temp->optimize(p.weightFixed, p.weightFree, p.BARRIER, samples, p.itMax);
+ p.SUCCESS = temp->optimize(p.weightFixed, p.weightFree, p.BARRIER, samples, p.itMax, p.minJac, p.maxJac);
temp->mesh.updateGEntityPositions();
(*itr)->tetrahedra = tets;
- temp->mesh.writeMSH("final.msh");
+ // temp->mesh.writeMSH("final.msh");
}
}
}
--
GitLab