Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 28038)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 28039)
@@ -286,5 +286,5 @@
 	IssmDouble  groundedice,bed,sealevel;
 
-	/*Depth average B for stress calculation*/
+	/*Depth average velocity for stress calculation*/
 	this->InputDepthAverageAtBase(VxEnum,VxAverageEnum);
 	this->InputDepthAverageAtBase(VyEnum,VyAverageEnum);
@@ -625,4 +625,95 @@
 		delete gauss;
 	}
+}
+/*}}}*/
+void       Penta::CalvingRateCalvingMIP(){/*{{{*/
+
+	if(!this->IsOnBase()) return;
+
+	IssmDouble  calvingrate[NUMVERTICES];
+	IssmDouble  calvingratex[NUMVERTICES];
+	IssmDouble  calvingratey[NUMVERTICES];
+	int			experiment = 1;  /* exp:1 by default */
+	int         dim, domaintype;
+	IssmDouble	vx, vy, vel, c, wrate;
+	IssmDouble  time, groundedice, yts;
+
+	/*Get problem dimension and whether there is moving front or not*/
+	this->FindParam(&domaintype,DomainTypeEnum);
+	this->FindParam(&time,TimeEnum);
+	this->FindParam(&yts,ConstantsYtsEnum);
+
+	switch(domaintype){
+		case Domain2DverticalEnum:   dim = 1; break;
+		case Domain2DhorizontalEnum: dim = 2; break;
+		case Domain3DEnum:           dim = 2; break;
+		default: _error_("mesh "<<EnumToStringx(domaintype)<<" not supported yet");
+	}
+	if(dim==1) _error_("not implemented in 1D...");
+
+	/*Depth average velocity for stress calculation*/
+	this->InputDepthAverageAtBase(VxEnum,VxAverageEnum);
+	this->InputDepthAverageAtBase(VyEnum,VyAverageEnum);
+
+	/*Retrieve all inputs and parameters we will need*/
+	Input *vx_input      = this->GetInput(VxAverageEnum);                                _assert_(vx_input);
+	Input *vy_input      = this->GetInput(VyAverageEnum);                                _assert_(vy_input);
+	Input *wrate_input   = this->GetInput(CalvingAblationrateEnum);               _assert_(wrate_input); 
+	Input* gr_input      = this->GetInput(MaskOceanLevelsetEnum);						_assert_(gr_input);
+
+	/* Use which experiment: use existing Enum */
+	this->FindParam(&experiment, CalvingUseParamEnum);
+
+	/* Start looping on the number of vertices: */
+	GaussPenta gauss;
+	for(int iv=0;iv<3;iv++){
+		gauss.GaussVertex(iv);
+
+		/*Get velocity components */
+		vx_input->GetInputValue(&vx,&gauss);
+		vy_input->GetInputValue(&vy,&gauss);
+		vel=sqrt(vx*vx+vy*vy)+1.e-14;
+
+		/* no calving for grounded ice in EXP4 */
+		gr_input->GetInputValue(&groundedice,&gauss);
+
+		switch (experiment) { 
+			case 1:
+			case 3:
+				/* Exp 1 and 3: set c=v-wrate, wrate=0, so that w=0 */
+				wrate = 0.0;
+				break;
+			case 2:
+				/* Exp 2: set c=v-wrate(given)*/
+				wrate_input->GetInputValue(&wrate,&gauss);
+				break;
+			case 4:
+				/* Exp 4: set c=v-wrate(given), for the first 500 years, then c=0 for the second 500 years*/
+				if((groundedice<0) && (time<=500.0*yts)) {
+				//	wrate_input->GetInputValue(&wrate,&gauss);
+					wrate = -750*sin(2.0*M_PI*time/yts/1000)/yts;  // m/a -> m/s
+				}
+				else {
+					/* no calving on the grounded ice*/
+					wrate = vel;
+				}
+				break;
+			default:
+				_error_("The experiment is not supported yet!");
+		}
+
+		calvingrate[iv] = vel - wrate;
+		calvingratex[iv] = vx - wrate*vx/vel;
+		calvingratey[iv] = vy - wrate*vy/vel;
+	}
+	/*Add input*/
+	this->AddBasalInput(CalvingCalvingrateEnum,&calvingrate[0],P1DGEnum);
+	this->AddBasalInput(CalvingratexEnum,&calvingratex[0],P1DGEnum);
+	this->AddBasalInput(CalvingrateyEnum,&calvingratey[0],P1DGEnum);
+
+	/*Extrude*/
+	this->InputExtrude(CalvingCalvingrateEnum,-1);
+	this->InputExtrude(CalvingratexEnum,-1);
+	this->InputExtrude(CalvingrateyEnum,-1);
 }
 /*}}}*/
@@ -2810,4 +2901,5 @@
 		case CalvingVonmisesEnum:
 		case CalvingLevermannEnum:
+		case CalvingCalvingMIPEnum:
 			calvingratex_input=this->GetInput(CalvingratexEnum); _assert_(calvingratex_input);
 			calvingratey_input=this->GetInput(CalvingrateyEnum); _assert_(calvingratey_input);
@@ -2847,4 +2939,5 @@
 			case CalvingVonmisesEnum:
 			case CalvingLevermannEnum:
+			case CalvingCalvingMIPEnum:
 				calvingratex_input->GetInputValue(&c[0],&gauss);
 				calvingratey_input->GetInputValue(&c[1],&gauss);
Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 28038)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 28039)
@@ -58,4 +58,5 @@
 		void           CalvingFluxLevelset();
 		void           CalvingMeltingFluxLevelset();
+		void           CalvingRateCalvingMIP();
 		IssmDouble     CharacteristicLength(void){_error_("not implemented yet");};
 		void           ComputeBasalStress(void);