Index: /issm/trunk/src/c/objects/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Penta.cpp	(revision 1195)
+++ /issm/trunk/src/c/objects/Penta.cpp	(revision 1196)
@@ -1222,5 +1222,4 @@
 void  Penta::GradjDrag(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type){
 	
-	
 	Tria* tria=NULL;
 	
@@ -1241,5 +1240,21 @@
 #define __FUNCT__ "Penta::GradjB"
 void  Penta::GradjB(Vec grad_g,void* inputs,int analysis_type,int sub_analysis_type){
-	throw ErrorException(__FUNCT__," not supported yet!");
+
+	Tria* tria=NULL;
+
+	if (collapse){
+		/*Bail out element if collapsed (2d) and not on bed*/
+		if (!onbed) return;
+
+		/*This element should be collapsed into a tria element at its base. Create this tria element, 
+		 * and compute gardj*/
+		tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria (lower face).
+		tria->GradjB(grad_g,inputs,analysis_type,sub_analysis_type);
+		delete tria;
+		return;
+	}
+	else{
+		throw ErrorException(__FUNCT__," not supported yet!");
+	}
 }
         
@@ -1888,5 +1903,5 @@
 
 #undef __FUNCT__ 
-#define __FUNCT__ "Penta::CreateKMatrix"
+#define __FUNCT__ "Penta::GetParameterValue"
 void Penta::GetParameterValue(double* pvalue, double* v_list,double* gauss_l1l2l3l4){
 	
@@ -2011,4 +2026,34 @@
 			}
 		} //if (strcmp(field_name,"velocity")==0)
+		else if (strcmp(field_name,"gradj")==0){
+
+			/* node data: */
+			int          dof1;
+			double       fieldel;
+
+			/*this penta is a collapsed macayeal. For each node on the base of this penta, 
+			 * we grab the field. Once we know the field, we follow the upper nodes, 
+			 * inserting the same field value into field, until we reach the surface: */
+			for(i=0;i<3;i++){
+
+				node=nodes[i]; //base nodes
+				dof1=node->GetDofList1();
+
+				/*get field for this base node: */
+				fieldel=field_serial[dof1];
+
+				//go throfieldn all nodes which sit on top of this node, until we reach the surface, 
+				//and plfield  field in field
+				for(;;){
+
+					dof1=node->GetDofList1();
+					VecSetValues(field,1,&dof1,&fieldel,INSERT_VALUES);
+
+					if (node->IsOnSurface())break;
+					/*get next node: */
+					node=node->GetUpperNode();
+				}
+			}
+		}
 		else if ( 
 				(strcmp(field_name,"thickness")==0) ||
@@ -2016,5 +2061,5 @@
 				(strcmp(field_name,"bed")==0)  ||
 				(strcmp(field_name,"slopex")==0)  ||
-				(strcmp(field_name,"slopey")==0) 
+				(strcmp(field_name,"slopey")==0)
 				){
 
Index: /issm/trunk/src/c/parallel/GradJCompute.cpp
===================================================================
--- /issm/trunk/src/c/parallel/GradJCompute.cpp	(revision 1195)
+++ /issm/trunk/src/c/parallel/GradJCompute.cpp	(revision 1196)
@@ -41,4 +41,7 @@
 	Vec grad_g=NULL;
 
+	/*flags: */
+	int debug=0;
+	int dim=-1;
 	
 	/*some parameters:*/
@@ -48,4 +51,6 @@
 	femmodel->parameters->FindParam((void*)&solverstring,"solverstring");
 	femmodel->parameters->FindParam((void*)&control_type,"control_type");
+	femmodel->parameters->FindParam((void*)&debug,"debug");
+	femmodel->parameters->FindParam((void*)&dim,"dim");
 
 	_printf_("%s\n","      recover solution for this stiffness and right hand side:");
@@ -75,4 +80,10 @@
 	Gradjx( &grad_g, numberofnodes,femmodel->elements,femmodel->nodes, femmodel->loads, femmodel->materials, 
 				inputs,analysis_type,sub_analysis_type,control_type);
+
+	if (dim==3){
+
+		_printf_("%s\n","      extruding gradient...");
+		FieldExtrudex( grad_g, femmodel->elements,femmodel->nodes,femmodel->loads,femmodel->materials,"gradj",1);
+	}
 	
 	/*Free ressources:*/