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

*** empty log message ***

parent 78e8c6a8
No related branches found
No related tags found
No related merge requests found
...@@ -15,6 +15,14 @@ ...@@ -15,6 +15,14 @@
GEdgeCompound::GEdgeCompound(GModel *m, int tag, std::vector<GEdge*> &compound) GEdgeCompound::GEdgeCompound(GModel *m, int tag, std::vector<GEdge*> &compound)
: GEdge(m, tag, 0 , 0), _compound(compound) : GEdge(m, tag, 0 , 0), _compound(compound)
{ {
for (std::vector<GEdge*>::iterator it = compound.begin(); it != compound.end(); ++it){
if (!(*it)) {
Msg::Error("Incorrect edge in compound line %d\n", tag);
Msg::Exit(1);
}
}
orderEdges (); orderEdges ();
int N = _compound.size(); int N = _compound.size();
v0 = _orientation[0] ? _compound[0]->getBeginVertex() : _compound[0]->getEndVertex(); v0 = _orientation[0] ? _compound[0]->getBeginVertex() : _compound[0]->getEndVertex();
...@@ -75,7 +83,7 @@ void GEdgeCompound::orderEdges() ...@@ -75,7 +83,7 @@ void GEdgeCompound::orderEdges()
} }
else{ else{
Msg::Error("EdgeCompound %d is wrong (it has %d end points)",tag(),tempv.size()); Msg::Error("EdgeCompound %d is wrong (it has %d end points)",tag(),tempv.size());
exit(1); Msg::Exit(1);
} }
//loop over all segments to order segments and store it in the list _c //loop over all segments to order segments and store it in the list _c
......
...@@ -174,14 +174,20 @@ GFaceCompound::GFaceCompound(GModel *m, int tag, std::list<GFace*> &compound, ...@@ -174,14 +174,20 @@ GFaceCompound::GFaceCompound(GModel *m, int tag, std::list<GFace*> &compound,
std::list<GEdge*> &U1, std::list<GEdge*> &V1) std::list<GEdge*> &U1, std::list<GEdge*> &V1)
: GFace(m, tag), _compound(compound), _U0(U0), _U1(U1), _V0(V0), _V1(V1), oct(0) : GFace(m, tag), _compound(compound), _U0(U0), _U1(U1), _V0(V0), _V1(V1), oct(0)
{ {
for (std::list<GFace*>::iterator it = _compound.begin(); it != _compound.end(); ++it){
if (!(*it)){
Msg::Error("Incorrect face in compound surface %d\n", tag);
Msg::Exit(1);
}
}
getBoundingEdges(); getBoundingEdges();
if (!_U0.size()) _type = UNITCIRCLE; if (!_U0.size()) _type = UNITCIRCLE;
else if (!_V1.size()) _type = UNITCIRCLE; else if (!_V1.size()) _type = UNITCIRCLE;
else _type = SQUARE; else _type = SQUARE;
for (std::list<GFace*>::iterator it = _compound.begin(); it != _compound.end(); ++it){
if (!(*it)) Msg::Error("Incorrect face in compound surface %d\n", tag);
}
} }
...@@ -307,16 +313,16 @@ SPoint2 GFaceCompound::getCoordinates(MVertex *v) const ...@@ -307,16 +313,16 @@ SPoint2 GFaceCompound::getCoordinates(MVertex *v) const
vR = ge->mesh_vertices[j]; vR = ge->mesh_vertices[j];
vR->getParameter(0,tR); vR->getParameter(0,tR);
if (!vR->getParameter(0,tR)) { if (!vR->getParameter(0,tR)) {
printf("vertex vr %p not medgevertex \n", vR); Msg::Error("vertex vr %p not medgevertex \n", vR);
exit(1); Msg::Exit(1);
} }
//printf("***tLoc=%g tL=%g, tR=%g L=%p (%g,%g,%g) et R= %p (%g,%g,%g)\n", tLoc, tL, tR, vL, vL->x(),vL->y(),vL->z(),vR, vR->x(), vR->y(), vR->z()); //printf("***tLoc=%g tL=%g, tR=%g L=%p (%g,%g,%g) et R= %p (%g,%g,%g)\n", tLoc, tL, tR, vL, vL->x(),vL->y(),vL->z(),vR, vR->x(), vR->y(), vR->z());
if (tLoc >= tL && tLoc <= tR){ if (tLoc >= tL && tLoc <= tR){
found = true; found = true;
itR = coordinates.find(vR); itR = coordinates.find(vR);
if (itR == coordinates.end()){ if (itR == coordinates.end()){
printf("vertex %p (%g %g %g) not found\n", vR, vR->x(), vR->y(), vR->z()); Msg::Error("vertex %p (%g %g %g) not found\n", vR, vR->x(), vR->y(), vR->z());
exit(1); Msg::Exit(1);
} }
break; break;
} }
......
...@@ -378,7 +378,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER, ...@@ -378,7 +378,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER,
std::list<GEdge*> edges = gf->edges(); std::list<GEdge*> edges = gf->edges();
// here, we will replace edges by their compounds // here, we will replace edges by their compounds
//printf("***** In meshGFace: \n"); printf("***** In meshGFace: \n");
if (gf->geomType() == GEntity::CompoundSurface){ if (gf->geomType() == GEntity::CompoundSurface){
//printf("replace edges by compound lines \n"); //printf("replace edges by compound lines \n");
std::set<GEdge*> mySet; std::set<GEdge*> mySet;
...@@ -448,7 +448,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER, ...@@ -448,7 +448,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER,
double *V_ = new double[all_vertices.size()]; double *V_ = new double[all_vertices.size()];
v_container::iterator itv = all_vertices.begin(); v_container::iterator itv = all_vertices.begin();
//printf("boundary vertices size = %d \n", all_vertices.size()); printf("boundary vertices size = %d \n", all_vertices.size());
int count = 0; int count = 0;
SBoundingBox3d bbox; SBoundingBox3d bbox;
...@@ -548,7 +548,6 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER, ...@@ -548,7 +548,6 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER,
m->add_triangle(V1->getIndex(), V2->getIndex(), V3->getIndex()); m->add_triangle(V1->getIndex(), V2->getIndex(), V3->getIndex());
} }
// Recover the boundary edges and compute characteristic lenghts // Recover the boundary edges and compute characteristic lenghts
// using mesh edge spacing // using mesh edge spacing
if( debug && RECUR_ITER == 0){ if( debug && RECUR_ITER == 0){
...@@ -632,11 +631,12 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER, ...@@ -632,11 +631,12 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER,
return gmsh2DMeshGenerator(gf, RECUR_ITER+1, repairSelfIntersecting1dMesh, debug); return gmsh2DMeshGenerator(gf, RECUR_ITER+1, repairSelfIntersecting1dMesh, debug);
return false; return false;
} }
if(RECUR_ITER > 0) if(RECUR_ITER > 0)
Msg::Warning(":-) Gmsh was able to recover all edges after %d iterations", Msg::Warning(":-) Gmsh was able to recover all edges after %d iterations", RECUR_ITER);
RECUR_ITER);
Msg::Info("Boundary Edges recovered for surface %d", gf->tag());
// Msg::Info("Boundary Edges recovered for surface %d", gf->tag());
// Look for an edge that is on the boundary for which one of the two // Look for an edge that is on the boundary for which one of the two
// neighbors has a negative number node. The other triangle is // neighbors has a negative number node. The other triangle is
// inside the domain and, because all edges were recovered, // inside the domain and, because all edges were recovered,
...@@ -669,7 +669,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER, ...@@ -669,7 +669,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER,
++it; ++it;
} }
// compute characteristic lengths at vertices // compute characteristic lengths at vertices
Msg::Debug("Computing mesh size field at mesh vertices", edgesToRecover.size()); Msg::Debug("Computing mesh size field at mesh vertices %d", edgesToRecover.size());
for(int i = 0; i < doc.numPoints; i++){ for(int i = 0; i < doc.numPoints; i++){
MVertex *here = (MVertex *)doc.points[i].data; MVertex *here = (MVertex *)doc.points[i].data;
int num = here->getIndex(); int num = here->getIndex();
...@@ -733,7 +733,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER, ...@@ -733,7 +733,7 @@ static bool gmsh2DMeshGenerator(GFace *gf, int RECUR_ITER,
//outputScalarField(m->triangles, "beforeswop.pos",1); //outputScalarField(m->triangles, "beforeswop.pos",1);
Msg::Debug("Delaunizing the initial mesh"); Msg::Debug("Delaunizing the initial mesh");
gmshDelaunayizeBDS(gf, *m, nb_swap); gmshDelaunayizeBDS(gf, *m, nb_swap);
// outputScalarField(m->triangles, "afterswop.pos",0) //outputScalarField(m->triangles, "afterswop.pos",0);
Msg::Debug("Starting to add internal points"); Msg::Debug("Starting to add internal points");
// start mesh generation // start mesh generation
......
...@@ -449,9 +449,9 @@ void splitEdgePass(GFace *gf, BDS_Mesh &m, double MAXE_, int &nb_split) ...@@ -449,9 +449,9 @@ void splitEdgePass(GFace *gf, BDS_Mesh &m, double MAXE_, int &nb_split)
if (!e->deleted){ if (!e->deleted){
const double coord = 0.5; const double coord = 0.5;
BDS_Point *mid ; BDS_Point *mid ;
mid = m.add_point(++m.MAXPOINTNUMBER, double U = coord * e->p1->u + (1 - coord) * e->p2->u;
coord * e->p1->u + (1 - coord) * e->p2->u, double V = coord * e->p1->v + (1 - coord) * e->p2->v;
coord * e->p1->v + (1 - coord) * e->p2->v,gf); mid = m.add_point(++m.MAXPOINTNUMBER, U , V, gf);
mid->lcBGM() = BGM_MeshSize mid->lcBGM() = BGM_MeshSize
(gf, (gf,
(coord * e->p1->u + (1 - coord) * e->p2->u)*m.scalingU, (coord * e->p1->u + (1 - coord) * e->p2->u)*m.scalingU,
...@@ -534,11 +534,11 @@ void smoothVertexPass(GFace *gf, BDS_Mesh &m, int &nb_smooth, bool q) ...@@ -534,11 +534,11 @@ void smoothVertexPass(GFace *gf, BDS_Mesh &m, int &nb_smooth, bool q)
void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT, void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT,
const bool computeNodalSizeField) const bool computeNodalSizeField)
{ {
int IT = 0; int IT = 0;
int MAXNP = m.MAXPOINTNUMBER; int MAXNP = m.MAXPOINTNUMBER;
// classify correctly the embedded vertices // classify correctly the embedded vertices
// use a negative model face number to avoid // use a negative model face number to avoid
// mesh motion // mesh motion
...@@ -554,6 +554,7 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT, ...@@ -554,6 +554,7 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT,
++itvx; ++itvx;
} }
// IF ASKED , compute nodal size field using 1D Mesh // IF ASKED , compute nodal size field using 1D Mesh
if (computeNodalSizeField){ if (computeNodalSizeField){
...@@ -597,6 +598,7 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT, ...@@ -597,6 +598,7 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT,
int NN1 = m.edges.size(); int NN1 = m.edges.size();
int NN2 = 0; int NN2 = 0;
std::list<BDS_Edge*>::iterator it = m.edges.begin(); std::list<BDS_Edge*>::iterator it = m.edges.begin();
while (1){ while (1){
if (NN2++ >= NN1)break; if (NN2++ >= NN1)break;
if (!(*it)->deleted){ if (!(*it)->deleted){
...@@ -608,6 +610,8 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT, ...@@ -608,6 +610,8 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT,
} }
++it; ++it;
} }
if ((minL > MINE_ && maxL < MAXE_) || IT > (abs(NIT))) break; if ((minL > MINE_ && maxL < MAXE_) || IT > (abs(NIT))) break;
double maxE = MAXE_; double maxE = MAXE_;
double minE = MINE_; double minE = MINE_;
...@@ -623,7 +627,7 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT, ...@@ -623,7 +627,7 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT,
double t5 = Cpu(); double t5 = Cpu();
smoothVertexPass(gf, m, nb_smooth, false); smoothVertexPass(gf, m, nb_smooth, false);
double t6 = Cpu(); double t6 = Cpu();
// swapEdgePass ( gf, m, nb_swap); swapEdgePass ( gf, m, nb_swap);
double t7 = Cpu(); double t7 = Cpu();
// clean up the mesh // clean up the mesh
t_spl += t2 - t1; t_spl += t2 - t1;
...@@ -633,13 +637,16 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT, ...@@ -633,13 +637,16 @@ void gmshRefineMeshBDS(GFace *gf, BDS_Mesh &m, const int NIT,
t_col += t4 - t3; t_col += t4 - t3;
t_sm += t6 - t5; t_sm += t6 - t5;
m.cleanup(); m.cleanup();
IT++; IT++;
Msg::Debug(" iter %3d minL %8.3f/%8.3f maxL %8.3f/%8.3f : " Msg::Debug(" iter %3d minL %8.3f/%8.3f maxL %8.3f/%8.3f : "
"%6d splits, %6d swaps, %6d collapses, %6d moves", "%6d splits, %6d swaps, %6d collapses, %6d moves",
IT, minL, minE, maxL, maxE, nb_split, nb_swap, nb_collaps, nb_smooth); IT, minL, minE, maxL, maxE, nb_split, nb_swap, nb_collaps, nb_smooth);
if (nb_split == 0 && nb_collaps == 0) break; if (nb_split == 0 && nb_collaps == 0) break;
} }
double t_total = t_spl + t_sw + t_col + t_sm; double t_total = t_spl + t_sw + t_col + t_sm;
if(!t_total) t_total = 1.e-6; if(!t_total) t_total = 1.e-6;
Msg::Debug(" ---------------------------------------"); Msg::Debug(" ---------------------------------------");
......
lc = 0.5;
Point(1) = {0, 0, 0, lc};
Point(2) = {0, 2, 0, lc};
Point(3) = {2, 0, 0, lc};
Point(4) = {0, -2, 0, lc};
Point(5) = {-2, 0, 0, lc};
Circle(1) = {2, 1, 5};
Circle(2) = {5, 1, 4};
Circle(3) = {4, 1, 3};
Circle(4) = {3, 1, 2};
Point(22) = {0, 1, 0, lc};
Point(33) = {1, 0, 0, lc};
Point(44) = {0, -1, 0, lc};
Point(55) = {-1, 0, 0, lc};
Circle(11) = {22, 1, 55};
Circle(22) = {55, 1, 44};
Circle(33) = {44, 1, 33};
Circle(44) = {33, 1, 22};
Line Loop(45) = {44, 11, 22, 33};
Line Loop(46) = {4, 1, 2, 3};
Plane Surface(47) = {45, 46};
Mesh.CharacteristicLengthFactor=0.8;
Merge "circle2.msh";
CreateTopology;
Compound Line (100) = {1,2,3,4};
Compound Line (200) = {11,22,33,44};
Compound Surface(1000)={47} Boundary {{1,2,3,4}, {11,22,33,44}};
...@@ -3,11 +3,6 @@ Mesh.CharacteristicLengthFactor=6; ...@@ -3,11 +3,6 @@ Mesh.CharacteristicLengthFactor=6;
Merge "square_CLASS.msh"; Merge "square_CLASS.msh";
CreateTopology; CreateTopology;
//Compound Line(100)={1};
//Compound Line(200)={2};
//Compound Line(300)={3};
//Compound Line(400)={4};
Compound Line(150)={1,2}; Compound Line(150)={1,2};
Compound Line(160)={3,4}; Compound Line(160)={3,4};
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment