Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 15495)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 15496)
@@ -1011,23 +1011,8 @@
 	_assert_(sidlist);
 	_assert_(nodes);
-	int numnodes = this->GetNumberOfNodes();
+	int numnodes = this->NumberOfNodes();
 
 	for(int i=0;i<numnodes;i++){
 		sidlist[i]=nodes[i]->Sid();
-	}
-}
-/*}}}*/
-/*FUNCTION Penta::GetNumberOfNodes{{{*/
-int Penta::GetNumberOfNodes(void){
-
-	if(this->nodes==NULL) return 0;
-
-	switch(this->element_type){
-		case P1Enum:
-			return 6;
-		case P1DGEnum:
-			return 6;
-		default:
-			_error_("Element type "<<EnumToStringx(this->element_type)<<" not supported yet");
 	}
 }
@@ -2632,5 +2617,5 @@
 
 	/*Intermediaries*/
-	const int numnodes = this->GetNumberOfNodes();
+	const int numnodes = this->NumberOfNodes();
 
 	/*Output */
Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 15495)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 15496)
@@ -89,5 +89,4 @@
 		int    GetNodeIndex(Node* node);
 		void   GetNodesSidList(int* sidlist);
-		int    GetNumberOfNodes(void);
 		void   GetSolutionFromInputs(Vector<IssmDouble>* solution);
 		IssmDouble GetZcoord(GaussPenta* gauss);
Index: /issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp	(revision 15495)
+++ /issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp	(revision 15496)
@@ -18,5 +18,5 @@
 #define NUMNODESP1    6
 #define NUMNODESP1_2d 3
-#define NUMNODESMINI  7
+#define NUMNODESP1b  7
 
 /*Object constructors and destructor*/
@@ -98,5 +98,5 @@
 
 	int i;
-	IssmDouble dbasismini[3][NUMNODESMINI];
+	IssmDouble dbasismini[3][NUMNODESP1b];
 	IssmDouble basis[NUMNODESP1];
 
@@ -106,5 +106,5 @@
 
 	/*Build B: */
