diff --git a/Common/GmshDefines.h b/Common/GmshDefines.h
index 265b7ad6651a013ceaa58dae12609d98217e7e86..4920f3bdb785a588fc807e2de420eb8da623d659 100644
--- a/Common/GmshDefines.h
+++ b/Common/GmshDefines.h
@@ -52,11 +52,11 @@
 
 // Geometric entities
 #define ENT_NONE     0
-#define ENT_POINT    1
-#define ENT_LINE     2
-#define ENT_SURFACE  3
-#define ENT_VOLUME   4
-#define ENT_ALL      5
+#define ENT_POINT    (1<<0)
+#define ENT_LINE     (1<<1)
+#define ENT_SURFACE  (1<<2)
+#define ENT_VOLUME   (1<<3)
+#define ENT_ALL      (ENT_POINT | ENT_LINE | ENT_SURFACE | ENT_VOLUME)
 
 #define ELEMENTARY 1
 #define PHYSICAL   2
diff --git a/Common/Options.cpp b/Common/Options.cpp
index 0f4ccad4e509da099dcbcbbe29f1514e0a75043a..cad50ef4f90cea7aaadb7e285a5b58a51c3d3b9f 100644
--- a/Common/Options.cpp
+++ b/Common/Options.cpp
@@ -1,4 +1,4 @@
-// $Id: Options.cpp,v 1.301 2006-08-19 18:48:06 geuzaine Exp $
+// $Id: Options.cpp,v 1.302 2006-08-22 01:58:32 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -138,7 +138,7 @@ void Init_Options(int num)
   CTX.post.force_num = 0;
   CTX.print.gl_fonts = 1;
   CTX.threads_lock = 0; // very primitive locking
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = 0;
   CTX.mesh.oldxtrude = CTX.mesh.oldxtrude_recombine = 0; // old extrusion mesh generator
   CTX.mesh.bgmesh_type = WITHPOINTS;
   CTX.post.combine_time = 0; // try to combine_time views at startup
@@ -4046,7 +4046,8 @@ double opt_mesh_tangents(OPT_ARGS_NUM)
 double opt_mesh_explode(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.explode != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.explode != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.explode = val;
   }
 #if defined(HAVE_FLTK)
