Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 18055)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 18056)
@@ -3176,14 +3176,38 @@
 
 	/*Intermediaries*/
-	int    vertexpidlist[NUMVERTICES];
-	IssmDouble lambda[NUMVERTICES];
-	IssmDouble gradient_g[NUMVERTICES];
-
-	/*Compute Gradient*/
+	IssmDouble lambda,potential,weight,Jdet;
+	IssmDouble grade_g[NUMVERTICES];
+	IssmDouble basis[NUMVERTICES];
+	IssmDouble xyz_list[NUMVERTICES][3];
+	int        vertexpidlist[NUMVERTICES];
+
+	/*Recover inputs*/
+	Input* lambda_input  = inputs->GetInput(AdjointEnum);                            _assert_(lambda_input);
+	Input* weights_input = inputs->GetInput(InversionCostFunctionsCoefficientsEnum); _assert_(weights_input);
+
+	/* Get node coordinates and dof list: */
+	::GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
 	GradientIndexing(&vertexpidlist[0],control_index);
-	GetInputListOnVertices(&lambda[0],AdjointEnum);
-	for(int i=0;i<NUMVERTICES;i++) gradient_g[i]=lambda[i];
-
-	gradient->SetValues(NUMVERTICES,vertexpidlist,gradient_g,INS_VAL);
+
+	/* Start  looping on the number of gaussian points: */
+	Gauss* gauss=new GaussTria(2);
+	for(int ig=gauss->begin();ig<gauss->end();ig++){
+
+		gauss->GaussPoint(ig);
+
+		GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss);
+		GetNodalFunctions(&basis[0],gauss);
+		weights_input->GetInputValue(&weight,gauss,Balancethickness2MisfitEnum);
+		lambda_input->GetInputValue(&lambda,gauss);
+
+		/*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
+		for(int i=0;i<NUMVERTICES;i++){
+			grade_g[i]+= - weight*Jdet*gauss->weight*basis[i]*lambda;
+		}
+	}
+	gradient->SetValues(NUMVERTICES,vertexpidlist,grade_g,ADD_VAL);
+
+	/*Clean up and return*/
+	delete gauss;
 }
 /*}}}*/
Index: /issm/trunk-jpl/src/c/classes/FemModel.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/FemModel.cpp	(revision 18055)
+++ /issm/trunk-jpl/src/c/classes/FemModel.cpp	(revision 18056)
@@ -1208,5 +1208,5 @@
 	IssmDouble J_sum;
 
-	IssmDouble  weight,thicknessobs,thickness;
+	IssmDouble  weight,thicknessobs,thickness,potential;
 	IssmDouble  Jdet;
 	IssmDouble* xyz_list = NULL;
@@ -1224,4 +1224,5 @@
 		Input* thickness_input   =element->GetInput(ThicknessEnum);                          _assert_(thickness_input);
 		Input* thicknessobs_input=element->GetInput(InversionThicknessObsEnum);              _assert_(thicknessobs_input);
+		Input* potential_input   =element->GetInput(PotentialEnum);                          _assert_(potential_input);
 
 		/* Start  looping on the number of gaussian points: */
@@ -1238,7 +1239,9 @@
 			thickness_input->GetInputValue(&thickness,gauss);
 			thicknessobs_input->GetInputValue(&thicknessobs,gauss);
+			potential_input->GetInputValue(&potential,gauss);
 			gauss->GaussPoint(ig);
 
-			J +=0.5*(thickness*thickness - thicknessobs*thicknessobs)*(thickness*thickness - thicknessobs*thicknessobs)*weight*Jdet*gauss->weight;
+			//J +=0.5*(thickness*thickness - thicknessobs*thicknessobs)*(thickness*thickness - thicknessobs*thicknessobs)*weight*Jdet*gauss->weight;
+			J +=.5*potential*potential*weight*Jdet*gauss->weight;
 		}