Index: /issm/trunk/src/c/DataSet/DataSet.cpp
===================================================================
--- /issm/trunk/src/c/DataSet/DataSet.cpp	(revision 3810)
+++ /issm/trunk/src/c/DataSet/DataSet.cpp	(revision 3811)
@@ -680,5 +680,5 @@
 
 			element=(Element*)(*object);
-			element->Configure(loads,nodes,materials,parameters);
+			element->Configure(elements,loads,nodes,materials,parameters);
 		}
 		if(EnumIsLoad((*object)->Enum())){
Index: /issm/trunk/src/c/ModelProcessorx/Balancedthickness/CreateElementsNodesAndMaterialsBalancedthickness.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/Balancedthickness/CreateElementsNodesAndMaterialsBalancedthickness.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/Balancedthickness/CreateElementsNodesAndMaterialsBalancedthickness.cpp	(revision 3811)
@@ -91,4 +91,6 @@
 		IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate");
 		IoModelFetchData(&iomodel->dhdt,NULL,NULL,iomodel_handle,"dhdt");
+		IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+		IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 	
 		for (i=0;i<iomodel->numberofelements;i++){
@@ -116,4 +118,7 @@
 		xfree((void**)&iomodel->accumulation_rate);
 		xfree((void**)&iomodel->dhdt);
+		xfree((void**)&iomodel->upperelements);
+		xfree((void**)&iomodel->lowerelements);
+
 
 	} //if (strcmp(meshtype,"2d")==0)
Index: /issm/trunk/src/c/ModelProcessorx/Balancedvelocities/CreateElementsNodesAndMaterialsBalancedvelocities.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/Balancedvelocities/CreateElementsNodesAndMaterialsBalancedvelocities.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/Balancedvelocities/CreateElementsNodesAndMaterialsBalancedvelocities.cpp	(revision 3811)
@@ -90,4 +90,6 @@
 		IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate");
 		IoModelFetchData(&iomodel->accumulation_rate,NULL,NULL,iomodel_handle,"accumulation_rate");
+		IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+		IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 
 		for (i=0;i<iomodel->numberofelements;i++){
@@ -116,4 +118,7 @@
 		xfree((void**)&iomodel->melting_rate);
 		xfree((void**)&iomodel->accumulation_rate);
+		xfree((void**)&iomodel->upperelements);
+		xfree((void**)&iomodel->lowerelements);
+
 
 	} //if (strcmp(meshtype,"2d")==0)
Index: /issm/trunk/src/c/ModelProcessorx/DiagnosticHoriz/CreateElementsNodesAndMaterialsDiagnosticHoriz.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/DiagnosticHoriz/CreateElementsNodesAndMaterialsDiagnosticHoriz.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/DiagnosticHoriz/CreateElementsNodesAndMaterialsDiagnosticHoriz.cpp	(revision 3811)
@@ -122,4 +122,7 @@
 			IoModelFetchData(&iomodel->weights,NULL,NULL,iomodel_handle,"weights");
 		}
+		IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+		IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
+
 	
 		for (i=0;i<iomodel->numberofelements;i++){
@@ -160,4 +163,7 @@
 			xfree((void**)&iomodel->weights);
 		}
+		xfree((void**)&iomodel->upperelements);
+		xfree((void**)&iomodel->lowerelements);
+
 
 
Index: /issm/trunk/src/c/ModelProcessorx/DiagnosticStokes/CreateElementsNodesAndMaterialsDiagnosticStokes.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/DiagnosticStokes/CreateElementsNodesAndMaterialsDiagnosticStokes.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/DiagnosticStokes/CreateElementsNodesAndMaterialsDiagnosticStokes.cpp	(revision 3811)
@@ -58,4 +58,6 @@
 	IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy");
 	IoModelFetchData(&iomodel->vz,NULL,NULL,iomodel_handle,"vz");
+	IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+	IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 	
 	if (iomodel->control_analysis){
@@ -105,4 +107,7 @@
 		xfree((void**)&iomodel->weights);
 	}
+	xfree((void**)&iomodel->upperelements);
+	xfree((void**)&iomodel->lowerelements);
+
 
 	/*Add new constrant material property to materials, at the end: */
