From a7069c3fdfba674066d5bcad01719cb922e4e553 Mon Sep 17 00:00:00 2001
From: Thomas Toulorge <thomas.toulorge@mines-paristech.fr>
Date: Tue, 29 Oct 2013 14:50:27 +0000
Subject: [PATCH] Added "adaptive" visualization of high-order meshes
 (subdivide only curved elements)

---
 Geo/GModelVertexArrays.cpp | 19 ++++++----
 Geo/MElement.h             |  8 ++---
 Geo/MElementCut.h          | 16 ++++-----
 Geo/MHexahedron.cpp        | 54 +++++++++++++++-------------
 Geo/MHexahedron.h          | 32 ++++++++---------
 Geo/MLine.cpp              | 50 ++++++++++++++------------
 Geo/MLine.h                | 16 ++++-----
 Geo/MPoint.h               |  8 ++---
 Geo/MPrism.cpp             | 54 +++++++++++++++-------------
 Geo/MPrism.h               | 32 ++++++++---------
 Geo/MPyramid.cpp           | 74 +++++++++++++++++++++-----------------
 Geo/MPyramid.h             | 16 ++++-----
 Geo/MQuadrangle.cpp        | 56 +++++++++++++++++++----------
 Geo/MQuadrangle.h          | 32 ++++++++---------
 Geo/MTetrahedron.cpp       | 38 +++++++++++++-------
 Geo/MTetrahedron.h         | 24 ++++++-------
 Geo/MTriangle.cpp          | 37 ++++++++++++-------
 Geo/MTriangle.h            | 24 ++++++-------
 18 files changed, 335 insertions(+), 255 deletions(-)

diff --git a/Geo/GModelVertexArrays.cpp b/Geo/GModelVertexArrays.cpp
index 233a90bb76..feb6246c74 100644
--- a/Geo/GModelVertexArrays.cpp
+++ b/Geo/GModelVertexArrays.cpp
@@ -19,6 +19,8 @@
 #include "VertexArray.h"
 #include "SmoothData.h"
 
