Index: /issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp	(revision 15637)
+++ /issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp	(revision 15638)
@@ -1367,27 +1367,27 @@
 		case P2Enum:
 			/*Nodal function 1*/
-			dbasis[NUMNODESP2*0+0 ] = .5*zeta*(1.-zeta)*(-2.*gauss->coord1 + 0.5);
-			dbasis[NUMNODESP2*1+0 ] = .5*zeta*(1.-zeta)*(-2.*SQRT3/3.*gauss->coord1 + SQRT3/6.);
-			dbasis[NUMNODESP2*2+0 ] = .5*(1.-2.*zeta)*gauss->coord1*(2.*gauss->coord1-1.);
+			dbasis[NUMNODESP2*0+0 ] = .5*zeta*(zeta-1.)*(-2.*gauss->coord1 + 0.5);
+			dbasis[NUMNODESP2*1+0 ] = .5*zeta*(zeta-1.)*(-2.*SQRT3/3.*gauss->coord1 + SQRT3/6.);
+			dbasis[NUMNODESP2*2+0 ] = .5*(2.*zeta-1.)*gauss->coord1*(2.*gauss->coord1-1.);
 			/*Nodal function 2*/
-			dbasis[NUMNODESP2*0+1 ] = .5*zeta*(1.-zeta)*(+2.*gauss->coord2 - 0.5);
-			dbasis[NUMNODESP2*1+1 ] = .5*zeta*(1.-zeta)*(-2.*SQRT3/3.*gauss->coord2 + SQRT3/6.);
-			dbasis[NUMNODESP2*2+1 ] = .5*(1.-2.*zeta)*gauss->coord2*(2.*gauss->coord2-1.);
+			dbasis[NUMNODESP2*0+1 ] = .5*zeta*(zeta-1.)*(+2.*gauss->coord2 - 0.5);
+			dbasis[NUMNODESP2*1+1 ] = .5*zeta*(zeta-1.)*(-2.*SQRT3/3.*gauss->coord2 + SQRT3/6.);
+			dbasis[NUMNODESP2*2+1 ] = .5*(2.*zeta-1.)*gauss->coord2*(2.*gauss->coord2-1.);
 			/*Nodal function 3*/
 			dbasis[NUMNODESP2*0+2 ] = 0.;
-			dbasis[NUMNODESP2*1+2 ] = .5*zeta*(1.-zeta)*(4.*SQRT3/3.*gauss->coord3 - SQRT3/3.);
-			dbasis[NUMNODESP2*2+2 ] = .5*(1.-2.*zeta)*gauss->coord3*(2.*gauss->coord3-1.);
+			dbasis[NUMNODESP2*1+2 ] = .5*zeta*(zeta-1.)*(4.*SQRT3/3.*gauss->coord3 - SQRT3/3.);
+			dbasis[NUMNODESP2*2+2 ] = .5*(2.*zeta-1.)*gauss->coord3*(2.*gauss->coord3-1.);
 			/*Nodal function 4*/
-			dbasis[NUMNODESP2*0+3 ] = .5*zeta*(1.+zeta)*(-2.*gauss->coord1 + 0.5);
-			dbasis[NUMNODESP2*1+3 ] = .5*zeta*(1.+zeta)*(-2.*SQRT3/3.*gauss->coord1 + SQRT3/6.);
-			dbasis[NUMNODESP2*2+3 ] = .5*(1.+2.*zeta)*gauss->coord1*(2.*gauss->coord1-1.);
+			dbasis[NUMNODESP2*0+3 ] = .5*zeta*(zeta+1.)*(-2.*gauss->coord1 + 0.5);
+			dbasis[NUMNODESP2*1+3 ] = .5*zeta*(zeta+1.)*(-2.*SQRT3/3.*gauss->coord1 + SQRT3/6.);
+			dbasis[NUMNODESP2*2+3 ] = .5*(2.*zeta+1.)*gauss->coord1*(2.*gauss->coord1-1.);
 			/*Nodal function 5*/
-			dbasis[NUMNODESP2*0+4 ] = .5*zeta*(1.+zeta)*(+2.*gauss->coord2 - 0.5);
-			dbasis[NUMNODESP2*1+4 ] = .5*zeta*(1.+zeta)*(-2.*SQRT3/3.*gauss->coord2 + SQRT3/6.);
-			dbasis[NUMNODESP2*2+4 ] = .5*(1.+2.*zeta)*gauss->coord2*(2.*gauss->coord2-1.);
+			dbasis[NUMNODESP2*0+4 ] = .5*zeta*(zeta+1.)*(+2.*gauss->coord2 - 0.5);
+			dbasis[NUMNODESP2*1+4 ] = .5*zeta*(zeta+1.)*(-2.*SQRT3/3.*gauss->coord2 + SQRT3/6.);
+			dbasis[NUMNODESP2*2+4 ] = .5*(2.*zeta+1.)*gauss->coord2*(2.*gauss->coord2-1.);
 			/*Nodal function 6*/
 			dbasis[NUMNODESP2*0+5 ] = 0.;
-			dbasis[NUMNODESP2*1+5 ] = .5*zeta*(1.+zeta)*(4.*SQRT3/3.*gauss->coord3 - SQRT3/3.);
-			dbasis[NUMNODESP2*2+5 ] = .5*(1.+2.*zeta)*gauss->coord3*(2.*gauss->coord3-1.);
+			dbasis[NUMNODESP2*1+5 ] = .5*zeta*(zeta+1.)*(4.*SQRT3/3.*gauss->coord3 - SQRT3/3.);
+			dbasis[NUMNODESP2*2+5 ] = .5*(2.*zeta+1.)*gauss->coord3*(2.*gauss->coord3-1.);
 
 			/*Nodal function 7*/
