fix bug with slab geometries

This commit is contained in:
Steve Plimpton
2022-11-30 11:55:19 -07:00
parent cc18528ea1
commit 7ce4b2eb68
2 changed files with 17 additions and 29 deletions

View File

@ -489,10 +489,11 @@ void Grid2d::ghost_grid()
// ghost cell layers needed in each dim/dir = max of two extension effects // ghost cell layers needed in each dim/dir = max of two extension effects
// OFFSET allows generation of negative indices with static_cast // OFFSET allows generation of negative indices with static_cast
// out xyz lo/hi = index range of owned + ghost cells // out xyz lo/hi = index range of owned + ghost cells
// if yextra, ny and effective prd[1] are both larger, so dyinv is the same
double dxinv = nx / prd[0]; double dxinv = nx / prd[0];
double dyinv = ny / prd[1]; double dyinv = ny / prd[1];
double dyinv_slab = ny / (prd[1] * yfactor); if (yextra) dyinv = ny / (prd[1] * yfactor);
int lo, hi; int lo, hi;
@ -501,17 +502,10 @@ void Grid2d::ghost_grid()
outxlo = MIN(lo - stencil_atom_lo, inxlo - stencil_grid_lo); outxlo = MIN(lo - stencil_atom_lo, inxlo - stencil_grid_lo);
outxhi = MAX(hi + stencil_atom_hi, inxhi + stencil_grid_hi); outxhi = MAX(hi + stencil_atom_hi, inxhi + stencil_grid_hi);
if (yextra == 0) { lo = static_cast<int>((sublo[1]-dist[1]-boxlo[1]) * dyinv + shift_atom_lo + OFFSET) - OFFSET;
lo = static_cast<int>((sublo[1]-dist[1]-boxlo[1]) * dyinv + shift_atom_lo + OFFSET) - OFFSET; hi = static_cast<int>((subhi[1]+dist[1]-boxlo[1]) * dyinv + shift_atom_hi + OFFSET) - OFFSET;
hi = static_cast<int>((subhi[1]+dist[1]-boxlo[1]) * dyinv + shift_atom_hi + OFFSET) - OFFSET; outylo = MIN(lo - stencil_atom_lo, inylo - stencil_grid_lo);
outylo = MIN(lo - stencil_atom_lo, inylo - stencil_grid_lo); outyhi = MAX(hi + stencil_atom_hi, inyhi + stencil_grid_hi);
outyhi = MAX(hi + stencil_atom_hi, inyhi + stencil_grid_hi);
} else {
lo = static_cast<int>((sublo[1]-dist[1]-boxlo[1]) * dyinv_slab + shift_atom_lo + OFFSET) - OFFSET;
hi = static_cast<int>((subhi[1]+dist[1]-boxlo[1]) * dyinv_slab + shift_atom_hi + OFFSET) - OFFSET;
outylo = MIN(lo - stencil_atom_lo, inylo - stencil_grid_lo);
outyhi = MAX(hi + stencil_atom_hi, inyhi + stencil_grid_hi);
}
// if yextra = 1: // if yextra = 1:
// adjust grid boundaries for processors at +y end, // adjust grid boundaries for processors at +y end,
@ -532,7 +526,7 @@ void Grid2d::ghost_grid()
} }
// limit out xyz lo/hi indices to global grid for non-periodic dims // limit out xyz lo/hi indices to global grid for non-periodic dims
// if yextra = 1, grid is still fully periodic // if yextra = 1, treat y-dimension as if periodic
int *periodicity = domain->periodicity; int *periodicity = domain->periodicity;
@ -541,7 +535,7 @@ void Grid2d::ghost_grid()
outxhi = MIN(nx-1,outxhi); outxhi = MIN(nx-1,outxhi);
} }
if (!periodicity[1]) { if (!periodicity[1] && !yextra) {
outylo = MAX(0,outylo); outylo = MAX(0,outylo);
outyhi = MIN(ny-1,outyhi); outyhi = MIN(ny-1,outyhi);
} }