@@ -4092,7 +4093,8 @@ double opt_mesh_rand_factor(OPT_ARGS_NUM)
 double opt_mesh_quality_type(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.quality_type != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.quality_type != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.quality_type = (int)val;
     if(CTX.mesh.quality_type < 0 || CTX.mesh.quality_type > 2)
       CTX.mesh.quality_type = 0;
@@ -4108,7 +4110,8 @@ double opt_mesh_quality_type(OPT_ARGS_NUM)
 double opt_mesh_quality_inf(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.quality_inf != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.quality_inf != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.quality_inf = val;
   }
 #if defined(HAVE_FLTK)
@@ -4121,7 +4124,8 @@ double opt_mesh_quality_inf(OPT_ARGS_NUM)
 double opt_mesh_quality_sup(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.quality_sup != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.quality_sup != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.quality_sup = val;
   }
 #if defined(HAVE_FLTK)
@@ -4134,7 +4138,8 @@ double opt_mesh_quality_sup(OPT_ARGS_NUM)
 double opt_mesh_radius_inf(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.radius_inf != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.radius_inf != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.radius_inf = val;
   }
 #if defined(HAVE_FLTK)
@@ -4147,7 +4152,8 @@ double opt_mesh_radius_inf(OPT_ARGS_NUM)
 double opt_mesh_radius_sup(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.radius_sup != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.radius_sup != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.radius_sup = val;
   }
 #if defined(HAVE_FLTK)
@@ -4187,7 +4193,8 @@ double opt_mesh_points(OPT_ARGS_NUM)
 double opt_mesh_lines(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.lines != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.lines != val) 
+      CTX.mesh.changed = ENT_LINE;
     CTX.mesh.lines = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4200,7 +4207,8 @@ double opt_mesh_lines(OPT_ARGS_NUM)
 double opt_mesh_triangles(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.triangles != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.triangles != val) 
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.mesh.triangles = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4217,7 +4225,8 @@ double opt_mesh_triangles(OPT_ARGS_NUM)
 double opt_mesh_quadrangles(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.quadrangles != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.quadrangles != val) 
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.mesh.quadrangles = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4234,7 +4243,8 @@ double opt_mesh_quadrangles(OPT_ARGS_NUM)
 double opt_mesh_tetrahedra(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.tetrahedra != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.tetrahedra != val) 
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.mesh.tetrahedra = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4251,7 +4261,8 @@ double opt_mesh_tetrahedra(OPT_ARGS_NUM)
 double opt_mesh_hexahedra(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.hexahedra != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.hexahedra != val) 
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.mesh.hexahedra = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4268,7 +4279,8 @@ double opt_mesh_hexahedra(OPT_ARGS_NUM)
 double opt_mesh_prisms(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.prisms != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.prisms != val) 
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.mesh.prisms = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4285,7 +4297,8 @@ double opt_mesh_prisms(OPT_ARGS_NUM)
 double opt_mesh_pyramids(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.pyramids != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.pyramids != val)
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.mesh.pyramids = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4302,7 +4315,8 @@ double opt_mesh_pyramids(OPT_ARGS_NUM)
 double opt_mesh_surfaces_edges(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.surfaces_edges != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.surfaces_edges != val) 
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.mesh.surfaces_edges = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4315,7 +4329,8 @@ double opt_mesh_surfaces_edges(OPT_ARGS_NUM)
 double opt_mesh_surfaces_faces(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.surfaces_faces != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.surfaces_faces != val)
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.mesh.surfaces_faces = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4328,7 +4343,8 @@ double opt_mesh_surfaces_faces(OPT_ARGS_NUM)
 double opt_mesh_volumes_edges(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.volumes_edges != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.volumes_edges != val)
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.mesh.volumes_edges = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4341,7 +4357,8 @@ double opt_mesh_volumes_edges(OPT_ARGS_NUM)
 double opt_mesh_volumes_faces(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.volumes_faces != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.volumes_faces != val)
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.mesh.volumes_faces = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4464,7 +4481,8 @@ double opt_mesh_line_type(OPT_ARGS_NUM)
 double opt_mesh_reverse_all_normals(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.reverse_all_normals != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.reverse_all_normals != val)
+      CTX.mesh.changed = ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.reverse_all_normals = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4477,7 +4495,8 @@ double opt_mesh_reverse_all_normals(OPT_ARGS_NUM)
 double opt_mesh_smooth_normals(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.smooth_normals != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.smooth_normals != val)
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.mesh.smooth_normals = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4490,7 +4509,8 @@ double opt_mesh_smooth_normals(OPT_ARGS_NUM)
 double opt_mesh_angle_smooth_normals(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET) {
-    if(CTX.mesh.angle_smooth_normals != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.angle_smooth_normals != val)
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.mesh.angle_smooth_normals = val;
   }
 #if defined(HAVE_FLTK)
@@ -4814,7 +4834,8 @@ double opt_mesh_initial_only(OPT_ARGS_NUM)
 double opt_mesh_use_cut_plane(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.use_cut_plane != (int)val) CTX.mesh.changed = 1;
+    if(CTX.mesh.use_cut_plane != (int)val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.use_cut_plane = (int)val;
 #if defined(HAVE_FLTK)
     if(WID){
@@ -4837,7 +4858,8 @@ double opt_mesh_use_cut_plane(OPT_ARGS_NUM)
 double opt_mesh_cut_plane_draw_intersect(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.cut_plane_draw_intersect != (int)val) CTX.mesh.changed = 1;
+    if(CTX.mesh.cut_plane_draw_intersect != (int)val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.cut_plane_draw_intersect = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4850,7 +4872,8 @@ double opt_mesh_cut_plane_draw_intersect(OPT_ARGS_NUM)
 double opt_mesh_cut_plane_only_volume(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.cut_plane_only_volume != (int)val) CTX.mesh.changed = 1;
+    if(CTX.mesh.cut_plane_only_volume != (int)val)
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.cut_plane_only_volume = (int)val;
   }
 #if defined(HAVE_FLTK)
@@ -4863,7 +4886,8 @@ double opt_mesh_cut_plane_only_volume(OPT_ARGS_NUM)
 double opt_mesh_cut_planea(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.cut_planea != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.cut_planea != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.cut_planea = val;
   }
 #if defined(HAVE_FLTK)
@@ -4876,7 +4900,8 @@ double opt_mesh_cut_planea(OPT_ARGS_NUM)
 double opt_mesh_cut_planeb(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.cut_planeb != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.cut_planeb != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.cut_planeb = val;
   }
 #if defined(HAVE_FLTK)
@@ -4889,7 +4914,8 @@ double opt_mesh_cut_planeb(OPT_ARGS_NUM)
 double opt_mesh_cut_planec(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.cut_planec != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.cut_planec != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.cut_planec = val;
   }
 #if defined(HAVE_FLTK)
@@ -4902,7 +4928,8 @@ double opt_mesh_cut_planec(OPT_ARGS_NUM)
 double opt_mesh_cut_planed(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET){
-    if(CTX.mesh.cut_planed != val) CTX.mesh.changed = 1;
+    if(CTX.mesh.cut_planed != val) 
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.cut_planed = val;
   }
 #if defined(HAVE_FLTK)
@@ -4932,7 +4959,7 @@ double opt_mesh_color_carousel(OPT_ARGS_NUM)
     // vertex arrays need to be regenerated only when we color by
     // element type or by partition
     if(CTX.mesh.color_carousel != (int)val && (val == 0. || val == 3.))
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.mesh.color_carousel = (int)val;
     if(CTX.mesh.color_carousel < 0 || CTX.mesh.color_carousel > 3)
       CTX.mesh.color_carousel = 0;
@@ -7196,7 +7223,7 @@ unsigned int opt_mesh_color_lines(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.line != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.color.mesh.line = val;
   }
 #if defined(HAVE_FLTK)
@@ -7211,7 +7238,7 @@ unsigned int opt_mesh_color_triangles(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.triangle != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.color.mesh.triangle = val;
   }
 #if defined(HAVE_FLTK)
@@ -7226,7 +7253,7 @@ unsigned int opt_mesh_color_quadrangles(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.quadrangle != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_SURFACE;
     CTX.color.mesh.quadrangle = val;
   }
 #if defined(HAVE_FLTK)
@@ -7241,7 +7268,7 @@ unsigned int opt_mesh_color_tetrahedra(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.tetrahedron != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.color.mesh.tetrahedron = val;
   }
 #if defined(HAVE_FLTK)
@@ -7256,7 +7283,7 @@ unsigned int opt_mesh_color_hexahedra(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.hexahedron != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.color.mesh.hexahedron = val;
   }
 #if defined(HAVE_FLTK)
@@ -7271,7 +7298,7 @@ unsigned int opt_mesh_color_prisms(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.prism != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.color.mesh.prism = val;
   }
 #if defined(HAVE_FLTK)
@@ -7286,7 +7313,7 @@ unsigned int opt_mesh_color_pyramid(OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // element type
     if(CTX.color.mesh.pyramid != val && CTX.mesh.color_carousel == 0)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_VOLUME;
     CTX.color.mesh.pyramid = val;
   }
 #if defined(HAVE_FLTK)
@@ -7323,7 +7350,7 @@ unsigned int opt_mesh_color_(int i, OPT_ARGS_COL)
     // vertex arrays need to be regenerated only when we color by
     // partition
     if(CTX.color.mesh.carousel[i] != val && CTX.mesh.color_carousel == 3)
-      CTX.mesh.changed = 1;
+      CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     CTX.color.mesh.carousel[i] = val;
   }
 #if defined(HAVE_FLTK)
diff --git a/Common/VertexArray.cpp b/Common/VertexArray.cpp
index 5bf9c69eadbad2f5bc964be13f9f69756fbacf82..2b610dcd875990d45a75437f9177a629cd136439 100644
--- a/Common/VertexArray.cpp
+++ b/Common/VertexArray.cpp
@@ -1,4 +1,4 @@
-// $Id: VertexArray.cpp,v 1.12 2006-08-16 05:25:22 geuzaine Exp $
+// $Id: VertexArray.cpp,v 1.13 2006-08-22 01:58:32 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -54,44 +54,62 @@ int VertexArray::n()
   return List_Nbr(vertices) / 3;
 }
 
+template<class T>
+inline void List_Add_3(List_T * liste, T *data1, T *data2, T *data3)
+{
+  liste->n += 3;
+  List_Realloc(liste, liste->n);
+  liste->isorder = 0;
+  T* dest1 = (T*)&liste->array[(liste->n - 3) * liste->size];
+  T* dest2 = (T*)&liste->array[(liste->n - 2) * liste->size];
+  T* dest3 = (T*)&liste->array[(liste->n - 1) * liste->size];
+  *dest1 = *data1;
+  *dest2 = *data2;
+  *dest3 = *data3;
+}
+
+template<class T>
+inline void List_Add_4(List_T * liste, T *data1, T *data2, T *data3, T *data4)
+{
+  liste->n += 4;
+  List_Realloc(liste, liste->n);
+  liste->isorder = 0;
+  T* dest1 = (T*)&liste->array[(liste->n - 4) * liste->size];
+  T* dest2 = (T*)&liste->array[(liste->n - 3) * liste->size];
+  T* dest3 = (T*)&liste->array[(liste->n - 2) * liste->size];
+  T* dest4 = (T*)&liste->array[(liste->n - 1) * liste->size];
+  *dest1 = *data1;
+  *dest2 = *data2;
+  *dest3 = *data3;
+  *dest4 = *data4;
+}
+
 void VertexArray::add(float x, float y, float z, 
 		      float n0, float n1, float n2, unsigned int col)
 {
-  List_Add(vertices, &x);
-  List_Add(vertices, &y);
-  List_Add(vertices, &z);
+  List_Add_3(vertices, &x, &y, &z);
 
   char c0 = float2char(n0);
   char c1 = float2char(n1);
   char c2 = float2char(n2);
-  List_Add(normals, &c0);
-  List_Add(normals, &c1);
-  List_Add(normals, &c2);
+  List_Add_3(normals, &c0, &c1, &c2);
 
   unsigned char r = CTX.UNPACK_RED(col);
   unsigned char g = CTX.UNPACK_GREEN(col);
   unsigned char b = CTX.UNPACK_BLUE(col);
   unsigned char a = CTX.UNPACK_ALPHA(col);
-  List_Add(colors, &r);
-  List_Add(colors, &g);
-  List_Add(colors, &b);
-  List_Add(colors, &a);
+  List_Add_4(colors, &r, &g, &b, &a);
 }
 
 void VertexArray::add(float x, float y, float z, unsigned int col)
 {
-  List_Add(vertices, &x);
-  List_Add(vertices, &y);
-  List_Add(vertices, &z);
+  List_Add_3(vertices, &x, &y, &z);
 
   unsigned char r = CTX.UNPACK_RED(col);
   unsigned char g = CTX.UNPACK_GREEN(col);
   unsigned char b = CTX.UNPACK_BLUE(col);
   unsigned char a = CTX.UNPACK_ALPHA(col);
-  List_Add(colors, &r);
-  List_Add(colors, &g);
-  List_Add(colors, &b);
-  List_Add(colors, &a);
+  List_Add_4(colors, &r, &g, &b, &a);
 }
 
 static double theeye[3];
diff --git a/Fltk/Callbacks.cpp b/Fltk/Callbacks.cpp
index 34bac55ba99ee99fce1230af4d51bface8b8c814..354c0895b74ddac792b4991a02d3633d033eeac4 100644
--- a/Fltk/Callbacks.cpp
+++ b/Fltk/Callbacks.cpp
@@ -1,4 +1,4 @@
-// $Id: Callbacks.cpp,v 1.446 2006-08-20 17:02:27 geuzaine Exp $
+// $Id: Callbacks.cpp,v 1.447 2006-08-22 01:58:32 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -1351,7 +1351,7 @@ void visibility_ok_cb(CALLBACK_ARGS)
   // if the browser is not empty, get the selections made in the
   // browser and apply them into the model
   if(VisibilityManager::instance()->getNumEntities()){
-    CTX.mesh.changed = 1;
+    CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
     VisibilityManager::instance()->setAllInvisible(WID->vis_type->value());
     for(int i = 0; i < VisibilityManager::instance()->getNumEntities(); i++)
       if(WID->vis_browser->selected(i + 1))
@@ -1439,7 +1439,7 @@ void visibility_sort_cb(CALLBACK_ARGS)
 
 void visibility_number_cb(CALLBACK_ARGS)
 {
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
 
   int type = (int)(long)data;
   bool val;
@@ -2793,23 +2793,23 @@ void mesh_3d_cb(CALLBACK_ARGS)
   Msg(STATUS2N, " ");
 }
 
-void mesh_stl_cb(CALLBACK_ARGS)
+void mesh_edit_cb(CALLBACK_ARGS)
 {
-  WID->set_context(menu_mesh_stl, 0);
+  WID->set_context(menu_mesh_edit, 0);
+}
+
+void mesh_reparam_cb(CALLBACK_ARGS)
+{
+  printf("LAUCH REPARAMETERIZATION BIGNOU HERE!\n");
 }
 
 void mesh_degree_cb(CALLBACK_ARGS)
 {
-  switch ((long)data) {
-  case 2: 
-    Degre2(GMODEL->getMeshStatus());
-    break;
-  case 1:
-  default:
+  if((long)data == 2)
+    Degre2();
+  else
     Degre1();
-    break;
-  }
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
   Draw();
   Msg(STATUS2N, " ");
 }
@@ -2824,7 +2824,7 @@ void mesh_optimize_cb(CALLBACK_ARGS)
   Optimize_Netgen();
   CTX.threads_lock = 0;
 
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_LINE | ENT_SURFACE | ENT_VOLUME;
   Draw();
   Msg(STATUS2N, " ");
 }
diff --git a/Fltk/Callbacks.h b/Fltk/Callbacks.h
index 1881904c76bcb696759f88a98db374ac44c621c3..3f8d04de5247eedcfef26536c6a3a2a7cc9a4d55 100644
--- a/Fltk/Callbacks.h
+++ b/Fltk/Callbacks.h
@@ -274,10 +274,11 @@ void mesh_define_cb(CALLBACK_ARGS);
 void mesh_1d_cb(CALLBACK_ARGS);
 void mesh_2d_cb(CALLBACK_ARGS); 
 void mesh_3d_cb(CALLBACK_ARGS); 
-void mesh_stl_cb(CALLBACK_ARGS);
+void mesh_edit_cb(CALLBACK_ARGS);
 void mesh_remesh_cb(CALLBACK_ARGS); 
 void mesh_update_edges_cb(CALLBACK_ARGS); 
 void mesh_update_more_edges_cb(CALLBACK_ARGS);
+void mesh_reparam_cb(CALLBACK_ARGS);
 void mesh_degree_cb(CALLBACK_ARGS); 
 void mesh_optimize_cb(CALLBACK_ARGS); 
 void mesh_define_length_cb (CALLBACK_ARGS);
diff --git a/Fltk/GUI.cpp b/Fltk/GUI.cpp
index da09ea2ef3636b71f288c7de44a529b0b4fad882..442c4065596e86e35f18a95e3a54d78f217ed9d0 100644
--- a/Fltk/GUI.cpp
+++ b/Fltk/GUI.cpp
@@ -1,4 +1,4 @@
-// $Id: GUI.cpp,v 1.533 2006-08-19 19:46:07 geuzaine Exp $
+// $Id: GUI.cpp,v 1.534 2006-08-22 01:58:33 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -321,24 +321,25 @@ Context_Item menu_geometry[] = {
 
 Context_Item menu_mesh[] = {
   {"1Mesh", NULL} ,
-  {"Define", (Fl_Callback *)mesh_define_cb} ,
-  {"1D",     (Fl_Callback *)mesh_1d_cb} ,
-  {"2D",     (Fl_Callback *)mesh_2d_cb} , 
-  {"3D",     (Fl_Callback *)mesh_3d_cb} , 
-  {"STL",    (Fl_Callback *)mesh_stl_cb} , 
+  {"Define",       (Fl_Callback *)mesh_define_cb} ,
+  {"Edit",         (Fl_Callback *)mesh_edit_cb} , 
+  {"1D",           (Fl_Callback *)mesh_1d_cb} ,
+  {"2D",           (Fl_Callback *)mesh_2d_cb} , 
+  {"3D",           (Fl_Callback *)mesh_3d_cb} , 
   {"First order",  (Fl_Callback *)mesh_degree_cb, (void*)1 } , 
   {"Second order", (Fl_Callback *)mesh_degree_cb, (void*)2 } , 
 #if defined(HAVE_NETGEN)
   {"Optimize quality", (Fl_Callback *)mesh_optimize_cb} , 
 #endif
-  {"Save",   (Fl_Callback *)mesh_save_cb} ,
+  {"Save",         (Fl_Callback *)mesh_save_cb} ,
   {0} 
 };  
-    Context_Item menu_mesh_stl[] = {
-      {"1Mesh>STL", NULL} ,
-      {"Update edges", (Fl_Callback *)mesh_update_edges_cb} , 
+    Context_Item menu_mesh_edit[] = {
+      {"1Mesh>Edit", NULL} ,
+      {"Update edges",   (Fl_Callback *)mesh_update_edges_cb} ,
+      {"Reparameterize", (Fl_Callback *)mesh_reparam_cb} ,
       // {"Manually add edges", (Fl_Callback *)mesh_update_more_edges_cb} , 
-      {"Remesh", (Fl_Callback *)mesh_remesh_cb} , 
+      {"Remesh",         (Fl_Callback *)mesh_remesh_cb} , 
       {0} 
     };  
     Context_Item menu_mesh_define[] = {
@@ -346,7 +347,7 @@ Context_Item menu_mesh[] = {
       {"Characteristic length", (Fl_Callback *)mesh_define_length_cb  } ,
       {"Recombine",   (Fl_Callback *)mesh_define_recombine_cb  } ,
       {"Transfinite", (Fl_Callback *)mesh_define_transfinite_cb  } , 
-      {"Elliptic", (Fl_Callback *)mesh_define_elliptic_surface_cb  } , 
+      {"Elliptic",    (Fl_Callback *)mesh_define_elliptic_surface_cb  } , 
       {0} 
     };  
         Context_Item menu_mesh_define_transfinite[] = {
diff --git a/Fltk/GUI.h b/Fltk/GUI.h
index 91ddc8ce7ec9f68ded164348f02eefc77323db4f..a500418801da3eea6e3a541806c80feec56eeafe 100644
--- a/Fltk/GUI.h
+++ b/Fltk/GUI.h
@@ -87,7 +87,7 @@ extern        Context_Item menu_geometry_elementary_delete[];
 extern    Context_Item menu_geometry_physical[]; 
 extern        Context_Item menu_geometry_physical_add[]; 
 extern Context_Item menu_mesh[]; 
-extern    Context_Item menu_mesh_stl[]; 
+extern    Context_Item menu_mesh_edit[]; 
 extern    Context_Item menu_mesh_define[]; 
 extern        Context_Item menu_mesh_define_transfinite[]; 
 extern Context_Item menu_solver[]; 
diff --git a/Geo/MElement.h b/Geo/MElement.h
index 8565992b7e1d7469679a811a09b54006d2550707..a2dce869dad819b1ee68ccbdfa64994978b2b9c1 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -149,6 +149,11 @@ class MLine : public MElement {
   {
     _v[0] = v0; _v[1] = v1;
   }
+  MLine(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 2; i++) _v[i] = v[i];
+  }
   ~MLine(){}
   virtual int getDim(){ return 1; }
   inline int getNumVertices(){ return 2; }
@@ -171,6 +176,11 @@ class MLine2 : public MLine {
   {
     _vs[0] = v2;
   }
+  MLine2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MLine(v, num, part)
+  {
+    _vs[0] = v[2];
+  }
   ~MLine2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 3; }
@@ -200,6 +210,11 @@ class MTriangle : public MElement {
   {
     _v[0] = v0; _v[1] = v1; _v[2] = v2;
   }
+  MTriangle(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 3; i++) _v[i] = v[i];
+  }
   ~MTriangle(){}
   virtual int getDim(){ return 2; }
   inline int getNumVertices(){ return 3; }
@@ -229,6 +244,11 @@ class MTriangle2 : public MTriangle {
   {
     _vs[0] = v3; _vs[1] = v4; _vs[2] = v5;
   }
+  MTriangle2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MTriangle(v, num, part)
+  {
+    for(int i = 0; i < 3; i++) _vs[i] = v[3 + i];
+  }
   ~MTriangle2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 6; }
@@ -276,6 +296,11 @@ class MQuadrangle : public MElement {
   {
     _v[0] = v0; _v[1] = v1; _v[2] = v2; _v[3] = v3;
   }
+  MQuadrangle(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 4; i++) _v[i] = v[i];
+  }
   ~MQuadrangle(){}
   virtual int getDim(){ return 2; }
   inline int getNumVertices(){ return 4; }
@@ -302,6 +327,11 @@ class MQuadrangle2 : public MQuadrangle {
   {
     _vs[0] = v4; _v[1] = v5; _v[2] = v6; _v[3] = v7; _v[4] = v8;
   }
+  MQuadrangle2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MQuadrangle(v, num, part)
+  {
+    for(int i = 0; i < 5; i++) _vs[i] = v[4 + i];
+  }
   ~MQuadrangle2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 9; }
@@ -323,6 +353,11 @@ class MTetrahedron : public MElement {
   {
     _v[0] = v0; _v[1] = v1; _v[2] = v2; _v[3] = v3;
   }
+  MTetrahedron(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 4; i++) _v[i] = v[i];
+  }
   ~MTetrahedron(){}
   virtual int getDim(){ return 3; }
   inline int getNumVertices(){ return 4; }
@@ -379,6 +414,11 @@ class MTetrahedron2 : public MTetrahedron {
   {
     _vs[0] = v4; _vs[1] = v5; _vs[2] = v6; _vs[3] = v7; _vs[4] = v8; _vs[5] = v9;
   }
+  MTetrahedron2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MTetrahedron(v, num, part)
+  {
+    for(int i = 0; i < 6; i++) _vs[i] = v[4 + i];
+  }
   ~MTetrahedron2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 10; }
@@ -410,6 +450,11 @@ class MHexahedron : public MElement {
     _v[0] = v0; _v[1] = v1; _v[2] = v2; _v[3] = v3;
     _v[4] = v4; _v[5] = v5; _v[6] = v6; _v[7] = v7;
   }
+  MHexahedron(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 8; i++) _v[i] = v[i];
+  }
   ~MHexahedron(){}
   virtual int getDim(){ return 3; }
   inline int getNumVertices(){ return 8; }
@@ -471,6 +516,11 @@ class MHexahedron2 : public MHexahedron {
     _vs[10] = v18; _vs[11] = v19; _vs[12] = v20; _vs[13] = v21; _vs[14] = v22;
     _vs[15] = v23; _vs[16] = v24; _vs[17] = v25; _vs[18] = v26;
   }
+  MHexahedron2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MHexahedron(v, num, part)
+  {
+    for(int i = 0; i < 19; i++) _vs[i] = v[8 + i];
+  }
   ~MHexahedron2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 27; }
@@ -509,6 +559,11 @@ class MPrism : public MElement {
     _v[0] = v0; _v[1] = v1; _v[2] = v2; _v[3] = v3;
     _v[4] = v4; _v[5] = v5; 
   }
+  MPrism(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 6; i++) _v[i] = v[i];
+  }
   ~MPrism(){}
   virtual int getDim(){ return 3; }
   inline int getNumVertices(){ return 6; }
@@ -572,6 +627,11 @@ class MPrism2 : public MPrism {
     _vs[5] = v11; _vs[6] = v12; _vs[7] = v13; _vs[8] = v14; _vs[9] = v15; 
     _vs[10] = v16; _vs[11] = v17; 
   }
+  MPrism2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MPrism(v, num, part)
+  {
+    for(int i = 0; i < 12; i++) _vs[i] = v[6 + i];
+  }
   ~MPrism2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 18; }
