Index: /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscMat.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscMat.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscMat.cpp	(revision 13595)
@@ -31,5 +31,5 @@
 PetscMat::PetscMat(int M,int N){
 
-	this->matrix=NewMat(M,N);
+	this->matrix=NewMat(M,N,IssmComm::GetComm());
 }
 /*}}}*/
@@ -128,5 +128,5 @@
 	_assert_(this->matrix);
 	_assert_(X->vector);
-	MatMultPatch(this->matrix,X->vector,AX->vector);
+	MatMultPatch(this->matrix,X->vector,AX->vector,IssmComm::GetComm());
 
 }
@@ -150,5 +150,5 @@
 	 IssmDouble* output=NULL;
 
-	 MatToSerial(&output,this->matrix);
+	 MatToSerial(&output,this->matrix,IssmComm::GetComm());
 	 return output;
 
Index: /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscVec.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscVec.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/objects/PetscVec.cpp	(revision 13595)
@@ -167,5 +167,5 @@
 	
 	IssmDouble* vec_serial=NULL;
-	VecToMPISerial(&vec_serial, this->vector);
+	VecToMPISerial(&vec_serial, this->vector,IssmComm::GetComm());
 	return vec_serial;
 
Index: /issm/trunk-jpl/src/c/toolkits/petsc/patches/GetOwnershipBoundariesFromRange.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/patches/GetOwnershipBoundariesFromRange.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/patches/GetOwnershipBoundariesFromRange.cpp	(revision 13595)
@@ -14,5 +14,5 @@
 #include "../../../shared/shared.h"
 
-void GetOwnershipBoundariesFromRange(int* plower_row,int* pupper_row,int range){
+void GetOwnershipBoundariesFromRange(int* plower_row,int* pupper_row,int range,COMM comm){
 
 	/*externals :*/
@@ -29,5 +29,5 @@
 	/*Gather all range values into allranges, for all nodes*/
 	allranges=xNew<int>(num_procs);
-	MPI_Allgather(&range,1,MPI_INT,allranges,1,MPI_INT,MPI_COMM_WORLD);
+	MPI_Allgather(&range,1,MPI_INT,allranges,1,MPI_INT,comm);
 
 	/*From all ranges, get lower row and upper row*/
Index: /issm/trunk-jpl/src/c/toolkits/petsc/patches/MatMultPatch.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/patches/MatMultPatch.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/patches/MatMultPatch.cpp	(revision 13595)
@@ -17,8 +17,8 @@
 
 /*Function prototypes: */
-int MatMultCompatible(Mat A,Vec x);
+int MatMultCompatible(Mat A,Vec x,COMM comm);
 void VecRelocalize(Vec* outvector,Vec vector,int m);
 
-void MatMultPatch(Mat A,Vec X, Vec AX){ //same prototype as MatMult in Petsc
+void MatMultPatch(Mat A,Vec X, Vec AX,COMM comm){ //same prototype as MatMult in Petsc
 
 	int m,n;
@@ -27,5 +27,5 @@
 	_assert_(A); _assert_(X);
 
-	if (MatMultCompatible(A,X)){
+	if (MatMultCompatible(A,X,comm)){
 		MatMult(A,X,AX); 
 	}
@@ -42,5 +42,5 @@
 }
 
-int MatMultCompatible(Mat A,Vec x){
+int MatMultCompatible(Mat A,Vec x,COMM comm){
 	
 	/*error management*/
@@ -58,6 +58,6 @@
 	
 	/*synchronize result: */
-	MPI_Reduce (&result,&sumresult,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD );
-	MPI_Bcast(&sumresult,1,MPI_INT,0,MPI_COMM_WORLD);                
+	MPI_Reduce (&result,&sumresult,1,MPI_INT,MPI_SUM,0,comm );
+	MPI_Bcast(&sumresult,1,MPI_INT,0,comm);                
 	if (sumresult!=num_procs){
 		result=0;
Index: /issm/trunk-jpl/src/c/toolkits/petsc/patches/MatToSerial.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/patches/MatToSerial.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/patches/MatToSerial.cpp	(revision 13595)
@@ -13,5 +13,5 @@
 #include "../../../shared/shared.h"
 
-void MatToSerial(double** poutmatrix,Mat matrix){
+void MatToSerial(double** poutmatrix,Mat matrix,COMM comm){
 
 	int i;
@@ -66,10 +66,10 @@
 			buffer[1]=lower_row;
 			buffer[2]=range;
-			MPI_Send(buffer,3,MPI_INT,0,1,MPI_COMM_WORLD);   
-			if (range)MPI_Send(local_matrix,N*range,MPI_DOUBLE,0,1,MPI_COMM_WORLD); 
+			MPI_Send(buffer,3,MPI_INT,0,1,comm);   
+			if (range)MPI_Send(local_matrix,N*range,MPI_DOUBLE,0,1,comm); 
 		}
 		if (my_rank==0){
-			MPI_Recv(buffer,3,MPI_INT,i,1,MPI_COMM_WORLD,&status); 
-			if (buffer[2])MPI_Recv(outmatrix+(buffer[1]*N),N*buffer[2],MPI_DOUBLE,i,1,MPI_COMM_WORLD,&status);
+			MPI_Recv(buffer,3,MPI_INT,i,1,comm,&status); 
+			if (buffer[2])MPI_Recv(outmatrix+(buffer[1]*N),N*buffer[2],MPI_DOUBLE,i,1,comm,&status);
 		}
 	} 
Index: /issm/trunk-jpl/src/c/toolkits/petsc/patches/NewMat.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/patches/NewMat.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/patches/NewMat.cpp	(revision 13595)
@@ -21,5 +21,5 @@
 
 /*NewMat(int M,int N){{{*/
-Mat NewMat(int M,int N){
+Mat NewMat(int M,int N,COMM comm){
 
 	/*output:*/
@@ -40,7 +40,7 @@
 
 	#if _PETSC_MAJOR_ == 3 && _PETSC_MINOR_ > 2
-	MatCreateAIJ(MPI_COMM_WORLD,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
+	MatCreateAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
 	#else
-	MatCreateMPIAIJ(MPI_COMM_WORLD,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
+	MatCreateMPIAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
 	#endif
 
@@ -49,5 +49,5 @@
 /*}}}*/
 /*NewMat(int M,int N,double sparsity){{{*/
-Mat NewMat(int M,int N,double sparsity){
+Mat NewMat(int M,int N,double sparsity,COMM comm){
 
 	/*output:*/
@@ -69,11 +69,11 @@
 	#if _PETSC_MAJOR_ == 3 && _PETSC_MINOR_ > 2
 	if(sparsity==1){
-		MatCreateDense(MPI_COMM_WORLD,m,n,M,N,NULL,&outmatrix); 
+		MatCreateDense(comm,m,n,M,N,NULL,&outmatrix); 
 	}
 	else{
-		MatCreateAIJ(MPI_COMM_WORLD,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
+		MatCreateAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
 	}
 	#else
-	MatCreateMPIAIJ(MPI_COMM_WORLD,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
+	MatCreateMPIAIJ(comm,m,n,M,N,d_nz,NULL,o_nz,NULL,&outmatrix); 
 	#endif
 
@@ -82,5 +82,5 @@
 /*}}}*/
 /*NewMat(int M,int N,int connectivity,int numberofdofspernode){{{*/
-Mat NewMat(int M,int N,int connectivity,int numberofdofspernode){
+Mat NewMat(int M,int N,int connectivity,int numberofdofspernode,COMM comm){
 
 	/*output:*/
@@ -106,5 +106,5 @@
 	o_nz=(int)connectivity*numberofdofspernode/2;
 
-	MatCreate(MPI_COMM_WORLD,&outmatrix);
+	MatCreate(comm,&outmatrix);
 	MatSetSizes(outmatrix,m,n,M,N);
 	MatSetFromOptions(outmatrix);
Index: /issm/trunk-jpl/src/c/toolkits/petsc/patches/VecToMPISerial.cpp
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/patches/VecToMPISerial.cpp	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/patches/VecToMPISerial.cpp	(revision 13595)
@@ -13,5 +13,5 @@
 #include "../../../shared/shared.h"
 
-int VecToMPISerial(double** pgathered_vector, Vec vector){
+int VecToMPISerial(double** pgathered_vector, Vec vector,COMM comm){
 	
 	int i;
@@ -66,10 +66,10 @@
 			buffer[1]=lower_row;
 			buffer[2]=range;
-			MPI_Send(buffer,3,MPI_INT,0,1,MPI_COMM_WORLD);  
-			if (range)MPI_Send(local_vector,range,MPI_DOUBLE,0,1,MPI_COMM_WORLD); 
+			MPI_Send(buffer,3,MPI_INT,0,1,comm);  
+			if (range)MPI_Send(local_vector,range,MPI_DOUBLE,0,1,comm); 
 		}
 		if (my_rank==0){
-			MPI_Recv(buffer,3,MPI_INT,i,1,MPI_COMM_WORLD,&status); 
-			if (buffer[2])MPI_Recv(gathered_vector+buffer[1],buffer[2],MPI_DOUBLE,i,1,MPI_COMM_WORLD,&status);
+			MPI_Recv(buffer,3,MPI_INT,i,1,comm,&status); 
+			if (buffer[2])MPI_Recv(gathered_vector+buffer[1],buffer[2],MPI_DOUBLE,i,1,comm,&status);
 		}
 	}
@@ -81,5 +81,5 @@
 
 	/*Now, broadcast gathered_vector from node 0 to other nodes: */
-	MPI_Bcast(gathered_vector,vector_size,MPI_DOUBLE,0,MPI_COMM_WORLD);
+	MPI_Bcast(gathered_vector,vector_size,MPI_DOUBLE,0,comm);
 
 	/*Assign output pointers: */
Index: /issm/trunk-jpl/src/c/toolkits/petsc/patches/petscpatches.h
===================================================================
--- /issm/trunk-jpl/src/c/toolkits/petsc/patches/petscpatches.h	(revision 13594)
+++ /issm/trunk-jpl/src/c/toolkits/petsc/patches/petscpatches.h	(revision 13595)
@@ -13,15 +13,16 @@
 #include "./SolverEnum.h"
 #include "../../toolkitsenums.h"
+#include "../../../include/types.h"
 
 class Parameters;
 
 Vec NewVec(int size,bool fromlocalsize=false);
-Mat NewMat(int M,int N);
-Mat NewMat(int M,int N,double sparsity);
-Mat NewMat(int M,int N,int connectivity,int numberofdofspernode);
+Mat NewMat(int M,int N,COMM comm);
+Mat NewMat(int M,int N,double sparsity,COMM comm);
+Mat NewMat(int M,int N,int connectivity,int numberofdofspernode, COMM comm);
 
 int VecTranspose(Vec* tvector,Vec vector);
-int VecToMPISerial(double** pgathered_vector, Vec vector);
-void GetOwnershipBoundariesFromRange(int* plower_row,int* pupper_row,int range);
+int VecToMPISerial(double** pgathered_vector, Vec vector,COMM comm);
+void GetOwnershipBoundariesFromRange(int* plower_row,int* pupper_row,int range,COMM comm);
 void MatFree(Mat* pmat);
 void ISFree(IS* pis);
@@ -35,6 +36,6 @@
 void PetscOptionsDetermineSolverType(int* psolver_type);
 void VecMerge(Vec A, Vec B, double* row_partition_vector,int row_partition_size);
-void MatMultPatch(Mat A,Vec X, Vec AX);
-void MatToSerial(double** poutmatrix,Mat matrix);
+void MatMultPatch(Mat A,Vec X, Vec AX,COMM comm);
+void MatToSerial(double** poutmatrix,Mat matrix,COMM comm);
 void VecDuplicatePatch(Vec* output, Vec input);
 Vec  SerialToVec(double* vector,int vector_size);