Index: /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 24138)
+++ /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 24139)
@@ -157,4 +157,6 @@
 	/*Initialize cumdeltalthickness input*/
 	InputUpdateFromConstantx(elements,0.,SealevelriseCumDeltathicknessEnum);
+	/*Initialize ThicknessResidual input*/
+	InputUpdateFromConstantx(elements,0.,ThicknessResidualEnum);
 
 	/*Get what we need for ocean-induced basal melting*/
@@ -851,17 +853,25 @@
 	int numnodes = element->GetNumberOfNodes();
 	IssmDouble* newthickness = xNew<IssmDouble>(numnodes);
+	IssmDouble* thicknessresidual = xNew<IssmDouble>(numnodes);
 
 	/*Use the dof list to index into the solution vector: */
-   IssmDouble minthickness = element->FindParam(MasstransportMinThicknessEnum);
+	IssmDouble minthickness = element->FindParam(MasstransportMinThicknessEnum);
 	for(int i=0;i<numnodes;i++){
 		newthickness[i]=solution[doflist[i]];
+		thicknessresidual[i]=0.;	
 		/*Check solution*/
 		if(xIsNan<IssmDouble>(newthickness[i])) _error_("NaN found in solution vector");
 		if(xIsInf<IssmDouble>(newthickness[i])) _error_("Inf found in solution vector");
-      if(newthickness[i]<minthickness) newthickness[i]=minthickness;
+		if(newthickness[i]<minthickness){
+			thicknessresidual[i]=minthickness-newthickness[i];
+			newthickness[i]=minthickness;
+		}
 	}
 	element->AddBasalInput(ThicknessEnum,newthickness,element->GetElementType());
+	element->AddBasalInput(ThicknessResidualEnum,thicknessresidual,element->GetElementType());
+	
 	xDelete<int>(doflist);
 	xDelete<IssmDouble>(newthickness);
+ 	xDelete<IssmDouble>(thicknessresidual);
 
 	/*Update bed and surface accordingly*/
@@ -869,12 +879,12 @@
 	/*Get basal element*/
 	int domaintype; element->FindParam(&domaintype,DomainTypeEnum);
-   Element* basalelement=element;
+	Element* basalelement=element;
 	if(domaintype!=Domain2DhorizontalEnum) basalelement = element->SpawnBasalElement();
 
 	/*Fetch number of nodes and dof for this finite element*/
-   int numvertices = basalelement->GetNumberOfVertices();
+	int numvertices = basalelement->GetNumberOfVertices();
 
 	/*Now, we need to do some "processing"*/
-   newthickness  = xNew<IssmDouble>(numvertices);
+	newthickness  = xNew<IssmDouble>(numvertices);
 	IssmDouble* oldthickness      = xNew<IssmDouble>(numvertices);
 	IssmDouble* cumdeltathickness = xNew<IssmDouble>(numvertices);
@@ -889,5 +899,5 @@
 
 	/*Get previous base, thickness, surfac and current sealevel and bed:*/
-   basalelement->GetInputListOnVertices(&newthickness[0],ThicknessEnum);
+	basalelement->GetInputListOnVertices(&newthickness[0],ThicknessEnum);
 	basalelement->GetInputListOnVertices(&oldthickness[0],ThicknessOldEnum);
 	basalelement->GetInputListOnVertices(&oldbase[0],BaseEnum);
@@ -897,5 +907,5 @@
 	basalelement->GetInputListOnVertices(&cumdeltathickness[0],SealevelriseCumDeltathicknessEnum);
 
-   /*Do we do grounding line migration?*/
+	/*Do we do grounding line migration?*/
 	bool isgroundingline;
 	element->FindParam(&isgroundingline,TransientIsgroundinglineEnum);
@@ -909,5 +919,5 @@
 
 	/*Find MasstransportHydrostaticAdjustment to figure out how to update the geometry:*/
-   int hydroadjustment;
+	int hydroadjustment;
 	basalelement->FindParam(&hydroadjustment,MasstransportHydrostaticAdjustmentEnum);
 	IssmDouble rho_ice   = basalelement->FindParam(MaterialsRhoIceEnum);
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 24138)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 24139)
@@ -779,4 +779,5 @@
 	ThicknessOldEnum,
 	ThicknessPositiveEnum,
