From a3c77b7ec80b82352efc7c2d00c2cf7346e5e093 Mon Sep 17 00:00:00 2001
From: Amaury Johnan <amjohnen@gmail.com>
Date: Fri, 18 Oct 2013 12:42:22 +0000
Subject: [PATCH] fix getFaceVertices for serendipity hex & tet

---
 Geo/MElement.h     | 12 +++++-------
 Geo/MHexahedron.h  | 19 ++++++++++++++-----
 Geo/MPrism.cpp     |  2 +-
 Geo/MPyramid.h     | 18 ++++++++++++++----
 Geo/MTetrahedron.h | 19 ++++++++++++++-----
 5 files changed, 48 insertions(+), 22 deletions(-)

diff --git a/Geo/MElement.h b/Geo/MElement.h
index 5e50155ca2..67c280f409 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -130,12 +130,6 @@ class MElement
   virtual int getNumEdgesRep() = 0;
   virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) = 0;
 
-  // get all the vertices on an edge
-  virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
-  {
-    v.resize(0);
-  }
-
   // get the faces
   virtual int getNumFaces() = 0;
   virtual MFace getFace(int num) = 0;
@@ -150,7 +144,11 @@ class MElement
   virtual int getNumFacesRep() = 0;
   virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n) = 0;
 
-  // get all the vertices on a face
+  // get all the vertices on a edge or a face
+  virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
+  {
+    v.resize(0);
+  }
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(0);
diff --git a/Geo/MHexahedron.h b/Geo/MHexahedron.h
index d9cd208432..fa43326557 100644
--- a/Geo/MHexahedron.h
+++ b/Geo/MHexahedron.h
@@ -508,10 +508,16 @@ class MHexahedronN : public MHexahedron {
   }
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
-    v.resize((_order+1)*(_order+1));
-    MHexahedron::_getFaceVertices(num, v);
+    if (ElementType::SerendipityFromTag(getTypeForMSH()) > 0) {
+      v.resize(4 * _order);
+    }
+    else {
+      v.resize((_order+1) * (_order+1));
+    }
 
+    MHexahedron::_getFaceVertices(num, v);
     int count = 3;
+
     int n = _order-1;
     for (int i = 0; i < 4; i++) {
       if(faces2edge_hexa(num, i) > 0)
@@ -525,9 +531,12 @@ class MHexahedronN : public MHexahedron {
         for (int j = n-1; j >= 0; j--) v[++count] = _vs[n*edge_num + j];
       }
     }
-    int start = 12 * n + num * n*n;
-    for (int i = 0; i < n*n; i++){
-      v[++count] = _vs[start + i];
+
+    if (v.size() > count + 1) {
+      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/MPrism.cpp b/Geo/MPrism.cpp
index 3f11e25a5e..66b8f121d6 100644
--- a/Geo/MPrism.cpp
+++ b/Geo/MPrism.cpp
@@ -479,7 +479,7 @@ void _addFaceNodes(int num, int order, const std::vector<MVertex*> &vs,
 // To be tested
 void MPrismN::getFaceVertices(const int num, std::vector<MVertex*> &v) const
 {
-
+  // FIXME serendipity case
   static const int edge[5][4] = {
       {1, 3, 0, -1},
       {6, 8, 7, -1},
diff --git a/Geo/MPyramid.h b/Geo/MPyramid.h
index b230682f28..b8a4270e75 100644
--- a/Geo/MPyramid.h
+++ b/Geo/MPyramid.h
@@ -266,15 +266,25 @@ class MPyramidN : public MPyramid {
   }
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
+    if (ElementType::SerendipityFromTag(getTypeForMSH()) > 0) {
+      num == 4 ? v.resize(4 * _order)
+               : v.resize(3 * _order);
+    }
+    else {
+      num == 4 ? v.resize((_order+1) * (_order+1))
+               : v.resize((_order+1) * (_order+2) / 2);
+    }
+
+    // FIXME continue fix serendipity
+
     int j = 3;
     if (num == 4) {
       j = 4;
-      v.resize(_order * _order);
-    }
-    else {
-      v.resize(3 + 3 * (_order - 1) + (_order-1) * (_order - 2) /2);
     }
+
     MPyramid::_getFaceVertices(num, v);
+    //int count = num == 4 ? 3 : 2;
+
     int nbVQ =  (_order-1)*(_order-1);
     int nbVT = (_order - 1) * (_order - 2) / 2;
     const int ie = (num == 4) ? 4*nbVT + nbVQ : (num+1)*nbVT;
diff --git a/Geo/MTetrahedron.h b/Geo/MTetrahedron.h
index b33033a13b..091df43fb6 100644
--- a/Geo/MTetrahedron.h
+++ b/Geo/MTetrahedron.h
@@ -349,10 +349,16 @@ class MTetrahedronN : public MTetrahedron {
   }
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
-    v.resize((_order+1) * (_order+2) / 2);
-    MTetrahedron::_getFaceVertices(num, v);
+    if (ElementType::SerendipityFromTag(getTypeForMSH()) > 0) {
+      v.resize(3 * _order);
+    }
+    else {
+      v.resize((_order+1) * (_order+2) / 2);
+    }
 
+    MTetrahedron::_getFaceVertices(num, v);
     int count = 2;
+
     int n = _order-1;
     for (int i = 0; i < 3; i++) {
       if(faces2edge_tetra(num, i) > 0)
@@ -366,9 +372,12 @@ class MTetrahedronN : public MTetrahedron {
         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];
+
+    if (v.size() > count + 1) {
+      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
-- 
GitLab