-	for(i=0;i<NUMNODESMINI;i++){
+	for(i=0;i<NUMNODESP1b;i++){
 		B[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+0] = dbasismini[0][i];
 		B[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+1] = 0.;
@@ -221,5 +221,5 @@
 
 	int    i;
-	IssmDouble dbasismini[3][NUMNODESMINI];
+	IssmDouble dbasismini[3][NUMNODESP1b];
 
 	/*Get dbasis in actual coordinate system: */
@@ -227,5 +227,5 @@
 
 	/*Build Bprime: */
-	for(i=0;i<NUMNODESMINI;i++){
+	for(i=0;i<NUMNODESP1b;i++){
 		Bprime[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+0] = 2.*dbasismini[0][i];
 		Bprime[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+1] = dbasismini[1][i];
@@ -266,5 +266,5 @@
 	int i;
 
-	IssmDouble dbasismini[3][NUMNODESMINI];
+	IssmDouble dbasismini[3][NUMNODESP1b];
 	IssmDouble basis[NUMNODESP1];
 
@@ -274,5 +274,5 @@
 
 	/*Build B: */
-	for(i=0;i<NUMNODESMINI;i++){
+	for(i=0;i<NUMNODESP1b;i++){
 		B[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+0] = dbasismini[0][i+0]; //B[0][NDOF4*i+0] = dbasis[0][i+0];
 		B[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+1] = 0.;
@@ -399,5 +399,5 @@
 
 	int i;
-	IssmDouble dbasismini[3][NUMNODESMINI];
+	IssmDouble dbasismini[3][NUMNODESP1b];
 	IssmDouble basis[NUMNODESP1];
 
@@ -407,5 +407,5 @@
 
 	/*B_primeuild B_prime: */
-	for(i=0;i<NUMNODESMINI;i++){
+	for(i=0;i<NUMNODESP1b;i++){
 		B_prime[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+0] = dbasismini[0][i];
 		B_prime[(NDOF4*NUMNODESP1+3)*0+NDOF4*i+1] = 0.;
@@ -1198,56 +1198,57 @@
 		case P1Enum: case P1DGEnum:
 			/*Nodal function 1*/
-			dbasis[NUMNODESP1*0+0]=-0.5*(1.0-zeta)/2.0;
-			dbasis[NUMNODESP1*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
-			dbasis[NUMNODESP1*2+0]=-0.5*gauss->coord1;
+			dbasis[NUMNODESP1*0+0]   = (zeta-1.)/4.;
+			dbasis[NUMNODESP1*1+0]   = SQRT3/12.*(zeta-1.);
+			dbasis[NUMNODESP1*2+0]   = -.5*gauss->coord1;
 			/*Nodal function 2*/
-			dbasis[NUMNODESP1*0+1]=0.5*(1.0-zeta)/2.0;
-			dbasis[NUMNODESP1*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
-			dbasis[NUMNODESP1*2+1]=-0.5*gauss->coord2;
+			dbasis[NUMNODESP1*0+1]   = (1.-zeta)/4.;
+			dbasis[NUMNODESP1*1+1]   = SQRT3/12.*(zeta-1);
+			dbasis[NUMNODESP1*2+1]   = -.5*gauss->coord2;
 			/*Nodal function 3*/
-			dbasis[NUMNODESP1*0+2]=0.;
-			dbasis[NUMNODESP1*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
-			dbasis[NUMNODESP1*2+2]=-0.5*gauss->coord3;
+			dbasis[NUMNODESP1*0+2]   = 0.;
+			dbasis[NUMNODESP1*1+2]   = SQRT3/6.*(1.-zeta);
+			dbasis[NUMNODESP1*2+2]   = -.5*gauss->coord3;
 			/*Nodal function 4*/
-			dbasis[NUMNODESP1*0+3]=-0.5*(1.0+zeta)/2.0;
-			dbasis[NUMNODESP1*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
-			dbasis[NUMNODESP1*2+3]=0.5*gauss->coord1;
+			dbasis[NUMNODESP1*0+3]   = -(1.+zeta)/4.;
+			dbasis[NUMNODESP1*1+3]   = -SQRT3/12.*(1.+zeta);
+			dbasis[NUMNODESP1*2+3]   = .5*gauss->coord1;
 			/*Nodal function 5*/
-			dbasis[NUMNODESP1*0+4]=0.5*(1.0+zeta)/2.0;
-			dbasis[NUMNODESP1*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
-			dbasis[NUMNODESP1*2+4]=0.5*gauss->coord2;
+			dbasis[NUMNODESP1*0+4]   = (1.+zeta)/4.;
+			dbasis[NUMNODESP1*1+4]   = -SQRT3/12.*(1.+zeta);
+			dbasis[NUMNODESP1*2+4]   = .5*gauss->coord2;
 			/*Nodal function 6*/
-			dbasis[NUMNODESP1*0+5]=0.;
-			dbasis[NUMNODESP1*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
-			dbasis[NUMNODESP1*2+5]=0.5*gauss->coord3;
-		case MINIEnum:
+			dbasis[NUMNODESP1*0+5]   = 0.;
+			dbasis[NUMNODESP1*1+5]   = SQRT3/6.*(1.+zeta);
+			dbasis[NUMNODESP1*2+5]   = .5*gauss->coord3;
+			return;
+		case P1bubbleEnum:
 			/*Nodal function 1*/
-			dbasis[NUMNODESMINI*0+0]=-0.5*(1.0-zeta)/2.0;
-			dbasis[NUMNODESMINI*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
-			dbasis[NUMNODESMINI*2+0]=-0.5*gauss->coord1;
+			dbasis[NUMNODESP1*0+0]   = (zeta-1.)/4.;
+			dbasis[NUMNODESP1*1+0]   = SQRT3/12.*(zeta-1.);
+			dbasis[NUMNODESP1*2+0]   = -.5*gauss->coord1;
 			/*Nodal function 2*/
-			dbasis[NUMNODESMINI*0+1]=0.5*(1.0-zeta)/2.0;
-			dbasis[NUMNODESMINI*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
-			dbasis[NUMNODESMINI*2+1]=-0.5*gauss->coord2;
+			dbasis[NUMNODESP1*0+1]   = (1.-zeta)/4.;
+			dbasis[NUMNODESP1*1+1]   = SQRT3/12.*(zeta-1);
+			dbasis[NUMNODESP1*2+1]   = -.5*gauss->coord2;
 			/*Nodal function 3*/
-			dbasis[NUMNODESMINI*0+2]=0.;
-			dbasis[NUMNODESMINI*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
-			dbasis[NUMNODESMINI*2+2]=-0.5*gauss->coord3;
+			dbasis[NUMNODESP1*0+2]   = 0.;
+			dbasis[NUMNODESP1*1+2]   = SQRT3/6.*(1.-zeta);
+			dbasis[NUMNODESP1*2+2]   = -.5*gauss->coord3;
 			/*Nodal function 4*/
-			dbasis[NUMNODESMINI*0+3]=-0.5*(1.0+zeta)/2.0;
-			dbasis[NUMNODESMINI*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
-			dbasis[NUMNODESMINI*2+3]=0.5*gauss->coord1;
+			dbasis[NUMNODESP1*0+3]   = -(1.+zeta)/4.;
+			dbasis[NUMNODESP1*1+3]   = -SQRT3/12.*(1.+zeta);
+			dbasis[NUMNODESP1*2+3]   = .5*gauss->coord1;
 			/*Nodal function 5*/
-			dbasis[NUMNODESMINI*0+4]=0.5*(1.0+zeta)/2.0;
-			dbasis[NUMNODESMINI*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
-			dbasis[NUMNODESMINI*2+4]=0.5*gauss->coord2;
+			dbasis[NUMNODESP1*0+4]   = (1.+zeta)/4.;
+			dbasis[NUMNODESP1*1+4]   = -SQRT3/12.*(1.+zeta);
+			dbasis[NUMNODESP1*2+4]   = .5*gauss->coord2;
 			/*Nodal function 6*/
-			dbasis[NUMNODESMINI*0+5]=0.;
-			dbasis[NUMNODESMINI*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
-			dbasis[NUMNODESMINI*2+5]=0.5*gauss->coord3;
+			dbasis[NUMNODESP1*0+5]   = 0.;
+			dbasis[NUMNODESP1*1+5]   = SQRT3/6.*(1.+zeta);
+			dbasis[NUMNODESP1*2+5]   = .5*gauss->coord3;
 			/*Nodal function 7*/
-			dbasis[NUMNODESMINI*0+6]=27.*(1.+zeta)*(1.-zeta)*(-.5*gauss->coord2*gauss->coord3 + .5*gauss->coord1*gauss->coord3);
-			dbasis[NUMNODESMINI*1+6]=27.*(1.+zeta)*(1.-zeta)*SQRT3*(-1./6.*gauss->coord2*gauss->coord3 - 1./6.*gauss->coord1*gauss->coord3 +1./3.*gauss->coord1*gauss->coord2);
-			dbasis[NUMNODESMINI*2+6]=27*gauss->coord1*gauss->coord2*gauss->coord3*(-2.0*zeta);
+			dbasis[NUMNODESP1b*0+6] = 27.*(1.+zeta)*(1.-zeta)*(-.5*gauss->coord2*gauss->coord3 + .5*gauss->coord1*gauss->coord3);
+			dbasis[NUMNODESP1b*1+6] = 27.*(1.+zeta)*(1.-zeta)*SQRT3*(-1./6.*gauss->coord2*gauss->coord3 - 1./6.*gauss->coord1*gauss->coord3 +1./3.*gauss->coord1*gauss->coord2);
+			dbasis[NUMNODESP1b*2+6] = -54*gauss->coord1*gauss->coord2*gauss->coord3*zeta;
 			return;
 		default:
@@ -1278,5 +1279,5 @@
 	 * actual coordinate system): */
 
-	IssmDouble    dbasismini_ref[3][NUMNODESMINI];
+	IssmDouble    dbasismini_ref[3][NUMNODESP1b];
 	IssmDouble    Jinv[3][3];
 
@@ -1294,8 +1295,8 @@
 	 */
 
-	for(int i=0;i<NUMNODESMINI;i++){
-		*(dbasismini+NUMNODESMINI*0+i)=Jinv[0][0]*dbasismini_ref[0][i]+Jinv[0][1]*dbasismini_ref[1][i]+Jinv[0][2]*dbasismini_ref[2][i];
-		*(dbasismini+NUMNODESMINI*1+i)=Jinv[1][0]*dbasismini_ref[0][i]+Jinv[1][1]*dbasismini_ref[1][i]+Jinv[1][2]*dbasismini_ref[2][i];
-		*(dbasismini+NUMNODESMINI*2+i)=Jinv[2][0]*dbasismini_ref[0][i]+Jinv[2][1]*dbasismini_ref[1][i]+Jinv[2][2]*dbasismini_ref[2][i];
+	for(int i=0;i<NUMNODESP1b;i++){
+		*(dbasismini+NUMNODESP1b*0+i)=Jinv[0][0]*dbasismini_ref[0][i]+Jinv[0][1]*dbasismini_ref[1][i]+Jinv[0][2]*dbasismini_ref[2][i];
+		*(dbasismini+NUMNODESP1b*1+i)=Jinv[1][0]*dbasismini_ref[0][i]+Jinv[1][1]*dbasismini_ref[1][i]+Jinv[1][2]*dbasismini_ref[2][i];
+		*(dbasismini+NUMNODESP1b*2+i)=Jinv[2][0]*dbasismini_ref[0][i]+Jinv[2][1]*dbasismini_ref[1][i]+Jinv[2][2]*dbasismini_ref[2][i];
 	}
 
@@ -1310,31 +1311,31 @@
 
 	/*Nodal function 1*/
-	dbasis[NUMNODESMINI*0+0]=-0.5*(1.0-zeta)/2.0;
-	dbasis[NUMNODESMINI*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
-	dbasis[NUMNODESMINI*2+0]=-0.5*gauss->coord1;
+	dbasis[NUMNODESP1b*0+0]=-0.5*(1.0-zeta)/2.0;
+	dbasis[NUMNODESP1b*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
+	dbasis[NUMNODESP1b*2+0]=-0.5*gauss->coord1;
 	/*Nodal function 2*/
-	dbasis[NUMNODESMINI*0+1]=0.5*(1.0-zeta)/2.0;
-	dbasis[NUMNODESMINI*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
-	dbasis[NUMNODESMINI*2+1]=-0.5*gauss->coord2;
+	dbasis[NUMNODESP1b*0+1]=0.5*(1.0-zeta)/2.0;
+	dbasis[NUMNODESP1b*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
+	dbasis[NUMNODESP1b*2+1]=-0.5*gauss->coord2;
 	/*Nodal function 3*/
-	dbasis[NUMNODESMINI*0+2]=0.;
-	dbasis[NUMNODESMINI*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
-	dbasis[NUMNODESMINI*2+2]=-0.5*gauss->coord3;
+	dbasis[NUMNODESP1b*0+2]=0.;
+	dbasis[NUMNODESP1b*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
+	dbasis[NUMNODESP1b*2+2]=-0.5*gauss->coord3;
 	/*Nodal function 4*/
-	dbasis[NUMNODESMINI*0+3]=-0.5*(1.0+zeta)/2.0;
-	dbasis[NUMNODESMINI*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
-	dbasis[NUMNODESMINI*2+3]=0.5*gauss->coord1;
+	dbasis[NUMNODESP1b*0+3]=-0.5*(1.0+zeta)/2.0;
+	dbasis[NUMNODESP1b*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
+	dbasis[NUMNODESP1b*2+3]=0.5*gauss->coord1;
 	/*Nodal function 5*/
-	dbasis[NUMNODESMINI*0+4]=0.5*(1.0+zeta)/2.0;
-	dbasis[NUMNODESMINI*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
-	dbasis[NUMNODESMINI*2+4]=0.5*gauss->coord2;
+	dbasis[NUMNODESP1b*0+4]=0.5*(1.0+zeta)/2.0;
+	dbasis[NUMNODESP1b*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
+	dbasis[NUMNODESP1b*2+4]=0.5*gauss->coord2;
 	/*Nodal function 6*/
-	dbasis[NUMNODESMINI*0+5]=0.;
-	dbasis[NUMNODESMINI*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
-	dbasis[NUMNODESMINI*2+5]=0.5*gauss->coord3;
+	dbasis[NUMNODESP1b*0+5]=0.;
+	dbasis[NUMNODESP1b*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
+	dbasis[NUMNODESP1b*2+5]=0.5*gauss->coord3;
 	/*Nodal function 7*/
-	dbasis[NUMNODESMINI*0+6]=27.*(1.+zeta)*(1.-zeta)*(-.5*gauss->coord2*gauss->coord3 + .5*gauss->coord1*gauss->coord3);
-	dbasis[NUMNODESMINI*1+6]=27.*(1.+zeta)*(1.-zeta)*SQRT3*(-1./6.*gauss->coord2*gauss->coord3 - 1./6.*gauss->coord1*gauss->coord3 +1./3.*gauss->coord1*gauss->coord2);
-	dbasis[NUMNODESMINI*2+6]=27*gauss->coord1*gauss->coord2*gauss->coord3*(-2.0*zeta);
+	dbasis[NUMNODESP1b*0+6]=27.*(1.+zeta)*(1.-zeta)*(-.5*gauss->coord2*gauss->coord3 + .5*gauss->coord1*gauss->coord3);
+	dbasis[NUMNODESP1b*1+6]=27.*(1.+zeta)*(1.-zeta)*SQRT3*(-1./6.*gauss->coord2*gauss->coord3 - 1./6.*gauss->coord1*gauss->coord3 +1./3.*gauss->coord1*gauss->coord2);
+	dbasis[NUMNODESP1b*2+6]=27*gauss->coord1*gauss->coord2*gauss->coord3*(-2.0*zeta);
 }
 /*}}}*/
@@ -1504,5 +1505,5 @@
 	switch(this->element_type){
 		case P1Enum:   return NUMNODESP1;
-		case MINIEnum: return NUMNODESMINI;
+		case MINIEnum: return NUMNODESP1b;
 		default:       _error_("Element type "<<EnumToStringx(this->element_type)<<" not supported yet");
 	}
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 15495)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 15496)
@@ -122,5 +122,5 @@
 
 	/*Intermediaries*/
-	const int numnodes = this->GetNumberOfNodes();
+	const int numnodes = this->NumberOfNodes();
 
 	/*Output */
@@ -896,25 +896,8 @@
 	_assert_(sidlist);
 	_assert_(nodes);
-	int numnodes = this->GetNumberOfNodes();
+	int numnodes = this->NumberOfNodes();
 
 	for(int i=0;i<numnodes;i++){
 		sidlist[i]=nodes[i]->Sid();
-	}
-}
-/*}}}*/
-/*FUNCTION Tria::GetNumberOfNodes{{{*/
-int Tria::GetNumberOfNodes(void){
-
-	if(this->nodes==NULL) return 0;
-
-	switch(this->element_type){
-		case P1Enum:
-			return 3;
-		case P1DGEnum:
-			return 3;
-		case P2Enum:
-			return 6;
-		default:
-			_error_("Element type "<<EnumToStringx(this->element_type)<<" not supported yet");
 	}
 }
@@ -5216,5 +5199,5 @@
 
 	/*Constants*/
-	const int numnodes = this->GetNumberOfNodes();
+	const int numnodes = this->NumberOfNodes();
 	const int numdof   = NDOF2*numnodes;
 
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 15495)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 15496)
@@ -83,5 +83,4 @@
 		void   Delta18oParameterization(void);
 		int    GetNodeIndex(Node* node);
-		int    GetNumberOfNodes(void);
 		void   GetNodesSidList(int* sidlist);
 		int    Sid();
Index: /issm/trunk-jpl/src/c/classes/Elements/TriaRef.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/TriaRef.cpp	(revision 15495)
+++ /issm/trunk-jpl/src/c/classes/Elements/TriaRef.cpp	(revision 15496)
@@ -530,11 +530,11 @@
 			/*Nodal function 1*/
 			dbasis[NUMNODESP1*0+0] = -0.5;
-			dbasis[NUMNODESP1*1+0] = -1.0/(2.0*SQRT3);
+			dbasis[NUMNODESP1*1+0] = -SQRT3/6.;
 			/*Nodal function 2*/
 			dbasis[NUMNODESP1*0+1] = 0.5;
-			dbasis[NUMNODESP1*1+1] = -1.0/(2.0*SQRT3);
+			dbasis[NUMNODESP1*1+1] = -SQRT3/6.;
 			/*Nodal function 3*/
 			dbasis[NUMNODESP1*0+2] = 0;
-			dbasis[NUMNODESP1*1+2] = 1.0/SQRT3;
+			dbasis[NUMNODESP1*1+2] = SQRT3/3.;
 			return;
 		case P2Enum: