Index: /issm/trunk-jpl/src/c/classes/Loads/Friction.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Loads/Friction.cpp	(revision 22846)
+++ /issm/trunk-jpl/src/c/classes/Loads/Friction.cpp	(revision 22847)
@@ -60,9 +60,7 @@
 
 	/*diverse: */
-	int         CoupledFlag;
 	IssmDouble  q_exp;
 	IssmDouble  C_param;
 	IssmDouble  As;
-	IssmDouble  Neff;
 	IssmDouble  n;
 	IssmDouble  alpha;
@@ -70,37 +68,13 @@
 	IssmDouble  vx,vy,vz,vmag;
 	IssmDouble  alpha_complement;
-	IssmDouble  base,sealevel,thickness;
 
 	/*Recover parameters: */
 	element->GetInputValue(&q_exp,FrictionQEnum);
 	element->GetInputValue(&C_param,FrictionCEnum);
-
 	element->GetInputValue(&As,gauss,FrictionAsEnum);
 	element->GetInputValue(&n,gauss,MaterialsRheologyNEnum);
-	element->GetInputValue(&thickness, gauss,ThicknessEnum);
-	element->GetInputValue(&base, gauss,BaseEnum);
-	element->GetInputValue(&sealevel, gauss,SealevelEnum);
-
-	IssmDouble rho_water   = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
-	IssmDouble rho_ice     = element->GetMaterialParameter(MaterialsRhoIceEnum);
-	IssmDouble gravity     = element->GetMaterialParameter(ConstantsGEnum);
-	element->parameters->FindParam(&CoupledFlag,FrictionCouplingEnum);
-
-	switch(CoupledFlag){
-		//This should be removed at some point and the Enum uniformalized
-		case 0:
-			Neff=gravity*(rho_ice*thickness+rho_water*(base-sealevel));
-			break;
-		case 1:
-			element->GetInputValue(&Neff,gauss,FrictionEffectivePressureEnum);
-			break;
-		case 2:
-			element->GetInputValue(&Neff,gauss,EffectivePressureTimeAverageEnum);
-			break;
-		default:
-			_error_("not supported");
-	}
-
-	if(Neff<0)Neff=0;
+
+	/*Get effective pressure*/
+	IssmDouble Neff = EffectivePressure(gauss);
 
 	//We need the velocity magnitude to evaluate the basal stress:
@@ -180,11 +154,8 @@
 
 	/*diverse: */
-	int         CoupledFlag;
 	IssmDouble  r,s;
 	IssmDouble  vx,vy,vz,vmag;
 	IssmDouble  drag_p,drag_q;
-	IssmDouble  Neff;
 	IssmDouble  drag_coefficient;
-	IssmDouble  base,thickness,sealevel;
 	IssmDouble  alpha_complement;
 
@@ -192,31 +163,12 @@
 	element->GetInputValue(&drag_p,FrictionPEnum);
 	element->GetInputValue(&drag_q,FrictionQEnum);
-	element->GetInputValue(&thickness, gauss,ThicknessEnum);
-	element->GetInputValue(&base, gauss,BaseEnum);
-	element->GetInputValue(&sealevel, gauss,SealevelEnum);
 	element->GetInputValue(&drag_coefficient, gauss,FrictionCoefficientEnum);
-	IssmDouble rho_water   = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
-	IssmDouble rho_ice     = element->GetMaterialParameter(MaterialsRhoIceEnum);
-	IssmDouble gravity     = element->GetMaterialParameter(ConstantsGEnum);
-	element->parameters->FindParam(&CoupledFlag,FrictionCouplingEnum);
-	//compute r and q coefficients: */
+
+	/*compute r and q coefficients: */
 	r=drag_q/drag_p;
 	s=1./drag_p;
 
-	//From base and thickness, compute effective pressure when drag is viscous, or get Neff from forcing or coupled to hydrologymodel:
-	switch(CoupledFlag){
-		case 0:
-			Neff=gravity*(rho_ice*thickness+rho_water*(base-sealevel));
-			break;
-		case 1:
-			element->GetInputValue(&Neff,gauss,FrictionEffectivePressureEnum);
-			break;
-		case 2:
-			element->GetInputValue(&Neff,gauss,EffectivePressureTimeAverageEnum);
-			break;
-		default:
-			_error_("not supported");
-	}
-	if(Neff<0)Neff=0;
+	/*Get effective pressure*/
+	IssmDouble Neff = EffectivePressure(gauss);
 
 	//We need the velocity magnitude to evaluate the basal stress:
