Index: issm/trunk-jpl/src/m/classes/calvingdev2.m
===================================================================
--- issm/trunk-jpl/src/m/classes/calvingdev2.m	(revision 26660)
+++ issm/trunk-jpl/src/m/classes/calvingdev2.m	(revision 26661)
@@ -42,6 +42,6 @@
 			if (~strcmp(solution,'TransientSolution') | md.transient.ismovingfront==0), return; end
 
-			md = checkfield(md,'fieldname','calving.stress_threshold_groundedice','>',0,'nan',1,'Inf',1);
-			md = checkfield(md,'fieldname','calving.stress_threshold_floatingice','>',0,'nan',1,'Inf',1);
+			md = checkfield(md,'fieldname','calving.stress_threshold_groundedice','>',0,'NaN',1,'Inf',1);
+			md = checkfield(md,'fieldname','calving.stress_threshold_floatingice','>',0,'NaN',1,'Inf',1);
 			md = checkfield(md,'fieldname','calving.height_above_floatation','<=',0);
 		end % }}}
Index: issm/trunk-jpl/src/m/classes/calvingvonmises.m
===================================================================
--- issm/trunk-jpl/src/m/classes/calvingvonmises.m	(revision 26660)
+++ issm/trunk-jpl/src/m/classes/calvingvonmises.m	(revision 26661)
@@ -44,6 +44,6 @@
 			if (~strcmp(solution,'TransientSolution') | md.transient.ismovingfront==0), return; end
 
-			md = checkfield(md,'fieldname','calving.stress_threshold_groundedice','>',0,'nan',1,'Inf',1);
-			md = checkfield(md,'fieldname','calving.stress_threshold_floatingice','>',0,'nan',1,'Inf',1);
+			md = checkfield(md,'fieldname','calving.stress_threshold_groundedice','>',0,'NaN',1,'Inf',1);
+			md = checkfield(md,'fieldname','calving.stress_threshold_floatingice','>',0,'NaN',1,'Inf',1);
 			md = checkfield(md,'fieldname','calving.min_thickness','>=',0,'NaN',1,'Inf',1,'numel',1);
 		end % }}}
Index: issm/trunk-jpl/src/m/coordsystems/laea.py
===================================================================
--- issm/trunk-jpl/src/m/coordsystems/laea.py	(revision 26660)
+++ issm/trunk-jpl/src/m/coordsystems/laea.py	(revision 26661)
@@ -1,5 +1,4 @@
 def laea(lat, long): #{{{
-    '''
-    LAEA - Lambert Azimuthal Equal Area projection at lat, long projection 
+    """LAEA - Lambert Azimuthal Equal Area projection at lat, long projection 
     center.
 
@@ -10,5 +9,5 @@
             string = laea(45, -90)
             return string = '+proj=laea +lat_0=45 +lon_0=-90 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs'
-    '''
+    """
 
     return '+proj=laea +lat_0={} +lon_0={} +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs'.format(lat, long)
Index: issm/trunk-jpl/src/m/exp/expwrite.py
===================================================================
--- issm/trunk-jpl/src/m/exp/expwrite.py	(revision 26660)
+++ issm/trunk-jpl/src/m/exp/expwrite.py	(revision 26661)
@@ -3,30 +3,29 @@
 
 def expwrite(contours, filename):
-    """
-    EXPWRITE - write an Argus file from a dictionary given in input
+    """EXPWRITE - write an Argus file from a dictionary given in input
 
-       This routine writes an Argus file from a dict containing the fields:
-       x and y of the coordinates of the points.
-       The first argument is the list containing the points coordinates
-       and the second one the file to be written.
+    This routine writes an Argus file from a dict containing the fields:
+    x and y of the coordinates of the points.
+    The first argument is the list containing the points coordinates and the 
+    second one the file to be written.
 
-       Usage:
-          expwrite(contours, filename)
+    Usage:
+        expwrite(contours, filename)
 
-       Example:
-          expwrite(coordstruct, 'domainoutline.exp')
+    Example:
+        expwrite(coordstruct, 'domainoutline.exp')
 
-       See also EXPDOC, EXPREAD, EXPWRITEASVERTICES
+    See also EXPDOC, EXPREAD, EXPWRITEASVERTICES
     """
 
     fid = open(filename, 'w')
     for x, y in zip(contours['x'], contours['y']):
