diff --git a/Common/Context.cpp b/Common/Context.cpp
index 3b36d9fc547071784eacded10dc20eb494d1c922..b7790e68c39d04e638ebae6c1cd864de2399fb5b 100644
--- a/Common/Context.cpp
+++ b/Common/Context.cpp
@@ -1,4 +1,4 @@
-// $Id: Context.cpp,v 1.63 2008-07-05 23:00:57 geuzaine Exp $
+// $Id: Context.cpp,v 1.64 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -30,9 +30,9 @@ void Context_T::buildRotationMatrix(void)
     setEulerAnglesFromRotationMatrix();
   }
   else {
-    double x = r[0] * Pi / 180.;
-    double y = r[1] * Pi / 180.;
-    double z = r[2] * Pi / 180.;
+    double x = r[0] * M_PI / 180.;
+    double y = r[1] * M_PI / 180.;
+    double z = r[2] * M_PI / 180.;
     double A = cos(x);
     double B = sin(x);
     double C = cos(y);
@@ -58,7 +58,7 @@ void Context_T::addQuaternion(double p1x, double p1y, double p2x, double p2y)
 
 void Context_T::addQuaternionFromAxisAndAngle(double axis[3], double angle)
 {
-  double a = angle * Pi / 180.;
+  double a = angle * M_PI / 180.;
   double quat[4];
   axis_to_quat(axis, a, quat);
   add_quats(quat, quaternion, quaternion);  
@@ -74,9 +74,9 @@ void Context_T::setQuaternion(double q0, double q1, double q2, double q3)
 
 void Context_T::setQuaternionFromEulerAngles()
 {
-  double x = r[0] * Pi / 180.;
-  double y = r[1] * Pi / 180.;
-  double z = r[2] * Pi / 180.;
+  double x = r[0] * M_PI / 180.;
+  double y = r[1] * M_PI / 180.;
+  double z = r[2] * M_PI / 180.;
   double xx[3] = {1.,0.,0.};
   double yy[3] = {0.,1.,0.};
   double zz[3] = {0.,0.,1.};
@@ -92,20 +92,20 @@ void Context_T::setEulerAnglesFromRotationMatrix()
 {
   r[1] = asin(rot[8]); // Calculate Y-axis angle
   double C =  cos(r[1]);
-  r[1] *=  180. / Pi;
+  r[1] *=  180. / M_PI;
   if(fabs(C) > 0.005){ // Gimball lock?
     double tmpx =  rot[10] / C; // No, so get X-axis angle
     double tmpy = -rot[9] / C;
-    r[0] = atan2(tmpy, tmpx) * 180. / Pi;
+    r[0] = atan2(tmpy, tmpx) * 180. / M_PI;
     tmpx =  rot[0] / C; // Get Z-axis angle
     tmpy = -rot[4] / C;
-    r[2] = atan2(tmpy, tmpx) * 180. / Pi;
+    r[2] = atan2(tmpy, tmpx) * 180. / M_PI;
   }
   else{ // Gimball lock has occurred
     r[0] = 0.; // Set X-axis angle to zero
     double tmpx = rot[5]; // And calculate Z-axis angle
     double tmpy = rot[1];
-    r[2] = atan2(tmpy, tmpx) * 180. / Pi;
+    r[2] = atan2(tmpy, tmpx) * 180. / M_PI;
   }
   // return only positive angles in [0,360]
   if(r[0] < 0.) r[0] += 360.;
diff --git a/Common/OpenFile.cpp b/Common/OpenFile.cpp
index ea092bb6524e3351c3ab41db9b4958282e78f277..0ef419ad4c8f49986c062417e9aa597dd42b2e24 100644
--- a/Common/OpenFile.cpp
+++ b/Common/OpenFile.cpp
@@ -1,4 +1,4 @@
-// $Id: OpenFile.cpp,v 1.2 2008-07-04 16:58:48 geuzaine Exp $
+// $Id: OpenFile.cpp,v 1.3 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -50,6 +50,8 @@ extern GUI *WID;
 
 extern Context_T CTX;
 
+#define SQU(a)      ((a)*(a))
+
 static void FinishUpBoundingBox()
 {
   double range[3];
@@ -83,18 +85,18 @@ static void FinishUpBoundingBox()
     CTX.min[1] -= CTX.lc; CTX.max[1] += CTX.lc;
   }
   else if(range[0] < CTX.geom.tolerance) {
-    CTX.lc = sqrt(DSQR(range[1]) + DSQR(range[2]));
+    CTX.lc = sqrt(SQU(range[1]) + SQU(range[2]));
     CTX.min[0] -= CTX.lc; CTX.max[0] += CTX.lc;
   }
   else if(range[1] < CTX.geom.tolerance) {
-    CTX.lc = sqrt(DSQR(range[0]) + DSQR(range[2]));
+    CTX.lc = sqrt(SQU(range[0]) + SQU(range[2]));
     CTX.min[1] -= CTX.lc; CTX.max[1] += CTX.lc;
   }
   else if(range[2] < CTX.geom.tolerance) {
-    CTX.lc = sqrt(DSQR(range[0]) + DSQR(range[1]));
+    CTX.lc = sqrt(SQU(range[0]) + SQU(range[1]));
   }
   else {
-    CTX.lc = sqrt(DSQR(range[0]) + DSQR(range[1]) + DSQR(range[2]));
+    CTX.lc = sqrt(SQU(range[0]) + SQU(range[1]) + SQU(range[2]));
   }
 }
 
diff --git a/Common/SmoothData.cpp b/Common/SmoothData.cpp
index eca379b201248d7f57dae62154fdf738d1b2a4bc..29878378894793a5ea6addf33b7ce58fbfa9e7cc 100644
--- a/Common/SmoothData.cpp
+++ b/Common/SmoothData.cpp
@@ -1,4 +1,4 @@
-// $Id: SmoothData.cpp,v 1.8 2008-07-03 17:06:01 geuzaine Exp $
+// $Id: SmoothData.cpp,v 1.9 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -157,7 +157,7 @@ float xyzn::angle(int i, char nx, char ny, char nz)
   prosca(a, b, &cosc);
   double sinc = sqrt(c[0] * c[0] + c[1] * c[1] + c[2] * c[2]);
   double angplan = myatan2(sinc, cosc);
-  return (float)(angplan * 180. / Pi);
+  return (float)(angplan * 180. / M_PI);
 }
 
 void xyzn::update(char nx, char ny, char nz, float tol)
diff --git a/Fltk/GUI.cpp b/Fltk/GUI.cpp
index d737812a9fc6712acb594a1e60ab7f781c20c61f..c5261f78feeb9312b931b0bcb96af2d8417d57aa 100644
--- a/Fltk/GUI.cpp
+++ b/Fltk/GUI.cpp
@@ -1,4 +1,4 @@
-// $Id: GUI.cpp,v 1.694 2008-07-05 23:00:57 geuzaine Exp $
+// $Id: GUI.cpp,v 1.695 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -3416,7 +3416,7 @@ void GUI::update_view_window(int num)
   PViewData *data = view->getData();
   PViewOptions *opt = view->getOptions();
 
-  double maxval = MAX(fabs(data->getMin()), fabs(data->getMax()));
+  double maxval = std::max(fabs(data->getMin()), fabs(data->getMax()));
   if(!maxval) maxval = 1.;
   double val1 = 10. * CTX.lc;
   double val2 = 2. * CTX.lc / maxval;
diff --git a/Fltk/Opengl_Window.cpp b/Fltk/Opengl_Window.cpp
index c5a1cddcb5567edaecffab4d2f9f8e910490e9f3..8eaedc1d80ffad7d9fd33a1acac37cfab54cd620 100644
--- a/Fltk/Opengl_Window.cpp
+++ b/Fltk/Opengl_Window.cpp
@@ -1,4 +1,4 @@
-// $Id: Opengl_Window.cpp,v 1.85 2008-05-06 21:11:47 geuzaine Exp $
+// $Id: Opengl_Window.cpp,v 1.86 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -68,7 +68,7 @@ void lasso_zoom(MousePosition &click1, MousePosition &click2)
 
   CTX.s[0] *= (double)CTX.viewport[2] / (click2.win[0] - click1.win[0]);
   CTX.s[1] *= (double)CTX.viewport[3] / (click2.win[1] - click1.win[1]);
-  CTX.s[2] = MIN(CTX.s[0], CTX.s[1]); // bof...
+  CTX.s[2] = std::min(CTX.s[0], CTX.s[1]); // bof...
   
   // recenter around the center of the lasso rectangle
   MousePosition tmp(click1);
diff --git a/Geo/GFace.cpp b/Geo/GFace.cpp
index a450c243c12badfa182135e28c3441631ef5a9e8..6cf735bea49dcbc8790d1afb0c6f79dd8436037e 100644
--- a/Geo/GFace.cpp
+++ b/Geo/GFace.cpp
@@ -1,4 +1,4 @@
-// $Id: GFace.cpp,v 1.68 2008-07-08 12:44:33 remacle Exp $
+// $Id: GFace.cpp,v 1.69 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -44,6 +44,8 @@ void dsvdcmp(double **a, int m, int n, double w[], double **v);
 
 extern Context_T CTX;
 
+#define SQU(a)      ((a)*(a))
+
 GFace::GFace(GModel *model, int tag)
   : GEntity(model, tag), r1(0), r2(0), va_geom_triangles(0)
 {
@@ -333,7 +335,7 @@ void GFace::computeMeanPlane(const std::vector<SPoint3> &points)
     prosca(res, res2, &cosc);
     sinc = sqrt(c[0] * c[0] + c[1] * c[1] + c[2] * c[2]);
     angplan = myatan2(sinc, cosc);
-    angplan = angle_02pi(angplan) * 180. / Pi;
+    angplan = angle_02pi(angplan) * 180. / M_PI;
     if((angplan > 70 && angplan < 110) || (angplan > 260 && angplan < 280)) {
       Msg::Info("SVD failed (angle=%g): using rough algo...", angplan);
       res[0] = res2[0];
@@ -512,7 +514,7 @@ void GFace::XYZtoUV(const double X, const double Y, const double Z,
       iter = 1;
 
       GPoint P = point(U, V);
-      err2 = sqrt(DSQR(X - P.x()) + DSQR(Y - P.y()) + DSQR(Z - P.z()));
+      err2 = sqrt(SQU(X - P.x()) + SQU(Y - P.y()) + SQU(Z - P.z()));
       if (err2 < 1.e-8 * CTX.lc) return;
 
       while(err > Precision && iter < MaxIter) {
@@ -538,8 +540,8 @@ void GFace::XYZtoUV(const double X, const double Y, const double Z,
 
         //if(Unew > umax || Vnew > vmax || Unew < umin || Vnew < vmin) break;
 
-        err = DSQR(Unew - U) + DSQR(Vnew - V);
-        err2 = sqrt(DSQR(X - P.x()) + DSQR(Y - P.y()) + DSQR(Z - P.z()));
+        err = SQU(Unew - U) + SQU(Vnew - V);
+        err2 = sqrt(SQU(X - P.x()) + SQU(Y - P.y()) + SQU(Z - P.z()));
         iter++;
         U = Unew;
         V = Vnew;
diff --git a/Geo/Geo.cpp b/Geo/Geo.cpp
index 40a16d32300de0a0d623a20dbc72f605f225d563..f670e52b979547dcac2640979600915409f86ffb 100644
--- a/Geo/Geo.cpp
+++ b/Geo/Geo.cpp
@@ -1,4 +1,4 @@
-// $Id: Geo.cpp,v 1.121 2008-07-04 14:58:31 geuzaine Exp $
+// $Id: Geo.cpp,v 1.122 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -29,6 +29,8 @@
 #include "Field.h"
 #include "Context.h"
 
+#define SQU(a)      ((a)*(a))
+
 extern Context_T CTX;
 
 static List_T *ListOfTransformedPoints = NULL;
@@ -123,7 +125,7 @@ Vertex *Create_Vertex(int Num, double X, double Y, double Z, double lc, double u
   pV->w = 1.0;
   pV->Num = Num;
   GModel::current()->getGEOInternals()->MaxPointNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxPointNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxPointNum, Num);
   pV->u = u;
   pV->geometry = 0;
   return pV;
@@ -136,7 +138,7 @@ Vertex *Create_Vertex(int Num, double u, double v, gmshSurface *surf, double lc)
   pV->w = 1.0;
   pV->Num = Num;
   GModel::current()->getGEOInternals()->MaxPointNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxPointNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxPointNum, Num);
   pV->u = u;
   pV->geometry = surf;
   pV->pntOnGeometry = SPoint2(u,v);
@@ -158,7 +160,7 @@ PhysicalGroup *Create_PhysicalGroup(int Num, int typ, List_T *intlist)
   p->Entities = List_Create(List_Nbr(intlist), 1, sizeof(int));
   p->Num = Num;
   GModel::current()->getGEOInternals()->MaxPhysicalNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxPhysicalNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxPhysicalNum, Num);
   p->Typ = typ;
   p->Visible = 1;
   for(int i = 0; i < List_Nbr(intlist); i++) {
@@ -185,7 +187,7 @@ EdgeLoop *Create_EdgeLoop(int Num, List_T *intlist)
   l->Curves = List_Create(List_Nbr(intlist), 1, sizeof(int));
   l->Num = Num;
   GModel::current()->getGEOInternals()->MaxLineLoopNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxLineLoopNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxLineLoopNum, Num);
   for(int i = 0; i < List_Nbr(intlist); i++) {
     int j;
     List_Read(intlist, i, &j);
@@ -210,7 +212,7 @@ SurfaceLoop *Create_SurfaceLoop(int Num, List_T *intlist)
   l->Surfaces = List_Create(List_Nbr(intlist), 1, sizeof(int));
   l->Num = Num;
   GModel::current()->getGEOInternals()->MaxSurfaceLoopNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxSurfaceLoopNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxSurfaceLoopNum, Num);
   for(int i = 0; i < List_Nbr(intlist); i++) {
     int j;
     List_Read(intlist, i, &j);
@@ -397,11 +399,11 @@ void End_Curve(Curve *c)
         // sur y1/f2 ou y3/f2, qui peuvent legerement etre hors de
         // [-1,1]
         if(x1 < 0)
-          A1 = -myasin(y1 / f2) + A4 + Pi;
+          A1 = -myasin(y1 / f2) + A4 + M_PI;
         else
           A1 = myasin(y1 / f2) + A4;
         if(x3 < 0)
-          A3 = -myasin(y3 / f2) + A4 + Pi;
+          A3 = -myasin(y3 / f2) + A4 + M_PI;
         else
           A3 = myasin(y3 / f2) + A4;
       }
@@ -416,7 +418,7 @@ void End_Curve(Curve *c)
     A1 = angle_02pi(A1);
     A3 = angle_02pi(A3);
     if(A1 >= A3)
-      A3 += 2 * Pi;
+      A3 += 2 * M_PI;
 
     c->Circle.t1 = A1;
     c->Circle.t2 = A3;
@@ -427,7 +429,7 @@ void End_Curve(Curve *c)
     for(int i = 0; i < 4; i++)
       c->Circle.v[i] = v[i];
 
-    if(!CTX.expert_mode && c->Num > 0 && A3 - A1 > 1.01 * Pi){
+    if(!CTX.expert_mode && c->Num > 0 && A3 - A1 > 1.01 * M_PI){
       Msg::Error("Circle or ellipse arc %d greater than Pi (angle=%g)", c->Num, A3-A1);
       Msg::Error("(If you understand what this implies, you can disable this error");
       Msg::Error("message by selecting `Enable expert mode' in the option dialog.");
@@ -479,7 +481,7 @@ Curve *Create_Curve(int Num, int Typ, int Order, List_T *Liste,
   pC->Typ = Typ;
   pC->Num = Num;
   GModel::current()->getGEOInternals()->MaxLineNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxLineNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxLineNum, Num);
   pC->Method = MESH_UNSTRUCTURED;
   pC->degre = Order;
   pC->Circle.n[0] = 0.0;
@@ -593,7 +595,7 @@ Surface *Create_Surface(int Num, int Typ)
   pS->RuledSurfaceOptions = 0;
 
   GModel::current()->getGEOInternals()->MaxSurfaceNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxSurfaceNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxSurfaceNum, Num);
   pS->Typ = Typ;
   pS->Method = MESH_UNSTRUCTURED;
   pS->Recombine = 0;
@@ -629,7 +631,7 @@ Volume *Create_Volume(int Num, int Typ)
   pV->Visible = 1;
   pV->Num = Num;
   GModel::current()->getGEOInternals()->MaxVolumeNum = 
-    IMAX(GModel::current()->getGEOInternals()->MaxVolumeNum, Num);
+    std::max(GModel::current()->getGEOInternals()->MaxVolumeNum, Num);
   pV->Typ = Typ;
   pV->Method = MESH_UNSTRUCTURED;
   pV->TrsfPoints = List_Create(6, 6, sizeof(Vertex *));
@@ -713,12 +715,12 @@ int NEWPHYSICAL(void)
 
 int NEWREG(void)
 {
-  return (IMAX(GModel::current()->getGEOInternals()->MaxLineNum,
-               IMAX(GModel::current()->getGEOInternals()->MaxLineLoopNum,
-                    IMAX(GModel::current()->getGEOInternals()->MaxSurfaceNum,
-                         IMAX(GModel::current()->getGEOInternals()->MaxSurfaceLoopNum,
-                              IMAX(GModel::current()->getGEOInternals()->MaxVolumeNum,
-                                   GModel::current()->getGEOInternals()->MaxPhysicalNum)))))
+  return (std::max(GModel::current()->getGEOInternals()->MaxLineNum,
+            std::max(GModel::current()->getGEOInternals()->MaxLineLoopNum,
+              std::max(GModel::current()->getGEOInternals()->MaxSurfaceNum,
+                std::max(GModel::current()->getGEOInternals()->MaxSurfaceLoopNum,
+                  std::max(GModel::current()->getGEOInternals()->MaxVolumeNum,
+                           GModel::current()->getGEOInternals()->MaxPhysicalNum)))))
           + 1);
 }
 
@@ -2615,21 +2617,21 @@ static void MaxNumPoint(void *a, void *b)
 {
   Vertex *v = *(Vertex **)a;
   GModel::current()->getGEOInternals()->MaxPointNum = 
-    MAX(GModel::current()->getGEOInternals()->MaxPointNum, v->Num);
+    std::max(GModel::current()->getGEOInternals()->MaxPointNum, v->Num);
 }
 
 static void MaxNumCurve(void *a, void *b)
 {
   Curve *c = *(Curve **)a;
   GModel::current()->getGEOInternals()->MaxLineNum = 
-    MAX(GModel::current()->getGEOInternals()->MaxLineNum, c->Num);
+    std::max(GModel::current()->getGEOInternals()->MaxLineNum, c->Num);
 }
 
 static void MaxNumSurface(void *a, void *b)
 {
   Surface *s = *(Surface **)a;
   GModel::current()->getGEOInternals()->MaxSurfaceNum =
-    MAX(GModel::current()->getGEOInternals()->MaxSurfaceNum, s->Num);
+    std::max(GModel::current()->getGEOInternals()->MaxSurfaceNum, s->Num);
 }
 
 static void ReplaceDuplicatePoints()
@@ -2915,9 +2917,9 @@ static double projectPC(double u)
   if(u < CURVE->ubeg)
     u = CURVE->ubeg;
   Vertex c = InterpolateCurve(CURVE, u, 0);
-  return sqrt(DSQR(c.Pos.X - VERTEX->Pos.X) +
-              DSQR(c.Pos.Y - VERTEX->Pos.Y) + 
-              DSQR(c.Pos.Z - VERTEX->Pos.Z));
+  return sqrt(SQU(c.Pos.X - VERTEX->Pos.X) +
+              SQU(c.Pos.Y - VERTEX->Pos.Y) + 
+              SQU(c.Pos.Z - VERTEX->Pos.Z));
 }
 
 bool ProjectPointOnCurve(Curve *c, Vertex *v, Vertex *RES, Vertex *DER)
diff --git a/Geo/GeoInterpolation.cpp b/Geo/GeoInterpolation.cpp
index 823ae6e93dfbec18fdff16bbf1445abb4317fbd0..a64de3ee4001f47415d89e272a7c448eac994e1a 100644
--- a/Geo/GeoInterpolation.cpp
+++ b/Geo/GeoInterpolation.cpp
@@ -1,4 +1,4 @@
-// $Id: GeoInterpolation.cpp,v 1.40 2008-07-03 17:06:01 geuzaine Exp $
+// $Id: GeoInterpolation.cpp,v 1.41 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -25,6 +25,8 @@
 #include "GeoStringInterface.h"
 #include "Numeric.h"
 
+#define SQU(a)      ((a)*(a))
+
 // Curves
 
 static Vertex InterpolateCubicSpline(Vertex *v[4], double t, double mat[4][4],
@@ -438,11 +440,11 @@ static Vertex TransfiniteTri(Vertex c1, Vertex c2, Vertex c3,
 
 static void TransfiniteSph(Vertex S, Vertex center, Vertex *T)
 {
-  double r = sqrt(DSQR(S.Pos.X - center.Pos.X) + DSQR(S.Pos.Y - center.Pos.Y)
-                  + DSQR(S.Pos.Z - center.Pos.Z));
+  double r = sqrt(SQU(S.Pos.X - center.Pos.X) + SQU(S.Pos.Y - center.Pos.Y)
+                  + SQU(S.Pos.Z - center.Pos.Z));
 
-  double s = sqrt(DSQR(T->Pos.X - center.Pos.X) + DSQR(T->Pos.Y - center.Pos.Y)
-                  + DSQR(T->Pos.Z - center.Pos.Z));
+  double s = sqrt(SQU(T->Pos.X - center.Pos.X) + SQU(T->Pos.Y - center.Pos.Y)
+                  + SQU(T->Pos.Z - center.Pos.Z));
 
   double dirx = (T->Pos.X - center.Pos.X) / s;
   double diry = (T->Pos.Y - center.Pos.Y) / s;
diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 8b7b4567821198ae0a76a06ae537a0d88ee85aed..6a373763a5986a3378ce0d65b47b3146d7afdf23 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -1,4 +1,4 @@
-// $Id: MElement.cpp,v 1.78 2008-07-08 12:44:33 remacle Exp $
+// $Id: MElement.cpp,v 1.79 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -36,6 +36,8 @@
 #  include "qualityMeasures.h"
 #endif
 
+#define SQU(a)      ((a)*(a))
+
 extern Context_T CTX;
 
 int MElement::_globalNum = 0;
@@ -260,9 +262,9 @@ double MElement::getJacobian(double u, double v, double w, double jac[3][3])
       prodve(a, b, c);
       jac[2][0] = c[0]; jac[2][1] = c[1]; jac[2][2] = c[2]; 
     }
-    return sqrt(SQR(jac[0][0] * jac[1][1] - jac[0][1] * jac[1][0]) +
-                SQR(jac[0][2] * jac[1][0] - jac[0][0] * jac[1][2]) +
-                SQR(jac[0][1] * jac[1][2] - jac[0][2] * jac[1][1]));
+    return sqrt(SQU(jac[0][0] * jac[1][1] - jac[0][1] * jac[1][0]) +
+                SQU(jac[0][2] * jac[1][0] - jac[0][0] * jac[1][2]) +
+                SQU(jac[0][1] * jac[1][2] - jac[0][2] * jac[1][1]));
   case 1:
     for(int i = 0; i < getNumVertices(); i++) {
       getGradShapeFunction(i, u, v, w, s);
@@ -285,7 +287,7 @@ double MElement::getJacobian(double u, double v, double w, double jac[3][3])
       jac[1][0] = b[0]; jac[1][1] = b[1]; jac[1][2] = b[2]; 
       jac[2][0] = c[0]; jac[2][1] = c[1]; jac[2][2] = c[2]; 
     }
-    return sqrt(SQR(jac[0][0]) + SQR(jac[0][1]) + SQR(jac[0][2]));
+    return sqrt(SQU(jac[0][0]) + SQU(jac[0][1]) + SQU(jac[0][2]));
   default:
     return 1.;
   }
@@ -320,7 +322,7 @@ void MElement::xyz2uvw(double xyz[3], double uvw[3])
       inv[0][1] * (xyz[0] - xn) + inv[1][1] * (xyz[1] - yn) + inv[2][1] * (xyz[2] - zn) ;
     double wn = uvw[2] +
       inv[0][2] * (xyz[0] - xn) + inv[1][2] * (xyz[1] - yn) + inv[2][2] * (xyz[2] - zn) ;
-    error = sqrt(SQR(un - uvw[0]) + SQR(vn - uvw[1]) + SQR(wn - uvw[2]));
+    error = sqrt(SQU(un - uvw[0]) + SQU(vn - uvw[1]) + SQU(wn - uvw[2]));
     uvw[0] = un;
     uvw[1] = vn;
     uvw[2] = wn;
diff --git a/Geo/gmshEdge.cpp b/Geo/gmshEdge.cpp
index 15b4142c99498491b1027695095536012d4680a0..a512a7837f3b3c4307d2fcdec77e79e75d425007 100644
--- a/Geo/gmshEdge.cpp
+++ b/Geo/gmshEdge.cpp
@@ -1,4 +1,4 @@
-// $Id: gmshEdge.cpp,v 1.50 2008-07-03 17:06:02 geuzaine Exp $
+// $Id: gmshEdge.cpp,v 1.51 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -111,7 +111,7 @@ int gmshEdge::minimumMeshSegments () const
     np = GEdge::minimumMeshSegments();
   else if(geomType() == Circle || geomType() == Ellipse)
     np = (int)(fabs(c->Circle.t1 - c->Circle.t2) *
-                 (double)CTX.mesh.min_circ_points / Pi) - 1;
+                 (double)CTX.mesh.min_circ_points / M_PI) - 1;
   else
     np = CTX.mesh.min_curv_points - 1;
   return std::max(np, meshAttributes.minimumMeshSegments);
diff --git a/Graphics/Draw.cpp b/Graphics/Draw.cpp
index f8d47f52da383bb6f7ce9e37ea8830ca84622028..589c6febdc0069b5454c2be4ca3757109b0cce7f 100644
--- a/Graphics/Draw.cpp
+++ b/Graphics/Draw.cpp
@@ -1,4 +1,4 @@
-// $Id: Draw.cpp,v 1.119 2008-05-04 08:31:14 geuzaine Exp $
+// $Id: Draw.cpp,v 1.120 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -168,7 +168,7 @@ void InitProjection(int xpick, int ypick, int wpick, int hpick)
                   (GLdouble)wpick, (GLdouble)hpick, (GLint *)CTX.viewport);
 
   double grad_z, grad_xy;
-  double zmax = MAX(fabs(CTX.min[2]), fabs(CTX.max[2]));
+  double zmax = std::max(fabs(CTX.min[2]), fabs(CTX.max[2]));
   if(zmax < CTX.lc) zmax = CTX.lc;
 
   if(CTX.ortho) {
@@ -233,7 +233,7 @@ void InitProjection(int xpick, int ypick, int wpick, int hpick)
     else{ // radial
       double cx = grad_xy * (CTX.vxmin + CTX.vxmax) / 2.;
       double cy = grad_xy * (CTX.vymin + CTX.vymax) / 2.;
-      double r = grad_xy * MAX(CTX.vxmax - CTX.vxmin, CTX.vymax - CTX.vymin) / 2.;
+      double r = grad_xy * std::max(CTX.vxmax - CTX.vxmin, CTX.vymax - CTX.vymin) / 2.;
       glBegin(GL_TRIANGLE_FAN);
       glColor4ubv((GLubyte *) & CTX.color.bg_grad);
       glVertex3d(cx, cy, 0.);
diff --git a/Makefile b/Makefile
index a5fbdbaef47308469fcef1ac45fcb137a4c299c3..201d3d9648217b4fba31e5027171fc90dd058598 100644
--- a/Makefile
+++ b/Makefile
@@ -1,4 +1,4 @@
-# $Id: Makefile,v 1.492 2008-07-04 22:45:54 geuzaine Exp $
+# $Id: Makefile,v 1.493 2008-07-10 13:29:23 geuzaine Exp $
 #
 # Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 #
@@ -136,7 +136,8 @@ doc-info:
 purge:
 	rm -f `find . -name "*~" -o -name "*~~" -o -name ".gmsh-errors"\
                -o -name "\#*" -o -name ".\#*" -o -name ".DS_Store"\
-               -o -name "gmon.out" -o -name ".gdb_history" -o -name "debug?.pos"`
+               -o -name "gmon.out" -o -name ".gdb_history" -o -name "debug?.pos"\
+               -o -name "*.bak"`
 
 clean:
 	for i in doc lib ${GMSH_DIRS}; do (cd $$i && ${MAKE} clean); done
diff --git a/Mesh/DivideAndConquer.cpp b/Mesh/DivideAndConquer.cpp
index b8afbea596d877c3a31a831106533a41c30849bf..5329bbda03ba97c75beb92e13d76f6aac0257d08 100644
--- a/Mesh/DivideAndConquer.cpp
+++ b/Mesh/DivideAndConquer.cpp
@@ -1,4 +1,4 @@
-// $Id: DivideAndConquer.cpp,v 1.20 2008-06-27 18:00:52 geuzaine Exp $
+// $Id: DivideAndConquer.cpp,v 1.21 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -386,14 +386,14 @@ static int DListInsert(DListRecord **dlist, MPoint center, PointNumero newPoint)
   xx = (double)(pPointArray[newPoint].where.h - center.h);
   beta = atan2(yy, xx) - alpha1;
   if(beta <= 0)
-    beta += Deux_Pi;
+    beta += 2. * M_PI;
 
   do {
     yy = (double)(pPointArray[Succ(p)->point_num].where.v - center.v);
     xx = (double)(pPointArray[Succ(p)->point_num].where.h - center.h);
     alpha = atan2(yy, xx) - alpha1;
     if(alpha <= 1.e-15)
-      alpha += Deux_Pi;
+      alpha += 2. * M_PI;
     if(alpha >= beta) {
       Succ(newp) = Succ(p);
       Succ(p) = newp;
diff --git a/Mesh/HighOrder.cpp b/Mesh/HighOrder.cpp
index 1363055c2d4457625fa54156c7ea22b76084767f..5d30f79b006cc2db7b2a467821faa62d62b1c515 100644
--- a/Mesh/HighOrder.cpp
+++ b/Mesh/HighOrder.cpp
@@ -1,4 +1,4 @@
-// $Id: HighOrder.cpp,v 1.33 2008-07-08 12:44:34 remacle Exp $
+// $Id: HighOrder.cpp,v 1.34 2008-07-10 13:29:24 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -29,6 +29,8 @@
 #include "GmshMatrix.h"
 #include "FunctionSpace.h"
 
+#define SQU(a)      ((a)*(a))
+
 extern Context_T CTX;
 
 // for each pair of vertices (an edge), we build a list of vertices
@@ -915,7 +917,7 @@ static double mesh_functional_distorsion(MTriangle *t, double u, double v)
   double v2[3] = {mat[1][0], mat[1][1], mat[1][2]};
   double normal1[3];
   prodve(v1, v2, normal1);
-  double nn = sqrt(DSQR(normal1[0]) + DSQR(normal1[1]) + DSQR(normal1[2]));
+  double nn = sqrt(SQU(normal1[0]) + SQU(normal1[1]) + SQU(normal1[2]));
   
   // compute uncurved element jacobian d_u x and d_v x
   t->jac(u, v, 0,mat);
@@ -942,7 +944,7 @@ void getMinMaxJac (MTriangle *t, double &minJ, double &maxJ)
   double v2[3] = {mat[1][0], mat[1][1], mat[1][2]};
   double normal1[3], normal[3];
   prodve(v1, v2, normal1);
-  double nn = sqrt(DSQR(normal1[0]) + DSQR(normal1[1]) + DSQR(normal1[2]));
+  double nn = sqrt(SQU(normal1[0]) + SQU(normal1[1]) + SQU(normal1[2]));
   for(int i = 0; i < n; i++){
     for(int k = 0; k < n - i; k++){
       t->jac((double)i / (n - 1), (double)k / (n - 1), 0,mat);
diff --git a/Mesh/meshGEdge.cpp b/Mesh/meshGEdge.cpp
index 41ac4d298f75595da197a24b76a8213e916438b6..cb954c18acdf7047aaf51a15e60d38e38c698a6e 100644
--- a/Mesh/meshGEdge.cpp
+++ b/Mesh/meshGEdge.cpp
@@ -1,4 +1,4 @@
-// $Id: meshGEdge.cpp,v 1.68 2008-07-08 12:43:25 geuzaine Exp $
+// $Id: meshGEdge.cpp,v 1.69 2008-07-10 13:29:25 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -30,6 +30,8 @@
 #include "Context.h"
 #include "GModel.h"
 
+#define SQU(a)      ((a)*(a))
+
 extern Context_T CTX;
 
 typedef struct {
@@ -186,7 +188,7 @@ static double F_Transfinite(GEdge *ge, double t)
             / ((double)nbpt * ge->length());
         }
         double b = -a * ge->length() * ge->length() / (4. * (coef - 1.));
-        val = d / (-a * DSQR(t * ge->length() - (ge->length()) * 0.5) + b);
+        val = d / (-a * SQU(t * ge->length() - (ge->length()) * 0.5) + b);
       }
       break;
       
diff --git a/Mesh/meshGFaceTransfinite.cpp b/Mesh/meshGFaceTransfinite.cpp
index c726601441962aeb875043527ef0a544809683fc..b77e41cb869ea99afc847515a180bb654c4f7056 100644
--- a/Mesh/meshGFaceTransfinite.cpp
+++ b/Mesh/meshGFaceTransfinite.cpp
@@ -1,4 +1,4 @@
-// $Id: meshGFaceTransfinite.cpp,v 1.29 2008-07-01 14:24:07 geuzaine Exp $
+// $Id: meshGFaceTransfinite.cpp,v 1.30 2008-07-10 13:29:25 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -30,6 +30,8 @@
 #include "Message.h"
 #include "Numeric.h"
 
+#define SQU(a)      ((a)*(a))
+
 extern Context_T CTX;
 
 /*
@@ -318,12 +320,12 @@ int MeshTransfiniteSurface(GFace *gf)
           MVertex *v31 = tab[i + 1][j - 1];
           MVertex *v32 = tab[i + 1][j    ];
           MVertex *v33 = tab[i + 1][j + 1];
-          double alpha = 0.25 * (DSQR(v23->x() - v21->x()) +
-                                 DSQR(v23->y() - v21->y()) +
-                                 DSQR(v23->z() - v21->z()));
-          double gamma = 0.25 * (DSQR(v32->x() - v12->x()) +
-                                 DSQR(v32->y() - v12->y()) +
-                                 DSQR(v32->z() - v12->z()));
+          double alpha = 0.25 * (SQU(v23->x() - v21->x()) +
+                                 SQU(v23->y() - v21->y()) +
+                                 SQU(v23->z() - v21->z()));
+          double gamma = 0.25 * (SQU(v32->x() - v12->x()) +
+                                 SQU(v32->y() - v12->y()) +
+                                 SQU(v32->z() - v12->z()));
           double beta = 0.0625 * ((v32->x() - v12->x()) * (v23->x() - v21->x()) +
                                   (v32->y() - v12->y()) * (v23->y() - v21->y()) +
                                   (v32->z() - v12->z()) * (v23->z() - v21->z()));
diff --git a/Numeric/NumericEmbedded.cpp b/Numeric/NumericEmbedded.cpp
index d325b35bd8bc55285e687b742a6ebc6f5e0755a3..a98b74a40d41d421d335e61105205ff358bca489 100644
--- a/Numeric/NumericEmbedded.cpp
+++ b/Numeric/NumericEmbedded.cpp
@@ -1,4 +1,4 @@
-// $Id: NumericEmbedded.cpp,v 1.5 2008-05-04 08:31:16 geuzaine Exp $
+// $Id: NumericEmbedded.cpp,v 1.6 2008-07-10 13:29:25 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -31,6 +31,8 @@
 #include "Message.h"
 #endif
 
+#define SQU(a)      ((a)*(a))
+
 double myatan2(double a, double b)
 {
   if(a == 0.0 && b == 0)
@@ -41,9 +43,9 @@ double myatan2(double a, double b)
 double myasin(double a)
 {
   if(a <= -1.)
-    return -Pi / 2.;
+    return -M_PI / 2.;
   else if(a >= 1.)
-    return Pi / 2.;
+    return M_PI / 2.;
   else
     return asin(a);
 }
@@ -51,7 +53,7 @@ double myasin(double a)
 double myacos(double a)
 {
   if(a <= -1.)
-    return Pi;
+    return M_PI;
   else if(a >= 1.)
     return 0.;
   else
@@ -81,7 +83,7 @@ int sys2x2(double mat[2][2], double b[2], double res[2])
   double det, ud, norm;
   int i;
 
-  norm = DSQR(mat[0][0]) + DSQR(mat[1][1]) + DSQR(mat[0][1]) + DSQR(mat[1][0]);
+  norm = SQU(mat[0][0]) + SQU(mat[1][1]) + SQU(mat[0][1]) + SQU(mat[1][0]);
   det = mat[0][0] * mat[1][1] - mat[1][0] * mat[0][1];
 
   // TOLERANCE ! WARNING WARNING
@@ -162,9 +164,9 @@ int sys3x3_with_tol(double mat[3][3], double b[3], double res[3], double *det)
 
   out = sys3x3(mat, b, res, det);
   norm =
-    DSQR(mat[0][0]) + DSQR(mat[0][1]) + DSQR(mat[0][2]) +
-    DSQR(mat[1][0]) + DSQR(mat[1][1]) + DSQR(mat[1][2]) +
-    DSQR(mat[2][0]) + DSQR(mat[2][1]) + DSQR(mat[2][2]);
+    SQU(mat[0][0]) + SQU(mat[0][1]) + SQU(mat[0][2]) +
+    SQU(mat[1][0]) + SQU(mat[1][1]) + SQU(mat[1][2]) +
+    SQU(mat[2][0]) + SQU(mat[2][1]) + SQU(mat[2][2]);
 
   // TOLERANCE ! WARNING WARNING
   if(norm == 0.0 || fabs(*det) / norm < 1.e-12) {
@@ -238,7 +240,7 @@ double inv3x3(double mat[3][3], double inv[3][3])
 
 double angle_02pi(double A3)
 {
-  double DP = 2 * Pi;
+  double DP = 2 * M_PI;
   while(A3 > DP || A3 < 0.) {
     if(A3 > 0)
       A3 -= DP;
diff --git a/Numeric/NumericEmbedded.h b/Numeric/NumericEmbedded.h
index 219aa263031a94cd041a8ec78567ee12afe1f470..2cc8d1fc207921c36e479914007d870ae446dfa8 100644
--- a/Numeric/NumericEmbedded.h
+++ b/Numeric/NumericEmbedded.h
@@ -22,30 +22,6 @@
 
 #include <math.h>
 
-#define Pi      3.1415926535897932
-#define Deux_Pi 6.2831853071795865
-
-#if !defined(M_PI)
-#define M_PI Pi
-#endif
-
-#define MIN(a,b) (((a)<(b))?(a):(b))
-#define MAX(a,b) (((a)<(b))?(b):(a))
-#define SQR(a)   ((a)*(a))
-
-#define IMIN MIN
-#define LMIN MIN
-#define FMIN MIN
-#define DMIN MIN
-
-#define IMAX MAX
-#define LMAX MAX
-#define FMAX MAX
-#define DMAX MAX
-
-#define DSQR SQR
-#define FSQR SQR
-
 #define myhypot(a,b) (sqrt((a)*(a)+(b)*(b)))
 #define sign(x)      (((x)>=0)?1:-1)
 
diff --git a/Numeric/gsl_brent.cpp b/Numeric/gsl_brent.cpp
index 2863f6fd43bfa9b851502600a02631029a398e58..3bbaa4021340b8471856329bd19569fa093eb588 100644
--- a/Numeric/gsl_brent.cpp
+++ b/Numeric/gsl_brent.cpp
@@ -1,4 +1,4 @@
-// $Id: gsl_brent.cpp,v 1.20 2008-05-04 08:31:16 geuzaine Exp $
+// $Id: gsl_brent.cpp,v 1.21 2008-07-10 13:29:25 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -137,7 +137,7 @@ void mnbrak(double *ax, double *bx, double *cx,
     r = (*bx - *ax) * (f_b - f_c);
     q = (*bx - *cx) * (f_b - f_a);
     u = (*bx) - ((*bx - *cx) * q - (*bx - *ax) * r) / 
-      (2.0 * SIGN(MAX(fabs(q - r), MYTINY_), q - r));
+      (2.0 * SIGN(std::max(fabs(q - r), MYTINY_), q - r));
     ulim = *bx + MYLIMIT_ * (*cx - *bx);
 
     if((*bx - u) * (u - *cx) > 0.0) {
diff --git a/Parser/Gmsh.l b/Parser/Gmsh.l
index a29f0a464d11c98451cf857c1cfb0d171503c631..b22f72b865fd91b6311149af3bd2ed5d5cfd2b4d 100644
--- a/Parser/Gmsh.l
+++ b/Parser/Gmsh.l
@@ -1,5 +1,5 @@
 %{
-// $Id: Gmsh.l,v 1.103 2008-05-04 08:31:16 geuzaine Exp $
+// $Id: Gmsh.l,v 1.104 2008-07-10 13:29:25 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -26,7 +26,6 @@
 #include <math.h>
 
 #include "Message.h"
-#include "Numeric.h"
 #include "Geo.h"
 #include "Gmsh.tab.hpp"
 
@@ -326,7 +325,7 @@ void skip_until(const char *skip, const char *until)
       if(skip && chars[0] == skip[0]) break;
     }
 
-    l = MAX(l_skip,l_until);
+    l = std::max(l_skip,l_until);
     if(l >= (int)sizeof(chars)){
       Msg::Error("Search pattern too long in skip_until");
       return;
diff --git a/Parser/Gmsh.tab.cpp b/Parser/Gmsh.tab.cpp
index be783358b69de26c59619953380c8918de204281..1250610d816ffd7bda5f2e666a8c35391977793f 100644
--- a/Parser/Gmsh.tab.cpp
+++ b/Parser/Gmsh.tab.cpp
@@ -324,7 +324,7 @@
 /* Copy the first part of user declarations.  */
 #line 1 "Gmsh.y"
 
-// $Id: Gmsh.tab.cpp,v 1.374 2008-07-06 11:25:35 geuzaine Exp $
+// $Id: Gmsh.tab.cpp,v 1.375 2008-07-10 13:29:25 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
diff --git a/Parser/Gmsh.yy.cpp b/Parser/Gmsh.yy.cpp
index 28adfb3b5914880ad553643a93de23d231a1c3de..75353325c457ce1a59bcaf717a69fcc977cc926d 100644
--- a/Parser/Gmsh.yy.cpp
+++ b/Parser/Gmsh.yy.cpp
@@ -835,7 +835,7 @@ int gmsh_yy_flex_debug = 0;
 char *gmsh_yytext;
 #line 1 "Gmsh.l"
 #line 2 "Gmsh.l"
-// $Id: Gmsh.yy.cpp,v 1.373 2008-07-06 11:25:37 geuzaine Exp $
+// $Id: Gmsh.yy.cpp,v 1.374 2008-07-10 13:29:29 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -862,7 +862,6 @@ char *gmsh_yytext;
 #include <math.h>
 
 #include "Message.h"
-#include "Numeric.h"
 #include "Geo.h"
 #include "Gmsh.tab.hpp"
 
@@ -895,7 +894,7 @@ void   skipline(void);
 #define YY_NO_UNISTD_H
 #endif
 
-#line 899 "Gmsh.yy.cpp"
+#line 898 "Gmsh.yy.cpp"
 
 #define INITIAL 0
 
@@ -1048,10 +1047,10 @@ YY_DECL
 	register char *yy_cp, *yy_bp;
 	register int yy_act;
     
-#line 71 "Gmsh.l"
+#line 70 "Gmsh.l"
 
 
-#line 1055 "Gmsh.yy.cpp"
+#line 1054 "Gmsh.yy.cpp"
 
 	if ( !(yy_init) )
 		{
@@ -1137,721 +1136,721 @@ do_action:	/* This label is used only to access EOF actions. */
 case 1:
 /* rule 1 can match eol */
 YY_RULE_SETUP
-#line 73 "Gmsh.l"
+#line 72 "Gmsh.l"
 /* none */;
 	YY_BREAK
 case 2:
 YY_RULE_SETUP
-#line 74 "Gmsh.l"
+#line 73 "Gmsh.l"
 return tEND;
 	YY_BREAK
 case 3:
 YY_RULE_SETUP
-#line 75 "Gmsh.l"
+#line 74 "Gmsh.l"
 skipcomments();
 	YY_BREAK
 case 4:
 YY_RULE_SETUP
-#line 76 "Gmsh.l"
+#line 75 "Gmsh.l"
 skipline();
 	YY_BREAK
 case 5:
 YY_RULE_SETUP
-#line 77 "Gmsh.l"
+#line 76 "Gmsh.l"
 { parsestring('\"'); return tBIGSTR; }
 	YY_BREAK
 case 6:
 YY_RULE_SETUP
-#line 78 "Gmsh.l"
+#line 77 "Gmsh.l"
 { parsestring('\''); return tBIGSTR; }
 	YY_BREAK
 case 7:
 YY_RULE_SETUP
-#line 79 "Gmsh.l"
+#line 78 "Gmsh.l"
 { gmsh_yylval.d = NEWREG(); return tDOUBLE; }
 	YY_BREAK
 case 8:
 YY_RULE_SETUP
-#line 80 "Gmsh.l"
+#line 79 "Gmsh.l"
 { gmsh_yylval.d = NEWPOINT(); return tDOUBLE; }
 	YY_BREAK
 case 9:
 YY_RULE_SETUP
-#line 81 "Gmsh.l"
+#line 80 "Gmsh.l"
 { gmsh_yylval.d = NEWLINE(); return tDOUBLE; }
 	YY_BREAK
 case 10:
 YY_RULE_SETUP
-#line 82 "Gmsh.l"
+#line 81 "Gmsh.l"
 { gmsh_yylval.d = NEWLINE(); return tDOUBLE; }
 	YY_BREAK
 case 11:
 YY_RULE_SETUP
-#line 83 "Gmsh.l"
+#line 82 "Gmsh.l"
 { gmsh_yylval.d = NEWLINELOOP(); return tDOUBLE; }
 	YY_BREAK
 case 12:
 YY_RULE_SETUP
-#line 84 "Gmsh.l"
+#line 83 "Gmsh.l"
 { gmsh_yylval.d = NEWSURFACE(); return tDOUBLE; }
 	YY_BREAK
 case 13:
 YY_RULE_SETUP
-#line 85 "Gmsh.l"
+#line 84 "Gmsh.l"
 { gmsh_yylval.d = NEWSURFACELOOP(); return tDOUBLE; }
 	YY_BREAK
 case 14:
 YY_RULE_SETUP
-#line 86 "Gmsh.l"
+#line 85 "Gmsh.l"
 { gmsh_yylval.d = NEWVOLUME(); return tDOUBLE; }
 	YY_BREAK
 case 15:
 YY_RULE_SETUP
-#line 87 "Gmsh.l"
+#line 86 "Gmsh.l"
 { gmsh_yylval.d = NEWFIELD(); return tDOUBLE; }
 	YY_BREAK
 case 16:
 YY_RULE_SETUP
-#line 88 "Gmsh.l"
+#line 87 "Gmsh.l"
 return tAFFECT;
 	YY_BREAK
 case 17:
 YY_RULE_SETUP
-#line 89 "Gmsh.l"
+#line 88 "Gmsh.l"
 return tAFFECTPLUS;
 	YY_BREAK
 case 18:
 YY_RULE_SETUP
-#line 90 "Gmsh.l"
+#line 89 "Gmsh.l"
 return tAFFECTMINUS;
 	YY_BREAK
 case 19:
 YY_RULE_SETUP
-#line 91 "Gmsh.l"
+#line 90 "Gmsh.l"
 return tAFFECTTIMES;
 	YY_BREAK
 case 20:
 YY_RULE_SETUP
-#line 92 "Gmsh.l"
+#line 91 "Gmsh.l"
 return tAFFECTDIVIDE;
 	YY_BREAK
 case 21:
 YY_RULE_SETUP
-#line 93 "Gmsh.l"
+#line 92 "Gmsh.l"
 return tDOTS;
 	YY_BREAK
 case 22:
 YY_RULE_SETUP
-#line 94 "Gmsh.l"
+#line 93 "Gmsh.l"
 return tDOTS;
 	YY_BREAK
 case 23:
 YY_RULE_SETUP
-#line 95 "Gmsh.l"
+#line 94 "Gmsh.l"
 return tOR;
 	YY_BREAK
 case 24:
 YY_RULE_SETUP
-#line 96 "Gmsh.l"
+#line 95 "Gmsh.l"
 return tAND;
 	YY_BREAK
 case 25:
 YY_RULE_SETUP
-#line 97 "Gmsh.l"
+#line 96 "Gmsh.l"
 return tPLUSPLUS;
 	YY_BREAK
 case 26:
 YY_RULE_SETUP
-#line 98 "Gmsh.l"
+#line 97 "Gmsh.l"
 return tMINUSMINUS;
 	YY_BREAK
 case 27:
 YY_RULE_SETUP
-#line 99 "Gmsh.l"
+#line 98 "Gmsh.l"
 return tEQUAL;
 	YY_BREAK
 case 28:
 YY_RULE_SETUP
-#line 100 "Gmsh.l"
+#line 99 "Gmsh.l"
 return tNOTEQUAL;
 	YY_BREAK
 case 29:
 YY_RULE_SETUP
-#line 101 "Gmsh.l"
+#line 100 "Gmsh.l"
 return tLESSOREQUAL;
 	YY_BREAK
 case 30:
 YY_RULE_SETUP
-#line 102 "Gmsh.l"
+#line 101 "Gmsh.l"
 return tGREATEROREQUAL;
 	YY_BREAK
 case 31:
 YY_RULE_SETUP
-#line 104 "Gmsh.l"
+#line 103 "Gmsh.l"
 return tAcos;
 	YY_BREAK
 case 32:
 YY_RULE_SETUP
-#line 105 "Gmsh.l"
+#line 104 "Gmsh.l"
 return tAcos;
 	YY_BREAK
 case 33:
 YY_RULE_SETUP
-#line 106 "Gmsh.l"
+#line 105 "Gmsh.l"
 return tAlias;
 	YY_BREAK
 case 34:
 YY_RULE_SETUP
-#line 107 "Gmsh.l"
+#line 106 "Gmsh.l"
 return tAliasWithOptions;
 	YY_BREAK
 case 35:
 YY_RULE_SETUP
-#line 108 "Gmsh.l"
+#line 107 "Gmsh.l"
 return tAsin;
 	YY_BREAK
 case 36:
 YY_RULE_SETUP
-#line 109 "Gmsh.l"
+#line 108 "Gmsh.l"
 return tAsin;
 	YY_BREAK
 case 37:
 YY_RULE_SETUP
-#line 110 "Gmsh.l"
+#line 109 "Gmsh.l"
 return tAtan;
 	YY_BREAK
 case 38:
 YY_RULE_SETUP
-#line 111 "Gmsh.l"
+#line 110 "Gmsh.l"
 return tAtan;
 	YY_BREAK
 case 39:
 YY_RULE_SETUP
-#line 112 "Gmsh.l"
+#line 111 "Gmsh.l"
 return tAtan2;
 	YY_BREAK
 case 40:
 YY_RULE_SETUP
-#line 113 "Gmsh.l"
+#line 112 "Gmsh.l"
 return tAtan2;
 	YY_BREAK
 case 41:
 YY_RULE_SETUP
-#line 115 "Gmsh.l"
+#line 114 "Gmsh.l"
 return tBezier;
 	YY_BREAK
 case 42:
 YY_RULE_SETUP
-#line 116 "Gmsh.l"
+#line 115 "Gmsh.l"
 return tBoundary;
 	YY_BREAK
 case 43:
 YY_RULE_SETUP
-#line 117 "Gmsh.l"
+#line 116 "Gmsh.l"
 return tBump;
 	YY_BREAK
 case 44:
 YY_RULE_SETUP
-#line 118 "Gmsh.l"
+#line 117 "Gmsh.l"
 return tBSpline;
 	YY_BREAK
 case 45:
 YY_RULE_SETUP
-#line 119 "Gmsh.l"
+#line 118 "Gmsh.l"
 return tBoundingBox;
 	YY_BREAK
 case 46:
 YY_RULE_SETUP
-#line 121 "Gmsh.l"
+#line 120 "Gmsh.l"
 return tCeil;
 	YY_BREAK
 case 47:
 YY_RULE_SETUP
-#line 122 "Gmsh.l"
+#line 121 "Gmsh.l"
 return tCircle;
 	YY_BREAK
 case 48:
 YY_RULE_SETUP
-#line 123 "Gmsh.l"
+#line 122 "Gmsh.l"
 return tCoherence;
 	YY_BREAK
 case 49:
 YY_RULE_SETUP
-#line 124 "Gmsh.l"
+#line 123 "Gmsh.l"
 return tCombine;
 	YY_BREAK
 case 50:
 YY_RULE_SETUP
-#line 125 "Gmsh.l"
+#line 124 "Gmsh.l"
 return tCosh;
 	YY_BREAK
 case 51:
 YY_RULE_SETUP
-#line 126 "Gmsh.l"
+#line 125 "Gmsh.l"
 return tCos;
 	YY_BREAK
 case 52:
 YY_RULE_SETUP
-#line 127 "Gmsh.l"
+#line 126 "Gmsh.l"
 return tCharacteristic;
 	YY_BREAK
 case 53:
 YY_RULE_SETUP
-#line 128 "Gmsh.l"
+#line 127 "Gmsh.l"
 return tComplex;
 	YY_BREAK
 case 54:
 YY_RULE_SETUP
-#line 129 "Gmsh.l"
+#line 128 "Gmsh.l"
 return tColor;
 	YY_BREAK
 case 55:
 YY_RULE_SETUP
-#line 130 "Gmsh.l"
+#line 129 "Gmsh.l"
 return tColorTable;
 	YY_BREAK
 case 56:
 YY_RULE_SETUP
-#line 131 "Gmsh.l"
+#line 130 "Gmsh.l"
 return tCoordinates;
 	YY_BREAK
 case 57:
 YY_RULE_SETUP
-#line 132 "Gmsh.l"
+#line 131 "Gmsh.l"
 return tSpline;
 	YY_BREAK
 case 58:
 YY_RULE_SETUP
-#line 133 "Gmsh.l"
+#line 132 "Gmsh.l"
 return tCall;
 	YY_BREAK
 case 59:
 YY_RULE_SETUP
-#line 135 "Gmsh.l"
+#line 134 "Gmsh.l"
 return tDelete;
 	YY_BREAK
 case 60:
 YY_RULE_SETUP
-#line 136 "Gmsh.l"
+#line 135 "Gmsh.l"
 return tDilate;
 	YY_BREAK
 case 61:
 YY_RULE_SETUP
-#line 137 "Gmsh.l"
+#line 136 "Gmsh.l"
 return tDuplicata;
 	YY_BREAK
 case 62:
 YY_RULE_SETUP
-#line 138 "Gmsh.l"
+#line 137 "Gmsh.l"
 return tDraw;
 	YY_BREAK
 case 63:
 YY_RULE_SETUP
-#line 140 "Gmsh.l"
+#line 139 "Gmsh.l"
 return tExp;
 	YY_BREAK
 case 64:
 YY_RULE_SETUP
-#line 141 "Gmsh.l"
+#line 140 "Gmsh.l"
 return tEllipse;
 	YY_BREAK
 case 65:
 YY_RULE_SETUP
-#line 142 "Gmsh.l"
+#line 141 "Gmsh.l"
 return tEllipse;
 	YY_BREAK
 case 66:
 YY_RULE_SETUP
-#line 143 "Gmsh.l"
+#line 142 "Gmsh.l"
 return tExtrude;
 	YY_BREAK
 case 67:
 YY_RULE_SETUP
-#line 144 "Gmsh.l"
+#line 143 "Gmsh.l"
 return tElliptic;
 	YY_BREAK
 case 68:
 YY_RULE_SETUP
-#line 145 "Gmsh.l"
+#line 144 "Gmsh.l"
 return tEndFor;
 	YY_BREAK
 case 69:
 YY_RULE_SETUP
-#line 146 "Gmsh.l"
+#line 145 "Gmsh.l"
 return tEndIf;
 	YY_BREAK
 case 70:
 YY_RULE_SETUP
-#line 147 "Gmsh.l"
+#line 146 "Gmsh.l"
 return tEuclidian;
 	YY_BREAK
 case 71:
 YY_RULE_SETUP
-#line 148 "Gmsh.l"
+#line 147 "Gmsh.l"
 return tExit;
 	YY_BREAK
 case 72:
 YY_RULE_SETUP
-#line 150 "Gmsh.l"
+#line 149 "Gmsh.l"
 return tFabs;
 	YY_BREAK
 case 73:
 YY_RULE_SETUP
-#line 151 "Gmsh.l"
+#line 150 "Gmsh.l"
 return tField;
 	YY_BREAK
 case 74:
 YY_RULE_SETUP
-#line 152 "Gmsh.l"
+#line 151 "Gmsh.l"
 return tFloor;
 	YY_BREAK
 case 75:
 YY_RULE_SETUP
-#line 153 "Gmsh.l"
+#line 152 "Gmsh.l"
 return tFmod;
 	YY_BREAK
 case 76:
 YY_RULE_SETUP
-#line 154 "Gmsh.l"
+#line 153 "Gmsh.l"
 return tFor;
 	YY_BREAK
 case 77:
 YY_RULE_SETUP
-#line 155 "Gmsh.l"
+#line 154 "Gmsh.l"
 return tFunction;
 	YY_BREAK
 case 78:
 YY_RULE_SETUP
-#line 157 "Gmsh.l"
+#line 156 "Gmsh.l"
 return tGetValue;
 	YY_BREAK
 case 79:
 YY_RULE_SETUP
-#line 158 "Gmsh.l"
+#line 157 "Gmsh.l"
 return tGMSH_MAJOR_VERSION;
 	YY_BREAK
 case 80:
 YY_RULE_SETUP
-#line 159 "Gmsh.l"
+#line 158 "Gmsh.l"
 return tGMSH_MINOR_VERSION;
 	YY_BREAK
 case 81:
 YY_RULE_SETUP
-#line 160 "Gmsh.l"
+#line 159 "Gmsh.l"
 return tGMSH_PATCH_VERSION;
 	YY_BREAK
 case 82:
 YY_RULE_SETUP
-#line 162 "Gmsh.l"
+#line 161 "Gmsh.l"
 return tHide;
 	YY_BREAK
 case 83:
 YY_RULE_SETUP
-#line 163 "Gmsh.l"
+#line 162 "Gmsh.l"
 return tHole;
 	YY_BREAK
 case 84:
 YY_RULE_SETUP
-#line 164 "Gmsh.l"
+#line 163 "Gmsh.l"
 return tHypot;
 	YY_BREAK
 case 85:
 YY_RULE_SETUP
-#line 166 "Gmsh.l"
+#line 165 "Gmsh.l"
 return tIn;
 	YY_BREAK
 case 86:
 YY_RULE_SETUP
-#line 167 "Gmsh.l"
+#line 166 "Gmsh.l"
 return tIf;
 	YY_BREAK
 case 87:
 YY_RULE_SETUP
-#line 168 "Gmsh.l"
+#line 167 "Gmsh.l"
 return tIntersect;
 	YY_BREAK
 case 88:
 YY_RULE_SETUP
-#line 170 "Gmsh.l"
+#line 169 "Gmsh.l"
 return tKnots;
 	YY_BREAK
 case 89:
 YY_RULE_SETUP
-#line 172 "Gmsh.l"
+#line 171 "Gmsh.l"
 return tLength;
 	YY_BREAK
 case 90:
 YY_RULE_SETUP
-#line 173 "Gmsh.l"
+#line 172 "Gmsh.l"
 return tLine;
 	YY_BREAK
 case 91:
 YY_RULE_SETUP
-#line 174 "Gmsh.l"
+#line 173 "Gmsh.l"
 return tLoop;
 	YY_BREAK
 case 92:
 YY_RULE_SETUP
-#line 175 "Gmsh.l"
+#line 174 "Gmsh.l"
 return tLog;
 	YY_BREAK
 case 93:
 YY_RULE_SETUP
-#line 176 "Gmsh.l"
+#line 175 "Gmsh.l"
 return tLog10;
 	YY_BREAK
 case 94:
 YY_RULE_SETUP
-#line 177 "Gmsh.l"
+#line 176 "Gmsh.l"
 return tLayers;
 	YY_BREAK
 case 95:
 YY_RULE_SETUP
-#line 179 "Gmsh.l"
+#line 178 "Gmsh.l"
 return tModulo;
 	YY_BREAK
 case 96:
 YY_RULE_SETUP
-#line 180 "Gmsh.l"
+#line 179 "Gmsh.l"
 return tMPI_Rank;
 	YY_BREAK
 case 97:
 YY_RULE_SETUP
-#line 181 "Gmsh.l"
+#line 180 "Gmsh.l"
 return tMPI_Size;
 	YY_BREAK
 case 98:
 YY_RULE_SETUP
-#line 183 "Gmsh.l"
+#line 182 "Gmsh.l"
 return tNurbs;
 	YY_BREAK
 case 99:
 YY_RULE_SETUP
-#line 185 "Gmsh.l"
+#line 184 "Gmsh.l"
 return tOrder;
 	YY_BREAK
 case 100:
 YY_RULE_SETUP
-#line 187 "Gmsh.l"
+#line 186 "Gmsh.l"
 return tPhysical;
 	YY_BREAK
 case 101:
 YY_RULE_SETUP
-#line 188 "Gmsh.l"
+#line 187 "Gmsh.l"
 return tPi;
 	YY_BREAK
 case 102:
 YY_RULE_SETUP
-#line 189 "Gmsh.l"
+#line 188 "Gmsh.l"
 return tPlane;
 	YY_BREAK
 case 103:
 YY_RULE_SETUP
-#line 190 "Gmsh.l"
+#line 189 "Gmsh.l"
 return tPoint;
 	YY_BREAK
 case 104:
 YY_RULE_SETUP
-#line 191 "Gmsh.l"
+#line 190 "Gmsh.l"
 return tProgression;
 	YY_BREAK
 case 105:
 YY_RULE_SETUP
-#line 192 "Gmsh.l"
+#line 191 "Gmsh.l"
 return tProgression;
 	YY_BREAK
 case 106:
 YY_RULE_SETUP
-#line 193 "Gmsh.l"
+#line 192 "Gmsh.l"
 return tParametric;
 	YY_BREAK
 case 107:
 YY_RULE_SETUP
-#line 194 "Gmsh.l"
+#line 193 "Gmsh.l"
 return tPolarSphere;
 	YY_BREAK
 case 108:
 YY_RULE_SETUP
-#line 195 "Gmsh.l"
+#line 194 "Gmsh.l"
 return tPrintf;
 	YY_BREAK
 case 109:
 YY_RULE_SETUP
-#line 196 "Gmsh.l"
+#line 195 "Gmsh.l"
 return tPlugin;
 	YY_BREAK
 case 110:
 YY_RULE_SETUP
-#line 198 "Gmsh.l"
+#line 197 "Gmsh.l"
 return tRecombine;
 	YY_BREAK
 case 111:
 YY_RULE_SETUP
-#line 199 "Gmsh.l"
+#line 198 "Gmsh.l"
 return tRotate;
 	YY_BREAK
 case 112:
 YY_RULE_SETUP
-#line 200 "Gmsh.l"
+#line 199 "Gmsh.l"
 return tRuled;
 	YY_BREAK
 case 113:
 YY_RULE_SETUP
-#line 201 "Gmsh.l"
+#line 200 "Gmsh.l"
 return tRand;
 	YY_BREAK
 case 114:
 YY_RULE_SETUP
-#line 202 "Gmsh.l"
+#line 201 "Gmsh.l"
 return tReturn;
 	YY_BREAK
 case 115:
 YY_RULE_SETUP
-#line 204 "Gmsh.l"
+#line 203 "Gmsh.l"
 return tSmoother;
 	YY_BREAK
 case 116:
 YY_RULE_SETUP
-#line 205 "Gmsh.l"
+#line 204 "Gmsh.l"
 return tSqrt;
 	YY_BREAK
 case 117:
 YY_RULE_SETUP
-#line 206 "Gmsh.l"
+#line 205 "Gmsh.l"
 return tSin;
 	YY_BREAK
 case 118:
 YY_RULE_SETUP
-#line 207 "Gmsh.l"
+#line 206 "Gmsh.l"
 return tSinh;
 	YY_BREAK
 case 119:
 YY_RULE_SETUP
-#line 208 "Gmsh.l"
+#line 207 "Gmsh.l"
 return tSphere;
 	YY_BREAK
 case 120:
 YY_RULE_SETUP
-#line 209 "Gmsh.l"
+#line 208 "Gmsh.l"
 return tSpline;
 	YY_BREAK
 case 121:
 YY_RULE_SETUP
-#line 210 "Gmsh.l"
+#line 209 "Gmsh.l"
 return tSurface;
 	YY_BREAK
 case 122:
 YY_RULE_SETUP
-#line 211 "Gmsh.l"
+#line 210 "Gmsh.l"
 return tSymmetry;
 	YY_BREAK
 case 123:
 YY_RULE_SETUP
-#line 212 "Gmsh.l"
+#line 211 "Gmsh.l"
 return tSprintf;
 	YY_BREAK
 case 124:
 YY_RULE_SETUP
-#line 213 "Gmsh.l"
+#line 212 "Gmsh.l"
 return tStrCat;
 	YY_BREAK
 case 125:
 YY_RULE_SETUP
-#line 214 "Gmsh.l"
+#line 213 "Gmsh.l"
 return tStrPrefix;
 	YY_BREAK
 case 126:
 YY_RULE_SETUP
-#line 215 "Gmsh.l"
+#line 214 "Gmsh.l"
 return tStrRelative;
 	YY_BREAK
 case 127:
 YY_RULE_SETUP
-#line 216 "Gmsh.l"
+#line 215 "Gmsh.l"
 return tShow;
 	YY_BREAK
 case 128:
 YY_RULE_SETUP
-#line 218 "Gmsh.l"
+#line 217 "Gmsh.l"
 return tTransfinite;
 	YY_BREAK
 case 129:
 YY_RULE_SETUP
-#line 219 "Gmsh.l"
+#line 218 "Gmsh.l"
 return tTranslate;
 	YY_BREAK
 case 130:
 YY_RULE_SETUP
-#line 220 "Gmsh.l"
+#line 219 "Gmsh.l"
 return tTanh;
 	YY_BREAK
 case 131:
 YY_RULE_SETUP
-#line 221 "Gmsh.l"
+#line 220 "Gmsh.l"
 return tTan;
 	YY_BREAK
 case 132:
 YY_RULE_SETUP
-#line 222 "Gmsh.l"
+#line 221 "Gmsh.l"
 return tToday;
 	YY_BREAK
 case 133:
 YY_RULE_SETUP
-#line 224 "Gmsh.l"
+#line 223 "Gmsh.l"
 return tUsing;
 	YY_BREAK
 case 134:
 YY_RULE_SETUP
-#line 226 "Gmsh.l"
+#line 225 "Gmsh.l"
 return tVolume;
 	YY_BREAK
 case 135:
 YY_RULE_SETUP
-#line 228 "Gmsh.l"
+#line 227 "Gmsh.l"
 return tText2D;
 	YY_BREAK
 case 136:
 YY_RULE_SETUP
-#line 229 "Gmsh.l"
+#line 228 "Gmsh.l"
 return tText3D;
 	YY_BREAK
 case 137:
 YY_RULE_SETUP
-#line 230 "Gmsh.l"
+#line 229 "Gmsh.l"
 return tInterpolationScheme;
 	YY_BREAK
 case 138:
 YY_RULE_SETUP
-#line 231 "Gmsh.l"
+#line 230 "Gmsh.l"
 return tTime;
 	YY_BREAK
 case 139:
-#line 234 "Gmsh.l"
+#line 233 "Gmsh.l"
 case 140:
-#line 235 "Gmsh.l"
+#line 234 "Gmsh.l"
 case 141:
-#line 236 "Gmsh.l"
+#line 235 "Gmsh.l"
 case 142:
 YY_RULE_SETUP
-#line 236 "Gmsh.l"
+#line 235 "Gmsh.l"
 { gmsh_yylval.d = atof((char *)gmsh_yytext); return tDOUBLE; }
 	YY_BREAK
 case 143:
 YY_RULE_SETUP
-#line 238 "Gmsh.l"
+#line 237 "Gmsh.l"
 { gmsh_yylval.c = strsave((char*)gmsh_yytext); return tSTRING; }
 	YY_BREAK
 case 144:
 YY_RULE_SETUP
-#line 240 "Gmsh.l"
+#line 239 "Gmsh.l"
 return gmsh_yytext[0];
 	YY_BREAK
 case 145:
 YY_RULE_SETUP
-#line 242 "Gmsh.l"
+#line 241 "Gmsh.l"
 ECHO;
 	YY_BREAK
-#line 1855 "Gmsh.yy.cpp"
+#line 1854 "Gmsh.yy.cpp"
 case YY_STATE_EOF(INITIAL):
 	yyterminate();
 
@@ -2837,7 +2836,7 @@ void gmsh_yyfree (void * ptr )
 
 #define YYTABLES_NAME "yytables"
 
-#line 242 "Gmsh.l"
+#line 241 "Gmsh.l"
 
 
 
@@ -2926,7 +2925,7 @@ void skip_until(const char *skip, const char *until)
       if(skip && chars[0] == skip[0]) break;
     }
 
-    l = MAX(l_skip,l_until);
+    l = std::max(l_skip,l_until);
     if(l >= (int)sizeof(chars)){
       Msg::Error("Search pattern too long in skip_until");
       return;
diff --git a/Plugin/ExtractElements.cpp b/Plugin/ExtractElements.cpp
index 908cb97446cd2872bf0c03f33cc297f7064e4cf0..cbc3759dc4f00424b36e09a2a1859abbb3962e75 100644
--- a/Plugin/ExtractElements.cpp
+++ b/Plugin/ExtractElements.cpp
@@ -1,4 +1,4 @@
-// $Id: ExtractElements.cpp,v 1.14 2008-05-04 08:31:23 geuzaine Exp $
+// $Id: ExtractElements.cpp,v 1.15 2008-07-10 13:29:30 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -22,6 +22,8 @@
 #include "ExtractElements.h"
 #include "Numeric.h"
 
+#define SQU(a)      ((a)*(a))
+
 StringXNumber ExtractElementsOptions_Number[] = {
   {GMSH_FULLRC, "MinVal", NULL, 0.},
   {GMSH_FULLRC, "MaxVal", NULL, 1.},
@@ -93,7 +95,7 @@ static void extract(List_T *inList, int inNb,
         d += v[0];
         break;
       case 3 : // vector
-        d += sqrt(DSQR(v[0]) + DSQR(v[1]) + DSQR(v[2]));
+        d += sqrt(SQU(v[0]) + SQU(v[1]) + SQU(v[2]));
         break;
       case 9 : // tensor
         d += ComputeVonMises(v);
diff --git a/Post/ColorTable.cpp b/Post/ColorTable.cpp
index b80f70783e2f4dda1d6aa55865c848d800655344..466ebde493a5f4400a257ed9e58c630e55cc0c3f 100644
--- a/Post/ColorTable.cpp
+++ b/Post/ColorTable.cpp
@@ -1,4 +1,4 @@
-// $Id: ColorTable.cpp,v 1.6 2008-05-04 08:31:23 geuzaine Exp $
+// $Id: ColorTable.cpp,v 1.7 2008-07-10 13:29:30 geuzaine Exp $
 //
 // Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
 //
@@ -316,9 +316,9 @@ void ColorTable_Recompute(GmshColorTable * ct)
       b = (int)(255. * 1.);
       break;
     case 19: // matlab "copper"
-      r = (int)(255. * DMIN(1., gray(s - bias) * 1.25));
-      g = (int)(255. * DMIN(1., gray(s - bias) * 0.7812));
-      b = (int)(255. * DMIN(1., gray(s - bias) * 0.4975));
+      r = (int)(255. * std::min(1., gray(s - bias) * 1.25));
+      g = (int)(255. * std::min(1., gray(s - bias) * 0.7812));
+      b = (int)(255. * std::min(1., gray(s - bias) * 0.4975));
       break;
     default:
       r = g = b = 0;
diff --git a/Post/shapeFunctions.h b/Post/shapeFunctions.h
index 6f41b9400519cb4f0994504154c21a97c04bfccc..8c60099ff84ec41cffbeb7d43ff3fd3109e6c0cf 100644
--- a/Post/shapeFunctions.h
+++ b/Post/shapeFunctions.h
@@ -23,8 +23,9 @@
 #include "Numeric.h"
 #include "Message.h"
 
-#define ONE  (1. + 1.e-6)
-#define ZERO (-1.e-6)
+#define ONE     (1. + 1.e-6)
+#define ZERO    (-1.e-6)
+#define SQU(a)  ((a)*(a))
 
 class element{
 protected:
@@ -72,9 +73,9 @@ public:
         prodve(a, b, c);
         jac[2][0] = c[0]; jac[2][1] = c[1]; jac[2][2] = c[2]; 
       }
-      return sqrt(SQR(jac[0][0] * jac[1][1] - jac[0][1] * jac[1][0]) +
-                  SQR(jac[0][2] * jac[1][0] - jac[0][0] * jac[1][2]) +
-                  SQR(jac[0][1] * jac[1][2] - jac[0][2] * jac[1][1]));
+      return sqrt(SQU(jac[0][0] * jac[1][1] - jac[0][1] * jac[1][0]) +
+                  SQU(jac[0][2] * jac[1][0] - jac[0][0] * jac[1][2]) +
+                  SQU(jac[0][1] * jac[1][2] - jac[0][2] * jac[1][1]));
     case 1:
       for(int i = 0; i < getNumNodes(); i++) {
         getGradShapeFunction(i, u, v, w, s);
@@ -94,7 +95,7 @@ public:
         jac[1][0] = b[0]; jac[1][1] = b[1]; jac[1][2] = b[2]; 
         jac[2][0] = c[0]; jac[2][1] = c[1]; jac[2][2] = c[2]; 
       }
-      return sqrt(SQR(jac[0][0])+SQR(jac[0][1])+SQR(jac[0][2]));
+      return sqrt(SQU(jac[0][0])+SQU(jac[0][1])+SQU(jac[0][2]));
     default:
       return 1.;
     }
@@ -225,7 +226,7 @@ public:
       double wn = uvw[2] +
         inv[0][2] * (xyz[0] - xn) + inv[1][2] * (xyz[1] - yn) + inv[2][2] * (xyz[2] - zn) ;
       
-      error = sqrt(SQR(un - uvw[0]) + SQR(vn - uvw[1]) + SQR(wn - uvw[2]));
+      error = sqrt(SQU(un - uvw[0]) + SQU(vn - uvw[1]) + SQU(wn - uvw[2]));
       uvw[0] = un;
       uvw[1] = vn;
       uvw[2] = wn;
@@ -240,7 +241,7 @@ public:
     for(int i = 0; i < getNumEdges(); i++){
       int n1, n2;
       getEdge(i, n1, n2);
-      double d = sqrt(SQR(_x[n1]-_x[n2]) + SQR(_y[n1]-_y[n2]) + SQR(_z[n1]-_z[n2]));
+      double d = sqrt(SQU(_x[n1]-_x[n2]) + SQU(_y[n1]-_y[n2]) + SQU(_z[n1]-_z[n2]));
       if(d > max) max = d;
     }
     return max;
@@ -890,16 +891,16 @@ public:
       switch(num) {
       case 0  : s[0] = 0.25 * ( -(1.-v) + v*w/(1.-w) ) ;
                 s[1] = 0.25 * ( -(1.-u) + u*w/(1.-w) ) ;
-                s[2] = 0.25 * ( -1.     + u*v/SQR(1.-w) ) ; break ;
+                s[2] = 0.25 * ( -1.     + u*v/SQU(1.-w) ) ; break ;
       case 1  : s[0] = 0.25 * (  (1.-v) + v*w/(1.-w) ) ;
                 s[1] = 0.25 * ( -(1.+u) + u*w/(1.-w) ) ;
-                s[2] = 0.25 * ( -1.     + u*v/SQR(1.-w) ) ; break ;
+                s[2] = 0.25 * ( -1.     + u*v/SQU(1.-w) ) ; break ;
       case 2  : s[0] = 0.25 * (  (1.+v) + v*w/(1.-w) ) ;
                 s[1] = 0.25 * (  (1.+u) + u*w/(1.-w) ) ;
-                s[2] = 0.25 * ( -1.     + u*v/SQR(1.-w) ) ; break ;
+                s[2] = 0.25 * ( -1.     + u*v/SQU(1.-w) ) ; break ;
       case 3  : s[0] = 0.25 * ( -(1.+v) + v*w/(1.-w) ) ;
                 s[1] = 0.25 * (  (1.-u) + u*w/(1.-w) ) ;
-                s[2] = 0.25 * ( -1.     + u*v/SQR(1.-w) ) ; break ;
+                s[2] = 0.25 * ( -1.     + u*v/SQU(1.-w) ) ; break ;
       case 4  : s[0] = 0. ; 
                 s[1] = 0. ;
                 s[2] = 1. ; break ;
@@ -939,5 +940,6 @@ class elementFactory{
 
 #undef ONE
 #undef ZERO
+#undef SQU
 
 #endif