@@ -605,6 +665,11 @@ class MPyramid : public MElement {
   {
     _v[0] = v0; _v[1] = v1; _v[2] = v2; _v[3] = v3; _v[4] = v4;
   }
+  MPyramid(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MElement(num, part)
+  {
+    for(int i = 0; i < 5; i++) _v[i] = v[i];
+  }
   ~MPyramid(){}
   virtual int getDim(){ return 3; }
   inline int getNumVertices(){ return 5; }
@@ -666,6 +731,11 @@ class MPyramid2 : public MPyramid {
     _vs[0] = v5; _vs[1] = v6; _vs[2] = v7; _vs[3] = v8; _vs[4] = v9; 
     _vs[5] = v10; _vs[6] = v11; _vs[7] = v12; _vs[8] = v13; 
   }
+  MPyramid2(std::vector<MVertex*> &v, int num=0, int part=0) 
+    : MPyramid(v, num, part)
+  {
+    for(int i = 0; i < 9; i++) _vs[i] = v[5 + i];
+  }
   ~MPyramid2(){}
   inline int getPolynomialOrder(){ return 2; }
   inline int getNumVertices(){ return 14; }
diff --git a/Graphics/Mesh.cpp b/Graphics/Mesh.cpp
index 68b7a0d5fd41e7907327e756100719de060156c2..b24ef9b33b544c8b3c8dfc0da2f4698043fbab93 100644
--- a/Graphics/Mesh.cpp
+++ b/Graphics/Mesh.cpp
@@ -1,4 +1,4 @@
-// $Id: Mesh.cpp,v 1.181 2006-08-20 03:44:15 geuzaine Exp $
+// $Id: Mesh.cpp,v 1.182 2006-08-22 01:58:34 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -818,16 +818,16 @@ void Draw_Mesh()
     int status = GMODEL->getMeshStatus();
     if(CTX.mesh.changed) {
       Msg(DEBUG, "Mesh has changed: reinitializing drawing data");
-      if(status >= 1)
+      if(status >= 1 && CTX.mesh.changed & ENT_LINE)
 	std::for_each(GMODEL->firstEdge(), GMODEL->lastEdge(), initMeshGEdge());
-      if(status >= 2){
+      if(status >= 2 && CTX.mesh.changed & ENT_SURFACE){
 	if(GMODEL->normals) delete GMODEL->normals;
 	GMODEL->normals = new smooth_normals(CTX.mesh.angle_smooth_normals);
 	if(CTX.mesh.smooth_normals)
 	  std::for_each(GMODEL->firstFace(), GMODEL->lastFace(), initSmoothNormalsGFace());
 	std::for_each(GMODEL->firstFace(), GMODEL->lastFace(), initMeshGFace());
       }
-      if(status >= 3)
+      if(status >= 3 && CTX.mesh.changed & ENT_VOLUME)
 	std::for_each(GMODEL->firstRegion(), GMODEL->lastRegion(), initMeshGRegion());
     }
     if(status >= 0)
diff --git a/Mesh/2D_Recombine.cpp b/Mesh/2D_Recombine.cpp
index 8b0c212f7964441ffe98dc7317c927c3b5c57e10..e9a6b66092881d22b440e9369d6ba35476ff8d53 100644
--- a/Mesh/2D_Recombine.cpp
+++ b/Mesh/2D_Recombine.cpp
@@ -1,4 +1,4 @@
-// $Id: 2D_Recombine.cpp,v 1.29 2006-08-05 10:05:45 geuzaine Exp $
+// $Id: 2D_Recombine.cpp,v 1.30 2006-08-22 01:58:34 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -181,7 +181,7 @@ int Recombine(Tree_T * Vertices, Tree_T * Simplexes, Tree_T * Quadrangles,
   -) Enhancements are performned on the quad mesh.
  */
 
-int Recombine_All (Mesh *M)
+int Recombine_All(Mesh *M)
 {
   if(!Tree_Nbr(M->Surfaces)) 
     return 0;
@@ -205,7 +205,7 @@ int Recombine_All (Mesh *M)
   }
   
   // add 2nd order nodes to all elements 
-  Degre2(2);
+  Degre2();
 
   Msg(STATUS2, "Splitting all elements");
 
diff --git a/Mesh/3D_Mesh_Netgen.cpp b/Mesh/3D_Mesh_Netgen.cpp
index 9a86926c3c89d086181010b655274c2aec5bad92..faa8717169613f17de496c331675b39df2aa211a 100644
--- a/Mesh/3D_Mesh_Netgen.cpp
+++ b/Mesh/3D_Mesh_Netgen.cpp
@@ -1,4 +1,4 @@
-// $Id: 3D_Mesh_Netgen.cpp,v 1.25 2006-08-05 13:31:28 geuzaine Exp $
+// $Id: 3D_Mesh_Netgen.cpp,v 1.26 2006-08-22 01:58:34 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -285,7 +285,7 @@ int Mesh_Netgen(Volume * v)
     return 0;
   }
 
-  Msg(STATUS2, "Meshing volume %d", v->Num);
+  Msg(INFO, "Meshing volume %d", v->Num);
   Netgen ng(v);
   ng.MeshVolume();
   ng.TransferVolumeMesh();
@@ -302,7 +302,7 @@ void Optimize_Netgen(Volume * v)
      Extrude_Mesh(v) || !Tree_Nbr(v->Simplexes))
     return;
 
