Index: issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp
===================================================================
--- issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp	(revision 17208)
+++ issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp	(revision 17210)
@@ -826,5 +826,5 @@
 	int solution_type, i;
 	bool computebasalmeltingrates=true;
-	bool isdrainage=true;
+	bool isdrainage=false;
 	bool updatebasalconstraints=true;
 
Index: issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp
===================================================================
--- issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp	(revision 17208)
+++ issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp	(revision 17210)
@@ -67,12 +67,19 @@
 	bool save_results;
 
+	/* extrapolate */
+	Analysis* analysis = new ExtrapolationAnalysis();
+	femmodel->parameters->SetParam(VxEnum,ExtrapolationVariableEnum);
+	analysis->Core(femmodel);
+	delete analysis;
+
 	/*activate formulation: */
 	femmodel->SetCurrentConfiguration(LevelsetAnalysisEnum);
 
+	/*recover parameters: */
+	femmodel->parameters->FindParam(&save_results,SaveResultsEnum);
+
 	if(VerboseSolution()) _printf0_("call computational core:\n");
 	solutionsequence_linear(femmodel);
 
-	/*recover parameters: */
-	femmodel->parameters->FindParam(&save_results,SaveResultsEnum);
 	if(save_results){
 		if(VerboseSolution()) _printf0_("   saving results\n");
@@ -317,56 +324,2 @@
 }/*}}}*/
 
-/* Update of constraints */
-void LevelsetAnalysis::UpdateNoIceConstraints(FemModel* femmodel){/*{{{*/
-
-	IssmDouble* mask_ice = GetMaskOfIce(femmodel->elements, femmodel->nodes); 
-
-	for (int i=0;i<femmodel->nodes->Size();i++){
-		Node* node=dynamic_cast<Node*>(femmodel->nodes->GetObjectByOffset(i));
-		if(node->InAnalysis(LevelsetAnalysisEnum)){
-			if(mask_ice[node->Sid()]==1.){//FIXME: what should be done with actual spcs to ice model?
-// 				node->DofInFSet(0); /*remove spc*/ 
-			}
-			else{
-				IssmDouble defval=0.;
-				node->ApplyConstraint(1,defval); /*apply spc*/ 
-			}
-		}
-	}
-}/*}}}*/
-IssmDouble* LevelsetAnalysis::GetMaskOfIce(Elements* elements, Nodes* nodes){/*{{{*/
-
-	int                 i;
-	IssmDouble*         mask_ice      = NULL;
-	Vector<IssmDouble>* vec_mask_ice  = NULL;
-	Element*            element       = NULL;
-
-	/*Initialize vector with number of vertices*/
-	IssmDouble numnodes=nodes->NumberOfNodes(LevelsetAnalysisEnum);
-	vec_mask_ice=new Vector<IssmDouble>(numnodes); //nodes at ice front that have ice at next time step
-	for(i=0;i<numnodes;i++)
-		vec_mask_ice[i]=0.;
-	/*Fill vector vertices that have no contact to ice: */
-	for(i=0;i<elements->Size();i++){
-		element=dynamic_cast<Element*>(elements->GetObjectByOffset(i));
-		SetMaskOfIceElement(vec_mask_ice, element);
-	}
-
-	/*Assemble vector and serialize */
-	vec_mask_ice->Assemble();
-	mask_ice=vec_mask_ice->ToMPISerial();
-	delete vec_mask_ice;
-	return mask_ice;
-
-}/*}}}*/
-void LevelsetAnalysis::SetMaskOfIceElement(Vector<IssmDouble>* vec_mask_ice, Element* element){/*{{{*/
-
-	/* Intermediaries */
-	int numnodes = element->GetNumberOfNodes();
-	
-	if(element->IsIceInElement()){
-		for(int i = 0;i<numnodes;i++){
-			vec_mask_ice->SetValue(element->nodes[i]->Sid(),1.,INS_VAL);
-		}
-	}
-}/*}}}*/
Index: issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.h
===================================================================
--- issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.h	(revision 17208)
+++ issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.h	(revision 17210)
@@ -31,8 +31,4 @@
 	void GetBprime(IssmDouble* Bprime,Element* element,IssmDouble* xyz_list,Gauss* gauss);
 
-	/* Updating constraints */
-	void UpdateNoIceConstraints(FemModel* femmodel);
-	IssmDouble* GetMaskOfIce(Elements* elements, Nodes* nodes);
-	void SetMaskOfIceElement(Vector<IssmDouble>* vec_mask_ice, Element* element);
 };
 #endif
Index: issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp
===================================================================
--- issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp	(revision 17208)
+++ issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp	(revision 17210)
@@ -1368,5 +1368,5 @@
 
 	/*If no front, return NULL*/
-	if(!element->IsIcefront()) return NULL;
+	if(!element->IsZeroLevelset(MaskIceLevelsetEnum)) return NULL;
 
 	/*Intermediaries*/
