diff --git a/Geo/MHexahedron.h b/Geo/MHexahedron.h
index 871881b8ee24c2098e9f2855fa1112236cb7cb7e..4bd51e2bec383ae453c3481cc6d22dabb0a15304 100644
--- a/Geo/MHexahedron.h
+++ b/Geo/MHexahedron.h
@@ -185,6 +185,20 @@ class MHexahedron : public MElement {
     };
     return f[face][vert];
   }
+  static int faces2edge_hexa(const int face, const int edge)
+  {
+    // return -iedge - 1 if edge is inverted
+    //         iedge + 1 otherwise
+    static const int e[6][4] = {
+      {2, -6,  -4,  -1},
+      {1,  5,  -9,  -3},
+      {3, 10,  -8,  -2},
+      {4,  7, -11,  -5},
+      {6,  8, -12,  -7},
+      {9, 11,  12, -10}
+    };
+    return e[face][edge];
+  }
 };
 
 /*
@@ -499,40 +513,24 @@ class MHexahedronN : public MHexahedron {
   {
     v.resize((_order+1)*(_order+1));
     MHexahedron::_getFaceVertices(num, v);
-//    static const int f[6][4] = {
-//      {0, 3, 2, 1},
-//      {0, 1, 5, 4},
-//      {0, 4, 7, 3},
-//      {1, 2, 6, 5},
-//      {2, 3, 7, 6},
-//      {4, 5, 6, 7}
-//    };
-    // this is the local edge number indexed from 1. A minus sign is used to indicate that the nodes of the edge must be inverted in order to obtain a MQuandrangleN
-    static const int f[6][4] = {
-      {2,-6,-4,-1},
-      {1,5,-9,-3},
-      {3,10,-8,-2},
-      {4,7,-11,-5},
-      {6,8,-12,-7},
-      {9,11,12,-10}
-    };
-    int count = 4;
-    for (int i = 0; i < 4; i++){
-      if(f[num][i]>0)
+
+    int count = 3;
+    int n = _order-1;
+    for (int i = 0; i < 4; i++) {
+      if(faces2edge_hexa(num, i) > 0)
       {
-        int edge_num = f[num][i]-1;      
-        for (int j = 0; j < _order - 1; j++) v[count++] = _vs[(_order-1)*edge_num+j];
+        int edge_num = faces2edge_hexa(num, i) - 1;
+        for (int j = 0; j < n; j++) v[++count] = _vs[n*edge_num + j];
       }
       else
       {
-        int edge_num = -f[num][i]-1;
-        for (int j = _order-2; j > - 1; j--) v[count++] = _vs[(_order-1)*edge_num+j];
+        int edge_num = -faces2edge_hexa(num, i) - 1;
+        for (int j = n-1; j >= 0; j--) v[++count] = _vs[n*edge_num + j];
       }
     }
-    int N = _order - 1;
-    int start = 12 * N + num * (_order - 1) * (_order - 1);// -8 as _vs has not the 8 first order nodes
-    for (int i = 0; i < (_order - 1) * (_order - 1); i++){
-      v[count++] = _vs[start + i];
+    int start = 12 * n + num * n*n;
+    for (int i = 0; i < n*n; i++){
+      v[++count] = _vs[start + i];
     }
   }
   virtual int getTypeForMSH() const
diff --git a/Geo/MTetrahedron.h b/Geo/MTetrahedron.h
index 8ac0f1e003857a535da5596423c869b48377744c..dbee7499b6b23eadf09c1a93da012686e3e0b753 100644
--- a/Geo/MTetrahedron.h
+++ b/Geo/MTetrahedron.h
@@ -176,6 +176,18 @@ class MTetrahedron : public MElement {
     };
     return f[face][vert];
   }
+  static int faces2edge_tetra(const int face, const int vert)
+  {
+    // return -iedge - 1 if edge is inverted
+    //         iedge + 1 otherwise
+    static const int e[4][3] = {
+      {-3, -2, -1},
+      { 1, -6,  4},
+      {-4,  5,  3},
+      { 6,  2, -5}
+    };
+    return e[face][vert];
+  }
 };
 
 /*
@@ -334,12 +346,27 @@ class MTetrahedronN : public MTetrahedron {
   }
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
-    v.resize(3 + 3 * (_order - 1) + (_order-1) * (_order - 2) /2);
+    v.resize((_order+1) * (_order+2) / 2);
     MTetrahedron::_getFaceVertices(num, v);
-    int j = 3;
-    int nbV = (_order - 1) * (_order - 2) / 2;
-    const int ie = (num+1)*nbV;
-    for(int i = num*nbV; i != ie; ++i) v[j++] = _vs[i];
+
+    int count = 2;
+    int n = _order-1;
+    for (int i = 0; i < 3; i++) {
+      if(faces2edge_tetra(num, i) > 0)
+      {
+        int edge_num = faces2edge_tetra(num, i) - 1;
+        for (int j = 0; j < n; j++) v[++count] = _vs[n*edge_num + j];
+      }
+      else
+      {
+        int edge_num = -faces2edge_tetra(num, i) - 1;
+        for (int j = n-1; j >= 0; j--) v[++count] = _vs[n*edge_num + j];
+      }
+    }
+    int start = 6 * n + num * (n-1)*n/2;
+    for (int i = 0; i < (n-1)*n/2; i++){
+      v[++count] = _vs[start + i];
+    }
   }
   virtual int getNumVolumeVertices() const
   {
diff --git a/Mesh/HighOrder.cpp b/Mesh/HighOrder.cpp
index 03ee51ccff0f548f64b423f1dbdedc6ca669eac5..c19b899daea10e647d3ba1f16609014436394985 100644
--- a/Mesh/HighOrder.cpp
+++ b/Mesh/HighOrder.cpp
@@ -1348,6 +1348,7 @@ void SetOrderN(GModel *m, int order, bool linear, bool incomplete, bool onlyVisi
   double worst;
   checkHighOrderTriangles("Surface mesh", m, bad, worst);
   checkHighOrderTetrahedron("Volume Mesh", m, bad, worst);
+  // FIXME : add other element check
 
   Msg::StatusBar(true, "Done meshing order %d (%g s)", order, t2 - t1);
 }