diff --git a/Geo/GEdgeLoop.cpp b/Geo/GEdgeLoop.cpp
index 0d77a5ae9b3ccf6135d0242590cf19fe9da4b050..5803463c09d6d63116d87a8dbf1dcdc68f9265a5 100644
--- a/Geo/GEdgeLoop.cpp
+++ b/Geo/GEdgeLoop.cpp
@@ -101,25 +101,43 @@ int GEdgeLoop::count(GEdge* ge) const
GEdgeLoop::GEdgeLoop(const std::list<GEdge*> &cwire)
{
- std::list<GEdge*> wire(cwire);
- GEdgeSigned *prevOne = 0;
// Sometimes OCC puts a nasty degenerated edge in the middle of the wire ...
// pushing it to front fixes the problem as it concerns gmsh
- for (std::list<GEdge*>::iterator it = wire.begin(); it != wire.end() ; ++it){ // TEST
+
+ std::list<GEdge*> wire;
+ std::vector<GEdge*> degenerated;
+ GEdge *degeneratedToInsert = 0;
+ for (std::list<GEdge*>::const_iterator it = cwire.begin(); it != cwire.end() ; ++it){ // TEST
GEdge *ed = *it;
- if (ed->degenerate(0)){
- wire.erase(it);
- wire.push_front(ed);
- break;
- }
+ if (ed->degenerate(0))degenerated.push_back(ed);
+ else wire.push_back(ed);
}
-
+ if (degenerated.size() == 1){
+ wire.push_front(degenerated[0]);
+ }
+ else if (degenerated.size() == 2){
+ degeneratedToInsert = degenerated[1];
+ wire.push_front(degenerated[0]);
+ }
+ else if (degenerated.size() > 2){
+ Msg::Error("More than two degenerated edges in one model face");
+ }
+
+ GEdgeSigned *prevOne = 0;
GEdgeSigned ges(0,0);
+
+ // printf("go !!!\n");
while(wire.size()){
- ges = nextOne(prevOne, wire);
+ if (prevOne && degeneratedToInsert &&
+ degeneratedToInsert->getBeginVertex () == prevOne->getEndVertex()){
+ ges = GEdgeSigned(1,degeneratedToInsert);
+ degeneratedToInsert = 0;
+ // printf("second degenerated edge inserted\n");
+ }
+ else ges = nextOne(prevOne, wire);
if(ges.getSign() == 0){ // oops
Msg::Error("Something wrong in edge loop of size=%d, no sign !", wire.size());
for (std::list<GEdge* >::iterator it = wire.begin(); it != wire.end(); it++){
diff --git a/Mesh/meshGFace.cpp b/Mesh/meshGFace.cpp
index a3d8790ba1283fbb246511a5f8c9d8286b39d4d8..5aa10ac797751627644eaf21108ef4b9efb16c71 100644
--- a/Mesh/meshGFace.cpp
+++ b/Mesh/meshGFace.cpp
@@ -618,8 +618,9 @@ void modifyInitialMeshForTakingIntoAccountBoundaryLayers (GFace *gf){
deMeshGFace kil_;
kil_(gf);
-
+ printf("remeshing after hoplaboum\n");
meshGenerator(gf, 0, 0, true , true, &hop);
+ printf("remeshing after hoplaboum\n");
gf->quadrangles = blQuads;
gf->triangles.insert(gf->triangles.begin(),blTris.begin(),blTris.end());
@@ -1060,32 +1061,31 @@ bool meshGenerator(GFace *gf, int RECUR_ITER,
gf->model()->add(ge);
}
- gf->triangles.clear();
- gf->quadrangles.clear();
-
{
int nb_swap;
Msg::Debug("Delaunizing the initial mesh");
delaunayizeBDS(gf, *m, nb_swap);
}
+ gf->triangles.clear();
+ gf->quadrangles.clear();
Msg::Debug("Starting to add internal points");
// start mesh generation
if(!algoDelaunay2D(gf) && !onlyInitialMesh){
- if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine) {
- printf("coucou here !!!\n");
- backgroundMesh::unset();
- buildBackGroundMesh (gf);
- }
+ // if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine || 1) {
+ // printf("coucou here !!!\n");
+ // backgroundMesh::unset();
+ // buildBackGroundMesh (gf);
+ // }
refineMeshBDS(gf, *m, CTX::instance()->mesh.refineSteps, true,
&recoverMapInv);
optimizeMeshBDS(gf, *m, 2);
refineMeshBDS(gf, *m, CTX::instance()->mesh.refineSteps, false,
&recoverMapInv);
optimizeMeshBDS(gf, *m, 2);
- if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine) {
- backgroundMesh::unset();
- }
+ // if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine || 1) {
+ // backgroundMesh::unset();
+ // }
}
/*
computeMeshSizeFieldAccuracy(gf, *m, gf->meshStatistics.efficiency_index,
@@ -1636,11 +1636,11 @@ static bool meshGeneratorPeriodic(GFace *gf, bool debug = true)
// start mesh generation for periodic face
if(!algoDelaunay2D(gf)){
// need for a BGM for cross field
- if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine) {
- // printf("coucou here !!!\n");
- // backgroundMesh::unset();
- // buildBackGroundMesh (gf);
- }
+ // if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine || 1) {
+ // printf("coucou here !!!\n");
+ // backgroundMesh::unset();
+ // buildBackGroundMesh (gf);
+ // }
refineMeshBDS(gf, *m, CTX::instance()->mesh.refineSteps, true);
optimizeMeshBDS(gf, *m, 2);
refineMeshBDS(gf, *m, -CTX::instance()->mesh.refineSteps, false);
@@ -1653,9 +1653,9 @@ static bool meshGeneratorPeriodic(GFace *gf, bool debug = true)
gf->meshStatistics.nbEdge,
gf->meshStatistics.nbGoodLength);*/
gf->meshStatistics.status = GFace::DONE;
- if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine) {
- // backgroundMesh::unset();
- }
+ // if(CTX::instance()->mesh.recombineAll || gf->meshAttributes.recombine || 1) {
+ // backgroundMesh::unset();
+ // }
}
// fill the small gmsh structures
diff --git a/Mesh/meshGFaceBDS.cpp b/Mesh/meshGFaceBDS.cpp
index 1173d253194928e0ed5221b902c79f3d7526d723..143880b1bf6425ac367ba8751ce32db70c5802cd 100644
--- a/Mesh/meshGFaceBDS.cpp
+++ b/Mesh/meshGFaceBDS.cpp
@@ -28,8 +28,6 @@ double computeEdgeLinearLength(BDS_Point *p1, BDS_Point *p2)
const double dx = p1->X - p2->X;
const double dy = p1->Y - p2->Y;
const double dz = p1->Z - p2->Z;
- // FIXME SUPER HACK
- //const double l = std::max(fabs(dx),std::max(fabs(dy),fabs(dz)));
const double l = sqrt(dx * dx + dy * dy + dz * dz);
return l;
}
@@ -41,14 +39,14 @@ inline double computeEdgeLinearLength(BDS_Point *p1, BDS_Point *p2, GFace *f,
double SCALINGU, double SCALINGV)
{
// FIXME SUPER HACK
- /*
- if (CTX::instance()->mesh.recombineAll || f->meshAttributes.recombine){
- double quadAngle = backgroundMesh::current()->getAngle (0.5 * (p1->u + p2->u) * SCALINGU, 0.5 * (p1->v + p2->v) * SCALINGV,0);
- const double a [2] = {p1->u,p1->v};
- const double b [2] = {p2->u,p2->v};
- return lengthInfniteNorm(a,b, quadAngle);
- }
- */
+
+ // if (CTX::instance()->mesh.recombineAll || f->meshAttributes.recombine || 1){
+ // double quadAngle = backgroundMesh::current()->getAngle (0.5 * (p1->u + p2->u) * SCALINGU, 0.5 * (p1->v + p2->v) * SCALINGV,0);
+ // const double a [2] = {p1->u,p1->v};
+ // const double b [2] = {p2->u,p2->v};
+ // return lengthInfniteNorm(a,b, quadAngle);
+ // }
+
GPoint GP = f->point(SPoint2(0.5 * (p1->u + p2->u) * SCALINGU,
0.5 * (p1->v + p2->v) * SCALINGV));
diff --git a/Mesh/meshGFaceDelaunayInsertion.cpp b/Mesh/meshGFaceDelaunayInsertion.cpp
index ffda6e0a222dc01432de3f4bd39e90dabcc13bd8..2cdaab0df3535bc0063c0e7674cd76cdebe352b9 100644
--- a/Mesh/meshGFaceDelaunayInsertion.cpp
+++ b/Mesh/meshGFaceDelaunayInsertion.cpp
@@ -1410,6 +1410,9 @@ void optimalPointFrontalQuadAniso (GFace *gf,
void buildBackGroundMesh (GFace *gf) {
+ // printf("build bak mesh\n");
+ quadsToTriangles(gf, 100000);
+
if (!backgroundMesh::current()) {
std::vector<MTriangle*> TR;
for(int i=0;i<gf->triangles.size();i++){
@@ -1427,7 +1430,9 @@ void buildBackGroundMesh (GFace *gf) {
// Re-enable curv control if asked
CTX::instance()->mesh.lcFromCurvature = CurvControl;
// apply this to the BGM
+ // printf("1 end build bak mesh\n");
backgroundMesh::set(gf);
+ // printf("2 end build bak mesh\n");
char name[256];
if (CTX::instance()->mesh.saveAll){
sprintf(name,"bgm-%d.pos",gf->tag());
@@ -1435,7 +1440,10 @@ void buildBackGroundMesh (GFace *gf) {
sprintf(name,"cross-%d.pos",gf->tag());
backgroundMesh::current()->print(name,gf,1);
}
+ gf->triangles = TR;
}
+ // printf("end build bak mesh\n");
+
}
void bowyerWatsonFrontalLayers(GFace *gf, bool quad)
diff --git a/Mesh/meshGFaceOptimize.cpp b/Mesh/meshGFaceOptimize.cpp
index 28d0ade5677006ecc7d883a962500743619573c4..22b3c2d86ef37d7f76bc2b4055af65e92405fe45 100644
--- a/Mesh/meshGFaceOptimize.cpp
+++ b/Mesh/meshGFaceOptimize.cpp
@@ -1689,7 +1689,7 @@ static int _recombineIntoQuads(GFace *gf, int recur_level, bool cubicGraph = 1)
if (pairs[i].n2->onWhat()->dim() < 2)NB++;
if (pairs[i].n3->onWhat()->dim() < 2)NB++;
if (pairs[i].n4->onWhat()->dim() < 2)NB++;
- if (elen[i] > 60 && NB > 2) {elen[i] = 1000;}
+ if (elen[i] > gf->meshAttributes.recombineAngle && NB > 2) {elen[i] = 1000;}
}
if (cubicGraph){
@@ -1827,6 +1827,7 @@ static int _recombineIntoQuads(GFace *gf, int recur_level, bool cubicGraph = 1)
while(itp != pairs.end()){
// recombine if difference between max quad angle and right
// angle is smaller than tol
+ // printf("%g %g\n",gf->meshAttributes.recombineAngle,itp->angle);
if(itp->angle < gf->meshAttributes.recombineAngle){
MElement *t1 = itp->t1;
MElement *t2 = itp->t2;
@@ -1932,6 +1933,7 @@ void recombineIntoQuads(GFace *gf,
void quadsToTriangles(GFace *gf, double minqual)
{
+ printf("COUCOU %g\n",minqual);
std::vector<MQuadrangle*> qds;
for (int i=0;i<gf->quadrangles.size();i++){
MQuadrangle *q = gf->quadrangles[i];
@@ -1972,7 +1974,7 @@ void quadsToTriangles(GFace *gf, double minqual)
void recombineIntoQuadsIterative(GFace *gf)
{
recombineIntoQuads(gf);
- quadsToTriangles(gf,0.03);
+ // quadsToTriangles(gf, 1. - gf->meshAttributes.recombineAngle/90. );
return;
int COUNT = 0;
while (1){
diff --git a/benchmarks/2d/Square-01.geo b/benchmarks/2d/Square-01.geo
index ba31ef8cd470803eb26601bb78c1069022749dae..0f127051e6acc9f3f0508c5ee15eaf6f25869b2b 100644
--- a/benchmarks/2d/Square-01.geo
+++ b/benchmarks/2d/Square-01.geo
@@ -2,7 +2,7 @@ fact = 100;
lc = .1 * fact;
Point(1) = {0.0,0.0,0,lc};
Point(2) = {1* fact,0.0,0,lc};
-Point(3) = {1* fact,1* fact,0,lc*.5};
+Point(3) = {1* fact,1* fact,0,lc};
Point(4) = {0,1* fact,0,lc};
Line(1) = {3,2};
Line(2) = {2,1};
diff --git a/benchmarks/2d/wing-splines.geo b/benchmarks/2d/wing-splines.geo
index 839e8b4f9995d6e19e915464f7f59cd4bab4767d..4b9329a8641c9fefb639e251168df71e13d0ea9c 100644
--- a/benchmarks/2d/wing-splines.geo
+++ b/benchmarks/2d/wing-splines.geo
@@ -605,7 +605,7 @@ Circle(408) = {4355,4351,4354};
Circle(409) = {4354,4351,4353};
Circle(410) = {4353,4351,4352};
*/
-Recombine Surface {406};
+Recombine Surface {406} = 90;
Mesh.RecombinationAlgorithm=1;
Field[2] = BoundaryLayer;