fix a bug in converting Budd friction coefficient to Weertman. By def…

…ault, Weertman law has m=1/3
ISSMteam · Jun 27, 2024 · d8083d1 · d8083d1
1 parent 2341c36
commit d8083d1
Show file tree

Hide file tree

Showing 2 changed files with 6 additions and 6 deletions.
diff --git a/pinnicle/modeldata/issm_data.py b/pinnicle/modeldata/issm_data.py
@@ -53,22 +53,22 @@ def load_data(self):
         self.data_dict['s'] = md['geometry']['surface']
         self.data_dict['a'] = (md['smb']['mass_balance'] - md['balancethickness']['thickening_rate'])/self.yts
         self.data_dict['H'] = md['geometry']['thickness']
+        self.data_dict['B'] = md['materials']['rheology_B']
+        self.data_dict['vel'] = np.sqrt(self.data_dict['u']**2.0+self.data_dict['v']**2.0)
         # check the friction law
         if 'C' in md['friction']:
             self.data_dict['C'] = md['friction']['C'] # Weertman
         else:
-            # convert Budd to Weertman type friction coefficient
+            # convert Budd to Weertman type friction coefficient (m=1/3 by default)
             C_b = md['friction']['coefficient'] # Budd
             rho_ice = md['materials']['rho_ice']
             rho_w = md['materials']['rho_water']
             g = md['constants']['g']
             base = md['geometry']['base']
             N = rho_ice*g*self.data_dict['H'] + rho_w*g*base
-            N[np.where(N <= 0, True, False)] = 1
-            self.data_dict['C'] = C_b*np.sqrt(N)
+            N[np.where(N <= 1.0, True, False)] = 1.0
+            self.data_dict['C'] = C_b*np.sqrt(N)*(self.data_dict['vel']**(1.0/3.0))
 
-        self.data_dict['B'] = md['materials']['rheology_B']
-        self.data_dict['vel'] = np.sqrt(self.data_dict['u']**2.0+self.data_dict['v']**2.0)
         # clean up is any of the keys are empty
         self.data_dict = {k:self.data_dict[k] for k in self.data_dict if self.data_dict[k].shape != ()}
         # ice mask

diff --git a/pinnicle/physics/constants.py b/pinnicle/physics/constants.py
@@ -10,7 +10,7 @@ def __init__(self, **kwargs):
 
         # Typical range of variables
         self.variable_lb = {'u': -1.0e4/self.yts, 'v':-1.0e4/self.yts, 's':-1.0e3, 'H':10.0,   'C':0.01,  'a': -5.0/self.yts}
-        self.variable_ub = {'u': 1.0e4/self.yts,  'v':1.0e4/self.yts,  's':2.5e3,  'H':2500.0, 'C':1.0e4, 'a':  5.0/self.yts}
+        self.variable_ub = {'u': 1.0e4/self.yts,  'v':1.0e4/self.yts,  's':3.6e3,  'H':3500.0, 'C':1.0e4, 'a':  5.0/self.yts}
         # add more if needed
         self.variable_lb['u_base'] = -1.0e4/self.yts
         self.variable_ub['u_base'] =  1.0e4/self.yts