Index: /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp	(revision 23386)
+++ /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp	(revision 23387)
@@ -607,10 +607,4 @@
 				/*centered*/
 				//thickness[i] = old_thickness[i]*(1.0+opening-closing)/(1.0-opening+closing);
-
-				//if(element->nodes[i]->Sid()==2299){
-				//	printf("for node %i \n",element->nodes[i]->Sid());
-				//	printf(" old thickness is  %g \n",old_thickness[i]);
-				//	printf(" new thickness is  %g \n",thickness[i]);
-				//}
 				/*Take care of otherthikening*/
 				if(thickness[i]>max_thick){
@@ -709,21 +703,20 @@
 		if(old_active[i]>0.){
 			vec_mask->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
-			/* If epl thickness gets under colapse thickness, close the layer */
-			if(epl_thickness[i]<colapse_thick){
-				//if(element->nodes[i]->Sid()==2299){
-					//printf("masking node %i with thickness %f\n",element->nodes[i]->Sid(),epl_thickness[i]);}
+			/* If epl thickness gets under colapse thickness, close the layer if there is no residual*/
+			if(epl_thickness[i]<colapse_thick && residual[i]<=0.){
 				vec_mask->SetValue(basalelement->nodes[i]->Sid(),0.,INS_VAL);
 				recurence->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
 			}
 		}
-		/*If node is now closed bring its thickness back to initial*/
 		if (old_active[i]==0.){
-			epl_thickness[i]=init_thick;
 			/*Activate if we have a residual from sediment*/
 			if(residual[i]>0.){
 				vec_mask->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
-				if(old_active[i]==0.){
-					recurence->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
-				}
+				recurence->SetValue(basalelement->nodes[i]->Sid(),1.,INS_VAL);
+			}
+			else{
+				/*If node is now closed bring its thickness back to initial*/
+				epl_thickness[i]=init_thick;
+				vec_mask->SetValue(basalelement->nodes[i]->Sid(),0.,INS_VAL);
 			}
 		}
Index: /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp	(revision 23386)
+++ /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp	(revision 23387)
@@ -28,5 +28,5 @@
 
 	/*retrieve some parameters: */
-	 bool   issmb;
+	bool   issmb;
 	iomodel->FindConstant(&issmb,"md.transient.issmb");
 	iomodel->FindConstant(&hydrology_model,"md.hydrology.model");
@@ -374,5 +374,5 @@
 	}
 	if(smb_model==SMBgradientscomponentsEnum){
-		surface_runoff_input = basalelement->GetInput(SmbRunoffEnum); _assert_(surface_runoff_input); 
+		surface_runoff_input = basalelement->GetInput(SmbRunoffEnum); _assert_(surface_runoff_input);
 	}
 
@@ -518,5 +518,4 @@
 	/*Use the dof list to index into the solution vector: */
 
-	/*need to introduce thawed test herre*/
 	for(int i=0;i<numnodes;i++){
 		values[i] =solution[doflist[i]];
@@ -546,5 +545,13 @@
 		kappa=kmax*pow(10.,penalty_factor);
 
+
+		Input* thawed_element_input = basalelement->GetInput(HydrologydcMaskThawedEltEnum); _assert_(thawed_element_input);
+		thawed_element_input->GetInputValue(&thawed_element);
+
 		for(int i=0;i<numnodes;i++){
+			/*frozen elements heads are set to base elevation*/
+			if(!thawed_element){
+				values[i]=base[i];
+			}
 			GetHydrologyDCInefficientHmax(&h_max,basalelement,basalelement->GetNode(i));
 			if(values[i]>h_max) {
@@ -555,7 +562,6 @@
 				residual[i] = 0.;
 			}
-			//adding base in min to take into account heads under bed wich don't change N
-			//pressure[i]=(rho_ice*g*thickness[i])-(rho_freshwater*g*(max((min(h_max,values[i])-base[i]),0.0)));
 			pressure[i]=(rho_ice*g*thickness[i])-(rho_freshwater*g*(values[i]-base[i]));
+
 		}
 		xDelete<IssmDouble>(thickness);
Index: /issm/trunk-jpl/src/c/solutionsequences/solutionsequence_hydro_nonlinear.cpp
===================================================================
--- /issm/trunk-jpl/src/c/solutionsequences/solutionsequence_hydro_nonlinear.cpp	(revision 23386)
+++ /issm/trunk-jpl/src/c/solutionsequences/solutionsequence_hydro_nonlinear.cpp	(revision 23387)
@@ -190,10 +190,12 @@
 				Reduceloadx(pf,Kfs,ys); delete Kfs;
 				delete uf_epl;
+
 				femmodel->profiler->Start(SOLVER);
 				Solverx(&uf_epl,Kff,pf,uf_epl_sub_iter,df,femmodel->parameters);
 				femmodel->profiler->Stop(SOLVER);
+
 				delete Kff; delete pf; delete df;
 				delete uf_epl_sub_iter;
-				uf_epl_sub_iter=uf_epl->Duplicate();
+				uf_epl_sub_iter=uf_epl->Duplicate();_assert_(uf_epl_sub_iter);
 				uf_epl->Copy(uf_epl_sub_iter);
 				delete ug_epl;
@@ -215,6 +217,4 @@
 					_error_("   maximum number of EPL iterations (" << hydro_maxiter << ") exceeded");
 				}
-				//If there is some colapse go through sediment again
-				/* if(ThickCount<L2Count)eplconverged=true; */
 				eplcount++;