Index: /issm/trunk-jpl/src/c/classes/Elements/Element.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 25485)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.cpp	(revision 25486)
@@ -737,4 +737,5 @@
 	Gauss* gauss=this->NewGauss(5);
 	while(gauss->next()){
+
 		this->JacobianDeterminant(&Jdet,xyz_list,gauss);
 
@@ -3394,5 +3395,5 @@
 }
 /*}}}*/
-void       Element::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int* flagsindices_counter,int set1_enum,int set2_enum){/*{{{*/
+void       Element::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int* flagsindices_counter,int set1_enum,int set2_enum, int analysis_type){/*{{{*/
 
 	/*Intermediaries*/
@@ -3405,8 +3406,6 @@
 	/*Loop over all nodes*/
 	for(int i=0;i<numnodes;i++){
-
 		int nodelid = this->nodes[i]->Lid();
 		if(!flags[nodelid]){
-
 			/*flag current node so that no other element processes it*/
 			flags[nodelid]=true;
@@ -3448,6 +3447,5 @@
 	/*Special case: 2d/3d coupling, the node of this element might be connected
 	 *to the basal element*/
-	int analysis_type,approximation,numlayers;
-	parameters->FindParam(&analysis_type,AnalysisTypeEnum);
+	int approximation,numlayers;
 	if(analysis_type==StressbalanceAnalysisEnum){
 		this->GetInputValue(&approximation,ApproximationEnum);
Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 25485)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 25486)
@@ -164,5 +164,5 @@
 		void               ResultToMatrix(IssmDouble* values,int ncols,int output_enum);
 		void               ResultToVector(Vector<IssmDouble>* vector,int output_enum);
-		void               SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int* flagsindices_counter,int set1_enum,int set2_enum);
+		void               SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int* flagsindices_counter,int set1_enum,int set2_enum, int analysis_type);
 		void               SetBoolInput(Inputs* inputs,int enum_in,bool value);
 
Index: /issm/trunk-jpl/src/c/modules/AllocateSystemMatricesx/AllocateSystemMatricesx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/AllocateSystemMatricesx/AllocateSystemMatricesx.cpp	(revision 25485)
+++ /issm/trunk-jpl/src/c/modules/AllocateSystemMatricesx/AllocateSystemMatricesx.cpp	(revision 25486)
@@ -36,5 +36,5 @@
 	numberofdofspernode=femmodel->nodes->MaxNumDofs(GsetEnum);
 
-	/*if our matrices are coming from issm, we don't do dynamic allocation like Petsc 
+	/*if our matrices are coming from issm, we don't do dynamic allocation like Petsc
 	 * does, and this routine is essentially useless. Force standard alloc in this case: */
 	toolkittype=ToolkitOptions::GetToolkitType();
@@ -118,4 +118,6 @@
 	int numnodesperload     = femmodel->loads->MaxNumNodes();
 
+	int elementssize        = femmodel->elements->Size();
+	int loadssize           = femmodel->loads->Size();
 	/*First, we are building chaining vectors so that we know what nodes are
 	 * connected to what elements. These vectors are such that:
@@ -124,9 +126,9 @@
 	 * id*/
 	head_e         = xNew<int>(localnumnodes); for(i=0;i<localnumnodes;i++) head_e[i]=-1;
-	next_e         = xNew<int>(femmodel->elements->Size()*numnodesperelement);
-	count2offset_e = xNew<int>(femmodel->elements->Size()*numnodesperelement);
+	next_e         = xNew<int>(elementssize*numnodesperelement);
+	count2offset_e = xNew<int>(elementssize*numnodesperelement);
 
 	k=0;
-	for(i=0;i<femmodel->elements->Size();i++){
+	for(i=0;i<elementssize;i++){
 		element = xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(i));
 		int elementnumnodes = element->GetNumberOfNodes();
@@ -149,8 +151,8 @@
 	/*Chain for loads*/
 	head_l         = xNew<int>(localnumnodes); for(i=0;i<localnumnodes;i++) head_l[i]=-1;
-	next_l         = xNew<int>(femmodel->loads->Size()*numnodesperload);
-	count2offset_l = xNew<int>(femmodel->loads->Size()*numnodesperload);
+	next_l         = xNew<int>(loadssize*numnodesperload);
+	count2offset_l = xNew<int>(loadssize*numnodesperload);
 	k=0;
-	for(i=0;i<femmodel->loads->Size();i++){
+	for(i=0;i<loadssize;i++){
 		load = xDynamicCast<Load*>(femmodel->loads->GetObjectByOffset(i));
 		int loadnumnodes = load->GetNumberOfNodes();
@@ -177,7 +179,9 @@
 	int  *o_connectivity       = xNewZeroInit<int>(localnumnodes);
 	int   flagsindices_counter;
+	int   analysis_type;
 
 	Vector<IssmDouble> *connectivity_clone= new Vector<IssmDouble>(localmasters,numnodes);
 
+	femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum);
 	/*Resetting flags to false at each iteration takes a lot of time, so we keep track of the flags
 	 * to reset in flagsindices, initialized with -1*/
@@ -190,5 +194,4 @@
 		int   lid = node->Lid();
 		int   pid = node->Pid();
-
 		/*Reinitialize flags to false*/
 		j=0;
@@ -208,5 +211,5 @@
 			offset=count2offset_e[j];
 			element=xDynamicCast<Element*>(femmodel->elements->GetObjectByOffset(offset));
-			element->SetwiseNodeConnectivity(&d_nz,&o_nz,node,flags,flagsindices,&flagsindices_counter,set1enum,set2enum);
+			element->SetwiseNodeConnectivity(&d_nz,&o_nz,node,flags,flagsindices,&flagsindices_counter,set1enum,set2enum,analysis_type);
 			if(node->IsClone()){
 				connectivity_clone->SetValue(pid,d_nz+o_nz,ADD_VAL);