Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 15484)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 15485)
@@ -241,31 +241,32 @@
 ElementMatrix* Tria::CreateMassMatrix(void){
 
-	/*constants: */
-	const int    numdof=NDOF1*NUMVERTICES;
-
 	/* Intermediaries */
 	IssmDouble  D,Jdet;
 	IssmDouble  xyz_list[NUMVERTICES][3];
-	IssmDouble  basis[numdof];
-	GaussTria  *gauss = NULL;
-
-	/*Initialize Element matrix*/
-	ElementMatrix* Ke=new ElementMatrix(nodes,NUMVERTICES,this->parameters,NoneApproximationEnum);
-
+
+	/*Fetch number of nodes and dof for this finite element*/
+	int numnodes = this->NumberofNodes();
+	int numdof   = numnodes*NDOF2;
+
+	/*Initialize Element matrix and vectors*/
+	ElementMatrix* Ke    = new ElementMatrix(nodes,numnodes,this->parameters,NoneApproximationEnum);
+	IssmDouble*    basis = xNew<IssmDouble>(numdof);
+
+	/*Retrieve all inputs and parameters*/
 	GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
 
 	/* Start looping on the number of gaussian points: */
-	gauss=new GaussTria(2);
+	GaussTria* gauss=new GaussTria(2);
 	for(int ig=gauss->begin();ig<gauss->end();ig++){
 
 		gauss->GaussPoint(ig);
 
-		GetNodalFunctions(&basis[0],gauss);
+		GetNodalFunctions(basis,gauss);
 		GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss);
 		D=gauss->weight*Jdet;
 
-		TripleMultiply(&basis[0],1,3,1,
+		TripleMultiply(basis,1,numnodes,1,
 					&D,1,1,0,
-					&basis[0],1,3,0,
+					basis,1,numnodes,0,
 					&Ke->values[0],1);
 	}
@@ -273,4 +274,5 @@
 	/*Clean up and return*/
 	delete gauss;
+	xDelete<IssmDouble>(basis);
 	return Ke;
 }
@@ -2781,5 +2783,5 @@
 	int numdof   = numnodes*NDOF2;
 
-	/*Initialize Element matrix, vectors and Gaussian points*/
+	/*Initialize Element matrix and vectors*/
 	ElementMatrix* Ke     = new ElementMatrix(nodes,numnodes,this->parameters,MacAyealApproximationEnum);
 	IssmDouble*    B      = xNew<IssmDouble>(3*numdof);
@@ -2797,5 +2799,5 @@
 
 	/* Start  looping on the number of gaussian points: */
-	GaussTria*     gauss  = new GaussTria(2);
+	GaussTria* gauss  = new GaussTria(2);
 	for(int ig=gauss->begin();ig<gauss->end();ig++){