-  Msg(STATUS2, "Optimizing volume %d", v->Num);
+  Msg(INFO, "Optimizing volume %d", v->Num);
   Netgen ng(v, 1);
   ng.OptimizeVolume();
   ng.TransferVolumeMesh();
@@ -310,7 +310,7 @@ void Optimize_Netgen(Volume * v)
 
 void Optimize_Netgen()
 {
-  Msg(STATUS1, "Mesh optimize 3D...");
+  Msg(STATUS1, "Mesh optimize...");
   double t1 = Cpu();
 
   // cleanup 2nd order vertices, if any
@@ -329,7 +329,7 @@ void Optimize_Netgen()
   List_Delete(list);
 
   double t2 = Cpu();
-  Msg(STATUS1, "Mesh optimize 3D complete (%g s)", t2 - t1);
+  Msg(INFO, "Mesh optimize complete (%g s)", t2 - t1);
 }
 
 #endif // !HAVE_NETGEN
diff --git a/Mesh/DiscreteSurface.cpp b/Mesh/DiscreteSurface.cpp
index 8b2504678c1da76e1877541c8f2cb8643794ada3..3652282f77610d02861f24588e39ac769cc940bd 100644
--- a/Mesh/DiscreteSurface.cpp
+++ b/Mesh/DiscreteSurface.cpp
@@ -1,4 +1,4 @@
-// $Id: DiscreteSurface.cpp,v 1.44 2006-08-12 16:16:30 geuzaine Exp $
+// $Id: DiscreteSurface.cpp,v 1.45 2006-08-22 01:58:34 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -366,7 +366,7 @@ void BDS_To_Mesh(Mesh *M)
     ++it;
   }
 
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_SURFACE;
 }
 
 
