Index: /issm/trunk/src/c/modules/ModelProcessorx/CreateParameters.cpp
===================================================================
--- /issm/trunk/src/c/modules/ModelProcessorx/CreateParameters.cpp	(revision 7407)
+++ /issm/trunk/src/c/modules/ModelProcessorx/CreateParameters.cpp	(revision 7408)
@@ -72,4 +72,5 @@
 	parameters->AddObject(new IntParam(StabilizeConstraintsEnum,iomodel->stabilize_constraints));
 	parameters->AddObject(new DoubleParam(StokesReconditioningEnum,iomodel->stokesreconditioning));
+	parameters->AddObject(new IntParam(ShelfDampeningEnum,iomodel->shelf_dampening));
 	parameters->AddObject(new DoubleParam(ViscosityOvershootEnum,iomodel->viscosity_overshoot));
 	parameters->AddObject(new BoolParam(WaitOnLockEnum,iomodel->waitonlock));
Index: /issm/trunk/src/c/objects/Elements/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 7407)
+++ /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 7408)
@@ -2625,9 +2625,11 @@
 	/*Intermediaries*/
 	int         i,j,ig;
-	int         approximation;
+	int         approximation,shelf_dampening;
 	double      gravity,rho_water,bed,water_pressure;
+	double      damper,normal_vel,vx,vy,vz,dt;
 	double		xyz_list_tria[NUMVERTICES2D][3];
 	double      xyz_list[NUMVERTICES][3];
 	double		bed_normal[3];
+	double      dz[3];
 	double      l1l6[6]; //for the six nodes of the penta
 	double      Jdet2d;
@@ -2637,4 +2639,5 @@
 	if(!IsOnBed() || !IsOnShelf()) return NULL;
 	inputs->GetParameterValue(&approximation,ApproximationEnum);
+	this->parameters->FindParam(&shelf_dampening,ShelfDampeningEnum);
 	if(approximation!=StokesApproximationEnum && approximation!=MacAyealStokesApproximationEnum && approximation!=PattynStokesApproximationEnum) return NULL;
 	ElementVector* pe=new ElementVector(nodes,NUMVERTICES,this->parameters,StokesApproximationEnum);
@@ -2645,4 +2648,7 @@
 	GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
 	Input* bed_input=inputs->GetInput(BedEnum); _assert_(bed_input);
+	Input* vx_input=inputs->GetInput(VxEnum);   _assert_(vx_input);
+	Input* vy_input=inputs->GetInput(VyEnum);   _assert_(vy_input);
+	Input* vz_input=inputs->GetInput(VzEnum);   _assert_(vz_input);
 
 	for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) xyz_list_tria[i][j]=xyz_list[i][j];
@@ -2657,9 +2663,20 @@
 		GetNodalFunctionsP1(l1l6, gauss);
 
+		BedNormal(&bed_normal[0],xyz_list_tria);
 		bed_input->GetParameterValue(&bed, gauss);
-		BedNormal(&bed_normal[0],xyz_list_tria);
+		if(shelf_dampening){ //add dampening to avoid too high vertical velocities when not in hydrostatic equilibrium
+			bed_input->GetParameterDerivativeValue(&dz[0],&xyz_list[0][0],gauss);
+			vx_input->GetParameterValue(&vx, gauss);
+			vy_input->GetParameterValue(&vy, gauss);
+			vz_input->GetParameterValue(&vz, gauss);
+			dt=0;
+			normal_vel=bed_normal[0]*vx+bed_normal[1]*vy+bed_normal[2]*vz;
+			damper=gravity*rho_water*pow(1+pow(dz[0],2)+pow(dz[1],2),0.5)*normal_vel*dt;
+		}
+		else damper=0;
+		vz_input->GetParameterValue(&vz, gauss);
 		water_pressure=gravity*rho_water*bed;
 
-		for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) pe->values[i*NDOF4+j]+=water_pressure*gauss->weight*Jdet2d*l1l6[i]*bed_normal[j];
+		for(i=0;i<NUMVERTICES2D;i++) for(j=0;j<3;j++) pe->values[i*NDOF4+j]+=(water_pressure+damper)*gauss->weight*Jdet2d*l1l6[i]*bed_normal[j];
 	}
 
Index: /issm/trunk/src/c/objects/IoModel.cpp
===================================================================
--- /issm/trunk/src/c/objects/IoModel.cpp	(revision 7407)
+++ /issm/trunk/src/c/objects/IoModel.cpp	(revision 7408)
@@ -191,4 +191,5 @@
 	IoModelFetchData(&this->prognostic_DG,iomodel_handle,"prognostic_DG");
 	IoModelFetchData(&this->stokesreconditioning,iomodel_handle,"stokesreconditioning");
+	IoModelFetchData(&this->shelf_dampening,iomodel_handle,"shelf_dampening");
 	IoModelFetchData(&this->waitonlock,iomodel_handle,"waitonlock");
 	IoModelFetchData(&this->kff,iomodel_handle,"kff");
@@ -353,4 +354,5 @@
 	this->prognostic_DG=0;
 	this->stokesreconditioning=0;
+	this->shelf_dampening=0;
 	this->waitonlock=0;
 
Index: /issm/trunk/src/c/objects/IoModel.h
===================================================================
--- /issm/trunk/src/c/objects/IoModel.h	(revision 7407)
+++ /issm/trunk/src/c/objects/IoModel.h	(revision 7408)
@@ -119,4 +119,5 @@
 		double  viscosity_overshoot;
 		double  stokesreconditioning;
+		int     shelf_dampening;
 		double  cm_noisedmp;
 		double* cm_min;