Index: /issm/trunk/src/c/objects/Elements/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 5847)
+++ /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 5848)
@@ -2189,9 +2189,13 @@
 
 	int approximation;
+	ElementMatrix* Ke=NULL;
+
 	inputs->GetParameterValue(&approximation,ApproximationEnum);
 
 	switch(approximation){
 		case MacAyealApproximationEnum:
-			CreateKMatrixDiagnosticMacAyeal2d( Kgg);
+			Ke=CreateKMatrixDiagnosticMacAyeal2d();
+			if(Ke) Ke->AddToGlobal(Kgg,NULL,NULL);
+			delete Ke;
 			break;
 		case PattynApproximationEnum:
@@ -2291,36 +2295,24 @@
 }
 /*FUNCTION Penta::CreateKMatrixDiagnosticMacAyeal2d{{{1*/
-void Penta::CreateKMatrixDiagnosticMacAyeal2d(Mat Kgg){
+ElementMatrix* Penta::CreateKMatrixDiagnosticMacAyeal2d(void){
 	
-	/*If on water, skip stiffness: */
-	if(IsOnWater()) return;
-
-	/*Figure out if this pentaelem is collapsed. If so, then bailout, except if it is at the 
+	/*Figure out if this penta is collapsed. If so, then bailout, except if it is at the 
 	  bedrock, in which case we spawn a tria element using the 3 first grids, and use it to build 
 	  the stiffness matrix. */
-
-	if (!IsOnBed()){
-		/*This element should be collapsed, but this element is not on the bedrock, therefore all its 
-		 * dofs have already been frozen! Do nothing: */
-		return;
-	}
-	else{
-
-		/*This element should be collapsed into a tria element at its base. Create this tria element, 
-		 *and use its CreateKMatrix functionality to fill the global stiffness matrix: */
-
-		/*Depth Averaging B*/
-		this->InputDepthAverageAtBase(RheologyBEnum,RheologyBbarEnum,MaterialsEnum);
-
-		/*Call Tria function*/
-		Tria* tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
-		tria->CreateKMatrix(Kgg,NULL,NULL);
-		delete tria->matice; delete tria;
-
-		/*Delete B averaged*/
-		this->matice->inputs->DeleteInput(RheologyBbarEnum);
-
-		return;
-	}
+	if (!IsOnBed()) return NULL;
+
+	/*Depth Averaging B*/
+	this->InputDepthAverageAtBase(RheologyBEnum,RheologyBbarEnum,MaterialsEnum);
+
+	/*Call Tria function*/
+	Tria* tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
+	ElementMatrix* Ke=tria->CreateKMatrixDiagnosticMacAyeal();
+	delete tria->matice; delete tria;
+
+	/*Delete B averaged*/
+	this->matice->inputs->DeleteInput(RheologyBbarEnum);
+
+	/*clean up and return*/
+	return Ke;
 }
 /*}}}*/
@@ -5408,6 +5400,5 @@
 	int i;
 
-	const int    numdofpervertex=1;
-	const int    numdof=numdofpervertex*NUMVERTICES;
+	const int    numdof=NDOF1*NUMVERTICES;
 	int*         doflist=NULL;
 	double       values[numdof];
Index: /issm/trunk/src/c/objects/Elements/Penta.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.h	(revision 5847)
+++ /issm/trunk/src/c/objects/Elements/Penta.h	(revision 5848)
@@ -127,5 +127,5 @@
 		void	  CreateKMatrixDiagnosticHoriz( Mat Kgg);
 		void	  CreateKMatrixDiagnosticHutter( Mat Kgg);
-		void	  CreateKMatrixDiagnosticMacAyeal2d(Mat Kgg);
+		ElementMatrix* CreateKMatrixDiagnosticMacAyeal2d(void);
 		void	  CreateKMatrixDiagnosticMacAyeal3d(Mat Kgg);
 		void	  CreateKMatrixDiagnosticPattyn( Mat Kgg);