static void nsq_prepare_comm(GhostCommunicator *comm, int data_parts) { int n; /* no need for comm for only 1 node */ if (n_nodes == 1) { prepare_comm(comm, data_parts, 0); return; } prepare_comm(comm, data_parts, n_nodes); /* every node has its dedicated comm step */ for(n = 0; n < n_nodes; n++) { comm->comm[n].part_lists = malloc(sizeof(ParticleList *)); comm->comm[n].part_lists[0] = &cells[n]; comm->comm[n].n_part_lists = 1; comm->comm[n].node = n; comm->comm[n].mpi_comm = MPI_COMM_WORLD; } }
/** Create communicators for cell structure domain decomposition. (see \ref GhostCommunicator) */ void dd_prepare_comm(GhostCommunicator *comm, int data_parts) { int dir,lr,i,cnt, num, n_comm_cells[3]; int lc[3],hc[3],done[3]={0,0,0}; /* calculate number of communications */ num = 0; for(dir=0; dir<3; dir++) { for(lr=0; lr<2; lr++) { #ifdef PARTIAL_PERIODIC /* No communication for border of non periodic direction */ if( PERIODIC(dir) || (boundary[2*dir+lr] == 0) ) #endif { if(node_grid[dir] == 1 ) num++; else num += 2; } } } /* prepare communicator */ CELL_TRACE(fprintf(stderr,"%d Create Communicator: prep_comm data_parts %d num %d\n",this_node,data_parts,num)); prepare_comm(comm, data_parts, num); /* number of cells to communicate in a direction */ n_comm_cells[0] = dd.cell_grid[1] * dd.cell_grid[2]; n_comm_cells[1] = dd.cell_grid[2] * dd.ghost_cell_grid[0]; n_comm_cells[2] = dd.ghost_cell_grid[0] * dd.ghost_cell_grid[1]; cnt=0; /* direction loop: x, y, z */ for(dir=0; dir<3; dir++) { lc[(dir+1)%3] = 1-done[(dir+1)%3]; lc[(dir+2)%3] = 1-done[(dir+2)%3]; hc[(dir+1)%3] = dd.cell_grid[(dir+1)%3]+done[(dir+1)%3]; hc[(dir+2)%3] = dd.cell_grid[(dir+2)%3]+done[(dir+2)%3]; /* lr loop: left right */ /* here we could in principle build in a one sided ghost communication, simply by taking the lr loop only over one value */ for(lr=0; lr<2; lr++) { if(node_grid[dir] == 1) { /* just copy cells on a single node */ #ifdef PARTIAL_PERIODIC if( PERIODIC(dir ) || (boundary[2*dir+lr] == 0) ) #endif { comm->comm[cnt].type = GHOST_LOCL; comm->comm[cnt].node = this_node; /* Buffer has to contain Send and Recv cells -> factor 2 */ comm->comm[cnt].part_lists = malloc(2*n_comm_cells[dir]*sizeof(ParticleList *)); comm->comm[cnt].n_part_lists = 2*n_comm_cells[dir]; /* prepare folding of ghost positions */ if((data_parts & GHOSTTRANS_POSSHFTD) && boundary[2*dir+lr] != 0) comm->comm[cnt].shift[dir] = boundary[2*dir+lr]*box_l[dir]; /* fill send comm cells */ lc[(dir+0)%3] = hc[(dir+0)%3] = 1+lr*(dd.cell_grid[(dir+0)%3]-1); dd_fill_comm_cell_lists(comm->comm[cnt].part_lists,lc,hc); CELL_TRACE(fprintf(stderr,"%d: prep_comm %d copy to grid (%d,%d,%d)-(%d,%d,%d)\n",this_node,cnt, lc[0],lc[1],lc[2],hc[0],hc[1],hc[2])); /* fill recv comm cells */ lc[(dir+0)%3] = hc[(dir+0)%3] = 0+(1-lr)*(dd.cell_grid[(dir+0)%3]+1); /* place recieve cells after send cells */ dd_fill_comm_cell_lists(&comm->comm[cnt].part_lists[n_comm_cells[dir]],lc,hc); CELL_TRACE(fprintf(stderr,"%d: prep_comm %d copy from grid (%d,%d,%d)-(%d,%d,%d)\n",this_node,cnt,lc[0],lc[1],lc[2],hc[0],hc[1],hc[2])); cnt++; } } else { /* i: send/recv loop */ for(i=0; i<2; i++) { #ifdef PARTIAL_PERIODIC if( PERIODIC(dir) || (boundary[2*dir+lr] == 0) ) #endif if((node_pos[dir]+i)%2==0) { comm->comm[cnt].type = GHOST_SEND; comm->comm[cnt].node = node_neighbors[2*dir+lr]; comm->comm[cnt].part_lists = malloc(n_comm_cells[dir]*sizeof(ParticleList *)); comm->comm[cnt].n_part_lists = n_comm_cells[dir]; /* prepare folding of ghost positions */ if((data_parts & GHOSTTRANS_POSSHFTD) && boundary[2*dir+lr] != 0) comm->comm[cnt].shift[dir] = boundary[2*dir+lr]*box_l[dir]; lc[(dir+0)%3] = hc[(dir+0)%3] = 1+lr*(dd.cell_grid[(dir+0)%3]-1); dd_fill_comm_cell_lists(comm->comm[cnt].part_lists,lc,hc); CELL_TRACE(fprintf(stderr,"%d: prep_comm %d send to node %d grid (%d,%d,%d)-(%d,%d,%d)\n",this_node,cnt, comm->comm[cnt].node,lc[0],lc[1],lc[2],hc[0],hc[1],hc[2])); cnt++; } #ifdef PARTIAL_PERIODIC if( PERIODIC(dir) || (boundary[2*dir+(1-lr)] == 0) ) #endif if((node_pos[dir]+(1-i))%2==0) { comm->comm[cnt].type = GHOST_RECV; comm->comm[cnt].node = node_neighbors[2*dir+(1-lr)]; comm->comm[cnt].part_lists = malloc(n_comm_cells[dir]*sizeof(ParticleList *)); comm->comm[cnt].n_part_lists = n_comm_cells[dir]; lc[(dir+0)%3] = hc[(dir+0)%3] = 0+(1-lr)*(dd.cell_grid[(dir+0)%3]+1); dd_fill_comm_cell_lists(comm->comm[cnt].part_lists,lc,hc); CELL_TRACE(fprintf(stderr,"%d: prep_comm %d recv from node %d grid (%d,%d,%d)-(%d,%d,%d)\n",this_node,cnt, comm->comm[cnt].node,lc[0],lc[1],lc[2],hc[0],hc[1],hc[2])); cnt++; } } } done[dir]=1; } } }