@@ -1391,10 +1391,6 @@
 	IssmDouble gravity     = element->GetMaterialParameter(ConstantsGEnum);
 	element->GetVerticesCoordinates(&xyz_list);
-
-// 	element->ZeroLevelsetCoordinates(&xyz_list_front,xyz_list,MaskIceLevelsetEnum);
- 	element->GetIcefrontCoordinates(&xyz_list_front,xyz_list,MaskIceLevelsetEnum);
-
+	element->ZeroLevelsetCoordinates(&xyz_list_front,xyz_list,MaskIceLevelsetEnum);
 	element->NormalSection(&normal[0],xyz_list_front);
-// 	element->GetNormalFromLSF(&normal[0]);
 
 	/*Start looping on Gaussian points*/
@@ -1785,5 +1781,5 @@
 
 	/*If no front, return NULL*/
-	if(!element->IsIcefront()) return NULL;
+	if(!element->IsZeroLevelset(MaskIceLevelsetEnum)) return NULL;
 
 	/*Intermediaries*/
@@ -2208,5 +2204,5 @@
 
 	/*If no front, return NULL*/
-	if(!element->IsIcefront()) return NULL;
+	if(!element->IsZeroLevelset(MaskIceLevelsetEnum)) return NULL;
 
 	/*Intermediaries*/
@@ -2974,5 +2970,5 @@
 
 	/*If no front, return NULL*/
-	if(!element->IsIcefront()) return NULL;
+	if(!element->IsZeroLevelset(MaskIceLevelsetEnum)) return NULL;
 
 	/*Intermediaries*/
Index: issm/trunk-jpl/src/c/classes/Elements/Element.h
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 17208)
+++ issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 17210)
@@ -220,5 +220,4 @@
 		virtual int    PressureInterpolation()=0;
 		virtual bool   IsZeroLevelset(int levelset_enum)=0;
-		virtual bool   IsIcefront(void)=0;
 		virtual void   ZeroLevelsetCoordinates(IssmDouble** pxyz_zero,IssmDouble* xyz_list,int levelsetenum)=0;
 		virtual void   GetIcefrontCoordinates(IssmDouble** pxyz_front,IssmDouble* xyz_list,int levelsetenum)=0;
Index: issm/trunk-jpl/src/c/classes/Elements/Penta.h
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 17208)
+++ issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 17210)
@@ -114,5 +114,4 @@
 		int    PressureInterpolation();
 		bool   IsZeroLevelset(int levelset_enum);
-		bool   IsIcefront(void){_error_("not implemented yet");};
 		void   ZeroLevelsetCoordinates(IssmDouble** pxyz_zero,IssmDouble* xyz_list,int levelsetenum);
 		void   GetIcefrontCoordinates(IssmDouble** pxyz_front,IssmDouble* xyz_list,int levelsetenum){_error_("not implemented yet");};
Index: issm/trunk-jpl/src/c/classes/Elements/Seg.h
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Seg.h	(revision 17208)
+++ issm/trunk-jpl/src/c/classes/Elements/Seg.h	(revision 17210)
@@ -159,5 +159,4 @@
 		void        ViscousHeating(IssmDouble* pphi,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not implemented yet");};
 		bool        IsZeroLevelset(int levelset_enum){_error_("not implemented");};
-		bool		IsIcefront(void){_error_("not implemented yet");};
 		void        ZeroLevelsetCoordinates(IssmDouble** pxyz_zero,IssmDouble* xyz_list,int levelsetenum){_error_("not implemented");};
 		void		GetIcefrontCoordinates(IssmDouble** pxyz_front,IssmDouble* xyz_list,int levelsetenum){_error_("not implemented yet");};
Index: issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 17208)
+++ issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 17210)
@@ -909,28 +909,12 @@
 		}
 	}
-	_assert_(nrfrontnodes==2);
-	
-	/* arrange order of frontnodes such that they are oriented counterclockwise */
-	if((NUMVERTICES+indicesfront[0]-indicesfront[1])%NUMVERTICES!=2){
-		int index=indicesfront[0];
-		indicesfront[0]=indicesfront[1];
-		indicesfront[1]=index;
-	}
+
 	IssmDouble* xyz_front = xNew<IssmDouble>(3*nrfrontnodes);
 	/* Return nodes */
 	for(i=0;i<nrfrontnodes;i++){
 		for(dir=0;dir<3;dir++){
-			xyz_front[3*i+dir]=xyz_list[3*indicesfront[i]+dir];
-		}
-	}
-
-// 	for(i=0;i<nrfrontnodes;i++){
-// 		_printf0_("coords frontnode " << i << " :[");
-// 		for(dir=0;dir<3;dir++){
-// 			xyz_front[3*i+dir]=xyz_list[indicesfront[i]+dir];
-// 			_printf0_(xyz_front[3*i+dir] << "; ");
-// 		}
-// 		_printf0_("]\n");
-// 	}
+			xyz_front[3*i+dir]=xyz_list[indicesfront[i]+dir];
+		}
+	}
 
 	*pxyz_front=xyz_front;
