Index: /issm/trunk-jpl/src/m/exp/contourlevelzero.py
===================================================================
--- /issm/trunk-jpl/src/m/exp/contourlevelzero.py	(revision 27949)
+++ /issm/trunk-jpl/src/m/exp/contourlevelzero.py	(revision 27950)
@@ -3,229 +3,222 @@
 from collections import OrderedDict
 
+
 def contourlevelzero(md,mask,level):
-    """CONTOURLEVELZERO - figure out the zero level (or offset thereof, specified by the level value)
-                       of a vectorial mask, and vectorialize it into an exp or shp compatible structure.
-    
-       Usage:
-          contours=contourlevelzero(md,mask,level)
-    
-       See also: PLOT_CONTOUR
+    """CONTOURLEVELZERO - figure out the zero level (or offset thereof, 
+    specified by the level value of a vectorial mask, and vectorialize it into 
+    an exp or shp compatible structure.
+
+    Usage:
+        contours=contourlevelzero(md,mask,level)
+
+        See also: PLOT_CONTOUR
     """
-    
-    #process data 
-    if md.mesh.dimension()==3:
+
+    # Process data 
+    if md.mesh.dimension() == 3:
         x = md.mesh.x2d
         y = md.mesh.y2d
-        z=md.mesh.z
-        index=md.mesh.elements2d-1
+        z = md.mesh.z
+        index = md.mesh.elements2d - 1
     else:
-        x=md.mesh.x
-        y=md.mesh.y
-        index=md.mesh.elements-1
-        z=np.zeros((md.mesh.numberofvertices,1))
-        
-    if len(mask)==0:
+        x = md.mesh.x
+        y = md.mesh.y
+        index = md.mesh.elements - 1
+        z = np.zeros((md.mesh.numberofvertices, 1))
+
+    if len(mask) == 0:
         raise OSError("mask provided is empty")
-    
-    if md.mesh.dimension()==3:
-        if len(mask)!=md.mesh.numberofvertices2d: 
+
+    if md.mesh.dimension() == 3:
+        if len(mask) != md.mesh.numberofvertices2d: 
             raise OSError("mask provided should be specified at the vertices of the mesh")
     else:
-        if len(mask)!=md.mesh.numberofvertices:
+        if len(mask) != md.mesh.numberofvertices:
             raise OSError("mask provided should be specified at the vertices of the mesh")
+
+    # Initialization of some variables
+    numberofelements = np.size(index, 0)
+    elementslist = np.c_[0:numberofelements]
+    c = []
+    h = []
+
+    # Get unique edges in mesh
+    # 1: list of edges
+    edges = np.vstack((np.vstack((index[:, (0, 1)], index[:, (1, 2)])), index[:, (2, 0)]))
+    # 2: find unique edges
+    [edges, J] = np.unique(np.sort(edges, 1), axis=0, return_inverse=True)
+    # 3: unique edge numbers
+    vec = J
+    # 4: unique edges numbers in each triangle (2 triangles sharing the same 
+    # edge will have the same edge number)
+    edges_tria = np.hstack((np.hstack((vec[elementslist], vec[elementslist + numberofelements])), vec[elementslist + 2 * numberofelements]))
+
+    # Segments [nodes1 nodes2]
+    Seg1 = index[:, (0, 1)]
+    Seg2 = index[:, (1, 2)]
+    Seg3 = index[:, (2, 0)]
+
+    # Segment numbers [1;4;6;...]
+    Seg1_num = edges_tria[:, 0]
+    Seg2_num = edges_tria[:, 1]
+    Seg3_num = edges_tria[:, 2]
+
+    #value of data on each tips of the segments
+    Data1 = mask[Seg1]
+    Data2 = mask[Seg2]
+    Data3 = mask[Seg3]
+
+    # Get the ranges for each segment
+    Range1 = np.sort(Data1, 1)
+    Range2 = np.sort(Data2, 1)
+    Range3 = np.sort(Data3, 1)
+
+    # Find the segments that contain this value
+    pos1 = (Range1[:, 0] < level) & (Range1[:, 1] >= level)
+    pos2 = (Range2[:, 0] < level) & (Range2[:, 1] >= level)
+    pos3 = (Range3[:, 0] < level) & (Range3[:, 1] >= level)
+
+    # Get elements
+    poselem12 = (pos1) & (pos2)
+    poselem13 = (pos1) & (pos3)
+    poselem23 = (pos2) & (pos3)
+    poselem = np.where((poselem12) | (poselem13) | (poselem23))
+    poselem = poselem[0]
+    numelems = len(poselem)
+
+    # If no element has been flagged, skip to the next level
+    if numelems == 0:
+        raise Exception('contourlevelzero warning message: no elements found with corresponding level value in mask')
+        contours = []
+        return contours
+
+    # Go through the elements and build the coordinates for each segment (1 by element)
+    x1 = np.zeros((numelems, 1))
+    x2 = np.zeros((numelems, 1))
+    y1 = np.zeros((numelems, 1))
+    y2 = np.zeros((numelems, 1))
+    z1 = np.zeros((numelems, 1))
+    z2 = np.zeros((numelems, 1))
+
+    edge_l = np.zeros((numelems, 2))
+
+    for j in range(0, numelems):
+        with np.errstate(divide='ignore', invalid='ignore'):
+            weight1 = np.divide(level - Data1[poselem[j], 0],Data1[poselem[j], 1] - Data1[poselem[j], 0])
+            weight2 = np.divide(level - Data2[poselem[j], 0],Data2[poselem[j], 1] - Data2[poselem[j], 0])
+            weight3 = np.divide(level - Data3[poselem[j], 0],Data3[poselem[j], 1] - Data3[poselem[j], 0])
+
+        if poselem12[poselem[j]] == True:
+            x1[j] = x[Seg1[poselem[j], 0]] + weight1 * [x[Seg1[poselem[j], 1]] - x[Seg1[poselem[j], 0]]]
+            x2[j] = x[Seg2[poselem[j], 0]] + weight2 * [x[Seg2[poselem[j], 1]] - x[Seg2[poselem[j], 0]]]
+            y1[j] = y[Seg1[poselem[j], 0]] + weight1 * [y[Seg1[poselem[j], 1]] - y[Seg1[poselem[j], 0]]]
+            y2[j] = y[Seg2[poselem[j], 0]] + weight2 * [y[Seg2[poselem[j], 1]] - y[Seg2[poselem[j], 0]]]
+            z1[j] = z[Seg1[poselem[j], 0]] + weight1 * [z[Seg1[poselem[j], 1]] - z[Seg1[poselem[j], 0]]]
+            z2[j] = z[Seg2[poselem[j], 0]] + weight2 * [z[Seg2[poselem[j], 1]] - z[Seg2[poselem[j], 0]]]
+
+            edge_l[j, 0] = Seg1_num[poselem[j]]
+            edge_l[j, 1] = Seg2_num[poselem[j]]
+        elif poselem13[poselem[j]] == True:
+            x1[j] = x[Seg1[poselem[j], 0]] + weight1 * [x[Seg1[poselem[j], 1]] - x[Seg1[poselem[j], 0]]]
+            x2[j] = x[Seg3[poselem[j], 0]] + weight3 * [x[Seg3[poselem[j], 1]] - x[Seg3[poselem[j], 0]]]
+            y1[j] = y[Seg1[poselem[j], 0]] + weight1 * [y[Seg1[poselem[j], 1]] - y[Seg1[poselem[j], 0]]]
+            y2[j] = y[Seg3[poselem[j], 0]] + weight3 * [y[Seg3[poselem[j], 1]] - y[Seg3[poselem[j], 0]]]
+            z1[j] = z[Seg1[poselem[j], 0]] + weight1 * [z[Seg1[poselem[j], 1]] - z[Seg1[poselem[j], 0]]]
+            z2[j] = z[Seg3[poselem[j], 0]] + weight3 * [z[Seg3[poselem[j], 1]] - z[Seg3[poselem[j], 0]]]
+
+            edge_l[j, 0] = Seg1_num[poselem[j]]
+            edge_l[j, 1] = Seg3_num[poselem[j]]
+        elif poselem23[poselem[j]] == True:
+            x1[j] = x[Seg2[poselem[j], 0]] + weight2 * [x[Seg2[poselem[j], 1]] - x[Seg2[poselem[j], 0]]]
+            x2[j] = x[Seg3[poselem[j], 0]] + weight3 * [x[Seg3[poselem[j], 1]] - x[Seg3[poselem[j], 0]]]
+            y1[j] = y[Seg2[poselem[j], 0]] + weight2 * [y[Seg2[poselem[j], 1]] - y[Seg2[poselem[j], 0]]]
+            y2[j] = y[Seg3[poselem[j], 0]] + weight3 * [y[Seg3[poselem[j], 1]] - y[Seg3[poselem[j], 0]]]
+            z1[j] = z[Seg2[poselem[j], 0]] + weight2 * [z[Seg2[poselem[j], 1]] - z[Seg2[poselem[j], 0]]]
+            z2[j] = z[Seg3[poselem[j], 0]] + weight3 * [z[Seg3[poselem[j], 1]] - z[Seg3[poselem[j], 0]]]
+
+            edge_l[j, 0] = Seg2_num[poselem[j]]
+            edge_l[j, 1] = Seg3_num[poselem[j]]
+        # else:
+        # Should never get here
+
+    # Now that we have the segments, we must try to connect them...
+
+    # Loop over the subcontours
+    contours = []
+
+    while len(edge_l) > 0:
+        # Take the right edge of the second segment and connect it to the next segments if any
+        e1 = edge_l[0, 0]
+        e2 = edge_l[0, 1]
+        xc = np.vstack((x1[0], x2[0]))
+        yc = np.vstack((y1[0], y2[0]))
+        zc = np.vstack((z1[0], z2[0]))
+        # Erase the lines corresponding to this edge
+        edge_l = np.delete(edge_l, 0, axis=0)
+        x1 = np.delete(x1, 0, axis=0)
+        x2 = np.delete(x2, 0, axis=0)
+        y1 = np.delete(y1, 0, axis=0)
+        y2 = np.delete(y2, 0, axis=0)
+        z1 = np.delete(z1, 0, axis=0)
+        z2 = np.delete(z2,0,axis=0)
+        pos1 = np.where(edge_l == e1)
         
-    #initialization of some variables
-    numberofelements=np.size(index,0)
-    elementslist=np.c_[0:numberofelements]
-    c=[]
-    h=[]
-    
-    #get unique edges in mesh
-    #1: list of edges
-    edges=np.vstack((np.vstack((index[:,(0,1)],index[:,(1,2)])),index[:,(2,0)]))
-
-    #2: find unique edges
-    [edges,J]=np.unique(np.sort(edges,1),axis=0,return_inverse=True)
-    #3: unique edge numbers
-    vec=J
-    #4: unique edges numbers in each triangle (2 triangles sharing the same edge will have
-    #   the same edge number)
-    edges_tria=np.hstack((np.hstack((vec[elementslist],vec[elementslist+numberofelements])),vec[elementslist+2*numberofelements]))
-    
-    #segments [nodes1 nodes2]
-    Seg1=index[:,(0,1)]
-    Seg2=index[:,(1,2)]
-    Seg3=index[:,(2,0)]
-    
-    #segment numbers [1;4;6;...]
-    Seg1_num=edges_tria[:,0]
-    Seg2_num=edges_tria[:,1]
-    Seg3_num=edges_tria[:,2]
-    
-    #value of data on each tips of the segments
-    Data1=mask[Seg1]
-    Data2=mask[Seg2]
-    Data3=mask[Seg3]
-    
-    #get the ranges for each segment
-    Range1=np.sort(Data1,1)
-    Range2=np.sort(Data2,1)
-    Range3=np.sort(Data3,1)
-    
-    #find the segments that contain this value
-    pos1=(Range1[:,0]<level) & (Range1[:,1]>=level)
-    pos2=(Range2[:,0]<level) & (Range2[:,1]>=level)
-    pos3=(Range3[:,0]<level) & (Range3[:,1]>=level)
-    
-    #get elements
-    poselem12=(pos1) & (pos2)
-    poselem13=(pos1) & (pos3)
-    poselem23=(pos2) & (pos3)
-    poselem=np.where((poselem12) | (poselem13) | (poselem23))
-    poselem=poselem[0]
-    numelems=len(poselem)
-    
-    #if no element has been flagged, skip to the next level
-    if numelems==0:
-        raise Exception('contourlevelzero warning message: no elements found with corresponding level value in mask')
-        contours=[]
-        return contours
-    
-    #go through the elements and build the coordinates for each segment (1 by element)
-    x1=np.zeros((numelems,1))
-    x2=np.zeros((numelems,1))
-    y1=np.zeros((numelems,1))
-    y2=np.zeros((numelems,1))
-    z1=np.zeros((numelems,1))
-    z2=np.zeros((numelems,1))
-    
-    edge_l=np.zeros((numelems,2))
-    
-    for j in range(0,numelems):
-        
-        with np.errstate(divide='ignore', invalid='ignore'):
-            weight1=np.divide(level-Data1[poselem[j],0],Data1[poselem[j],1]-Data1[poselem[j],0])
-            weight2=np.divide(level-Data2[poselem[j],0],Data2[poselem[j],1]-Data2[poselem[j],0])
-            weight3=np.divide(level-Data3[poselem[j],0],Data3[poselem[j],1]-Data3[poselem[j],0])
-        
-        if poselem12[poselem[j]]==True:
-            
-            x1[j]=x[Seg1[poselem[j],0]]+weight1*[x[Seg1[poselem[j],1]]-x[Seg1[poselem[j],0]]]
-            x2[j]=x[Seg2[poselem[j],0]]+weight2*[x[Seg2[poselem[j],1]]-x[Seg2[poselem[j],0]]]
-            y1[j]=y[Seg1[poselem[j],0]]+weight1*[y[Seg1[poselem[j],1]]-y[Seg1[poselem[j],0]]]
-            y2[j]=y[Seg2[poselem[j],0]]+weight2*[y[Seg2[poselem[j],1]]-y[Seg2[poselem[j],0]]]
-            z1[j]=z[Seg1[poselem[j],0]]+weight1*[z[Seg1[poselem[j],1]]-z[Seg1[poselem[j],0]]]
-            z2[j]=z[Seg2[poselem[j],0]]+weight2*[z[Seg2[poselem[j],1]]-z[Seg2[poselem[j],0]]]
-            
-            edge_l[j,0]=Seg1_num[poselem[j]]
-            edge_l[j,1]=Seg2_num[poselem[j]]
-        elif poselem13[poselem[j]]==True:
-            
-            x1[j]=x[Seg1[poselem[j],0]]+weight1*[x[Seg1[poselem[j],1]]-x[Seg1[poselem[j],0]]]
-            x2[j]=x[Seg3[poselem[j],0]]+weight3*[x[Seg3[poselem[j],1]]-x[Seg3[poselem[j],0]]]
-            y1[j]=y[Seg1[poselem[j],0]]+weight1*[y[Seg1[poselem[j],1]]-y[Seg1[poselem[j],0]]]
-            y2[j]=y[Seg3[poselem[j],0]]+weight3*[y[Seg3[poselem[j],1]]-y[Seg3[poselem[j],0]]]
-            z1[j]=z[Seg1[poselem[j],0]]+weight1*[z[Seg1[poselem[j],1]]-z[Seg1[poselem[j],0]]]
-            z2[j]=z[Seg3[poselem[j],0]]+weight3*[z[Seg3[poselem[j],1]]-z[Seg3[poselem[j],0]]]
-            
-            edge_l[j,0]=Seg1_num[poselem[j]]
-            edge_l[j,1]=Seg3_num[poselem[j]]
-        elif poselem23[poselem[j]]==True:
-            
-            x1[j]=x[Seg2[poselem[j],0]]+weight2*[x[Seg2[poselem[j],1]]-x[Seg2[poselem[j],0]]]
-            x2[j]=x[Seg3[poselem[j],0]]+weight3*[x[Seg3[poselem[j],1]]-x[Seg3[poselem[j],0]]]
-            y1[j]=y[Seg2[poselem[j],0]]+weight2*[y[Seg2[poselem[j],1]]-y[Seg2[poselem[j],0]]]
-            y2[j]=y[Seg3[poselem[j],0]]+weight3*[y[Seg3[poselem[j],1]]-y[Seg3[poselem[j],0]]]
-            z1[j]=z[Seg2[poselem[j],0]]+weight2*[z[Seg2[poselem[j],1]]-z[Seg2[poselem[j],0]]]
-            z2[j]=z[Seg3[poselem[j],0]]+weight3*[z[Seg3[poselem[j],1]]-z[Seg3[poselem[j],0]]]
-
-            edge_l[j,0]=Seg2_num[poselem[j]]
-            edge_l[j,1]=Seg3_num[poselem[j]]
-
-        #else:
-	    #it shoud not go here
-            
-    #now that we have the segments, we must try to connect them...
-    
-    #loop over the subcontours
-    contours=[]
-    
-    while len(edge_l)>0:
-        
-        #take the right edge of the second segment and connect it to the next segments if any
-        e1=edge_l[0,0]
-        e2=edge_l[0,1]
-        xc=np.vstack((x1[0],x2[0]))
-        yc=np.vstack((y1[0],y2[0]))
-        zc=np.vstack((z1[0],z2[0]))
-        #erase the lines corresponding to this edge
-        edge_l=np.delete(edge_l,0,axis=0)
-        x1=np.delete(x1,0,axis=0)
-        x2=np.delete(x2,0,axis=0)
-        y1=np.delete(y1,0,axis=0)
-        y2=np.delete(y2,0,axis=0)
-        z1=np.delete(z1,0,axis=0)
-        z2=np.delete(z2,0,axis=0)
-        pos1=np.where(edge_l==e1)
-        
-        while len(pos1[0])>0:
-            
-            if np.all(pos1[1]==0):
-                xc=np.vstack((x2[pos1[0]],xc))
-                yc=np.vstack((y2[pos1[0]],yc))
-                zc=np.vstack((z2[pos1[0]],zc))
-                #next edge:
-                e1=edge_l[pos1[0],1]
+        while len(pos1[0]) > 0:
+            if np.all(pos1[1] == 0):
+                xc = np.vstack((x2[pos1[0]], xc))
+                yc = np.vstack((y2[pos1[0]], yc))
+                zc = np.vstack((z2[pos1[0]], zc))
+                # Next edge:
+                e1 = edge_l[pos1[0], 1]
             else:
-                xc=np.vstack((x1[pos1[0]],xc))
-                yc=np.vstack((y1[pos1[0]],yc))
-                zc=np.vstack((z1[pos1[0]],zc))
-                #next edge:
-                e1=edge_l[pos1[0],0]
-                
-            #erase the lines of this
-            edge_l=np.delete(edge_l,pos1[0],axis=0)
-            x1=np.delete(x1,pos1[0],axis=0)
-            x2=np.delete(x2,pos1[0],axis=0)
-            y1=np.delete(y1,pos1[0],axis=0)
-            y2=np.delete(y2,pos1[0],axis=0)
-            z1=np.delete(z1,pos1[0],axis=0)
-            z2=np.delete(z2,pos1[0],axis=0)
-            #next connection
-            pos1=np.where(edge_l==e1)
-            
-        #same thing the other way (to the right)
-        pos2=np.where(edge_l==e2)
-
-        while len(pos2[0])>0:
-            
-            if np.all(pos2[1]==0):
-                xc=np.vstack((xc,x2[pos2[0]]))
-                yc=np.vstack((yc,y2[pos2[0]]))
-                zc=np.vstack((zc,z2[pos2[0]]))
-                #next edge:
-                e2=edge_l[pos2[0],1]
+                xc = np.vstack((x1[pos1[0]], xc))
+                yc = np.vstack((y1[pos1[0]], yc))
+                zc = np.vstack((z1[pos1[0]], zc))
+                # Next edge:
+                e1 = edge_l[pos1[0], 0]
+
+            # Erase the lines of this
+            edge_l = np.delete(edge_l, pos1[0], axis=0)
+            x1 = np.delete(x1, pos1[0], axis=0)
+            x2 = np.delete(x2, pos1[0], axis=0)
+            y1 = np.delete(y1, pos1[0], axis=0)
+            y2 = np.delete(y2, pos1[0], axis=0)
+            z1 = np.delete(z1, pos1[0], axis=0)
+            z2 = np.delete(z2, pos1[0], axis=0)
+            # Next connection
+            pos1 = np.where(edge_l == e1)
+
+        # Same thing the other way (to the right)
+        pos2 = np.where(edge_l == e2)
+
+        while len(pos2[0]) > 0:
+            if np.all(pos2[1] == 0):
+                xc = np.vstack((xc, x2[pos2[0]]))
+                yc = np.vstack((yc, y2[pos2[0]]))
+                zc = np.vstack((zc, z2[pos2[0]]))
+                # Next edge:
+                e2 = edge_l[pos2[0], 1]
             else:
-                xc=np.vstack((xc,x1[pos2[0]]))
-                yc=np.vstack((yc,y1[pos2[0]]))
-                zc=np.vstack((zc,z1[pos2[0]]))
-                #next edge:
-                e2=edge_l[pos2[0],0]
-                
-            #erase the lines of this
-            edge_l=np.delete(edge_l,pos2[0],axis=0)
-            x1=np.delete(x1,pos2[0],axis=0)
-            x2=np.delete(x2,pos2[0],axis=0)
-            y1=np.delete(y1,pos2[0],axis=0)
-            y2=np.delete(y2,pos2[0],axis=0)
-            z1=np.delete(z1,pos2[0],axis=0)
-            z2=np.delete(z2,pos2[0],axis=0)
-            #next connection
-            pos2=np.where(edge_l==e2)
-            
-        #save xc,yc contour: 
+                xc = np.vstack((xc, x1[pos2[0]]))
+                yc = np.vstack((yc, y1[pos2[0]]))
+                zc = np.vstack((zc, z1[pos2[0]]))
+                # Next edge:
+                e2 = edge_l[pos2[0], 0]
+
+            # Erase the lines of this
+            edge_l = np.delete(edge_l, pos2[0], axis=0)
+            x1 = np.delete(x1, pos2[0], axis=0)
+            x2 = np.delete(x2, pos2[0], axis=0)
+            y1 = np.delete(y1, pos2[0], axis=0)
+            y2 = np.delete(y2, pos2[0], axis=0)
+            z1 = np.delete(z1, pos2[0], axis=0)
+            z2 = np.delete(z2, pos2[0], axis=0)
+            # Next connection
+            pos2 = np.where(edge_l == e2)
+
+        # Save xc, yc contour:
         newcontour = OrderedDict()
         newcontour['nods'] = np.size(xc)
-        newcontour['density'] = 1 
+        newcontour['density'] = 1
         newcontour['closed'] = 0
         newcontour['x'] = np.ma.filled(xc.astype(float), np.nan)
@@ -234,4 +227,4 @@
         newcontour['name'] = ''
         contours.append(newcontour)
-        
+
     return contours
Index: /issm/trunk-jpl/src/m/exp/exptool.m
===================================================================
--- /issm/trunk-jpl/src/m/exp/exptool.m	(revision 27949)
+++ /issm/trunk-jpl/src/m/exp/exptool.m	(revision 27950)
@@ -16,5 +16,5 @@
 %      - markersize (default=7)
 %      - markeredgecolor (default='r')
-%      - nofigurecopy (default=0) do not copy current figure, this is needed on some plateform to avoid an offset in the figure
+%      - nofigurecopy (default=0) do not copy current figure, this is needed on some platform to avoid an offset in the figure
 %
 %   Usage:
@@ -139,5 +139,5 @@
 disableDefaultInteractivity(gca); %disables the built-in interactions for the specified axes
 
-%Build backup structre for do and redo
+%Build backup structure for do and redo
 backup=cell(1,3);
 backup{1,1}=A;
Index: /issm/trunk-jpl/src/m/solve/loadresultsfromcluster.py
===================================================================
--- /issm/trunk-jpl/src/m/solve/loadresultsfromcluster.py	(revision 27949)
+++ /issm/trunk-jpl/src/m/solve/loadresultsfromcluster.py	(revision 27950)
@@ -17,5 +17,5 @@
 
 def loadresultsfromcluster(md, *args):  #{{{
-    """LOADRESULTSFROMCLUSTER - load results of solution sequence from cluster
+    """loadresultsfromcluster - load results of solution sequence from cluster
 
     Usage:
Index: /issm/trunk-jpl/src/m/solve/marshall.py
===================================================================
--- /issm/trunk-jpl/src/m/solve/marshall.py	(revision 27949)
+++ /issm/trunk-jpl/src/m/solve/marshall.py	(revision 27950)
@@ -5,5 +5,5 @@
 
 def marshall(md):
-    """MARSHALL - outputs a compatible binary file from @model md, for certain solution type.
+    """marshall - outputs a compatible binary file from @model md, for certain solution type.
 
     The routine creates a compatible binary file from @model md
@@ -15,5 +15,5 @@
 
     if md.verbose.solution:
-        print("marshalling file {}.bin".format(md.miscellaneous.name))
+        print('marshalling file \'{}\'.bin'.format(md.miscellaneous.name))
 
     # Open file for binary writing
@@ -21,5 +21,5 @@
         fid = open(md.miscellaneous.name + '.bin', 'wb')
     except IOError as e:
-        raise IOError("marshall error message: could not open '%s.bin' file for binary writing. Due to: ".format(md.miscellaneous.name), e)
+        print('marshall error message: could not open \'{}.bin\' file for binary writing due to: {}'.format(md.miscellaneous.name, e))
 
     fields = md.properties()
@@ -32,5 +32,5 @@
         # Check that current field is an object
         if not hasattr(getattr(md, field), 'marshall'):
-            raise TypeError("field '{}' is not an object.".format(field))
+            raise TypeError('field \'{}\' is not an object.'.format(field))
 
         # Marshall current object
@@ -46,5 +46,5 @@
 
     except IOError as e:
-        print('marshall error message: could not close file \'{}.bin\' due to:'.format(md.miscellaneous.name), e)
+        print('marshall error message: could not close \'{}.bin\' file for binary writing due to: {}'.format(md.miscellaneous.name, e))
 
     # Uncomment the following to make a copy of the binary input file for 
Index: /issm/trunk-jpl/src/m/solve/solve.py
===================================================================
--- /issm/trunk-jpl/src/m/solve/solve.py	(revision 27949)
+++ /issm/trunk-jpl/src/m/solve/solve.py	(revision 27950)
@@ -11,5 +11,5 @@
 
 def solve(md, solutionstring, *args):
-    """SOLVE - apply solution sequence for this model
+    """solve - apply solution sequence for this model
 
     Usage:
@@ -98,5 +98,5 @@
     if options.getfieldvalue('checkconsistency', 'yes') == 'yes':
         if md.verbose.solution:
-            print("checking model consistency")
+            print('checking model consistency')
         ismodelselfconsistent(md)
 
@@ -112,5 +112,5 @@
             if options.getfieldvalue('runtimename', True):
                 c = datetime.now()
-                md.private.runtimename = "%s-%02i-%02i-%04i-%02i-%02i-%02i-%i" % (md.miscellaneous.name, c.month, c.day, c.year, c.hour, c.minute, c.second, os.getpid())
+                md.private.runtimename = '%s-%02i-%02i-%04i-%02i-%02i-%02i-%i' % (md.miscellaneous.name, c.month, c.day, c.year, c.hour, c.minute, c.second, os.getpid())
             else:
                 md.private.runtimename = md.miscellaneous.name
Index: /issm/trunk-jpl/src/m/solve/solveslm.py
===================================================================
--- /issm/trunk-jpl/src/m/solve/solveslm.py	(revision 27949)
+++ /issm/trunk-jpl/src/m/solve/solveslm.py	(revision 27950)
@@ -10,8 +10,8 @@
 
 def solveslm(slm, solutionstringi, *args):
-    """SOLVESLM - apply solution sequence for this sealevel model
+    """solveslm - apply solution sequence for this sealevel model
 
     Usage:
-        slm=solveslm(slm,solutionstring,varargin)
+        slm = solveslm(slm,solutionstring,varargin)
         where varargin is a lit of paired arguments of string OR enums