From 59cfe84d26ad5f5c468a6850b5d69a179d2258ed Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Sun, 7 Dec 2003 00:23:08 +0000
Subject: [PATCH] Added support to visualize mesh partitions. The partition
 info is currently set as the physical entity number in the MSH file, but it
 would be trivial to change it to something else (e.g., as soon as we
 generalize the MSH format to take arbitrary number of tags per element).

---
 Box/Main.cpp            |   3 +-
 Common/DefaultOptions.h |   2 +-
 Common/Options.cpp      |   9 ++-
 Common/Visibility.cpp   |  49 +++++++++++++--
 Common/Visibility.h     |   1 +
 Fltk/Callbacks.cpp      |  16 +++--
 Fltk/GUI.cpp            |  21 ++++---
 Fltk/Main.cpp           |   3 +-
 Geo/Geo.h               |   1 +
 Graphics/Mesh.cpp       |  57 +++++++++++------
 Mesh/Create.cpp         |  58 +++++++++++++++++-
 Mesh/Create.h           |  13 +++-
 Mesh/Generator.cpp      |   8 ++-
 Mesh/Mesh.h             | 131 +++++++++++++++++++++-------------------
 Mesh/Read_Mesh.cpp      |  22 +++++--
 Mesh/Simplex.cpp        |   4 +-
 Mesh/Simplex.h          |  21 ++++---
 17 files changed, 297 insertions(+), 122 deletions(-)

diff --git a/Box/Main.cpp b/Box/Main.cpp
index 4cac5ffefb..f03db52f50 100644
--- a/Box/Main.cpp
+++ b/Box/Main.cpp
@@ -1,4 +1,4 @@
-// $Id: Main.cpp,v 1.34 2003-12-02 01:31:22 geuzaine Exp $
+// $Id: Main.cpp,v 1.35 2003-12-07 00:23:06 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -106,6 +106,7 @@ int main(int argc, char *argv[])
   M.Surfaces = NULL;
   M.Volumes = NULL;
   M.PhysicalGroups = NULL;
+  M.Partitions = NULL;
   M.Metric = NULL;
 
   signal(SIGINT, Signal);
diff --git a/Common/DefaultOptions.h b/Common/DefaultOptions.h
index c0db1eda14..891a6d366a 100644
--- a/Common/DefaultOptions.h
+++ b/Common/DefaultOptions.h
@@ -696,7 +696,7 @@ StringXNumber MeshOptions_Number[] = {
   { F|O, "CharacteristicLengthFactor" , opt_mesh_lc_factor , 1.0 ,
     "Factor applied to all characteristic lengths (and background meshes)" },
   { F|O, "ColorCarousel" , opt_mesh_color_carousel , 1. ,
-    "Use a `color by region number' coloring scheme" },
+    "Use a `color by region' coloring scheme (0=no, 1=by geometrical entity, 2=by partition)" },
   { F|O, "ColorScheme" , opt_mesh_color_scheme , 0. , 
     "Default mesh color scheme (0, 1 or 2)" },
   { F|O, "ConstrainedBackgroundMesh" , opt_mesh_constrained_bgmesh, 0. ,
diff --git a/Common/Options.cpp b/Common/Options.cpp
index ffdedf9274..3095d2e9bb 100644
--- a/Common/Options.cpp
+++ b/Common/Options.cpp
@@ -1,4 +1,4 @@
-// $Id: Options.cpp,v 1.126 2003-12-04 02:10:30 geuzaine Exp $
+// $Id: Options.cpp,v 1.127 2003-12-07 00:23:06 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -3412,8 +3412,11 @@ double opt_mesh_color_carousel(OPT_ARGS_NUM)
     CTX.mesh.changed = 1;
   }
 #if defined(HAVE_FLTK)
-  if(WID && (action & GMSH_GUI))
-    WID->mesh_butt[17]->value(CTX.mesh.color_carousel);
+  if(WID && (action & GMSH_GUI)){
+    WID->mesh_butt[17]->value(CTX.mesh.color_carousel==0);
+    WID->mesh_butt[18]->value(CTX.mesh.color_carousel==1);
+    WID->mesh_butt[19]->value(CTX.mesh.color_carousel==2);
+  }
 #endif
   return CTX.mesh.color_carousel;
 }
diff --git a/Common/Visibility.cpp b/Common/Visibility.cpp
index 9758ed52df..7787b77381 100644
--- a/Common/Visibility.cpp
+++ b/Common/Visibility.cpp
@@ -1,4 +1,4 @@
-// $Id: Visibility.cpp,v 1.2 2003-12-01 23:57:17 geuzaine Exp $
+// $Id: Visibility.cpp,v 1.3 2003-12-07 00:23:06 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -29,7 +29,7 @@
 
 extern Mesh *THEM;
 
-static List_T *ElementaryEntities = NULL, *PhysicalEntities = NULL;
+static List_T *ElementaryEntities = NULL, *PhysicalEntities = NULL, *Partitions = NULL;
 static Tree_T *VisibleThroughPhysical[4] = { NULL, NULL, NULL, NULL };
 static List_T *NumxSymb = NULL;
 static int Sort = 1;
@@ -50,6 +50,8 @@ int Entity::Num()
     return data.surface ? data.surface->Num : 0;
   case 8:
     return data.volume ? data.volume->Num : 0;
+  case 9:
+    return data.partition ? data.partition->Num : 0;
   default:
     return 0;
   }
@@ -70,6 +72,8 @@ char *Entity::Type()
   case 4:
   case 8:
     return "Volume";
+  case 9:
+    return "Partition";
   default:
     return "Unknown";
   }
@@ -108,6 +112,8 @@ int Entity::Visible()
     return data.surface ? data.surface->Visible : 0;
   case 8:
     return data.volume ? data.volume->Visible : 0;
+  case 9:
+    return data.partition ? data.partition->Visible : 0;
   default:
     return 0;
   }