@@ -938,5 +922,5 @@
 	xDelete<int>(indicesfront);
 }/*}}}*/
-void  Tria::GetNormalFromLSF(IssmDouble *normal){/*{{{*/
+void  Tria::GetNormalFromLSF(IssmDouble *pnormal){/*{{{*/
 
 	/* Intermediaries */
@@ -945,4 +929,5 @@
 	IssmDouble* xyz_list = NULL;
 	IssmDouble  dlevelset[dim], norm_dlevelset;
+	IssmDouble  normal[dim]={0.};
 
 	/*Retrieve all inputs and parameters*/
@@ -951,6 +936,4 @@
 	
 	counter=0;
-	for(i=0;i<dim;i++) normal[i]=0.;
-
 	Gauss* gauss = this->NewGauss(2);
 	for(int ig=gauss->begin();ig<gauss->end();ig++){
@@ -964,5 +947,7 @@
 	}
 	_assert_(counter>0);
-	for(i=0;i<dim;i++) normal[i]/=IssmDouble(counter);
+	for(i=0;i<dim;i++) normal[i]/counter;
+	
+	pnormal=&normal[0];
 
 	delete gauss;
@@ -2601,5 +2586,5 @@
 	GetInputListOnVertices(&ls[0],levelset_enum);
 
-	/* If levelset function changes sign, or is zero on at least two vertices, there is a zero level set here */
+	/*If the level set is awlays <0, there is no ice front here*/
 	iszerols= false;
 	if(IsIceInElement()){
@@ -2612,26 +2597,4 @@
 }
 /*}}}*/
-/*FUNCTION Tria::IsIcefront{{{*/
-bool Tria::IsIcefront(void){
-
-	bool isicefront;
-	int i,nrice;
-	IssmDouble ls[NUMVERTICES];
-
-	/*Retrieve all inputs and parameters*/
-	GetInputListOnVertices(&ls[0],MaskIceLevelsetEnum);
-
-	/* If only one vertex has ice, there is an ice front here */
-	isicefront= false;
-	if(IsIceInElement()){
-		nrice=0;	
-		for(i=0;i<NUMVERTICES;i++)
-			if(ls[i]<0.) nrice++;
-		if(nrice==1) isicefront= true;
-	}
-	return isicefront;
-}
-/*}}}*/
-
 
 #ifdef _HAVE_RESPONSES_
Index: issm/trunk-jpl/src/c/classes/Elements/Tria.h
===================================================================
--- issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 17208)
+++ issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 17210)
@@ -132,7 +132,6 @@
 		void        ZeroLevelsetCoordinates(IssmDouble** pxyz_zero,IssmDouble* xyz_list,int levelsetenum);
 		void	    GetIcefrontCoordinates(IssmDouble** pxyz_front,IssmDouble* xyz_list,int levelsetenum);
-		void        GetNormalFromLSF(IssmDouble *normal);
+	    void        GetNormalFromLSF(IssmDouble *pnormal);
 		bool        IsZeroLevelset(int levelset_enum);
-		bool        IsIcefront(void);
 
 		#ifdef _HAVE_RESPONSES_
Index: issm/trunk-jpl/src/c/classes/FemModel.cpp
===================================================================
--- issm/trunk-jpl/src/c/classes/FemModel.cpp	(revision 17208)
+++ issm/trunk-jpl/src/c/classes/FemModel.cpp	(revision 17210)
@@ -407,4 +407,6 @@
 	if(VerboseModule()) _printf0_("   Updating constraints for time: " << time << "\n");
 
+	// analysis->UpdateConstraints();
+	
 	/*Second, constraints might be time dependent: */
 	SpcNodesx(nodes,constraints,parameters,analysis_type); 
Index: issm/trunk-jpl/src/c/cores/transient_core.cpp
===================================================================
--- issm/trunk-jpl/src/c/cores/transient_core.cpp	(revision 17208)
+++ issm/trunk-jpl/src/c/cores/transient_core.cpp	(revision 17210)
@@ -162,20 +162,7 @@
 		if(islevelset){
 			if(VerboseSolution()) _printf0_("   computing movement of ice boundaries\n");
-
-			/* extrapolate */
-			Analysis* extanalysis = new ExtrapolationAnalysis();
-			int vars[2] = {VxEnum, VyEnum};
-			for(int iv=0; i<2;i++){
-				femmodel->parameters->SetParam(vars[i],ExtrapolationVariableEnum);
-				extanalysis->Core(femmodel);
-			}
-			delete extanalysis;
-
-			LevelsetAnalysis* lsanalysis = new LevelsetAnalysis();
-			/* solve level-set equation */
-			lsanalysis->Core(femmodel);
-			/* update spc boundary conditions for ice model */
-// 			lsanalysis->UpdateNoIceConstraints(femmodel);
-			delete lsanalysis;
+			analysis = new LevelsetAnalysis();
+			analysis->Core(femmodel);
+			delete analysis;
 		}
 
