Skip to content
Snippets Groups Projects
Commit ad08dd6a authored by Jean-François Remacle's avatar Jean-François Remacle
Browse files

removed bugs in BL

parent 6e28aae0
No related branches found
No related tags found
No related merge requests found
......@@ -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);
}
......
......@@ -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) {
......
......@@ -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;
......@@ -337,7 +338,7 @@ BoundaryLayerColumns* buidAdditionalPoints2D (GFace *gf)
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);
}
}
......
......@@ -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;
......
......@@ -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);
......
......@@ -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");
}
}
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment