Index: /issm/trunk-jpl/src/c/classes/Elements/Element.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 19153)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 19154)
@@ -1367,4 +1367,22 @@
 IssmDouble Element::PureIceEnthalpy(IssmDouble pressure){/*{{{*/
 	return this->matpar->PureIceEnthalpy(pressure);
+}/*}}}*/
+void       Element::ResetIceLevelset(void){/*{{{*/
+
+	/*Get Ice Levelset on vertices*/
+	int         numvertices = this->GetNumberOfVertices();
+	IssmDouble *values      = xNew<IssmDouble>(numvertices);
+	this->GetInputListOnVertices(values,MaskIceLevelsetEnum);
+
+	/*Project on {-1,+1}*/
+	for(int i=0;i<numvertices;i++){
+		if(values[i]>0) values[i]=+1.;
+		else            values[i]=-1.;
+	}
+	this->AddInput(MaskIceLevelsetEnum,values,P1Enum);
+
+	/*Cleanup and return*/
+	xDelete<IssmDouble>(values);
+
 }/*}}}*/
 void       Element::ResultInterpolation(int* pinterpolation,int* pnodesperelement,int output_enum){/*{{{*/
Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 19153)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 19154)
@@ -123,4 +123,5 @@
 		ElementVector*     NewElementVector(int approximation_enum=NoneApproximationEnum);
 		IssmDouble         PureIceEnthalpy(IssmDouble pressure);
+		void               ResetIceLevelset(void);
 		void               ResultInterpolation(int* pinterpolation,int*nodesperelement,int output_enum);
 		void               ResultToPatch(IssmDouble* values,int nodesperelement,int output_enum);
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 19153)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 19154)
@@ -336,4 +336,5 @@
 	IssmDouble  calvingrate[NUMVERTICES];
 	IssmDouble  lambda1,lambda2,ex,ey,vx,vy,vel;
+	IssmDouble  sigma_vm,sigma_max;
 
 	/* Get node coordinates and dof list: */
@@ -341,6 +342,8 @@
 
 	/*Retrieve all inputs and parameters we will need*/
-	Input* vx_input=inputs->GetInput(VxEnum);        _assert_(vx_input);
-	Input* vy_input=inputs->GetInput(VyEnum);        _assert_(vy_input);
+	Input* vx_input = inputs->GetInput(VxEnum); _assert_(vx_input);
+	Input* vy_input = inputs->GetInput(VyEnum); _assert_(vy_input);
+	IssmDouble  B   = this->GetMaterialParameter(MaterialsRheologyBbarEnum);
+	IssmDouble  n   = this->GetMaterialParameter(MaterialsRheologyNEnum);
 
 	/* Start looping on the number of vertices: */
@@ -366,7 +369,11 @@
 		lambda2 = max(lambda2,0.);
 
+		/*Calculate sigma_vm*/
+		sigma_vm  = sqrt(3./2.) * B * pow(lambda1*lambda1 + lambda2*lambda2,1./(2.*n));
+		sigma_max = 350.e+3;
+
 		/*Assign values*/
-		calvingratex[iv]=vx*pow(lambda1 + lambda2,1./3.)*3.e+2;
-		calvingratey[iv]=vy*pow(lambda1 + lambda2,1./3.)*3.e+2;
+		calvingratex[iv]=vx*sigma_vm/sigma_max;
+		calvingratey[iv]=vy*sigma_vm/sigma_max;
 		calvingrate[iv]=sqrt(calvingratex[iv]*calvingratex[iv] + calvingratey[iv]*calvingratey[iv]);
 	}