diff --git a/Common/GmshDefines.h b/Common/GmshDefines.h
index adc4bae0686bb4acfd4dac6c7b2b693423fe1169..b954314d4d5d13465f3a3aa5f45e93638207069b 100644
--- a/Common/GmshDefines.h
+++ b/Common/GmshDefines.h
@@ -142,6 +142,12 @@
#define MSH_TET_164 82
#define MSH_TET_202 83
#define MSH_NUM_TYPE 84
+#define MSH_LIN_1 85
+#define MSH_TRI_1 86
+#define MSH_QUA_1 87
+#define MSH_TET_1 88
+#define MSH_HEX_1 89
+#define MSH_PRI_1 90
// Geometric entities
#define ENT_NONE 0
diff --git a/Geo/MLine.cpp b/Geo/MLine.cpp
index 4771a12613fd9a27bea33fa5972bd3fdf60a32f0..66e84489c29464304f97e2ddb5acd3518f715669 100644
--- a/Geo/MLine.cpp
+++ b/Geo/MLine.cpp
@@ -14,6 +14,7 @@ const polynomialBasis* MLine::getFunctionSpace(int o) const
int order = (o == -1) ? getPolynomialOrder() : o;
switch (order) {
+ case 0: return polynomialBases::find(MSH_LIN_1);
case 1: return polynomialBases::find(MSH_LIN_2);
case 2: return polynomialBases::find(MSH_LIN_3);
case 3: return polynomialBases::find(MSH_LIN_4);
diff --git a/Geo/MPrism.cpp b/Geo/MPrism.cpp
index 22586f844df6f32916430fca77f0eb810bd185e2..ea385278fb5ed115ecaee1339f81bf4e3a5c6cad 100644
--- a/Geo/MPrism.cpp
+++ b/Geo/MPrism.cpp
@@ -38,6 +38,7 @@ const polynomialBasis* MPrism::getFunctionSpace(int o) const
if ((nv == 0) && (o == -1)) {
switch (order) {
+ case 0: return polynomialBases::find(MSH_PRI_1);
case 1: return polynomialBases::find(MSH_PRI_6);
case 2: return polynomialBases::find(MSH_PRI_18);
default: Msg::Error("Order %d prism function space not implemented", order);
@@ -45,6 +46,7 @@ const polynomialBasis* MPrism::getFunctionSpace(int o) const
}
else {
switch (order) {
+ case 0: return polynomialBases::find(MSH_PRI_1);
case 1: return polynomialBases::find(MSH_PRI_6);
case 2: return polynomialBases::find(MSH_PRI_18);
default: Msg::Error("Order %d prism function space not implemented", order);
diff --git a/Geo/MQuadrangle.cpp b/Geo/MQuadrangle.cpp
index 13eb65e2cc4487ec3d87ed35b64dbcf6689beb1c..0563ae6d29f0fecd4fd95368bedf569f771d0258 100644
--- a/Geo/MQuadrangle.cpp
+++ b/Geo/MQuadrangle.cpp
@@ -24,6 +24,7 @@ const polynomialBasis* MQuadrangle::getFunctionSpace(int o) const
if ((nf == 0) && (o == -1)) {
switch (order) {
+ case 0: return polynomialBases::find(MSH_QUA_1);
case 1: return polynomialBases::find(MSH_QUA_4);
case 2: return polynomialBases::find(MSH_QUA_8);
case 3: return polynomialBases::find(MSH_QUA_12);
@@ -37,6 +38,7 @@ const polynomialBasis* MQuadrangle::getFunctionSpace(int o) const
}
}
switch (order) {
+ case 0: return polynomialBases::find(MSH_QUA_1);
case 1: return polynomialBases::find(MSH_QUA_4);
case 2: return polynomialBases::find(MSH_QUA_9);
case 3: return polynomialBases::find(MSH_QUA_16);
diff --git a/Geo/MTetrahedron.cpp b/Geo/MTetrahedron.cpp
index aea4336b8d0f60d29fdbbdbca6119483fa19f8de..83785d50e3a53ea0d4c8e754cadf94af27165ab6 100644
--- a/Geo/MTetrahedron.cpp
+++ b/Geo/MTetrahedron.cpp
@@ -125,6 +125,7 @@ const polynomialBasis* MTetrahedron::getFunctionSpace(int o) const
if ((nv == 0) && (o == -1)) {
switch (order) {
+ case 0: return polynomialBases::find(MSH_TET_1);
case 1: return polynomialBases::find(MSH_TET_4);
case 2: return polynomialBases::find(MSH_TET_10);
case 3: return polynomialBases::find(MSH_TET_20);
@@ -140,6 +141,7 @@ const polynomialBasis* MTetrahedron::getFunctionSpace(int o) const
}
else {
switch (order) {
+ case 0: return polynomialBases::find(MSH_TET_1);
case 1: return polynomialBases::find(MSH_TET_4);
case 2: return polynomialBases::find(MSH_TET_10);
case 3: return polynomialBases::find(MSH_TET_20);
diff --git a/Geo/MTriangle.cpp b/Geo/MTriangle.cpp
index 721fa4de2268d4eae79e250f639abf4311cddfd1..0c4f6322bad7be7d23a6c2145752a039f632f3a2 100644
--- a/Geo/MTriangle.cpp
+++ b/Geo/MTriangle.cpp
@@ -84,6 +84,7 @@ const polynomialBasis* MTriangle::getFunctionSpace(int o) const
if ((nf == 0) && (o == -1)) {
switch (order) {
+ case 0: return polynomialBases::find(MSH_TRI_1);
case 1: return polynomialBases::find(MSH_TRI_3);
case 2: return polynomialBases::find(MSH_TRI_6);
case 3: return polynomialBases::find(MSH_TRI_9);
@@ -99,6 +100,7 @@ const polynomialBasis* MTriangle::getFunctionSpace(int o) const
}
else {
switch (order) {
+ case 0: return polynomialBases::find(MSH_TRI_1);
case 1: return polynomialBases::find(MSH_TRI_3);
case 2: return polynomialBases::find(MSH_TRI_6);
case 3: return polynomialBases::find(MSH_TRI_10);
@@ -117,12 +119,14 @@ const polynomialBasis* MTriangle::getFunctionSpace(int o) const
const JacobianBasis* MTriangle::getJacobianFuncSpace(int o) const
{
+
int order = (o == -1) ? getPolynomialOrder() : o;
int nf = getNumFaceVertices();
if ((nf == 0) && (o == -1)) {
switch (order) {
+ case 0: return JacobianBases::find(MSH_TRI_1);
case 1: return JacobianBases::find(MSH_TRI_3);
case 2: return JacobianBases::find(MSH_TRI_6);
case 3: return JacobianBases::find(MSH_TRI_9);
diff --git a/Numeric/polynomialBasis.cpp b/Numeric/polynomialBasis.cpp
index 79d31aca83abeb8d822118abd92556f42adb1065..da15f93cde44b302de0595cfd0882b05de1fac0c 100644
--- a/Numeric/polynomialBasis.cpp
+++ b/Numeric/polynomialBasis.cpp
@@ -14,6 +14,7 @@
static int getTriangleType (int order) {
switch(order) {
+ case 0 : return MSH_TRI_1;
case 1 : return MSH_TRI_3;
case 2 : return MSH_TRI_6;
case 3 : return MSH_TRI_10;
@@ -29,6 +30,7 @@ static int getTriangleType (int order) {
}
static int getQuadType (int order) {
switch(order) {
+ case 0 : return MSH_QUA_1;
case 1 : return MSH_QUA_4;
case 2 : return MSH_QUA_9;
case 3 : return MSH_QUA_16;
@@ -44,6 +46,7 @@ static int getQuadType (int order) {
}
static int getLineType (int order) {
switch(order) {
+ case 0 : return MSH_LIN_1;
case 1 : return MSH_LIN_2;
case 2 : return MSH_LIN_3;
case 3 : return MSH_LIN_4;
@@ -1018,6 +1021,16 @@ static void generate1dVertexClosure(polynomialBasis::clCont &closure)
closure[1].type = MSH_PNT;
}
+static void generateClosureOrder0(polynomialBasis::clCont &closure, int nb)
+{
+ closure.clear();
+ closure.resize(nb);
+ for (int i=0; i < nb; i++) {
+ closure[i].push_back(0);
+ closure[i].type = MSH_PNT;
+ }
+}
+
std::map<int, polynomialBasis> polynomialBases::fs;
const polynomialBasis *polynomialBases::find(int tag)
@@ -1034,6 +1047,13 @@ const polynomialBasis *polynomialBases::find(int tag)
F.points = generate1DPoints(0);
F.parentType = TYPE_PNT;
break;
+ case MSH_LIN_1 :
+ F.numFaces = 2;
+ F.monomials = generate1DMonomials(0);
+ F.points = generate1DPoints(0);
+ generateClosureOrder0(F.closures,2);
+ F.parentType = TYPE_LIN;
+ break;
case MSH_LIN_2 :
F.numFaces = 2;
F.monomials = generate1DMonomials(1);
@@ -1104,6 +1124,13 @@ const polynomialBasis *polynomialBases::find(int tag)
F.parentType = TYPE_LIN;
generate1dVertexClosure(F.closures);
break;
+ case MSH_TRI_1 :
+ F.numFaces = 3;
+ F.monomials = generatePascalTriangle(0);
+ F.points = gmshGeneratePointsTriangle(0, false);
+ F.parentType = TYPE_TRI;
+ generateClosureOrder0(F.closures,6);
+ break;
case MSH_TRI_3 :
F.numFaces = 3;
F.monomials = generatePascalTriangle(1);
@@ -1230,6 +1257,13 @@ const polynomialBasis *polynomialBases::find(int tag)
F.parentType = TYPE_TRI;
generate2dEdgeClosure(F.closures, 10);
break;
+ case MSH_TET_1 :
+ F.numFaces = 4;
+ F.monomials = generatePascalTetrahedron(0);
+ F.points = gmshGeneratePointsTetrahedron(0, false);
+ F.parentType = TYPE_TET;
+ generateClosureOrder0(F.closures,24);
+ break;
case MSH_TET_4 :
F.numFaces = 4;
F.monomials = generatePascalTetrahedron(1);
@@ -1349,6 +1383,13 @@ const polynomialBasis *polynomialBases::find(int tag)
F.parentType = TYPE_TET;
generateFaceClosureTet(F.closures, 10);
break;
+ case MSH_QUA_1 :
+ F.numFaces = 4;
+ F.monomials = generatePascalQuad(0);
+ F.points = gmshGeneratePointsQuad(0,false);
+ F.parentType = TYPE_QUA;
+ generateClosureOrder0(F.closures,8);
+ break;
case MSH_QUA_4 :
F.numFaces = 4;
F.monomials = generatePascalQuad(1);
@@ -1482,6 +1523,13 @@ const polynomialBasis *polynomialBases::find(int tag)
F.parentType = TYPE_QUA;
generate2dEdgeClosure(F.closures, 10, 4);
break;
+ case MSH_PRI_1 :
+ F.numFaces = 5;
+ F.monomials = generatePascalPrism(0);
+ F.points = gmshGeneratePointsPrism(1, false);
+ F.parentType = TYPE_PRI;
+ generateClosureOrder0(F.closures,48);
+ break;
case MSH_PRI_6 :
F.numFaces = 5;
F.monomials = generatePascalPrism(1);