+static const double curvedRepTol = 1.e-5;
+
 unsigned int getColorByEntity(GEntity *e)
 {
   if(e->getSelection()){ // selection
@@ -161,12 +163,14 @@ static void addSmoothNormals(GEntity *e, std::vector<T*> &elements)
 {
   for(unsigned int i = 0; i < elements.size(); i++){
     MElement *ele = elements[i];
+    const bool curved = (ele->getPolynomialOrder() > 1) &&
+                        (ele->maxDistToStraight() > curvedRepTol*ele->getInnerRadius());
     SPoint3 pc(0., 0., 0.);
     if(CTX::instance()->mesh.explode != 1.) pc = ele->barycenter();
-    for(int j = 0; j < ele->getNumFacesRep(); j++){
+    for(int j = 0; j < ele->getNumFacesRep(curved); j++){
       double x[3], y[3], z[3];
       SVector3 n[3];
-      ele->getFaceRep(j, x, y, z, n);
+      ele->getFaceRep(curved, j, x, y, z, n);
       for(int k = 0; k < 3; k++){
         if(CTX::instance()->mesh.explode != 1.){
           x[k] = pc[0] + CTX::instance()->mesh.explode * (x[k] - pc[0]);
@@ -191,15 +195,18 @@ static void addElementsInArrays(GEntity *e, std::vector<T*> &elements,
     unsigned int c = getColorByElement(ele);
     unsigned int col[4] = {c, c, c, c};
 
+    const bool curved = (ele->getPolynomialOrder() > 1) &&
+                        (ele->maxDistToStraight() > curvedRepTol*ele->getInnerRadius());
+
     SPoint3 pc(0., 0., 0.);
     if(CTX::instance()->mesh.explode != 1.) pc = ele->barycenter();
 
     if(edges){
       bool unique = e->dim() > 1 && !CTX::instance()->pickElements;
-      for(int j = 0; j < ele->getNumEdgesRep(); j++){
+      for(int j = 0; j < ele->getNumEdgesRep(curved); j++){
         double x[2], y[2], z[2];
         SVector3 n[2];
-        ele->getEdgeRep(j, x, y, z, n);
+        ele->getEdgeRep(curved, j, x, y, z, n);
         if(CTX::instance()->mesh.explode != 1.){
           for(int k = 0; k < 2; k++){
             x[k] = pc[0] + CTX::instance()->mesh.explode * (x[k] - pc[0]);
@@ -217,10 +224,10 @@ static void addElementsInArrays(GEntity *e, std::vector<T*> &elements,
     if(faces){
       bool unique = e->dim() > 2 && !CTX::instance()->pickElements;
       bool skin = e->dim() > 2 && CTX::instance()->mesh.drawSkinOnly;
-      for(int j = 0; j < ele->getNumFacesRep(); j++){
+      for(int j = 0; j < ele->getNumFacesRep(curved); j++){
         double x[3], y[3], z[3];
         SVector3 n[3];
-        ele->getFaceRep(j, x, y, z, n);
+        ele->getFaceRep(curved, j, x, y, z, n);
         if(CTX::instance()->mesh.explode != 1.){
           for(int k = 0; k < 3; k++){
             x[k] = pc[0] + CTX::instance()->mesh.explode * (x[k] - pc[0]);
diff --git a/Geo/MElement.h b/Geo/MElement.h
index 12e362dc63..cc19534964 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -127,8 +127,8 @@ class MElement
   }
 
   // get an edge representation for drawing
-  virtual int getNumEdgesRep() = 0;
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) = 0;
+  virtual int getNumEdgesRep(bool curved) = 0;
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) = 0;
 
   // get the faces
   virtual int getNumFaces() = 0;
@@ -141,8 +141,8 @@ class MElement
   }
 
   // get a face representation for drawing
-  virtual int getNumFacesRep() = 0;
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n) = 0;
+  virtual int getNumFacesRep(bool curved) = 0;
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) = 0;
 
   // get all the vertices on a edge or a face
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
diff --git a/Geo/MElementCut.h b/Geo/MElementCut.h
index ba3206ebfc..500c3d23a7 100644
--- a/Geo/MElementCut.h
+++ b/Geo/MElementCut.h
@@ -72,8 +72,8 @@ class MPolyhedron : public MElement {
   }
   virtual int getNumEdges() { return _edges.size(); }
   virtual MEdge getEdge(int num) const{ return _edges[num]; }
-  virtual int getNumEdgesRep() { return _edges.size(); }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved) { return _edges.size(); }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     MEdge e(getEdge(num));
     for(unsigned int i = 0; i < _faces.size(); i++)
@@ -89,8 +89,8 @@ class MPolyhedron : public MElement {
   }
   virtual int getNumFaces() { return _faces.size(); }
   virtual MFace getFace(int num) { return _faces[num]; }
-  virtual int getNumFacesRep() { return _faces.size(); }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved) { return _faces.size(); }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     _getFaceRep(_faces[num].getVertex(0), _faces[num].getVertex(1),
                 _faces[num].getVertex(2), x, y, z, n);
@@ -229,8 +229,8 @@ class MPolygon : public MElement {
   }
   virtual int getNumEdges() { return _edges.size(); }
   virtual MEdge getEdge(int num) const{ return _edges[num]; }
-  virtual int getNumEdgesRep() { return getNumEdges(); }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved) { return getNumEdges(); }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     MEdge e(getEdge(num));
     _getEdgeRep(e.getVertex(0), e.getVertex(1), x, y, z, n, 0);
@@ -243,8 +243,8 @@ class MPolygon : public MElement {
   }
   virtual int getNumFaces() { return 1; }
   virtual MFace getFace(int num) { return MFace(_vertices); }
-  virtual int getNumFacesRep() { return _parts.size(); }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved) { return _parts.size(); }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     _getFaceRep(_parts[num]->getVertex(0), _parts[num]->getVertex(1),
                 _parts[num]->getVertex(2), x, y, z, n);
diff --git a/Geo/MHexahedron.cpp b/Geo/MHexahedron.cpp
index 413efd7720..7e4a20ae73 100644
--- a/Geo/MHexahedron.cpp
+++ b/Geo/MHexahedron.cpp
@@ -163,31 +163,34 @@ static void _myGetEdgeRep(MHexahedron *hex, int num, double *x, double *y, doubl
   n[0] = n[1] = 1 ;
 }
 
-void MHexahedron20::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+void MHexahedron20::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) {
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MHexahedron::getEdgeRep(false, num, x, y, z, n);
 }
 
-void MHexahedron27::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+void MHexahedron27::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) {
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MHexahedron::getEdgeRep(false, num, x, y, z, n);
 }
 
-void MHexahedronN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+void MHexahedronN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) {
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MHexahedron::getEdgeRep(false, num, x, y, z, n);
 }
 
-int MHexahedron20::getNumEdgesRep()
+int MHexahedron20::getNumEdgesRep(bool curved)
 {
-  return 12 * CTX::instance()->mesh.numSubEdges;
+  return curved ? 12 * CTX::instance()->mesh.numSubEdges : 12;
 }
 
-int MHexahedron27::getNumEdgesRep()
+int MHexahedron27::getNumEdgesRep(bool curved)
 {
-  return 12 * CTX::instance()->mesh.numSubEdges;
+  return curved ? 12 * CTX::instance()->mesh.numSubEdges : 12;
 }
 
-int MHexahedronN::getNumEdgesRep()
+int MHexahedronN::getNumEdgesRep(bool curved)
 {
-  return 12 * CTX::instance()->mesh.numSubEdges;
+  return curved ? 12 * CTX::instance()->mesh.numSubEdges : 12;
 }
 
 const nodalBasis* MHexahedron::getFunctionSpace(int order) const
@@ -373,33 +376,36 @@ void _myGetFaceRep(MHexahedron *hex, int num, double *x, double *y, double *z,
 }
 
 
-void MHexahedron20::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MHexahedron20::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MHexahedron::getFaceRep(false, num, x, y, z, n);
 }
 
-void MHexahedron27::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MHexahedron27::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MHexahedron::getFaceRep(false, num, x, y, z, n);
 }
 
-void MHexahedronN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MHexahedronN::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MHexahedron::getFaceRep(false, num, x, y, z, n);
 }
 
-int MHexahedron20::getNumFacesRep()
+int MHexahedron20::getNumFacesRep(bool curved)
 {
-  return 6 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
+  return curved ? 6 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2) : 12;
 }
 
 
-int MHexahedron27::getNumFacesRep()
+int MHexahedron27::getNumFacesRep(bool curved)
 {
-  return 6 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
+  return curved ? 6 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2) : 12;
 }
 
-int MHexahedronN::getNumFacesRep()
+int MHexahedronN::getNumFacesRep(bool curved)
 {
-  return 6 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
+  return curved ? 6 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2) : 12;
 }
diff --git a/Geo/MHexahedron.h b/Geo/MHexahedron.h
index fe7e3016f3..17596c1a7f 100644
--- a/Geo/MHexahedron.h
+++ b/Geo/MHexahedron.h
@@ -71,8 +71,8 @@ class MHexahedron : public MElement {
   {
     return MEdge(_v[edges_hexa(num, 0)], _v[edges_hexa(num, 1)]);
   }
-  virtual int getNumEdgesRep(){ return 12; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 12; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[12] = {0, 0, 1, 0, 1, 0, 3, 2, 1, 2, 3, 4};
     MEdge e(getEdge(num));
@@ -94,8 +94,8 @@ class MHexahedron : public MElement {
   virtual double getInnerRadius();
   virtual double angleShapeMeasure();
   virtual void getFaceInfo (const MFace & face, int &ithFace, int &sign, int &rot)const;
-  virtual int getNumFacesRep(){ return 12; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved){ return 12; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[12][3] = {
       {0, 3, 2}, {0, 2, 1},
@@ -269,16 +269,16 @@ class MHexahedron20 : public MHexahedron {
     return getVertex(map[num]);
   }
   virtual int getNumEdgeVertices() const { return 12; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MHexahedron::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(8);
@@ -376,16 +376,16 @@ class MHexahedron27 : public MHexahedron {
   virtual int getNumEdgeVertices() const { return 12; }
   virtual int getNumFaceVertices() const { return 6; }
   virtual int getNumVolumeVertices() const { return 1; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MHexahedron::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(9);
@@ -498,8 +498,8 @@ class MHexahedronN : public MHexahedron {
     else
       return (_order - 1) * (_order - 1) * (_order - 1);
   }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(_order+1);
@@ -562,8 +562,8 @@ class MHexahedronN : public MHexahedron {
     Msg::Error("no tag matches a p%d hexahedron with %d vertices", _order, 8+_vs.size());
     return 0;
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void reverse()
   {
     Msg::Error("Reverse not implemented yet for MHexahedronN");
diff --git a/Geo/MLine.cpp b/Geo/MLine.cpp
index 4fd16e3f44..fa296bf2ff 100644
--- a/Geo/MLine.cpp
+++ b/Geo/MLine.cpp
@@ -48,36 +48,42 @@ double MLine::getVolume()
   return getLength();
 }
 
-int MLine3::getNumEdgesRep()
+int MLine3::getNumEdgesRep(bool curved)
 {
-  return  CTX::instance()->mesh.numSubEdges;
+  return curved ? CTX::instance()->mesh.numSubEdges : 1;
 }
 
-void MLine3::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MLine3::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  int numSubEdges = CTX::instance()->mesh.numSubEdges;
-  SPoint3 pnt1, pnt2;
-  pnt(-1. + 2 * (double)num / numSubEdges, 0., 0., pnt1);
-  pnt(-1. + 2 * (double)(num + 1) / numSubEdges, 0., 0, pnt2);
-  x[0] = pnt1.x(); x[1] = pnt2.x();
-  y[0] = pnt1.y(); y[1] = pnt2.y();
-  z[0] = pnt1.z(); z[1] = pnt2.z();
-  n[0] = n[1] = MEdge(_v[0], _v[1]).normal();
+  if (curved) {
+    int numSubEdges = CTX::instance()->mesh.numSubEdges;
+    SPoint3 pnt1, pnt2;
+    pnt(-1. + 2 * (double)num / numSubEdges, 0., 0., pnt1);
+    pnt(-1. + 2 * (double)(num + 1) / numSubEdges, 0., 0, pnt2);
+    x[0] = pnt1.x(); x[1] = pnt2.x();
+    y[0] = pnt1.y(); y[1] = pnt2.y();
+    z[0] = pnt1.z(); z[1] = pnt2.z();
+    n[0] = n[1] = MEdge(_v[0], _v[1]).normal();
+  }
+  else MLine::getEdgeRep(false, num, x, y, z, n);
 }
 
-int MLineN::getNumEdgesRep()
+int MLineN::getNumEdgesRep(bool curved)
 {
-  return  CTX::instance()->mesh.numSubEdges;
+  return curved ? CTX::instance()->mesh.numSubEdges : 1;
 }
 
-void MLineN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MLineN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  int numSubEdges = CTX::instance()->mesh.numSubEdges;
-  SPoint3 pnt1, pnt2;
-  pnt(-1. + 2 * (double)num / numSubEdges, 0., 0., pnt1);
-  pnt(-1. + 2 * (double)(num + 1) / numSubEdges, 0., 0, pnt2);
-  x[0] = pnt1.x(); x[1] = pnt2.x();
-  y[0] = pnt1.y(); y[1] = pnt2.y();
-  z[0] = pnt1.z(); z[1] = pnt2.z();
-  n[0] = n[1] = MEdge(_v[0], _v[1]).normal();
+  if (curved) {
+    int numSubEdges = CTX::instance()->mesh.numSubEdges;
+    SPoint3 pnt1, pnt2;
+    pnt(-1. + 2 * (double)num / numSubEdges, 0., 0., pnt1);
+    pnt(-1. + 2 * (double)(num + 1) / numSubEdges, 0., 0, pnt2);
+    x[0] = pnt1.x(); x[1] = pnt2.x();
+    y[0] = pnt1.y(); y[1] = pnt2.y();
+    z[0] = pnt1.z(); z[1] = pnt2.z();
+    n[0] = n[1] = MEdge(_v[0], _v[1]).normal();
+  }
+  else MLine::getEdgeRep(false, num, x, y, z, n);
 }
diff --git a/Geo/MLine.h b/Geo/MLine.h
index a4b9883bce..12c192d20b 100644
--- a/Geo/MLine.h
+++ b/Geo/MLine.h
@@ -49,8 +49,8 @@ class MLine : public MElement {
   }
   virtual int getNumEdges(){ return 1; }
   virtual MEdge getEdge(int num) const{ return MEdge(_v[0], _v[1]); }
-  virtual int getNumEdgesRep(){ return 1; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 1; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     _getEdgeRep(_v[0], _v[1], x, y, z, n);
   }
@@ -61,8 +61,8 @@ class MLine : public MElement {
   }
   virtual int getNumFaces(){ return 0; }
   virtual MFace getFace(int num){ return MFace(); }
-  virtual int getNumFacesRep(){ return 0; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n){}
+  virtual int getNumFacesRep(bool curved){ return 0; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n){}
   virtual int getType() const { return TYPE_LIN; }
   virtual int getTypeForMSH() const { return MSH_LIN_2; }
   virtual int getTypeForUNV() const { return 21; } // linear beam
@@ -133,8 +133,8 @@ class MLine3 : public MLine {
   }
   virtual MVertex *getVertexINP(int num){ return getVertexUNV(num); }
   virtual int getNumEdgeVertices() const { return 1; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
@@ -182,8 +182,8 @@ class MLineN : public MLine {
   virtual MVertex *getVertex(int num){ return num < 2 ? _v[num] : _vs[num - 2]; }
   virtual const MVertex *getVertex(int num) const{ return num < 2 ? _v[num] : _vs[num - 2]; }
   virtual int getNumEdgeVertices() const { return _vs.size(); }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(2 + _vs.size());
diff --git a/Geo/MPoint.h b/Geo/MPoint.h
index 08ccc7fbc1..06e89b0481 100644
--- a/Geo/MPoint.h
+++ b/Geo/MPoint.h
@@ -35,12 +35,12 @@ class MPoint : public MElement {
   virtual void setVertex(int num,  MVertex *v){ _v[0] = v; }
   virtual int getNumEdges(){ return 0; }
   virtual MEdge getEdge(int num) const{ return MEdge(); }
-  virtual int getNumEdgesRep(){ return 0; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n){}
+  virtual int getNumEdgesRep(bool curved){ return 0; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n){}
   virtual int getNumFaces(){ return 0; }
   virtual MFace getFace(int num){ return MFace(); }
-  virtual int getNumFacesRep(){ return 0; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n){}
+  virtual int getNumFacesRep(bool curved){ return 0; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n){}
   virtual int getType() const { return TYPE_PNT; }
   virtual int getTypeForMSH() const { return MSH_PNT; }
   virtual int getTypeForVTK() const { return 1; }
diff --git a/Geo/MPrism.cpp b/Geo/MPrism.cpp
index 66b8f121d6..35cf9d96b2 100644
--- a/Geo/MPrism.cpp
+++ b/Geo/MPrism.cpp
@@ -168,28 +168,31 @@ static void _myGetEdgeRep(MPrism *pri, int num, double *x, double *y, double *z,
 }
 
 
-void MPrism15::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+void MPrism15::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) {
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPrism::getEdgeRep(false, num, x, y, z, n);
 }
 
-void MPrism18::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+void MPrism18::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) {
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPrism::getEdgeRep(false, num, x, y, z, n);
 }
 
-void MPrismN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n) {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+void MPrismN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n) {
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPrism::getEdgeRep(false, num, x, y, z, n);
 }
 
-int MPrism15::getNumEdgesRep() {
-  return 9 * CTX::instance()->mesh.numSubEdges;
+int MPrism15::getNumEdgesRep(bool curved) {
+  return curved ? 9 * CTX::instance()->mesh.numSubEdges : 9;
 }
 
-int MPrism18::getNumEdgesRep() {
-  return 9 * CTX::instance()->mesh.numSubEdges;
+int MPrism18::getNumEdgesRep(bool curved) {
+  return curved ? 9 * CTX::instance()->mesh.numSubEdges : 9;
 }
 
-int MPrismN::getNumEdgesRep() {
-  return 9 * CTX::instance()->mesh.numSubEdges;
+int MPrismN::getNumEdgesRep(bool curved) {
+  return curved ? 9 * CTX::instance()->mesh.numSubEdges : 9;
 }
 
 
@@ -411,31 +414,34 @@ void _myGetFaceRep(MPrism *pri, int num, double *x, double *y, double *z,
   n[2] = n[0];
 }
 
-void MPrism15::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MPrism15::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPrism::getFaceRep(false, num, x, y, z, n);
 }
 
-void MPrism18::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MPrism18::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPrism::getFaceRep(false, num, x, y, z, n);
 }
 
-void MPrismN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MPrismN::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPrism::getFaceRep(false, num, x, y, z, n);
 }
 
-int MPrism15::getNumFacesRep() {
-  return 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
+int MPrism15::getNumFacesRep(bool curved) {
+  return curved ? 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2) : 8;
 }
 
-int MPrism18::getNumFacesRep() {
-  return 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
+int MPrism18::getNumFacesRep(bool curved) {
+  return curved ? 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2) : 8;
 }
 
-int MPrismN::getNumFacesRep() {
-  return 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2);
+int MPrismN::getNumFacesRep(bool curved) {
+  return curved ? 4 * (CTX::instance()->mesh.numSubEdges * CTX::instance()->mesh.numSubEdges * 2) : 8;
 }
 
 namespace {
diff --git a/Geo/MPrism.h b/Geo/MPrism.h
index dc9a8b9bff..3bf01b3e84 100644
--- a/Geo/MPrism.h
+++ b/Geo/MPrism.h
@@ -71,8 +71,8 @@ class MPrism : public MElement {
   {
     return MEdge(_v[edges_prism(num, 0)], _v[edges_prism(num, 1)]);
   }
-  virtual int getNumEdgesRep(){ return 9; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 9; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[9] = {0, 1, 2, 0, 2, 3, 1, 1, 1};
     MEdge e(getEdge(num));
@@ -97,8 +97,8 @@ class MPrism : public MElement {
                    _v[faces_prism(num, 2)],
                    _v[faces_prism(num, 3)]);
   }
-  virtual int getNumFacesRep(){ return 8; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved){ return 8; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[8][3] = {
       {0, 2, 1},
@@ -241,16 +241,16 @@ class MPrism15 : public MPrism {
   }
   virtual MVertex *getVertexINP(int num){ return getVertexBDF(num); }
   virtual int getNumEdgeVertices() const { return 9; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MPrism::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize((num < 2) ? 6 : 8);
@@ -335,16 +335,16 @@ class MPrism18 : public MPrism {
   virtual const MVertex *getVertex(int num) const{ return num < 6 ? _v[num] : _vs[num - 6]; }
   virtual int getNumEdgeVertices() const { return 9; }
   virtual int getNumFaceVertices() const { return 3; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MPrism::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize((num < 2) ? 6 : 9);
@@ -425,8 +425,8 @@ class MPrismN : public MPrism {
     else
       {int n = _order-1; return n * (n * (n+1) / 2);}
   }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(_order+1);
@@ -434,8 +434,8 @@ class MPrismN : public MPrism {
     const int n = _order-1;
     for (int i=0; i<n; i++) v[2+i] = _vs[num*n+i];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const;
   virtual int getTypeForMSH() const {
     switch (_order) {
diff --git a/Geo/MPyramid.cpp b/Geo/MPyramid.cpp
index 6dc0d132b6..1f496faccb 100644
--- a/Geo/MPyramid.cpp
+++ b/Geo/MPyramid.cpp
@@ -58,42 +58,49 @@ const JacobianBasis* MPyramid::getJacobianFuncSpace(int o) const
 
 MPyramidN::~MPyramidN() {}
 
-int MPyramidN::getNumEdgesRep(){ return 8 * CTX::instance()->mesh.numSubEdges; }
+int MPyramidN::getNumEdgesRep(bool curved) {
+  return curved ? 8 * CTX::instance()->mesh.numSubEdges : 8;
+}
 
-void MPyramidN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MPyramidN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  int numSubEdges = CTX::instance()->mesh.numSubEdges;
-  static double pp[5][3] = {{-1,-1,0},{1,-1,0},{1,1,0},{-1,1,0},{0,0,1}};
-  static int ed [8][2] = {{0,1},{0,3},{0,4},{1,2},{1,4},{2,3},{2,4},{3,4}};
-  int iEdge = num / numSubEdges;
-  int iSubEdge = num % numSubEdges;
-
-  int iVertex1 = ed [iEdge][0];
-  int iVertex2 = ed [iEdge][1];
-  double t1 = (double) iSubEdge / (double) numSubEdges;
-  double u1 = pp[iVertex1][0] * (1.-t1) + pp[iVertex2][0] * t1;
-  double v1 = pp[iVertex1][1] * (1.-t1) + pp[iVertex2][1] * t1;
-  double w1 = pp[iVertex1][2] * (1.-t1) + pp[iVertex2][2] * t1;
-
-  double t2 = (double) (iSubEdge+1) / (double) numSubEdges;
-  double u2 = pp[iVertex1][0] * (1.-t2) + pp[iVertex2][0] * t2;
-  double v2 = pp[iVertex1][1] * (1.-t2) + pp[iVertex2][1] * t2;
-  double w2 = pp[iVertex1][2] * (1.-t2) + pp[iVertex2][2] * t2;
-
-  SPoint3 pnt1, pnt2;
-  pnt(u1,v1,w1,pnt1);
-  pnt(u2,v2,w2,pnt2);
-  x[0] = pnt1.x(); x[1] = pnt2.x();
-  y[0] = pnt1.y(); y[1] = pnt2.y();
-  z[0] = pnt1.z(); z[1] = pnt2.z();
-
-  // not great, but better than nothing
-  static const int f[8] = {0, 1, 0, 2, 2, 3, 2, 3};
-  n[0] = n[1] = getFace(f[iEdge]).normal();
+  if (curved) {
+    int numSubEdges = CTX::instance()->mesh.numSubEdges;
+    static double pp[5][3] = {{-1,-1,0},{1,-1,0},{1,1,0},{-1,1,0},{0,0,1}};
+    static int ed [8][2] = {{0,1},{0,3},{0,4},{1,2},{1,4},{2,3},{2,4},{3,4}};
+    int iEdge = num / numSubEdges;
+    int iSubEdge = num % numSubEdges;
+
+    int iVertex1 = ed [iEdge][0];
+    int iVertex2 = ed [iEdge][1];
+    double t1 = (double) iSubEdge / (double) numSubEdges;
+    double u1 = pp[iVertex1][0] * (1.-t1) + pp[iVertex2][0] * t1;
+    double v1 = pp[iVertex1][1] * (1.-t1) + pp[iVertex2][1] * t1;
+    double w1 = pp[iVertex1][2] * (1.-t1) + pp[iVertex2][2] * t1;
+
+    double t2 = (double) (iSubEdge+1) / (double) numSubEdges;
+    double u2 = pp[iVertex1][0] * (1.-t2) + pp[iVertex2][0] * t2;
+    double v2 = pp[iVertex1][1] * (1.-t2) + pp[iVertex2][1] * t2;
+    double w2 = pp[iVertex1][2] * (1.-t2) + pp[iVertex2][2] * t2;
+
+    SPoint3 pnt1, pnt2;
+    pnt(u1,v1,w1,pnt1);
+    pnt(u2,v2,w2,pnt2);
+    x[0] = pnt1.x(); x[1] = pnt2.x();
+    y[0] = pnt1.y(); y[1] = pnt2.y();
+    z[0] = pnt1.z(); z[1] = pnt2.z();
+
+    // not great, but better than nothing
+    static const int f[8] = {0, 1, 0, 2, 2, 3, 2, 3};
+    n[0] = n[1] = getFace(f[iEdge]).normal();
+  }
+  else MPyramid::getEdgeRep(false, num, x, y, z, n);
 }
 
 
-int MPyramidN::getNumFacesRep(){ return 6 * SQU(CTX::instance()->mesh.numSubEdges); }
+int MPyramidN::getNumFacesRep(bool curved) {
+  return curved ? 6 * SQU(CTX::instance()->mesh.numSubEdges) : 6;
+}
 
 static void _myGetFaceRep(MPyramid *pyr, int num, double *x, double *y, double *z,
                           SVector3 *n, int numSubEdges)
@@ -331,7 +338,8 @@ static void _myGetFaceRep(MPyramid *pyr, int num, double *x, double *y, double *
   }
 }
 
-void MPyramidN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MPyramidN::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MPyramid::getFaceRep(false, num, x, y, z, n);
 }
diff --git a/Geo/MPyramid.h b/Geo/MPyramid.h
index b8a4270e75..73c71a9d80 100644
--- a/Geo/MPyramid.h
+++ b/Geo/MPyramid.h
@@ -76,8 +76,8 @@ class MPyramid : public MElement {
   {
     return MEdge(_v[edges_pyramid(num, 0)], _v[edges_pyramid(num, 1)]);
   }
-  virtual int getNumEdgesRep(){ return 8; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 8; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[8] = {0, 1, 1, 2, 0, 3, 2, 3};
     MEdge e(getEdge(num));
@@ -98,8 +98,8 @@ class MPyramid : public MElement {
     else
       return MFace(_v[0], _v[3], _v[2], _v[1]);
   }
-  virtual int getNumFacesRep(){ return 6; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved){ return 6; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[6][3] = {
       {0, 1, 4},
@@ -333,10 +333,10 @@ class MPyramidN : public MPyramid {
     Msg::Error("Reverse not implemented yet for MPyramidN");
 
   }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
-  virtual int getNumEdgesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
-  virtual int getNumFacesRep();
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
   virtual void getNode(int num, double &u, double &v, double &w) const
   {
     num < 5 ? MPyramid::getNode(num, u, v, w) : MElement::getNode(num, u, v, w);
diff --git a/Geo/MQuadrangle.cpp b/Geo/MQuadrangle.cpp
index 8f968870e5..2e921ac533 100644
--- a/Geo/MQuadrangle.cpp
+++ b/Geo/MQuadrangle.cpp
@@ -33,9 +33,17 @@ const JacobianBasis* MQuadrangle::getJacobianFuncSpace(int order) const
   return tag ? BasisFactory::getJacobianBasis(tag) : NULL;
 }
 
-int MQuadrangleN::getNumEdgesRep(){ return 4 * CTX::instance()->mesh.numSubEdges; }
-int MQuadrangle8::getNumEdgesRep(){ return 4 * CTX::instance()->mesh.numSubEdges; }
-int MQuadrangle9::getNumEdgesRep(){ return 4 * CTX::instance()->mesh.numSubEdges; }
+int MQuadrangleN::getNumEdgesRep(bool curved){
+  return curved ? 4 * CTX::instance()->mesh.numSubEdges : 4;
+}
+
+int MQuadrangle8::getNumEdgesRep(bool curved){
+  return curved ? 4 * CTX::instance()->mesh.numSubEdges : 4;
+}
+
+int MQuadrangle9::getNumEdgesRep(bool curved){
+  return curved ? 4 * CTX::instance()->mesh.numSubEdges : 4;
+}
 
 double MQuadrangle::getVolume()
 {
@@ -84,23 +92,32 @@ static void _myGetEdgeRep(MQuadrangle *q, int num, double *x, double *y, double
   z[0] = pnt1.z(); z[1] = pnt2.z();
 }
 
-void MQuadrangleN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MQuadrangleN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MQuadrangle::getEdgeRep(false, num, x, y, z, n);
 }
-void MQuadrangle8::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MQuadrangle8::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MQuadrangle::getEdgeRep(false, num, x, y, z, n);
 }
-void MQuadrangle9::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MQuadrangle9::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MQuadrangle::getEdgeRep(false, num, x, y, z, n);
 }
 
 
-int MQuadrangleN::getNumFacesRep(){ return 2*SQU(CTX::instance()->mesh.numSubEdges); }
-int MQuadrangle8::getNumFacesRep(){ return 2*SQU(CTX::instance()->mesh.numSubEdges); }
-int MQuadrangle9::getNumFacesRep(){ return 2*SQU(CTX::instance()->mesh.numSubEdges); }
+int MQuadrangleN::getNumFacesRep(bool curved){
+  return curved ? 2*SQU(CTX::instance()->mesh.numSubEdges) : 2;
+}
+int MQuadrangle8::getNumFacesRep(bool curved){
+  return curved ? 2*SQU(CTX::instance()->mesh.numSubEdges) : 2;
+}
+int MQuadrangle9::getNumFacesRep(bool curved){
+  return curved ? 2*SQU(CTX::instance()->mesh.numSubEdges) : 2;
+}
 
 static void _myGetFaceRep(MQuadrangle *t, int num, double *x, double *y, double *z,
                           SVector3 *n, int numSubEdges)
@@ -154,17 +171,20 @@ static void _myGetFaceRep(MQuadrangle *t, int num, double *x, double *y, double
   z[0] = pnt1.z(); z[1] = pnt2.z(); z[2] = pnt3.z();
 }
 
-void MQuadrangleN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MQuadrangleN::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MQuadrangle::getFaceRep(false, num, x, y, z, n);
 }
-void MQuadrangle8::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MQuadrangle8::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MQuadrangle::getFaceRep(false, num, x, y, z, n);
 }
-void MQuadrangle9::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MQuadrangle9::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MQuadrangle::getFaceRep(false, num, x, y, z, n);
 }
 
 void MQuadrangle::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
diff --git a/Geo/MQuadrangle.h b/Geo/MQuadrangle.h
index 90fca18b9c..f2e696e0ec 100644
--- a/Geo/MQuadrangle.h
+++ b/Geo/MQuadrangle.h
@@ -82,8 +82,8 @@ class MQuadrangle : public MElement {
     }
     Msg::Error("Could not get edge information for quadranglee %d", getNum());
   }
-  virtual int getNumEdgesRep(){ return 4; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 4; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     MEdge e(getEdge(num));
     _getEdgeRep(e.getVertex(0), e.getVertex(1), x, y, z, n, 0);
@@ -95,8 +95,8 @@ class MQuadrangle : public MElement {
   }
   virtual int getNumFaces(){ return 1; }
   virtual MFace getFace(int num){ return MFace(_v[0], _v[1], _v[2], _v[3]); }
-  virtual int getNumFacesRep(){ return 2; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved){ return 2; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[2][3] = {
       {0, 1, 2}, {0, 2, 3}
@@ -212,16 +212,16 @@ class MQuadrangle8 : public MQuadrangle {
     return getVertex(map[num]);
   }
   virtual int getNumEdgeVertices() const { return 4; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MQuadrangle::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(8);
@@ -295,16 +295,16 @@ class MQuadrangle9 : public MQuadrangle {
   }
   virtual int getNumEdgeVertices() const { return 4; }
   virtual int getNumFaceVertices() const { return 1; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MQuadrangle::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(9);
@@ -381,9 +381,9 @@ class MQuadrangleN : public MQuadrangle {
       return  (_order - 1) * (_order - 1);
   }
   virtual int getNumEdgeVertices() const { return 4 * (_order - 1); }
-  virtual int getNumEdgesRep();
-  virtual int getNumFacesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(_order + 1);
@@ -393,7 +393,7 @@ class MQuadrangleN : public MQuadrangle {
     for(int i = num * (_order-1); i != ie; ++i)
       v[j++] = _vs[i];
   }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(4 + _vs.size());
diff --git a/Geo/MTetrahedron.cpp b/Geo/MTetrahedron.cpp
index 4a542fad8c..680ee8f8bc 100644
--- a/Geo/MTetrahedron.cpp
+++ b/Geo/MTetrahedron.cpp
@@ -116,8 +116,13 @@ const JacobianBasis* MTetrahedron::getJacobianFuncSpace(int order) const
   return tag ? BasisFactory::getJacobianBasis(tag) : NULL;
 }
 
-int MTetrahedron10::getNumEdgesRep(){ return 6 * CTX::instance()->mesh.numSubEdges; }
-int MTetrahedronN::getNumEdgesRep(){ return 6 * CTX::instance()->mesh.numSubEdges; }
+int MTetrahedron10::getNumEdgesRep(bool curved){
+  return curved ? 6 * CTX::instance()->mesh.numSubEdges : 6;
+}
+
+int MTetrahedronN::getNumEdgesRep(bool curved){
+  return curved ? 6 * CTX::instance()->mesh.numSubEdges : 6;
+}
 
 static void _myGetEdgeRep(MTetrahedron *tet, int num, double *x, double *y, double *z,
                           SVector3 *n, int numSubEdges)
@@ -151,18 +156,25 @@ static void _myGetEdgeRep(MTetrahedron *tet, int num, double *x, double *y, doub
   n[0] = n[1] = tet->getFace(f[iEdge]).normal();
 }
 
-void MTetrahedron10::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTetrahedron10::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTetrahedron::getEdgeRep(false, num, x, y, z, n);
 }
 
-void MTetrahedronN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTetrahedronN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTetrahedron::getEdgeRep(false, num, x, y, z, n);
 }
 
-int MTetrahedronN::getNumFacesRep(){ return 4 * SQU(CTX::instance()->mesh.numSubEdges); }
-int MTetrahedron10::getNumFacesRep(){ return 4 * SQU(CTX::instance()->mesh.numSubEdges); }
+int MTetrahedronN::getNumFacesRep(bool curved) {
+  return curved ? 4 * SQU(CTX::instance()->mesh.numSubEdges) : 4;
+}
+
+int MTetrahedron10::getNumFacesRep(bool curved) {
+  return curved ? 4 * SQU(CTX::instance()->mesh.numSubEdges) : 4;
+}
 
 static void _myGetFaceRep(MTetrahedron *tet, int num, double *x, double *y, double *z,
                           SVector3 *n, int numSubEdges)
@@ -243,14 +255,16 @@ static void _myGetFaceRep(MTetrahedron *tet, int num, double *x, double *y, doub
   n[2] = n[0];
 }
 
-void MTetrahedronN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTetrahedronN::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTetrahedron::getFaceRep(false, num, x, y, z, n);
 }
 
-void MTetrahedron10::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTetrahedron10::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTetrahedron::getFaceRep(false, num, x, y, z, n);
 }
 
 void MTetrahedron::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
diff --git a/Geo/MTetrahedron.h b/Geo/MTetrahedron.h
index e8f32c5590..de26df6f3e 100644
--- a/Geo/MTetrahedron.h
+++ b/Geo/MTetrahedron.h
@@ -67,8 +67,8 @@ class MTetrahedron : public MElement {
   {
     return MEdge(_v[edges_tetra(num, 0)], _v[edges_tetra(num, 1)]);
   }
-  virtual int getNumEdgesRep(){ return 6; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 6; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     static const int f[6] = {0, 0, 0, 1, 2, 3};
     MEdge e(getEdge(num));
@@ -87,8 +87,8 @@ class MTetrahedron : public MElement {
                  _v[faces_tetra(num, 2)]);
   }
   virtual void getFaceInfo(const MFace & face, int &ithFace, int &sign, int &rot) const;
-  virtual int getNumFacesRep(){ return 4; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved){ return 4; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     MFace f(getFace(num));
     _getFaceRep(f.getVertex(0), f.getVertex(1), f.getVertex(2), x, y, z, n);
@@ -243,10 +243,10 @@ class MTetrahedron10 : public MTetrahedron {
   virtual MVertex *getVertexDIFF(int num){ return getVertexBDF(num); }
   virtual MVertex *getVertexINP(int num){ return getVertexBDF(num); }
   virtual int getNumEdgeVertices() const { return 6; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
-  virtual int getNumEdgesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
-  virtual int getNumFacesRep();
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
@@ -420,10 +420,10 @@ class MTetrahedronN : public MTetrahedron {
       inv[i] = _vs[reverseIndices[i + 4] - 4];
     _vs = inv;
   }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
-  virtual int getNumEdgesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
-  virtual int getNumFacesRep();
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
   virtual void getNode(int num, double &u, double &v, double &w) const
   {
     num < 4 ? MTetrahedron::getNode(num, u, v, w) : MElement::getNode(num, u, v, w);
diff --git a/Geo/MTriangle.cpp b/Geo/MTriangle.cpp
index 939ab6ccf5..e7f2a8b37f 100644
--- a/Geo/MTriangle.cpp
+++ b/Geo/MTriangle.cpp
@@ -132,8 +132,13 @@ const JacobianBasis* MTriangle::getJacobianFuncSpace(int order) const
   return tag ? BasisFactory::getJacobianBasis(tag) : NULL;
 }
 
-int MTriangleN::getNumEdgesRep(){ return 3 * CTX::instance()->mesh.numSubEdges; }
-int MTriangle6::getNumEdgesRep(){ return 3 * CTX::instance()->mesh.numSubEdges; }
+int MTriangleN::getNumEdgesRep(bool curved) {
+  return curved ? 3 * CTX::instance()->mesh.numSubEdges : 3;
+}
+
+int MTriangle6::getNumEdgesRep(bool curved) {
+  return curved ? 3 * CTX::instance()->mesh.numSubEdges : 3;
+}
 
 static void _myGetEdgeRep(MTriangle *t, int num, double *x, double *y, double *z,
                           SVector3 *n, int numSubEdges)
@@ -170,18 +175,24 @@ static void _myGetEdgeRep(MTriangle *t, int num, double *x, double *y, double *z
   }
 }
 
-void MTriangleN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTriangleN::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTriangle::getEdgeRep(false, num, x, y, z, n);
 }
 
-void MTriangle6::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTriangle6::getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetEdgeRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTriangle::getEdgeRep(false, num, x, y, z, n);
 }
 
-int MTriangle6::getNumFacesRep(){ return SQU(CTX::instance()->mesh.numSubEdges); }
-int MTriangleN::getNumFacesRep(){ return SQU(CTX::instance()->mesh.numSubEdges); }
+int MTriangle6::getNumFacesRep(bool curved) {
+  return curved ? SQU(CTX::instance()->mesh.numSubEdges) : 1;
+}
+int MTriangleN::getNumFacesRep(bool curved) {
+  return curved ? SQU(CTX::instance()->mesh.numSubEdges) : 1;
+}
 
 static void _myGetFaceRep(MTriangle *t, int num, double *x, double *y, double *z,
                           SVector3 *n, int numSubEdges)
@@ -246,13 +257,15 @@ static void _myGetFaceRep(MTriangle *t, int num, double *x, double *y, double *z
   z[0] = pnt1.z(); z[1] = pnt2.z(); z[2] = pnt3.z();
 }
 
-void MTriangleN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTriangleN::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTriangle::getFaceRep(false, num, x, y, z, n);
 }
-void MTriangle6::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+void MTriangle6::getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
 {
-  _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  if (curved) _myGetFaceRep(this, num, x, y, z, n, CTX::instance()->mesh.numSubEdges);
+  else MTriangle::getFaceRep(false, num, x, y, z, n);
 }
 
 void MTriangle::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
diff --git a/Geo/MTriangle.h b/Geo/MTriangle.h
index d4eac411c3..db9c31c19b 100644
--- a/Geo/MTriangle.h
+++ b/Geo/MTriangle.h
@@ -86,8 +86,8 @@ class MTriangle : public MElement {
     }
     Msg::Error("Could not get edge information for triangle %d", getNum());
   }
-  virtual int getNumEdgesRep(){ return 3; }
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumEdgesRep(bool curved){ return 3; }
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     MEdge e(getEdge(num));
     _getEdgeRep(e.getVertex(0), e.getVertex(1), x, y, z, n, 0);
@@ -102,8 +102,8 @@ class MTriangle : public MElement {
   {
     return MFace(_v[0], _v[1], _v[2]);
   }
-  virtual int getNumFacesRep(){ return 1; }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+  virtual int getNumFacesRep(bool curved){ return 1; }
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n)
   {
     _getFaceRep(_v[0], _v[1], _v[2], x, y, z, n);
   }
@@ -202,16 +202,16 @@ class MTriangle6 : public MTriangle {
   }
   virtual void xyz2uvw(double xyz[3], double uvw[3]) const{ MElement::xyz2uvw(xyz, uvw); }
   virtual int getNumEdgeVertices() const { return 3; }
-  virtual int getNumEdgesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3);
     MTriangle::_getEdgeVertices(num, v);
     v[2] = _vs[num];
   }
-  virtual int getNumFacesRep();
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(6);
@@ -286,9 +286,9 @@ class MTriangleN : public MTriangle {
   }
   virtual void xyz2uvw(double xyz[3], double uvw[3]) const { MElement::xyz2uvw(xyz, uvw); }
   virtual int getNumEdgeVertices() const { return 3 * (_order - 1); }
-  virtual int getNumEdgesRep();
-  virtual int getNumFacesRep();
-  virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual int getNumEdgesRep(bool curved);
+  virtual int getNumFacesRep(bool curved);
+  virtual void getEdgeRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(_order + 1);
@@ -297,7 +297,7 @@ class MTriangleN : public MTriangle {
     const int ie = (num + 1) * (_order - 1);
     for(int i = num * (_order-1); i != ie; ++i) v[j++] = _vs[i];
   }
-  virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
+  virtual void getFaceRep(bool curved, int num, double *x, double *y, double *z, SVector3 *n);
   virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
   {
     v.resize(3 + _vs.size());
-- 
GitLab