git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@12240 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp 2014-08-01 21:06:10 +00:00
parent 5ac9e9a785
commit 4fbebcc31d
5 changed files with 169 additions and 147 deletions

View File

@ -504,11 +504,12 @@ void CommBrick::forward_comm(int dummy)
for (int iswap = 0; iswap < nswap; iswap++) { for (int iswap = 0; iswap < nswap; iswap++) {
if (sendproc[iswap] != me) { if (sendproc[iswap] != me) {
if (comm_x_only) { if (comm_x_only) {
if (size_forward_recv[iswap]) buf = x[firstrecv[iswap]]; if (size_forward_recv[iswap]) {
else buf = NULL; if (size_forward_recv[iswap]) buf = x[firstrecv[iswap]];
if (size_forward_recv[iswap]) else buf = NULL;
MPI_Irecv(buf,size_forward_recv[iswap],MPI_DOUBLE, MPI_Irecv(buf,size_forward_recv[iswap],MPI_DOUBLE,
recvproc[iswap],0,world,&request); recvproc[iswap],0,world,&request);
}
n = avec->pack_comm(sendnum[iswap],sendlist[iswap], n = avec->pack_comm(sendnum[iswap],sendlist[iswap],
buf_send,pbc_flag[iswap],pbc[iswap]); buf_send,pbc_flag[iswap],pbc[iswap]);
if (n) MPI_Send(buf_send,n,MPI_DOUBLE,sendproc[iswap],0,world); if (n) MPI_Send(buf_send,n,MPI_DOUBLE,sendproc[iswap],0,world);
@ -575,11 +576,12 @@ void CommBrick::reverse_comm()
if (size_reverse_recv[iswap]) if (size_reverse_recv[iswap])
MPI_Irecv(buf_recv,size_reverse_recv[iswap],MPI_DOUBLE, MPI_Irecv(buf_recv,size_reverse_recv[iswap],MPI_DOUBLE,
sendproc[iswap],0,world,&request); sendproc[iswap],0,world,&request);
if (size_reverse_send[iswap]) buf = f[firstrecv[iswap]]; if (size_reverse_send[iswap]) {
else buf = NULL; if (size_reverse_send[iswap]) buf = f[firstrecv[iswap]];
if (size_reverse_send[iswap]) else buf = NULL;
MPI_Send(buf,size_reverse_send[iswap],MPI_DOUBLE, MPI_Send(buf,size_reverse_send[iswap],MPI_DOUBLE,
recvproc[iswap],0,world); recvproc[iswap],0,world);
}
if (size_reverse_recv[iswap]) MPI_Wait(&request,&status); if (size_reverse_recv[iswap]) MPI_Wait(&request,&status);
} else { } else {
if (size_reverse_recv[iswap]) if (size_reverse_recv[iswap])
@ -620,7 +622,7 @@ void CommBrick::exchange()
int i,m,nsend,nrecv,nrecv1,nrecv2,nlocal; int i,m,nsend,nrecv,nrecv1,nrecv2,nlocal;
double lo,hi,value; double lo,hi,value;
double **x; double **x;
double *sublo,*subhi,*buf; double *sublo,*subhi;
MPI_Request request; MPI_Request request;
MPI_Status status; MPI_Status status;
AtomVec *avec = atom->avec; AtomVec *avec = atom->avec;
@ -656,7 +658,9 @@ void CommBrick::exchange()
// loop over dimensions // loop over dimensions
for (int dim = 0; dim < 3; dim++) { int dimension = domain->dimension;
for (int dim = 0; dim < dimension; dim++) {
// fill buffer with atoms leaving my box, using < and >= // fill buffer with atoms leaving my box, using < and >=
// when atom is deleted, fill it in with last atom // when atom is deleted, fill it in with last atom
@ -679,16 +683,13 @@ void CommBrick::exchange()
// send/recv atoms in both directions // send/recv atoms in both directions
// send size of message first so receiver can realloc buf_recv if needed // send size of message first so receiver can realloc buf_recv if needed
// if 1 proc in dimension, no send/recv, set recv buf to send buf // if 1 proc in dimension, no send/recv
// set nrecv = 0 so buf_send atoms will be lost
// if 2 procs in dimension, single send/recv // if 2 procs in dimension, single send/recv
// if more than 2 procs in dimension, send/recv to both neighbors // if more than 2 procs in dimension, send/recv to both neighbors
if (procgrid[dim] == 1) { if (procgrid[dim] == 1) nrecv = 0;
nrecv = nsend; // NOTE: could just be nrecv = 0 ?? else {
// no buf, just buf_recv
buf = buf_send;
} else {
MPI_Sendrecv(&nsend,1,MPI_INT,procneigh[dim][0],0, MPI_Sendrecv(&nsend,1,MPI_INT,procneigh[dim][0],0,
&nrecv1,1,MPI_INT,procneigh[dim][1],0,world,&status); &nrecv1,1,MPI_INT,procneigh[dim][1],0,world,&status);
nrecv = nrecv1; nrecv = nrecv1;
@ -710,19 +711,18 @@ void CommBrick::exchange()
MPI_Send(buf_send,nsend,MPI_DOUBLE,procneigh[dim][1],0,world); MPI_Send(buf_send,nsend,MPI_DOUBLE,procneigh[dim][1],0,world);
MPI_Wait(&request,&status); MPI_Wait(&request,&status);
} }
buf = buf_recv;
} }
// check incoming atoms to see if they are in my box // check incoming atoms to see if they are in my box
// if so, add to my list // if so, add to my list
// check is only for this dimension, may be passed to another proc // check is only for this dimension,
// may be passed to another proc on later dims
m = 0; m = 0;
while (m < nrecv) { while (m < nrecv) {
value = buf[m+dim+1]; value = buf_recv[m+dim+1];
if (value >= lo && value < hi) m += avec->unpack_exchange(&buf[m]); if (value >= lo && value < hi) m += avec->unpack_exchange(&buf_recv[m]);
else m += static_cast<int> (buf[m]); else m += static_cast<int> (buf_recv[m]);
} }
} }

View File

@ -239,9 +239,9 @@ void CommTiled::setup()
double lo1[3],hi1[3],lo2[3],hi2[3]; double lo1[3],hi1[3],lo2[3],hi2[3];
int one,two; int one,two;
nswap = 0; int iswap = 0;
for (int idim = 0; idim < dimension; idim++) { for (int idim = 0; idim < dimension; idim++) {
for (int iswap = 0; iswap < 2; iswap++) { for (int idir = 0; idir < 2; idir++) {
// one = first ghost box in same periodic image // one = first ghost box in same periodic image
// two = second ghost box wrapped across periodic boundary // two = second ghost box wrapped across periodic boundary
@ -250,7 +250,7 @@ void CommTiled::setup()
one = 1; one = 1;
lo1[0] = sublo[0]; lo1[1] = sublo[1]; lo1[2] = sublo[2]; lo1[0] = sublo[0]; lo1[1] = sublo[1]; lo1[2] = sublo[2];
hi1[0] = subhi[0]; hi1[1] = subhi[1]; hi1[2] = subhi[2]; hi1[0] = subhi[0]; hi1[1] = subhi[1]; hi1[2] = subhi[2];
if (iswap == 0) { if (idir == 0) {
lo1[idim] = sublo[idim] - cut; lo1[idim] = sublo[idim] - cut;
hi1[idim] = sublo[idim]; hi1[idim] = sublo[idim];
} else { } else {
@ -259,13 +259,13 @@ void CommTiled::setup()
} }
two = 0; two = 0;
if (iswap == 0 && periodicity[idim] && lo1[idim] < boxlo[idim]) two = 1; if (idir == 0 && periodicity[idim] && lo1[idim] < boxlo[idim]) two = 1;
if (iswap == 1 && periodicity[idim] && hi1[idim] > boxhi[idim]) two = 1; if (idir == 1 && periodicity[idim] && hi1[idim] > boxhi[idim]) two = 1;
if (two) { if (two) {
lo2[0] = sublo[0]; lo2[1] = sublo[1]; lo2[2] = sublo[2]; lo2[0] = sublo[0]; lo2[1] = sublo[1]; lo2[2] = sublo[2];
hi2[0] = subhi[0]; hi2[1] = subhi[1]; hi2[2] = subhi[2]; hi2[0] = subhi[0]; hi2[1] = subhi[1]; hi2[2] = subhi[2];
if (iswap == 0) { if (idir == 0) {
lo2[idim] = lo1[idim] + prd[idim]; lo2[idim] = lo1[idim] + prd[idim];
hi2[idim] = hi1[idim] + prd[idim]; hi2[idim] = hi1[idim] + prd[idim];
if (sublo[idim] == boxlo[idim]) { if (sublo[idim] == boxlo[idim]) {
@ -301,44 +301,44 @@ void CommTiled::setup()
} }
// reallocate 2nd dimensions of all send/recv arrays, based on noverlap // reallocate 2nd dimensions of all send/recv arrays, based on noverlap
// # of sends of this swap = # of recvs of nswap +/- 1 // # of sends of this swap = # of recvs of iswap +/- 1
if (noverlap > nprocmax[nswap]) { if (noverlap > nprocmax[iswap]) {
int oldmax = nprocmax[nswap]; int oldmax = nprocmax[iswap];
while (nprocmax[nswap] < noverlap) nprocmax[nswap] += DELTA_PROCS; while (nprocmax[iswap] < noverlap) nprocmax[iswap] += DELTA_PROCS;
grow_swap_send(nswap,nprocmax[nswap],oldmax); grow_swap_send(iswap,nprocmax[iswap],oldmax);
if (iswap == 0) grow_swap_recv(nswap+1,nprocmax[nswap]); if (idir == 0) grow_swap_recv(iswap+1,nprocmax[iswap]);
else grow_swap_recv(nswap-1,nprocmax[nswap]); else grow_swap_recv(iswap-1,nprocmax[iswap]);
} }
// overlap how has list of noverlap procs // overlap how has list of noverlap procs
// includes PBC effects // includes PBC effects
if (overlap[noverlap-1] == me) sendself[nswap] = 1; if (overlap[noverlap-1] == me) sendself[iswap] = 1;
else sendself[nswap] = 0; else sendself[iswap] = 0;
if (noverlap-sendself[nswap]) sendother[nswap] = 1; if (noverlap-sendself[iswap]) sendother[iswap] = 1;
else sendother[nswap] = 0; else sendother[iswap] = 0;
//MPI_Barrier(world); //MPI_Barrier(world);
//printf("AAA nswap %d me %d: noverlap %d: %g %g: %g %g\n", //printf("AAA idir %d me %d: noverlap %d: %g %g: %g %g\n",
// nswap,me,noverlap,sublo[0],sublo[1],subhi[0],subhi[1]); // idir,me,noverlap,sublo[0],sublo[1],subhi[0],subhi[1]);
//if (nswap == 0) error->all(FLERR,"ALL DONE"); //if (idir == 0) error->all(FLERR,"ALL DONE");
nsendproc[nswap] = noverlap; nsendproc[iswap] = noverlap;
for (i = 0; i < noverlap; i++) sendproc[nswap][i] = overlap[i]; for (i = 0; i < noverlap; i++) sendproc[iswap][i] = overlap[i];
if (iswap == 0) { if (idir == 0) {
nrecvproc[nswap+1] = noverlap; nrecvproc[iswap+1] = noverlap;
for (i = 0; i < noverlap; i++) recvproc[nswap+1][i] = overlap[i]; for (i = 0; i < noverlap; i++) recvproc[iswap+1][i] = overlap[i];
} else { } else {
nrecvproc[nswap-1] = noverlap; nrecvproc[iswap-1] = noverlap;
for (i = 0; i < noverlap; i++) recvproc[nswap-1][i] = overlap[i]; for (i = 0; i < noverlap; i++) recvproc[iswap-1][i] = overlap[i];
} }
// compute sendbox for each of my sends // compute sendbox for each of my sends
// obox = intersection of ghostbox with other proc's sub-domain // obox = intersection of ghostbox with other proc's sub-domain
// sbox = what I need to send to other proc // sbox = what I need to send to other proc
// = sublo to MIN(sublo+cut,subhi) in idim, for iswap = 0 // = sublo to MIN(sublo+cut,subhi) in idim, for idir = 0
// = MIN(subhi-cut,sublo) to subhi in idim, for iswap = 1 // = MIN(subhi-cut,sublo) to subhi in idim, for idir = 1
// = obox in other 2 dims // = obox in other 2 dims
// if sbox touches sub-box boundaries in lower dims, // if sbox touches sub-box boundaries in lower dims,
// extend sbox in those lower dims to include ghost atoms // extend sbox in those lower dims to include ghost atoms
@ -346,11 +346,11 @@ void CommTiled::setup()
double oboxlo[3],oboxhi[3],sbox[6]; double oboxlo[3],oboxhi[3],sbox[6];
for (i = 0; i < noverlap; i++) { for (i = 0; i < noverlap; i++) {
pbc_flag[nswap][i] = 0; pbc_flag[iswap][i] = 0;
pbc[nswap][i][0] = pbc[nswap][i][1] = pbc[nswap][i][2] = pbc[iswap][i][0] = pbc[iswap][i][1] = pbc[iswap][i][2] =
pbc[nswap][i][3] = pbc[nswap][i][4] = pbc[nswap][i][5] = 0; pbc[iswap][i][3] = pbc[iswap][i][4] = pbc[iswap][i][5] = 0;
(this->*box_other)(idim,iswap,overlap[i],oboxlo,oboxhi); (this->*box_other)(idim,idir,overlap[i],oboxlo,oboxhi);
if (i < noverlap1) { if (i < noverlap1) {
sbox[0] = MAX(oboxlo[0],lo1[0]); sbox[0] = MAX(oboxlo[0],lo1[0]);
@ -360,9 +360,9 @@ void CommTiled::setup()
sbox[4] = MIN(oboxhi[1],hi1[1]); sbox[4] = MIN(oboxhi[1],hi1[1]);
sbox[5] = MIN(oboxhi[2],hi1[2]); sbox[5] = MIN(oboxhi[2],hi1[2]);
} else { } else {
pbc_flag[nswap][i] = 1; pbc_flag[iswap][i] = 1;
if (iswap == 0) pbc[nswap][i][idim] = 1; if (idir == 0) pbc[iswap][i][idim] = 1;
else pbc[nswap][i][idim] = -1; else pbc[iswap][i][idim] = -1;
sbox[0] = MAX(oboxlo[0],lo2[0]); sbox[0] = MAX(oboxlo[0],lo2[0]);
sbox[1] = MAX(oboxlo[1],lo2[1]); sbox[1] = MAX(oboxlo[1],lo2[1]);
sbox[2] = MAX(oboxlo[2],lo2[2]); sbox[2] = MAX(oboxlo[2],lo2[2]);
@ -371,7 +371,7 @@ void CommTiled::setup()
sbox[5] = MIN(oboxhi[2],hi2[2]); sbox[5] = MIN(oboxhi[2],hi2[2]);
} }
if (iswap == 0) { if (idir == 0) {
sbox[idim] = sublo[idim]; sbox[idim] = sublo[idim];
if (i < noverlap1) sbox[3+idim] = MIN(sbox[3+idim]+cut,subhi[idim]); if (i < noverlap1) sbox[3+idim] = MIN(sbox[3+idim]+cut,subhi[idim]);
else sbox[3+idim] = MIN(sbox[3+idim]-prd[idim]+cut,subhi[idim]); else sbox[3+idim] = MIN(sbox[3+idim]-prd[idim]+cut,subhi[idim]);
@ -390,10 +390,10 @@ void CommTiled::setup()
if (sbox[4] == subhi[1]) sbox[4] += cut; if (sbox[4] == subhi[1]) sbox[4] += cut;
} }
memcpy(sendbox[nswap][i],sbox,6*sizeof(double)); memcpy(sendbox[iswap][i],sbox,6*sizeof(double));
} }
nswap++; iswap++;
} }
} }
@ -404,56 +404,51 @@ void CommTiled::setup()
// sets nesendproc, nerecvproc, esendproc, erecvproc // sets nesendproc, nerecvproc, esendproc, erecvproc
// resets neprocmax // resets neprocmax
// NOTE: should there be a unique neprocmax? iswap = 0;
printf("SUBBOX %d: %g %g: %g %g\n",me,sublo[0],sublo[1],subhi[0],subhi[1]);
MPI_Barrier(world);
nswap = 0;
for (int idim = 0; idim < dimension; idim++) { for (int idim = 0; idim < dimension; idim++) {
for (int iswap = 0; iswap < 2; iswap++) { for (int idir = 0; idir < 2; idir++) {
noverlap = 0; noverlap = 0;
n = nsendproc[nswap]; n = nsendproc[iswap];
for (i = 0; i < n; i++) { for (i = 0; i < n; i++) {
if (sendproc[nswap][i] == me) continue; if (sendproc[iswap][i] == me) continue;
if ((this->*box_touch)(sendproc[nswap][i],idim,iswap)) if ((this->*box_touch)(sendproc[iswap][i],idim,idir))
overlap[noverlap++] = sendproc[nswap][i]; overlap[noverlap++] = sendproc[iswap][i];
} }
// reallocate esendproc or erecvproc if needed based on novlerap // reallocate esendproc or erecvproc if needed based on novlerap
if (noverlap > neprocmax[nswap]) { if (noverlap > neprocmax[iswap]) {
while (neprocmax[nswap] < noverlap) neprocmax[nswap] += DELTA_PROCS; while (neprocmax[iswap] < noverlap) neprocmax[iswap] += DELTA_PROCS;
n = neprocmax[nswap]; n = neprocmax[iswap];
delete [] esendproc[nswap]; delete [] esendproc[iswap];
esendproc[nswap] = new int[n]; esendproc[iswap] = new int[n];
if (iswap == 0) { if (idir == 0) {
delete [] erecvproc[nswap+1]; delete [] erecvproc[iswap+1];
erecvproc[nswap+1] = new int[n]; erecvproc[iswap+1] = new int[n];
} else { } else {
delete [] erecvproc[nswap-1]; delete [] erecvproc[iswap-1];
erecvproc[nswap-1] = new int[n]; erecvproc[iswap-1] = new int[n];
} }
} }
nesendproc[nswap] = noverlap; nesendproc[iswap] = noverlap;
for (i = 0; i < noverlap; i++) esendproc[nswap][i] = overlap[i]; for (i = 0; i < noverlap; i++) esendproc[iswap][i] = overlap[i];
if (iswap == 0) { if (idir == 0) {
nerecvproc[nswap+1] = noverlap; nerecvproc[iswap+1] = noverlap;
for (i = 0; i < noverlap; i++) erecvproc[nswap+1][i] = overlap[i]; for (i = 0; i < noverlap; i++) erecvproc[iswap+1][i] = overlap[i];
} else { } else {
nerecvproc[nswap-1] = noverlap; nerecvproc[iswap-1] = noverlap;
for (i = 0; i < noverlap; i++) erecvproc[nswap-1][i] = overlap[i]; for (i = 0; i < noverlap; i++) erecvproc[iswap-1][i] = overlap[i];
} }
nswap++; iswap++;
} }
} }
// reallocate MPI Requests and Statuses as needed // reallocate MPI Requests and Statuses as needed
int nmax = 0; int nmax = 0;
for (i = 0; i < nswap; i++) nmax = MAX(nmax,nprocmax[i]); for (iswap = 0; iswap < nswap; iswap++) nmax = MAX(nmax,nprocmax[iswap]);
if (nmax > maxreqstat) { if (nmax > maxreqstat) {
maxreqstat = nmax; maxreqstat = nmax;
delete [] requests; delete [] requests;
@ -465,6 +460,9 @@ void CommTiled::setup()
// DEBUG // DEBUG
/* /*
printf("SUBBOX %d: %g %g: %g %g\n",me,sublo[0],sublo[1],subhi[0],subhi[1]);
MPI_Barrier(world);
for (i = 0; i < nswap; i++) { for (i = 0; i < nswap; i++) {
if (nsendproc[i] == 1) if (nsendproc[i] == 1)
printf("SETUP SEND %d %d: nsend %d self %d sproc0 %d: " printf("SETUP SEND %d %d: nsend %d self %d sproc0 %d: "
@ -703,15 +701,11 @@ void CommTiled::reverse_comm()
void CommTiled::exchange() void CommTiled::exchange()
{ {
int i,m,nsend,nrecv,nsendsize,nrecvsize,nlocal,dim,offset; int i,m,nsend,nrecv,nsendsize,nrecvsize,nlocal,proc,offset;
double lo,hi,value; double lo,hi,value;
double **x; double **x;
AtomVec *avec = atom->avec; AtomVec *avec = atom->avec;
MPI_Barrier(world);
printf("PREEXCH %d %d\n",me,atom->nlocal);
MPI_Barrier(world);
// clear global->local map for owned and ghost atoms // clear global->local map for owned and ghost atoms
// b/c atoms migrate to new procs in exchange() and // b/c atoms migrate to new procs in exchange() and
// new ghosts are created in borders() // new ghosts are created in borders()
@ -731,8 +725,13 @@ void CommTiled::exchange()
if (bufextra > bufextra_old) if (bufextra > bufextra_old)
memory->grow(buf_send,maxsend+bufextra,"comm:buf_send"); memory->grow(buf_send,maxsend+bufextra,"comm:buf_send");
// domain properties used in exchange method and methods it calls
// subbox bounds for orthogonal or triclinic // subbox bounds for orthogonal or triclinic
prd = domain->prd;
boxlo = domain->boxlo;
boxhi = domain->boxhi;
if (triclinic == 0) { if (triclinic == 0) {
sublo = domain->sublo; sublo = domain->sublo;
subhi = domain->subhi; subhi = domain->subhi;
@ -741,10 +740,11 @@ void CommTiled::exchange()
subhi = domain->subhi_lamda; subhi = domain->subhi_lamda;
} }
// loop over swaps in all dimensions // loop over dimensions
for (int iswap = 0; iswap < nswap; iswap++) { int dimension = domain->dimension;
dim = iswap/2;
for (int dim = 0; dim < dimension; dim++) {
// fill buffer with atoms leaving my box, using < and >= // fill buffer with atoms leaving my box, using < and >=
// when atom is deleted, fill it in with last atom // when atom is deleted, fill it in with last atom
@ -755,28 +755,17 @@ void CommTiled::exchange()
nlocal = atom->nlocal; nlocal = atom->nlocal;
i = nsendsize = 0; i = nsendsize = 0;
if (iswap % 2 == 0) { while (i < nlocal) {
while (i < nlocal) { if (x[i][dim] < lo || x[i][dim] >= hi) {
if (x[i][dim] < lo) { if (nsendsize > maxsend) grow_send(nsendsize,1);
printf("SEND1 from me %d on swap %d: %d: %24.18g %24.18g\n", proc = (this->*point_drop)(dim,x[i]);
me,iswap,atom->tag[i],x[i][dim],lo); if (proc != me) {
if (nsendsize > maxsend) grow_send(nsendsize,1); buf_send[nsendsize++] = proc;
nsendsize += avec->pack_exchange(i,&buf_send[nsendsize]); nsendsize += avec->pack_exchange(i,&buf_send[nsendsize]);
avec->copy(nlocal-1,i,1); avec->copy(nlocal-1,i,1);
nlocal--; nlocal--;
} else i++; } else i++;
} } else i++;
} else {
while (i < nlocal) {
if (x[i][dim] >= hi) {
printf("SEND2 from me %d on swap %d: %d: %24.18g %24.18g\n",
me,iswap,atom->tag[i],x[i][dim],hi);
if (nsendsize > maxsend) grow_send(nsendsize,1);
nsendsize += avec->pack_exchange(i,&buf_send[nsendsize]);
avec->copy(nlocal-1,i,1);
nlocal--;
} else i++;
}
} }
atom->nlocal = nlocal; atom->nlocal = nlocal;
@ -784,52 +773,51 @@ void CommTiled::exchange()
// send and recv atoms from neighbor procs that touch my sub-box in dim // send and recv atoms from neighbor procs that touch my sub-box in dim
// no send/recv with self // no send/recv with self
// send size of message first // send size of message first
// receiver may receive multiple messages, can realloc buf_recv if needed // receiver may receive multiple messages, realloc buf_recv if needed
nsend = nesendproc[iswap]; nsend = nesendproc[dim];
nrecv = nerecvproc[iswap]; nrecv = nerecvproc[dim];
for (m = 0; m < nrecv; m++) for (m = 0; m < nrecv; m++)
MPI_Irecv(&recvnum[iswap][m],1,MPI_INT, MPI_Irecv(&recvnum[dim][m],1,MPI_INT,
erecvproc[iswap][m],0,world,&requests[m]); erecvproc[dim][m],0,world,&requests[m]);
for (m = 0; m < nsend; m++) for (m = 0; m < nsend; m++)
MPI_Send(&nsendsize,1,MPI_INT,esendproc[iswap][m],0,world); MPI_Send(&nsendsize,1,MPI_INT,esendproc[dim][m],0,world);
MPI_Waitall(nrecv,requests,statuses); MPI_Waitall(nrecv,requests,statuses);
nrecvsize = 0; nrecvsize = 0;
for (m = 0; m < nrecv; m++) nrecvsize += recvnum[iswap][m]; for (m = 0; m < nrecv; m++) nrecvsize += recvnum[dim][m];
if (nrecvsize > maxrecv) grow_recv(nrecvsize); if (nrecvsize > maxrecv) grow_recv(nrecvsize);
offset = 0; offset = 0;
for (m = 0; m < nrecv; m++) { for (m = 0; m < nrecv; m++) {
MPI_Irecv(&buf_recv[offset],recvnum[iswap][m], MPI_Irecv(&buf_recv[offset],recvnum[dim][m],
MPI_DOUBLE,erecvproc[iswap][m],0,world,&requests[m]); MPI_DOUBLE,erecvproc[dim][m],0,world,&requests[m]);
offset += recvnum[iswap][m]; offset += recvnum[dim][m];
} }
for (m = 0; m < nsend; m++) for (m = 0; m < nsend; m++)
MPI_Send(buf_send,nsendsize,MPI_DOUBLE,esendproc[iswap][m],0,world); MPI_Send(buf_send,nsendsize,MPI_DOUBLE,esendproc[dim][m],0,world);
MPI_Waitall(nrecv,requests,statuses); MPI_Waitall(nrecv,requests,statuses);
// check incoming atoms to see if they are in my box // check incoming atoms to see if I own it and they are in my box
// if so, add to my list // if so, add to my list
// check is only for this dimension, may be passed to another proc // box check is only for this dimension,
// atom may be passed to another proc in later dims
m = 0; m = 0;
while (m < nrecvsize) { while (m < nrecvsize) {
value = buf_recv[m+dim+1]; proc = static_cast<int> (buf_recv[m++]);
if (value >= lo && value < hi) { if (proc == me) {
m += avec->unpack_exchange(&buf_recv[m]); value = buf_recv[m+dim+1];
printf("RECV from me %d on swap %d: %d\n",me,iswap, if (value >= lo && value < hi) {
atom->tag[atom->nlocal-1]); m += avec->unpack_exchange(&buf_recv[m]);
continue;
}
} }
else m += static_cast<int> (buf_recv[m]); m += static_cast<int> (buf_recv[m]);
} }
} }
MPI_Barrier(world);
printf("POSTEXCH %d %d\n",me,atom->nlocal);
MPI_Barrier(world);
if (atom->firstgroupname) atom->first_reorder(); if (atom->firstgroupname) atom->first_reorder();
} }
@ -1494,7 +1482,7 @@ void CommTiled::box_drop_tiled_recurse(double *lo, double *hi,
return other box owned by proc as lo/hi corner pts return other box owned by proc as lo/hi corner pts
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void CommTiled::box_other_brick(int idim, int iswap, void CommTiled::box_other_brick(int idim, int idir,
int proc, double *lo, double *hi) int proc, double *lo, double *hi)
{ {
lo[0] = sublo[0]; lo[1] = sublo[1]; lo[2] = sublo[2]; lo[0] = sublo[0]; lo[1] = sublo[1]; lo[2] = sublo[2];
@ -1515,7 +1503,7 @@ void CommTiled::box_other_brick(int idim, int iswap,
} }
int dir = -1; int dir = -1;
if (iswap) dir = 1; if (idir) dir = 1;
int index = myloc[idim]; int index = myloc[idim];
int n = procgrid[idim]; int n = procgrid[idim];
@ -1542,7 +1530,7 @@ void CommTiled::box_other_brick(int idim, int iswap,
return other box owned by proc as lo/hi corner pts return other box owned by proc as lo/hi corner pts
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void CommTiled::box_other_tiled(int idim, int iswap, void CommTiled::box_other_tiled(int idim, int idir,
int proc, double *lo, double *hi) int proc, double *lo, double *hi)
{ {
double (*split)[2] = rcbinfo[proc].mysplit; double (*split)[2] = rcbinfo[proc].mysplit;
@ -1565,9 +1553,9 @@ void CommTiled::box_other_tiled(int idim, int iswap,
procneigh stores 6 procs that touch me procneigh stores 6 procs that touch me
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
int CommTiled::box_touch_brick(int proc, int idim, int iswap) int CommTiled::box_touch_brick(int proc, int idim, int idir)
{ {
if (procneigh[idim][iswap] == proc) return 1; if (procneigh[idim][idir] == proc) return 1;
return 0; return 0;
} }
@ -1575,12 +1563,12 @@ int CommTiled::box_touch_brick(int proc, int idim, int iswap)
return 1 if proc's box touches me, else 0 return 1 if proc's box touches me, else 0
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
int CommTiled::box_touch_tiled(int proc, int idim, int iswap) int CommTiled::box_touch_tiled(int proc, int idim, int idir)
{ {
// sending to left // sending to left
// only touches if proc hi = my lo, or if proc hi = boxhi and my lo = boxlo // only touches if proc hi = my lo, or if proc hi = boxhi and my lo = boxlo
if (iswap == 0) { if (idir == 0) {
if (rcbinfo[proc].mysplit[idim][1] == rcbinfo[me].mysplit[idim][0]) if (rcbinfo[proc].mysplit[idim][1] == rcbinfo[me].mysplit[idim][0])
return 1; return 1;
else if (rcbinfo[proc].mysplit[idim][1] == 1.0 && else if (rcbinfo[proc].mysplit[idim][1] == 1.0 &&
@ -1601,6 +1589,33 @@ int CommTiled::box_touch_tiled(int proc, int idim, int iswap)
return 0; return 0;
} }
/* ----------------------------------------------------------------------
------------------------------------------------------------------------- */
int CommTiled::point_drop_brick(int idim, double *x)
{
double deltalo,deltahi;
if (sublo[idim] == boxlo[idim])
deltalo = fabs(x[idim]-prd[idim] - sublo[idim]);
else deltalo = fabs(x[idim] - sublo[idim]);
if (subhi[idim] == boxhi[idim])
deltahi = fabs(x[idim]+prd[idim] - subhi[idim]);
else deltahi = fabs(x[idim] - subhi[idim]);
if (deltalo < deltahi) return procneigh[idim][0];
return procneigh[idim][1];
}
/* ----------------------------------------------------------------------
------------------------------------------------------------------------- */
int CommTiled::point_drop_tiled(int idim, double *x)
{
return 0;
}
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
realloc the size of the send buffer as needed with BUFFACTOR and bufextra realloc the size of the send buffer as needed with BUFFACTOR and bufextra
if flag = 1, realloc if flag = 1, realloc

View File

@ -133,6 +133,11 @@ class CommTiled : public Comm {
int box_touch_brick(int, int, int); int box_touch_brick(int, int, int);
int box_touch_tiled(int, int, int); int box_touch_tiled(int, int, int);
typedef int (CommTiled::*PointDropPtr)(int, double *);
PointDropPtr point_drop;
int point_drop_brick(int, double *);
int point_drop_tiled(int, double *);
void grow_send(int, int); // reallocate send buffer void grow_send(int, int); // reallocate send buffer
void grow_recv(int); // free/allocate recv buffer void grow_recv(int); // free/allocate recv buffer
void grow_list(int, int, int); // reallocate sendlist for one swap/proc void grow_list(int, int, int); // reallocate sendlist for one swap/proc

View File

@ -172,6 +172,9 @@ void FixBalance::setup(int vflag)
void FixBalance::setup_pre_exchange() void FixBalance::setup_pre_exchange()
{ {
// insure atoms are in current box & update box via shrink-wrap // insure atoms are in current box & update box via shrink-wrap
// has to be be done before rebalance() invokes Irregular::migrate_atoms()
// since it requires atoms be inside simulation box
// even though pbc() will be done again in Verlet::run()
// no exchange() since doesn't matter if atoms are assigned to correct procs // no exchange() since doesn't matter if atoms are assigned to correct procs
if (domain->triclinic) domain->x2lamda(atom->nlocal); if (domain->triclinic) domain->x2lamda(atom->nlocal);

View File

@ -102,7 +102,6 @@ RCB::~RCB()
all proc particles will be inside or on surface of 3-d box all proc particles will be inside or on surface of 3-d box
defined by final lo/hi defined by final lo/hi
// NOTE: worry about re-use of data structs for fix balance? // NOTE: worry about re-use of data structs for fix balance?
// NOTE: could get rid of wt all together, will it be used?
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void RCB::compute(int dimension, int n, double **x, double *wt, void RCB::compute(int dimension, int n, double **x, double *wt,