-        #if np.size(contour['x']) != np.size(contour['y']):
         if len(x) != len(y):
-            raise RuntimeError("contours x and y coordinates must be of identical size")
+            raise RuntimeError('contours x and y coordinates must be of identical size')
+
         if 'name' in contours:
-            fid.write("%s%s\n" % ('  # Name:', contours['name']))
+            fid.write('{}{}\n'.format('# Name:', contours['name']))
         else:
-            fid.write("%s%s\n" % ('  # Name:', filename))
+            fid.write('{}{}\n'.format('# Name:', filename))
 
         #Add density if it's not there FIXME what is this ever used for?
@@ -35,13 +34,17 @@
         density = 1
 
-        fid.write("%s\n" % '  # Icon:0')
-        fid.write("%s\n" % '  # Points Count Value')
-    #fid.write("%i %f\n" % (np.size(contour['x']), contour['density']))
-        fid.write("%i %f\n" % (np.size(x), density))
-        fid.write("%s\n" % '  # X pos Y pos')
-    #for x, y in zip(contour['x'], contour['y']):
+        fid.write('{}\n'.format('## Icon:0'))
+        fid.write('{}\n'.format('# Points Count Value')
+        if 'density' in contours:
+            if isinstance(contours['density'], int):
+                fid.write('{} {}\n'.format(np.size(x), density))
+            else:
+                fid.write('{} {}\n'.format(np.size(x), 1.))
+        else:
+            fid.write('{} {}\n'.format(np.size(x), 1.))
+        fid.write('{}\n'.format('# X pos Y pos'))
         for xi, yi in zip(x, y):
-            fid.write("%10.10f %10.10f\n" % (xi, yi))
-        fid.write("\n")
+            fid.write('%10.10f %10.10f\n' % (xi, yi))
+        fid.write('\n')
 
     fid.close()
Index: issm/trunk-jpl/src/m/modules/ExpToLevelSet.py
===================================================================
--- issm/trunk-jpl/src/m/modules/ExpToLevelSet.py	(revision 26660)
+++ issm/trunk-jpl/src/m/modules/ExpToLevelSet.py	(revision 26661)
@@ -30,5 +30,5 @@
             contourname = shpread(contourname)
 
-    # NOTE: This library does not currently exist! See TODO list in function 
+    # NOTE: This module does not currently exist! See TODO list in function 
     #       header.
     distance = ExpToLevelSet_python(x, y, contourname)
Index: issm/trunk-jpl/src/m/shp/shp2exp.m
===================================================================
--- issm/trunk-jpl/src/m/shp/shp2exp.m	(revision 26660)
+++ issm/trunk-jpl/src/m/shp/shp2exp.m	(revision 26661)
@@ -1,4 +1,4 @@
 function shp2exp(shpfilename,expfilename)
-%SHP2EXP- transform shape file to Argus .exp file
+%SHP2EXP - transform shape file to Argus .exp file
 %
 %   Usage:
Index: issm/trunk-jpl/src/m/shp/shp2exp.py
===================================================================
--- issm/trunk-jpl/src/m/shp/shp2exp.py	(revision 26660)
+++ issm/trunk-jpl/src/m/shp/shp2exp.py	(revision 26661)
@@ -1,12 +1,12 @@
+import os
 import shapefile
-import os
 from expwrite import expwrite
 
 
-def shp2exp(shapefilename, * expfilename):
-    '''
-    Convert a shapefile to an .exp file.  Optionally, expfilename can be
-    specified to give a name for the .exp file to be created, otherwise the
-    .exp file will have the same prefix as the .shp file.
+def shp2exp(shapefilename, *expfilename):
+    """SHP2EXP - Convert a shapefile to an Argus .exp file. Optionally, 
+    expfilename can be specified to give a name for the .exp file to be 
+    created, otherwise the .exp file will have the same prefix as the .shp 
+    file.
 
     Usage:
@@ -15,10 +15,12 @@
 
     Examples:
-        shp2exp('Domain.shp') % creates Domain.exp
-        shp2exp('Domain.shp', 'DomainForISSM.exp')
-    '''
+        shp2exp('Domain.shp') # Creates Domain.exp
+        shp2exp('Domain.shp', 'Domain.exp')
+
+    See also EXPMASTER, EXPDOC
+    """
 
     if not os.path.exists(shapefilename):
-        raise IOError("shp2exp error message: file '%s' not found!" % shapefilename)
+        raise IOError('shp2exp error message: file {} not found!'.format(shapefilename))
     if not len(expfilename):
         expfile = os.path.splitext(shapefilename)[0] + '.exp'
@@ -33,5 +35,5 @@
     for i in range(len(shp.shapes())):
         geom = shp.shapes()[i].shapeType
-        if geom == 5:  # polygon
+        if geom == 5: # polygon
             expdict['closed'] = 1
             tmpx = [p[0] for p in shp.shapes()[i].points]
@@ -39,5 +41,5 @@
             x.append(tmpx)
             y.append(tmpy)
-        elif geom == 3:  # line
+        elif geom == 3: # line
             expdict['closed'] = 0
             tmpx = [p[0] for p in shp.shapes()[i].points]
@@ -45,5 +47,5 @@
             x.append(tmpx)
             y.append(tmpy)
-        elif geom == 1:  # point
+        elif geom == 1: # point
             expdict['closed'] = 0
             x.append(shp.shapes()[i].points[0][0])
Index: issm/trunk-jpl/src/m/solve/loadresultsfromcluster.py
===================================================================
--- issm/trunk-jpl/src/m/solve/loadresultsfromcluster.py	(revision 26660)
+++ issm/trunk-jpl/src/m/solve/loadresultsfromcluster.py	(revision 26661)
@@ -22,4 +22,9 @@
         md = loadresultsfromcluster(md)
         md = loadresultsfromcluster(md, 'runtimename', runtimename)
+
+        Options include: 'runtimename', 'nolog'
+
+    Example:
+        md = loadresultsfromcluster(md, 'runtimename', 'test101-06-15-2021-13-24-18-4883')
     """
 
Index: issm/trunk-jpl/src/m/solvers/asmoptions.m
===================================================================
--- issm/trunk-jpl/src/m/solvers/asmoptions.m	(revision 26660)
+++ issm/trunk-jpl/src/m/solvers/asmoptions.m	(revision 26661)
@@ -1,4 +1,4 @@
 function asm=asmoptions(varargin)
-%ASMOPTIONS - return Additive Schwartz Method petsc options
+%ASMOPTIONS - return Additive Schwartz Method PETSc options
 %
 %   Usage:
Index: issm/trunk-jpl/src/m/solvers/asmoptions.py
===================================================================
--- issm/trunk-jpl/src/m/solvers/asmoptions.py	(revision 26660)
+++ issm/trunk-jpl/src/m/solvers/asmoptions.py	(revision 26661)
@@ -3,22 +3,25 @@
 
 def asmoptions(*args):
-    #ASMOPTIONS - return ASM petsc options
-    #
-    #   Usage:
-    #      options = asmoptions
+    """ASMOPTIONS - Return Additive Schwartz Method PETSc options
 
-    #retrieve options provided in *args
+    Usage:
+        options = asmoptions
+    """
+
+    # Retrieve options provided in *args
     arguments = pairoptions(*args)
 
-    options = [['toolkit', 'petsc'],
-               ['mat_type', 'mpiaij'],
-               ['ksp_type', 'gmres'],
-               ['pc_type', 'asm'],
-               ['sub_pc_type', 'lu'],
-               ['pc_asm_overlap', 3],
-               ['ksp_max_it', 100],
-               ['ksp_rtol', 1e-30]]
+    options = [
+        ['toolkit', 'petsc'],
+        ['mat_type', 'mpiaij'],
+        ['ksp_type', 'gmres'],
+        ['pc_type', 'asm'],
+        ['sub_pc_type', 'lu'],
+        ['pc_asm_overlap', 3],
+        ['ksp_max_it', 100],
+        ['ksp_rtol', 1e-30]
+    ]
 
-    #now, go through our arguments, and write over default options.
+    # Now, go through our arguments, and write over default options
     for i in range(len(arguments.list)):
         arg1 = arguments.list[i][0]
@@ -33,5 +36,5 @@
                 break
         if not found:
-            #this option did not exist, add it:
+            # This option did not exist; add it
             options.append([arg1, arg2])
 
Index: issm/trunk-jpl/src/m/solvers/asmstokesoptions.m
===================================================================
--- issm/trunk-jpl/src/m/solvers/asmstokesoptions.m	(revision 26660)
+++ issm/trunk-jpl/src/m/solvers/asmstokesoptions.m	(revision 26661)
@@ -1,7 +1,7 @@
-function asm=asmoptions(varargin)
-%ASMOPTIONS - return Additive Schwartz Method petsc options
+function asm=asmstokesoptions(varargin)
+%ASMSTOKESOPTIONS - return Additive Schwartz Method Stokes PETSc options
 %
 %   Usage:
-%      options=asmoptions;
+%      options=asmstokesoptions;
 
 %retrieve options provided in varargin
Index: issm/trunk-jpl/src/m/solvers/asmstokesoptions.py
===================================================================
--- issm/trunk-jpl/src/m/solvers/asmstokesoptions.py	(revision 26661)
+++ issm/trunk-jpl/src/m/solvers/asmstokesoptions.py	(revision 26661)
@@ -0,0 +1,46 @@
+from collections import OrderedDict
+from pairoptions import pairoptions
+
+def asmstokesoptions(*args):
+    """ASMSTOKESOPTIONS - Return Additive Schwartz Method Stokes PETSc options
+
+    Usage:
+        options = asmstokesoptions
+    """
+
+    # Retrieve options provided in *args
+    arguments = pairoptions(*args)
+
+    options = [
+        ['toolkit', 'petsc'],
+        ['mat_type', 'mpiaij'],
+        ['ksp_type', 'gmres'],
+        ['pc_type', 'asm'],
+        ['sub_pc_type', 'lu'],
+        ['pc_asm_overlap', 1], # COMSOL's default
+        ['ksp_max_it', 100],
+        ['ksp_rtol', 1e-7], # Tuned for best performance and to fit ISMIP-HOM-C 5km with MUMPS
+        ['ksp_atol', 1e-10] # Tuned for best performance and to fit ISMIP-HOM-C 5km with MUMPS
+    ]
+
+    # Now, go through our arguments, and write over default options
+    for i in range(len(arguments.list)):
+        arg1 = arguments.list[i][0]
+        arg2 = arguments.list[i][1]
+        found = 0
+        for j in range(len(options)):
+            joption = options[j][0]
+            if joption == arg1:
+                joption[1] = arg2
+                options[j] = joption
+                found = 1
+                break
+        if not found:
+            # This option did not exist; add it
+            options.append([arg1, arg2])
+
+    asmoptions = OrderedDict()
+    for j in range(len(options)):
+        asmoptions[options[j][0]]=options[j][1]
+
+    return asmoptions