@@ -1405,27 +1405,27 @@
 
 			/*Nodal function 10*/
-			dbasis[NUMNODESP2*0+9 ] = zeta*(1.-zeta)*gauss->coord3;
-			dbasis[NUMNODESP2*1+9 ] = .5*zeta*(1.-zeta)*(+4.*SQRT3/3.*gauss->coord2 - 2.*SQRT3/3.*gauss->coord3);
-			dbasis[NUMNODESP2*2+9 ] = 2.*gauss->coord3*gauss->coord2*(1.-2.*zeta);
+			dbasis[NUMNODESP2*0+9 ] = zeta*(zeta-1.)*gauss->coord3;
+			dbasis[NUMNODESP2*1+9 ] = .5*zeta*(zeta-1.)*(+4.*SQRT3/3.*gauss->coord2 - 2.*SQRT3/3.*gauss->coord3);
+			dbasis[NUMNODESP2*2+9 ] = 2.*gauss->coord3*gauss->coord2*(2.*zeta-1.);
 			/*Nodal function 11*/
-			dbasis[NUMNODESP2*0+10] = -zeta*(1.-zeta)*gauss->coord3;
-			dbasis[NUMNODESP2*1+10] = .5*zeta*(1.-zeta)*(+4.*SQRT3/3.*gauss->coord1 - 2.*SQRT3/3.*gauss->coord3);
-			dbasis[NUMNODESP2*2+10] = 2.*gauss->coord3*gauss->coord1*(1.-2.*zeta);
+			dbasis[NUMNODESP2*0+10] = -zeta*(zeta-1.)*gauss->coord3;
+			dbasis[NUMNODESP2*1+10] = .5*zeta*(zeta-1.)*(+4.*SQRT3/3.*gauss->coord1 - 2.*SQRT3/3.*gauss->coord3);
+			dbasis[NUMNODESP2*2+10] = 2.*gauss->coord3*gauss->coord1*(2.*zeta-1.);
 			/*Nodal function 12*/
-			dbasis[NUMNODESP2*0+11] = zeta*(1.-zeta)*(gauss->coord1-gauss->coord2);
-			dbasis[NUMNODESP2*1+11] = .5*zeta*(1.-zeta)*(-2.*SQRT3/3.*(gauss->coord1+gauss->coord2));
-			dbasis[NUMNODESP2*2+11] = 2.*gauss->coord1*gauss->coord2*(1.-2.*zeta);
+			dbasis[NUMNODESP2*0+11] = zeta*(zeta-1.)*(gauss->coord1-gauss->coord2);
+			dbasis[NUMNODESP2*1+11] = .5*zeta*(zeta-1.)*(-2.*SQRT3/3.*(gauss->coord1+gauss->coord2));
+			dbasis[NUMNODESP2*2+11] = 2.*gauss->coord1*gauss->coord2*(2.*zeta-1.);
 			/*Nodal function 13*/
-			dbasis[NUMNODESP2*0+12] = zeta*(1.+zeta)*gauss->coord3;
-			dbasis[NUMNODESP2*1+12] = .5*zeta*(1.+zeta)*(+4.*SQRT3/3.*gauss->coord2 - 2.*SQRT3/3.*gauss->coord3);
-			dbasis[NUMNODESP2*2+12] = 2.*gauss->coord3*gauss->coord2*(1.+2.*zeta);
+			dbasis[NUMNODESP2*0+12] = zeta*(zeta+1.)*gauss->coord3;
+			dbasis[NUMNODESP2*1+12] = .5*zeta*(zeta+1.)*(+4.*SQRT3/3.*gauss->coord2 - 2.*SQRT3/3.*gauss->coord3);
+			dbasis[NUMNODESP2*2+12] = 2.*gauss->coord3*gauss->coord2*(2.*zeta+1.);
 			/*Nodal function 14*/
-			dbasis[NUMNODESP2*0+13] = -zeta*(1.+zeta)*gauss->coord3;
-			dbasis[NUMNODESP2*1+13] = .5*zeta*(1.+zeta)*(+4.*SQRT3/3.*gauss->coord1 - 2.*SQRT3/3.*gauss->coord3);
-			dbasis[NUMNODESP2*2+13] = 2.*gauss->coord3*gauss->coord1*(1.+2.*zeta);
+			dbasis[NUMNODESP2*0+13] = -zeta*(zeta+1.)*gauss->coord3;
+			dbasis[NUMNODESP2*1+13] = .5*zeta*(zeta+1.)*(+4.*SQRT3/3.*gauss->coord1 - 2.*SQRT3/3.*gauss->coord3);
+			dbasis[NUMNODESP2*2+13] = 2.*gauss->coord3*gauss->coord1*(2.*zeta+1.);
 			/*Nodal function 15*/
-			dbasis[NUMNODESP2*0+14] = zeta*(1.+zeta)*(gauss->coord1-gauss->coord2);
-			dbasis[NUMNODESP2*1+14] = .5*zeta*(1.+zeta)*(-2.*SQRT3/3.*(gauss->coord1+gauss->coord2));
-			dbasis[NUMNODESP2*2+14] = 2.*gauss->coord1*gauss->coord2*(1.+2.*zeta);
+			dbasis[NUMNODESP2*0+14] = zeta*(zeta+1.)*(gauss->coord1-gauss->coord2);
+			dbasis[NUMNODESP2*1+14] = .5*zeta*(zeta+1.)*(-2.*SQRT3/3.*(gauss->coord1+gauss->coord2));
+			dbasis[NUMNODESP2*2+14] = 2.*gauss->coord1*gauss->coord2*(2.*zeta+1.);
 			return;
 		default: