diff --git a/Common/CommandLine.cpp b/Common/CommandLine.cpp
index 68d229c24558295b0e126495d9e59bef3fe09103..4dfc8ca79aeb83460d201ae35f80d0ca5a690ac9 100644
--- a/Common/CommandLine.cpp
+++ b/Common/CommandLine.cpp
@@ -1,4 +1,4 @@
-// $Id: CommandLine.cpp,v 1.44 2004-06-26 17:58:14 geuzaine Exp $
+// $Id: CommandLine.cpp,v 1.45 2004-06-30 17:49:51 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -74,6 +74,9 @@ void Print_Usage(char *name){
   Msg(DIRECT, "  -format msh|unv|gref  set output mesh format (default: msh)");
   Msg(DIRECT, "  -algo iso|tri|aniso   select mesh algorithm (default: iso)");
   Msg(DIRECT, "  -smooth int           set mesh smoothing (default: 0)");
+#if defined(HAVE_NETGEN)
+  Msg(DIRECT, "  -optimize             optimize quality of tetrahedral elements");
+#endif
   Msg(DIRECT, "  -order int            set mesh order (default: 1)");
   Msg(DIRECT, "  -scale float          set global scaling factor (default: 1.0)");
   Msg(DIRECT, "  -meshscale float      set mesh scaling factor (default: 1.0)");
@@ -231,6 +234,10 @@ void Get_Options(int argc, char *argv[], int *nbfiles)
         CTX.mesh.histogram = 1;
         i++;
       }
+      else if(!strcmp(argv[i] + 1, "optimize")) {
+        CTX.mesh.optimize = 1;
+        i++;
+      }
       else if(!strcmp(argv[i] + 1, "option")) {
         i++;
         if(argv[i] != NULL)
diff --git a/Common/Context.h b/Common/Context.h
index b22d843841992b9070cd7ae409b71ea02dd8456b..ba835e87c2b0e6d20f417d9d10360537f58b378f 100644
--- a/Common/Context.h
+++ b/Common/Context.h
@@ -162,6 +162,7 @@ public :
     int points_num, lines_num, surfaces_num, volumes_num;
     int point_type, line_type; // flat or 3D
     double point_size, line_width;
+    int optimize;
     double quality;
     double gamma_inf, gamma_sup, radius_inf, radius_sup;
     double scaling_factor, lc_factor, rand_factor;
diff --git a/Common/DefaultOptions.h b/Common/DefaultOptions.h
index 36697b2ce643be828db8ac53f76847eb4fcb7c6c..e1ea0cd62bd80d337db46b6f0e7b767d4dccdd95 100644
--- a/Common/DefaultOptions.h
+++ b/Common/DefaultOptions.h
@@ -783,6 +783,9 @@ StringXNumber MeshOptions_Number[] = {
   { F|O, "Normals" , opt_mesh_normals , 0.0 ,
     "Display size of normal vectors (in pixels)" }, 
 
+  { F|O, "Optimize" , opt_mesh_optimize , 0. , 
+    "Optimize the mesh using Netgen to improve the quality of tetrahedral elements" },
+
   { F|O, "Points" , opt_mesh_points , 1. , 
     "Display mesh vertices (nodes)?" },
   { F|O, "PointInsertion" , opt_mesh_point_insertion, CENTER_CIRCCIRC ,
@@ -794,8 +797,8 @@ StringXNumber MeshOptions_Number[] = {
   { F|O, "PointType" , opt_mesh_point_type , 0. , 
     "Display mesh vertices as solid color dots (0) or 3D spheres (1)" },
 
-  { F|O, "Quality" , opt_mesh_quality , 0.0 ,
-    "Target quality for tetrahedral elements (not fully functional)" },
+  //{ F|O, "Quality" , opt_mesh_quality , 0.0 ,
+  //  "Target quality for tetrahedral elements (not fully functional)" },
 
   { F|O, "RadiusInf" , opt_mesh_radius_inf , 0.0 , 
     "Only display elements whose Radius is greater than RadiusInf" },
diff --git a/Common/Options.cpp b/Common/Options.cpp
index 7f2dd4cb7700cad1579410377a0602c7f8503f43..e26058d2399b35ff3b7d87e4ea2ccf45c6a7abe0 100644
--- a/Common/Options.cpp
+++ b/Common/Options.cpp
@@ -1,4 +1,4 @@
-// $Id: Options.cpp,v 1.168 2004-06-20 23:25:31 geuzaine Exp $
+// $Id: Options.cpp,v 1.169 2004-06-30 17:49:51 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -3107,6 +3107,17 @@ double opt_geometry_stl_create_physical(OPT_ARGS_NUM)
   return CTX.geom.stl_create_physical;
 }
 
+double opt_mesh_optimize(OPT_ARGS_NUM)
+{
+  if(action & GMSH_SET)
+    CTX.mesh.optimize =(int) val;
+#if defined(HAVE_FLTK)
+  if(WID && (action & GMSH_GUI))
+    WID->mesh_butt[2]->value(CTX.mesh.optimize);
+#endif
+  return CTX.mesh.optimize;
+}
+
 double opt_mesh_quality(OPT_ARGS_NUM)
 {
   if(action & GMSH_SET)
diff --git a/Common/Options.h b/Common/Options.h
index 2268f396cfb21e75741cc73aebc6dd5fb98ea191..910ee4a52db2ca2ea933463da8b94b66f68c49ac 100644
--- a/Common/Options.h
+++ b/Common/Options.h
@@ -336,6 +336,7 @@ double opt_geometry_scaling_factor(OPT_ARGS_NUM);
 double opt_geometry_color_scheme(OPT_ARGS_NUM);
 double opt_geometry_stl_create_elementary(OPT_ARGS_NUM);
 double opt_geometry_stl_create_physical(OPT_ARGS_NUM);
+double opt_mesh_optimize(OPT_ARGS_NUM);
 double opt_mesh_quality(OPT_ARGS_NUM);
 double opt_mesh_normals(OPT_ARGS_NUM);
 double opt_mesh_tangents(OPT_ARGS_NUM);
diff --git a/Fltk/Callbacks.cpp b/Fltk/Callbacks.cpp
index 1cd05cd267a748f92da0190151549de6261b4330..71d72f12d2fb191600ba2bbd729b95a7a76aa8ee 100644
--- a/Fltk/Callbacks.cpp
+++ b/Fltk/Callbacks.cpp
@@ -1,4 +1,4 @@
-// $Id: Callbacks.cpp,v 1.251 2004-06-30 07:51:07 geuzaine Exp $
+// $Id: Callbacks.cpp,v 1.252 2004-06-30 17:49:51 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -972,6 +972,7 @@ void mesh_options_color_scheme_cb(CALLBACK_ARGS)
 
 void mesh_options_ok_cb(CALLBACK_ARGS)
 {
+  opt_mesh_optimize(0, GMSH_SET, WID->mesh_butt[2]->value());
   opt_mesh_order(0, GMSH_SET, WID->mesh_butt[3]->value()? 2 : 1);
   opt_mesh_interactive(0, GMSH_SET, WID->mesh_butt[4]->value());
   opt_mesh_constrained_bgmesh(0, GMSH_SET, WID->mesh_butt[5]->value());
@@ -2490,13 +2491,14 @@ void mesh_degree_cb(CALLBACK_ARGS)
 
 void mesh_optimize_cb(CALLBACK_ARGS)
 {
-  List_T *list = Tree2List(THEM->Volumes);
-  for(int i = 0; i < List_Nbr(list); i++){
-    Volume *v;
-    List_Read(list, i, &v);
-    Optimize_Netgen(v);
+  if(CTX.threads_lock) {
+    Msg(INFO, "I'm busy! Ask me that later...");
+    return;
   }
-  List_Delete(list);
+  CTX.threads_lock = 1;
+  Optimize_Netgen(THEM);
+  CTX.threads_lock = 0;
+
   CTX.mesh.changed = 1;
   Draw();
   Msg(STATUS3N, "Ready");
diff --git a/Fltk/GUI.cpp b/Fltk/GUI.cpp
index 1e56f3503a7ac2aa82e00544556555beeb98324c..c9afd06d439f31005fa5cab86db56ddfba80b809 100644
--- a/Fltk/GUI.cpp
+++ b/Fltk/GUI.cpp
@@ -1,4 +1,4 @@
-// $Id: GUI.cpp,v 1.320 2004-06-30 16:38:58 geuzaine Exp $
+// $Id: GUI.cpp,v 1.321 2004-06-30 17:49:51 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -285,11 +285,11 @@ Context_Item menu_mesh[] = {
   { "1D",     (Fl_Callback *)mesh_1d_cb } ,
   { "2D",     (Fl_Callback *)mesh_2d_cb } , 
   { "3D",     (Fl_Callback *)mesh_3d_cb } , 
-#if defined(HAVE_NETGEN)
-  { "Optimize 3D", (Fl_Callback *)mesh_optimize_cb } , 
-#endif
   { "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 } ,
   { NULL } 
 };  
@@ -1788,12 +1788,20 @@ void GUI::create_option_window()
         mesh_value[i]->align(FL_ALIGN_RIGHT);
       }
 
-      mesh_butt[3] = new Fl_Check_Button(2 * WB, 2 * WB + 8 * BH, BW, BH, "Generate second order elements");
+      mesh_butt[2] = new Fl_Check_Button(2 * WB, 2 * WB + 8 * BH, BW, BH, "Optimize quality of tetrahedral elements");
+      mesh_butt[2]->type(FL_TOGGLE_BUTTON);
+      mesh_butt[2]->down_box(TOGGLE_BOX);
+      mesh_butt[2]->selection_color(TOGGLE_COLOR);
+#if !defined(HAVE_NETGEN)
+      mesh_butt[2]->deactivate();
+#endif
+
+      mesh_butt[3] = new Fl_Check_Button(2 * WB, 2 * WB + 9 * BH, BW, BH, "Generate second order elements");
       mesh_butt[3]->type(FL_TOGGLE_BUTTON);
       mesh_butt[3]->down_box(TOGGLE_BOX);
       mesh_butt[3]->selection_color(TOGGLE_COLOR);
 
-      mesh_butt[5] = new Fl_Check_Button(2 * WB, 2 * WB + 9 * BH, BW, BH, "Constrain background mesh");
+      mesh_butt[5] = new Fl_Check_Button(2 * WB, 2 * WB + 10 * BH, BW, BH, "Constrain background mesh");
       mesh_butt[5]->type(FL_TOGGLE_BUTTON);
       mesh_butt[5]->down_box(TOGGLE_BOX);
       mesh_butt[5]->selection_color(TOGGLE_COLOR);
diff --git a/Mesh/3D_Mesh_Netgen.cpp b/Mesh/3D_Mesh_Netgen.cpp
index 149e121d4c289bb7088832d43f6887e0cb79b524..b7d07b09e88781c28732ce027fe5104826d287da 100644
--- a/Mesh/3D_Mesh_Netgen.cpp
+++ b/Mesh/3D_Mesh_Netgen.cpp
@@ -1,4 +1,4 @@
-// $Id: 3D_Mesh_Netgen.cpp,v 1.8 2004-06-30 16:38:58 geuzaine Exp $
+// $Id: 3D_Mesh_Netgen.cpp,v 1.9 2004-06-30 17:49:51 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -46,6 +46,11 @@ void Optimize_Netgen(Volume * v)
   Msg(GERROR, "Netgen is not compiled in this version of Gmsh");
 }
 
+void Optimize_Netgen(Mesh * m)
+{
+  Msg(GERROR, "Netgen is not compiled in this version of Gmsh");
+}
+
 #else
 
 #include "nglib.h"
@@ -53,12 +58,12 @@ void Optimize_Netgen(Volume * v)
 
 class Netgen{
  private:
-  List_T *_vertices, *_volverts;
+  List_T *_surverts, *_volverts;
   Volume *_vol;
   Ng_Mesh *_ngmesh;
  public:
   Netgen(Volume *vol, int importVolumeMesh = 0);
-  Netgen(Surface *s, int importSurfaceMesh = 0);
+  Netgen(Surface *sur, int importSurfaceMesh = 0);
   ~Netgen();
   void MeshVolume();
   void TransferVolumeMesh();
@@ -81,8 +86,8 @@ Netgen::Netgen(Volume *vol, int importVolumeMesh)
     List_Read(_vol->Surfaces, i, &s);
     Tree_Unit(tree, s->Vertices);
   }
-  _vertices = Tree2List(tree);
-  List_Sort(_vertices, compareVertex);
+  _surverts = Tree2List(tree);
+  List_Sort(_surverts, compareVertex);
 
   if(importVolumeMesh){
     Tree_T *tree2 = Tree_Soustraction(_vol->Vertices, tree);
@@ -94,9 +99,9 @@ Netgen::Netgen(Volume *vol, int importVolumeMesh)
   Tree_Delete(tree);
   
   // Transfer the vertices
-  for(int i = 0; i < List_Nbr(_vertices); i++){
+  for(int i = 0; i < List_Nbr(_surverts); i++){
     Vertex *v;
-    List_Read(_vertices, i, &v);
+    List_Read(_surverts, i, &v);
     double tmp[3];
     tmp[0] = v->Pos.X;
     tmp[1] = v->Pos.Y;
@@ -134,9 +139,9 @@ Netgen::Netgen(Volume *vol, int importVolumeMesh)
 	index[1] = 2;
 	index[2] = 1;
       }
-      tmp[0] = 1 + List_ISearch(_vertices, &simp->V[index[0]], compareVertex);
-      tmp[1] = 1 + List_ISearch(_vertices, &simp->V[index[1]], compareVertex);
-      tmp[2] = 1 + List_ISearch(_vertices, &simp->V[index[2]], compareVertex);
+      tmp[0] = 1 + List_ISearch(_surverts, &simp->V[index[0]], compareVertex);
+      tmp[1] = 1 + List_ISearch(_surverts, &simp->V[index[1]], compareVertex);
+      tmp[2] = 1 + List_ISearch(_surverts, &simp->V[index[2]], compareVertex);
       Ng_AddSurfaceElement(_ngmesh, NG_TRIG, tmp);
     }
     List_Delete(simplist);
@@ -154,11 +159,11 @@ Netgen::Netgen(Volume *vol, int importVolumeMesh)
 	simp->V[1] = temp;
       }
       int tmp[4];
-      tmp[0] = 1 + List_ISearch(_vertices, &simp->V[0], compareVertex);
-      tmp[1] = 1 + List_ISearch(_vertices, &simp->V[1], compareVertex);
-      tmp[2] = 1 + List_ISearch(_vertices, &simp->V[2], compareVertex);
-      tmp[3] = 1 + List_ISearch(_vertices, &simp->V[3], compareVertex);
-      int n = List_Nbr(_vertices) + 1;
+      tmp[0] = 1 + List_ISearch(_surverts, &simp->V[0], compareVertex);
+      tmp[1] = 1 + List_ISearch(_surverts, &simp->V[1], compareVertex);
+      tmp[2] = 1 + List_ISearch(_surverts, &simp->V[2], compareVertex);
+      tmp[3] = 1 + List_ISearch(_surverts, &simp->V[3], compareVertex);
+      int n = List_Nbr(_surverts) + 1;
       if(!tmp[0]) tmp[0] = n + List_ISearch(_volverts, &simp->V[0], compareVertex);
       if(!tmp[1]) tmp[1] = n + List_ISearch(_volverts, &simp->V[1], compareVertex);
       if(!tmp[2]) tmp[2] = n + List_ISearch(_volverts, &simp->V[2], compareVertex);
@@ -176,7 +181,7 @@ Netgen::Netgen(Surface *sur, int importSurfaceMesh)
 
 Netgen::~Netgen()
 {
-  List_Delete(_vertices);
+  List_Delete(_surverts);
   List_Delete(_volverts);
   Ng_DeleteMesh(_ngmesh);
   Ng_Exit();
@@ -201,14 +206,14 @@ void Netgen::TransferVolumeMesh()
   Vertex **vtable = (Vertex **)Malloc(nbv * sizeof(Vertex*));
   
   // Get existing unmodified surface vertices
-  for(int i = 0; i < List_Nbr(_vertices); i++){
-    List_Read(_vertices, i, &vtable[i]);
+  for(int i = 0; i < List_Nbr(_surverts); i++){
+    List_Read(_surverts, i, &vtable[i]);
     Tree_Insert(_vol->Vertices, &vtable[i]);
     Tree_Insert(THEM->Vertices, &vtable[i]);
   }
 
   // Create new volume vertices
-  for(int i = List_Nbr(_vertices); i < nbv; i++) {
+  for(int i = List_Nbr(_surverts); i < nbv; i++) {
     double tmp[3];
     Ng_GetPoint(_ngmesh, i+1, tmp);
     vtable[i] = Create_Vertex(++(THEM->MaxPointNum), tmp[0], tmp[1], tmp[2], 1., 0);
@@ -292,4 +297,21 @@ void Optimize_Netgen(Volume * v)
   ng.TransferVolumeMesh();
 }
 
+void Optimize_Netgen(Mesh *m)
+{
+  Msg(STATUS2, "Optimize volume mesh...");
+  double t1 = Cpu();
+
+  List_T *list = Tree2List(m->Volumes);
+  for(int i = 0; i < List_Nbr(list); i++){
+    Volume *v;
+    List_Read(list, i, &v);
+    Optimize_Netgen(v);
+  }
+  List_Delete(list);
+
+  double t2 = Cpu();
+  Msg(STATUS2, "Optimize volume mesh complete (%g s)", t2 - t1);
+}
+
 #endif // !HAVE_NETGEN
diff --git a/Mesh/Generator.cpp b/Mesh/Generator.cpp
index d807fa068a628811a0a8a34f49ea5d90da13f52b..f53534ef553d8448bfc728ac765f45c93e374f2d 100644
--- a/Mesh/Generator.cpp
+++ b/Mesh/Generator.cpp
@@ -1,4 +1,4 @@
-// $Id: Generator.cpp,v 1.58 2004-06-30 07:27:20 geuzaine Exp $
+// $Id: Generator.cpp,v 1.59 2004-06-30 17:49:51 geuzaine Exp $
 //
 // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
 //
@@ -295,14 +295,6 @@ void Maillage_Dimension_3(Mesh * M)
   Tree_Suppress(M->Volumes, &v);
   Free_Volume_But_Not_Elements(&v, NULL);
 
-  // optimize quality
-  if(CTX.mesh.quality) {
-    for(int i = 0; i < List_Nbr(list); i++){
-      List_Read(list, i, &v);
-      Optimize_Netgen(v);
-    }
-  }
-
   List_Delete(list);
 
   t2 = Cpu();
@@ -465,6 +457,11 @@ void mai3d(Mesh * M, int Asked)
     M->status = 3;
   }
 
+  // Optimize quality
+
+  if(M->status == 3 && CTX.mesh.optimize)
+    Optimize_Netgen(M);
+
   // Second order elements
 
   if(M->status && CTX.mesh.order == 2)
diff --git a/Mesh/Mesh.h b/Mesh/Mesh.h
index da8a0fbd1d224fd334e99774c5bea5f237d987a5..503126b30a172f033a55f846f651c9a1c94adcec 100644
--- a/Mesh/Mesh.h
+++ b/Mesh/Mesh.h
@@ -470,6 +470,7 @@ void Maillage_Automatique_VieuxCode(Surface * pS, Mesh * m, int ori);
 int Mesh_Shewchuk(Surface *s);
 int Mesh_Netgen(Volume * v);
 void Optimize_Netgen(Volume * v);
+void Optimize_Netgen(Mesh * m);
 
 int Calcule_Contours(Surface * s);
 void Link_Simplexes(List_T * Sim, Tree_T * Tim);
diff --git a/doc/texinfo/command_line.texi b/doc/texinfo/command_line.texi
index 9eefdd874db6ebfd706d582fe0d21aba7cde22c0..4c075fc9f38622d5a8ebde55d9c1ed4d7d064e06 100644
--- a/doc/texinfo/command_line.texi
+++ b/doc/texinfo/command_line.texi
@@ -22,6 +22,8 @@ set output mesh format (default: msh)
 select 2D mesh algorithm (default: iso)
 @item -smooth int
 set mesh smoothing (default: 0)
+@item -optimize
+optimize quality of tetrahedral elements
 @item -order int
 set the order of the generated elements (default: 1)
 @item -scale float