@@ -410,7 +410,7 @@ int ReMesh()
 
   MeshDiscreteSurface((Surface *) 0);
   CreateVolumeWithAllSurfaces(THEM);
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_SURFACE;
   return 1;
 }
 
@@ -425,7 +425,7 @@ int MeshDiscreteSurface(Surface * s)
                       THEM->bds->LC,
                       CTX.mesh.beta_smooth_metric, CTX.mesh.nb_elem_per_rc);
     if(!THEM->bds_mesh) {
-      Msg(STATUS1, "Remesh 2D...");      
+      Msg(STATUS1, "Remeshing 2D...");      
       double t1 = Cpu();
 
       THEM->bds_mesh = new BDS_Mesh(*(THEM->bds));
@@ -440,7 +440,7 @@ int MeshDiscreteSurface(Surface * s)
       // THEM->bds_mesh->save_gmsh_format("3.msh");
 
       double t2 = Cpu();
-      Msg(STATUS1, "Remesh 2D complete (%g s)", t2 - t1);
+      Msg(INFO, "Remesh 2D complete (%g s)", t2 - t1);
       //      NITER++;
       return 1;
     }
diff --git a/Mesh/Generator.cpp b/Mesh/Generator.cpp
index 7231ecee198ec68ac4a8e137cb08c08046ee02f4..550ca47179e2d63d72d84d410be1391a46fe9d68 100644
--- a/Mesh/Generator.cpp
+++ b/Mesh/Generator.cpp
@@ -1,4 +1,4 @@
-// $Id: Generator.cpp,v 1.93 2006-08-12 17:44:25 geuzaine Exp $
+// $Id: Generator.cpp,v 1.94 2006-08-22 01:58:35 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -458,7 +458,7 @@ void Init_Mesh()
   THEM->status = 0;
 
   CTX.mesh.bgmesh_type = WITHPOINTS;
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_ALL;
 }
 
 void mai3d(int ask)