@@ -214,6 +220,11 @@ void Entity::Visible(int mode)
       break;
     data.volume->Visible = mode;
     break;
+  case 9:
+    if(!data.partition)
+      break;
+    data.partition->Visible = mode;
+    break;
   }
 }
 
@@ -326,6 +337,8 @@ void Entity::RecurVisible()
         break;
       Recur(data.volume, mode);
       break;
+    case 9:
+      break;
     }
   }
 }
@@ -401,6 +414,16 @@ static void AddPhysical(void *a, void *b)
   List_Add(PhysicalEntities, &e);
 }
 
+static void AddPartition(void *a, void *b)
+{
+  MeshPartition *p = *(MeshPartition **) a;
+  Entity e;
+  e.type = 9;
+  e.data.partition = p;
+  e.str = NULL;
+  List_Add(Partitions, &e);
+}
+
 static void AddVertex(void *a, void *b)
 {
   Vertex *v = *(Vertex **) a;
@@ -467,6 +490,11 @@ List_T *GetVisibilityList(int type)
   else
     List_Reset(PhysicalEntities);
 
+  if(!Partitions)
+    Partitions = List_Create(100, 100, sizeof(Entity));
+  else
+    List_Reset(Partitions);
+
   if(!NumxSymb)
     NumxSymb = List_Create(100, 100, sizeof(NxS));
   else
@@ -475,6 +503,7 @@ List_T *GetVisibilityList(int type)
   Tree_Action(Symbol_T, addInNumxSymb);
 
   List_Action(THEM->PhysicalGroups, AddPhysical);
+  List_Action(THEM->Partitions, AddPartition);
 
   Tree_Action(THEM->Points, AddVertex);
   Tree_Action(THEM->Curves, AddCurve);
@@ -483,14 +512,16 @@ List_T *GetVisibilityList(int type)
 
   List_Sort(ElementaryEntities, CompareEntity);
   List_Sort(PhysicalEntities, CompareEntity);
+  List_Sort(Partitions, CompareEntity);
 
   switch (type) {
   case ELEMENTARY:
     return ElementaryEntities;
   case PHYSICAL:
     return PhysicalEntities;
+  case PARTITION:
   default:
-    return NULL;
+    return Partitions;
   }
 }
 
