Index: /issm/trunk-jpl/src/c/cores/movingfront_core.cpp
===================================================================
--- /issm/trunk-jpl/src/c/cores/movingfront_core.cpp	(revision 24035)
+++ /issm/trunk-jpl/src/c/cores/movingfront_core.cpp	(revision 24036)
@@ -88,4 +88,8 @@
 	delete analysis;
 
+	/*Kill ice berg to avoid free body motion*/
+	if(VerboseSolution()) _printf0_("   looking for icebergs to kill\n");
+	KillIcebergsx(femmodel);
+
 	/*Reset levelset if needed*/
 	if(reinit_frequency && (step%reinit_frequency==0)){
Index: /issm/trunk-jpl/src/c/modules/KillIcebergsx/KillIcebergsx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/KillIcebergsx/KillIcebergsx.cpp	(revision 24035)
+++ /issm/trunk-jpl/src/c/modules/KillIcebergsx/KillIcebergsx.cpp	(revision 24036)
@@ -7,30 +7,19 @@
 #include "../../toolkits/toolkits.h"
 
+#include "../InputUpdateFromVectorx/InputUpdateFromVectorx.h"
+
 void KillIcebergsx(FemModel* femmodel){
 
 	/*Intermediaries*/
-	IssmDouble* local_mask = NULL;
-	const int MAXVERTICES = 6;
-	bool      found1;
-	int       sidlist[MAXVERTICES];
-	int       lidlist[MAXVERTICES];
+	int lid;
 
-	/*retrieve vertex info*/
-	int nbv_global = femmodel->vertices->NumberOfVertices();
-	int nbv_local  = femmodel->vertices->NumberOfVerticesLocal();
-	if(nbv_global==0)  return;
-	Vector<IssmDouble>* vec_connected_to_land=new Vector<IssmDouble>(nbv_local,nbv_global);
+	/*retrieve vertex info and prepare element flag to speed up process*/
+	int         nbv_local    = femmodel->vertices->Size();
+	IssmDouble *local_mask   = xNewZeroInit<IssmDouble>(nbv_local);
+	bool       *element_flag = xNewZeroInit<bool>(femmodel->elements->Size());
 
-	/*Prepare element flag to speed up process*/
-	bool* element_flag = xNewZeroInit<bool>(femmodel->elements->Size());
-
-	/*Fill vector with 1 once for all*/
-	IssmDouble eflags[MAXVERTICES];
-	for(int i=0;i<MAXVERTICES;i++) eflags[i] = 1.;
-
-	/*Step 1, go through all elements and put 1 in vec_connected_to_land if the element is grounded*/
+	/*Step 1, go through all elements and put 1 in local_mask if the element is grounded*/
 	for(int i=0;i<femmodel->elements->Size();i++){
 		Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
-
 		if(!element->IsIceInElement()){
 			/*Nothing to do, just flag element to speed up the computation*/
@@ -39,29 +28,26 @@
 		else{
 			if(element->IsGrounded()){
-				int  numvertices = element->GetNumberOfVertices();
-				if(numvertices>MAXVERTICES) _error_("need to increase MAXVERTICES");
-				element->GetVerticesSidList(&sidlist[0]);
-				vec_connected_to_land->SetValues(numvertices,&sidlist[0],&eflags[0],ADD_VAL);
+				int numvertices = element->GetNumberOfVertices();
+				for(int v=0;v<numvertices;v++) local_mask[element->vertices[v]->Lid()] = 1.;
 			}
 		}
 	}
-	vec_connected_to_land->Assemble();
 
-	/*Now we have 2 loops, one across cpus, and one for each cpus: we are going to propagate the mask if an element
-	 * is connected to a positive mask already.
-	 * We then communicate to the other partitions. We stop when the mask stops changing*/
+	/*Now we have 2 loops, one across cpus, and one for each cpus: we are going
+	 * to propagate the mask if an element is connected to a positive mask
+	 * already.  We then communicate to the other partitions. We stop when the
+	 * mask stops changing*/
 	bool keepsyncing = true;
 	while(keepsyncing){
 
 		/*Get local mask from parallel vector*/
-		if(local_mask) xDelete<IssmDouble>(local_mask);
-		femmodel->GetLocalVectorWithClonesVertices(&local_mask,vec_connected_to_land);
+		femmodel->SyncLocalVectorWithClonesVertices(local_mask);
 
 		/*Local iterations on partition*/
-		bool keepgoing = true;
+		bool keepgoing    = true;
 		int  didsomething = 0;
-		int  iter      = 1;
+		int  iter         = 1;
 		while(keepgoing){
-			_printf0_("   -- Kill icebergs: iteration "<<iter<<"\n");
+			_printf0_("   -- Kill icebergs: local iteration "<<iter<<"\n");
 
 			keepgoing    = false;
@@ -72,8 +58,8 @@
 				if(!element_flag[i]){
 					int numvertices = element->GetNumberOfVertices();
-					element->GetVerticesLidList(&lidlist[0]);
 					bool found1 = false;
 					for(int j=0;j<numvertices;j++){
-						if(local_mask[lidlist[j]]>0.){
+						lid = element->vertices[j]->Lid();
+						if(local_mask[lid]>0.){
 							found1 = true;
 							break;
@@ -83,6 +69,7 @@
 						element_flag[i] = true;
 						for(int j=0;j<numvertices;j++){
-							if(local_mask[lidlist[j]]==0.){
-								local_mask[lidlist[j]]=1.;
+							lid = element->vertices[j]->Lid();
+							if(local_mask[lid]==0.){
+								local_mask[lid]=1.;
 								keepgoing = true;
 								didsomething = 1;
@@ -96,29 +83,16 @@
 
 		/*Check how many iterations all cpus did*/
-		int didsomething_max;
-		ISSM_MPI_Reduce(&didsomething,&didsomething_max,1,ISSM_MPI_INT,ISSM_MPI_MAX,0,IssmComm::GetComm());
-		ISSM_MPI_Bcast(&didsomething_max,1,ISSM_MPI_INT,0,IssmComm::GetComm());
-		if(didsomething_max==0){
+		int iter_max;
+		ISSM_MPI_Reduce(&iter,&iter_max,1,ISSM_MPI_INT,ISSM_MPI_MAX,0,IssmComm::GetComm());
+		ISSM_MPI_Bcast(&iter_max,1,ISSM_MPI_INT,0,IssmComm::GetComm());
+		if(iter_max==2){
+			/*If iter is only 2, nothing else was changed in the while loop above (iter is initialized as 1 and then ++)*/
 			keepsyncing = false;
-		}
-		else{
-			for(int i=0;i<femmodel->elements->Size();i++){
-				Element* element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
-
-				if(element->IsIceInElement()){
-						int  numvertices = element->GetNumberOfVertices();
-						element->GetVerticesSidList(&sidlist[0]);
-						element->GetVerticesLidList(&lidlist[0]);
-						for(int j=0;j<numvertices;j++) eflags[j] = local_mask[lidlist[j]];
-						vec_connected_to_land->SetValues(numvertices,&sidlist[0],&eflags[0],ADD_VAL);
-					}
-			}
-			vec_connected_to_land->Assemble();
 		}
 	}
 
+	InputUpdateFromVectorx(femmodel,local_mask,PressureEnum,VertexLIdEnum);
 	/*Cleanup*/
 	xDelete<bool>(element_flag);
-	delete vec_connected_to_land;
 
 	/*OK, now deactivate iceberg and count the number of deactivated vertices*/
@@ -128,8 +102,8 @@
 		if(element->IsIceInElement()){
 			int  numvertices = element->GetNumberOfVertices();
-			element->GetVerticesLidList(&lidlist[0]);
 			bool deactivate = false;
 			for(int j=0;j<numvertices;j++){
-				if(local_mask[lidlist[j]]==0.){
+				lid = element->vertices[j]->Lid();
+				if(local_mask[lid]==0.){
 					deactivate = true;
 					break;