Index: /issm/trunk-jpl/src/m/classes/thermal.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/thermal.py	(revision 24566)
+++ /issm/trunk-jpl/src/m/classes/thermal.py	(revision 24567)
@@ -15,5 +15,5 @@
 
     def __init__(self):  # {{{
-        self.spctemperature = float('NaN')
+        self.spctemperature = np.nan
         self.penalty_threshold = 0
         self.stabilization = 0
@@ -107,5 +107,5 @@
             md = checkfield(md, 'fieldname', 'thermal.watercolumn_upperlimit', '>=', 0)
 
-            TEMP = md.thermal.spctemperature[:-1].flatten(-1)
+            TEMP = md.thermal.spctemperature[:-1].flatten()
             pos = np.where(~np.isnan(TEMP))
             try:
@@ -114,7 +114,7 @@
                 spccol = 1
 
-            replicate = np.tile(md.geometry.surface - md.mesh.z, (spccol)).flatten(-1)
+            replicate = np.tile(md.geometry.surface - md.mesh.z, (spccol)).flatten()
             control = md.materials.meltingpoint - md.materials.beta * md.materials.rho_ice * md.constants.g * replicate + 1.0e-5
-            md = checkfield(md, 'fieldname', 'thermal.spctemperature', 'field', md.thermal.spctemperature.flatten(-1)[pos], '<=', control[pos], 'message', "spctemperature should be below the adjusted melting point")
+            md = checkfield(md, 'fieldname', 'thermal.spctemperature', 'field', md.thermal.spctemperature.flatten()[pos], '<=', control[pos], 'message', "spctemperature should be below the adjusted melting point")
             md = checkfield(md, 'fieldname', 'thermal.isenthalpy', 'numel', [1], 'values', [0, 1])
             md = checkfield(md, 'fieldname', 'thermal.isdynamicbasalspc', 'numel', [1], 'values', [0, 1])
Index: /issm/trunk-jpl/src/m/mesh/roundmesh.py
===================================================================
--- /issm/trunk-jpl/src/m/mesh/roundmesh.py	(revision 24566)
+++ /issm/trunk-jpl/src/m/mesh/roundmesh.py	(revision 24567)
@@ -17,11 +17,9 @@
           md = roundmesh(md, radius, resolution)
     """
+    # First we have to create the domain outline
+    # Get number of points on the circle
+    pointsonedge = int(np.floor((2. * np.pi * radius) / resolution) + 1)  # + 1 to close the outline
 
-    #First we have to create the domain outline
-
-    #Get number of points on the circle
-    pointsonedge = np.floor((2. * np.pi * radius) / resolution) + 1  # + 1 to close the outline
-
-    #Calculate the cartesians coordinates of the points
+    # Calculate the cartesians coordinates of the points
     theta = np.linspace(0., 2. * np.pi, pointsonedge)
     x_list = roundsigfig(radius * np.cos(theta), 12)
@@ -33,14 +31,14 @@
     expwrite(A, 'RoundDomainOutline.exp')
 
-    #Call Bamg
+    # Call Bamg
     md = triangle(md, 'RoundDomainOutline.exp', resolution)
-    #md = bamg(md, 'domain', 'RoundDomainOutline.exp', 'hmin', resolution)
+    # md = bamg(md, 'domain', 'RoundDomainOutline.exp', 'hmin', resolution)
 
-    #move the closest node to the center
+    # move the closest node to the center
     pos = np.argmin(md.mesh.x**2 + md.mesh.y**2)
     md.mesh.x[pos] = 0.
     md.mesh.y[pos] = 0.
 
-    #delete domain
+    # delete domain
     os.remove('RoundDomainOutline.exp')