Added ScaLBL_Gradient_Unpack function

common/ScaLBL.cpp

@@ -3582,6 +3582,91 @@ void ScaLBL_Communicator::RecvD3Q19AA(double *dist){
 
 }
 
+void ScaLBL_Communicator::RecvGrad(double *grad){
+
+	// Recieves halo and incorporates into D3Q19 based stencil gradient computation
+	//...................................................................................
+	// Wait for completion of D3Q19 communication
+	MPI_Waitall(18,req1,stat1);
+	MPI_Waitall(18,req2,stat2);
+	ScaLBL_DeviceBarrier();
+
+	//...................................................................................
+	// NOTE: AA Routine writes to opposite 
+	// Unpack the gradributions on the device
+	//...................................................................................
+	//...Unpacking for x face(2,8,10,12,14)................................
+	ScaLBL_Gradient_Unpack(2,dvcRecvDist_x,0,recvCount_x,recvbuf_x,grad,N);
+	ScaLBL_Gradient_Unpack(8,dvcRecvDist_x,recvCount_x,recvCount_x,recvbuf_x,grad,N);
+	ScaLBL_Gradient_Unpack(10,dvcRecvDist_x,2*recvCount_x,recvCount_x,recvbuf_x,grad,N);
+	ScaLBL_Gradient_Unpack(12,dvcRecvDist_x,3*recvCount_x,recvCount_x,recvbuf_x,grad,N);
+	ScaLBL_Gradient_Unpack(14,dvcRecvDist_x,4*recvCount_x,recvCount_x,recvbuf_x,grad,N);
+	//...................................................................................
+	//...Packing for X face(1,7,9,11,13)................................
+	ScaLBL_Gradient_Unpack(1,dvcRecvDist_X,0,recvCount_X,recvbuf_X,grad,N);
+	ScaLBL_Gradient_Unpack(7,dvcRecvDist_X,recvCount_X,recvCount_X,recvbuf_X,grad,N);
+	ScaLBL_Gradient_Unpack(9,dvcRecvDist_X,2*recvCount_X,recvCount_X,recvbuf_X,grad,N);
+	ScaLBL_Gradient_Unpack(11,dvcRecvDist_X,3*recvCount_X,recvCount_X,recvbuf_X,grad,N);
+	ScaLBL_Gradient_Unpack(13,dvcRecvDist_X,4*recvCount_X,recvCount_X,recvbuf_X,grad,N);
+	//...................................................................................
+	//...Packing for y face(4,8,9,16,18).................................
+	ScaLBL_Gradient_Unpack(4,dvcRecvDist_y,0,recvCount_y,recvbuf_y,grad,N);
+	ScaLBL_Gradient_Unpack(8,dvcRecvDist_y,recvCount_y,recvCount_y,recvbuf_y,grad,N);
+	ScaLBL_Gradient_Unpack(9,dvcRecvDist_y,2*recvCount_y,recvCount_y,recvbuf_y,grad,N);
+	ScaLBL_Gradient_Unpack(16,dvcRecvDist_y,3*recvCount_y,recvCount_y,recvbuf_y,grad,N);
+	ScaLBL_Gradient_Unpack(18,dvcRecvDist_y,4*recvCount_y,recvCount_y,recvbuf_y,grad,N);
+	//...................................................................................
+	//...Packing for Y face(3,7,10,15,17).................................
+	ScaLBL_Gradient_Unpack(3,dvcRecvDist_Y,0,recvCount_Y,recvbuf_Y,grad,N);
+	ScaLBL_Gradient_Unpack(7,dvcRecvDist_Y,recvCount_Y,recvCount_Y,recvbuf_Y,grad,N);
+	ScaLBL_Gradient_Unpack(10,dvcRecvDist_Y,2*recvCount_Y,recvCount_Y,recvbuf_Y,grad,N);
+	ScaLBL_Gradient_Unpack(15,dvcRecvDist_Y,3*recvCount_Y,recvCount_Y,recvbuf_Y,grad,N);
+	ScaLBL_Gradient_Unpack(17,dvcRecvDist_Y,4*recvCount_Y,recvCount_Y,recvbuf_Y,grad,N);
+	//...................................................................................
+	//...Packing for z face(6,12,13,16,17)................................
+	ScaLBL_Gradient_Unpack(6,dvcRecvDist_z,0,recvCount_z,recvbuf_z,grad,N);
+	ScaLBL_Gradient_Unpack(12,dvcRecvDist_z,recvCount_z,recvCount_z,recvbuf_z,grad,N);
+	ScaLBL_Gradient_Unpack(13,dvcRecvDist_z,2*recvCount_z,recvCount_z,recvbuf_z,grad,N);
+	ScaLBL_Gradient_Unpack(16,dvcRecvDist_z,3*recvCount_z,recvCount_z,recvbuf_z,grad,N);
+	ScaLBL_Gradient_Unpack(17,dvcRecvDist_z,4*recvCount_z,recvCount_z,recvbuf_z,grad,N);
+	//...Packing for Z face(5,11,14,15,18)................................
+	ScaLBL_Gradient_Unpack(5,dvcRecvDist_Z,0,recvCount_Z,recvbuf_Z,grad,N);
+	ScaLBL_Gradient_Unpack(11,dvcRecvDist_Z,recvCount_Z,recvCount_Z,recvbuf_Z,grad,N);
+	ScaLBL_Gradient_Unpack(14,dvcRecvDist_Z,2*recvCount_Z,recvCount_Z,recvbuf_Z,grad,N);
+	ScaLBL_Gradient_Unpack(15,dvcRecvDist_Z,3*recvCount_Z,recvCount_Z,recvbuf_Z,grad,N);
+	ScaLBL_Gradient_Unpack(18,dvcRecvDist_Z,4*recvCount_Z,recvCount_Z,recvbuf_Z,grad,N);
+	//..................................................................................
+	//...Pack the xy edge (8)................................
+	ScaLBL_Gradient_Unpack(8,dvcRecvDist_xy,0,recvCount_xy,recvbuf_xy,grad,N);
+	//...Pack the Xy edge (9)................................
+	ScaLBL_Gradient_Unpack(9,dvcRecvDist_Xy,0,recvCount_Xy,recvbuf_Xy,grad,N);
+	//...Pack the xY edge (10)................................
+	ScaLBL_Gradient_Unpack(10,dvcRecvDist_xY,0,recvCount_xY,recvbuf_xY,grad,N);
+	//...Pack the XY edge (7)................................
+	ScaLBL_Gradient_Unpack(7,dvcRecvDist_XY,0,recvCount_XY,recvbuf_XY,grad,N);
+	//...Pack the xz edge (12)................................
+	ScaLBL_Gradient_Unpack(12,dvcRecvDist_xz,0,recvCount_xz,recvbuf_xz,grad,N);
+	//...Pack the xZ edge (14)................................
+	ScaLBL_Gradient_Unpack(14,dvcRecvDist_xZ,0,recvCount_xZ,recvbuf_xZ,grad,N);
+	//...Pack the Xz edge (13)................................
+	ScaLBL_Gradient_Unpack(13,dvcRecvDist_Xz,0,recvCount_Xz,recvbuf_Xz,grad,N);
+	//...Pack the XZ edge (11)................................
+	ScaLBL_Gradient_Unpack(11,dvcRecvDist_XZ,0,recvCount_XZ,recvbuf_XZ,grad,N);
+	//...Pack the yz edge (16)................................
+	ScaLBL_Gradient_Unpack(16,dvcRecvDist_yz,0,recvCount_yz,recvbuf_yz,grad,N);
+	//...Pack the yZ edge (18)................................
+	ScaLBL_Gradient_Unpack(18,dvcRecvDist_yZ,0,recvCount_yZ,recvbuf_yZ,grad,N);
+	//...Pack the Yz edge (17)................................
+	ScaLBL_Gradient_Unpack(17,dvcRecvDist_Yz,0,recvCount_Yz,recvbuf_Yz,grad,N);
+	//...Pack the YZ edge (15)................................
+	ScaLBL_Gradient_Unpack(15,dvcRecvDist_YZ,0,recvCount_YZ,recvbuf_YZ,grad,N);
+	//...................................................................................
+	Lock=false; // unlock the communicator after communications complete
+	//...................................................................................
+
+}
+
+
 
 void ScaLBL_Communicator::TestSendD3Q19(double *f_even, double *f_odd){
 
@@ -4025,8 +4110,6 @@ void ScaLBL_Communicator::SendHalo(double *data){
 	// Send / Recv all the phase indcator field values
 	//...................................................................................
 
-
-
 	MPI_Isend(sendbuf_x, sendCount_x,MPI_DOUBLE,rank_x,sendtag,MPI_COMM_SCALBL,&req1[0]);
 	MPI_Irecv(recvbuf_X, recvCount_X,MPI_DOUBLE,rank_X,recvtag,MPI_COMM_SCALBL,&req2[0]);
 	MPI_Isend(sendbuf_X, sendCount_X,MPI_DOUBLE,rank_X,sendtag,MPI_COMM_SCALBL,&req1[1]);

common/ScaLBL.h

@@ -37,6 +37,9 @@ extern "C" void ScaLBL_Scalar_Pack(int *list, int count, double *sendbuf, double
 
 extern "C" void ScaLBL_Scalar_Unpack(int *list, int count, double *recvbuf, double *Data, int N);
 
+extern "C" void ScaLBL_Gradient_Unpack(double weight, double Cqx, double Cqy, double Cqz, 
+		int *list, int start, int count, double *recvbuf, double *grad, int N);
+
 extern "C" void ScaLBL_PackDenD3Q7(int *list, int count, double *sendbuf, int number, double *Data, int N);
 
 extern "C" void ScaLBL_UnpackDenD3Q7(int *list, int count, double *recvbuf, int number, double *Data, int N);
@@ -294,6 +297,7 @@ public:
 	void TriRecvD3Q7AA(double *Aq, double *Bq, double *Cq);
 	void SendHalo(double *data);
 	void RecvHalo(double *data);
+	void RecvGrad(double *Gradient);
 	void RegularLayout(IntArray map, double *data, DoubleArray &regdata);
 
 	// Routines to set boundary conditions

cpu/dfh.cpp

@@ -1,6 +1,28 @@
 #include <math.h>
 #include <stdio.h>
 
+extern "C" void ScaLBL_Gradient_Unpack(double weight, double Cqx, double Cqy, double Cqz, 
+		int *list, int start, int count, double *recvbuf, double *grad, int N){
+	//....................................................................................
+	// unpack halo and incorporate into D3Q19 based gradient 
+	// Distribution q matche Cqx, Cqy, Cqz
+	//....................................................................................
+	int n,idx;
+	double value;
+	for (idx=0; idx<count; idx++){
+		// Get the index from the list
+		n = list[start+idx];
+		// unpack the distribution to the proper location
+		value=weight*recvbuf[idx];
+		if (!(n<0)){
+			// PARALLEL UPDATE MUST BE DONE ATOMICALLY
+			grad[n] += Cqx*value;
+			grad[N+n] += Cqy*value;
+			grad[2*N+n] += Cqz*value;
+		}
+	}
+}
+
 extern "C" void ScaLBL_DFH_Init(double *Phi, double *Den, double *Aq, double *Bq, int start, int finish, int Np){
 	int idx,n;
 	double phi,nA,nB;

tests/TestColorGradDFH.cpp

@@ -192,18 +192,17 @@ int main(int argc, char **argv)
 				}
 			}
 		}
-		Np+=32;
+		int Npad=(Np/16 + 2)*16;
 
-		int neighborSize=18*Np*sizeof(int);
 		int *neighborList;
 		IntArray Map(Nx,Ny,Nz);
-		neighborList= new int[18*Np];
-
+		neighborList= new int[18*Npad];
 		Np = ScaLBL_Comm.MemoryOptimizedLayoutAA(Map,neighborList,Dm.id,Np);
 		MPI_Barrier(comm);
 
 		//......................device distributions.................................
 		int dist_mem_size = Np*sizeof(double);
+		int neighborSize=18*Np*sizeof(int);
 		if (rank==0)	printf ("Allocating distributions \n");
 
 		int *NeighborList;