@@ -498,7 +529,13 @@ void ClearVisibilityList(int type)
 {
   int i;
   Entity *e;
-  List_T *list = (type == ELEMENTARY) ? ElementaryEntities : PhysicalEntities;
+  List_T *list;
+  switch(type){
+  case ELEMENTARY : list = ElementaryEntities; break;
+  case PHYSICAL : list = PhysicalEntities; break;
+  case PARTITION : list = Partitions; break;
+  default: return;
+  }
   for(i = 0; i < List_Nbr(list); i++) {
     e = (Entity *) List_Pointer(list, i);
     e->Visible(0);
@@ -648,6 +685,10 @@ static void vis_pyr(void *a, void *b)
 {
   (*(Pyramid **) a)->Visible = vmode;
 }
+static void vis_part(void *a, void *b)
+{
+  (*(MeshPartition **) a)->Visible = vmode;
+}
 static void vis_cur(void *a, void *b)
 {
   (*(Curve **) a)->Visible = vmode;
diff --git a/Common/Visibility.h b/Common/Visibility.h
index b90d51f294..49f1f81db4 100644
--- a/Common/Visibility.h
+++ b/Common/Visibility.h
@@ -37,6 +37,7 @@ public :
     Surface *surface;
     Volume *volume;
     PhysicalGroup *physical;
+    MeshPartition *partition;
   } data;
   char *str;
 
diff --git a/Fltk/Callbacks.cpp b/Fltk/Callbacks.cpp
index 58a66f91a8..1422f978f9 100644
--- a/Fltk/Callbacks.cpp
+++ b/Fltk/Callbacks.cpp
@@ -1,4 +1,4 @@
-// $Id: Callbacks.cpp,v 1.195 2003-12-04 02:10:31 geuzaine Exp $
+// $Id: Callbacks.cpp,v 1.196 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -1038,7 +1038,9 @@ void mesh_options_ok_cb(CALLBACK_ARGS)
   opt_mesh_aspect(0, GMSH_SET,
                   WID->mesh_butt[14]->value()? 0 :
                   WID->mesh_butt[15]->value()? 1 : 2);
-  opt_mesh_color_carousel(0, GMSH_SET, WID->mesh_butt[17]->value());
+  opt_mesh_color_carousel(0, GMSH_SET,
+			  WID->mesh_butt[17]->value()? 0 :
+			  WID->mesh_butt[18]->value()? 1 : 2);
 
   opt_mesh_nb_smoothing(0, GMSH_SET, WID->mesh_value[0]->value());
   opt_mesh_scaling_factor(0, GMSH_SET, WID->mesh_value[1]->value());
@@ -1233,9 +1235,12 @@ void visibility_cb(CALLBACK_ARGS)
   case 0:
     type = ELEMENTARY;
     break;
-  default:
+  case 1:
     type = PHYSICAL;
     break;
+  default :
+    type = PARTITION;
+    break;
   }
   switch (WID->vis_browser_mode->value()) {
   case 0:
@@ -1269,9 +1274,12 @@ void visibility_ok_cb(CALLBACK_ARGS)
   case 0:
     ClearVisibilityList(PHYSICAL);
     break;
-  default:
+  case 1:
     ClearVisibilityList(ELEMENTARY);
     break;
+  default:
+    // partitions: do nothing
+    break;
   }
   switch (WID->vis_browser_mode->value()) {
   case 0:
diff --git a/Fltk/GUI.cpp b/Fltk/GUI.cpp
index ff4fcf228f..b99b63ccb7 100644
--- a/Fltk/GUI.cpp
+++ b/Fltk/GUI.cpp
@@ -1,4 +1,4 @@
-// $Id: GUI.cpp,v 1.264 2003-12-04 02:10:31 geuzaine Exp $
+// $Id: GUI.cpp,v 1.265 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -1824,22 +1824,26 @@ void GUI::create_option_window()
     {
       Fl_Group *o = new Fl_Group(WB, WB + BH, width - 2 * WB, height - 2 * WB - BH, "Colors");
       o->hide();
-      mesh_butt[17] = new Fl_Check_Button(2 * WB, 2 * WB + 1 * BH, BW, BH, "Switch color by entity");
-      mesh_butt[17]->type(FL_TOGGLE_BUTTON);
-      mesh_butt[17]->down_box(TOGGLE_BOX);
-      mesh_butt[17]->selection_color(TOGGLE_COLOR);
+      mesh_butt[17] = new Fl_Check_Button(2 * WB, 2 * WB + 1 * BH, BW, BH, "Color by element type");
+      mesh_butt[18] = new Fl_Check_Button(2 * WB, 2 * WB + 2 * BH, BW, BH, "Color by geometrical entity");
+      mesh_butt[19] = new Fl_Check_Button(2 * WB, 2 * WB + 3 * BH, BW, BH, "Color by partition");
+      for(i = 17; i < 20; i++) {
+        mesh_butt[i]->type(FL_RADIO_BUTTON);
+        mesh_butt[i]->down_box(RADIO_BOX);
+        mesh_butt[i]->selection_color(RADIO_COLOR);
+      }
 
-      mesh_value[12] = new Fl_Value_Input(2 * WB, 2 * WB + 2 * BH, IW, BH, "Predefined color scheme");
+      mesh_value[12] = new Fl_Value_Input(2 * WB, 2 * WB + 4 * BH, IW, BH, "Predefined color scheme");
       mesh_value[12]->minimum(0);
       mesh_value[12]->maximum(2);
       mesh_value[12]->step(1);
       mesh_value[12]->align(FL_ALIGN_RIGHT);
       mesh_value[12]->callback(mesh_options_color_scheme_cb);
 
-      Fl_Scroll *s = new Fl_Scroll(2 * WB, 3 * WB + 3 * BH, IW + 20, height - 5 * WB - 3 * BH);
+      Fl_Scroll *s = new Fl_Scroll(2 * WB, 3 * WB + 5 * BH, IW + 20, height - 5 * WB - 5 * BH);
       i = 0;
       while(MeshOptions_Color[i].str) {
-        mesh_col[i] = new Fl_Button(2 * WB, 3 * WB + (3 + i) * BH, IW, BH, MeshOptions_Color[i].str);
+        mesh_col[i] = new Fl_Button(2 * WB, 3 * WB + (5 + i) * BH, IW, BH, MeshOptions_Color[i].str);
         mesh_col[i]->callback(color_cb, (void *)MeshOptions_Color[i].function);
         i++;
       }
@@ -2823,6 +2827,7 @@ void GUI::create_visibility_window()
   static Fl_Menu_Item type_table[] = {
     {"Elementary", 0, (Fl_Callback *) visibility_cb},
     {"Physical", 0, (Fl_Callback *) visibility_cb},
+    {"Partitions", 0, (Fl_Callback *) visibility_cb},
     {0}
   };
   static Fl_Menu_Item browser_mode_table[] = {
diff --git a/Fltk/Main.cpp b/Fltk/Main.cpp
index 21a80eb352..16a24f5467 100644
--- a/Fltk/Main.cpp
+++ b/Fltk/Main.cpp
@@ -1,4 +1,4 @@
-// $Id: Main.cpp,v 1.53 2003-11-29 19:29:26 geuzaine Exp $
+// $Id: Main.cpp,v 1.54 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -114,6 +114,7 @@ int main(int argc, char *argv[])
   M.Surfaces = NULL;
   M.Volumes = NULL;
   M.PhysicalGroups = NULL;
+  M.Partitions = NULL;
   M.Metric = NULL;
 
   // Signal handling
diff --git a/Geo/Geo.h b/Geo/Geo.h
index c332e6cedf..d04154b1fc 100644
--- a/Geo/Geo.h
+++ b/Geo/Geo.h
@@ -24,6 +24,7 @@
 
 #define ELEMENTARY 1
 #define PHYSICAL   2
+#define PARTITION  3
 
 #define INFILE     1
 #define INSTRING   2
diff --git a/Graphics/Mesh.cpp b/Graphics/Mesh.cpp
index c5c42b8675..6712687824 100644
--- a/Graphics/Mesh.cpp
+++ b/Graphics/Mesh.cpp
@@ -1,4 +1,4 @@
-// $Id: Mesh.cpp,v 1.62 2003-06-14 04:37:42 geuzaine Exp $
+// $Id: Mesh.cpp,v 1.63 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -59,12 +59,11 @@ void draw_polygon_2d(double r, double g, double b, int n,
 
 }
 
-static int iColor;
 static DrawingColor theColor;
 
 void ColorSwitch(int i)
 {
-  glColor4ubv((GLubyte *) & CTX.color.mesh.carousel[i % 10]);
+  glColor4ubv((GLubyte *) & CTX.color.mesh.carousel[abs(i % 10)]);
 }
 
 void Draw_Mesh(Mesh * M)
@@ -123,8 +122,6 @@ void Draw_Mesh(Mesh * M)
     glEnd();
   }
 
-  iColor = 0;
-
   if(CTX.mesh.draw && CTX.render_mode != GMSH_SELECT) {
 
     static int first = 1, listnum;
@@ -197,7 +194,6 @@ void Draw_Mesh_Volumes(void *a, void *b)
 {
   Volume *v;
   v = *(Volume **) a;
-  iColor++;
   theColor = v->Color;
   // FIXME: this is the correct method, but will only work when a
   // coherent datastruct exists for volumes
@@ -222,7 +218,6 @@ void Draw_Mesh_Surfaces(void *a, void *b)
 {
   Surface *s;
   s = *(Surface **) a;
-  iColor++;
   theColor = s->Color;
   if(!(s->Visible & VIS_MESH))
     return;
@@ -250,7 +245,6 @@ void Draw_Mesh_Curves(void *a, void *b)
   c = *(Curve **) a;
   if(c->Num < 0)
     return;
-  iColor++;
   theColor = c->Color;
   if(!(c->Visible & VIS_MESH))
     return;
@@ -317,6 +311,10 @@ void Draw_Simplex_Volume(void *a, void *b)
   if(!s->V[3] || !(s->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, s->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   // FIXME: remove as soon as a coherent structure exists for volumes
   Volume *V;
   if((V = FindVolume(s->iEnt, THEM)) && !(V->Visible & VIS_MESH))
@@ -349,8 +347,9 @@ void Draw_Simplex_Volume(void *a, void *b)
     fulldraw = 1;
   }
 
-  if(CTX.mesh.color_carousel && !fulldraw)
-    ColorSwitch(s->iEnt);
+  if(CTX.mesh.color_carousel && !fulldraw){
+    ColorSwitch((CTX.mesh.color_carousel == 2) ? s->iPart : s->iEnt);
+  }
   else if(fulldraw)
     glColor4ubv((GLubyte *) & CTX.color.mesh.line);
   else
@@ -430,7 +429,7 @@ void Draw_Simplex_Volume(void *a, void *b)
   double n[4], x1x0, y1y0, z1z0, x2x0, y2y0, z2z0;
 
   if(CTX.mesh.color_carousel)
-    ColorSwitch(s->iEnt);
+    ColorSwitch((CTX.mesh.color_carousel == 2) ? s->iPart : s->iEnt);
   else
     glColor4ubv((GLubyte *) & CTX.color.mesh.tetrahedron);
 
@@ -546,7 +545,7 @@ void Draw_Simplex_Surface_Common(Simplex * s, double *pX, double *pY,
       if(theColor.type)
         glColor4ubv((GLubyte *) & theColor.mesh);
       else
-        ColorSwitch(iColor);
+        ColorSwitch((CTX.mesh.color_carousel == 2) ? s->iPart : s->iEnt);
     }
     else {
       glColor4ubv((GLubyte *) & CTX.color.mesh.line);
@@ -589,7 +588,7 @@ void Draw_Simplex_Surface_Common(Simplex * s, double *pX, double *pY,
     if(theColor.type)
       glColor4ubv((GLubyte *) & theColor.mesh);
     else
-      ColorSwitch(iColor);
+      ColorSwitch((CTX.mesh.color_carousel == 2) ? s->iPart : s->iEnt);
   }
   else {
     if(K == 3)
@@ -663,6 +662,10 @@ void Draw_Simplex_Surface_Simple(void *a, void *b)
   if(!s->V[2] || !(s->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, s->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   Draw_Simplex_Surface_Common(s, NULL, NULL, NULL, n);
 }
 
@@ -679,6 +682,10 @@ void Draw_Simplex_Surface(void *a, void *b)
   if(!s->V[2] || !(s->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, s->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   L = (s->VSUP) ? 1 : 0;
   K = (s->V[3]) ? 4 : 3;
 
@@ -769,6 +776,10 @@ void Draw_Simplex_Curves(void *a, void *b)
   if(!(s->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, s->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   Xc = 0.5 * (s->V[0]->Pos.X + s->V[1]->Pos.X);
   Yc = 0.5 * (s->V[0]->Pos.Y + s->V[1]->Pos.Y);
   Zc = 0.5 * (s->V[0]->Pos.Z + s->V[1]->Pos.Z);
@@ -795,7 +806,7 @@ void Draw_Simplex_Curves(void *a, void *b)
     if(theColor.type)
       glColor4ubv((GLubyte *) & theColor.mesh);
     else
-      ColorSwitch(iColor);
+      ColorSwitch((CTX.mesh.color_carousel == 2) ? s->iPart : s->iEnt);
   }
   else
     glColor4ubv((GLubyte *) & CTX.color.mesh.line);
@@ -849,6 +860,10 @@ void Draw_Hexahedron_Volume(void *a, void *b)
   if(!(h->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, h->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   // FIXME: remove as soon as a coherent structure exists for volumes
   Volume *V;
   if((V = FindVolume(h->iEnt, THEM)) && !(V->Visible & VIS_MESH))
@@ -869,7 +884,7 @@ void Draw_Hexahedron_Volume(void *a, void *b)
   }
 
   if(CTX.mesh.color_carousel)
-    ColorSwitch(h->iEnt);
+    ColorSwitch((CTX.mesh.color_carousel == 2) ? h->iPart : h->iEnt);
   else
     glColor4ubv((GLubyte *) & CTX.color.mesh.hexahedron);
 
@@ -990,6 +1005,10 @@ void Draw_Prism_Volume(void *a, void *b)
   if(!(p->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, p->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   // FIXME: remove as soon as a coherent structure exists for volumes
   Volume *V;
   if((V = FindVolume(p->iEnt, THEM)) && !(V->Visible & VIS_MESH))
@@ -1010,7 +1029,7 @@ void Draw_Prism_Volume(void *a, void *b)
   }
 
   if(CTX.mesh.color_carousel)
-    ColorSwitch(p->iEnt);
+    ColorSwitch((CTX.mesh.color_carousel == 2) ? p->iPart : p->iEnt);
   else
     glColor4ubv((GLubyte *) & CTX.color.mesh.prism);
 
@@ -1112,6 +1131,10 @@ void Draw_Pyramid_Volume(void *a, void *b)
   if(!(p->Visible & VIS_MESH))
     return;
 
+  MeshPartition **part = (MeshPartition**)List_Pointer_Test(THEM->Partitions, p->iPart);
+  if(part && !(*part)->Visible)
+    return;
+
   // FIXME: remove as soon as a coherent structure exists for volumes
   Volume *V;
   if((V = FindVolume(p->iEnt, THEM)) && !(V->Visible & VIS_MESH))
@@ -1132,7 +1155,7 @@ void Draw_Pyramid_Volume(void *a, void *b)
   }
 
   if(CTX.mesh.color_carousel)
-    ColorSwitch(p->iEnt);
+    ColorSwitch((CTX.mesh.color_carousel == 2) ? p->iPart : p->iEnt);
   else
     glColor4ubv((GLubyte *) & CTX.color.mesh.pyramid);
 
diff --git a/Mesh/Create.cpp b/Mesh/Create.cpp
index 137ce51eae..1e50f8846e 100644
--- a/Mesh/Create.cpp
+++ b/Mesh/Create.cpp
@@ -1,4 +1,4 @@
-// $Id: Create.cpp,v 1.45 2003-05-29 14:36:56 geuzaine Exp $
+// $Id: Create.cpp,v 1.46 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -202,6 +202,59 @@ void Add_PhysicalGroup(int Num, int typ, List_T * intlist, Mesh * M)
   List_Add(M->PhysicalGroups, &p);
 }
 
+void Free_PhysicalGroup(void *a, void *b)
+{
+  PhysicalGroup *p = *(PhysicalGroup **) a;
+  if(p) {
+    List_Delete(p->Entities);
+    Free(p);
+    p = NULL;
+  }
+}
+
+int Add_MeshPartition(int Num, Mesh * M)
+{
+  MeshPartition P, *p, **pp;
+  p = &P;
+  p->Num = Num;
+  if((pp = (MeshPartition**)List_PQuery(M->Partitions, &p, compareMeshPartitionNum))){
+    return (*pp)->Index;
+  }
+  else{
+    p = (MeshPartition*)Malloc(sizeof(MeshPartition));
+    p->Num = Num;
+    p->Visible = VIS_GEOM | VIS_MESH;
+    p->Index = List_Nbr(M->Partitions);
+    List_Add(M->Partitions, &p);
+    return p->Index;
+  }
+}
+
+void Free_MeshPartition(void *a, void *b)
+{
+  MeshPartition *p = *(MeshPartition **) a;
+  if(p) {
+    Free(p);
+    p = NULL;
+  }
+}
+
+int compareMeshPartitionNum(const void *a, const void *b)
+{
+  MeshPartition *q, *w;
+  q = *(MeshPartition **) a;
+  w = *(MeshPartition **) b;
+  return (q->Num - w->Num);
+}
+
+int compareMeshPartitionIndex(const void *a, const void *b)
+{
+  MeshPartition *q, *w;
+  q = *(MeshPartition **) a;
+  w = *(MeshPartition **) b;
+  return (q->Index - w->Index);
+}
+
 void Add_EdgeLoop(int Num, List_T * intlist, Mesh * M)
 {
   EdgeLoop *pEL;
@@ -692,6 +745,7 @@ Hexahedron *Create_Hexahedron(Vertex * v1, Vertex * v2, Vertex * v3,
 
   h = (Hexahedron *) Malloc(sizeof(Hexahedron));
   h->iEnt = -1;
+  h->iPart = -1;
   h->Num = ++THEM->MaxSimplexNum;
   h->Visible = VIS_MESH;
   h->V[0] = v1;
@@ -723,6 +777,7 @@ Prism *Create_Prism(Vertex * v1, Vertex * v2, Vertex * v3,
 
   p = (Prism *) Malloc(sizeof(Prism));
   p->iEnt = -1;
+  p->iPart = -1;
   p->Num = ++THEM->MaxSimplexNum;
   p->Visible = VIS_MESH;
   p->V[0] = v1;
@@ -752,6 +807,7 @@ Pyramid *Create_Pyramid(Vertex * v1, Vertex * v2, Vertex * v3,
 
   p = (Pyramid *) Malloc(sizeof(Pyramid));
   p->iEnt = -1;
+  p->iPart = -1;
   p->Num = ++THEM->MaxSimplexNum;
   p->Visible = VIS_MESH;
   p->V[0] = v1;
diff --git a/Mesh/Create.h b/Mesh/Create.h
index b5db504342..f61f63d388 100644
--- a/Mesh/Create.h
+++ b/Mesh/Create.h
@@ -37,21 +37,28 @@ int compareAttractor (const void *a, const void *b);
 int compareSurface (const void *a, const void *b);
 int compareVolume (const void *a, const void *b);
 int compareSxF (const void *a, const void *b);
+int compareMeshPartitionNum(const void *a, const void *b);
+int compareMeshPartitionIndex(const void *a, const void *b);
 
 Attractor * Create_Attractor (int Num, double lc1, double lc2, double Radius,
                               Vertex * v, Curve * c, Surface * s);
+
+void Add_EdgeLoop (int Num, List_T * intlist, Mesh * M);
 void Add_SurfaceLoop (int Num, List_T * intlist, Mesh * M);
+
 void Add_PhysicalGroup (int Num, int typ, List_T * intlist, Mesh * M);
-void Add_EdgeLoop (int Num, List_T * intlist, Mesh * M);
+void Free_PhysicalGroup(void *a, void *b);
 
-void End_Curve (Curve * c);
-void End_Surface (Surface * s);
+int  Add_MeshPartition(int Num, Mesh * M);
+void Free_MeshPartition(void *a, void *b);
 
 Curve *Create_Curve (int Num, int Typ, int Order, List_T * Liste,
                      List_T * Knots, int p1, int p2, double u1, double u2);
+void End_Curve (Curve * c);
 void Free_Curve(void *a, void *b);
 
 Surface * Create_Surface (int Num, int Typ);
+void End_Surface (Surface * s);
 void Free_Surface(void *a, void *b);
 
 Volume * Create_Volume (int Num, int Typ);
diff --git a/Mesh/Generator.cpp b/Mesh/Generator.cpp
index 642e32ef54..47052c8a16 100644
--- a/Mesh/Generator.cpp
+++ b/Mesh/Generator.cpp
@@ -1,4 +1,4 @@
-// $Id: Generator.cpp,v 1.45 2003-06-14 04:37:42 geuzaine Exp $
+// $Id: Generator.cpp,v 1.46 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -212,9 +212,12 @@ void Init_Mesh(Mesh * M, int all)
   Tree_Action(M->Volumes, Free_Volume);
   Tree_Delete(M->Volumes);
 
-  //Tree_Action (M->PhysicalGroups, Free_PhysicalGroup); // todo
+  List_Action(M->PhysicalGroups, Free_PhysicalGroup);
   List_Delete(M->PhysicalGroups);
 
+  List_Action(M->Partitions, Free_MeshPartition);
+  List_Delete(M->Partitions);
+
   if(M->Metric) {
     delete M->Metric;
   }
@@ -228,6 +231,7 @@ void Init_Mesh(Mesh * M, int all)
   M->Surfaces = Tree_Create(sizeof(Surface *), compareSurface);
   M->Volumes = Tree_Create(sizeof(Volume *), compareVolume);
   M->PhysicalGroups = List_Create(5, 5, sizeof(PhysicalGroup *));
+  M->Partitions = List_Create(5, 5, sizeof(MeshPartition *));
   M->Metric = new GMSHMetric;
   M->BGM.bgm = NULL;
 
diff --git a/Mesh/Mesh.h b/Mesh/Mesh.h
index 8a729fa9c1..b7e7276b0b 100644
--- a/Mesh/Mesh.h
+++ b/Mesh/Mesh.h
@@ -99,11 +99,11 @@ typedef struct _DELAUNAY Delaunay, *delpeek;
 typedef int PointNumero;
 
 struct _DOC{
-  PointRecord *points;  /* points a trianguler */
-  List_T *hotpoints;    /* hotpoints */
-  int numPoints;        /* nombre de points */
-  int numTriangles;     /* nombre de triangles */
-  Delaunay *delaunay;   /* resultats 2D */
+  PointRecord *points;  // points to triangulate
+  List_T *hotpoints;    // hotpoints
+  int numPoints;        // number of points
+  int numTriangles;     // number of triangles
+  Delaunay *delaunay;   // 2D results
 };
 
 typedef struct{
@@ -124,7 +124,7 @@ typedef struct{
 }IntPoint;
 
 struct _CDLIST{
-  PointNumero point_num; /* numero du point */
+  PointNumero point_num;
   DListPeek next, prev;
 };
 
@@ -197,27 +197,30 @@ class NXE{
 };
 
 typedef struct{
-  int Num;              /* Numero                                       */
-  int iEnt;             /* Entite geometrique                           */
-  char Visible;         /* Visualization flag                           */
-  Vertex *V[8];         /* 8 noeuds                                     */
-  Vertex **VSUP;        /* noeuds supplem pour les elts de degre eleves */
+  int Num; 
+  int iEnt; // parent geometrical entity
+  int iPart; // mesh partition index
+  char Visible;
+  Vertex *V[8];
+  Vertex **VSUP;
 }Hexahedron;
 
 typedef struct{
-  int Num;              /* Numero                                       */
-  int iEnt;             /* Entite geometrique                           */
-  char Visible;         /* Visualization flag                           */
-  Vertex *V[6];         /* 6 noeuds                                     */
-  Vertex **VSUP;        /* noeuds supplem pour les elts de degre eleves */
+  int Num;
+  int iEnt; // parent geometrical entity
+  int iPart; // mesh partition index
+  char Visible;
+  Vertex *V[6];
+  Vertex **VSUP;
 }Prism;
 
 typedef struct{
-  int Num;              /* Numero                                       */
-  int iEnt;             /* Entite geometrique                           */
-  char Visible;         /* Visualization flag                           */
-  Vertex *V[5];         /* 5 noeuds                                     */
-  Vertex **VSUP;        /* noeuds supplem pour les elts de degre eleves */
+  int Num;
+  int iEnt; // parent geometrical entity
+  int iPart; // mesh partition index
+  char Visible;
+  Vertex *V[5];
+  Vertex **VSUP;
 }Pyramid;
 
 typedef struct{
@@ -254,8 +257,8 @@ struct _Surf{
   double RecombineAngle;
   int ipar[5];
   int Nu, Nv;
-  List_T *Generatrices; /* Surface reglee    */
-  List_T *Control_Points; /* Patchs bicubiques */
+  List_T *Generatrices;
+  List_T *Control_Points;
   double plan[3][3];
   double invplan[3][3];
   double a, b, c, d;
@@ -270,14 +273,13 @@ struct _Surf{
   float *ku, *kv, *cp;
   struct _Surf *Support;
   CylParam Cyl;
-  Grid_T Grid;          /* Grille de recherches rapides */
+  Grid_T Grid;  // fast search grid
   ExtrudeParams *Extrude;
   STL_Data *STL; // stl representation of the surface
   POLY_rep *thePolyRep;
-  int Dirty; //flag to prevent any meshing
+  int Dirty; // flag to prevent any meshing
   DrawingColor Color;
-  /// a pointer to a solid model entity
-  void * aSolidModelEntity;
+  void * aSolidModelEntity; // pointer to a solid model entity
 };
 
 typedef struct _Surf Surface;
@@ -299,6 +301,12 @@ typedef struct{
   List_T *Entities;
 }PhysicalGroup;
 
+typedef struct{
+  int Index;
+  int Num;
+  char Visible;
+}MeshPartition;
+
 typedef struct{
   Face F;
   Face Sorted;
@@ -318,7 +326,7 @@ typedef struct {
   Tree_T *Edges;
   Tree_T *Faces;
   Tree_T *Simplexes;
-  Tree_T *Simp_Surf;//for old extrusion mesh generator
+  Tree_T *Simp_Surf; // for old extrusion mesh generator
   Tree_T *Hexahedra;
   Tree_T *Prisms;
   Tree_T *Pyramids;
@@ -331,23 +339,23 @@ typedef struct {
   int iFac;
 }exf_T;
 
-/* Structure intersection arete - Simplexe */
+// Edge-Simplex intersections
 
 typedef struct{
-  int NbIntersect;      /* nombre total d'intersections                   */
-  Edge *e;              /* arete                                          */
-  Simplex *s;           /* simplexe                                       */
-  Face *f;              /* face                                           */
-  int NbVertex;         /* nombre de noeuds du simplexe que coupe l'arete */
-  Vertex *V[12];        /* noeuds du simplexe que coupe l'arete           */
-  int iV[12];           /* noeuds du simplexe que coupe l'arete           */
-  int NbEdge;           /* nombre d'intersections arete-arete             */
-  int E[12];            /* aretes                                         */
-  Vertex *VE[12];       /* noeuds d'intersection                          */
-  int NbFace;           /* nombre d'intersections face-arete              */
-  Face *F[12];          /* faces                                          */
-  int iF[12];           /* faces                                          */
-  Vertex *VF[12];       /* position des points d'intersections face-arete */
+  int NbIntersect;      // nombre total d'intersections
+  Edge *e;              // arete
+  Simplex *s;           // simplexe
+  Face *f;              // face
+  int NbVertex;         // nombre de noeuds du simplexe que coupe l'arete
+  Vertex *V[12];        // noeuds du simplexe que coupe l'arete
+  int iV[12];           // noeuds du simplexe que coupe l'arete
+  int NbEdge;           // nombre d'intersections arete-arete
+  int E[12];            // aretes
+  Vertex *VE[12];       // noeuds d'intersection
+  int NbFace;           // nombre d'intersections face-arete
+  Face *F[12];          // faces
+  int iF[12];           // faces
+  Vertex *VF[12];       // position des points d'intersections face-arete
 }Intersection;
 
 typedef struct _Mesh Mesh;
@@ -386,7 +394,7 @@ typedef struct{
   int degre;
   CircParam Circle;
   char functu[256], functv[256], functw[256];
-  int Dirty; //flag to prevent any meshing
+  int Dirty; // flag to prevent any meshing
   DrawingColor Color;
 }Curve;
 
@@ -421,22 +429,23 @@ class MeshParameters{
 };
 
 struct _Mesh{
-  char name[256];               /* Nom du probleme                       */
-  int status;                   /* Etat actuel du maillage               */
-  Tree_T *Points;               /* Points de controle                    */
-  Tree_T *Vertices;             /* Noeuds du maillage                    */
-  Tree_T *Simplexes;            /* Simplexes                             */
-  Tree_T *Curves;               /* Courbes                               */
-  Tree_T *Surfaces;             /* Surfaces                              */
-  Tree_T *Volumes;              /* Volumes                               */
-  Tree_T *SurfaceLoops;         /* Surface Loops                         */
-  Tree_T *EdgeLoops;            /* Edge Loops                            */
-  List_T *PhysicalGroups;       /* Physical Groups                       */
-  Grid_T Grid;                  /* Grille de recherches rapides          */
-  LcField BGM;                  /* Background mesh                       */
-  double Statistics[50];        /* Mesh statistics                       */
-  int Histogram[3][NB_HISTOGRAM]; /* Quality histograms                 */
-  GMSHMetric *Metric;           /* Metric                                */
+  char name[256];
+  int status; // current state of the mesh
+  Tree_T *Points;
+  Tree_T *Vertices;
+  Tree_T *Simplexes;
+  Tree_T *Curves;
+  Tree_T *Surfaces;
+  Tree_T *Volumes;
+  Tree_T *SurfaceLoops;
+  Tree_T *EdgeLoops;
+  List_T *PhysicalGroups;
+  List_T *Partitions;
+  Grid_T Grid; // fast search grid
+  LcField BGM; // background mesh
+  double Statistics[50]; // mesh statistics
+  int Histogram[3][NB_HISTOGRAM]; // quality histograms
+  GMSHMetric *Metric;
   MeshParameters MeshParams;
   int MaxPointNum, MaxLineNum, MaxLineLoopNum, MaxSurfaceNum;
   int MaxSurfaceLoopNum, MaxVolumeNum, MaxPhysicalNum;
@@ -455,7 +464,7 @@ struct Map{
 };
 
 
-/* public functions */
+// public functions
 
 void mai3d (Mesh * M, int Asked);
 
diff --git a/Mesh/Read_Mesh.cpp b/Mesh/Read_Mesh.cpp
index fd3c802205..50af7084ec 100644
--- a/Mesh/Read_Mesh.cpp
+++ b/Mesh/Read_Mesh.cpp
@@ -1,4 +1,4 @@
-// $Id: Read_Mesh.cpp,v 1.59 2003-11-27 02:33:31 geuzaine Exp $
+// $Id: Read_Mesh.cpp,v 1.60 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -172,7 +172,6 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
 
         fscanf(fp, "%d %d %d %d %d",
 	       &Num, &Type, &Physical, &Elementary, &Nbr_Nodes);
-	//&Num, &Type, &Elementary, &Physical, &Nbr_Nodes) ;
 
         for(j = 0; j < Nbr_Nodes; j++)
           fscanf(fp, "%d", &verts[j].Num);
@@ -266,6 +265,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
           simp = Create_Simplex(vertsp[0], vertsp[1], NULL, NULL);
           simp->Num = Num;
           simp->iEnt = Elementary;
+          simp->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == LGN2){
 	    simp->VSUP = (Vertex **) Malloc(1 * sizeof(Vertex *));
 	    simp->VSUP[0] = vertsp[2];
@@ -281,6 +281,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
           simp = Create_Simplex(vertsp[0], vertsp[1], vertsp[2], NULL);
           simp->Num = Num;
           simp->iEnt = Elementary;
+          simp->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == TRI2){
 	    simp->VSUP = (Vertex **) Malloc(3 * sizeof(Vertex *));
 	    for(i = 0; i < 3; i++){
@@ -300,6 +301,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
           simp = Create_Quadrangle(vertsp[0], vertsp[1], vertsp[2], vertsp[3]);
           simp->Num = Num;
           simp->iEnt = Elementary;
+          simp->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == QUA2){
 	    simp->VSUP = (Vertex **) Malloc(4 * sizeof(Vertex *));
 	    for(i = 0; i < 4; i++){
@@ -320,6 +322,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
           simp = Create_Simplex(vertsp[0], vertsp[1], vertsp[2], vertsp[3]);
           simp->Num = Num;
           simp->iEnt = Elementary;
+          simp->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == TET2){
 	    simp->VSUP = (Vertex **) Malloc(6 * sizeof(Vertex *));
 	    for(i = 0; i < 6; i++){
@@ -340,6 +343,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
                                   vertsp[4], vertsp[5], vertsp[6], vertsp[7]);
           hex->Num = Num;
           hex->iEnt = Elementary;
+          hex->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == HEX2){
 	    hex->VSUP = (Vertex **) Malloc(12 * sizeof(Vertex *));
 	    for(i = 0; i < 12; i++){
@@ -360,6 +364,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
                              vertsp[3], vertsp[4], vertsp[5]);
           pri->Num = Num;
           pri->iEnt = Elementary;
+          pri->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == PRI2){
 	    pri->VSUP = (Vertex **) Malloc(9 * sizeof(Vertex *));
 	    for(i = 0; i < 9; i++){
@@ -380,6 +385,7 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
                                vertsp[3], vertsp[4]);
           pyr->Num = Num;
           pyr->iEnt = Elementary;
+          pyr->iPart = Add_MeshPartition(Physical, M);
 	  if(Type == PYR2){
 	    pyr->VSUP = (Vertex **) Malloc(8 * sizeof(Vertex *));
 	    for(i = 0; i < 8; i++){
@@ -427,20 +433,26 @@ void Read_Mesh_MSH(Mesh * M, FILE * fp)
 
   if(Tree_Nbr(M->Volumes)) {
     M->status = 3;
-    M->Statistics[6] = Tree_Nbr(M->Vertices);   //incorrect, mais...
+    M->Statistics[6] = Tree_Nbr(M->Vertices);   // wrong, but...
   }
   else if(Tree_Nbr(M->Surfaces)) {
     M->status = 2;
-    M->Statistics[5] = Tree_Nbr(M->Vertices);   //incorrect, mais...
+    M->Statistics[5] = Tree_Nbr(M->Vertices);   // wrong, but...
   }
   else if(Tree_Nbr(M->Curves)) {
     M->status = 1;
-    M->Statistics[4] = Tree_Nbr(M->Vertices);   //incorrect, mais...
+    M->Statistics[4] = Tree_Nbr(M->Vertices);   // wrong, but...
   }
   else if(Tree_Nbr(M->Points))
     M->status = 0;
   else
     M->status = -1;
+
+  // 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
+  // access through Llist_Pointer & co.
+  List_Sort(M->Partitions, compareMeshPartitionIndex);
 }
 
 
diff --git a/Mesh/Simplex.cpp b/Mesh/Simplex.cpp
index 668ded107f..37d72412ac 100644
--- a/Mesh/Simplex.cpp
+++ b/Mesh/Simplex.cpp
@@ -1,4 +1,4 @@
-// $Id: Simplex.cpp,v 1.27 2003-03-21 00:52:42 geuzaine Exp $
+// $Id: Simplex.cpp,v 1.28 2003-12-07 00:23:07 geuzaine Exp $
 //
 // Copyright (C) 1997-2003 C. Geuzaine, J.-F. Remacle
 //
@@ -44,6 +44,7 @@ Simplex::Simplex()
   V[0] = V[1] = V[2] = V[3] = NULL;
   S[0] = S[1] = S[2] = S[3] = NULL;
   iEnt = -1;
+  iPart = -1;
   Quality = 0.;
   Num = TotalNumber;
   Visible = VIS_MESH;
@@ -60,6 +61,7 @@ Simplex::Simplex(Vertex * v1, Vertex * v2, Vertex * v3, Vertex * v4)
   Num = TotalNumber;
   THEM->MaxSimplexNum = IMAX(THEM->MaxSimplexNum, Num);
   iEnt = -1;
+  iPart = -1;
   Visible = VIS_MESH;
 }
 
diff --git a/Mesh/Simplex.h b/Mesh/Simplex.h
index f95a3c2bfd..b667ea14bf 100644
--- a/Mesh/Simplex.h
+++ b/Mesh/Simplex.h
@@ -30,16 +30,17 @@ typedef struct {
 class Simplex{
 
 public:
-  int     Num;           /* Numero                                       */
-  int     iEnt;          /* Entite geometrique                           */
-  char    Visible;       /* Visualization flag                           */
-  Face    F[4];          /* 4 faces                                      */
-  Vertex  **VSUP;        /* noeuds supplem pour les elts de degre eleves */
-  Vertex  *V[4];         /* 4 noeuds                                     */
-  double  Quality;       /* Qualite du simplexe                          */
-  Coord   Center;        /* centre du CC                                 */
-  double  Radius;        /* Rayon du CC                                  */
-  Simplex *S[4];         /* 4 Voisins                                    */
+  int     Num;           // Number
+  int     iEnt;          // Elementary geometrical entity
+  int     iPart;         // Mesh partition index
+  char    Visible;       // Visualization flag
+  Face    F[4];          // 4 faces
+  Vertex  **VSUP;        // suppl. nodes for higher order elts
+  Vertex  *V[4];         // 4 nodes
+  double  Quality;       // simplex quality
+  Coord   Center;        // CC center
+  double  Radius;        // CC radius
+  Simplex *S[4];         // 4 neighbours
   static  int TotalNumber;
   static  int TotalAllocated;
   Simplex();
-- 
GitLab