Index: /issm/trunk/src/c/objects/Elements/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 8413)
+++ /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 8414)
@@ -2827,12 +2827,63 @@
 ElementVector* Penta::CreatePVectorDiagnosticVertBase(void){
 
+
+	/*Constants*/
+	const int    numdof=NDOF1*NUMVERTICES;
+
+	/*Intermediaries */
+	int        i,j,ig;
+	int        approximation;
+	double     xyz_list[NUMVERTICES][3];
+	double     xyz_list_tria[NUMVERTICES2D][3];
+	double     Jdet2d;
+	double     vx,vy,vz,dbdx,dbdy,basalmeltingvalue;
+	double     slope[2];
+	double     l1l6[NUMVERTICES];
+	GaussPenta* gauss=NULL;
+
 	if (!IsOnBed()) return NULL;
 
-	/*Call Tria function*/
-	Tria* tria=(Tria*)SpawnTria(0,1,2); //nodes 0, 1 and 2 make the new tria.
-	ElementVector* pe=tria->CreatePVectorDiagnosticBaseVert();
-	delete tria->matice; delete tria;
+	/*Initialize Element vector*/
+	ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters);
+
+	/*Retrieve all inputs and parameters*/
+	GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
+	for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) xyz_list_tria[i][j]=xyz_list[i][j];
+	inputs->GetParameterValue(&approximation,ApproximationEnum);
+	Input* bed_input=inputs->GetInput(BedEnum);             _assert_(bed_input);
+	Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input);
+	Input* vx_input=inputs->GetInput(VxEnum);               _assert_(vx_input);
+	Input* vy_input=inputs->GetInput(VyEnum);               _assert_(vy_input);
+	Input* vzstokes_input=NULL;
+	if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
+		vzstokes_input=inputs->GetInput(VzStokesEnum);       _assert_(vzstokes_input);
+	}
+
+	/* Start  looping on the number of gaussian points: */
+	gauss=new GaussPenta(0,1,2,2);
+	for(ig=gauss->begin();ig<gauss->end();ig++){
+
+		gauss->GaussPoint(ig);
+
+		basal_melting_input->GetParameterValue(&basalmeltingvalue, gauss);
+		bed_input->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],gauss);
+		vx_input->GetParameterValue(&vx, gauss);
+		vy_input->GetParameterValue(&vy, gauss);
+		if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
+			vzstokes_input->GetParameterValue(&vz, gauss);
+		}
+		else vz=0;
+
+		dbdx=slope[0];
+		dbdy=slope[1];
+
+		GetTriaJacobianDeterminand(&Jdet2d, &xyz_list_tria[0][0],gauss);
+		GetNodalFunctionsP1(&l1l6[0], gauss);
+
+		for(i=0;i<numdof;i++) pe->values[i]+=-Jdet2d*gauss->weight*(vx*dbdx+vy*dbdy-vz-basalmeltingvalue)*l1l6[i];
+	}
 
 	/*Clean up and return*/
+	delete gauss;
 	return pe;
 }
Index: /issm/trunk/src/c/objects/Elements/Tria.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Tria.cpp	(revision 8413)
+++ /issm/trunk/src/c/objects/Elements/Tria.cpp	(revision 8414)
@@ -1717,64 +1717,4 @@
 
 		for(i=0;i<numdof;i++) pe->values[i]+=Jdettria*gauss->weight*(surface_mass_balance_g-basal_melting_g)*L[i];
-	}
-
-	/*Clean up and return*/
-	delete gauss;
-	return pe;
-}
-/*}}}*/
-/*FUNCTION Tria::CreatePVectorDiagnosticBaseVert {{{1*/
-ElementVector* Tria::CreatePVectorDiagnosticBaseVert(void){
-
-	/*Constants*/
-	const int    numdof=NDOF1*NUMVERTICES;
-
-	/*Intermediaries */
-	int        i,j,ig;
-	int        approximation;
-	double     xyz_list[NUMVERTICES][3];
-	double     Jdet;
-	double     vx,vy,vz,dbdx,dbdy,basalmeltingvalue;
-	double     slope[2];
-	double     L[NUMVERTICES];
-	GaussTria* gauss=NULL;
-
-	/*Initialize Element vector*/
-	ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters);
-
-	/*Retrieve all inputs and parameters*/
-	GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
-	inputs->GetParameterValue(&approximation,ApproximationEnum);
-	Input* bed_input=inputs->GetInput(BedEnum);             _assert_(bed_input);
-	Input* basal_melting_input=inputs->GetInput(BasalMeltingRateEnum); _assert_(basal_melting_input);
-	Input* vx_input=inputs->GetInput(VxEnum);               _assert_(vx_input);
-	Input* vy_input=inputs->GetInput(VyEnum);               _assert_(vy_input);
-	Input* vzstokes_input=NULL;
-	if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
-		vzstokes_input=inputs->GetInput(VzStokesEnum);       _assert_(vzstokes_input);
-	}
-
-	/* Start  looping on the number of gaussian points: */
-	gauss=new GaussTria(2);
-	for(ig=gauss->begin();ig<gauss->end();ig++){
-
-		gauss->GaussPoint(ig);
-
-		basal_melting_input->GetParameterValue(&basalmeltingvalue, gauss);
-		bed_input->GetParameterDerivativeValue(&slope[0],&xyz_list[0][0],gauss);
-		vx_input->GetParameterValue(&vx, gauss);
-		vy_input->GetParameterValue(&vy, gauss);
-		if(approximation==PattynStokesApproximationEnum || approximation==MacAyealStokesApproximationEnum){
-			vzstokes_input->GetParameterValue(&vz, gauss);
-		}
-		else vz=0;
-
-		dbdx=slope[0];
-		dbdy=slope[1];
-
-		GetJacobianDeterminant3d(&Jdet, &xyz_list[0][0],gauss);
-		GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);
-
-		for(i=0;i<numdof;i++) pe->values[i]+=-Jdet*gauss->weight*(vx*dbdx+vy*dbdy-vz-basalmeltingvalue)*L[i];
 	}
 
Index: /issm/trunk/src/c/objects/Elements/Tria.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Tria.h	(revision 8413)
+++ /issm/trunk/src/c/objects/Elements/Tria.h	(revision 8414)
@@ -164,5 +164,4 @@
 		ElementVector* CreatePVectorBalancethickness_CG(void);
 		ElementVector* CreatePVectorBalancevelocities(void);
-		ElementVector* CreatePVectorDiagnosticBaseVert(void);
 		ElementVector* CreatePVectorDiagnosticMacAyeal(void);
 		ElementVector* CreatePVectorAdjointHoriz(void);