@@ -482,30 +482,30 @@ void mai3d(int ask)
 
   // 1D mesh
   if((ask > old && ask > 0 && old < 1) || (ask < old && ask > 0)) {
-    Msg(STATUS1, "Mesh 1D...");
+    Msg(STATUS1, "Meshing 1D...");
     if(GMODEL->getMeshStatus() > 1){
       OpenProblem(CTX.filename);
     }
     Maillage_Dimension_1();
-    Msg(STATUS1, "Mesh 1D complete (%g s)", CTX.mesh_timer[0]);
+    Msg(INFO, "Mesh 1D complete (%g s)", CTX.mesh_timer[0]);
   }
 
   // 2D mesh
   if((ask > old && ask > 1 && old < 2) || (ask < old && ask > 1)) {
-    Msg(STATUS1, "Mesh 2D...");
+    Msg(STATUS1, "Meshing 2D...");
     if(GMODEL->getMeshStatus() > 2) {
       OpenProblem(CTX.filename);
       Maillage_Dimension_1();
     }
     Maillage_Dimension_2();
-    Msg(STATUS1, "Mesh 2D complete (%g s)", CTX.mesh_timer[1]);
+    Msg(INFO, "Mesh 2D complete (%g s)", CTX.mesh_timer[1]);
   }
 
   // 3D mesh
   if((ask > old && ask > 2 && old < 3) || (ask < old && ask > 2)) {
-    Msg(STATUS1, "Mesh 3D...");
+    Msg(STATUS1, "Meshing 3D...");
     Maillage_Dimension_3();
-    Msg(STATUS1, "Mesh 3D complete (%g s)", CTX.mesh_timer[2]);
+    Msg(INFO, "Mesh 3D complete (%g s)", CTX.mesh_timer[2]);
   }
 
   // Optimize quality
@@ -513,13 +513,14 @@ void mai3d(int ask)
     Optimize_Netgen();
 
   // Create second order elements
-  if(GMODEL->getMeshStatus() && CTX.mesh.order == 2)
-    Degre2(GMODEL->getMeshStatus());
+  if(GMODEL->getMeshStatus() && CTX.mesh.order == 2) 
+    Degre2();
 
   // Partition
   if(GMODEL->getMeshStatus() > 1 && CTX.mesh.nbPartitions != 1)
     PartitionMesh(THEM, CTX.mesh.nbPartitions);
 
+  Msg(STATUS1, "Mesh");
   CTX.threads_lock = 0;
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_ALL;
 }
diff --git a/Mesh/Mesh.h b/Mesh/Mesh.h
index 2af1cbb557908fc6e961404e58495421c1f0cc75..8ae5cf31210a389b0cb3bb8c38efceb7e05b58c9 100644
--- a/Mesh/Mesh.h
+++ b/Mesh/Mesh.h
@@ -397,6 +397,6 @@ int Recombine_All(Mesh *M);
 void ApplyLcFactor();
 
 void Degre1();
-void Degre2(int dim);
+void Degre2();
 
 #endif