@@ -290,8 +242,6 @@
 
 	/*diverse: */
-	int         CoupledFlag;
 	IssmDouble  r,s;
 	IssmDouble  drag_p, drag_q;
-	IssmDouble  Neff;
 	IssmDouble  thickness,base,floatation_thickness,sealevel;
 	IssmDouble  vx,vy,vz,vmag;
@@ -307,8 +257,7 @@
 	element->GetInputValue(&drag_coefficient, gauss,FrictionCoefficientEnum);
 	element->GetInputValue(&drag_coefficient_coulomb, gauss,FrictionCoefficientcoulombEnum);
-	IssmDouble rho_water        = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
-	IssmDouble rho_ice          = element->GetMaterialParameter(MaterialsRhoIceEnum);
-	IssmDouble gravity          = element->GetMaterialParameter(ConstantsGEnum);
-	element->parameters->FindParam(&CoupledFlag,FrictionCouplingEnum);
+	IssmDouble rho_water = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
+	IssmDouble rho_ice   = element->GetMaterialParameter(MaterialsRhoIceEnum);
+	IssmDouble gravity   = element->GetMaterialParameter(ConstantsGEnum);
 
 	//compute r and q coefficients: */
@@ -316,21 +265,9 @@
 	s=1./drag_p;
 
-	//From base and thickness, compute effective pressure when drag is viscous:
-	switch(CoupledFlag){
-		case 0:
-			Neff=gravity*(rho_ice*thickness+rho_water*(base-sealevel));
-			break;
-		case 1:
-			element->GetInputValue(&Neff,gauss,FrictionEffectivePressureEnum);
-			break;
-		case 2:
-			element->GetInputValue(&Neff,gauss,EffectivePressureTimeAverageEnum);
-			break;
-		default:
-			_error_("not supported");
-	}
-
-	if(Neff<0)Neff=0;
-
+
+	/*Get effective pressure*/
+	IssmDouble Neff = EffectivePressure(gauss);
+
+	/*Compute velocity magnitude*/
 	switch(dim){
 		case 1:
@@ -379,10 +316,7 @@
 
 	/*diverse: */
-	int         CoupledFlag;
 	IssmDouble  q_exp;
 	IssmDouble  C_param;
 	IssmDouble  As;
-
-	IssmDouble  Neff;
 	IssmDouble  n;
 
@@ -392,5 +326,4 @@
 	IssmDouble  vx,vy,vz,vmag;
 	IssmDouble  alpha2;
-	IssmDouble  base,thickness,sealevel;
 
 	/*Recover parameters: */
@@ -399,28 +332,9 @@
 	element->GetInputValue(&As,gauss,FrictionAsEnum);
 	element->GetInputValue(&n,gauss,MaterialsRheologyNEnum);
-	element->GetInputValue(&thickness, gauss,ThicknessEnum);
-	element->GetInputValue(&base, gauss,BaseEnum);
-	element->GetInputValue(&sealevel, gauss,SealevelEnum);
-	IssmDouble rho_water   = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
-	IssmDouble rho_ice     = element->GetMaterialParameter(MaterialsRhoIceEnum);
-	IssmDouble gravity     = element->GetMaterialParameter(ConstantsGEnum);
-	element->parameters->FindParam(&CoupledFlag,FrictionCouplingEnum);
-
-	switch(CoupledFlag){
-		case 0:
-			Neff=gravity*(rho_ice*thickness+rho_water*(base-sealevel));
-			break;
-		case 1:
-			element->GetInputValue(&Neff,gauss,FrictionEffectivePressureEnum);
-			break;
-		case 2:
-			element->GetInputValue(&Neff,gauss,EffectivePressureTimeAverageEnum);
-			break;
-		default:
-			_error_("not supported");
-	}
-
-	if(Neff<0)Neff=0;
-
+
+	/*Get effective pressure*/
+	IssmDouble Neff = EffectivePressure(gauss);
+
+	/*Compute velocity magnitude*/
 	switch(dim){
 		case 1:
@@ -575,9 +489,6 @@
 
 	/*diverse: */
-	int         CoupledFlag;
 	IssmDouble  r,s;
 	IssmDouble  drag_p, drag_q;
-	IssmDouble  Neff;
-	IssmDouble  thickness,base,sealevel;
 	IssmDouble  vx,vy,vz,vmag;
 	IssmDouble  drag_coefficient;
@@ -587,32 +498,14 @@
 	element->GetInputValue(&drag_p,FrictionPEnum);
 	element->GetInputValue(&drag_q,FrictionQEnum);
-	element->GetInputValue(&thickness, gauss,ThicknessEnum);
-	element->GetInputValue(&base, gauss,BaseEnum);
-	element->GetInputValue(&sealevel, gauss,SealevelEnum);
 	element->GetInputValue(&drag_coefficient, gauss,FrictionCoefficientEnum);
-	IssmDouble rho_water   = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
-	IssmDouble rho_ice     = element->GetMaterialParameter(MaterialsRhoIceEnum);
-	IssmDouble gravity     = element->GetMaterialParameter(ConstantsGEnum);
-	element->parameters->FindParam(&CoupledFlag,FrictionCouplingEnum);
-	//compute r and q coefficients: */
+
+	/*compute r and q coefficients: */
 	r=drag_q/drag_p;
 	s=1./drag_p;
 
-	//From base and thickness, compute effective pressure when drag is viscous, or get Neff from forcing:
-	switch(CoupledFlag){
-		case 0:
-			Neff=gravity*(rho_ice*thickness+rho_water*(base-sealevel));
-			break;
-		case 1:
-			element->GetInputValue(&Neff,gauss,FrictionEffectivePressureEnum);
-			break;
-		case 2:
-			element->GetInputValue(&Neff,gauss,EffectivePressureTimeAverageEnum);
-			break;
-		default:
-			_error_("not supported");
-	}
-	if(Neff<0)Neff=0;
-
+	/*Get effective pressure*/
+	IssmDouble Neff = EffectivePressure(gauss);
+
+	/*Compute velocity magnitude*/
 	switch(dim){
 		case 1:
@@ -777,2 +670,47 @@
 	*palpha2=alpha2;
 }/*}}}*/
+IssmDouble Friction::EffectivePressure(Gauss* gauss){/*{{{*/
+	/*Get effective pressure as a function of 
+	 *  - coupling
+	 *  - type
+	 * Neff=rho_ice*g*thickness+rho_ice*g*base
+	 */
+
+
+	/*diverse: */
+	int         coupled_flag;
+	IssmDouble  Neff;
+
+	/*Recover parameters: */
+	element->parameters->FindParam(&coupled_flag,FrictionCouplingEnum);
+
+	/*From base and thickness, compute effective pressure when drag is viscous, or get Neff from forcing:*/
+	switch(coupled_flag){
+		case 0:{
+		   IssmDouble  thickness,base,sealevel;
+			element->GetInputValue(&thickness, gauss,ThicknessEnum);
+			element->GetInputValue(&base, gauss,BaseEnum);
+			element->GetInputValue(&sealevel, gauss,SealevelEnum);
+			IssmDouble rho_water = element->GetMaterialParameter(MaterialsRhoSeawaterEnum);
+			IssmDouble rho_ice   = element->GetMaterialParameter(MaterialsRhoIceEnum);
+			IssmDouble gravity   = element->GetMaterialParameter(ConstantsGEnum);
+			Neff=gravity*(rho_ice*thickness+rho_water*(base-sealevel));
+				 }
+			break;
+		case 1:
+			element->GetInputValue(&Neff,gauss,FrictionEffectivePressureEnum);
+			break;
+		case 2:
+			element->GetInputValue(&Neff,gauss,EffectivePressureTimeAverageEnum);
+			break;
+		default:
+			_error_("not supported");
+	}
+
+	/*Make sure Neff is positive*/
+	if(Neff<0.) Neff=0.;
+
+	/*Return effective pressure*/
+	return Neff;
+
+}/*}}}*/
Index: /issm/trunk-jpl/src/c/classes/Loads/Friction.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Loads/Friction.h	(revision 22846)
+++ /issm/trunk-jpl/src/c/classes/Loads/Friction.h	(revision 22847)
@@ -43,4 +43,6 @@
 		void  GetAlpha2Weertman(IssmDouble* palpha2,Gauss* gauss);
 		void  GetAlpha2WeertmanTemp(IssmDouble* palpha2,Gauss* gauss);
+
+		IssmDouble EffectivePressure(Gauss* gauss);
 };