Index: /issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp	(revision 21207)
+++ /issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp	(revision 21208)
@@ -120,5 +120,7 @@
 	iomodel->FetchDataToInput(elements,"md.initialization.vx",VxEnum);
 	iomodel->FetchDataToInput(elements,"md.initialization.vy",VyEnum);
-
+   iomodel->FetchDataToInput(elements,"md.hydrology.head",HydrologyHeadOldEnum); 
+	/*iomodel->FetchDataToInput(elements,"md.hydrology.eff_pressure",EffectivePressureEnum);
+*/ 
 	iomodel->FindConstant(&frictionlaw,"md.friction.law");
 	/*Friction law variables*/
@@ -210,5 +212,5 @@
 	/*Intermediaries */
 	IssmDouble  Jdet,meltrate,G,dh[2],B,A,n;
-	IssmDouble  gap,bed,thickness,head,ieb;
+	IssmDouble  gap,bed,thickness,head,ieb,head_old;
 	IssmDouble  lr,br,vx,vy,beta;
 	IssmDouble  alpha2,frictionheat;
@@ -241,4 +243,5 @@
 	Input* lr_input             = element->GetInput(HydrologyBumpSpacingEnum);       _assert_(lr_input);
 	Input* br_input             = element->GetInput(HydrologyBumpHeightEnum);        _assert_(br_input);
+   Input* headold_input        = element->GetInput(HydrologyHeadOldEnum);           _assert_(headold_input);
 
 	/*Get conductivity from inputs*/
@@ -270,4 +273,5 @@
 		vx_input->GetInputValue(&vx,gauss);
 		vy_input->GetInputValue(&vy,gauss);
+      headold_input->GetInputValue(&head_old,gauss);
 
 		/*Get ice A parameter*/
@@ -293,4 +297,8 @@
 		if(pressure_water>pressure_ice) pressure_water = pressure_ice;
 
+		/*Get water pressure from previous time step to use in lagged creep term*/
+		IssmDouble pressure_water_old = rho_water*g*(head_old-bed);
+		if(pressure_water_old>pressure_ice) pressure_water_old = pressure_ice;
+
 		/*Compute change in sensible heat due to changes in pressure melting point*/
    	dpressure_water[0] = rho_water*g*(dh[0] - dbed[0]);
@@ -304,10 +312,9 @@
 		 (
 		  meltrate*(1/rho_water-1/rho_ice)
-		  +A*pow(fabs(pressure_ice-pressure_water),n-1)*(pressure_ice-pressure_water)*gap
+		  +A*pow(fabs(pressure_ice - pressure_water_old),n-1)*(pressure_ice - pressure_water_old)*gap
 		  -beta*sqrt(vx*vx+vy*vy)
 		  +ieb
 		  )*basis[i];     	
-      	}
-
+	}
 	/*Clean up and return*/
 	xDelete<IssmDouble>(xyz_list);
@@ -379,7 +386,4 @@
 	head_input->GetInputDerivativeAverageValue(&dh[0],xyz_list);
 	IssmDouble conductivity = GetConductivity(element);
-//if(element->Id()==1){
-//	printf("Conductivity in UpdateInputFromSolution: %g \n",conductivity);
-//}
 
 	IssmDouble reynolds = conductivity*sqrt(dh[0]*dh[0]+dh[1]*dh[1])/(2.*NU);
@@ -535,5 +539,5 @@
 	/*Divide by connectivity*/
 	newgap = newgap/totalweights;
-	IssmDouble mingap = 0.001;
+	IssmDouble mingap = 1e-5;
 	if(newgap<mingap) newgap=mingap;
 
Index: /issm/trunk-jpl/src/c/classes/Elements/Element.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 21207)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 21208)
@@ -1524,4 +1524,6 @@
 				name==HydrologydcEplThicknessEnum ||
 				name==HydrologydcMaskEplactiveNodeEnum ||
+				name==HydrologyHeadEnum ||
+	         name==HydrologyHeadOldEnum ||		
 				name==StressbalanceConvergenceNumStepsEnum || 
 				name==MeshVertexonbaseEnum 
Index: /issm/trunk-jpl/src/c/cores/hydrology_core.cpp
===================================================================
--- /issm/trunk-jpl/src/c/cores/hydrology_core.cpp	(revision 21207)
+++ /issm/trunk-jpl/src/c/cores/hydrology_core.cpp	(revision 21208)
@@ -78,11 +78,12 @@
 	}
 
-	else if (hydrology_model==HydrologysommersEnum){
-		femmodel->SetCurrentConfiguration(HydrologySommersAnalysisEnum);
-		solutionsequence_nonlinear(femmodel,modify_loads);
+	else if (hydrology_model==HydrologysommersEnum){	
+		femmodel->SetCurrentConfiguration(HydrologySommersAnalysisEnum);	
+      InputDuplicatex(femmodel,HydrologyHeadEnum,HydrologyHeadOldEnum);	
+		solutionsequence_nonlinear(femmodel,modify_loads); 
 		if(VerboseSolution()) _printf0_("   updating gap height\n");
 		HydrologySommersAnalysis* analysis = new HydrologySommersAnalysis();
 		analysis->UpdateGapHeight(femmodel);
-		delete analysis;	
+		delete analysis;
 		
 		if(save_results){
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 21207)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 21208)
@@ -166,4 +166,5 @@
 	HydrologysommersEnum,
 	HydrologyHeadEnum,
+   HydrologyHeadOldEnum,
 	HydrologyGapHeightEnum,
 	HydrologyBumpSpacingEnum,
@@ -544,5 +545,5 @@
 	GiaCrossSectionShapeEnum,
 	GiadWdtEnum,
-	GiaWEnum,
+	GiaWEnum, 
 	/*}}}*/
 	/*Results{{{*/