Index: /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.cpp	(revision 26649)
+++ /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.cpp	(revision 26650)
@@ -154,4 +154,25 @@
    parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.relaxation",HydrologyRelaxationEnum));
 	parameters->AddObject(iomodel->CopyConstantObject("md.hydrology.storage",HydrologyStorageEnum));
+
+	/*Deal with friction parameters*/
+	int frictionlaw;
+	iomodel->FindConstant(&frictionlaw,"md.friction.law");
+	if(frictionlaw==6){
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.gamma",FrictionGammaEnum));
+	}
+	if(frictionlaw==4){
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.gamma",FrictionGammaEnum));
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.coupling",FrictionCouplingEnum));
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.effective_pressure_limit",FrictionEffectivePressureLimitEnum));
+	}
+	if(frictionlaw==1 || frictionlaw==3 || frictionlaw==7){
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.coupling",FrictionCouplingEnum));
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.effective_pressure_limit",FrictionEffectivePressureLimitEnum));
+	}
+	if(frictionlaw==9){
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.gamma",FrictionGammaEnum));
+		parameters->AddObject(iomodel->CopyConstantObject("md.friction.effective_pressure_limit",FrictionEffectivePressureLimitEnum));
+		parameters->AddObject(new IntParam(FrictionCouplingEnum,0));
+	}
 
   /*Requested outputs*/
@@ -465,5 +486,4 @@
 
 	/*Get thickness and base on nodes to apply cap on water head*/
-   IssmDouble* eff_pressure = xNew<IssmDouble>(numnodes);
 	IssmDouble* thickness = xNew<IssmDouble>(numnodes);
 	IssmDouble* bed       = xNew<IssmDouble>(numnodes);
@@ -496,7 +516,4 @@
 	   values[i] = head_old[i] - relaxation*(head_old[i]-values[i]);
 
-		/*Calculate effective pressure*/
-		eff_pressure[i] = rho_ice*g*thickness[i] - rho_water*g*(values[i]-bed[i]);
-
 		if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
 		if(xIsInf<IssmDouble>(values[i])) _error_("Inf found in solution vector");
@@ -505,5 +522,4 @@
 	/*Add input to the element: */
 	element->AddInput(HydrologyHeadEnum,values,element->GetElementType());
-   element->AddInput(EffectivePressureEnum,eff_pressure,P1Enum);
 
 	/*Update reynolds number according to new solution*/
@@ -519,4 +535,7 @@
 	IssmDouble reynolds = conductivity*sqrt(dh[0]*dh[0]+dh[1]*dh[1])/NU;
 	element->AddInput(HydrologyReynoldsEnum,&reynolds,P0Enum);
+
+   /*Compute new effective pressure*/
+   this->UpdateEffectivePressure(element);
 
 	/*Free resources:*/
@@ -526,5 +545,4 @@
 	xDelete<IssmDouble>(xyz_list);
 	xDelete<int>(doflist);
-	xDelete<IssmDouble>(eff_pressure);
    xDelete<IssmDouble>(head_old);
 }/*}}}*/
@@ -724,2 +742,49 @@
 	delete gauss;
 }/*}}}*/
+void HydrologyShaktiAnalysis::UpdateEffectivePressure(FemModel* femmodel){/*{{{*/
+
+	for(Object* & object : femmodel->elements->objects){
+		Element* element=xDynamicCast<Element*>(object);
+		UpdateEffectivePressure(element);
+	}
+
+}/*}}}*/
+void HydrologyShaktiAnalysis::UpdateEffectivePressure(Element* element){/*{{{*/
+
+	/*Skip if water or ice shelf element*/
+	if(element->IsAllFloating()) return;
+
+	/*Intermediaries*/
+	IssmDouble bed,thickness,head;
+
+	/* Fetch number of nodes and allocate output*/
+   int numnodes = element->GetNumberOfNodes();
+   IssmDouble* N = xNew<IssmDouble>(numnodes);
+
+	/*Retrieve all inputs and parameters*/
+	IssmDouble  g          = element->FindParam(ConstantsGEnum);
+	IssmDouble  rho_ice    = element->FindParam(MaterialsRhoIceEnum);
+	IssmDouble  rho_water  = element->FindParam(MaterialsRhoFreshwaterEnum);
+	Input* head_input      = element->GetInput(HydrologyHeadEnum); _assert_(head_input);
+	Input* thickness_input = element->GetInput(ThicknessEnum);     _assert_(thickness_input);
+	Input* base_input      = element->GetInput(BaseEnum);          _assert_(base_input);
+
+
+   Gauss* gauss=element->NewGauss();
+   for (int i=0;i<numnodes;i++){
+      gauss->GaussNode(element->GetElementType(),i);
+
+		base_input->GetInputValue(&bed,gauss);
+		thickness_input->GetInputValue(&thickness,gauss);
+		head_input->GetInputValue(&head,gauss);
+
+		N[i] = rho_ice*g*thickness - rho_water*g*(head-bed);
+	}
+
+	/*Add new gap as an input*/
+	element->AddInput(EffectivePressureEnum,N,element->GetElementType());
+
+	/*Clean up and return*/
+   xDelete<IssmDouble>(N);
+	delete gauss;
+}/*}}}*/
Index: /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.h
===================================================================
--- /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.h	(revision 26649)
+++ /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.h	(revision 26650)
@@ -36,4 +36,6 @@
 		void UpdateGapHeight(FemModel* femmodel);
 		void UpdateGapHeight(Element* element);
+		void UpdateEffectivePressure(FemModel* femmodel);
+		void UpdateEffectivePressure(Element* element);
 };
 #endif