Index: /issm/trunk-jpl/src/c/modules/IoModelToConstraintsx/IoModelToConstraintsx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/IoModelToConstraintsx/IoModelToConstraintsx.cpp	(revision 15629)
+++ /issm/trunk-jpl/src/c/modules/IoModelToConstraintsx/IoModelToConstraintsx.cpp	(revision 15630)
@@ -57,6 +57,6 @@
 			for(i=0;i<iomodel->numberofedges;i++){
 				if(my_edges[i]){
-					v1 = iomodel->edges[2*i+0]-1;
-					v2 = iomodel->edges[2*i+1]-1;
+					v1 = iomodel->edges[3*i+0]-1;
+					v2 = iomodel->edges[3*i+1]-1;
 					if(!xIsNan<IssmDouble>(spcdata[v1]) && !xIsNan<IssmDouble>(spcdata[v2])){
 						constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
@@ -99,6 +99,6 @@
 			for(i=0;i<iomodel->numberofedges;i++){
 				if(my_edges[i]){
-					v1 = iomodel->edges[2*i+0]-1;
-					v2 = iomodel->edges[2*i+1]-1;
+					v1 = iomodel->edges[3*i+0]-1;
+					v2 = iomodel->edges[3*i+1]-1;
 					values=xNew<IssmDouble>(N);
 					spcpresent=false;
Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateEdges.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateEdges.cpp	(revision 15629)
+++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateEdges.cpp	(revision 15630)
@@ -19,5 +19,6 @@
 	int  i,j,v1,v2,v3;
 	int  maxnbe,nbe,elementnbe,elementnbv;
-	int *elementedges = NULL;
+	int *elementedges         = NULL;
+	int *elementedges_markers = NULL;
 
 	/*Mesh dependent variables*/
@@ -25,22 +26,24 @@
 		elementnbv = 3;
 		elementnbe = 3;
-		elementedges=xNew<int>(elementnbe*2);
-		elementedges[2*0+0] = 1; elementedges[2*0+1] = 2;
-		elementedges[2*1+0] = 2; elementedges[2*1+1] = 0;
-		elementedges[2*2+0] = 0; elementedges[2*2+1] = 1;
+		elementedges         = xNew<int>(elementnbe*2);
+		elementedges_markers = xNew<int>(elementnbe);
+		elementedges[2*0+0] = 1; elementedges[2*0+1] = 2; elementedges_markers[0] = 1;
+		elementedges[2*1+0] = 2; elementedges[2*1+1] = 0; elementedges_markers[1] = 1;
+		elementedges[2*2+0] = 0; elementedges[2*2+1] = 1; elementedges_markers[2] = 1;
 	}
 	else if(iomodel->dim==3){
 		elementnbv = 6;
 		elementnbe = 9;
-		elementedges=xNew<int>(elementnbe*2);
-		elementedges[2*0+0] = 0; elementedges[2*0+1] = 3;
-		elementedges[2*1+0] = 1; elementedges[2*1+1] = 4;
-		elementedges[2*2+0] = 2; elementedges[2*2+1] = 5;
-		elementedges[2*3+0] = 1; elementedges[2*3+1] = 2;
-		elementedges[2*4+0] = 2; elementedges[2*4+1] = 0;
-		elementedges[2*5+0] = 0; elementedges[2*5+1] = 1;
-		elementedges[2*6+0] = 4; elementedges[2*6+1] = 5;
-		elementedges[2*7+0] = 5; elementedges[2*7+1] = 3;
-		elementedges[2*8+0] = 3; elementedges[2*8+1] = 4;
+		elementedges         = xNew<int>(elementnbe*2);
+		elementedges_markers = xNew<int>(elementnbe);
+		elementedges[2*0+0] = 0; elementedges[2*0+1] = 3; elementedges_markers[0] = 2;
+		elementedges[2*1+0] = 1; elementedges[2*1+1] = 4; elementedges_markers[1] = 2;
+		elementedges[2*2+0] = 2; elementedges[2*2+1] = 5; elementedges_markers[2] = 2;
+		elementedges[2*3+0] = 1; elementedges[2*3+1] = 2; elementedges_markers[3] = 1;
+		elementedges[2*4+0] = 2; elementedges[2*4+1] = 0; elementedges_markers[4] = 1;
+		elementedges[2*5+0] = 0; elementedges[2*5+1] = 1; elementedges_markers[5] = 1;
+		elementedges[2*6+0] = 4; elementedges[2*6+1] = 5; elementedges_markers[6] = 3;
+		elementedges[2*7+0] = 5; elementedges[2*7+1] = 3; elementedges_markers[7] = 3;
+		elementedges[2*8+0] = 3; elementedges[2*8+1] = 4; elementedges_markers[8] = 3;
 	}
 	else{
@@ -52,5 +55,5 @@
 
 	/*Initialize intermediaries*/
-	int *edgestemp                 = xNew<int>(maxnbe*2);                             /*format: [vertex1 vertex2]       */
+	int *edgestemp                 = xNew<int>(maxnbe*3);                             /*format: [vertex1 vertex2 marker]       */
 	int *element_edge_connectivity = xNew<int>(iomodel->numberofelements*elementnbe); /*format: [edge1 edge2 ... edgen] */
 
@@ -80,5 +83,5 @@
 			/*Go through all processed edges connected to v1 and check whether we have seen this edge yet*/
 			for(int e=head_minv[v1]; e!=-1; e=next_edge[e]){
-				if(edgestemp[e*2+1]==v2+1){
+				if(edgestemp[e*3+1]==v2+1){
 					exist = true;
 					element_edge_connectivity[i*elementnbe+j]=e;
@@ -92,6 +95,7 @@
 
 				/*Update edges*/
-				edgestemp[nbe*2+0] = v1+1;
-				edgestemp[nbe*2+1] = v2+1;
+				edgestemp[nbe*3+0] = v1+1;
+				edgestemp[nbe*3+1] = v2+1;
+				edgestemp[nbe*3+2] = elementedges_markers[j];
 
 				/*Update Connectivity*/
@@ -113,6 +117,6 @@
 
 	/*Create final edges*/
-	int* edges = xNew<int>(nbe*2); /*vertex1 vertex2*/
-	for(int i=0;i<2*nbe;i++) edges[i] = edgestemp[i];
+	int* edges = xNew<int>(nbe*2); /*format: [vertex1 vertex2 marker]*/
+	for(int i=0;i<3*nbe;i++) edges[i] = edgestemp[i];
 
 	/*Clean up*/