View File

@ -536,11 +536,12 @@ void Grid3d::ghost_grid()
// ghost cell layers needed in each dim/dir = max of two extension effects // ghost cell layers needed in each dim/dir = max of two extension effects
// OFFSET allows generation of negative indices with static_cast // OFFSET allows generation of negative indices with static_cast
// out xyz lo/hi = index range of owned + ghost cells // out xyz lo/hi = index range of owned + ghost cells
// if zextra, nz and effective prd[2] are both larger, so dzinv is the same
double dxinv = nx / prd[0]; double dxinv = nx / prd[0];
double dyinv = ny / prd[1]; double dyinv = ny / prd[1];
double dzinv = nz / prd[2]; double dzinv = nz / prd[2];
double dzinv_slab = nz / (prd[2] * zfactor); if (zextra) dzinv = nz / (prd[2] * zfactor);
int lo, hi; int lo, hi;
@ -554,17 +555,10 @@ void Grid3d::ghost_grid()
outylo = MIN(lo - stencil_atom_lo, inylo - stencil_grid_lo); outylo = MIN(lo - stencil_atom_lo, inylo - stencil_grid_lo);
outyhi = MAX(hi + stencil_atom_hi, inyhi + stencil_grid_hi); outyhi = MAX(hi + stencil_atom_hi, inyhi + stencil_grid_hi);
if (zextra == 0) { lo = static_cast<int>((sublo[2]-dist[2]-boxlo[2]) * dzinv + shift_atom_lo + OFFSET) - OFFSET;
lo = static_cast<int>((sublo[2]-dist[2]-boxlo[2]) * dzinv + shift_atom_lo + OFFSET) - OFFSET; hi = static_cast<int>((subhi[2]+dist[2]-boxlo[2]) * dzinv + shift_atom_hi + OFFSET) - OFFSET;
hi = static_cast<int>((subhi[2]+dist[2]-boxlo[2]) * dzinv + shift_atom_hi + OFFSET) - OFFSET; outzlo = MIN(lo - stencil_atom_lo, inzlo - stencil_grid_lo);
outzlo = MIN(lo - stencil_atom_lo, inzlo - stencil_grid_lo); outzhi = MAX(hi + stencil_atom_hi, inzhi + stencil_grid_hi);
outzhi = MAX(hi + stencil_atom_hi, inzhi + stencil_grid_hi);
} else {
lo = static_cast<int>((sublo[2]-dist[2]-boxlo[2]) * dzinv_slab + shift_atom_lo + OFFSET) - OFFSET;
hi = static_cast<int>((subhi[2]+dist[2]-boxlo[2]) * dzinv_slab + shift_atom_hi + OFFSET) - OFFSET;
outzlo = MIN(lo - stencil_atom_lo, inzlo - stencil_grid_lo);
outzhi = MAX(hi + stencil_atom_hi, inzhi + stencil_grid_hi);
}
// if zextra = 1: // if zextra = 1:
// adjust grid boundaries for processors at +z end, // adjust grid boundaries for processors at +z end,
@ -585,7 +579,7 @@ void Grid3d::ghost_grid()
} }
// limit out xyz lo/hi indices to global grid for non-periodic dims // limit out xyz lo/hi indices to global grid for non-periodic dims
// if zextra = 1 (e.g. PPPM), grid is still fully periodic // if zextra = 1 (e.g. PPPM), treat z-dimension as if periodic
int *periodicity = domain->periodicity; int *periodicity = domain->periodicity;
@ -599,7 +593,7 @@ void Grid3d::ghost_grid()
outyhi = MIN(ny-1,outyhi); outyhi = MIN(ny-1,outyhi);
} }
if (!periodicity[2]) { if (!periodicity[2] && !zextra) {
outzlo = MAX(0,outzlo); outzlo = MAX(0,outzlo);
outzhi = MIN(nz-1,outzhi); outzhi = MIN(nz-1,outzhi);
} }