diff --git a/Geo/GeoUtils.cpp b/Geo/GeoUtils.cpp
index 671bd2ffe99d30266d613760279f9c0d7afa2bae..34ad2d10086d6922c00338e784c7aa860002490c 100644
--- a/Geo/GeoUtils.cpp
+++ b/Geo/GeoUtils.cpp
@@ -1,4 +1,4 @@
-// $Id: GeoUtils.cpp,v 1.12 2006-04-16 03:27:30 geuzaine Exp $
+// $Id: GeoUtils.cpp,v 1.13 2006-04-16 18:55:18 geuzaine Exp $
//
// Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
//
@@ -51,22 +51,7 @@ void sortEdgesInLoop(int num, List_T *edges)
for(int k = List_Nbr(rc->Vertices) - 1; k >= 0; k--)
List_Add(c->Vertices, List_Pointer(rc->Vertices, k));
}
- // if the loop is composed of a single discrete curve we check
- // if the first vertex and the last vertex are the same. If
- // not, we assume that they should, and we add the missing end
- // vertex. (This situation can arise due to a limitation in
- // Read_Mesh, where we do a List_Insert() to add the vertices
- // in discrete curves--so that we never get the last vertex
- // for single, closed curves.)
- if(nbEdges == 1 && List_Nbr(c->Vertices)){
- Vertex *first = *(Vertex**)List_Pointer(c->Vertices, 0);
- Vertex *last = *(Vertex**)List_Pointer(c->Vertices, List_Nbr(c->Vertices) - 1);
- if(first != last){
- Msg(INFO, "Adding last vertex in closed discrete curve %d", c->Num);
- List_Add(c->Vertices, &first);
- }
- }
- // setting end points for discrete curves
+ // set end points for discrete curves
if(!c->beg && !c->end && List_Nbr(c->Vertices)){
Vertex *first = *(Vertex**)List_Pointer(c->Vertices, 0);
Vertex *last = *(Vertex**)List_Pointer(c->Vertices, List_Nbr(c->Vertices) - 1);
@@ -80,7 +65,6 @@ void sortEdgesInLoop(int num, List_T *edges)
}
List_Reset(edges);
-
int j = 0, k = 0;
c0 = c1 = *(Curve **) List_Pointer(temp, 0);
List_Add(edges, &c1->Num);
diff --git a/Mesh/Read_Mesh.cpp b/Mesh/Read_Mesh.cpp
index 42adc7324352cffb4d38309d8365f31aa413fe28..1007662f2c170b105235049687ccb889daba020e 100644
--- a/Mesh/Read_Mesh.cpp
+++ b/Mesh/Read_Mesh.cpp
@@ -1,4 +1,4 @@
-// $Id: Read_Mesh.cpp,v 1.103 2006-04-04 04:32:31 geuzaine Exp $
+// $Id: Read_Mesh.cpp,v 1.104 2006-04-16 18:55:18 geuzaine Exp $
//
// Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
//
@@ -20,6 +20,7 @@
// Please report all bugs and problems to <gmsh@geuz.org>.
#include <map>
+#include <vector>
#include "Gmsh.h"
#include "Geo.h"
#include "CAD.h"
@@ -66,6 +67,80 @@ Curve *addElementaryCurve(Mesh * M, int Num)
return c;
}
+class nodex{
+ public:
+ Vertex *v;
+ int left, right; // index left/right segment in simplex list
+ nodex(Vertex *vv) : v(vv), left(-1), right(-1) {}
+};
+
+void addVerticesInCurve(void *a, void *b)
+{
+ Curve *c = *(Curve**)a;
+
+ // check if already done or if nothing to do
+ if(List_Nbr(c->Vertices) || !Tree_Nbr(c->SimplexesBase))
+ return;
+
+ List_T *elm = Tree2List(c->SimplexesBase);
+
+ // construct node x segment connectivity
+ std::map<int, nodex*> nodes;
+ for(int i = 0; i < List_Nbr(elm); i++){
+ SimplexBase *s = *(SimplexBase**)List_Pointer(elm, i);
+ for(int j = 0; j < 2; j++){
+ if(!nodes.count(s->V[j]->Num))
+ nodes[s->V[j]->Num] = new nodex(s->V[j]);
+ if(!j)
+ nodes[s->V[j]->Num]->right = i;
+ else
+ nodes[s->V[j]->Num]->left = i;
+ }
+ }
+
+ // find starting element (or use the first one for closed curve)
+ std::map<int, nodex*>::const_iterator beg = nodes.begin();
+ std::map<int, nodex*>::const_iterator end = nodes.end();
+ int start = 0;
+ while(beg != end){
+ if((*beg).second->left < 0){ // no element to the left
+ start = (*beg).second->right;
+ break;
+ }
+ beg++;
+ }
+
+ // add nodes by following the line segments
+ if(start < 0){
+ Msg(GERROR, "Something is wrong in mesh of curve %d", c->Num);
+ }
+ else{
+ int N = List_Nbr(elm);
+ SimplexBase *s = *(SimplexBase**)List_Pointer(elm, start);
+ List_Add(c->Vertices, &s->V[0]);
+ while(N > 0){
+ if(s->VSUP) List_Add(c->Vertices, &s->VSUP[0]);
+ List_Add(c->Vertices, &s->V[1]);
+ nodex *n = nodes[s->V[1]->Num];
+ if(n->left != start){
+ Msg(GERROR, "Wrong orientation of element %d in curve %d", s->Num, c->Num);
+ break;
+ }
+ start = n->right;
+ if(start >= 0){
+ s = *(SimplexBase**)List_Pointer(elm, start);
+ }
+ else if(N != 1){
+ Msg(GERROR, "Something is wrong in mesh of curve %d", c->Num);
+ break;
+ }
+ N--;
+ }
+ }
+
+ List_Delete(elm);
+}
+
Surface *addElementarySurface(Mesh * M, int Num)
{
Surface *s;
@@ -430,11 +505,11 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
Msg(GERROR, "Line element %d already exists", simp->Num);
Free_SimplexBase(&simp, 0);
}
- else{
- // this can be quite slow if there are many nodes on the curve...
- for(i = 0; i < Nbr_Nodes; i++)
- List_Insert(c->Vertices, &vertsp[i], compareVertex);
- }
+ // we don't insert the vertices in the list of vertices at
+ // this point, since we need the list to ordered (and
+ // consistent!), and simply doing a List_Insert and sorting
+ // according to node numbers is not enough (in addition to
+ // being slow): see addVerticesInCurve() below.
break;
case TRI1:
case TRI2:
@@ -641,6 +716,9 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
else
M->status = -1;
+ // add vertices in curves
+ Tree_Action(M->Curves, addVerticesInCurve);
+
// For efficiency reasons, we store the partition index (and not the
// partition number) in the various mesh elements. We need to
// re-sort the list according to these indices to allow direct