diff --git a/Mesh/SecondOrder.cpp b/Mesh/SecondOrder.cpp
index f737b56fcd6ec5cab69fec2804061851b95956f2..ad20d1ee22b6e08fe0a92b727c7cee59fd78b96f 100644
--- a/Mesh/SecondOrder.cpp
+++ b/Mesh/SecondOrder.cpp
@@ -1,4 +1,4 @@
-// $Id: SecondOrder.cpp,v 1.42 2006-08-20 17:02:28 geuzaine Exp $
+// $Id: SecondOrder.cpp,v 1.43 2006-08-22 01:58:35 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -133,7 +133,6 @@ void getFaceVertices(GFace *gf, MElement *ele, std::vector<MVertex*> &vf,
 		     std::map<std::vector<MVertex*>, MVertex* > &faceVertices,
 		     bool linear)
 {
-  vf.clear();
   for(int i = 0; i < ele->getNumFaces(); i++){
     MFace face = ele->getFace(i);
     if(face.getNumVertices() != 4) continue;
@@ -169,7 +168,6 @@ void getFaceVertices(GRegion *gr, MElement *ele, std::vector<MVertex*> &vf,
 		     std::map<std::vector<MVertex*>, MVertex* > &faceVertices,
 		     bool linear)
 {
-  vf.clear();
   for(int i = 0; i < ele->getNumFaces(); i++){
     MFace face = ele->getFace(i);
     if(face.getNumVertices() != 4) continue;
@@ -192,16 +190,17 @@ void setSecondOrder(GEdge *ge,
 		    std::map<std::pair<MVertex*,MVertex*>, MVertex* > &edgeVertices,
 		    bool linear)
 {
-  std::vector<MLine*> l2;
+  std::vector<MLine*> lines2;
   for(unsigned int i = 0; i < ge->lines.size(); i++){
     MLine *l = ge->lines[i];
     std::vector<MVertex*> ve;
     getEdgeVertices(ge, l, ve, edgeVertices, linear);
-    l2.push_back(new MLine2(l->getVertex(0), l->getVertex(1), ve[0]));
+    lines2.push_back(new MLine2(l->getVertex(0), l->getVertex(1), ve[0]));
     delete l;
   }
-  ge->lines = l2;
-  ge->meshRep->destroy();
+  ge->lines = lines2;
+
+  if(ge->meshRep) ge->meshRep->destroy();
 }
 
 void setSecondOrder(GFace *gf,
@@ -209,30 +208,32 @@ void setSecondOrder(GFace *gf,
 		    std::map<std::vector<MVertex*>, MVertex* > &faceVertices,
 		    bool linear)
 {
-  std::vector<MTriangle*> t2;
+  std::vector<MTriangle*> triangles2;
   for(unsigned int i = 0; i < gf->triangles.size(); i++){
     MTriangle *t = gf->triangles[i];
     std::vector<MVertex*> ve;
     getEdgeVertices(gf, t, ve, edgeVertices, linear);
-    t2.push_back(new MTriangle2(t->getVertex(0), t->getVertex(1), t->getVertex(2),
-				ve[0], ve[1], ve[2]));
+    triangles2.push_back
+      (new MTriangle2(t->getVertex(0), t->getVertex(1), t->getVertex(2),
+		      ve[0], ve[1], ve[2]));
     delete t;
   }
-  gf->triangles = t2;
+  gf->triangles = triangles2;
 
-  std::vector<MQuadrangle*> q2;
+  std::vector<MQuadrangle*> quadrangles2;
   for(unsigned int i = 0; i < gf->quadrangles.size(); i++){
     MQuadrangle *q = gf->quadrangles[i];
     std::vector<MVertex*> ve, vf;
     getEdgeVertices(gf, q, ve, edgeVertices, linear);
     getFaceVertices(gf, q, vf, faceVertices, linear);
-    q2.push_back(new MQuadrangle2(q->getVertex(0), q->getVertex(1), q->getVertex(2),
-				  q->getVertex(3), ve[0], ve[1], ve[2], ve[3], vf[0]));
+    quadrangles2.push_back
+      (new MQuadrangle2(q->getVertex(0), q->getVertex(1), q->getVertex(2),
+			q->getVertex(3), ve[0], ve[1], ve[2], ve[3], vf[0]));
     delete q;
   }
-  gf->quadrangles = q2;
+  gf->quadrangles = quadrangles2;
   
-  gf->meshRep->destroy();
+  if(gf->meshRep) gf->meshRep->destroy();
 }
 
 void setSecondOrder(GRegion *gr,
@@ -240,45 +241,147 @@ void setSecondOrder(GRegion *gr,
 		    std::map<std::vector<MVertex*>, MVertex* > &faceVertices,
 		    bool linear)
 {
-  std::vector<MTetrahedron*> t2;
+  std::vector<MTetrahedron*> tetrahedra2;
   for(unsigned int i = 0; i < gr->tetrahedra.size(); i++){
     MTetrahedron *t = gr->tetrahedra[i];
     std::vector<MVertex*> ve;
     getEdgeVertices(gr, t, ve, edgeVertices, linear);
-    t2.push_back(new MTetrahedron2(t->getVertex(0), t->getVertex(1), t->getVertex(2), 
-				   t->getVertex(3), 
-				   ve[0], ve[1], ve[2], ve[3], ve[4], ve[5]));
+    tetrahedra2.push_back
+      (new MTetrahedron2(t->getVertex(0), t->getVertex(1), t->getVertex(2), 
+			 t->getVertex(3), ve[0], ve[1], ve[2], ve[3], ve[4], ve[5]));
     delete t;
   }
-  gr->tetrahedra = t2;
+  gr->tetrahedra = tetrahedra2;
 
-  // hexa prisms pyr
+  std::vector<MHexahedron*> hexahedra2;
+  for(unsigned int i = 0; i < gr->hexahedra.size(); i++){
+    MHexahedron *h = gr->hexahedra[i];
+    std::vector<MVertex*> ve, vf;
+    getEdgeVertices(gr, h, ve, edgeVertices, linear);
+    getFaceVertices(gr, h, vf, faceVertices, linear);
+    SPoint3 pc = h->barycenter();
+    MVertex *v = new MVertex(pc.x(), pc.y(), pc.z(), gr);
+    gr->mesh_vertices.push_back(v);
+    hexahedra2.push_back
+      (new MHexahedron2(h->getVertex(0), h->getVertex(1), h->getVertex(2), 
+			h->getVertex(3), h->getVertex(4), h->getVertex(5), 
+			h->getVertex(6), h->getVertex(7), ve[0], ve[1], ve[2], 
+			ve[3], ve[4], ve[5], ve[6], ve[7], ve[8], ve[9], ve[10], 
+			ve[11], vf[0], vf[1], vf[2], vf[3], vf[4], vf[5], v));
+    delete h;
+  }
+  gr->hexahedra = hexahedra2;
 
-  gr->meshRep->destroy();
+  std::vector<MPrism*> prisms2;
+  for(unsigned int i = 0; i < gr->prisms.size(); i++){
+    MPrism *p = gr->prisms[i];
+    std::vector<MVertex*> ve, vf;
+    getEdgeVertices(gr, p, ve, edgeVertices, linear);
+    getFaceVertices(gr, p, vf, faceVertices, linear);
+    prisms2.push_back
+      (new MPrism2(p->getVertex(0), p->getVertex(1), p->getVertex(2), 
+		   p->getVertex(3), p->getVertex(4), p->getVertex(5), 
+		   ve[0], ve[1], ve[2], ve[3], ve[4], ve[5], ve[6], ve[7],
+		   ve[8], vf[0], vf[1], vf[2]));
+    delete p;
+  }
+  gr->prisms = prisms2;
+
+  std::vector<MPyramid*> pyramids2;
+  for(unsigned int i = 0; i < gr->pyramids.size(); i++){
+    MPyramid *p = gr->pyramids[i];
+    std::vector<MVertex*> ve, vf;
+    getEdgeVertices(gr, p, ve, edgeVertices, linear);
+    getFaceVertices(gr, p, vf, faceVertices, linear);
+    pyramids2.push_back
+      (new MPyramid2(p->getVertex(0), p->getVertex(1), p->getVertex(2), 
+		     p->getVertex(3), p->getVertex(4), ve[0], ve[1], ve[2], 
+		     ve[3], ve[4], ve[5], ve[6], ve[7], vf[0]));
+    delete p;
+  }
+  gr->pyramids = pyramids2;
+
+  if(gr->meshRep) gr->meshRep->destroy();
 }
 
-// Main routines
+template<class T>
+void setFirstOrder(GEntity *e, std::vector<T*> &elements)
+{
+  std::vector<T*> elements1;
+  for(unsigned int i = 0; i < elements.size(); i++){
+    T *ele = elements[i];
+    int n = ele->getNumVertices() - ele->getNumEdgeVertices() - 
+      ele->getNumFaceVertices() - ele->getNumVolumeVertices();
+    std::vector<MVertex*> v1;
+    for(int j = 0; j < n; j++)
+      v1.push_back(ele->getVertex(j));
+    for(int j = n; j < ele->getNumVertices(); j++)
+      ele->getVertex(j)->setVisibility(-1);
+    elements1.push_back(new T(v1));
+    delete ele;
+  }
+  elements = elements1;
+  
+  if(e->meshRep) e->meshRep->destroy();
+}
+
+void removeSecondOrderVertices(GEntity *e)
+{
+  std::vector<MVertex*> v1;
+  for(unsigned int i = 0; i < e->mesh_vertices.size(); i++){
+    if(e->mesh_vertices[i]->getVisibility() < 0)
+      delete e->mesh_vertices[i];
+    else
+      v1.push_back(e->mesh_vertices[i]);
+  }
+  e->mesh_vertices = v1;
+}
 
 void Degre1()
 {
-  // loop on all elements 
-  // - if polynomialOrder() == 2
-  // - get their edge/face/volume vertices mark them with setVisibility(-1);
-  // - generate a first order element for each element
-  // - swap lists at the end
-  // loop on all vertices and delete all vertices marked (-1)
+  // replace all elements with first order elements and mark all
+  // unused nodes with a -1 visibility flag
+  for(GModel::eiter it = GMODEL->firstEdge(); it != GMODEL->lastEdge(); ++it){
+    setFirstOrder(*it, (*it)->lines);
+  }
+  for(GModel::fiter it = GMODEL->firstFace(); it != GMODEL->lastFace(); ++it){
+    setFirstOrder(*it, (*it)->triangles);
+    setFirstOrder(*it, (*it)->quadrangles);
+  }
+  for(GModel::riter it = GMODEL->firstRegion(); it != GMODEL->lastRegion(); ++it){
+    setFirstOrder(*it, (*it)->tetrahedra);
+    setFirstOrder(*it, (*it)->hexahedra);
+    setFirstOrder(*it, (*it)->prisms);
+    setFirstOrder(*it, (*it)->pyramids);
+  }
+
+  // remove all nodes with a -1 visibility flag
+  for(GModel::eiter it = GMODEL->firstEdge(); it != GMODEL->lastEdge(); ++it)
+    removeSecondOrderVertices(*it);
+  for(GModel::fiter it = GMODEL->firstFace(); it != GMODEL->lastFace(); ++it)
+    removeSecondOrderVertices(*it);
+  for(GModel::riter it = GMODEL->firstRegion(); it != GMODEL->lastRegion(); ++it)
+    removeSecondOrderVertices(*it);
 }
 
-void Degre2(int dim)
+void Degre2()
 {
-  Msg(STATUS1, "Mesh second order...");
+  Msg(STATUS1, "Meshing second order...");
   double t1 = Cpu();
 
   bool linear = true;
   //bool linear = false;
+
   std::map<std::pair<MVertex*,MVertex*>, MVertex* > edgeVertices;
   std::map<std::vector<MVertex*>, MVertex* > faceVertices;
 
+  // replace all elements with second order elements by creating
+  // unique nodes on the edges/faces of the mesh. (To generate nodes
+  // on the exact geometrical edges/faces this assumes that the
+  // geometrical edges/faces are discretized with 1D/2D
+  // elements. I.e., if there are only 3D elements in the mesh then
+  // any new nodes on the boundary will simply be added by linear
+  // interpolation.)
   for(GModel::eiter it = GMODEL->firstEdge(); it != GMODEL->lastEdge(); ++it)
     setSecondOrder(*it, edgeVertices, linear);
   for(GModel::fiter it = GMODEL->firstFace(); it != GMODEL->lastFace(); ++it)
@@ -287,5 +390,6 @@ void Degre2(int dim)
     setSecondOrder(*it, edgeVertices, faceVertices, linear);
 
   double t2 = Cpu();
-  Msg(STATUS1, "Mesh second order complete (%g s)", t2 - t1);
+  Msg(INFO, "Mesh second order complete (%g s)", t2 - t1);
+  Msg(STATUS1, "Mesh");
 }
diff --git a/Parser/OpenFile.cpp b/Parser/OpenFile.cpp
index e3f48e1c50503f63cc54aa9af6056e6768ea4d4d..5a4e69e333cb2f99c4dd6b85758556c0f8608a2b 100644
--- a/Parser/OpenFile.cpp
+++ b/Parser/OpenFile.cpp
@@ -1,4 +1,4 @@
-// $Id: OpenFile.cpp,v 1.115 2006-08-20 14:12:42 geuzaine Exp $
+// $Id: OpenFile.cpp,v 1.116 2006-08-22 01:58:35 geuzaine Exp $
 //
 // Copyright (C) 1997-2006 C. Geuzaine, J.-F. Remacle
 //
@@ -332,7 +332,7 @@ int MergeProblem(char *name, int warn_if_missing)
   }
 
   SetBoundingBox();
-  CTX.mesh.changed = 1;
+  CTX.mesh.changed = ENT_ALL;
   Msg(STATUS2, "Read '%s'", name);
   return status;
 }