Index: /issm/trunk-jpl/src/m/classes/model.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/model.m	(revision 17990)
+++ /issm/trunk-jpl/src/m/classes/model.m	(revision 17991)
@@ -208,5 +208,5 @@
 			md.stressbalance.loadingforce=project2d(md,md.stressbalance.loadingforce,md.mesh.numberoflayers);
 			md.masstransport.spcthickness=project2d(md,md.masstransport.spcthickness,md.mesh.numberoflayers);
-			md.damage.spcdamage=project2d(md,md.damage.spcdamage,md.mesh.numberoflayers)
+			if ~isnan(md.damage.spcdamage), md.damage.spcdamage=project2d(md,md.damage.spcdamage,md.mesh.numberoflayers); end
 			md.thermal.spctemperature=project2d(md,md.thermal.spctemperature,md.mesh.numberoflayers);
 
@@ -236,4 +236,5 @@
 			md.geometry.base=project2d(md,md.geometry.base,1);
 			md.geometry.bed=project2d(md,md.geometry.bed,1);
+
 			md.mask.groundedice_levelset=project2d(md,md.mask.groundedice_levelset,1);
 			md.mask.ice_levelset=project2d(md,md.mask.ice_levelset,1);
@@ -250,6 +251,6 @@
 			mesh.numberofelements=md.mesh.numberofelements2d;
 			mesh.elements=md.mesh.elements2d;
-			mesh.vertexonboundary=project2d(md,md.mesh.vertexonboundary,1);
-			mesh.elementconnectivity=project2d(md,md.mesh.elementconnectivity,1);
+			if ~isnan(md.mesh.vertexonboundary), mesh.vertexonboundary=project2d(md,md.mesh.vertexonboundary,1); end
+			if ~isnan(md.mesh.elementconnectivity), mesh.elementconnectivity=project2d(md,md.mesh.elementconnectivity,1); end
 			md.mesh=mesh;
 
@@ -760,11 +761,13 @@
 
 			%connectivity
-			md.mesh.elementconnectivity=repmat(md.mesh.elementconnectivity,numlayers-1,1);
-			md.mesh.elementconnectivity(find(md.mesh.elementconnectivity==0))=NaN;
-			for i=2:numlayers-1,
-				md.mesh.elementconnectivity((i-1)*md.mesh.numberofelements2d+1:(i)*md.mesh.numberofelements2d,:)...
-					=md.mesh.elementconnectivity((i-1)*md.mesh.numberofelements2d+1:(i)*md.mesh.numberofelements2d,:)+md.mesh.numberofelements2d;
-			end
-			md.mesh.elementconnectivity(find(isnan(md.mesh.elementconnectivity)))=0;
+			if ~isnan(md.mesh.elementconnectivity)
+				md.mesh.elementconnectivity=repmat(md.mesh.elementconnectivity,numlayers-1,1);
+				md.mesh.elementconnectivity(find(md.mesh.elementconnectivity==0))=NaN;
+				for i=2:numlayers-1,
+					md.mesh.elementconnectivity((i-1)*md.mesh.numberofelements2d+1:(i)*md.mesh.numberofelements2d,:)...
+						=md.mesh.elementconnectivity((i-1)*md.mesh.numberofelements2d+1:(i)*md.mesh.numberofelements2d,:)+md.mesh.numberofelements2d;
+				end
+				md.mesh.elementconnectivity(find(isnan(md.mesh.elementconnectivity)))=0;
+			end
 
 			%materials
Index: /issm/trunk-jpl/src/m/classes/model.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/model.py	(revision 17990)
+++ /issm/trunk-jpl/src/m/classes/model.py	(revision 17991)
@@ -657,8 +657,9 @@
 		md.mesh.elementconnectivity=numpy.tile(md.mesh.elementconnectivity,(numlayers-1,1))
 		md.mesh.elementconnectivity[numpy.nonzero(md.mesh.elementconnectivity==0)]=-sys.maxint-1
-		for i in xrange(1,numlayers-1):
-			md.mesh.elementconnectivity[i*md.mesh.numberofelements2d:(i+1)*md.mesh.numberofelements2d,:] \
-				=md.mesh.elementconnectivity[i*md.mesh.numberofelements2d:(i+1)*md.mesh.numberofelements2d,:]+md.mesh.numberofelements2d
-		md.mesh.elementconnectivity[numpy.nonzero(md.mesh.elementconnectivity<0)]=0
+		if not numpy.isnan(md.mesh.lementconnectivity).all():
+			for i in xrange(1,numlayers-1):
+				md.mesh.elementconnectivity[i*md.mesh.numberofelements2d:(i+1)*md.mesh.numberofelements2d,:] \
+						=md.mesh.elementconnectivity[i*md.mesh.numberofelements2d:(i+1)*md.mesh.numberofelements2d,:]+md.mesh.numberofelements2d
+				md.mesh.elementconnectivity[numpy.nonzero(md.mesh.elementconnectivity<0)]=0
 
 		#materials
@@ -764,5 +765,5 @@
 		md.stressbalance.loadingforce=project2d(md,md.stressbalance.loadingforce,md.mesh.numberoflayers)
 		md.masstransport.spcthickness=project2d(md,md.masstransport.spcthickness,md.mesh.numberoflayers)
-		md.damage.spcdamage=project2d(md,md.damage.spcdamage,md.mesh.numberoflayers-1)
+		if not numpy.isnan(md.damage.spcdamage).all(): md.damage.spcdamage=project2d(md,md.damage.spcdamage,md.mesh.numberoflayers-1)
 		md.thermal.spctemperature=project2d(md,md.thermal.spctemperature,md.mesh.numberoflayers-1)
 
@@ -805,6 +806,6 @@
 		mesh.numberofelements=md.mesh.numberofelements2d
 		mesh.elements=md.mesh.elements2d
-		mesh.vertexonboundary=project2d(md,md.mesh.vertexonboundary,1)
-		mesh.elementconnectivity=project2d(md,md.mesh.elementconnectivity,1)
+		if not numpy.isnan(md.mesh.vertexonboundary).all(): mesh.vertexonboundary=project2d(md,md.mesh.vertexonboundary,1)
+		if not numpy.isnan(md.mesh.elementconnectivity).all(): mesh.elementconnectivity=project2d(md,md.mesh.elementconnectivity,1)
 		md.mesh=mesh