int tMPI_Thread_barrier_wait(tMPI_Thread_barrier_t *barrier) { int cycle; BOOL rc=FALSE; int ret=0; /*tMPI_Thread_pthread_barrier_t *p;*/ /* check whether the barrier is initialized */ if (tMPI_Atomic_get( &(barrier->initialized) ) == 0) { tMPI_Thread_barrier_init_once(barrier,barrier->threshold); } #if 0 EnterCriticalSection( &(barrier->barrierp->cs) ); #else tMPI_Thread_mutex_lock( &(barrier->barrierp->cs) ); #endif cycle = barrier->cycle; /* Decrement the count atomically and check if it is zero. * This will only be true for the last thread calling us. */ if( --(barrier->count) <= 0 ) { barrier->cycle = !barrier->cycle; barrier->count = barrier->threshold; #if 0 WakeAllConditionVariable( &(barrier->barrierp->cv) ); #else tMPI_Thread_cond_broadcast( &(barrier->barrierp->cv) ); #endif } else { while(cycle == barrier->cycle) { #if 0 rc=SleepConditionVariableCS (&(barrier->barrierp->cv), &(barrier->barrierp->cs), INFINITE); if(!rc) { ret=-1; break; } #else rc = tMPI_Thread_cond_wait(&barrier->barrierp->cv, &barrier->barrierp->cs); if(rc != 0) break; #endif } } #if 0 LeaveCriticalSection( &(barrier->barrierp->cs) ); #else tMPI_Thread_mutex_unlock( &(barrier->barrierp->cs) ); #endif return ret; }
/* this is the main comm creation function. All other functions that create comms use this*/ int tMPI_Comm_split(tMPI_Comm comm, int color, int key, tMPI_Comm *newcomm) { int i, j; int N = tMPI_Comm_N(comm); volatile tMPI_Comm *newcomm_list; volatile int colors[MAX_PREALLOC_THREADS]; /* array with the colors of each thread */ volatile int keys[MAX_PREALLOC_THREADS]; /* same for keys (only one of the threads actually suplies these arrays to the comm structure) */ tmpi_bool i_am_first = FALSE; int myrank = tMPI_Comm_seek_rank(comm, tMPI_Get_current()); struct tmpi_split *spl; int ret; #ifdef TMPI_TRACE tMPI_Trace_print("tMPI_Comm_split(%p, %d, %d, %p)", comm, color, key, newcomm); #endif if (!comm) { *newcomm = NULL; return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_COMM); } ret = tMPI_Thread_mutex_lock(&(comm->comm_create_lock)); if (ret != 0) { return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_IO); } /* first get the colors */ if (!comm->new_comm) { /* i am apparently first */ comm->split = (struct tmpi_split*)tMPI_Malloc(sizeof(struct tmpi_split)); comm->new_comm = (tMPI_Comm*)tMPI_Malloc(N*sizeof(tMPI_Comm)); if (N <= MAX_PREALLOC_THREADS) { comm->split->colors = colors; comm->split->keys = keys; } else { comm->split->colors = (int*)tMPI_Malloc(N*sizeof(int)); comm->split->keys = (int*)tMPI_Malloc(N*sizeof(int)); } comm->split->Ncol_init = tMPI_Comm_N(comm); comm->split->can_finish = FALSE; i_am_first = TRUE; /* the main communicator contains a list the size of grp.N */ } newcomm_list = comm->new_comm; /* we copy it to the local stacks because we can later erase comm->new_comm safely */ spl = comm->split; /* we do the same for spl */ spl->colors[myrank] = color; spl->keys[myrank] = key; spl->Ncol_init--; if (spl->Ncol_init == 0) { ret = tMPI_Thread_cond_signal(&(comm->comm_create_prep)); if (ret != 0) { return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_IO); } } if (!i_am_first) { /* all other threads can just wait until the creator thread is finished */ while (!spl->can_finish) { ret = tMPI_Thread_cond_wait(&(comm->comm_create_finish), &(comm->comm_create_lock) ); if (ret != 0) { return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_IO); } } } else { int Ncomms = 0; int comm_color_[MAX_PREALLOC_THREADS]; int comm_N_[MAX_PREALLOC_THREADS]; int *comm_color = comm_color_; /* there can't be more comms than N*/ int *comm_N = comm_N_; /* the number of procs in a group */ int *comm_groups; /* the groups */ tMPI_Comm *comms; /* the communicators */ /* wait for the colors to be done */ /*if (N>1)*/ while (spl->Ncol_init > 0) { ret = tMPI_Thread_cond_wait(&(comm->comm_create_prep), &(comm->comm_create_lock)); if (ret != 0) { return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_IO); } } /* reset the state so that a new comm creating function can run */ spl->Ncol_destroy = N; comm->new_comm = 0; comm->split = 0; comm_groups = (int*)tMPI_Malloc(N*N*sizeof(int)); if (N > MAX_PREALLOC_THREADS) { comm_color = (int*)tMPI_Malloc(N*sizeof(int)); comm_N = (int*)tMPI_Malloc(N*sizeof(int)); } /* count colors, allocate and split up communicators */ tMPI_Split_colors(N, (int*)spl->colors, (int*)spl->keys, &Ncomms, comm_N, comm_color, comm_groups); /* allocate a bunch of communicators */ comms = (tMPI_Comm*)tMPI_Malloc(Ncomms*sizeof(tMPI_Comm)); for (i = 0; i < Ncomms; i++) { ret = tMPI_Comm_alloc(&(comms[i]), comm, comm_N[i]); if (ret != TMPI_SUCCESS) { return ret; } } /* now distribute the comms */ for (i = 0; i < Ncomms; i++) { comms[i]->grp.N = comm_N[i]; for (j = 0; j < comm_N[i]; j++) { comms[i]->grp.peers[j] = comm->grp.peers[comm_groups[i*comm->grp.N + j]]; } } /* and put them into the newcomm_list */ for (i = 0; i < N; i++) { newcomm_list[i] = TMPI_COMM_NULL; for (j = 0; j < Ncomms; j++) { if (spl->colors[i] == comm_color[j]) { newcomm_list[i] = comms[j]; break; } } } #ifdef TMPI_DEBUG /* output */ for (i = 0; i < Ncomms; i++) { printf("Group %d (color %d) has %d members: ", i, comm_color[i], comm_N[i]); for (j = 0; j < comm_N[i]; j++) { printf(" %d ", comm_groups[comm->grp.N*i + j]); } printf(" rank: "); for (j = 0; j < comm_N[i]; j++) { printf(" %d ", spl->keys[comm_groups[N*i + j]]); } printf(" color: "); for (j = 0; j < comm_N[i]; j++) { printf(" %d ", spl->colors[comm_groups[N*i + j]]); } printf("\n"); } #endif if (N > MAX_PREALLOC_THREADS) { free((int*)spl->colors); free((int*)spl->keys); free(comm_color); free(comm_N); } free(comm_groups); free(comms); spl->can_finish = TRUE; /* tell the waiting threads that there's a comm ready */ ret = tMPI_Thread_cond_broadcast(&(comm->comm_create_finish)); if (ret != 0) { return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_IO); } } /* here the individual threads get their comm object */ *newcomm = newcomm_list[myrank]; /* free when we have assigned them all, so we can reuse the object*/ spl->Ncol_destroy--; if (spl->Ncol_destroy == 0) { free((void*)newcomm_list); free(spl); } ret = tMPI_Thread_mutex_unlock(&(comm->comm_create_lock)); if (ret != 0) { return tMPI_Error(TMPI_COMM_WORLD, TMPI_ERR_IO); } return TMPI_SUCCESS; }