Changeset 23841

Timestamp:

04/11/19 15:51:27 (6 years ago)

Author:

Mathieu Morlighem

Message:

BUG: fixing alpha2 <0 in coulomb friction law

Location:

issm/trunk-jpl

Files:

• src/c/classes/Loads/Friction.cpp (modified) (1 diff)
• test/Archives/Archive408.arch (modified) ( previous)
• test/Archives/Archive443.arch (modified) ( previous)

issm/trunk-jpl/src/c/classes/Loads/Friction.cpp

@@ -241 +241 @@
 
         /*Check to prevent dividing by zero if vmag==0*/
-        if(vmag==0. && (s-1.)<0.) alpha2=0.;
-        else alpha2=drag_coefficient*drag_coefficient*pow(Neff,r)*pow(vmag,(s-1.));
+        if(vmag==0. && (s-1.)<0.){
+                alpha2=0.;
+        }
+        else{
+                alpha2=drag_coefficient*drag_coefficient*pow(Neff,r)*pow(vmag,(s-1.));
+        }
 
         floatation_thickness=0;
         if(base<0) floatation_thickness=-(rho_water/rho_ice)*base;
-        if(vmag==0.) alpha2_coulomb=0.;
-        else alpha2_coulomb=drag_coefficient_coulomb*drag_coefficient_coulomb*rho_ice*gravity*(thickness-floatation_thickness)/vmag;
+        if(vmag==0.){
+                alpha2_coulomb=0.;
+        }
+        else{
+                alpha2_coulomb=drag_coefficient_coulomb*drag_coefficient_coulomb*rho_ice*gravity*max(0.,thickness-floatation_thickness)/vmag;
+        }
 
         if(alpha2_coulomb<alpha2) alpha2=alpha2_coulomb;