+	ThicknessResidualEnum,
 	VelEnum,
 	VxAverageEnum,
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 24138)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 24139)
@@ -785,4 +785,5 @@
 		case ThicknessOldEnum : return "ThicknessOld";
 		case ThicknessPositiveEnum : return "ThicknessPositive";
+		case ThicknessResidualEnum : return "ThicknessResidual";
 		case VelEnum : return "Vel";
 		case VxAverageEnum : return "VxAverage";
Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 24138)
+++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 24139)
@@ -803,4 +803,5 @@
 	      else if (strcmp(name,"ThicknessOld")==0) return ThicknessOldEnum;
 	      else if (strcmp(name,"ThicknessPositive")==0) return ThicknessPositiveEnum;
+	      else if (strcmp(name,"ThicknessResidual")==0) return ThicknessResidualEnum;
 	      else if (strcmp(name,"Vel")==0) return VelEnum;
 	      else if (strcmp(name,"VxAverage")==0) return VxAverageEnum;
@@ -874,9 +875,9 @@
 	      else if (strcmp(name,"Outputdefinition54")==0) return Outputdefinition54Enum;
 	      else if (strcmp(name,"Outputdefinition55")==0) return Outputdefinition55Enum;
-	      else if (strcmp(name,"Outputdefinition56")==0) return Outputdefinition56Enum;
          else stage=8;
    }
    if(stage==8){
-	      if (strcmp(name,"Outputdefinition57")==0) return Outputdefinition57Enum;
+	      if (strcmp(name,"Outputdefinition56")==0) return Outputdefinition56Enum;
+	      else if (strcmp(name,"Outputdefinition57")==0) return Outputdefinition57Enum;
 	      else if (strcmp(name,"Outputdefinition58")==0) return Outputdefinition58Enum;
 	      else if (strcmp(name,"Outputdefinition59")==0) return Outputdefinition59Enum;
@@ -997,9 +998,9 @@
 	      else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
 	      else if (strcmp(name,"DoubleInput")==0) return DoubleInputEnum;
-	      else if (strcmp(name,"DoubleMatArrayParam")==0) return DoubleMatArrayParamEnum;
          else stage=9;
    }
    if(stage==9){
-	      if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
+	      if (strcmp(name,"DoubleMatArrayParam")==0) return DoubleMatArrayParamEnum;
+	      else if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
 	      else if (strcmp(name,"DoubleMatParam")==0) return DoubleMatParamEnum;
 	      else if (strcmp(name,"DoubleParam")==0) return DoubleParamEnum;
@@ -1120,9 +1121,9 @@
 	      else if (strcmp(name,"Materials")==0) return MaterialsEnum;
 	      else if (strcmp(name,"Matestar")==0) return MatestarEnum;
-	      else if (strcmp(name,"Matice")==0) return MaticeEnum;
          else stage=10;
    }
    if(stage==10){
-	      if (strcmp(name,"Matlitho")==0) return MatlithoEnum;
+	      if (strcmp(name,"Matice")==0) return MaticeEnum;
+	      else if (strcmp(name,"Matlitho")==0) return MatlithoEnum;
 	      else if (strcmp(name,"MatrixParam")==0) return MatrixParamEnum;
 	      else if (strcmp(name,"MaxAbsVx")==0) return MaxAbsVxEnum;
@@ -1243,9 +1244,9 @@
 	      else if (strcmp(name,"StringArrayParam")==0) return StringArrayParamEnum;
 	      else if (strcmp(name,"StringExternalResult")==0) return StringExternalResultEnum;
-	      else if (strcmp(name,"StringParam")==0) return StringParamEnum;
          else stage=11;
    }
    if(stage==11){
-	      if (strcmp(name,"SubelementFriction1")==0) return SubelementFriction1Enum;
+	      if (strcmp(name,"StringParam")==0) return StringParamEnum;
+	      else if (strcmp(name,"SubelementFriction1")==0) return SubelementFriction1Enum;
 	      else if (strcmp(name,"SubelementFriction2")==0) return SubelementFriction2Enum;
 	      else if (strcmp(name,"SubelementMelt1")==0) return SubelementMelt1Enum;