Index: /issm/trunk/src/c/ModelProcessorx/DiagnosticVert/CreateElementsNodesAndMaterialsDiagnosticVert.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/DiagnosticVert/CreateElementsNodesAndMaterialsDiagnosticVert.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/DiagnosticVert/CreateElementsNodesAndMaterialsDiagnosticVert.cpp	(revision 3811)
@@ -51,4 +51,6 @@
 	IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx");
 	IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy");
+	IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+	IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 	
 	for (i=0;i<iomodel->numberofelements;i++){
@@ -80,4 +82,6 @@
 	xfree((void**)&iomodel->vx);
 	xfree((void**)&iomodel->vy);
+	xfree((void**)&iomodel->upperelements);
+	xfree((void**)&iomodel->lowerelements);
 
 	/*Add new constrant material property to materials, at the end: */
Index: /issm/trunk/src/c/ModelProcessorx/Melting/CreateElementsNodesAndMaterialsMelting.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/Melting/CreateElementsNodesAndMaterialsMelting.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/Melting/CreateElementsNodesAndMaterialsMelting.cpp	(revision 3811)
@@ -50,4 +50,6 @@
 	IoModelFetchData(&iomodel->melting_rate,NULL,NULL,iomodel_handle,"melting_rate");
 	IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater");
+	IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+	IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 	
 	for (i=0;i<iomodel->numberofelements;i++){
@@ -81,4 +83,6 @@
 	xfree((void**)&iomodel->melting_rate);
 	xfree((void**)&iomodel->elementonwater);
+	xfree((void**)&iomodel->upperelements);
+	xfree((void**)&iomodel->lowerelements);
 
 	/*Add new constrant material property tgo materials, at the end: */
Index: /issm/trunk/src/c/ModelProcessorx/Prognostic/CreateElementsNodesAndMaterialsPrognostic.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/Prognostic/CreateElementsNodesAndMaterialsPrognostic.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/Prognostic/CreateElementsNodesAndMaterialsPrognostic.cpp	(revision 3811)
@@ -91,4 +91,6 @@
 		IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx");
 		IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy");
+		IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+		IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 
 		for (i=0;i<iomodel->numberofelements;i++){
@@ -117,4 +119,7 @@
 		xfree((void**)&iomodel->vx);
 		xfree((void**)&iomodel->vy);
+		xfree((void**)&iomodel->upperelements);
+		xfree((void**)&iomodel->lowerelements);
+
 
 	} //if (strcmp(meshtype,"2d")==0)
Index: /issm/trunk/src/c/ModelProcessorx/SlopeCompute/CreateElementsNodesAndMaterialsSlopeCompute.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/SlopeCompute/CreateElementsNodesAndMaterialsSlopeCompute.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/SlopeCompute/CreateElementsNodesAndMaterialsSlopeCompute.cpp	(revision 3811)
@@ -71,4 +71,6 @@
 		IoModelFetchData(&iomodel->elementonbed,NULL,NULL,iomodel_handle,"elementonbed");
 		IoModelFetchData(&iomodel->elementonwater,NULL,NULL,iomodel_handle,"elementonwater");
+		IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+		IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 	
 		for (i=0;i<iomodel->numberofelements;i++){
@@ -89,4 +91,6 @@
 		xfree((void**)&iomodel->elementonbed);
 		xfree((void**)&iomodel->elementonwater);
+		xfree((void**)&iomodel->upperelements);
+		xfree((void**)&iomodel->lowerelements);
 
 	} //if (strcmp(meshtype,"2d")==0)
Index: /issm/trunk/src/c/ModelProcessorx/Thermal/CreateElementsNodesAndMaterialsThermal.cpp
===================================================================
--- /issm/trunk/src/c/ModelProcessorx/Thermal/CreateElementsNodesAndMaterialsThermal.cpp	(revision 3810)
+++ /issm/trunk/src/c/ModelProcessorx/Thermal/CreateElementsNodesAndMaterialsThermal.cpp	(revision 3811)
@@ -53,4 +53,6 @@
 	IoModelFetchData(&iomodel->vx,NULL,NULL,iomodel_handle,"vx");
 	IoModelFetchData(&iomodel->vy,NULL,NULL,iomodel_handle,"vy");
+	IoModelFetchData(&iomodel->upperelements,NULL,NULL,iomodel_handle,"upperelements");
+	IoModelFetchData(&iomodel->lowerelements,NULL,NULL,iomodel_handle,"lowerelements");
 	
 	for (i=0;i<iomodel->numberofelements;i++){
@@ -85,4 +87,6 @@
 	xfree((void**)&iomodel->vx);
 	xfree((void**)&iomodel->vy);
+	xfree((void**)&iomodel->upperelements);
+	xfree((void**)&iomodel->lowerelements);
 
 	/*Add new constrant material property tgo materials, at the end: */
Index: /issm/trunk/src/c/objects/Elements/Beam.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Beam.cpp	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Beam.cpp	(revision 3811)
@@ -164,5 +164,5 @@
 /*Object management*/
 /*FUNCTION Beam::Configure{{{1*/
-void  Beam::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
+void  Beam::Configure(DataSet* elementsin,DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
 
 	int i;
Index: /issm/trunk/src/c/objects/Elements/Beam.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Beam.h	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Beam.h	(revision 3811)
@@ -49,5 +49,5 @@
 		int   Id(); 
 		int   MyRank();
-		void  Configure(DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
+		void  Configure(DataSet* elements,DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
 		Object* copy();
 		void  SetClone(int* minranks);
Index: /issm/trunk/src/c/objects/Elements/Element.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Element.h	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Element.h	(revision 3811)
@@ -24,5 +24,5 @@
 		
 		virtual        ~Element(){};
-		virtual void   Configure(DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters)=0;
+		virtual void   Configure(DataSet* elements,DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters)=0;
 		
 		virtual void   CreateKMatrix(Mat Kgg,int analysis_type,int sub_analysis_type)=0;
Index: /issm/trunk/src/c/objects/Elements/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 3811)
@@ -26,4 +26,5 @@
 	this->matice=NULL;
 	this->matpar=NULL;
+	this->neighboors=NULL;
 	
 	this->inputs=NULL;
@@ -44,4 +45,5 @@
 	int penta_matice_id;
 	int penta_matpar_id;
+	int penta_element_ids[2];
 	double nodeinputs[6];
 	bool collapse;
@@ -56,8 +58,11 @@
 	penta_matice_id=index+1; //refers to the corresponding ice material object
 	penta_matpar_id=iomodel->numberofelements+1; //refers to the constant material parameters object
+	penta_elements_ids[0]=(int)*(iomodel->upperelements[index]);
+	penta_elements_ids[1]=(int)*(iomodel->lowerelements[index]);
 
 	this->InitHookNodes(penta_node_ids); this->nodes=NULL;
 	this->InitHookMatice(penta_matice_id);this->matice=NULL;
 	this->InitHookMatpar(penta_matpar_id);this->matpar=NULL;
+	this->InitHookNeighboors(penta_element_ids);this->neighboors=NULL;
 
 	//intialize inputs, and add as many inputs per element as requested: 
@@ -197,4 +202,5 @@
 	penta->hmatice.copy(&this->hmatice);
 	penta->hmatpar.copy(&this->hmatpar);
+	penta->hneighboors.copy(&this->hneighboors);
 
 	/*recover objects: */
@@ -202,4 +208,5 @@
 	penta->matice=(Matice*)penta->hmatice.delivers();
 	penta->matpar=(Matpar*)penta->hmatpar.delivers();
+	penta->neighboors=(Penta**)penta->hneighboors.deliverp();
 
 	return penta;
@@ -209,5 +216,5 @@
 /*Object management: */
 /*FUNCTION Penta::Configure {{{1*/
-void  Penta::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
+void  Penta::Configure(DataSet* elementsin, DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
 
 	/*Take care of hooking up all objects for this element, ie links the objects in the hooks to their respective 
@@ -216,4 +223,5 @@
 	this->hmatice.configure(materialsin);
 	this->hmatpar.configure(materialsin);
+	this->hneighboors.configure(elementsin);
 
 	/*Now, go pick up the objects inside the hooks: */
@@ -221,4 +229,5 @@
 	this->matice=(Matice*)this->hmatice.delivers();
 	this->matpar=(Matpar*)this->hmatpar.delivers();
+	this->neighboors=(Penta**)this->hneighboors.deliverp();
 
 	/*point parameters to real dataset: */
@@ -245,4 +254,5 @@
 	hmatice.Demarshall(&marshalled_dataset);
 	hmatpar.Demarshall(&marshalled_dataset);
+	hneighboors.Demarshall(&marshalled_dataset);
 
 	/*pointers are garbabe, until configuration is carried out: */
@@ -250,4 +260,5 @@
 	matice=NULL;
 	matpar=NULL;
+	neighboors=NULL;
 	
 	/*demarshall inputs: */
@@ -275,4 +286,5 @@
 	matice->DeepEcho();
 	matpar->DeepEcho();
+	printf("   neigboor ids: %i-%i\n",neigboors[0]->Id(),neigboors[1]->Id());
 	printf("   parameters\n");
 	parameters->DeepEcho();
@@ -320,4 +332,5 @@
 	hmatice.Marshall(&marshalled_dataset);
 	hmatpar.Marshall(&marshalled_dataset);
+	hneighboors.Marshall(&marshalled_dataset);
 
 	/*Marshall inputs: */
@@ -342,4 +355,5 @@
 		+hmatice.MarshallSize()
 		+hmatpar.MarshallSize()
+		+hneighboors.MarshallSize()
 		+inputs->MarshallSize()
 		+sizeof(int); //sizeof(int) for enum type
Index: /issm/trunk/src/c/objects/Elements/Penta.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.h	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Penta.h	(revision 3811)
@@ -27,12 +27,13 @@
 	public:
 
-		int id;
-		
-		Node** nodes; // 6 nodes
-		Matice* matice;  // 1 material ice
-		Matpar* matpar;  // 1 material parameter
+		int          id;
 
-		Parameters* parameters; //pointer to solution parameters
-		Inputs* inputs;
+		Node       **nodes;        // 6 nodes
+		Matice      *matice;       // 1 material ice
+		Matpar      *matpar;       // 1 material parameter
+		Penta       *neighboors;   // 2 neighboors: first one under, second one above
+
+		Parameters  *parameters;   //pointer to solution parameters
+		Inputs      *inputs;
 
 		/*FUNCTION constructors, destructors {{{1*/
@@ -42,5 +43,5 @@
 		/*}}}*/
 		/*FUNCTION object management {{{1*/
-		void  Configure(DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
+		void  Configure(DataSet* elements,DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
 		Object* copy();
 		void  DeepEcho();
Index: /issm/trunk/src/c/objects/Elements/PentaHook.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/PentaHook.cpp	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/PentaHook.cpp	(revision 3811)
@@ -45,2 +45,8 @@
 }
 
+/*FUNCTION PentaHook::InitHookNeighboors(int* element_ids){{{1*/
+void PentaHook::InitHookNeighboors(int* element_ids){
+	this->hneighboors.Init(element_ids,2);
+
+}
+/*}}}*/
Index: /issm/trunk/src/c/objects/Elements/PentaHook.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/PentaHook.h	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/PentaHook.h	(revision 3811)
@@ -14,4 +14,5 @@
 		Hook hmatice; // 1 ice material
 		Hook hmatpar; // 1 material parameter
+		Hook hneighboors // 2 elements, first down, second up
 
 
@@ -22,4 +23,5 @@
 		void InitHookMatice(int matice_id);
 		void InitHookMatpar(int matpar_id);
+		void InitHookNeighboors(int* element_ids);
 		/*}}}*/
 
Index: /issm/trunk/src/c/objects/Elements/Sing.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Sing.cpp	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Sing.cpp	(revision 3811)
@@ -106,5 +106,5 @@
 /*Object management*/
 /*FUNCTION Sing::Configure {{{1*/
-void  Sing::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
+void  Sing::Configure(DataSet* elementsin,DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
 
 	int i;
Index: /issm/trunk/src/c/objects/Elements/Sing.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Sing.h	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Sing.h	(revision 3811)
@@ -41,5 +41,5 @@
 		/*}}}*/
 		/*object management: {{{1*/
-		void  Configure(DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
+		void  Configure(DataSet* elements,DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
 		Object* copy();
 		void  DeepEcho();
Index: /issm/trunk/src/c/objects/Elements/Tria.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Tria.cpp	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Tria.cpp	(revision 3811)
@@ -211,5 +211,5 @@
 /*Object management: */
 /*FUNCTION Tria::Configure {{{1*/
-void  Tria::Configure(DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
+void  Tria::Configure(DataSet* elementsin, DataSet* loadsin, DataSet* nodesin, DataSet* materialsin, Parameters* parametersin){
 
 	/*Take care of hooking up all objects for this element, ie links the objects in the hooks to their respective 
Index: /issm/trunk/src/c/objects/Elements/Tria.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Tria.h	(revision 3810)
+++ /issm/trunk/src/c/objects/Elements/Tria.h	(revision 3811)
@@ -40,5 +40,5 @@
 		/*}}}*/
 		/*FUNCTION object management {{{1*/
-		void  Configure(DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
+		void  Configure(DataSet* elements,DataSet* loads,DataSet* nodes,DataSet* materials,Parameters* parameters);
 		Object* copy();
 		void  DeepEcho();
Index: /issm/trunk/src/c/objects/IoModel.cpp
===================================================================
--- /issm/trunk/src/c/objects/IoModel.cpp	(revision 3810)
+++ /issm/trunk/src/c/objects/IoModel.cpp	(revision 3811)
@@ -45,4 +45,6 @@
 		xfree((void**)&this->deadgrids);
 		xfree((void**)&this->uppernodes);
+		xfree((void**)&this->upperelements);
+		xfree((void**)&this->lowerelements);
 		xfree((void**)&this->gridonpattyn);
 	}
@@ -248,4 +250,6 @@
 	this->numlayers=0;
 	this->uppernodes=NULL;
+	this->upperelements=NULL;
+	this->lowerelements=NULL;
 	this->gridonhutter=NULL;
 	this->gridonmacayeal=NULL;
Index: /issm/trunk/src/c/objects/IoModel.h
===================================================================
--- /issm/trunk/src/c/objects/IoModel.h	(revision 3810)
+++ /issm/trunk/src/c/objects/IoModel.h	(revision 3811)
@@ -42,4 +42,5 @@
 		int     numlayers;
 		double* uppernodes;
+		double* upperelements;
 
 		/*elements type: */
Index: /issm/trunk/src/m/classes/@model/model.m
===================================================================
--- /issm/trunk/src/m/classes/@model/model.m	(revision 3810)
+++ /issm/trunk/src/m/classes/@model/model.m	(revision 3811)
@@ -88,4 +88,6 @@
 	%Projections
 	md.uppergrids=NaN;
+	md.upperelements=NaN;
+	md.lowerelements=NaN;
 	md.lowergrids=NaN;
 	md.deadgrids=NaN;
Index: /issm/trunk/src/m/classes/public/extrude.m
===================================================================
--- /issm/trunk/src/m/classes/public/extrude.m	(revision 3810)
+++ /issm/trunk/src/m/classes/public/extrude.m	(revision 3811)
@@ -103,4 +103,14 @@
 md.uppergrids=uppergrids;
 
+%same for lower and upper elements
+lowerelements=NaN*ones(number_el3d,1);
+upperelements=NaN*ones(number_el3d,1);
+lowerelements(md.numberofelements+1:end)=1:(numlayers-2)*md.numberofelements;
+upperelements(1:(numlayers-2)*md.numberofelements)=md.numberofelements+1:(numlayers-1)*md.numberofelements;
+md.lowerelements=lowerelements;
+md.upperelements=upperelements;
+
+
+
 %Save old mesh 
 md.x2d=md.x;
Index: /issm/trunk/src/m/classes/public/marshall.m
===================================================================
--- /issm/trunk/src/m/classes/public/marshall.m	(revision 3810)
+++ /issm/trunk/src/m/classes/public/marshall.m	(revision 3811)
@@ -44,4 +44,6 @@
 end
 WriteData(fid,md.uppergrids,'Mat','uppergrids');
+WriteData(fid,md.upperelements,'Mat','upperelements');
+WriteData(fid,md.lowerelements,'Mat','lowerelements');
 WriteData(fid,md.elementonbed,'Mat','elementonbed');
 WriteData(fid,md.elementonsurface,'Mat','elementonsurface');