int Zoltan_Comm_Do( ZOLTAN_COMM_OBJ * plan, /* communication data structure */ int tag, /* message tag for communicating */ char *send_data, /* array of data I currently own */ int nbytes, /* multiplier for sizes */ char *recv_data) /* array of data I'll own after comm */ { int status = ZOLTAN_OK; if (!plan->maxed_recvs){ status = Zoltan_Comm_Do_Post (plan, tag, send_data, nbytes, recv_data); if (status == ZOLTAN_OK) status = Zoltan_Comm_Do_Wait (plan, tag, send_data, nbytes, recv_data); } else{ status = Zoltan_Comm_Do_AlltoAll(plan, send_data, nbytes, recv_data); } return status; }
int Zoltan_Comm_Do_Reverse( ZOLTAN_COMM_OBJ *plan, /* communication data structure */ int tag, /* message tag for communicating */ char *send_data, /* array of data I currently own */ int nbytes, /* # bytes per data item */ int *sizes, /* variable size of objects (if not NULL) */ char *recv_data) /* array of data I'll own after reverse comm */ { int status; /* create plan->plan_reverse */ status = create_reverse_plan(plan, tag, sizes); if (status == ZOLTAN_OK || status == ZOLTAN_WARN){ if (plan->plan_reverse->maxed_recvs){ /* use MPI_Alltoallv to implement plan->plan_reverse, because comm_do_post * would post more receives that allowed on this machine */ status = Zoltan_Comm_Do_AlltoAll(plan->plan_reverse, send_data, nbytes, recv_data); } else{ /* use post/wait which is faster when each sends to few */ status = Zoltan_Comm_Do_Post(plan->plan_reverse, tag, send_data, nbytes, recv_data); if (status == ZOLTAN_OK || status == ZOLTAN_WARN){ status = Zoltan_Comm_Do_Wait (plan->plan_reverse, tag, send_data, nbytes, recv_data); } } } free_reverse_plan(plan); return status; }
int Zoltan_Comm_Do_Post( ZOLTAN_COMM_OBJ * plan, /* communication data structure */ int tag, /* message tag for communicating */ char *send_data, /* array of data I currently own */ int nbytes, /* multiplier for sizes */ char *recv_data) /* array of data I'll own after comm */ { char *send_buff; /* space to buffer outgoing data */ int my_proc; /* processor ID */ unsigned int self_recv_address = 0;/* where in recv_data self info starts */ int self_num=0; /* where in send list my_proc appears */ int offset; /* offset into array I'm copying into */ int self_index = 0; /* send offset for data I'm keeping */ int out_of_mem; /* am I out of memory? */ int nblocks; /* number of procs who need my data */ int proc_index; /* loop counter over procs to send to */ int i, j, k, jj; /* loop counters */ static char *yo = "Zoltan_Comm_Do_Post"; /* Check input parameters */ if (!plan) { MPI_Comm_rank(MPI_COMM_WORLD, &my_proc); ZOLTAN_COMM_ERROR("Communication plan = NULL", yo, my_proc); return ZOLTAN_FATAL; } /* If not point to point, currently we do synchroneous communications */ if (plan->maxed_recvs){ int status; status = Zoltan_Comm_Do_AlltoAll(plan, send_data, nbytes, recv_data); return (status); } MPI_Comm_rank(plan->comm, &my_proc); if ((plan->nsends + plan->self_msg) && !send_data) { int sum = 0; if (plan->sizes_to) /* Not an error if all sizes_to == 0 */ for (i = 0; i < (plan->nsends + plan->self_msg); i++) sum += plan->sizes_to[i]; if (!plan->sizes_to || (plan->sizes_to && sum)) { ZOLTAN_COMM_ERROR("nsends not zero, but send_data = NULL", yo, my_proc); return ZOLTAN_FATAL; } } if ((plan->nrecvs + plan->self_msg) && !recv_data) { int sum = 0; if (plan->sizes_from) /* Not an error if all sizes_from == 0 */ for (i = 0; i < (plan->nrecvs + plan->self_msg); i++) sum += plan->sizes_from[i]; if (!plan->sizes_from || (plan->sizes_from && sum)) { ZOLTAN_COMM_ERROR("nrecvs not zero, but recv_data = NULL", yo, my_proc); return ZOLTAN_FATAL; } } if (nbytes < 0) { ZOLTAN_COMM_ERROR("Scale factor nbytes is negative", yo, my_proc); return ZOLTAN_FATAL; } /* Post irecvs */ out_of_mem = 0; if (plan->indices_from == NULL) { /* Data can go directly into user space. */ plan->recv_buff = recv_data; } else { /* Need to buffer receive to reorder */ plan->recv_buff = (char *) ZOLTAN_MALLOC(plan->total_recv_size * nbytes); if (plan->recv_buff == NULL && plan->total_recv_size * nbytes != 0) out_of_mem = 1; } if (!out_of_mem) { if (plan->sizes == NULL) { /* All data the same size */ k = 0; for (i = 0; i < plan->nrecvs + plan->self_msg; i++) { if (plan->procs_from[i] != my_proc) { MPI_Irecv((void *) & plan->recv_buff[plan->starts_from[i] * nbytes], plan->lengths_from[i] * nbytes, (MPI_Datatype) MPI_BYTE, plan->procs_from[i], tag, plan->comm, &plan->request[k]); k++; } else { self_recv_address = plan->starts_from[i] * nbytes; } } } else { /* Data of varying sizes */ k = 0; for (i = 0; i < plan->nrecvs + plan->self_msg; i++) { if (plan->procs_from[i] != my_proc) { if (plan->sizes_from[i]) MPI_Irecv((void *) &plan->recv_buff[plan->starts_from_ptr[i] * nbytes], plan->sizes_from[i] * nbytes, (MPI_Datatype) MPI_BYTE, plan->procs_from[i], tag, plan->comm, &plan->request[k]); else plan->request[k] = MPI_REQUEST_NULL; k++; } else { self_recv_address = plan->starts_from_ptr[i] * nbytes; } } } } /* Do remaining allocation to check for any mem problems. */ if (plan->indices_to != NULL) { /* can't sent straight from input */ send_buff = (char *) ZOLTAN_MALLOC(plan->max_send_size * nbytes); if (send_buff == 0 && plan->max_send_size * nbytes != 0) out_of_mem = 1; } else send_buff = NULL; /* Barrier to ensure irecvs are posted before doing any sends. */ /* Simultaneously see if anyone out of memory */ MPI_Allreduce(&out_of_mem, &j, 1, MPI_INT, MPI_SUM, plan->comm); if (j > 0) { /* Some proc is out of memory -> Punt */ ZOLTAN_FREE(&send_buff); if (plan->indices_from != NULL) ZOLTAN_FREE(&plan->recv_buff); return (ZOLTAN_MEMERR); } /* Send out data */ /* Scan through procs_to list to start w/ higher numbered procs */ /* This should balance message traffic. */ nblocks = plan->nsends + plan->self_msg; proc_index = 0; while (proc_index < nblocks && plan->procs_to[proc_index] < my_proc) proc_index++; if (proc_index == nblocks) proc_index = 0; if (plan->sizes == NULL) { /* Data all of same size */ if (plan->indices_to == NULL) { /* data already blocked by processor. */ for (i = proc_index, j = 0; j < nblocks; j++) { if (plan->procs_to[i] != my_proc) { MPI_Rsend((void *) &send_data[plan->starts_to[i] * nbytes], plan->lengths_to[i] * nbytes, (MPI_Datatype) MPI_BYTE, plan->procs_to[i], tag, plan->comm); } else self_num = i; if (++i == nblocks) i = 0; } if (plan->self_msg) { /* Copy data to self. */ /* I use array+offset instead of &(array[offset]) because of a bug with PGI v9 */ /* I use memmove because I'm not sure that the pointer are not overlapped. */ memmove(plan->recv_buff+self_recv_address, send_data+plan->starts_to[self_num] * nbytes, plan->lengths_to[self_num]*nbytes); } } else { /* Not blocked by processor. Need to buffer. */ for (i = proc_index, jj = 0; jj < nblocks; jj++) { if (plan->procs_to[i] != my_proc) { /* Need to pack message first. */ offset = 0; j = plan->starts_to[i]; for (k = 0; k < plan->lengths_to[i]; k++) { memcpy(&send_buff[offset], &send_data[plan->indices_to[j++] * nbytes], nbytes); offset += nbytes; } MPI_Rsend((void *) send_buff, plan->lengths_to[i] * nbytes, (MPI_Datatype) MPI_BYTE, plan->procs_to[i], tag, plan->comm); } else { self_num = i; self_index = plan->starts_to[i]; } if (++i == nblocks) i = 0; } if (plan->self_msg) { /* Copy data to self. */ for (k = 0; k < plan->lengths_to[self_num]; k++) { memcpy(&plan->recv_buff[self_recv_address], &send_data[plan->indices_to[self_index++] * nbytes], nbytes); self_recv_address += nbytes; } } ZOLTAN_FREE(&send_buff); } } else { /* Data of differing sizes */ if (plan->indices_to == NULL) { /* data already blocked by processor. */ for (i = proc_index, j = 0; j < nblocks; j++) { if (plan->procs_to[i] != my_proc) { if (plan->sizes_to[i]) { MPI_Rsend((void *) &send_data[plan->starts_to_ptr[i] * nbytes], plan->sizes_to[i] * nbytes, (MPI_Datatype) MPI_BYTE, plan->procs_to[i], tag, plan->comm); } } else self_num = i; if (++i == nblocks) i = 0; } if (plan->self_msg) { /* Copy data to self. */ if (plan->sizes_to[self_num]) { char* lrecv = &plan->recv_buff[self_recv_address]; char* lsend = &send_data[plan->starts_to_ptr[self_num] * nbytes]; int sindex = plan->sizes_to[self_num], idx; for (idx=0; idx<nbytes; idx++) { memcpy(lrecv, lsend, sindex); lrecv += sindex; lsend += sindex; } } } } else { /* Not blocked by processor. Need to buffer. */ for (i = proc_index, jj = 0; jj < nblocks; jj++) { if (plan->procs_to[i] != my_proc) { /* Need to pack message first. */ offset = 0; j = plan->starts_to[i]; for (k = 0; k < plan->lengths_to[i]; k++) { if (plan->sizes[plan->indices_to[j]]) { memcpy(&send_buff[offset], &send_data[plan->indices_to_ptr[j] * nbytes], plan->sizes[plan->indices_to[j]] * nbytes); offset += plan->sizes[plan->indices_to[j]] * nbytes; } j++; } if (plan->sizes_to[i]) { MPI_Rsend((void *) send_buff, plan->sizes_to[i] * nbytes, (MPI_Datatype) MPI_BYTE, plan->procs_to[i], tag, plan->comm); } } else self_num = i; if (++i == nblocks) i = 0; } if (plan->self_msg) { /* Copy data to self. */ if (plan->sizes_to[self_num]) { j = plan->starts_to[self_num]; for (k = 0; k < plan->lengths_to[self_num]; k++) { int kk = plan->indices_to_ptr[j]; char* lrecv = &plan->recv_buff[self_recv_address]; unsigned int send_idx = kk * nbytes; char* lsend = &send_data[send_idx]; int sindex = plan->sizes[plan->indices_to[j]], idx; for (idx=0; idx<nbytes; idx++) { memcpy(lrecv, lsend, sindex); lrecv += sindex; lsend += sindex; } self_recv_address += plan->sizes[plan->indices_to[j]] * nbytes; j++; } } } ZOLTAN_FREE(&send_buff); } } return (ZOLTAN_OK); }