diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 125af1ae6d66d4cc6dc37ccf726c39cfc22dadaa..e500de751cb2fbff1b0439d5dbb00f9622b2240d 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -719,7 +719,72 @@ int MElement::getInfoMSH(const int typeMSH, const char **const name)
if(name) *name = "Unknown";
return 0;
}
-}
+}
+
+void MTriangle::getShapeFunction(int num,double uu, double vv,double ww,double& s) {
+
+#if !defined(HAVE_GMSH_EMBEDDED)
+
+ int nf = getNumFaceVertices();
+ double fv[256];
+
+ if (!nf){
+ switch(getPolynomialOrder()){
+ case 1: gmshFunctionSpaces::find(MSH_TRI_3).f(uu, vv, fv); break;
+ case 2: gmshFunctionSpaces::find(MSH_TRI_6).f(uu, vv, fv); break;
+ case 3: gmshFunctionSpaces::find(MSH_TRI_9).f(uu, vv, fv); break;
+ case 4: gmshFunctionSpaces::find(MSH_TRI_12).f(uu, vv, fv); break;
+ case 5: gmshFunctionSpaces::find(MSH_TRI_15I).f(uu, vv, fv); break;
+ default: Msg::Error("Order %d triangle shape function not implemented", getPolynomialOrder()); break;
+ }
+ }
+ else{
+ switch(getPolynomialOrder()){
+ case 1: gmshFunctionSpaces::find(MSH_TRI_3).f(uu, vv, fv); break;
+ case 2: gmshFunctionSpaces::find(MSH_TRI_6).f(uu, vv, fv); break;
+ case 3: gmshFunctionSpaces::find(MSH_TRI_10).f(uu, vv, fv); break;
+ case 4: gmshFunctionSpaces::find(MSH_TRI_15).f(uu, vv, fv); break;
+ case 5: gmshFunctionSpaces::find(MSH_TRI_21).f(uu, vv, fv); break;
+ default: Msg::Error("Order %d triangle jac not implemented", getPolynomialOrder()); break;
+ }
+ }
+ s = fv[num];
+#endif
+
+}
+
+void MTriangle::getGradShapeFunction(int num,double uu,double vv,double ww,double s[3]) {
+
+
+#if !defined(HAVE_GMSH_EMBEDDED)
+ double grads[256][2];
+ int nf = getNumFaceVertices();
+
+ if (!nf){
+ switch(getPolynomialOrder()){
+ case 1: gmshFunctionSpaces::find(MSH_TRI_3).df(uu, vv, ww, grads); break;
+ case 2: gmshFunctionSpaces::find(MSH_TRI_6).df(uu, vv, ww, grads); break;
+ case 3: gmshFunctionSpaces::find(MSH_TRI_9).df(uu, vv, ww, grads); break;
+ case 4: gmshFunctionSpaces::find(MSH_TRI_12).df(uu, vv, ww, grads); break;
+ case 5: gmshFunctionSpaces::find(MSH_TRI_15I).df(uu, vv, ww, grads); break;
+ default: Msg::Error("Order %d triangle jac not implemented", getPolynomialOrder()); break;
+ }
+ }
+ else{
+ switch(getPolynomialOrder()){
+ case 1: gmshFunctionSpaces::find(MSH_TRI_3).df(uu, vv, ww, grads); break;
+ case 2: gmshFunctionSpaces::find(MSH_TRI_6).df(uu, vv, ww, grads); break;
+ case 3: gmshFunctionSpaces::find(MSH_TRI_10).df(uu, vv, ww, grads); break;
+ case 4: gmshFunctionSpaces::find(MSH_TRI_15).df(uu, vv, ww, grads); break;
+ case 5: gmshFunctionSpaces::find(MSH_TRI_21).df(uu, vv, ww, grads); break;
+ default: Msg::Error("Order %d triangle jac not implemented", getPolynomialOrder()); break;
+ }
+ }
+
+ for (int i=0;i<3;i++) s[i] = grads[num][i];
+
+#endif
+}
void MTriangle::jac(int ord, MVertex *vs[], double uu, double vv, double ww,
double j[2][3])
diff --git a/Geo/MElement.h b/Geo/MElement.h
index fc96c6362864157706d2a99318498af66165563f..806cba53e7316e0a81f50fa03655a9c802db4a43 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -566,24 +566,28 @@ class MTriangle : public MElement {
{
MVertex *tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
}
- virtual void getShapeFunction(int num, double u, double v, double w, double &s)
- {
- switch(num){
- case 0 : s = 1. - u - v; break;
- case 1 : s = u ; break;
- case 2 : s = v; break;
- default : s = 0.; break;
- }
- }
- virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3])
- {
- switch(num) {
- case 0 : s[0] = -1.; s[1] = -1.; s[2] = 0.; break;
- case 1 : s[0] = 1.; s[1] = 0.; s[2] = 0.; break;
- case 2 : s[0] = 0.; s[1] = 1.; s[2] = 0.; break;
- default : s[0] = s[1] = s[2] = 0.; break;
- }
- }
+
+ virtual void getShapeFunction(int num, double u, double v, double w, double &s);
+ virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]);
+
+/* virtual void getShapeFunction(int num, double u, double v, double w, double &s) */
+/* { */
+/* switch(num){ */
+/* case 0 : s = 1. - u - v; break; */
+/* case 1 : s = u ; break; */
+/* case 2 : s = v; break; */
+/* default : s = 0.; break; */
+/* } */
+/* } */
+/* virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]) */
+/* { */
+/* switch(num) { */
+/* case 0 : s[0] = -1.; s[1] = -1.; s[2] = 0.; break; */
+/* case 1 : s[0] = 1.; s[1] = 0.; s[2] = 0.; break; */
+/* case 2 : s[0] = 0.; s[1] = 1.; s[2] = 0.; break; */
+/* default : s[0] = s[1] = s[2] = 0.; break; */
+/* } */
+/* } */
virtual bool isInside(double u, double v, double w, double tol=1.e-8)
{
if(u < (-tol) || v < (-tol) || u > ((1. + tol) - v))
@@ -704,6 +708,7 @@ class MTriangle6 : public MTriangle {
tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
tmp = _vs[0]; _vs[0] = _vs[2]; _vs[2] = tmp;
}
+
virtual void jac(double u, double v, double w, double j[2][3])
{
MTriangle::jac(2, _vs, u, v, w, j);
diff --git a/Parser/Gmsh.tab.cpp b/Parser/Gmsh.tab.cpp
index 3d2a661fe48223ec6b6cda5500c0ea6f2feaa541..9a22a17e6bdcf22d8a18336d2f8483a1cec83e70 100644
--- a/Parser/Gmsh.tab.cpp
+++ b/Parser/Gmsh.tab.cpp
@@ -376,7 +376,7 @@ static PViewDataList *ViewData;
static ExtrudeParams extr;
static gmshSurface *myGmshSurface = 0;
static List_T *ViewValueList = 0;
-static double ViewCoord[100];
+static double ViewCoord[105];
static int *ViewNumList = 0, ViewCoordIdx = 0;
#define MAX_RECUR_LOOPS 100
static int ImbricatedLoop = 0;
diff --git a/Parser/Gmsh.y b/Parser/Gmsh.y
index cb013f861f072c55230095834102b8790e12629c..3806432de133b25b07f7a95cbd1fcec21792ec5a 100644
--- a/Parser/Gmsh.y
+++ b/Parser/Gmsh.y
@@ -51,7 +51,7 @@ static PViewDataList *ViewData;
static ExtrudeParams extr;
static gmshSurface *myGmshSurface = 0;
static List_T *ViewValueList = 0;
-static double ViewCoord[100];
+static double ViewCoord[105]; // KH: support up to order 4 mappings
static int *ViewNumList = 0, ViewCoordIdx = 0;
#define MAX_RECUR_LOOPS 100
static int ImbricatedLoop = 0;
@@ -453,7 +453,7 @@ Element :
#if !defined(HAVE_NO_POST)
if(ViewValueList){
for(int i = 0; i < 3; i++)
- for(int j = 0; j < ViewCoordIdx / 3; j++)
+ for(int j = 0; j < ViewCoordIdx / 3; j++)
List_Add(ViewValueList, &ViewCoord[3 * j + i]);
}
#endif