1
0
mirror of https://github.com/paboyle/Grid.git synced 2025-04-09 21:50:45 +01:00

add space in stencil output

This commit is contained in:
nmeyer-ur 2020-05-29 17:11:17 +02:00
parent 5cb3530c34
commit e947b563ea

View File

@ -1,6 +1,6 @@
/************************************************************************************* /*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/Stencil.h Source file: ./lib/Stencil.h
@ -41,13 +41,13 @@
// Stencil based code will exchange haloes and use a table lookup for neighbours. // Stencil based code will exchange haloes and use a table lookup for neighbours.
// This will be done with generality to allow easier efficient implementations. // This will be done with generality to allow easier efficient implementations.
// Overlap of comms and compute is enabled by tabulating off-node connected, // Overlap of comms and compute is enabled by tabulating off-node connected,
// //
// Generic services // Generic services
// 0) Prebuild neighbour tables // 0) Prebuild neighbour tables
// 1) Compute sizes of all haloes/comms buffers; allocate them. // 1) Compute sizes of all haloes/comms buffers; allocate them.
// 2) Gather all faces, and communicate. // 2) Gather all faces, and communicate.
// 3) Loop over result sites, giving nbr index/offnode info for each // 3) Loop over result sites, giving nbr index/offnode info for each
// //
////////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////////
NAMESPACE_BEGIN(Grid); NAMESPACE_BEGIN(Grid);
@ -59,10 +59,10 @@ NAMESPACE_BEGIN(Grid);
void Gather_plane_table_compute (GridBase *grid,int dimension,int plane,int cbmask, void Gather_plane_table_compute (GridBase *grid,int dimension,int plane,int cbmask,
int off,Vector<std::pair<int,int> > & table); int off,Vector<std::pair<int,int> > & table);
template<class vobj,class cobj,class compressor> template<class vobj,class cobj,class compressor>
void Gather_plane_simple_table (Vector<std::pair<int,int> >& table,const Lattice<vobj> &rhs,cobj *buffer,compressor &compress, int off,int so) __attribute__((noinline)); void Gather_plane_simple_table (Vector<std::pair<int,int> >& table,const Lattice<vobj> &rhs,cobj *buffer,compressor &compress, int off,int so) __attribute__((noinline));
template<class vobj,class cobj,class compressor> template<class vobj,class cobj,class compressor>
void Gather_plane_simple_table (Vector<std::pair<int,int> >& table,const Lattice<vobj> &rhs,cobj *buffer,compressor &compress, int off,int so) void Gather_plane_simple_table (Vector<std::pair<int,int> >& table,const Lattice<vobj> &rhs,cobj *buffer,compressor &compress, int off,int so)
{ {
int num=table.size(); int num=table.size();
@ -92,13 +92,13 @@ void Gather_plane_exchange_table(Vector<std::pair<int,int> >& table,const Lattic
{ {
assert( (table.size()&0x1)==0); assert( (table.size()&0x1)==0);
int num=table.size()/2; int num=table.size()/2;
int so = plane*rhs.Grid()->_ostride[dimension]; // base offset for start of plane int so = plane*rhs.Grid()->_ostride[dimension]; // base offset for start of plane
auto rhs_v = rhs.View(); auto rhs_v = rhs.View();
auto p0=&pointers[0][0]; auto p0=&pointers[0][0];
auto p1=&pointers[1][0]; auto p1=&pointers[1][0];
auto tp=&table[0]; auto tp=&table[0];
accelerator_forNB(j, num, 1, { accelerator_forNB(j, num, 1, {
compress.CompressExchange(p0,p1, &rhs_v[0], j, compress.CompressExchange(p0,p1, &rhs_v[0], j,
so+tp[2*j ].second, so+tp[2*j ].second,
so+tp[2*j+1].second, so+tp[2*j+1].second,
@ -106,20 +106,20 @@ void Gather_plane_exchange_table(Vector<std::pair<int,int> >& table,const Lattic
}); });
} }
struct StencilEntry { struct StencilEntry {
#ifdef GRID_NVCC #ifdef GRID_NVCC
uint64_t _byte_offset; // 8 bytes uint64_t _byte_offset; // 8 bytes
uint32_t _offset; // 4 bytes uint32_t _offset; // 4 bytes
#else #else
uint64_t _byte_offset; // 8 bytes uint64_t _byte_offset; // 8 bytes
uint64_t _offset; // 8 bytes (8 ever required?) uint64_t _offset; // 8 bytes (8 ever required?)
#endif #endif
uint8_t _is_local; // 1 bytes uint8_t _is_local; // 1 bytes
uint8_t _permute; // 1 bytes uint8_t _permute; // 1 bytes
uint8_t _around_the_world; // 1 bytes uint8_t _around_the_world; // 1 bytes
uint8_t _pad; // 1 bytes uint8_t _pad; // 1 bytes
}; };
// Could pack to 8 + 4 + 4 = 128 bit and use // Could pack to 8 + 4 + 4 = 128 bit and use
template<class vobj,class cobj,class Parameters> template<class vobj,class cobj,class Parameters>
class CartesianStencilView { class CartesianStencilView {
@ -145,18 +145,18 @@ class CartesianStencilView {
accelerator_inline cobj *CommBuf(void) { return u_recv_buf_p; } accelerator_inline cobj *CommBuf(void) { return u_recv_buf_p; }
accelerator_inline int GetNodeLocal(int osite,int point) { accelerator_inline int GetNodeLocal(int osite,int point) {
return this->_entries_p[point+this->_npoints*osite]._is_local; return this->_entries_p[point+this->_npoints*osite]._is_local;
} }
accelerator_inline StencilEntry * GetEntry(int &ptype,int point,int osite) { accelerator_inline StencilEntry * GetEntry(int &ptype,int point,int osite) {
ptype = this->_permute_type[point]; return & this->_entries_p[point+this->_npoints*osite]; ptype = this->_permute_type[point]; return & this->_entries_p[point+this->_npoints*osite];
} }
accelerator_inline uint64_t GetInfo(int &ptype,int &local,int &perm,int point,int ent,uint64_t base) { accelerator_inline uint64_t GetInfo(int &ptype,int &local,int &perm,int point,int ent,uint64_t base) {
uint64_t cbase = (uint64_t)&u_recv_buf_p[0]; uint64_t cbase = (uint64_t)&u_recv_buf_p[0];
local = this->_entries_p[ent]._is_local; local = this->_entries_p[ent]._is_local;
perm = this->_entries_p[ent]._permute; perm = this->_entries_p[ent]._permute;
if (perm) ptype = this->_permute_type[point]; if (perm) ptype = this->_permute_type[point];
if (local) { if (local) {
return base + this->_entries_p[ent]._byte_offset; return base + this->_entries_p[ent]._byte_offset;
} else { } else {
@ -171,7 +171,7 @@ class CartesianStencilView {
else return cbase + this->_entries_p[ent]._byte_offset; else return cbase + this->_entries_p[ent]._byte_offset;
} }
accelerator_inline void iCoorFromIindex(Coordinate &coor,int lane) accelerator_inline void iCoorFromIindex(Coordinate &coor,int lane)
{ {
Lexicographic::CoorFromIndex(coor,lane,this->_simd_layout); Lexicographic::CoorFromIndex(coor,lane,this->_simd_layout);
} }
@ -211,12 +211,12 @@ public:
cobj * mpi_p; cobj * mpi_p;
Integer buffer_size; Integer buffer_size;
}; };
protected: protected:
GridBase * _grid; GridBase * _grid;
public: public:
GridBase *Grid(void) const { return _grid; } GridBase *Grid(void) const { return _grid; }
//////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////
@ -230,7 +230,7 @@ public:
View_type accessor(*( (View_type *) this)); View_type accessor(*( (View_type *) this));
return accessor; return accessor;
} }
int face_table_computed; int face_table_computed;
std::vector<Vector<std::pair<int,int> > > face_table ; std::vector<Vector<std::pair<int,int> > > face_table ;
Vector<int> surface_list; Vector<int> surface_list;
@ -280,7 +280,7 @@ public:
//////////////////////////////////////// ////////////////////////////////////////
// Stencil query // Stencil query
//////////////////////////////////////// ////////////////////////////////////////
inline int SameNode(int point) { inline int SameNode(int point) {
int dimension = this->_directions[point]; int dimension = this->_directions[point];
int displacement = this->_distances[point]; int displacement = this->_distances[point];
@ -304,7 +304,7 @@ public:
// FIXME this logic needs to be sorted for three link term // FIXME this logic needs to be sorted for three link term
// assert( (displacement==1) || (displacement==-1)); // assert( (displacement==1) || (displacement==-1));
// Present hack only works for >= 4^4 subvol per node // Present hack only works for >= 4^4 subvol per node
_grid->ShiftedRanks(dimension,nbr_proc,xmit_to_rank,recv_from_rank); _grid->ShiftedRanks(dimension,nbr_proc,xmit_to_rank,recv_from_rank);
void *shm = (void *) _grid->ShmBufferTranslate(recv_from_rank,this->u_recv_buf_p); void *shm = (void *) _grid->ShmBufferTranslate(recv_from_rank,this->u_recv_buf_p);
@ -344,7 +344,7 @@ public:
comm_time_thr[mythread] += comm_leave_thr[mythread] - comm_enter_thr[mythread]; comm_time_thr[mythread] += comm_leave_thr[mythread] - comm_enter_thr[mythread];
} }
} }
void CollateThreads(void) void CollateThreads(void)
{ {
int nthreads = CartesianCommunicator::nCommThreads; int nthreads = CartesianCommunicator::nCommThreads;
@ -368,7 +368,7 @@ public:
if ( (t0 > 0.0) && ( t0 < first ) ) first = t0; // min time seen if ( (t0 > 0.0) && ( t0 < first ) ) first = t0; // min time seen
if ( t1 > last ) last = t1; // max time seen if ( t1 > last ) last = t1; // max time seen
} }
commtime+= last-first; commtime+= last-first;
} }
@ -430,30 +430,30 @@ public:
this->CommunicateBegin(reqs); this->CommunicateBegin(reqs);
this->CommunicateComplete(reqs); this->CommunicateComplete(reqs);
} }
} }
template<class compressor> void HaloExchange(const Lattice<vobj> &source,compressor &compress) template<class compressor> void HaloExchange(const Lattice<vobj> &source,compressor &compress)
{ {
Prepare(); Prepare();
HaloGather(source,compress); HaloGather(source,compress);
Communicate(); Communicate();
CommsMergeSHM(compress); CommsMergeSHM(compress);
CommsMerge(compress); CommsMerge(compress);
} }
template<class compressor> int HaloGatherDir(const Lattice<vobj> &source,compressor &compress,int point,int & face_idx) template<class compressor> int HaloGatherDir(const Lattice<vobj> &source,compressor &compress,int point,int & face_idx)
{ {
int dimension = this->_directions[point]; int dimension = this->_directions[point];
int displacement = this->_distances[point]; int displacement = this->_distances[point];
int fd = _grid->_fdimensions[dimension]; int fd = _grid->_fdimensions[dimension];
int rd = _grid->_rdimensions[dimension]; int rd = _grid->_rdimensions[dimension];
// Map to always positive shift modulo global full dimension. // Map to always positive shift modulo global full dimension.
int shift = (displacement+fd)%fd; int shift = (displacement+fd)%fd;
assert (source.Checkerboard()== this->_checkerboard); assert (source.Checkerboard()== this->_checkerboard);
// the permute type // the permute type
int simd_layout = _grid->_simd_layout[dimension]; int simd_layout = _grid->_simd_layout[dimension];
int comm_dim = _grid->_processors[dimension] >1 ; int comm_dim = _grid->_processors[dimension] >1 ;
@ -471,7 +471,7 @@ public:
auto tmp = GatherSimd(source,dimension,shift,0x3,compress,face_idx); auto tmp = GatherSimd(source,dimension,shift,0x3,compress,face_idx);
is_same_node = is_same_node && tmp; is_same_node = is_same_node && tmp;
splicetime+=usecond(); splicetime+=usecond();
} else { } else {
nosplicetime-=usecond(); nosplicetime-=usecond();
auto tmp = Gather(source,dimension,shift,0x3,compress,face_idx); auto tmp = Gather(source,dimension,shift,0x3,compress,face_idx);
is_same_node = is_same_node && tmp; is_same_node = is_same_node && tmp;
@ -497,7 +497,7 @@ public:
} }
return is_same_node; return is_same_node;
} }
template<class compressor> template<class compressor>
void HaloGather(const Lattice<vobj> &source,compressor &compress) void HaloGather(const Lattice<vobj> &source,compressor &compress)
{ {
@ -508,9 +508,9 @@ public:
// conformable(source.Grid(),_grid); // conformable(source.Grid(),_grid);
assert(source.Grid()==_grid); assert(source.Grid()==_grid);
halogtime-=usecond(); halogtime-=usecond();
u_comm_offset=0; u_comm_offset=0;
// Gather all comms buffers // Gather all comms buffers
int face_idx=0; int face_idx=0;
for(int point = 0 ; point < this->_npoints; point++) { for(int point = 0 ; point < this->_npoints; point++) {
@ -523,16 +523,16 @@ public:
accelerator_barrier(); accelerator_barrier();
halogtime+=usecond(); halogtime+=usecond();
} }
///////////////////////// /////////////////////////
// Implementation // Implementation
///////////////////////// /////////////////////////
void Prepare(void) void Prepare(void)
{ {
Decompressions.resize(0); Decompressions.resize(0);
DecompressionsSHM.resize(0); DecompressionsSHM.resize(0);
Mergers.resize(0); Mergers.resize(0);
MergersSHM.resize(0); MergersSHM.resize(0);
Packets.resize(0); Packets.resize(0);
calls++; calls++;
} }
@ -561,22 +561,22 @@ public:
mv.push_back(m); mv.push_back(m);
} }
template<class decompressor> void CommsMerge(decompressor decompress) { template<class decompressor> void CommsMerge(decompressor decompress) {
CommsMerge(decompress,Mergers,Decompressions); CommsMerge(decompress,Mergers,Decompressions);
} }
template<class decompressor> void CommsMergeSHM(decompressor decompress) { template<class decompressor> void CommsMergeSHM(decompressor decompress) {
mpi3synctime-=usecond(); mpi3synctime-=usecond();
_grid->StencilBarrier();// Synch shared memory on a single nodes _grid->StencilBarrier();// Synch shared memory on a single nodes
mpi3synctime+=usecond(); mpi3synctime+=usecond();
shmmergetime-=usecond(); shmmergetime-=usecond();
CommsMerge(decompress,MergersSHM,DecompressionsSHM); CommsMerge(decompress,MergersSHM,DecompressionsSHM);
shmmergetime+=usecond(); shmmergetime+=usecond();
} }
template<class decompressor> template<class decompressor>
void CommsMerge(decompressor decompress,std::vector<Merge> &mm,std::vector<Decompress> &dd) { void CommsMerge(decompressor decompress,std::vector<Merge> &mm,std::vector<Decompress> &dd) {
mergetime-=usecond(); mergetime-=usecond();
for(int i=0;i<mm.size();i++){ for(int i=0;i<mm.size();i++){
auto mp = &mm[i].mpointer[0]; auto mp = &mm[i].mpointer[0];
auto vp0= &mm[i].vpointers[0][0]; auto vp0= &mm[i].vpointers[0][0];
auto vp1= &mm[i].vpointers[1][0]; auto vp1= &mm[i].vpointers[1][0];
@ -588,7 +588,7 @@ public:
mergetime+=usecond(); mergetime+=usecond();
decompresstime-=usecond(); decompresstime-=usecond();
for(int i=0;i<dd.size();i++){ for(int i=0;i<dd.size();i++){
auto kp = dd[i].kernel_p; auto kp = dd[i].kernel_p;
auto mp = dd[i].mpi_p; auto mp = dd[i].mpi_p;
accelerator_forNB(o,dd[i].buffer_size,1,{ accelerator_forNB(o,dd[i].buffer_size,1,{
@ -604,7 +604,7 @@ public:
for(int i=0;i<_entries.size();i++){ for(int i=0;i<_entries.size();i++){
if( _entries[i]._is_local ) { if( _entries[i]._is_local ) {
_entries[i]._byte_offset = _entries[i]._offset*sizeof(vobj); _entries[i]._byte_offset = _entries[i]._offset*sizeof(vobj);
} else { } else {
_entries[i]._byte_offset = _entries[i]._offset*sizeof(cobj); _entries[i]._byte_offset = _entries[i]._offset*sizeof(cobj);
} }
} }
@ -619,15 +619,15 @@ public:
for(int point=0;point<this->_npoints;point++){ for(int point=0;point<this->_npoints;point++){
this->same_node[point] = this->SameNode(point); this->same_node[point] = this->SameNode(point);
} }
for(int site = 0 ;site< vol4;site++){ for(int site = 0 ;site< vol4;site++){
int local = 1; int local = 1;
for(int point=0;point<this->_npoints;point++){ for(int point=0;point<this->_npoints;point++){
if( (!this->GetNodeLocal(site*Ls,point)) && (!this->same_node[point]) ){ if( (!this->GetNodeLocal(site*Ls,point)) && (!this->same_node[point]) ){
local = 0; local = 0;
} }
} }
if(local == 0) { if(local == 0) {
surface_list.push_back(site); surface_list.push_back(site);
} }
} }
@ -638,11 +638,11 @@ public:
int checkerboard, int checkerboard,
const std::vector<int> &directions, const std::vector<int> &directions,
const std::vector<int> &distances, const std::vector<int> &distances,
Parameters p) Parameters p)
: shm_bytes_thr(npoints), : shm_bytes_thr(npoints),
comm_bytes_thr(npoints), comm_bytes_thr(npoints),
comm_enter_thr(npoints), comm_enter_thr(npoints),
comm_leave_thr(npoints), comm_leave_thr(npoints),
comm_time_thr(npoints) comm_time_thr(npoints)
{ {
face_table_computed=0; face_table_computed=0;
@ -653,7 +653,7 @@ public:
///////////////////////////////////// /////////////////////////////////////
this->_npoints = npoints; this->_npoints = npoints;
this->_comm_buf_size.resize(npoints), this->_comm_buf_size.resize(npoints),
this->_permute_type.resize(npoints), this->_permute_type.resize(npoints),
this->_simd_layout = _grid->_simd_layout; // copy simd_layout to give access to Accelerator Kernels this->_simd_layout = _grid->_simd_layout; // copy simd_layout to give access to Accelerator Kernels
this->_directions = StencilVector(directions); this->_directions = StencilVector(directions);
this->_distances = StencilVector(distances); this->_distances = StencilVector(distances);
@ -663,24 +663,24 @@ public:
surface_list.resize(0); surface_list.resize(0);
int osites = _grid->oSites(); int osites = _grid->oSites();
_entries.resize(this->_npoints* osites); _entries.resize(this->_npoints* osites);
this->_entries_p = &_entries[0]; this->_entries_p = &_entries[0];
for(int ii=0;ii<npoints;ii++){ for(int ii=0;ii<npoints;ii++){
int i = ii; // reverse direction to get SIMD comms done first int i = ii; // reverse direction to get SIMD comms done first
int point = i; int point = i;
int dimension = directions[i]; int dimension = directions[i];
int displacement = distances[i]; int displacement = distances[i];
int shift = displacement; int shift = displacement;
int fd = _grid->_fdimensions[dimension]; int fd = _grid->_fdimensions[dimension];
int rd = _grid->_rdimensions[dimension]; int rd = _grid->_rdimensions[dimension];
this->_permute_type[point]=_grid->PermuteType(dimension); this->_permute_type[point]=_grid->PermuteType(dimension);
this->_checkerboard = checkerboard; this->_checkerboard = checkerboard;
////////////////////////// //////////////////////////
// the permute type // the permute type
////////////////////////// //////////////////////////
@ -690,25 +690,25 @@ public:
int rotate_dim = _grid->_simd_layout[dimension]>2; int rotate_dim = _grid->_simd_layout[dimension]>2;
assert ( (rotate_dim && comm_dim) == false) ; // Do not think spread out is supported assert ( (rotate_dim && comm_dim) == false) ; // Do not think spread out is supported
int sshift[2]; int sshift[2];
////////////////////////// //////////////////////////
// Underlying approach. For each local site build // Underlying approach. For each local site build
// up a table containing the npoint "neighbours" and whether they // up a table containing the npoint "neighbours" and whether they
// live in lattice or a comms buffer. // live in lattice or a comms buffer.
////////////////////////// //////////////////////////
if ( !comm_dim ) { if ( !comm_dim ) {
sshift[0] = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,Even); sshift[0] = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,Even);
sshift[1] = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,Odd); sshift[1] = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,Odd);
if ( sshift[0] == sshift[1] ) { if ( sshift[0] == sshift[1] ) {
Local(point,dimension,shift,0x3); Local(point,dimension,shift,0x3);
} else { } else {
Local(point,dimension,shift,0x1);// if checkerboard is unfavourable take two passes Local(point,dimension,shift,0x1);// if checkerboard is unfavourable take two passes
Local(point,dimension,shift,0x2);// both with block stride loop iteration Local(point,dimension,shift,0x2);// both with block stride loop iteration
} }
} else { } else {
// All permute extract done in comms phase prior to Stencil application // All permute extract done in comms phase prior to Stencil application
// So tables are the same whether comm_dim or splice_dim // So tables are the same whether comm_dim or splice_dim
sshift[0] = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,Even); sshift[0] = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,Even);
@ -750,23 +750,23 @@ public:
int ld = _grid->_ldimensions[dimension]; int ld = _grid->_ldimensions[dimension];
int gd = _grid->_gdimensions[dimension]; int gd = _grid->_gdimensions[dimension];
int ly = _grid->_simd_layout[dimension]; int ly = _grid->_simd_layout[dimension];
// Map to always positive shift modulo global full dimension. // Map to always positive shift modulo global full dimension.
int shift = (shiftpm+fd)%fd; int shift = (shiftpm+fd)%fd;
// the permute type // the permute type
int permute_dim =_grid->PermuteDim(dimension); int permute_dim =_grid->PermuteDim(dimension);
for(int x=0;x<rd;x++){ for(int x=0;x<rd;x++){
// int o = 0; // int o = 0;
int bo = x * _grid->_ostride[dimension]; int bo = x * _grid->_ostride[dimension];
int cb= (cbmask==0x2)? Odd : Even; int cb= (cbmask==0x2)? Odd : Even;
int sshift = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,cb); int sshift = _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,cb);
int sx = (x+sshift)%rd; int sx = (x+sshift)%rd;
int wraparound=0; int wraparound=0;
if ( (shiftpm==-1) && (sx>x) ) { if ( (shiftpm==-1) && (sx>x) ) {
wraparound = 1; wraparound = 1;
@ -774,7 +774,7 @@ public:
if ( (shiftpm== 1) && (sx<x) ) { if ( (shiftpm== 1) && (sx<x) ) {
wraparound = 1; wraparound = 1;
} }
int permute_slice=0; int permute_slice=0;
if(permute_dim){ if(permute_dim){
int wrap = sshift/rd; wrap=wrap % ly; // but it is local anyway int wrap = sshift/rd; wrap=wrap % ly; // but it is local anyway
@ -782,66 +782,66 @@ public:
if ( x< rd-num ) permute_slice=wrap; if ( x< rd-num ) permute_slice=wrap;
else permute_slice = (wrap+1)%ly; else permute_slice = (wrap+1)%ly;
} }
CopyPlane(point,dimension,x,sx,cbmask,permute_slice,wraparound); CopyPlane(point,dimension,x,sx,cbmask,permute_slice,wraparound);
} }
} }
void Comms (int point,int dimension,int shiftpm,int cbmask) void Comms (int point,int dimension,int shiftpm,int cbmask)
{ {
GridBase *grid=_grid; GridBase *grid=_grid;
const int Nsimd = grid->Nsimd(); const int Nsimd = grid->Nsimd();
int fd = _grid->_fdimensions[dimension]; int fd = _grid->_fdimensions[dimension];
int ld = _grid->_ldimensions[dimension]; int ld = _grid->_ldimensions[dimension];
int rd = _grid->_rdimensions[dimension]; int rd = _grid->_rdimensions[dimension];
int pd = _grid->_processors[dimension]; int pd = _grid->_processors[dimension];
int simd_layout = _grid->_simd_layout[dimension]; int simd_layout = _grid->_simd_layout[dimension];
int comm_dim = _grid->_processors[dimension] >1 ; int comm_dim = _grid->_processors[dimension] >1 ;
assert(comm_dim==1); assert(comm_dim==1);
int shift = (shiftpm + fd) %fd; int shift = (shiftpm + fd) %fd;
assert(shift>=0); assert(shift>=0);
assert(shift<fd); assert(shift<fd);
// done in reduced dims, so SIMD factored // done in reduced dims, so SIMD factored
int buffer_size = _grid->_slice_nblock[dimension]*_grid->_slice_block[dimension]; int buffer_size = _grid->_slice_nblock[dimension]*_grid->_slice_block[dimension];
this->_comm_buf_size[point] = buffer_size; // Size of _one_ plane. Multiple planes may be gathered and this->_comm_buf_size[point] = buffer_size; // Size of _one_ plane. Multiple planes may be gathered and
// send to one or more remote nodes. // send to one or more remote nodes.
int cb= (cbmask==0x2)? Odd : Even; int cb= (cbmask==0x2)? Odd : Even;
int sshift= _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,cb); int sshift= _grid->CheckerBoardShiftForCB(this->_checkerboard,dimension,shift,cb);
for(int x=0;x<rd;x++){ for(int x=0;x<rd;x++){
int permute_type=grid->PermuteType(dimension); int permute_type=grid->PermuteType(dimension);
int sx = (x+sshift)%rd; int sx = (x+sshift)%rd;
int offnode = 0; int offnode = 0;
if ( simd_layout > 1 ) { if ( simd_layout > 1 ) {
for(int i=0;i<Nsimd;i++){ for(int i=0;i<Nsimd;i++){
int inner_bit = (Nsimd>>(permute_type+1)); int inner_bit = (Nsimd>>(permute_type+1));
int ic= (i&inner_bit)? 1:0; int ic= (i&inner_bit)? 1:0;
int my_coor = rd*ic + x; int my_coor = rd*ic + x;
int nbr_coor = my_coor+sshift; int nbr_coor = my_coor+sshift;
int nbr_proc = ((nbr_coor)/ld) % pd;// relative shift in processors int nbr_proc = ((nbr_coor)/ld) % pd;// relative shift in processors
if ( nbr_proc ) { if ( nbr_proc ) {
offnode =1; offnode =1;
} }
} }
} else { } else {
int comm_proc = ((x+sshift)/rd)%pd; int comm_proc = ((x+sshift)/rd)%pd;
offnode = (comm_proc!= 0); offnode = (comm_proc!= 0);
} }
int wraparound=0; int wraparound=0;
if ( (shiftpm==-1) && (sx>x) && (grid->_processor_coor[dimension]==0) ) { if ( (shiftpm==-1) && (sx>x) && (grid->_processor_coor[dimension]==0) ) {
wraparound = 1; wraparound = 1;
@ -850,24 +850,24 @@ public:
wraparound = 1; wraparound = 1;
} }
if (!offnode) { if (!offnode) {
int permute_slice=0; int permute_slice=0;
CopyPlane(point,dimension,x,sx,cbmask,permute_slice,wraparound); CopyPlane(point,dimension,x,sx,cbmask,permute_slice,wraparound);
} else { } else {
int words = buffer_size; int words = buffer_size;
if (cbmask != 0x3) words=words>>1; if (cbmask != 0x3) words=words>>1;
// int rank = grid->_processor; // int rank = grid->_processor;
// int recv_from_rank; // int recv_from_rank;
// int xmit_to_rank; // int xmit_to_rank;
int unified_buffer_offset = _unified_buffer_size; int unified_buffer_offset = _unified_buffer_size;
_unified_buffer_size += words; _unified_buffer_size += words;
ScatterPlane(point,dimension,x,cbmask,unified_buffer_offset,wraparound); // permute/extract/merge is done in comms phase ScatterPlane(point,dimension,x,cbmask,unified_buffer_offset,wraparound); // permute/extract/merge is done in comms phase
} }
} }
} }
@ -875,13 +875,13 @@ public:
void CopyPlane(int point, int dimension,int lplane,int rplane,int cbmask,int permute,int wrap) void CopyPlane(int point, int dimension,int lplane,int rplane,int cbmask,int permute,int wrap)
{ {
int rd = _grid->_rdimensions[dimension]; int rd = _grid->_rdimensions[dimension];
if ( !_grid->CheckerBoarded(dimension) ) { if ( !_grid->CheckerBoarded(dimension) ) {
int o = 0; // relative offset to base within plane int o = 0; // relative offset to base within plane
int ro = rplane*_grid->_ostride[dimension]; // base offset for start of plane int ro = rplane*_grid->_ostride[dimension]; // base offset for start of plane
int lo = lplane*_grid->_ostride[dimension]; // offset in buffer int lo = lplane*_grid->_ostride[dimension]; // offset in buffer
// Simple block stride gather of SIMD objects // Simple block stride gather of SIMD objects
for(int n=0;n<_grid->_slice_nblock[dimension];n++){ for(int n=0;n<_grid->_slice_nblock[dimension];n++){
for(int b=0;b<_grid->_slice_block[dimension];b++){ for(int b=0;b<_grid->_slice_block[dimension];b++){
@ -893,18 +893,18 @@ public:
} }
o +=_grid->_slice_stride[dimension]; o +=_grid->_slice_stride[dimension];
} }
} else { } else {
int ro = rplane*_grid->_ostride[dimension]; // base offset for start of plane int ro = rplane*_grid->_ostride[dimension]; // base offset for start of plane
int lo = lplane*_grid->_ostride[dimension]; // base offset for start of plane int lo = lplane*_grid->_ostride[dimension]; // base offset for start of plane
int o = 0; // relative offset to base within plane int o = 0; // relative offset to base within plane
for(int n=0;n<_grid->_slice_nblock[dimension];n++){ for(int n=0;n<_grid->_slice_nblock[dimension];n++){
for(int b=0;b<_grid->_slice_block[dimension];b++){ for(int b=0;b<_grid->_slice_block[dimension];b++){
int ocb=1<<_grid->CheckerBoardFromOindex(o+b); int ocb=1<<_grid->CheckerBoardFromOindex(o+b);
if ( ocb&cbmask ) { if ( ocb&cbmask ) {
int idx = point+(lo+o+b)*this->_npoints; int idx = point+(lo+o+b)*this->_npoints;
_entries[idx]._offset =ro+o+b; _entries[idx]._offset =ro+o+b;
@ -912,24 +912,24 @@ public:
_entries[idx]._permute=permute; _entries[idx]._permute=permute;
_entries[idx]._around_the_world=wrap; _entries[idx]._around_the_world=wrap;
} }
} }
o +=_grid->_slice_stride[dimension]; o +=_grid->_slice_stride[dimension];
} }
} }
} }
// Routine builds up integer table for each site in _offsets, _is_local, _permute // Routine builds up integer table for each site in _offsets, _is_local, _permute
void ScatterPlane (int point,int dimension,int plane,int cbmask,int offset, int wrap) void ScatterPlane (int point,int dimension,int plane,int cbmask,int offset, int wrap)
{ {
int rd = _grid->_rdimensions[dimension]; int rd = _grid->_rdimensions[dimension];
if ( !_grid->CheckerBoarded(dimension) ) { if ( !_grid->CheckerBoarded(dimension) ) {
int so = plane*_grid->_ostride[dimension]; // base offset for start of plane int so = plane*_grid->_ostride[dimension]; // base offset for start of plane
int o = 0; // relative offset to base within plane int o = 0; // relative offset to base within plane
int bo = 0; // offset in buffer int bo = 0; // offset in buffer
// Simple block stride gather of SIMD objects // Simple block stride gather of SIMD objects
for(int n=0;n<_grid->_slice_nblock[dimension];n++){ for(int n=0;n<_grid->_slice_nblock[dimension];n++){
for(int b=0;b<_grid->_slice_block[dimension];b++){ for(int b=0;b<_grid->_slice_block[dimension];b++){
@ -941,16 +941,16 @@ public:
} }
o +=_grid->_slice_stride[dimension]; o +=_grid->_slice_stride[dimension];
} }
} else { } else {
int so = plane*_grid->_ostride[dimension]; // base offset for start of plane int so = plane*_grid->_ostride[dimension]; // base offset for start of plane
int o = 0; // relative offset to base within plane int o = 0; // relative offset to base within plane
int bo = 0; // offset in buffer int bo = 0; // offset in buffer
for(int n=0;n<_grid->_slice_nblock[dimension];n++){ for(int n=0;n<_grid->_slice_nblock[dimension];n++){
for(int b=0;b<_grid->_slice_block[dimension];b++){ for(int b=0;b<_grid->_slice_block[dimension];b++){
int ocb=1<<_grid->CheckerBoardFromOindex(o+b);// Could easily be a table lookup int ocb=1<<_grid->CheckerBoardFromOindex(o+b);// Could easily be a table lookup
if ( ocb & cbmask ) { if ( ocb & cbmask ) {
int idx = point+(so+o+b)*this->_npoints; int idx = point+(so+o+b)*this->_npoints;
@ -964,16 +964,16 @@ public:
} }
} }
} }
template<class compressor> template<class compressor>
int Gather(const Lattice<vobj> &rhs,int dimension,int shift,int cbmask,compressor & compress,int &face_idx) int Gather(const Lattice<vobj> &rhs,int dimension,int shift,int cbmask,compressor & compress,int &face_idx)
{ {
typedef typename cobj::vector_type vector_type; typedef typename cobj::vector_type vector_type;
typedef typename cobj::scalar_type scalar_type; typedef typename cobj::scalar_type scalar_type;
assert(rhs.Grid()==_grid); assert(rhs.Grid()==_grid);
// conformable(_grid,rhs.Grid()); // conformable(_grid,rhs.Grid());
int fd = _grid->_fdimensions[dimension]; int fd = _grid->_fdimensions[dimension];
int rd = _grid->_rdimensions[dimension]; int rd = _grid->_rdimensions[dimension];
int pd = _grid->_processors[dimension]; int pd = _grid->_processors[dimension];
@ -985,37 +985,37 @@ public:
assert(shift<fd); assert(shift<fd);
int buffer_size = _grid->_slice_nblock[dimension]*_grid->_slice_block[dimension]; int buffer_size = _grid->_slice_nblock[dimension]*_grid->_slice_block[dimension];
int cb= (cbmask==0x2)? Odd : Even; int cb= (cbmask==0x2)? Odd : Even;
int sshift= _grid->CheckerBoardShiftForCB(rhs.Checkerboard(),dimension,shift,cb); int sshift= _grid->CheckerBoardShiftForCB(rhs.Checkerboard(),dimension,shift,cb);
int shm_receive_only = 1; int shm_receive_only = 1;
for(int x=0;x<rd;x++){ for(int x=0;x<rd;x++){
int sx = (x+sshift)%rd; int sx = (x+sshift)%rd;
int comm_proc = ((x+sshift)/rd)%pd; int comm_proc = ((x+sshift)/rd)%pd;
if (comm_proc) { if (comm_proc) {
int words = buffer_size; int words = buffer_size;
if (cbmask != 0x3) words=words>>1; if (cbmask != 0x3) words=words>>1;
int bytes = words * compress.CommDatumSize(); int bytes = words * compress.CommDatumSize();
int so = sx*rhs.Grid()->_ostride[dimension]; // base offset for start of plane int so = sx*rhs.Grid()->_ostride[dimension]; // base offset for start of plane
if ( !face_table_computed ) { if ( !face_table_computed ) {
face_table.resize(face_idx+1); face_table.resize(face_idx+1);
Gather_plane_table_compute ((GridBase *)_grid,dimension,sx,cbmask,u_comm_offset,face_table[face_idx]); Gather_plane_table_compute ((GridBase *)_grid,dimension,sx,cbmask,u_comm_offset,face_table[face_idx]);
} }
// int rank = _grid->_processor; // int rank = _grid->_processor;
int recv_from_rank; int recv_from_rank;
int xmit_to_rank; int xmit_to_rank;
_grid->ShiftedRanks(dimension,comm_proc,xmit_to_rank,recv_from_rank); _grid->ShiftedRanks(dimension,comm_proc,xmit_to_rank,recv_from_rank);
assert (xmit_to_rank != _grid->ThisRank()); assert (xmit_to_rank != _grid->ThisRank());
assert (recv_from_rank != _grid->ThisRank()); assert (recv_from_rank != _grid->ThisRank());
///////////////////////////////////////////////////////// /////////////////////////////////////////////////////////
// try the direct copy if possible // try the direct copy if possible
///////////////////////////////////////////////////////// /////////////////////////////////////////////////////////
@ -1028,13 +1028,13 @@ public:
} }
send_buf = (cobj *)_grid->ShmBufferTranslate(xmit_to_rank,recv_buf); send_buf = (cobj *)_grid->ShmBufferTranslate(xmit_to_rank,recv_buf);
if ( send_buf==NULL ) { if ( send_buf==NULL ) {
send_buf = this->u_send_buf_p; send_buf = this->u_send_buf_p;
} }
// Find out if we get the direct copy. // Find out if we get the direct copy.
void *success = (void *) _grid->ShmBufferTranslate(recv_from_rank,this->u_send_buf_p); void *success = (void *) _grid->ShmBufferTranslate(recv_from_rank,this->u_send_buf_p);
if (success==NULL) { if (success==NULL) {
// we found a packet that comes from MPI and contributes to this leg of stencil // we found a packet that comes from MPI and contributes to this leg of stencil
shm_receive_only = 0; shm_receive_only = 0;
} }
@ -1043,9 +1043,9 @@ public:
assert(send_buf!=NULL); assert(send_buf!=NULL);
Gather_plane_simple_table(face_table[face_idx],rhs,send_buf,compress,u_comm_offset,so); face_idx++; Gather_plane_simple_table(face_table[face_idx],rhs,send_buf,compress,u_comm_offset,so); face_idx++;
gathertime+=usecond(); gathertime+=usecond();
if ( compress.DecompressionStep() ) { if ( compress.DecompressionStep() ) {
if ( shm_receive_only ) { // Early decompress before MPI is finished is possible if ( shm_receive_only ) { // Early decompress before MPI is finished is possible
AddDecompress(&this->u_recv_buf_p[u_comm_offset], AddDecompress(&this->u_recv_buf_p[u_comm_offset],
&recv_buf[u_comm_offset], &recv_buf[u_comm_offset],
@ -1074,7 +1074,7 @@ public:
} }
return shm_receive_only; return shm_receive_only;
} }
template<class compressor> template<class compressor>
int GatherSimd(const Lattice<vobj> &rhs,int dimension,int shift,int cbmask,compressor &compress,int & face_idx) int GatherSimd(const Lattice<vobj> &rhs,int dimension,int shift,int cbmask,compressor &compress,int & face_idx)
{ {
@ -1102,7 +1102,7 @@ public:
/////////////////////////////////////////////// ///////////////////////////////////////////////
int buffer_size = _grid->_slice_nblock[dimension]*_grid->_slice_block[dimension]; int buffer_size = _grid->_slice_nblock[dimension]*_grid->_slice_block[dimension];
// int words = sizeof(cobj)/sizeof(vector_type); // int words = sizeof(cobj)/sizeof(vector_type);
assert(cbmask==0x3); // Fixme think there is a latent bug if not true assert(cbmask==0x3); // Fixme think there is a latent bug if not true
// This assert will trap it if ever hit. Not hit normally so far // This assert will trap it if ever hit. Not hit normally so far
int reduced_buffer_size = buffer_size; int reduced_buffer_size = buffer_size;
@ -1118,22 +1118,22 @@ public:
/////////////////////////////////////////// ///////////////////////////////////////////
// Work out what to send where // Work out what to send where
/////////////////////////////////////////// ///////////////////////////////////////////
int cb = (cbmask==0x2)? Odd : Even; int cb = (cbmask==0x2)? Odd : Even;
int sshift= _grid->CheckerBoardShiftForCB(rhs.Checkerboard(),dimension,shift,cb); int sshift= _grid->CheckerBoardShiftForCB(rhs.Checkerboard(),dimension,shift,cb);
// loop over outer coord planes orthog to dim // loop over outer coord planes orthog to dim
int shm_receive_only = 1; int shm_receive_only = 1;
for(int x=0;x<rd;x++){ for(int x=0;x<rd;x++){
int any_offnode = ( ((x+sshift)%fd) >= rd ); int any_offnode = ( ((x+sshift)%fd) >= rd );
if ( any_offnode ) { if ( any_offnode ) {
for(int i=0;i<maxl;i++){ for(int i=0;i<maxl;i++){
spointers[i] = (cobj *) &u_simd_send_buf[i][u_comm_offset]; spointers[i] = (cobj *) &u_simd_send_buf[i][u_comm_offset];
} }
int sx = (x+sshift)%rd; int sx = (x+sshift)%rd;
if ( !face_table_computed ) { if ( !face_table_computed ) {
@ -1168,13 +1168,13 @@ public:
int recv_from_rank; int recv_from_rank;
int xmit_to_rank; int xmit_to_rank;
_grid->ShiftedRanks(dimension,nbr_proc,xmit_to_rank,recv_from_rank); _grid->ShiftedRanks(dimension,nbr_proc,xmit_to_rank,recv_from_rank);
// shm == receive pointer if offnode // shm == receive pointer if offnode
// shm == Translate[send pointer] if on node -- my view of his send pointer // shm == Translate[send pointer] if on node -- my view of his send pointer
cobj *shm = (cobj *) _grid->ShmBufferTranslate(recv_from_rank,sp); cobj *shm = (cobj *) _grid->ShmBufferTranslate(recv_from_rank,sp);
if (shm==NULL) { if (shm==NULL) {
shm = rp; shm = rp;
// we found a packet that comes from MPI and contributes to this shift. // we found a packet that comes from MPI and contributes to this shift.
// is_same_node is only used in the WilsonStencil, and gets set for this point in the stencil. // is_same_node is only used in the WilsonStencil, and gets set for this point in the stencil.
@ -1188,15 +1188,15 @@ public:
AddPacket((void *)sp,(void *)rp,xmit_to_rank,recv_from_rank,bytes); AddPacket((void *)sp,(void *)rp,xmit_to_rank,recv_from_rank,bytes);
} else { } else {
rpointers[i] = sp; rpointers[i] = sp;
} }
} }
if ( shm_receive_only ) { if ( shm_receive_only ) {
AddMerge(&this->u_recv_buf_p[u_comm_offset],rpointers,reduced_buffer_size,permute_type,MergersSHM); AddMerge(&this->u_recv_buf_p[u_comm_offset],rpointers,reduced_buffer_size,permute_type,MergersSHM);
} else { } else {
AddMerge(&this->u_recv_buf_p[u_comm_offset],rpointers,reduced_buffer_size,permute_type,Mergers); AddMerge(&this->u_recv_buf_p[u_comm_offset],rpointers,reduced_buffer_size,permute_type,Mergers);
@ -1231,9 +1231,9 @@ public:
shm_bytes = 0.; shm_bytes = 0.;
calls = 0.; calls = 0.;
}; };
void Report(void) { void Report(void) {
#define AVERAGE(A) #define AVERAGE(A)
#define PRINTIT(A) AVERAGE(A); std::cout << GridLogMessage << " Stencil " << #A << " "<< A/calls<<std::endl; #define PRINTIT(A) AVERAGE(A); std::cout << GridLogMessage << " Stencil " << #A << " "<< A/calls<<std::endl;
RealD NP = _grid->_Nprocessors; RealD NP = _grid->_Nprocessors;
RealD NN = _grid->NodeCount(); RealD NN = _grid->NodeCount();
@ -1250,7 +1250,7 @@ public:
} }
} }
if (threaded) commtime += t; if (threaded) commtime += t;
_grid->GlobalSum(commtime); commtime/=NP; _grid->GlobalSum(commtime); commtime/=NP;
if ( calls > 0. ) { if ( calls > 0. ) {
std::cout << GridLogMessage << " Stencil calls "<<calls<<std::endl; std::cout << GridLogMessage << " Stencil calls "<<calls<<std::endl;
@ -1273,8 +1273,8 @@ public:
std::cout << GridLogMessage << " Stencil SHM " << (shm_bytes)/gatheralltime/1000.*NP/NN << " GB/s per node"<<std::endl; std::cout << GridLogMessage << " Stencil SHM " << (shm_bytes)/gatheralltime/1000.*NP/NN << " GB/s per node"<<std::endl;
auto all_bytes = comms_bytes+shm_bytes; auto all_bytes = comms_bytes+shm_bytes;
std::cout << GridLogMessage << " Stencil SHM all" << (all_bytes)/gatheralltime/1000. << " GB/s per rank"<<std::endl; std::cout << GridLogMessage << " Stencil SHM all " << (all_bytes)/gatheralltime/1000. << " GB/s per rank"<<std::endl;
std::cout << GridLogMessage << " Stencil SHM all" << (all_bytes)/gatheralltime/1000.*NP/NN << " GB/s per node"<<std::endl; std::cout << GridLogMessage << " Stencil SHM all " << (all_bytes)/gatheralltime/1000.*NP/NN << " GB/s per node"<<std::endl;
auto membytes = (shm_bytes + comms_bytes/2) // read/write auto membytes = (shm_bytes + comms_bytes/2) // read/write
+ (shm_bytes+comms_bytes)/2 * sizeof(vobj)/sizeof(cobj); + (shm_bytes+comms_bytes)/2 * sizeof(vobj)/sizeof(cobj);
@ -1292,7 +1292,7 @@ public:
#undef PRINTIT #undef PRINTIT
#undef AVERAGE #undef AVERAGE
}; };
}; };
NAMESPACE_END(Grid); NAMESPACE_END(Grid);