From 791cb050c8e79bd0994e474044da5f055a9c05e8 Mon Sep 17 00:00:00 2001 From: paboyle Date: Tue, 1 Nov 2016 11:35:43 +0000 Subject: [PATCH 01/14] Comms improvements --- benchmarks/Benchmark_comms.cc | 86 ++- bootstrap.sh | 6 - configure.ac | 4 + lib/Cshift.h | 4 + lib/Makefile.am | 5 + lib/communicator/Communicator_base.cc | 15 +- lib/communicator/Communicator_base.h | 43 +- lib/communicator/Communicator_mpi.cc | 12 +- lib/communicator/Communicator_mpi3.cc | 12 +- lib/communicator/Communicator_mpi3_leader.cc | 733 +++++++++++++++++++ lib/communicator/Communicator_none.cc | 8 +- lib/communicator/Communicator_shmem.cc | 31 +- scripts/arm_configure.experimental | 1 - scripts/arm_configure.experimental_cortex57 | 3 - scripts/bench_wilson.sh | 9 - scripts/build-all | 46 -- scripts/configure-all | 11 - scripts/configure-commands | 89 --- scripts/configure-cray | 10 - scripts/configure-mic | 10 - scripts/copyright | 16 +- scripts/cray-modules | 2 - scripts/reconfigure_script | 4 - scripts/update_fftw.sh | 18 - scripts/wilson.gnu | 7 - 25 files changed, 898 insertions(+), 287 deletions(-) create mode 100644 lib/communicator/Communicator_mpi3_leader.cc delete mode 100644 scripts/arm_configure.experimental delete mode 100644 scripts/arm_configure.experimental_cortex57 delete mode 100755 scripts/bench_wilson.sh delete mode 100755 scripts/build-all delete mode 100755 scripts/configure-all delete mode 100755 scripts/configure-commands delete mode 100755 scripts/configure-cray delete mode 100755 scripts/configure-mic delete mode 100644 scripts/cray-modules delete mode 100755 scripts/reconfigure_script delete mode 100755 scripts/update_fftw.sh delete mode 100644 scripts/wilson.gnu diff --git a/benchmarks/Benchmark_comms.cc b/benchmarks/Benchmark_comms.cc index d0cd96c7..234d2fbb 100644 --- a/benchmarks/Benchmark_comms.cc +++ b/benchmarks/Benchmark_comms.cc @@ -42,14 +42,13 @@ int main (int argc, char ** argv) int Nloop=10; int nmu=0; - for(int mu=0;mu<4;mu++) if (mpi_layout[mu]>1) nmu++; + for(int mu=0;mu1) nmu++; + std::cout< latt_size ({lat*mpi_layout[0], + lat*mpi_layout[1], + lat*mpi_layout[2], + lat*mpi_layout[3]}); + + GridCartesian Grid(latt_size,simd_layout,mpi_layout); + + std::vector xbuf(8); + std::vector rbuf(8); + Grid.ShmBufferFreeAll(); + for(int d=0;d<8;d++){ + xbuf[d] = (HalfSpinColourVectorD *)Grid.ShmBufferMalloc(lat*lat*lat*Ls*sizeof(HalfSpinColourVectorD)); + rbuf[d] = (HalfSpinColourVectorD *)Grid.ShmBufferMalloc(lat*lat*lat*Ls*sizeof(HalfSpinColourVectorD)); + } + + int ncomm; + int bytes=lat*lat*lat*Ls*sizeof(HalfSpinColourVectorD); + + double start=usecond(); + for(int i=0;i requests; + + ncomm=0; + for(int mu=0;mu<4;mu++){ + + if (mpi_layout[mu]>1 ) { + + ncomm++; + int comm_proc=1; + int xmit_to_rank; + int recv_from_rank; + + Grid.ShiftedRanks(mu,comm_proc,xmit_to_rank,recv_from_rank); + Grid.StencilSendToRecvFromBegin(requests, + (void *)&xbuf[mu][0], + xmit_to_rank, + (void *)&rbuf[mu][0], + recv_from_rank, + bytes); + + comm_proc = mpi_layout[mu]-1; + + Grid.ShiftedRanks(mu,comm_proc,xmit_to_rank,recv_from_rank); + Grid.StencilSendToRecvFromBegin(requests, + (void *)&xbuf[mu+4][0], + xmit_to_rank, + (void *)&rbuf[mu+4][0], + recv_from_rank, + bytes); + + } + } + Grid.StencilSendToRecvFromComplete(requests); + Grid.Barrier(); + + } + double stop=usecond(); + + double dbytes = bytes; + double xbytes = Nloop*dbytes*2.0*ncomm; + double rbytes = xbytes; + double bidibytes = xbytes+rbytes; + + double time = stop-start; // microseconds + + std::cout< #include #endif +#ifdef GRID_COMMS_MPI3L +#include +#endif + #ifdef GRID_COMMS_SHMEM #include // uses same implementation of communicator #endif diff --git a/lib/Makefile.am b/lib/Makefile.am index 3285a545..a779135f 100644 --- a/lib/Makefile.am +++ b/lib/Makefile.am @@ -9,6 +9,11 @@ if BUILD_COMMS_MPI3 extra_sources+=communicator/Communicator_base.cc endif +if BUILD_COMMS_MPI3L + extra_sources+=communicator/Communicator_mpi3_leader.cc + extra_sources+=communicator/Communicator_base.cc +endif + if BUILD_COMMS_SHMEM extra_sources+=communicator/Communicator_shmem.cc extra_sources+=communicator/Communicator_base.cc diff --git a/lib/communicator/Communicator_base.cc b/lib/communicator/Communicator_base.cc index 9810987d..77b4f800 100644 --- a/lib/communicator/Communicator_base.cc +++ b/lib/communicator/Communicator_base.cc @@ -31,13 +31,6 @@ namespace Grid { /////////////////////////////////////////////////////////////// // Info that is setup once and indept of cartesian layout /////////////////////////////////////////////////////////////// -int CartesianCommunicator::ShmRank; -int CartesianCommunicator::ShmSize; -int CartesianCommunicator::GroupRank; -int CartesianCommunicator::GroupSize; -int CartesianCommunicator::WorldRank; -int CartesianCommunicator::WorldSize; -int CartesianCommunicator::Slave; void * CartesianCommunicator::ShmCommBuf; ///////////////////////////////// @@ -70,12 +63,6 @@ int CartesianCommunicator::ProcessorCount(void) { return //////////////////////////////////////////////////////////////////////////////// // very VERY rarely (Log, serial RNG) we need world without a grid //////////////////////////////////////////////////////////////////////////////// -int CartesianCommunicator::RankWorld(void){ return WorldRank; }; -int CartesianCommunicator::Ranks (void) { return WorldSize; }; -int CartesianCommunicator::Nodes (void) { return GroupSize; }; -int CartesianCommunicator::Cores (void) { return ShmSize; }; -int CartesianCommunicator::NodeRank (void) { return GroupRank; }; -int CartesianCommunicator::CoreRank (void) { return ShmRank; }; void CartesianCommunicator::GlobalSum(ComplexF &c) { @@ -94,7 +81,7 @@ void CartesianCommunicator::GlobalSumVector(ComplexD *c,int N) GlobalSumVector((double *)c,2*N); } -#ifndef GRID_COMMS_MPI3 +#if !defined( GRID_COMMS_MPI3) && !defined (GRID_COMMS_MPI3L) void CartesianCommunicator::StencilSendToRecvFromBegin(std::vector &list, void *xmit, diff --git a/lib/communicator/Communicator_base.h b/lib/communicator/Communicator_base.h index 9a15d2b0..d4bb50d9 100644 --- a/lib/communicator/Communicator_base.h +++ b/lib/communicator/Communicator_base.h @@ -1,3 +1,4 @@ + /************************************************************************************* Grid physics library, www.github.com/paboyle/Grid @@ -37,6 +38,9 @@ Author: Peter Boyle #ifdef GRID_COMMS_MPI3 #include #endif +#ifdef GRID_COMMS_MPI3L +#include +#endif #ifdef GRID_COMMS_SHMEM #include #endif @@ -60,9 +64,9 @@ class CartesianCommunicator { std::vector _processor_coor; // linear processor coordinate unsigned long _ndimension; -#if defined (GRID_COMMS_MPI) || defined (GRID_COMMS_MPI3) - MPI_Comm communicator; +#if defined (GRID_COMMS_MPI) || defined (GRID_COMMS_MPI3) || defined (GRID_COMMS_MPI3L) static MPI_Comm communicator_world; + MPI_Comm communicator; typedef MPI_Request CommsRequest_t; #else typedef int CommsRequest_t; @@ -75,7 +79,15 @@ class CartesianCommunicator { // cartesian communicator on a subset of ranks, slave ranks controlled // by group leader with data xfer via shared memory //////////////////////////////////////////////////////////////////// -#ifdef GRID_COMMS_MPI3 +#ifdef GRID_COMMS_MPI3 + + static int ShmRank; + static int ShmSize; + static int GroupRank; + static int GroupSize; + static int WorldRank; + static int WorldSize; + std::vector WorldDims; std::vector GroupDims; std::vector ShmDims; @@ -83,7 +95,7 @@ class CartesianCommunicator { std::vector GroupCoor; std::vector ShmCoor; std::vector WorldCoor; - + static std::vector GroupRanks; static std::vector MyGroup; static int ShmSetup; @@ -93,13 +105,20 @@ class CartesianCommunicator { std::vector LexicographicToWorldRank; static std::vector ShmCommBufs; + #else static void ShmInitGeneric(void); static commVector ShmBufStorageVector; #endif + + ///////////////////////////////// + // Grid information and queries + // Implemented in Communicator_base.C + ///////////////////////////////// static void * ShmCommBuf; size_t heap_top; size_t heap_bytes; + void *ShmBufferSelf(void); void *ShmBuffer(int rank); void *ShmBufferTranslate(int rank,void * local_p); @@ -123,28 +142,12 @@ class CartesianCommunicator { int RankFromProcessorCoor(std::vector &coor); void ProcessorCoorFromRank(int rank,std::vector &coor); - ///////////////////////////////// - // Grid information and queries - ///////////////////////////////// - static int ShmRank; - static int ShmSize; - static int GroupSize; - static int GroupRank; - static int WorldRank; - static int WorldSize; - static int Slave; - int IsBoss(void) ; int BossRank(void) ; int ThisRank(void) ; const std::vector & ThisProcessorCoor(void) ; const std::vector & ProcessorGrid(void) ; int ProcessorCount(void) ; - static int Ranks (void); - static int Nodes (void); - static int Cores (void); - static int NodeRank (void); - static int CoreRank (void); //////////////////////////////////////////////////////////////////////////////// // very VERY rarely (Log, serial RNG) we need world without a grid diff --git a/lib/communicator/Communicator_mpi.cc b/lib/communicator/Communicator_mpi.cc index a638eebb..65ced9c7 100644 --- a/lib/communicator/Communicator_mpi.cc +++ b/lib/communicator/Communicator_mpi.cc @@ -44,13 +44,6 @@ void CartesianCommunicator::Init(int *argc, char ***argv) { MPI_Init(argc,argv); } MPI_Comm_dup (MPI_COMM_WORLD,&communicator_world); - MPI_Comm_rank(communicator_world,&WorldRank); - MPI_Comm_size(communicator_world,&WorldSize); - ShmRank=0; - ShmSize=1; - GroupRank=WorldRank; - GroupSize=WorldSize; - Slave =0; ShmInitGeneric(); } @@ -198,6 +191,11 @@ void CartesianCommunicator::Broadcast(int root,void* data, int bytes) // Should only be used prior to Grid Init finished. // Check for this? /////////////////////////////////////////////////////// +int CartesianCommunicator::RankWorld(void){ + int r; + MPI_Comm_rank(communicator_world,&r); + return r; +} void CartesianCommunicator::BroadcastWorld(int root,void* data, int bytes) { int ierr= MPI_Bcast(data, diff --git a/lib/communicator/Communicator_mpi3.cc b/lib/communicator/Communicator_mpi3.cc index 14a152ab..c707ec1f 100644 --- a/lib/communicator/Communicator_mpi3.cc +++ b/lib/communicator/Communicator_mpi3.cc @@ -30,12 +30,18 @@ Author: Peter Boyle namespace Grid { - /////////////////////////////////////////////////////////////////////////////////////////////////// // Info that is setup once and indept of cartesian layout /////////////////////////////////////////////////////////////////////////////////////////////////// int CartesianCommunicator::ShmSetup = 0; +int CartesianCommunicator::ShmRank; +int CartesianCommunicator::ShmSize; +int CartesianCommunicator::GroupRank; +int CartesianCommunicator::GroupSize; +int CartesianCommunicator::WorldRank; +int CartesianCommunicator::WorldSize; + MPI_Comm CartesianCommunicator::communicator_world; MPI_Comm CartesianCommunicator::ShmComm; MPI_Win CartesianCommunicator::ShmWindow; @@ -97,15 +103,15 @@ void CartesianCommunicator::Init(int *argc, char ***argv) { std::vector world_ranks(WorldSize); GroupRanks.resize(WorldSize); - MyGroup.resize(ShmSize); for(int r=0;r + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License along + with this program; if not, write to the Free Software Foundation, Inc., + 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + + See the full license in the file "LICENSE" in the top level distribution directory + *************************************************************************************/ + /* END LEGAL */ +#include "Grid.h" +#include + +namespace Grid { + +enum { COMMAND_ISEND, COMMAND_IRECV, COMMAND_WAITALL }; + +struct Descriptor { + uint64_t buf; + size_t bytes; + int rank; + int tag; + int command; + MPI_Request request; +}; + +const int pool = 48; + +class SlaveState { +public: + volatile int head; + volatile int start; + volatile int tail; + volatile Descriptor Descrs[pool]; +}; + +class Slave { +public: + SlaveState *state; + MPI_Comm squadron; + uint64_t base; + //////////////////////////////////////////////////////////// + // Descriptor circular pointers + //////////////////////////////////////////////////////////// + Slave() {}; + + void Init(SlaveState * _state,MPI_Comm _squadron); + + void EventLoop (void) { + std::cerr<< " Entering even loop "<tail != state->head ); + std::cerr<< " FIFO drained "<< state->tail < Slaves; + + static int ShmSetup; + + static int UniverseRank; + static int UniverseSize; + + static MPI_Comm communicator_universe; + static MPI_Comm communicator_cached; + + static MPI_Comm HorizontalComm; + static int HorizontalRank; + static int HorizontalSize; + + static MPI_Comm VerticalComm; + static MPI_Win VerticalWindow; + static int VerticalSize; + static int VerticalRank; + + static std::vector VerticalShmBufs; + static std::vector > UniverseRanks; + static std::vector UserCommunicatorToWorldRanks; + + static MPI_Group WorldGroup, CachedGroup; + + static void CommunicatorInit (MPI_Comm &communicator_world, + MPI_Comm &ShmComm, + void * &ShmCommBuf); + + static void MapCommRankToWorldRank(int &hashtag, int & comm_world_peer,int tag, MPI_Comm comm,int rank); + + ///////////////////////////////////////////////////////// + // routines for master proc must handle any communicator + ///////////////////////////////////////////////////////// + + static uint64_t QueueSend(int slave,void *buf, int bytes, int tag, MPI_Comm comm,int rank) { + std::cerr<< " Queueing send "<< bytes<= units ) { + mywork = myoff = 0; + } else { + mywork = (nwork+me)/units; + myoff = basework * me; + if ( me > backfill ) + myoff+= (me-backfill); + } + return; + }; + + static void QueueMultiplexedSend(void *buf, int bytes, int tag, MPI_Comm comm,int rank) { + uint8_t * cbuf = (uint8_t *) buf; + int mywork, myoff, procs; + procs = VerticalSize-1; + for(int s=0;s MPIoffloadEngine::Slaves; + +int MPIoffloadEngine::UniverseRank; +int MPIoffloadEngine::UniverseSize; + +MPI_Comm MPIoffloadEngine::communicator_universe; +MPI_Comm MPIoffloadEngine::communicator_cached; +MPI_Group MPIoffloadEngine::WorldGroup; +MPI_Group MPIoffloadEngine::CachedGroup; + +MPI_Comm MPIoffloadEngine::HorizontalComm; +int MPIoffloadEngine::HorizontalRank; +int MPIoffloadEngine::HorizontalSize; + +MPI_Comm MPIoffloadEngine::VerticalComm; +MPI_Win MPIoffloadEngine::VerticalWindow; +int MPIoffloadEngine::VerticalSize; +int MPIoffloadEngine::VerticalRank; + +std::vector MPIoffloadEngine::VerticalShmBufs; +std::vector > MPIoffloadEngine::UniverseRanks; +std::vector MPIoffloadEngine::UserCommunicatorToWorldRanks; + +int MPIoffloadEngine::ShmSetup = 0; + +void MPIoffloadEngine::CommunicatorInit (MPI_Comm &communicator_world, + MPI_Comm &ShmComm, + void * &ShmCommBuf) +{ + int flag; + assert(ShmSetup==0); + + ////////////////////////////////////////////////////////////////////// + // Universe is all nodes prior to squadron grouping + ////////////////////////////////////////////////////////////////////// + MPI_Comm_dup (MPI_COMM_WORLD,&communicator_universe); + MPI_Comm_rank(communicator_universe,&UniverseRank); + MPI_Comm_size(communicator_universe,&UniverseSize); + + ///////////////////////////////////////////////////////////////////// + // Split into groups that can share memory (Verticals) + ///////////////////////////////////////////////////////////////////// + // MPI_Comm_split_type(communicator_universe, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL,&VerticalComm); + MPI_Comm_split(communicator_universe,(UniverseRank&0x1),UniverseRank,&VerticalComm); + MPI_Comm_rank(VerticalComm ,&VerticalRank); + MPI_Comm_size(VerticalComm ,&VerticalSize); + + ////////////////////////////////////////////////////////////////////// + // Split into horizontal groups by rank in squadron + ////////////////////////////////////////////////////////////////////// + MPI_Comm_split(communicator_universe,VerticalRank,UniverseRank,&HorizontalComm); + MPI_Comm_rank(HorizontalComm,&HorizontalRank); + MPI_Comm_size(HorizontalComm,&HorizontalSize); + assert(HorizontalSize*VerticalSize==UniverseSize); + + //////////////////////////////////////////////////////////////////////////////// + // What is my place in the world + //////////////////////////////////////////////////////////////////////////////// + int WorldRank=0; + if(VerticalRank==0) WorldRank = HorizontalRank; + int ierr=MPI_Allreduce(MPI_IN_PLACE,&WorldRank,1,MPI_INT,MPI_SUM,VerticalComm); + assert(ierr==0); + + //////////////////////////////////////////////////////////////////////////////// + // Where is the world in the universe? + //////////////////////////////////////////////////////////////////////////////// + UniverseRanks = std::vector >(HorizontalSize,std::vector(VerticalSize,0)); + UniverseRanks[WorldRank][VerticalRank] = UniverseRank; + for(int w=0;w cached_ranks(size); + + for(int r=0;r>0 )&0xFFFF)^((icomm>>16)&0xFFFF) + ^ ((icomm>>32)&0xFFFF)^((icomm>>48)&0xFFFF); + + hashtag = (comm_hash<<15) | tag; + +}; + +void Slave::Init(SlaveState * _state,MPI_Comm _squadron) +{ + squadron=_squadron; + state =_state; + state->head = state->tail = state->start = 0; + MPI_Barrier(squadron); + base = (uint64_t)MPIoffloadEngine::VerticalShmBufs[0]; + int rank; MPI_Comm_rank(_squadron,&rank); +} +#define PERI_PLUS(A) ( (A+1)%pool ) +void Slave::Event (void) { + + static int tail_last; + static int head_last; + static int start_last; + int ierr; + + if ( (state->tail != tail_last) + ||(state->head != head_last) + ||(state->start != start_last) + ) { + std::cerr<< " Event loop "<< state->tail << " "<< state->start<< " "<< state->head <tail; + if ( t != state->start ) { + int flag=0; + + std::cerr<< " Testing tail "<< t<<" "<< (void *)&state->Descrs[t].request + << " "<Descrs[t].request<Descrs[t].request,&flag,MPI_STATUS_IGNORE); + // ierr=MPI_Test((MPI_Request *)&state->Descrs[t].request,&flag,MPI_STATUS_IGNORE); + assert(ierr==0); + if ( flag ) { + state->tail = PERI_PLUS(t); + std::cerr<< " Tail advanced from "<< t<start; + if ( s != state->head ) { + switch ( state->Descrs[s].command ) { + case COMMAND_ISEND: + std::cerr<< " Send "<Descrs[s].buf<< "["<Descrs[s].bytes<<"]" + << " to " << state->Descrs[s].rank<< " tag" << state->Descrs[s].tag + << " Comm " << MPIoffloadEngine::communicator_universe<< std::endl; + + std::cerr<< " Request was "<Descrs[s].request<Descrs[s].buf+base), + state->Descrs[s].bytes, + MPI_CHAR, + state->Descrs[s].rank, + state->Descrs[s].tag, + MPIoffloadEngine::communicator_universe, + (MPI_Request *)&state->Descrs[s].request); + std::cerr<< " Request is "<Descrs[s].request<Descrs[0].request<start = PERI_PLUS(s); + break; + + case COMMAND_IRECV: + std::cerr<< " Recv "<Descrs[s].buf<< "["<Descrs[s].bytes<<"]" + << " to " << state->Descrs[s].rank<< " tag" << state->Descrs[s].tag + << " Comm " << MPIoffloadEngine::communicator_universe<< std::endl; + + std::cerr<< " Request was "<Descrs[s].request<Descrs[s].buf+base), + state->Descrs[s].bytes, + MPI_CHAR, + state->Descrs[s].rank, + state->Descrs[s].tag, + MPIoffloadEngine::communicator_universe, + (MPI_Request *)&state->Descrs[s].request); + std::cerr<< " Request is "<Descrs[s].request<Descrs[0].request<start = PERI_PLUS(s); + break; + + case COMMAND_WAITALL: + std::cerr<< " Wait all "<tail;t!=s; t=PERI_PLUS(t) ){ + std::cerr<< " Wait ["<Descrs[t].request <Descrs[0].request<Descrs[t].request,MPI_STATUS_IGNORE); + }; + s=PERI_PLUS(s); + state->start = s; + state->tail = s; + break; + + default: + assert(0); + break; + } + return; + } +} + ////////////////////////////////////////////////////////////////////////////// + // External interaction with the queue + ////////////////////////////////////////////////////////////////////////////// + +uint64_t Slave::QueueCommand(int command,void *buf, int bytes, int tag, MPI_Comm comm,int commrank) +{ + ///////////////////////////////////////// + // Spin; if FIFO is full until not full + ///////////////////////////////////////// + int head =state->head; + int next = PERI_PLUS(head); + + // Set up descriptor + int worldrank; + int hashtag; + MPI_Comm communicator; + MPI_Request request; + + MPIoffloadEngine::MapCommRankToWorldRank(hashtag,commrank,tag,comm,worldrank); + + int VerticalRank = MPIoffloadEngine::VerticalRank; + uint64_t relative= (uint64_t)buf - base; + state->Descrs[head].buf = relative; + state->Descrs[head].bytes = bytes; + state->Descrs[head].rank = MPIoffloadEngine::UniverseRanks[worldrank][VerticalRank]; + state->Descrs[head].tag = hashtag; + state->Descrs[head].command= command; + std::cerr<< " QueueCommand "<tail==next ); + + // Msync on weak order architectures + // Advance pointer + state->head = next; + + return 0; +} + + +/////////////////////////////////////////////////////////////////////////////////////////////////// +// Info that is setup once and indept of cartesian layout +/////////////////////////////////////////////////////////////////////////////////////////////////// + +MPI_Comm CartesianCommunicator::communicator_world; + +void CartesianCommunicator::Init(int *argc, char ***argv) +{ + int flag; + MPI_Initialized(&flag); // needed to coexist with other libs apparently + if ( !flag ) { + MPI_Init(argc,argv); + } + communicator_world = MPI_COMM_WORLD; + MPI_Comm ShmComm; + MPIoffloadEngine::CommunicatorInit (communicator_world,ShmComm,ShmCommBuf); +} +void CartesianCommunicator::ShiftedRanks(int dim,int shift,int &source,int &dest) +{ + int ierr=MPI_Cart_shift(communicator,dim,shift,&source,&dest); + assert(ierr==0); +} +int CartesianCommunicator::RankFromProcessorCoor(std::vector &coor) +{ + int rank; + int ierr=MPI_Cart_rank (communicator, &coor[0], &rank); + assert(ierr==0); + return rank; +} +void CartesianCommunicator::ProcessorCoorFromRank(int rank, std::vector &coor) +{ + coor.resize(_ndimension); + int ierr=MPI_Cart_coords (communicator, rank, _ndimension,&coor[0]); + assert(ierr==0); +} + +CartesianCommunicator::CartesianCommunicator(const std::vector &processors) +{ + _ndimension = processors.size(); + std::vector periodic(_ndimension,1); + + _Nprocessors=1; + _processors = processors; + + for(int i=0;i<_ndimension;i++){ + _Nprocessors*=_processors[i]; + } + + int Size; + MPI_Comm_size(communicator_world,&Size); + assert(Size==_Nprocessors); + + _processor_coor.resize(_ndimension); + MPI_Cart_create(communicator_world, _ndimension,&_processors[0],&periodic[0],1,&communicator); + MPI_Comm_rank (communicator,&_processor); + MPI_Cart_coords(communicator,_processor,_ndimension,&_processor_coor[0]); +}; + +void CartesianCommunicator::GlobalSum(uint32_t &u){ + int ierr=MPI_Allreduce(MPI_IN_PLACE,&u,1,MPI_UINT32_T,MPI_SUM,communicator); + assert(ierr==0); +} +void CartesianCommunicator::GlobalSum(uint64_t &u){ + int ierr=MPI_Allreduce(MPI_IN_PLACE,&u,1,MPI_UINT64_T,MPI_SUM,communicator); + assert(ierr==0); +} +void CartesianCommunicator::GlobalSum(float &f){ + int ierr=MPI_Allreduce(MPI_IN_PLACE,&f,1,MPI_FLOAT,MPI_SUM,communicator); + assert(ierr==0); +} +void CartesianCommunicator::GlobalSumVector(float *f,int N) +{ + int ierr=MPI_Allreduce(MPI_IN_PLACE,f,N,MPI_FLOAT,MPI_SUM,communicator); + assert(ierr==0); +} +void CartesianCommunicator::GlobalSum(double &d) +{ + int ierr = MPI_Allreduce(MPI_IN_PLACE,&d,1,MPI_DOUBLE,MPI_SUM,communicator); + assert(ierr==0); +} +void CartesianCommunicator::GlobalSumVector(double *d,int N) +{ + int ierr = MPI_Allreduce(MPI_IN_PLACE,d,N,MPI_DOUBLE,MPI_SUM,communicator); + assert(ierr==0); +} + +// Basic Halo comms primitive +void CartesianCommunicator::SendToRecvFrom(void *xmit, + int dest, + void *recv, + int from, + int bytes) +{ + std::vector reqs(0); + SendToRecvFromBegin(reqs,xmit,dest,recv,from,bytes); + SendToRecvFromComplete(reqs); +} + +void CartesianCommunicator::SendRecvPacket(void *xmit, + void *recv, + int sender, + int receiver, + int bytes) +{ + MPI_Status stat; + assert(sender != receiver); + int tag = sender; + if ( _processor == sender ) { + MPI_Send(xmit, bytes, MPI_CHAR,receiver,tag,communicator); + } + if ( _processor == receiver ) { + MPI_Recv(recv, bytes, MPI_CHAR,sender,tag,communicator,&stat); + } +} + +// Basic Halo comms primitive +void CartesianCommunicator::SendToRecvFromBegin(std::vector &list, + void *xmit, + int dest, + void *recv, + int from, + int bytes) +{ + MPI_Request xrq; + MPI_Request rrq; + int rank = _processor; + int ierr; + ierr =MPI_Isend(xmit, bytes, MPI_CHAR,dest,_processor,communicator,&xrq); + ierr|=MPI_Irecv(recv, bytes, MPI_CHAR,from,from,communicator,&rrq); + + assert(ierr==0); + + list.push_back(xrq); + list.push_back(rrq); +} + +void CartesianCommunicator::StencilSendToRecvFromBegin(std::vector &list, + void *xmit, + int dest, + void *recv, + int from, + int bytes) +{ + uint64_t xmit_i = (uint64_t) xmit; + uint64_t recv_i = (uint64_t) recv; + uint64_t shm = (uint64_t) ShmCommBuf; + // assert xmit and recv lie in shared memory region + assert( (xmit_i >= shm) && (xmit_i+bytes <= shm+MAX_MPI_SHM_BYTES) ); + assert( (recv_i >= shm) && (recv_i+bytes <= shm+MAX_MPI_SHM_BYTES) ); + MPIoffloadEngine::QueueMultiplexedSend(xmit,bytes,_processor,communicator,dest); + MPIoffloadEngine::QueueMultiplexedRecv(recv,bytes,from,communicator,from); +} + + +void CartesianCommunicator::StencilSendToRecvFromComplete(std::vector &list) +{ + MPIoffloadEngine::WaitAll(); +} + +void CartesianCommunicator::StencilBarrier(void) +{ +} + +void CartesianCommunicator::SendToRecvFromComplete(std::vector &list) +{ + int nreq=list.size(); + std::vector status(nreq); + int ierr = MPI_Waitall(nreq,&list[0],&status[0]); + assert(ierr==0); +} + +void CartesianCommunicator::Barrier(void) +{ + int ierr = MPI_Barrier(communicator); + assert(ierr==0); +} + +void CartesianCommunicator::Broadcast(int root,void* data, int bytes) +{ + int ierr=MPI_Bcast(data, + bytes, + MPI_BYTE, + root, + communicator); + assert(ierr==0); +} + +void CartesianCommunicator::BroadcastWorld(int root,void* data, int bytes) +{ + int ierr= MPI_Bcast(data, + bytes, + MPI_BYTE, + root, + communicator_world); + assert(ierr==0); +} + +void *CartesianCommunicator::ShmBufferSelf(void) { return ShmCommBuf; } + +void *CartesianCommunicator::ShmBuffer(int rank) { + return NULL; +} +void *CartesianCommunicator::ShmBufferTranslate(int rank,void * local_p) { + return NULL; +} + + +}; + diff --git a/lib/communicator/Communicator_none.cc b/lib/communicator/Communicator_none.cc index 198c1add..0f43f1f5 100644 --- a/lib/communicator/Communicator_none.cc +++ b/lib/communicator/Communicator_none.cc @@ -34,13 +34,6 @@ namespace Grid { void CartesianCommunicator::Init(int *argc, char *** arv) { - WorldRank = 0; - WorldSize = 1; - ShmRank=0; - ShmSize=1; - GroupRank=WorldRank; - GroupSize=WorldSize; - Slave =0; ShmInitGeneric(); } @@ -99,6 +92,7 @@ void CartesianCommunicator::SendToRecvFromComplete(std::vector & assert(0); } +int CartesianCommunicator::RankWorld(void){return 0;} void CartesianCommunicator::Barrier(void){} void CartesianCommunicator::Broadcast(int root,void* data, int bytes) {} void CartesianCommunicator::BroadcastWorld(int root,void* data, int bytes) { } diff --git a/lib/communicator/Communicator_shmem.cc b/lib/communicator/Communicator_shmem.cc index d5f3ed76..56e03224 100644 --- a/lib/communicator/Communicator_shmem.cc +++ b/lib/communicator/Communicator_shmem.cc @@ -50,11 +50,16 @@ typedef struct HandShake_t { uint64_t seq_remote; } HandShake; +std::array make_psync_init(void) { + array ret; + ret.fill(SHMEM_SYNC_VALUE); + return ret; +} +static std::array psync_init = make_psync_init(); static Vector< HandShake > XConnections; static Vector< HandShake > RConnections; - void CartesianCommunicator::Init(int *argc, char ***argv) { shmem_init(); XConnections.resize(shmem_n_pes()); @@ -65,13 +70,6 @@ void CartesianCommunicator::Init(int *argc, char ***argv) { RConnections[pe].seq_local = 0; RConnections[pe].seq_remote= 0; } - WorldSize = shmem_n_pes(); - WorldRank = shmem_my_pe(); - ShmRank=0; - ShmSize=1; - GroupRank=WorldRank; - GroupSize=WorldSize; - Slave =0; shmem_barrier_all(); ShmInitGeneric(); } @@ -103,7 +101,7 @@ void CartesianCommunicator::GlobalSum(uint32_t &u){ static long long source ; static long long dest ; static long long llwrk[_SHMEM_REDUCE_MIN_WRKDATA_SIZE]; - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; // int nreduce=1; // int pestart=0; @@ -119,7 +117,7 @@ void CartesianCommunicator::GlobalSum(uint64_t &u){ static long long source ; static long long dest ; static long long llwrk[_SHMEM_REDUCE_MIN_WRKDATA_SIZE]; - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; // int nreduce=1; // int pestart=0; @@ -135,7 +133,7 @@ void CartesianCommunicator::GlobalSum(float &f){ static float source ; static float dest ; static float llwrk[_SHMEM_REDUCE_MIN_WRKDATA_SIZE]; - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; source = f; dest =0.0; @@ -147,7 +145,7 @@ void CartesianCommunicator::GlobalSumVector(float *f,int N) static float source ; static float dest = 0 ; static float llwrk[_SHMEM_REDUCE_MIN_WRKDATA_SIZE]; - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; if ( shmem_addr_accessible(f,_processor) ){ shmem_float_sum_to_all(f,f,N,0,0,_Nprocessors,llwrk,psync); @@ -166,7 +164,7 @@ void CartesianCommunicator::GlobalSum(double &d) static double source; static double dest ; static double llwrk[_SHMEM_REDUCE_MIN_WRKDATA_SIZE]; - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; source = d; dest = 0; @@ -178,7 +176,8 @@ void CartesianCommunicator::GlobalSumVector(double *d,int N) static double source ; static double dest ; static double llwrk[_SHMEM_REDUCE_MIN_WRKDATA_SIZE]; - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; + if ( shmem_addr_accessible(d,_processor) ){ shmem_double_sum_to_all(d,d,N,0,0,_Nprocessors,llwrk,psync); @@ -295,7 +294,7 @@ void CartesianCommunicator::Barrier(void) } void CartesianCommunicator::Broadcast(int root,void* data, int bytes) { - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; static uint32_t word; uint32_t *array = (uint32_t *) data; assert( (bytes % 4)==0); @@ -318,7 +317,7 @@ void CartesianCommunicator::Broadcast(int root,void* data, int bytes) } void CartesianCommunicator::BroadcastWorld(int root,void* data, int bytes) { - static long psync[_SHMEM_REDUCE_SYNC_SIZE]; + static std::array psync = psync_init; static uint32_t word; uint32_t *array = (uint32_t *) data; assert( (bytes % 4)==0); diff --git a/scripts/arm_configure.experimental b/scripts/arm_configure.experimental deleted file mode 100644 index fbad84a5..00000000 --- a/scripts/arm_configure.experimental +++ /dev/null @@ -1 +0,0 @@ -./configure --host=arm-linux-gnueabihf CXX=clang++-3.5 CXXFLAGS='-std=c++11 -O3 -target arm-linux-gnueabihf -I/usr/arm-linux-gnueabihf/include/ -I/home/neo/Codes/gmp6.0/gmp-arm/include/ -I/usr/arm-linux-gnueabihf/include/c++/4.8.2/arm-linux-gnueabihf/ -L/home/neo/Codes/gmp6.0/gmp-arm/lib/ -I/home/neo/Codes/mpfr3.1.2/mpfr-arm/include/ -L/home/neo/Codes/mpfr3.1.2/mpfr-arm/lib/ -static -mcpu=cortex-a7' --enable-simd=NEONv7 diff --git a/scripts/arm_configure.experimental_cortex57 b/scripts/arm_configure.experimental_cortex57 deleted file mode 100644 index d229763e..00000000 --- a/scripts/arm_configure.experimental_cortex57 +++ /dev/null @@ -1,3 +0,0 @@ -#./configure --host=arm-linux-gnueabihf CXX=clang++-3.5 CXXFLAGS='-std=c++11 -O3 -target arm-linux-gnueabihf -I/usr/arm-linux-gnueabihf/include/ -I/home/neo/Codes/gmp6.0/gmp-arm/include/ -I/usr/lib/llvm-3.5/lib/clang/3.5.0/include/ -L/home/neo/Codes/gmp6.0/gmp-arm/lib/ -I/home/neo/Codes/mpfr3.1.2/mpfr-arm/include/ -L/home/neo/Codes/mpfr3.1.2/mpfr-arm/lib/ -static -mcpu=cortex-a57' --enable-simd=NEONv7 - -./configure --host=aarch64-linux-gnu CXX=clang++-3.5 CXXFLAGS='-std=c++11 -O3 -target aarch64-linux-gnu -static -I/home/neo/Codes/gmp6.0/gmp-armv8/include/ -L/home/neo/Codes/gmp6.0/gmp-armv8/lib/ -I/home/neo/Codes/mpfr3.1.2/mpfr-armv8/include/ -L/home/neo/Codes/mpfr3.1.2/mpfr-armv8/lib/ -I/usr/aarch64-linux-gnu/include/ -I/usr/aarch64-linux-gnu/include/c++/4.8.2/aarch64-linux-gnu/' --enable-simd=NEONv7 diff --git a/scripts/bench_wilson.sh b/scripts/bench_wilson.sh deleted file mode 100755 index af73d591..00000000 --- a/scripts/bench_wilson.sh +++ /dev/null @@ -1,9 +0,0 @@ -for omp in 1 2 4 -do -echo > wilson.t$omp -for vol in 4.4.4.4 4.4.4.8 4.4.8.8 4.8.8.8 8.8.8.8 8.8.8.16 8.8.16.16 8.16.16.16 -do -perf=` ./benchmarks/Grid_wilson --grid $vol --omp $omp | grep mflop | awk '{print $3}'` -echo $vol $perf >> wilson.t$omp -done -done \ No newline at end of file diff --git a/scripts/build-all b/scripts/build-all deleted file mode 100755 index b97dca19..00000000 --- a/scripts/build-all +++ /dev/null @@ -1,46 +0,0 @@ -#!/bin/bash -e - -DIRS="clang-avx clang-avx-openmp clang-avx-openmp-mpi clang-avx-mpi clang-avx2 clang-avx2-openmp clang-avx2-openmp-mpi clang-avx2-mpi clang-sse" -EXTRADIRS="g++-avx g++-sse4 icpc-avx icpc-avx2 icpc-avx512" -BLACK="\033[30m" -RED="\033[31m" -GREEN="\033[32m" -YELLOW="\033[33m" -BLUE="\033[34m" -PINK="\033[35m" -CYAN="\033[36m" -WHITE="\033[37m" -NORMAL="\033[0;39m" - -for D in $DIRS -do - -echo -echo -e $RED ============================== -echo -e $GREEN $D -echo -e $RED ============================== -echo -e $BLUE - - cd builds/$D - make clean all -j 8 - cd ../../ -echo -e $NORMAL -done - -if [ "X$1" == "Xextra" ] -then -for D in $EXTRADIRS -do - -echo -echo -e $RED ============================== -echo -e $RED $D -echo -e $RED ============================== -echo -e $BLUE - - cd builds/$D - make clean all -j 8 - cd ../../ -echo -e $NORMAL -done -fi \ No newline at end of file diff --git a/scripts/configure-all b/scripts/configure-all deleted file mode 100755 index ad91034d..00000000 --- a/scripts/configure-all +++ /dev/null @@ -1,11 +0,0 @@ -#!/bin/bash - -DIRS="clang-avx clang-avx-openmp clang-avx-openmp-mpi clang-avx-mpi clang-avx2 clang-avx2-openmp clang-avx2-openmp-mpi clang-avx2-mpi icpc-avx icpc-avx2 icpc-avx512 g++-sse4 g++-avx clang-sse icpc-avx-openmp-mpi icpc-avx-openmp" - -for D in $DIRS -do - mkdir -p builds/$D - cd builds/$D - ../../scripts/configure-commands $D - cd ../.. -done diff --git a/scripts/configure-commands b/scripts/configure-commands deleted file mode 100755 index a3599d1f..00000000 --- a/scripts/configure-commands +++ /dev/null @@ -1,89 +0,0 @@ -#!/bin/bash -WD=$1 -BLACK="\033[30m" -RED="\033[31m" -GREEN="\033[32m" -YELLOW="\033[33m" -BLUE="\033[34m" -PINK="\033[35m" -CYAN="\033[36m" -WHITE="\033[37m" -NORMAL="\033[0;39m" -echo -echo -e $RED ============================== -echo -e $GREEN $WD -echo -e $RED ============================== -echo -e $YELLOW - -case $WD in -g++-avx) - CXX=g++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -g++-avx-openmp) - CXX=g++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -fopenmp -O3 -std=c++11" LIBS="-fopenmp -lgmp -lmpfr" --enable-comms=none - ;; -g++5-sse4) - CXX=g++-5 ../../configure --enable-simd=SSE4 CXXFLAGS="-msse4 -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -g++5-avx) - CXX=g++-5 ../../configure --enable-simd=AVX CXXFLAGS="-mavx -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -icpc-avx) - CXX=icpc ../../configure --enable-simd=AVX CXXFLAGS="-mavx -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -icpc-avx-openmp-mpi) -CXX=icpc ../../configure --enable-simd=AVX CXXFLAGS="-mavx -fopenmp -O3 -I/opt/local/include/openmpi-mp/ -std=c++11" LDFLAGS=-L/opt/local/lib/openmpi-mp/ LIBS="-lmpi -lmpi_cxx -fopenmp -lgmp -lmpfr" --enable-comms=mpi -;; -icpc-avx-openmp) -CXX=icpc ../../configure --enable-precision=single --enable-simd=AVX CXXFLAGS="-mavx -fopenmp -O3 -std=c++11" LIBS="-fopenmp -lgmp -lmpfr" --enable-comms=mpi -;; -icpc-avx2) - CXX=icpc ../../configure --enable-simd=AVX2 CXXFLAGS="-march=core-avx2 -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -icpc-avx512) - CXX=icpc ../../configure --enable-simd=AVX512 CXXFLAGS="-xCOMMON-AVX512 -O3 -std=c++11" --host=none LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -icpc-mic) - CXX=icpc ../../configure --host=none --enable-simd=IMCI CXXFLAGS="-mmic -O3 -std=c++11" LDFLAGS=-mmic LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -icpc-mic-avx512) - CXX=icpc ../../configure --host=none --enable-simd=IMCI CXXFLAGS="-xCOMMON_AVX512 -O3 -std=c++11" LDFLAGS=-xCOMMON_AVX512 LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-sse) -CXX=clang++ ../../configure --enable-precision=single --enable-simd=SSE4 CXXFLAGS="-msse4 -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-avx) -CXX=clang++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-avx2) -CXX=clang++ ../../configure --enable-simd=AVX2 CXXFLAGS="-mavx2 -mfma -O3 -std=c++11" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-avx-openmp) -CXX=clang-omp++ ../../configure --enable-precision=double --enable-simd=AVX CXXFLAGS="-mavx -fopenmp -O3 -std=c++11" LDFLAGS="-fopenmp" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-xc30) -CXX=$HOME/Clang/install/bin/clang++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -O3 -std=c++11 -I/opt/gcc/4.9.2/snos/include/g++/x86_64-suse-linux/ -I/opt/gcc/4.9.2/snos/include/g++/ " LDFLAGS="" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-xc30-openmp) -CXX=$HOME/Clang/install/bin/clang++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -fopenmp -O3 -std=c++11 -I/opt/gcc/4.9.2/snos/include/g++/x86_64-suse-linux/ -I/opt/gcc/4.9.2/snos/include/g++/ " LDFLAGS="-fopenmp" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-avx2-openmp) -CXX=clang-omp++ ../../configure --enable-simd=AVX2 CXXFLAGS="-mavx2 -mfma -fopenmp -O3 -std=c++11" LDFLAGS="-fopenmp" LIBS="-lgmp -lmpfr" --enable-comms=none - ;; -clang-avx-openmp-mpi) -CXX=clang-omp++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -fopenmp -O3 -I/opt/local/include/openmpi-mp/ -std=c++11" LDFLAGS=-L/opt/local/lib/openmpi-mp/ LIBS="-lmpi -lmpi_cxx -fopenmp -lgmp -lmpfr" --enable-comms=mpi -;; -clang-avx2-openmp-mpi) -CXX=clang-omp++ ../../configure --enable-simd=AVX2 CXXFLAGS="-mavx2 -mfma -fopenmp -O3 -I/opt/local/include/openmpi-mp/ -std=c++11" LDFLAGS=-L/opt/local/lib/openmpi-mp/ LIBS="-lmpi -lmpi_cxx -fopenmp -lgmp -lmpfr" --enable-comms=mpi -;; -clang-avx-mpi) -CXX=clang++ ../../configure --enable-simd=AVX CXXFLAGS="-mavx -O3 -I/opt/local/include/openmpi-mp/ -std=c++11" LDFLAGS=-L/opt/local/lib/openmpi-mp/ LIBS="-lmpi -lmpi_cxx -lgmp -lmpfr" --enable-comms=mpi -;; -clang-avx2-mpi) -CXX=clang++ ../../configure --enable-simd=AVX2 CXXFLAGS="-mavx2 -mfma -O3 -I/opt/local/include/openmpi-mp/ -std=c++11" LDFLAGS=-L/opt/local/lib/openmpi-mp/ LIBS="-lmpi -lmpi_cxx -lgmp -lmpfr" --enable-comms=mpi -;; -clang-avx2) -CXX=clang++ ../../configure --enable-simd=AVX2 CXXFLAGS="-mavx2 -mfma -O3 -std=c++11" LDFLAGS="-L/usr/local/lib/" LIBS="-lgmp -lmpfr" --enable-comms=none -;; -esac -echo -e $NORMAL diff --git a/scripts/configure-cray b/scripts/configure-cray deleted file mode 100755 index db581493..00000000 --- a/scripts/configure-cray +++ /dev/null @@ -1,10 +0,0 @@ -#!/bin/bash -DIRS="g++-avx-openmp g++-avx clang-xc30 clang-xc30-openmp" - -for D in $DIRS -do - mkdir -p builds/$D - cd builds/$D - ../../scripts/configure-commands $D - cd ../.. -done diff --git a/scripts/configure-mic b/scripts/configure-mic deleted file mode 100755 index 668845fe..00000000 --- a/scripts/configure-mic +++ /dev/null @@ -1,10 +0,0 @@ -#!/bin/bash -DIRS="build-icpc-mic" - -for D in $DIRS -do - mkdir -p $D - cd $D - ../configure-commands - cd .. -done diff --git a/scripts/copyright b/scripts/copyright index 92772f16..cc9ed6e5 100755 --- a/scripts/copyright +++ b/scripts/copyright @@ -12,6 +12,7 @@ Grid physics library, www.github.com/paboyle/Grid Source file: $1 Copyright (C) 2015 +Copyright (C) 2016 EOF @@ -38,8 +39,21 @@ See the full license in the file "LICENSE" in the top level distribution directo /* END LEGAL */ EOF + cat message > tmp.fil -cat $1 >> tmp.fil + +NOTICE=`grep -n "END LEGAL" $1 | awk '{ print $1 }' ` + +if [ "X$NOTICE" != "X" ] +then + echo "found notice ending on line $NOTICE" + awk 'BEGIN { P=0 } { if ( P ) print } /END LEGAL/{P=1} ' $1 >> tmp.fil +else + cat $1 >> tmp.fil + +fi + + cp tmp.fil $1 shift diff --git a/scripts/cray-modules b/scripts/cray-modules deleted file mode 100644 index 95de2b0b..00000000 --- a/scripts/cray-modules +++ /dev/null @@ -1,2 +0,0 @@ -module swap PrgEnv-cray PrgEnv-intel -module swap intel/14.0.4.211 intel/15.0.2.164 diff --git a/scripts/reconfigure_script b/scripts/reconfigure_script deleted file mode 100755 index d8d7212d..00000000 --- a/scripts/reconfigure_script +++ /dev/null @@ -1,4 +0,0 @@ -aclocal -I m4 -autoheader -f -automake -f --add-missing -autoconf -f diff --git a/scripts/update_fftw.sh b/scripts/update_fftw.sh deleted file mode 100755 index 20e42080..00000000 --- a/scripts/update_fftw.sh +++ /dev/null @@ -1,18 +0,0 @@ -#!/usr/bin/env bash - -if (( $# != 1 )); then - echo "usage: `basename $0` " 1>&2 - exit 1 -fi -ARC=$1 - -INITDIR=`pwd` -rm -rf lib/fftw -mkdir lib/fftw - -ARCDIR=`tar -tf ${ARC} | head -n1 | sed -e 's@/.*@@'` -tar -xf ${ARC} -cp ${ARCDIR}/api/fftw3.h lib/fftw/ - -cd ${INITDIR} -rm -rf ${ARCDIR} diff --git a/scripts/wilson.gnu b/scripts/wilson.gnu deleted file mode 100644 index 69bca5b5..00000000 --- a/scripts/wilson.gnu +++ /dev/null @@ -1,7 +0,0 @@ -plot 'wilson.t1' u 2 w l t "AVX1-OMP=1" -replot 'wilson.t2' u 2 w l t "AVX1-OMP=2" -replot 'wilson.t4' u 2 w l t "AVX1-OMP=4" -set terminal 'pdf' -set output 'wilson_clang.pdf' -replot -quit From bb94ddd0ebca02fb7c82b47e1befd411fafe0283 Mon Sep 17 00:00:00 2001 From: paboyle Date: Wed, 2 Nov 2016 08:07:09 +0000 Subject: [PATCH 02/14] Tidy up of mpi3; also some cleaning of the dslash controls. --- benchmarks/Benchmark_comms.cc | 134 +------------- benchmarks/Benchmark_dwf.cc | 63 +++++-- benchmarks/Benchmark_dwf_ntpf.cc | 153 ---------------- benchmarks/Benchmark_dwf_sweep.cc | 16 +- benchmarks/Benchmark_wilson_sweep.cc | 13 ++ benchmarks/Benchmark_zmm.cc | 175 ------------------- lib/Init.cc | 156 +++++++++++------ lib/communicator/Communicator_base.cc | 9 +- lib/communicator/Communicator_base.h | 2 +- lib/communicator/Communicator_mpi3_leader.cc | 147 ++++++++-------- lib/qcd/action/fermion/WilsonKernels.cc | 3 +- lib/qcd/action/fermion/WilsonKernels.h | 80 ++++++--- lib/stencil/Lebesgue.cc | 2 +- 13 files changed, 321 insertions(+), 632 deletions(-) delete mode 100644 benchmarks/Benchmark_dwf_ntpf.cc delete mode 100644 benchmarks/Benchmark_zmm.cc diff --git a/benchmarks/Benchmark_comms.cc b/benchmarks/Benchmark_comms.cc index 234d2fbb..eab9b9b9 100644 --- a/benchmarks/Benchmark_comms.cc +++ b/benchmarks/Benchmark_comms.cc @@ -48,8 +48,8 @@ int main (int argc, char ** argv) std::cout< latt_size ({lat*mpi_layout[0], @@ -124,7 +124,7 @@ int main (int argc, char ** argv) std::cout< latt_size ({lat,lat,lat,lat}); @@ -199,7 +199,7 @@ int main (int argc, char ** argv) std::cout< latt_size ({lat*mpi_layout[0], @@ -271,131 +271,5 @@ int main (int argc, char ** argv) std::cout< latt_size ({lat,lat,lat,lat}); - - GridCartesian Grid(latt_size,simd_layout,mpi_layout); - - std::vector > xbuf(8,std::vector(lat*lat*lat*Ls)); - std::vector > rbuf(8,std::vector(lat*lat*lat*Ls)); - - - int ncomm; - int bytes=lat*lat*lat*Ls*sizeof(HalfSpinColourVectorD); - - - std::vector empty; - std::vector > requests_fwd(Nd,empty); - std::vector > requests_bwd(Nd,empty); - - for(int mu=0;mu<4;mu++){ - ncomm=0; - if (mpi_layout[mu]>1 ) { - ncomm++; - - int comm_proc; - int xmit_to_rank; - int recv_from_rank; - - comm_proc=1; - Grid.ShiftedRanks(mu,comm_proc,xmit_to_rank,recv_from_rank); - Grid.SendToRecvFromInit(requests_fwd[mu], - (void *)&xbuf[mu][0], - xmit_to_rank, - (void *)&rbuf[mu][0], - recv_from_rank, - bytes); - - comm_proc = mpi_layout[mu]-1; - Grid.ShiftedRanks(mu,comm_proc,xmit_to_rank,recv_from_rank); - Grid.SendToRecvFromInit(requests_bwd[mu], - (void *)&xbuf[mu+4][0], - xmit_to_rank, - (void *)&rbuf[mu+4][0], - recv_from_rank, - bytes); - - } - } - - { - double start=usecond(); - for(int i=0;i1 ) { - - Grid.SendToRecvFromBegin(requests_fwd[mu]); - Grid.SendToRecvFromComplete(requests_fwd[mu]); - Grid.SendToRecvFromBegin(requests_bwd[mu]); - Grid.SendToRecvFromComplete(requests_bwd[mu]); - } - } - Grid.Barrier(); - } - - double stop=usecond(); - - double dbytes = bytes; - double xbytes = Nloop*dbytes*2.0*ncomm; - double rbytes = xbytes; - double bidibytes = xbytes+rbytes; - - double time = stop-start; - - std::cout<1 ) { - - Grid.SendToRecvFromBegin(requests_fwd[mu]); - Grid.SendToRecvFromBegin(requests_bwd[mu]); - Grid.SendToRecvFromComplete(requests_fwd[mu]); - Grid.SendToRecvFromComplete(requests_bwd[mu]); - } - } - Grid.Barrier(); - } - - double stop=usecond(); - - double dbytes = bytes; - double xbytes = Nloop*dbytes*2.0*ncomm; - double rbytes = xbytes; - double bidibytes = xbytes+rbytes; - - double time = stop-start; - - std::cout< WilsonFermion5DR; typedef WilsonFermion5D WilsonFermion5DF; typedef WilsonFermion5D WilsonFermion5DD; @@ -54,10 +53,6 @@ int main (int argc, char ** argv) { Grid_init(&argc,&argv); - if( GridCmdOptionExists(argv,argv+argc,"--asynch") ){ - overlapComms = true; - } - int threads = GridThread::GetThreads(); std::cout<_Nprocessors; - for(int doasm=1;doasm<2;doasm++){ - - QCD::WilsonKernelsStatic::AsmOpt=doasm; - DomainWallFermionR Dw(Umu,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5); - - std::cout<::Dhop "< WilsonFermion5DR; LatticeFermion ssrc(sFGrid); LatticeFermion sref(sFGrid); @@ -248,6 +261,16 @@ int main (int argc, char ** argv) sr_e = zero; sr_o = zero; + std::cout << GridLogMessage<< "*********************************************************" <::DhopEO "< -Author: paboyle - - This program is free software; you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation; either version 2 of the License, or - (at your option) any later version. - - This program is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License along - with this program; if not, write to the Free Software Foundation, Inc., - 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. - - See the full license in the file "LICENSE" in the top level distribution directory - *************************************************************************************/ - /* END LEGAL */ -#include - -using namespace std; -using namespace Grid; -using namespace Grid::QCD; - -template -struct scal { - d internal; -}; - - Gamma::GammaMatrix Gmu [] = { - Gamma::GammaX, - Gamma::GammaY, - Gamma::GammaZ, - Gamma::GammaT - }; - -bool overlapComms = false; - - -int main (int argc, char ** argv) -{ - Grid_init(&argc,&argv); - - if( GridCmdOptionExists(argv,argv+argc,"--asynch") ){ - overlapComms = true; - } - - int threads = GridThread::GetThreads(); - std::cout< latt4 = GridDefaultLatt(); - const int Ls=16; - GridCartesian * UGrid = SpaceTimeGrid::makeFourDimGrid(GridDefaultLatt(), GridDefaultSimd(Nd,vComplex::Nsimd()),GridDefaultMpi()); - GridRedBlackCartesian * UrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid); - GridCartesian * FGrid = SpaceTimeGrid::makeFiveDimGrid(Ls,UGrid); - GridRedBlackCartesian * FrbGrid = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls,UGrid); - - std::vector seeds4({1,2,3,4}); - std::vector seeds5({5,6,7,8}); - - GridParallelRNG RNG4(UGrid); RNG4.SeedFixedIntegers(seeds4); - GridParallelRNG RNG5(FGrid); RNG5.SeedFixedIntegers(seeds5); - - LatticeFermion src (FGrid); random(RNG5,src); - LatticeFermion result(FGrid); result=zero; - LatticeFermion ref(FGrid); ref=zero; - LatticeFermion tmp(FGrid); - LatticeFermion err(FGrid); - - ColourMatrix cm = Complex(1.0,0.0); - - LatticeGaugeField Umu(UGrid); - random(RNG4,Umu); - - LatticeGaugeField Umu5d(FGrid); - - // replicate across fifth dimension - for(int ss=0;ssoSites();ss++){ - for(int s=0;s U(4,FGrid); - for(int mu=0;mu(Umu5d,mu); - } - - if (1) - { - ref = zero; - for(int mu=0;mu_Nprocessors; - - - QCD::WilsonKernelsStatic::AsmOpt=1; - - DomainWallFermionR Dw(Umu,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5,params); - - std::cout< seeds({1,2,3,4}); RealD mass = 0.1; + std::cout << GridLogMessage<< "*****************************************************************" < - - This program is free software; you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation; either version 2 of the License, or - (at your option) any later version. - - This program is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License along - with this program; if not, write to the Free Software Foundation, Inc., - 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. - - See the full license in the file "LICENSE" in the top level distribution directory - *************************************************************************************/ - /* END LEGAL */ -#include - - -using namespace Grid; -using namespace Grid::QCD; - - -int bench(std::ofstream &os, std::vector &latt4,int Ls); - -int main(int argc,char **argv) -{ - Grid_init(&argc,&argv); - std::ofstream os("zmm.dat"); - - os << "#V Ls Lxy Lzt C++ Asm OMP L1 " < grid({L,L,m*L,m*L}); - std::cout << GridLogMessage <<"\t"; - for(int i=0;i<4;i++) { - std::cout << grid[i]<<"x"; - } - std::cout << Ls<<"\t\t"; - bench(os,grid,Ls); - } - } - } -} - -int bench(std::ofstream &os, std::vector &latt4,int Ls) -{ - - GridCartesian * UGrid = SpaceTimeGrid::makeFourDimGrid(latt4, GridDefaultSimd(Nd,vComplex::Nsimd()),GridDefaultMpi()); - GridRedBlackCartesian * UrbGrid = SpaceTimeGrid::makeFourDimRedBlackGrid(UGrid); - GridCartesian * FGrid = SpaceTimeGrid::makeFiveDimGrid(Ls,UGrid); - GridRedBlackCartesian * FrbGrid = SpaceTimeGrid::makeFiveDimRedBlackGrid(Ls,UGrid); - - std::vector simd_layout = GridDefaultSimd(Nd,vComplex::Nsimd()); - std::vector mpi_layout = GridDefaultMpi(); - int threads = GridThread::GetThreads(); - - std::vector seeds4({1,2,3,4}); - std::vector seeds5({5,6,7,8}); - - GridSerialRNG sRNG; sRNG.SeedFixedIntegers(seeds4); - - LatticeFermion src (FGrid); - LatticeFermion tmp (FGrid); - LatticeFermion srce(FrbGrid); - - LatticeFermion resulto(FrbGrid); resulto=zero; - LatticeFermion resulta(FrbGrid); resulta=zero; - LatticeFermion junk(FrbGrid); junk=zero; - LatticeFermion diff(FrbGrid); - LatticeGaugeField Umu(UGrid); - - double mfc, mfa, mfo, mfl1; - - GridParallelRNG RNG4(UGrid); RNG4.SeedFixedIntegers(seeds4); - GridParallelRNG RNG5(FGrid); RNG5.SeedFixedIntegers(seeds5); - random(RNG5,src); -#if 1 - random(RNG4,Umu); -#else - int mmu=2; - std::vector U(4,UGrid); - for(int mu=0;mu(Umu,mu); - if ( mu!=mmu ) U[mu] = zero; - if ( mu==mmu ) U[mu] = 1.0; - PokeIndex(Umu,U[mu],mu); - } -#endif - pickCheckerboard(Even,srce,src); - - RealD mass=0.1; - RealD M5 =1.8; - DomainWallFermionR Dw(Umu,*FGrid,*FrbGrid,*UGrid,*UrbGrid,mass,M5); - - int ncall=50; - double t0=usecond(); - for(int i=0;i & vec) return; } +void GridCmdOptionInt(std::string &str,int & val) +{ + std::stringstream ss(str); + ss>>val; + return; +} + void GridParseLayout(char **argv,int argc, std::vector &latt, @@ -153,14 +160,12 @@ void GridParseLayout(char **argv,int argc, assert(ompthreads.size()==1); GridThread::SetThreads(ompthreads[0]); } - if( GridCmdOptionExists(argv,argv+argc,"--cores") ){ - std::vector cores(0); + int cores; arg= GridCmdOptionPayload(argv,argv+argc,"--cores"); - GridCmdOptionIntVector(arg,cores); - GridThread::SetCores(cores[0]); + GridCmdOptionInt(arg,cores); + GridThread::SetCores(cores); } - } std::string GridCmdVectorIntToString(const std::vector & vec){ @@ -169,7 +174,7 @@ std::string GridCmdVectorIntToString(const std::vector & vec){ return oss.str(); } ///////////////////////////////////////////////////////// -// +// Reinit guard ///////////////////////////////////////////////////////// static int Grid_is_initialised = 0; @@ -178,27 +183,31 @@ void Grid_init(int *argc,char ***argv) { GridLogger::StopWatch.Start(); + std::string arg; + + //////////////////////////////////// + // Shared memory block size + //////////////////////////////////// + if( GridCmdOptionExists(*argv,*argv+*argc,"--shm") ){ + int MB; + arg= GridCmdOptionPayload(*argv,*argv+*argc,"--shm"); + GridCmdOptionInt(arg,MB); + CartesianCommunicator::MAX_MPI_SHM_BYTES = MB*1024*1024; + } + CartesianCommunicator::Init(argc,argv); - // Parse command line args. + //////////////////////////////////// + // Logging + //////////////////////////////////// - std::string arg; std::vector logstreams; std::string defaultLog("Error,Warning,Message,Performance"); - GridCmdOptionCSL(defaultLog,logstreams); GridLogConfigure(logstreams); - if( GridCmdOptionExists(*argv,*argv+*argc,"--help") ){ - std::cout<= MAX_MPI_SHM_BYTES) { + std::cout<< " ShmBufferMalloc exceeded shared heap size -- try increasing with --shm flag" < /* END LEGAL */ #include "Grid.h" #include +#include namespace Grid { @@ -45,6 +46,7 @@ const int pool = 48; class SlaveState { public: + sem_t sem; volatile int head; volatile int start; volatile int tail; @@ -56,29 +58,32 @@ public: SlaveState *state; MPI_Comm squadron; uint64_t base; + int universe_rank; + int vertical_rank; //////////////////////////////////////////////////////////// // Descriptor circular pointers //////////////////////////////////////////////////////////// Slave() {}; - void Init(SlaveState * _state,MPI_Comm _squadron); + void Init(SlaveState * _state,MPI_Comm _squadron,int _universe_rank,int _vertical_rank); void EventLoop (void) { - std::cerr<< " Entering even loop "<tail != state->head ); - std::cerr<< " FIFO drained "<< state->tail <"<"<"<>0 )&0xFFFF)^((icomm>>16)&0xFFFF) ^ ((icomm>>32)&0xFFFF)^((icomm>>48)&0xFFFF); - hashtag = (comm_hash<<15) | tag; + // hashtag = (comm_hash<<15) | tag; + hashtag = tag; }; -void Slave::Init(SlaveState * _state,MPI_Comm _squadron) +void Slave::Init(SlaveState * _state,MPI_Comm _squadron,int _universe_rank,int _vertical_rank) { squadron=_squadron; + universe_rank=_universe_rank; + vertical_rank=_vertical_rank; state =_state; + std::cout << "state "<<_state<<" comm "<<_squadron<<" universe_rank"<head = state->tail = state->start = 0; - MPI_Barrier(squadron); base = (uint64_t)MPIoffloadEngine::VerticalShmBufs[0]; int rank; MPI_Comm_rank(_squadron,&rank); } #define PERI_PLUS(A) ( (A+1)%pool ) -void Slave::Event (void) { +int Slave::Event (void) { static int tail_last; static int head_last; static int start_last; int ierr; - if ( (state->tail != tail_last) - ||(state->head != head_last) - ||(state->start != start_last) - ) { - std::cerr<< " Event loop "<< state->tail << " "<< state->start<< " "<< state->head <tail; - if ( t != state->start ) { - int flag=0; - - std::cerr<< " Testing tail "<< t<<" "<< (void *)&state->Descrs[t].request - << " "<Descrs[t].request<Descrs[t].request,&flag,MPI_STATUS_IGNORE); - // ierr=MPI_Test((MPI_Request *)&state->Descrs[t].request,&flag,MPI_STATUS_IGNORE); - assert(ierr==0); - if ( flag ) { - state->tail = PERI_PLUS(t); - std::cerr<< " Tail advanced from "<< t<head ) { switch ( state->Descrs[s].command ) { case COMMAND_ISEND: - std::cerr<< " Send "<Descrs[s].buf<< "["<Descrs[s].bytes<<"]" - << " to " << state->Descrs[s].rank<< " tag" << state->Descrs[s].tag - << " Comm " << MPIoffloadEngine::communicator_universe<< std::endl; - - std::cerr<< " Request was "<Descrs[s].request<Descrs[s].buf<< "["<Descrs[s].bytes<<"]" + << " to " << state->Descrs[s].rank<< " tag" << state->Descrs[s].tag + << " Comm " << MPIoffloadEngine::communicator_universe<< " me " <Descrs[s].buf+base), state->Descrs[s].bytes, MPI_CHAR, @@ -438,18 +430,17 @@ void Slave::Event (void) { state->Descrs[s].tag, MPIoffloadEngine::communicator_universe, (MPI_Request *)&state->Descrs[s].request); - std::cerr<< " Request is "<Descrs[s].request<Descrs[0].request<start = PERI_PLUS(s); + return 1; break; case COMMAND_IRECV: - std::cerr<< " Recv "<Descrs[s].buf<< "["<Descrs[s].bytes<<"]" - << " to " << state->Descrs[s].rank<< " tag" << state->Descrs[s].tag - << " Comm " << MPIoffloadEngine::communicator_universe<< std::endl; - - std::cerr<< " Request was "<Descrs[s].request<Descrs[s].buf<< "["<Descrs[s].bytes<<"]" + << " from " << state->Descrs[s].rank<< " tag" << state->Descrs[s].tag + << " Comm " << MPIoffloadEngine::communicator_universe<< " me "<< universe_rank<< std::endl; + */ ierr=MPI_Irecv((void *)(state->Descrs[s].buf+base), state->Descrs[s].bytes, MPI_CHAR, @@ -457,30 +448,32 @@ void Slave::Event (void) { state->Descrs[s].tag, MPIoffloadEngine::communicator_universe, (MPI_Request *)&state->Descrs[s].request); - std::cerr<< " Request is "<Descrs[s].request<Descrs[0].request<Descrs[s].request<Descrs[0].request<start = PERI_PLUS(s); + return 1; break; case COMMAND_WAITALL: - std::cerr<< " Wait all "<tail;t!=s; t=PERI_PLUS(t) ){ - std::cerr<< " Wait ["<Descrs[t].request <Descrs[0].request<Descrs[t].request,MPI_STATUS_IGNORE); }; s=PERI_PLUS(s); state->start = s; state->tail = s; + MPI_Barrier(squadron); + return 1; break; default: assert(0); break; } - return; } + return 0; } ////////////////////////////////////////////////////////////////////////////// // External interaction with the queue @@ -500,17 +493,29 @@ uint64_t Slave::QueueCommand(int command,void *buf, int bytes, int tag, MPI_Comm MPI_Comm communicator; MPI_Request request; - MPIoffloadEngine::MapCommRankToWorldRank(hashtag,commrank,tag,comm,worldrank); + MPIoffloadEngine::MapCommRankToWorldRank(hashtag,worldrank,tag,comm,commrank); - int VerticalRank = MPIoffloadEngine::VerticalRank; uint64_t relative= (uint64_t)buf - base; state->Descrs[head].buf = relative; state->Descrs[head].bytes = bytes; - state->Descrs[head].rank = MPIoffloadEngine::UniverseRanks[worldrank][VerticalRank]; + state->Descrs[head].rank = MPIoffloadEngine::UniverseRanks[worldrank][vertical_rank]; state->Descrs[head].tag = hashtag; state->Descrs[head].command= command; - std::cerr<< " QueueCommand "<tail==next ); @@ -671,6 +676,8 @@ void CartesianCommunicator::StencilSendToRecvFromBegin(std::vector= shm) && (xmit_i+bytes <= shm+MAX_MPI_SHM_BYTES) ); assert( (recv_i >= shm) && (recv_i+bytes <= shm+MAX_MPI_SHM_BYTES) ); + assert(from!=_processor); + assert(dest!=_processor); MPIoffloadEngine::QueueMultiplexedSend(xmit,bytes,_processor,communicator,dest); MPIoffloadEngine::QueueMultiplexedRecv(recv,bytes,from,communicator,from); } diff --git a/lib/qcd/action/fermion/WilsonKernels.cc b/lib/qcd/action/fermion/WilsonKernels.cc index 52ee8704..43776c86 100644 --- a/lib/qcd/action/fermion/WilsonKernels.cc +++ b/lib/qcd/action/fermion/WilsonKernels.cc @@ -32,8 +32,7 @@ directory namespace Grid { namespace QCD { -int WilsonKernelsStatic::HandOpt; -int WilsonKernelsStatic::AsmOpt; +int WilsonKernelsStatic::Opt; template WilsonKernels::WilsonKernels(const ImplParams &p) : Base(p){}; diff --git a/lib/qcd/action/fermion/WilsonKernels.h b/lib/qcd/action/fermion/WilsonKernels.h index 669ee4be..919f7540 100644 --- a/lib/qcd/action/fermion/WilsonKernels.h +++ b/lib/qcd/action/fermion/WilsonKernels.h @@ -40,9 +40,9 @@ namespace QCD { //////////////////////////////////////////////////////////////////////////////////////////////////////////////// class WilsonKernelsStatic { public: + enum { OptGeneric, OptHandUnroll, OptInlineAsm }; // S-direction is INNERMOST and takes no part in the parity. - static int AsmOpt; // these are a temporary hack - static int HandOpt; // these are a temporary hack + static int Opt; // these are a temporary hack }; template class WilsonKernels : public FermionOperator , public WilsonKernelsStatic { @@ -56,24 +56,40 @@ public: template typename std::enable_if::type DiracOptDhopSite(StencilImpl &st, LebesgueOrder &lo, DoubledGaugeField &U, SiteHalfSpinor * buf, - int sF, int sU, int Ls, int Ns, const FermionField &in, FermionField &out) { + int sF, int sU, int Ls, int Ns, const FermionField &in, FermionField &out) + { + switch(Opt) { #ifdef AVX512 - if (AsmOpt) { - WilsonKernels::DiracOptAsmDhopSite(st,lo,U,buf,sF,sU,Ls,Ns,in,out); - } else { -#else - { -#endif + case OptInlineAsm: for (int site = 0; site < Ns; site++) { for (int s = 0; s < Ls; s++) { - if (HandOpt) - WilsonKernels::DiracOptHandDhopSite(st,lo,U,buf,sF,sU,in,out); - else - WilsonKernels::DiracOptGenericDhopSite(st,lo,U,buf,sF,sU,in,out); + WilsonKernels::DiracOptAsmDhopSite(st,lo,U,buf,sF,sU,Ls,Ns,in,out); sF++; } sU++; } + break; +#endif + case OptHandUnroll: + for (int site = 0; site < Ns; site++) { + for (int s = 0; s < Ls; s++) { + WilsonKernels::DiracOptHandDhopSite(st,lo,U,buf,sF,sU,in,out); + sF++; + } + sU++; + } + break; + case OptGeneric: + for (int site = 0; site < Ns; site++) { + for (int s = 0; s < Ls; s++) { + WilsonKernels::DiracOptGenericDhopSite(st,lo,U,buf,sF,sU,in,out); + sF++; + } + sU++; + } + break; + default: + assert(0); } } @@ -81,7 +97,7 @@ public: typename std::enable_if<(Impl::Dimension != 3 || (Impl::Dimension == 3 && Nc != 3)) && EnableBool, void>::type DiracOptDhopSite(StencilImpl &st, LebesgueOrder &lo, DoubledGaugeField &U, SiteHalfSpinor * buf, int sF, int sU, int Ls, int Ns, const FermionField &in, FermionField &out) { - + // no kernel choice for (int site = 0; site < Ns; site++) { for (int s = 0; s < Ls; s++) { WilsonKernels::DiracOptGenericDhopSite(st, lo, U, buf, sF, sU, in, out); @@ -95,23 +111,39 @@ public: typename std::enable_if::type DiracOptDhopSiteDag(StencilImpl &st, LebesgueOrder &lo, DoubledGaugeField &U, SiteHalfSpinor * buf, int sF, int sU, int Ls, int Ns, const FermionField &in, FermionField &out) { + + switch(Opt) { #ifdef AVX512 - if (AsmOpt) { - WilsonKernels::DiracOptAsmDhopSiteDag(st,lo,U,buf,sF,sU,Ls,Ns,in,out); - } else { -#else - { -#endif + case OptInlineAsm: for (int site = 0; site < Ns; site++) { for (int s = 0; s < Ls; s++) { - if (HandOpt) - WilsonKernels::DiracOptHandDhopSiteDag(st,lo,U,buf,sF,sU,in,out); - else - WilsonKernels::DiracOptGenericDhopSiteDag(st,lo,U,buf,sF,sU,in,out); + WilsonKernels::DiracOptAsmDhopSiteDag(st,lo,U,buf,sF,sU,Ls,Ns,in,out); sF++; } sU++; } + break; +#endif + case OptHandUnroll: + for (int site = 0; site < Ns; site++) { + for (int s = 0; s < Ls; s++) { + WilsonKernels::DiracOptHandDhopSiteDag(st,lo,U,buf,sF,sU,in,out); + sF++; + } + sU++; + } + break; + case OptGeneric: + for (int site = 0; site < Ns; site++) { + for (int s = 0; s < Ls; s++) { + WilsonKernels::DiracOptGenericDhopSiteDag(st,lo,U,buf,sF,sU,in,out); + sF++; + } + sU++; + } + break; + default: + assert(0); } } diff --git a/lib/stencil/Lebesgue.cc b/lib/stencil/Lebesgue.cc index c34b5c96..c30c14c7 100644 --- a/lib/stencil/Lebesgue.cc +++ b/lib/stencil/Lebesgue.cc @@ -32,7 +32,7 @@ Author: paboyle namespace Grid { int LebesgueOrder::UseLebesgueOrder; -std::vector LebesgueOrder::Block({2,2,2,2}); +std::vector LebesgueOrder::Block({8,2,2,2}); LebesgueOrder::IndexInteger LebesgueOrder::alignup(IndexInteger n){ n--; // 1000 0011 --> 1000 0010 From 32375aca657c87674c321ce5e8fc2a82b48dcec2 Mon Sep 17 00:00:00 2001 From: paboyle Date: Wed, 2 Nov 2016 09:27:20 +0000 Subject: [PATCH 03/14] Semaphore sleep/wake up on remote processes. --- lib/communicator/Communicator_mpi3_leader.cc | 84 +++++++++++++++----- 1 file changed, 65 insertions(+), 19 deletions(-) diff --git a/lib/communicator/Communicator_mpi3_leader.cc b/lib/communicator/Communicator_mpi3_leader.cc index db0638da..f160c681 100644 --- a/lib/communicator/Communicator_mpi3_leader.cc +++ b/lib/communicator/Communicator_mpi3_leader.cc @@ -27,7 +27,24 @@ Author: Peter Boyle /* END LEGAL */ #include "Grid.h" #include + +/// bloody mac os doesn't implement unnamed semaphores since it is "optional" posix. +/// darwin dispatch semaphores don't seem to be multiprocess. +//#ifdef __APPLE__ #include +typedef sem_t *Grid_semaphore; +#define SEM_INIT(S) S = sem_open(sem_name,0,00777,0); +#define SEM_INIT_EXCL(S) sem_unlink(sem_name); S = sem_open(sem_name,O_CREAT|O_EXCL,00777,0); +#define SEM_POST(S) sem_post(S); +#define SEM_WAIT(S) sem_wait(S); +//#else +//#include +//typedef sem_t Grid_semaphore; +//#define SEM_INIT(S) +//#define SEM_INIT_EXCL(S) sem_init(&S,0); +//#define SEM_POST(S) sem_post(&S); +//#define SEM_WAIT(S) sem_wait(&S); +//#endif namespace Grid { @@ -46,7 +63,8 @@ const int pool = 48; class SlaveState { public: - sem_t sem; + Grid_semaphore sem_head; + Grid_semaphore sem_tail; volatile int head; volatile int start; volatile int tail; @@ -60,18 +78,43 @@ public: uint64_t base; int universe_rank; int vertical_rank; + char sem_name [NAME_MAX]; //////////////////////////////////////////////////////////// // Descriptor circular pointers //////////////////////////////////////////////////////////// Slave() {}; void Init(SlaveState * _state,MPI_Comm _squadron,int _universe_rank,int _vertical_rank); - + + void SemInit(void) { + sprintf(sem_name,"/Grid_mpi3_head_%d",universe_rank); + SEM_INIT(state->sem_head); + sprintf(sem_name,"/Grid_mpi3_tail_%d",universe_rank); + SEM_INIT(state->sem_tail); + } + void SemInitExcl(void) { + sprintf(sem_name,"/Grid_mpi3_head_%d",universe_rank); + SEM_INIT_EXCL(state->sem_head); + sprintf(sem_name,"/Grid_mpi3_tail_%d",universe_rank); + SEM_INIT_EXCL(state->sem_tail); + } + void WakeUpDMA(void) { + SEM_POST(state->sem_head); + }; + void WakeUpCompute(void) { + SEM_POST(state->sem_tail); + }; + void WaitForCommand(void) { + SEM_WAIT(state->sem_head); + }; + void WaitForComplete(void) { + SEM_WAIT(state->sem_tail); + }; void EventLoop (void) { std::cout<< " Entering event loop "<tail != state->head ); + QueueCommand(COMMAND_WAITALL,0,0,0,squadron,0); + WakeUpDMA(); + WaitForComplete(); + assert ( state->tail == state->head ); } }; @@ -130,25 +173,22 @@ public: // routines for master proc must handle any communicator ///////////////////////////////////////////////////////// - static uint64_t QueueSend(int slave,void *buf, int bytes, int tag, MPI_Comm comm,int rank) { + static void QueueSend(int slave,void *buf, int bytes, int tag, MPI_Comm comm,int rank) { // std::cout<< " Queueing send "<< bytes<< " slave "<< slave << " to comm "<start = s; state->tail = s; - MPI_Barrier(squadron); + + WakeUpCompute(); + return 1; break; From 757a928f9a70c0f2448cc6844f1112bab9cdf7d9 Mon Sep 17 00:00:00 2001 From: paboyle Date: Wed, 2 Nov 2016 12:37:46 +0000 Subject: [PATCH 04/14] Improvement to use own SHM_OPEN call to avoid openmpi bug. --- benchmarks/Benchmark_comms.cc | 4 +- lib/communicator/Communicator_mpi3_leader.cc | 159 ++++++++++++++----- 2 files changed, 124 insertions(+), 39 deletions(-) diff --git a/benchmarks/Benchmark_comms.cc b/benchmarks/Benchmark_comms.cc index eab9b9b9..de73bc81 100644 --- a/benchmarks/Benchmark_comms.cc +++ b/benchmarks/Benchmark_comms.cc @@ -48,7 +48,7 @@ int main (int argc, char ** argv) std::cout< #include "Grid.h" #include -/// bloody mac os doesn't implement unnamed semaphores since it is "optional" posix. -/// darwin dispatch semaphores don't seem to be multiprocess. -//#ifdef __APPLE__ +//////////////////////////////////////////////////////////////////////////////////////////////////////////////// +/// Workarounds: +/// i) bloody mac os doesn't implement unnamed semaphores since it is "optional" posix. +/// darwin dispatch semaphores don't seem to be multiprocess. +/// +/// ii) openmpi under --mca shmem posix works with two squadrons per node; +/// openmpi under default mca settings (I think --mca shmem mmap) on MacOS makes two squadrons map the SAME +/// memory as each other, despite their living on different communicators. This appears to be a bug in OpenMPI. +/// +//////////////////////////////////////////////////////////////////////////////////////////////////////////////// #include typedef sem_t *Grid_semaphore; -#define SEM_INIT(S) S = sem_open(sem_name,0,00777,0); -#define SEM_INIT_EXCL(S) sem_unlink(sem_name); S = sem_open(sem_name,O_CREAT|O_EXCL,00777,0); -#define SEM_POST(S) sem_post(S); -#define SEM_WAIT(S) sem_wait(S); -//#else -//#include -//typedef sem_t Grid_semaphore; -//#define SEM_INIT(S) -//#define SEM_INIT_EXCL(S) sem_init(&S,0); -//#define SEM_POST(S) sem_post(&S); -//#define SEM_WAIT(S) sem_wait(&S); -//#endif + +#define SEM_INIT(S) S = sem_open(sem_name,0,0600,0); assert ( S != SEM_FAILED ); +#define SEM_INIT_EXCL(S) sem_unlink(sem_name); S = sem_open(sem_name,O_CREAT|O_EXCL,0600,0); assert ( S != SEM_FAILED ); +#define SEM_POST(S) assert ( sem_post(S) == 0 ); +#define SEM_WAIT(S) assert ( sem_wait(S) == 0 ); + +#include + namespace Grid { @@ -63,8 +66,6 @@ const int pool = 48; class SlaveState { public: - Grid_semaphore sem_head; - Grid_semaphore sem_tail; volatile int head; volatile int start; volatile int tail; @@ -73,6 +74,8 @@ public: class Slave { public: + Grid_semaphore sem_head; + Grid_semaphore sem_tail; SlaveState *state; MPI_Comm squadron; uint64_t base; @@ -87,33 +90,38 @@ public: void Init(SlaveState * _state,MPI_Comm _squadron,int _universe_rank,int _vertical_rank); void SemInit(void) { - sprintf(sem_name,"/Grid_mpi3_head_%d",universe_rank); - SEM_INIT(state->sem_head); - sprintf(sem_name,"/Grid_mpi3_tail_%d",universe_rank); - SEM_INIT(state->sem_tail); + sprintf(sem_name,"/Grid_mpi3_sem_head_%d",universe_rank); + printf("SEM_NAME: %s \n",sem_name); + SEM_INIT(sem_head); + sprintf(sem_name,"/Grid_mpi3_sem_tail_%d",universe_rank); + printf("SEM_NAME: %s \n",sem_name); + SEM_INIT(sem_tail); } void SemInitExcl(void) { - sprintf(sem_name,"/Grid_mpi3_head_%d",universe_rank); - SEM_INIT_EXCL(state->sem_head); - sprintf(sem_name,"/Grid_mpi3_tail_%d",universe_rank); - SEM_INIT_EXCL(state->sem_tail); + sprintf(sem_name,"/Grid_mpi3_sem_head_%d",universe_rank); + printf("SEM_INIT_EXCL: %s \n",sem_name); + SEM_INIT_EXCL(sem_head); + sprintf(sem_name,"/Grid_mpi3_sem_tail_%d",universe_rank); + printf("SEM_INIT_EXCL: %s \n",sem_name); + SEM_INIT_EXCL(sem_tail); } void WakeUpDMA(void) { - SEM_POST(state->sem_head); + SEM_POST(sem_head); }; void WakeUpCompute(void) { - SEM_POST(state->sem_tail); + SEM_POST(sem_tail); }; void WaitForCommand(void) { - SEM_WAIT(state->sem_head); + SEM_WAIT(sem_head); }; void WaitForComplete(void) { - SEM_WAIT(state->sem_tail); + SEM_WAIT(sem_tail); }; void EventLoop (void) { std::cout<< " Entering event loop "<tail == state->head ); } }; @@ -176,19 +188,25 @@ public: static void QueueSend(int slave,void *buf, int bytes, int tag, MPI_Comm comm,int rank) { // std::cout<< " Queueing send "<< bytes<< " slave "<< slave << " to comm "< MPIoffloadEngine::VerticalShmBufs; std::vector > MPIoffloadEngine::UniverseRanks; std::vector MPIoffloadEngine::UserCommunicatorToWorldRanks; @@ -274,8 +291,12 @@ void MPIoffloadEngine::CommunicatorInit (MPI_Comm &communicator_world, ///////////////////////////////////////////////////////////////////// // Split into groups that can share memory (Verticals) ///////////////////////////////////////////////////////////////////// - // MPI_Comm_split_type(communicator_universe, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL,&VerticalComm); - MPI_Comm_split(communicator_universe,(UniverseRank/2),UniverseRank,&VerticalComm); +#define MPI_SHARED_MEM_DEBUG +#ifdef MPI_SHARED_MEM_DEBUG + MPI_Comm_split(communicator_universe,(UniverseRank/4),UniverseRank,&VerticalComm); +#else + MPI_Comm_split_type(communicator_universe, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL,&VerticalComm); +#endif MPI_Comm_rank(VerticalComm ,&VerticalRank); MPI_Comm_size(VerticalComm ,&VerticalSize); @@ -308,19 +329,83 @@ void MPIoffloadEngine::CommunicatorInit (MPI_Comm &communicator_world, ////////////////////////////////////////////////////////////////////////////////////////////////////////// // allocate the shared window for our group, pass back Shm info to CartesianCommunicator ////////////////////////////////////////////////////////////////////////////////////////////////////////// + VerticalShmBufs.resize(VerticalSize); + +#undef MPI_SHARED_MEM +#ifdef MPI_SHARED_MEM ierr = MPI_Win_allocate_shared(CartesianCommunicator::MAX_MPI_SHM_BYTES,1,MPI_INFO_NULL,VerticalComm,&ShmCommBuf,&VerticalWindow); ierr|= MPI_Win_lock_all (MPI_MODE_NOCHECK, VerticalWindow); assert(ierr==0); - - std::cout<<"SHM "<0 ) size = sizeof(SlaveState); + + sprintf(shm_name,"/Grid_mpi3_shm_%d_%d",WorldRank,r); + + shm_unlink(shm_name); + + int fd=shm_open(shm_name,O_RDWR|O_CREAT,0600); + if ( fd < 0 ) { + perror("failed shm_open"); + assert(0); + } + + ftruncate(fd, size); + + VerticalShmBufs[r] = mmap(NULL,size, PROT_READ | PROT_WRITE, MAP_SHARED, fd, 0); + + if ( VerticalShmBufs[r] == MAP_FAILED ) { + perror("failed mmap"); + assert(0); + } + + uint64_t * check = (uint64_t *) VerticalShmBufs[r]; + check[0] = WorldRank; + check[1] = r; + + // std::cout<<"SHM "<0 ) size = sizeof(SlaveState); + + sprintf(shm_name,"/Grid_mpi3_shm_%d_%d",WorldRank,r); + + int fd=shm_open(shm_name,O_RDWR|O_CREAT,0600); + if ( fd<0 ) { + perror("failed shm_open"); + assert(0); + } + VerticalShmBufs[r] = mmap(NULL,size, PROT_READ | PROT_WRITE, MAP_SHARED, fd, 0); + + uint64_t * check = (uint64_t *) VerticalShmBufs[r]; + assert(check[0]== WorldRank); + assert(check[1]== r); + std::cerr<<"SHM "< Date: Wed, 2 Nov 2016 19:54:03 +0000 Subject: [PATCH 05/14] Decrease mpi3l verbose --- lib/communicator/Communicator_mpi3_leader.cc | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/lib/communicator/Communicator_mpi3_leader.cc b/lib/communicator/Communicator_mpi3_leader.cc index b6995a44..49a3e948 100644 --- a/lib/communicator/Communicator_mpi3_leader.cc +++ b/lib/communicator/Communicator_mpi3_leader.cc @@ -48,7 +48,6 @@ typedef sem_t *Grid_semaphore; #include - namespace Grid { enum { COMMAND_ISEND, COMMAND_IRECV, COMMAND_WAITALL }; @@ -91,18 +90,18 @@ public: void SemInit(void) { sprintf(sem_name,"/Grid_mpi3_sem_head_%d",universe_rank); - printf("SEM_NAME: %s \n",sem_name); + // printf("SEM_NAME: %s \n",sem_name); SEM_INIT(sem_head); sprintf(sem_name,"/Grid_mpi3_sem_tail_%d",universe_rank); - printf("SEM_NAME: %s \n",sem_name); + // printf("SEM_NAME: %s \n",sem_name); SEM_INIT(sem_tail); } void SemInitExcl(void) { sprintf(sem_name,"/Grid_mpi3_sem_head_%d",universe_rank); - printf("SEM_INIT_EXCL: %s \n",sem_name); + // printf("SEM_INIT_EXCL: %s \n",sem_name); SEM_INIT_EXCL(sem_head); sprintf(sem_name,"/Grid_mpi3_sem_tail_%d",universe_rank); - printf("SEM_INIT_EXCL: %s \n",sem_name); + // printf("SEM_INIT_EXCL: %s \n",sem_name); SEM_INIT_EXCL(sem_tail); } void WakeUpDMA(void) { @@ -118,7 +117,7 @@ public: SEM_WAIT(sem_tail); }; void EventLoop (void) { - std::cout<< " Entering event loop "<head = state->tail = state->start = 0; base = (uint64_t)MPIoffloadEngine::VerticalShmBufs[0]; int rank; MPI_Comm_rank(_squadron,&rank); From 111bfbc6bc60b4edcb696f21c2edff4606fdb5f0 Mon Sep 17 00:00:00 2001 From: paboyle Date: Thu, 3 Nov 2016 11:40:26 +0000 Subject: [PATCH 06/14] notimestamp by default --- lib/Init.cc | 6 ++++-- tests/forces/Test_wilson_force.cc | 18 ++++++++---------- 2 files changed, 12 insertions(+), 12 deletions(-) diff --git a/lib/Init.cc b/lib/Init.cc index d15e4bd1..695bb01b 100644 --- a/lib/Init.cc +++ b/lib/Init.cc @@ -234,7 +234,7 @@ void Grid_init(int *argc,char ***argv) std::cout<(mom,mommu,mu); // fourth order exponential approx - parallel_for(auto i=mom.begin();i Date: Thu, 3 Nov 2016 13:48:07 +0000 Subject: [PATCH 07/14] MPI auto configure fix --- configure.ac | 31 ++++++++++++++++++------------- 1 file changed, 18 insertions(+), 13 deletions(-) diff --git a/configure.ac b/configure.ac index a8ad40dd..9178a11d 100644 --- a/configure.ac +++ b/configure.ac @@ -253,18 +253,23 @@ AC_ARG_ENABLE([comms],[AC_HELP_STRING([--enable-comms=none|mpi|mpi-auto|mpi3|mpi case ${ac_COMMS} in none) AC_DEFINE([GRID_COMMS_NONE],[1],[GRID_COMMS_NONE] ) + comms_type='none' ;; - mpi|mpi-auto) - AC_DEFINE([GRID_COMMS_MPI],[1],[GRID_COMMS_MPI] ) - ;; - mpi3|mpi3-auto) - AC_DEFINE([GRID_COMMS_MPI3],[1],[GRID_COMMS_MPI3] ) - ;; - mpi3l) + mpi3l*) AC_DEFINE([GRID_COMMS_MPI3L],[1],[GRID_COMMS_MPI3L] ) + comms_type='mpi3l' + ;; + mpi3*) + AC_DEFINE([GRID_COMMS_MPI3],[1],[GRID_COMMS_MPI3] ) + comms_type='mpi3' + ;; + mpi*) + AC_DEFINE([GRID_COMMS_MPI],[1],[GRID_COMMS_MPI] ) + comms_type='mpi' ;; shmem) AC_DEFINE([GRID_COMMS_SHMEM],[1],[GRID_COMMS_SHMEM] ) + comms_type='shmem' ;; *) AC_MSG_ERROR([${ac_COMMS} unsupported --enable-comms option]); @@ -282,11 +287,11 @@ case ${ac_COMMS} in ;; esac -AM_CONDITIONAL(BUILD_COMMS_SHMEM,[ test "X${ac_COMMS}X" == "XshmemX" ]) -AM_CONDITIONAL(BUILD_COMMS_MPI,[ test "X${ac_COMMS}X" == "XmpiX" || test "X${ac_COMMS}X" == "Xmpi-autoX" ]) -AM_CONDITIONAL(BUILD_COMMS_MPI3,[ test "X${ac_COMMS}X" == "Xmpi3X"] ) -AM_CONDITIONAL(BUILD_COMMS_MPI3L,[ test "X${ac_COMMS}X" == "Xmpi3lX"] ) -AM_CONDITIONAL(BUILD_COMMS_NONE,[ test "X${ac_COMMS}X" == "XnoneX" ]) +AM_CONDITIONAL(BUILD_COMMS_SHMEM, [ test "${comms_type}X" == "shmemX" ]) +AM_CONDITIONAL(BUILD_COMMS_MPI, [ test "${comms_type}X" == "mpiX" ]) +AM_CONDITIONAL(BUILD_COMMS_MPI3, [ test "${comms_type}X" == "mpi3X" ] ) +AM_CONDITIONAL(BUILD_COMMS_MPI3L, [ test "${comms_type}X" == "mpi3lX" ] ) +AM_CONDITIONAL(BUILD_COMMS_NONE, [ test "${comms_type}X" == "noneX" ]) ############### RNG selection AC_ARG_ENABLE([rng],[AC_HELP_STRING([--enable-rng=ranlux48|mt19937],\ @@ -379,7 +384,7 @@ compiler version : ${ax_cv_gxx_version} ----- BUILD OPTIONS ----------------------------------- SIMD : ${ac_SIMD} Threading : ${ac_openmp} -Communications type : ${ac_COMMS} +Communications type : ${comms_type} Default precision : ${ac_PRECISION} RNG choice : ${ac_RNG} GMP : `if test "x$have_gmp" = xtrue; then echo yes; else echo no; fi` From c65d23935a4e30d8c3a501772c5afd270a27fccc Mon Sep 17 00:00:00 2001 From: Antonin Portelli Date: Thu, 3 Nov 2016 13:48:20 +0000 Subject: [PATCH 08/14] README update --- README.md | 30 +++++++++++++++--------------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/README.md b/README.md index faf86faf..c7461368 100644 --- a/README.md +++ b/README.md @@ -20,7 +20,7 @@ License: GPL v2. Last update Nov 2016. -_Please send all pull requests to the `develop` branch._ +_Please do not send pull requests to the `master` branch which is reserved for releases._ ### Bug report @@ -29,7 +29,7 @@ _To help us tracking and solving more efficiently issues with Grid, please repor When you file an issue, please go though the following checklist: 1. Check that the code is pointing to the `HEAD` of `develop` or any commit in `master` which is tagged with a version number. -2. Give a description of the target platform (CPU, network, compiler). +2. Give a description of the target platform (CPU, network, compiler). Please give the full CPU part description, using for example `cat /proc/cpuinfo | grep 'model name' | uniq` (Linux) or `sysctl machdep.cpu.brand_string` (macOS) and the full output the `--version` option of your compiler. 3. Give the exact `configure` command used. 4. Attach `config.log`. 5. Attach `config.summary`. @@ -45,7 +45,7 @@ are provided, similar to HPF and cmfortran, and user control is given over the m array indices to both MPI tasks and SIMD processing elements. * Identically shaped arrays then be processed with perfect data parallelisation. -* Such identically shapped arrays are called conformable arrays. +* Such identically shaped arrays are called conformable arrays. The transformation is based on the observation that Cartesian array processing involves identical processing to be performed on different regions of the Cartesian array. @@ -127,14 +127,15 @@ make -C tests/ tests The following options can be use with the `--enable-simd=` option to target different communication interfaces: -| `` | Description | -| ------------- | -------------------------------------------- | -| `none` | no communications | -| `mpi[-auto]` | MPI communications | -| `mpi3[-auto]` | MPI communications using MPI 3 shared memory | -| `shmem ` | Cray SHMEM communications | +| `` | Description | +| -------------- | ------------------------------------------------------------- | +| `none` | no communications | +| `mpi[-auto]` | MPI communications | +| `mpi3[-auto]` | MPI communications using MPI 3 shared memory | +| `mpi3l[-auto]` | MPI communications using MPI 3 shared memory and leader model | +| `shmem ` | Cray SHMEM communications | -For `mpi` and `mpi3` the optional `-auto` suffix instructs the `configure` scripts to determine all the necessary compilation and linking flags. This is done by extracting the informations from the MPI wrapper specified in the environment variable `MPICXX` (if not specified `configure` will scan though a list of default names). +For the MPI interfaces the optional `-auto` suffix instructs the `configure` scripts to determine all the necessary compilation and linking flags. This is done by extracting the informations from the MPI wrapper specified in the environment variable `MPICXX` (if not specified `configure` will scan though a list of default names). ### Possible SIMD types @@ -160,7 +161,7 @@ Alternatively, some CPU codenames can be directly used: | `BGQ` | Blue Gene/Q | #### Notes: -- We currently support AVX512 only for the Intel compiler. Support for GCC and clang will appear in future versions. +- We currently support AVX512 only for the Intel compiler. Support for GCC and clang will appear in future versions of Grid when the AVX512 support within GCC and clang will be more advanced. - For BG/Q only [bgclang](http://trac.alcf.anl.gov/projects/llvm-bgq) is supported. We do not presently plan to support more compilers for this platform. - BG/Q performances are currently rather poor. This is being investigated for future versions. @@ -171,7 +172,7 @@ The following configuration is recommended for the Intel Knights Landing platfor ``` bash ../configure --enable-precision=double\ --enable-simd=KNL \ - --enable-comms=mpi3-auto \ + --enable-comms=mpi-auto \ --with-gmp= \ --with-mpfr= \ --enable-mkl \ @@ -183,10 +184,9 @@ where `` is the UNIX prefix where GMP and MPFR are installed. If you are w ``` bash ../configure --enable-precision=double\ --enable-simd=KNL \ - --enable-comms=mpi3 \ + --enable-comms=mpi \ --with-gmp= \ --with-mpfr= \ --enable-mkl \ CXX=CC CC=cc -``` - +``` \ No newline at end of file From c30d96ea5097ab1760a3f0a6ea1aed8ac1e6142b Mon Sep 17 00:00:00 2001 From: James Harrison Date: Wed, 9 Nov 2016 11:06:20 +0000 Subject: [PATCH 09/14] QedFVol: x86intrin.h namespace fix --- lib/PerfCount.h | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/lib/PerfCount.h b/lib/PerfCount.h index 9ac58883..5ab07c02 100644 --- a/lib/PerfCount.h +++ b/lib/PerfCount.h @@ -43,6 +43,9 @@ Author: paboyle #else #include #endif +#ifdef __x86_64__ +#include +#endif namespace Grid { @@ -86,7 +89,6 @@ inline uint64_t cyclecount(void){ return tmp; } #elif defined __x86_64__ -#include inline uint64_t cyclecount(void){ return __rdtsc(); // unsigned int dummy; From cf167d0cd1c561bed3557eaf89350b8d8eb8d9b1 Mon Sep 17 00:00:00 2001 From: James Harrison Date: Mon, 14 Nov 2016 17:02:29 +0000 Subject: [PATCH 10/14] QedFVol: implement exponentiation of photon field --- programs/qed-fvol/WilsonLoops.h | 19 ++++++++++--------- programs/qed-fvol/qed-fvol.cc | 32 +++++++++++++++++++++++++------- 2 files changed, 35 insertions(+), 16 deletions(-) diff --git a/programs/qed-fvol/WilsonLoops.h b/programs/qed-fvol/WilsonLoops.h index 610fdc7b..c40fbaf3 100644 --- a/programs/qed-fvol/WilsonLoops.h +++ b/programs/qed-fvol/WilsonLoops.h @@ -5,7 +5,7 @@ BEGIN_QEDFVOL_NAMESPACE -template class WilsonLoops : public Gimpl { +template class NewWilsonLoops : public Gimpl { public: INHERIT_GIMPL_TYPES(Gimpl); @@ -55,7 +55,7 @@ public: ////////////////////////////////////////////////// // sum over all x,y,z,t and over all planes of plaquette ////////////////////////////////////////////////// - static RealD sumPlaquette(const GaugeLorentz &Umu) { + static Real sumPlaquette(const GaugeLorentz &Umu) { std::vector U(4, Umu._grid); for (int mu = 0; mu < Nd; mu++) { @@ -73,8 +73,8 @@ public: ////////////////////////////////////////////////// // average over all x,y,z,t and over all planes of plaquette ////////////////////////////////////////////////// - static RealD avgPlaquette(const GaugeLorentz &Umu) { - RealD sumplaq = sumPlaquette(Umu); + static Real avgPlaquette(const GaugeLorentz &Umu) { + Real sumplaq = sumPlaquette(Umu); double vol = Umu._grid->gSites(); double faces = (1.0 * Nd * (Nd - 1)) / 2.0; return sumplaq / vol / faces / Nc; // Nd , Nc dependent... FIXME @@ -112,14 +112,14 @@ public: const int Rmu, const int Rnu, const int mu, const int nu) { GaugeMat sp(U[0]._grid); - WilsonLoop(sp, U, Rmu, Rnu, mu, nu); + wilsonLoop(sp, U, Rmu, Rnu, mu, nu); wl = trace(sp); } ////////////////////////////////////////////////// // sum over all planes of Wilson loop ////////////////////////////////////////////////// static void siteWilsonLoop(LatticeComplex &Wl, - const std::vector &U + const std::vector &U, const int R1, const int R2) { LatticeComplex siteWl(U[0]._grid); Wl = zero; @@ -135,7 +135,7 @@ public: ////////////////////////////////////////////////// // sum over all x,y,z,t and over all planes of Wilson loop ////////////////////////////////////////////////// - static RealD sumWilsonLoop(const GaugeLorentz &Umu, + static Real sumWilsonLoop(const GaugeLorentz &Umu, const int R1, const int R2) { std::vector U(4, Umu._grid); @@ -154,13 +154,14 @@ public: ////////////////////////////////////////////////// // average over all x,y,z,t and over all planes of Wilson loop ////////////////////////////////////////////////// - static RealD avgPlaquette(const GaugeLorentz &Umu, + static Real avgWilsonLoop(const GaugeLorentz &Umu, const int R1, const int R2) { - RealD sumWl = sumWilsonLoop(Umu); + Real sumWl = sumWilsonLoop(Umu, R1, R2); double vol = Umu._grid->gSites(); double faces = 1.0 * Nd * (Nd - 1); return sumWl / vol / faces / Nc; // Nd , Nc dependent... FIXME } +}; END_QEDFVOL_NAMESPACE diff --git a/programs/qed-fvol/qed-fvol.cc b/programs/qed-fvol/qed-fvol.cc index fd780edf..68705b8f 100644 --- a/programs/qed-fvol/qed-fvol.cc +++ b/programs/qed-fvol/qed-fvol.cc @@ -1,4 +1,5 @@ #include +#include using namespace Grid; using namespace QCD; @@ -24,10 +25,11 @@ public: typedef Lattice GaugeField; }; -typedef QedGimpl QedGimplR; -typedef Photon PhotonR; -typedef PhotonR::GaugeField EmField; -typedef PhotonR::GaugeLinkField EmComp; +typedef QedGimpl QedGimplR; +typedef PeriodicGaugeImpl QedPeriodicGimplR; +typedef Photon PhotonR; +typedef PhotonR::GaugeField EmField; +typedef PhotonR::GaugeLinkField EmComp; int main(int argc, char *argv[]) { @@ -60,11 +62,18 @@ int main(int argc, char *argv[]) PhotonR photon(PhotonR::Gauge::Feynman, PhotonR::ZmScheme::QedL); EmField a(&grid); + EmField expA(&grid); + + Real avgPlaqAexp, avgWl2x2Aexp; pRNG.SeedRandomDevice(); photon.StochasticField(a, pRNG); - // Calculate log of plaquette + // Exponentiate photon field + Complex imag_unit(0, 1); + expA = exp(imag_unit*0.5*(a+conjugate(a))); + + // Calculate plaquette from photon field EmComp plaqA(&grid); EmComp wlA(&grid); EmComp tmp(&grid); @@ -105,8 +114,17 @@ int main(int argc, char *argv[]) peekSite(tplaqsite, plaqtrace, site0); Complex plaqsite = TensorRemove(tplaqsite); - LOG(Message) << "Plaquette average: " << avgPlaqA << std::endl; - LOG(Message) << "2x2 Wilson Loop average: " << avgWlA << std::endl; + // Calculate plaquette from exponentiated photon field + avgPlaqAexp = NewWilsonLoops::avgPlaquette(expA); + avgWl2x2Aexp = NewWilsonLoops::avgWilsonLoop(expA, 2, 2); + + avgPlaqAexp = avgPlaqAexp*3; + avgWl2x2Aexp = avgWl2x2Aexp*3; + + LOG(Message) << "Plaquette average (from A): " << avgPlaqA << std::endl; + LOG(Message) << "Plaquette average (from exp(A)): " << avgPlaqAexp << std::endl; + LOG(Message) << "2x2 Wilson Loop average (from A): " << avgWlA << std::endl; + LOG(Message) << "2x2 Wilson Loop average (from exp(A)): " << avgWl2x2Aexp << std::endl; LOG(Message) << "Plaquette (one site): " << plaqsite / faces << std::endl; // epilogue From f4ebea3381046026276864f3f908cb10b114d6a5 Mon Sep 17 00:00:00 2001 From: James Harrison Date: Mon, 14 Nov 2016 17:51:53 +0000 Subject: [PATCH 11/14] QedFVol: add functions for computing spatial and timelike Wilson loops --- programs/qed-fvol/WilsonLoops.h | 117 +++++++++++++++++++++++++++++--- programs/qed-fvol/qed-fvol.cc | 8 ++- 2 files changed, 114 insertions(+), 11 deletions(-) diff --git a/programs/qed-fvol/WilsonLoops.h b/programs/qed-fvol/WilsonLoops.h index c40fbaf3..98db6b7a 100644 --- a/programs/qed-fvol/WilsonLoops.h +++ b/programs/qed-fvol/WilsonLoops.h @@ -45,7 +45,7 @@ public: const std::vector &U) { LatticeComplex sitePlaq(U[0]._grid); Plaq = zero; - for (int mu = 1; mu < Nd; mu++) { + for (int mu = 1; mu < U[0]._grid->_ndimension; mu++) { for (int nu = 0; nu < mu; nu++) { traceDirPlaquette(sitePlaq, U, mu, nu); Plaq = Plaq + sitePlaq; @@ -58,7 +58,7 @@ public: static Real sumPlaquette(const GaugeLorentz &Umu) { std::vector U(4, Umu._grid); - for (int mu = 0; mu < Nd; mu++) { + for (int mu = 0; mu < Umu._grid->_ndimension; mu++) { U[mu] = PeekIndex(Umu, mu); } @@ -74,10 +74,11 @@ public: // average over all x,y,z,t and over all planes of plaquette ////////////////////////////////////////////////// static Real avgPlaquette(const GaugeLorentz &Umu) { + int ndim = Umu._grid->_ndimension; Real sumplaq = sumPlaquette(Umu); - double vol = Umu._grid->gSites(); - double faces = (1.0 * Nd * (Nd - 1)) / 2.0; - return sumplaq / vol / faces / Nc; // Nd , Nc dependent... FIXME + Real vol = Umu._grid->gSites(); + Real faces = (1.0 * ndim * (ndim - 1)) / 2.0; + return sumplaq / vol / faces / Nc; // Nc dependent... FIXME } ////////////////////////////////////////////////// @@ -123,7 +124,42 @@ public: const int R1, const int R2) { LatticeComplex siteWl(U[0]._grid); Wl = zero; - for (int mu = 1; mu < Nd; mu++) { + for (int mu = 1; mu < U[0]._grid->_ndimension; mu++) { + for (int nu = 0; nu < mu; nu++) { + traceWilsonLoop(siteWl, U, R1, R2, mu, nu); + Wl = Wl + siteWl; + traceWilsonLoop(siteWl, U, R2, R1, mu, nu); + Wl = Wl + siteWl; + } + } + } + ////////////////////////////////////////////////// + // sum over planes of Wilson loop with length R1 + // in the time direction + ////////////////////////////////////////////////// + static void siteTimelikeWilsonLoop(LatticeComplex &Wl, + const std::vector &U, + const int R1, const int R2) { + LatticeComplex siteWl(U[0]._grid); + + int ndim = U[0]._grid->_ndimension; + + Wl = zero; + for (int nu = 0; nu < ndim - 1; nu++) { + traceWilsonLoop(siteWl, U, R1, R2, ndim-1, nu); + Wl = Wl + siteWl; + } + } + ////////////////////////////////////////////////// + // sum Wilson loop over all planes orthogonal to the time direction + ////////////////////////////////////////////////// + static void siteSpatialWilsonLoop(LatticeComplex &Wl, + const std::vector &U, + const int R1, const int R2) { + LatticeComplex siteWl(U[0]._grid); + + Wl = zero; + for (int mu = 1; mu < U[0]._grid->_ndimension - 1; mu++) { for (int nu = 0; nu < mu; nu++) { traceWilsonLoop(siteWl, U, R1, R2, mu, nu); Wl = Wl + siteWl; @@ -139,7 +175,7 @@ public: const int R1, const int R2) { std::vector U(4, Umu._grid); - for (int mu = 0; mu < Nd; mu++) { + for (int mu = 0; mu < Umu._grid->_ndimension; mu++) { U[mu] = PeekIndex(Umu, mu); } @@ -152,14 +188,75 @@ public: return p.real(); } ////////////////////////////////////////////////// + // sum over all x,y,z,t and over all planes of timelike Wilson loop + ////////////////////////////////////////////////// + static Real sumTimelikeWilsonLoop(const GaugeLorentz &Umu, + const int R1, const int R2) { + std::vector U(4, Umu._grid); + + for (int mu = 0; mu < Umu._grid->_ndimension; mu++) { + U[mu] = PeekIndex(Umu, mu); + } + + LatticeComplex Wl(Umu._grid); + + siteTimelikeWilsonLoop(Wl, U, R1, R2); + + TComplex Tp = sum(Wl); + Complex p = TensorRemove(Tp); + return p.real(); + } + ////////////////////////////////////////////////// + // sum over all x,y,z,t and over all planes of spatial Wilson loop + ////////////////////////////////////////////////// + static Real sumSpatialWilsonLoop(const GaugeLorentz &Umu, + const int R1, const int R2) { + std::vector U(4, Umu._grid); + + for (int mu = 0; mu < Umu._grid->_ndimension; mu++) { + U[mu] = PeekIndex(Umu, mu); + } + + LatticeComplex Wl(Umu._grid); + + siteSpatialWilsonLoop(Wl, U, R1, R2); + + TComplex Tp = sum(Wl); + Complex p = TensorRemove(Tp); + return p.real(); + } + ////////////////////////////////////////////////// // average over all x,y,z,t and over all planes of Wilson loop ////////////////////////////////////////////////// static Real avgWilsonLoop(const GaugeLorentz &Umu, const int R1, const int R2) { + int ndim = Umu._grid->_ndimension; Real sumWl = sumWilsonLoop(Umu, R1, R2); - double vol = Umu._grid->gSites(); - double faces = 1.0 * Nd * (Nd - 1); - return sumWl / vol / faces / Nc; // Nd , Nc dependent... FIXME + Real vol = Umu._grid->gSites(); + Real faces = 1.0 * ndim * (ndim - 1); + return sumWl / vol / faces / Nc; // Nc dependent... FIXME + } + ////////////////////////////////////////////////// + // average over all x,y,z,t and over all planes of timelike Wilson loop + ////////////////////////////////////////////////// + static Real avgTimelikeWilsonLoop(const GaugeLorentz &Umu, + const int R1, const int R2) { + int ndim = Umu._grid->_ndimension; + Real sumWl = sumTimelikeWilsonLoop(Umu, R1, R2); + Real vol = Umu._grid->gSites(); + Real faces = 1.0 * (ndim - 1); + return sumWl / vol / faces / Nc; // Nc dependent... FIXME + } + ////////////////////////////////////////////////// + // average over all x,y,z,t and over all planes of spatial Wilson loop + ////////////////////////////////////////////////// + static Real avgSpatialWilsonLoop(const GaugeLorentz &Umu, + const int R1, const int R2) { + int ndim = Umu._grid->_ndimension; + Real sumWl = sumSpatialWilsonLoop(Umu, R1, R2); + Real vol = Umu._grid->gSites(); + Real faces = 1.0 * (ndim - 1) * (ndim - 2); + return sumWl / vol / faces / Nc; // Nc dependent... FIXME } }; diff --git a/programs/qed-fvol/qed-fvol.cc b/programs/qed-fvol/qed-fvol.cc index 68705b8f..d026057e 100644 --- a/programs/qed-fvol/qed-fvol.cc +++ b/programs/qed-fvol/qed-fvol.cc @@ -64,7 +64,7 @@ int main(int argc, char *argv[]) EmField a(&grid); EmField expA(&grid); - Real avgPlaqAexp, avgWl2x2Aexp; + Real avgPlaqAexp, avgWl2x2Aexp, avgWl2x2Aexp_time, avgWl2x2Aexp_space; pRNG.SeedRandomDevice(); photon.StochasticField(a, pRNG); @@ -117,14 +117,20 @@ int main(int argc, char *argv[]) // Calculate plaquette from exponentiated photon field avgPlaqAexp = NewWilsonLoops::avgPlaquette(expA); avgWl2x2Aexp = NewWilsonLoops::avgWilsonLoop(expA, 2, 2); + avgWl2x2Aexp_time = NewWilsonLoops::avgTimelikeWilsonLoop(expA, 2, 2); + avgWl2x2Aexp_space = NewWilsonLoops::avgSpatialWilsonLoop(expA, 2, 2); avgPlaqAexp = avgPlaqAexp*3; avgWl2x2Aexp = avgWl2x2Aexp*3; + avgWl2x2Aexp_time = avgWl2x2Aexp_time*3; + avgWl2x2Aexp_space = avgWl2x2Aexp_space*3; LOG(Message) << "Plaquette average (from A): " << avgPlaqA << std::endl; LOG(Message) << "Plaquette average (from exp(A)): " << avgPlaqAexp << std::endl; LOG(Message) << "2x2 Wilson Loop average (from A): " << avgWlA << std::endl; LOG(Message) << "2x2 Wilson Loop average (from exp(A)): " << avgWl2x2Aexp << std::endl; + LOG(Message) << "2x2 Wilson Loop timelike average (from exp(A)): " << avgWl2x2Aexp_time << std::endl; + LOG(Message) << "2x2 Wilson Loop spatial average (from exp(A)): " << avgWl2x2Aexp_space << std::endl; LOG(Message) << "Plaquette (one site): " << plaqsite / faces << std::endl; // epilogue From 92ec3404f8a404e7d6420ebfa0f113af5eb6ec6d Mon Sep 17 00:00:00 2001 From: James Harrison Date: Mon, 14 Nov 2016 17:59:02 +0000 Subject: [PATCH 12/14] Set imaginary part of stochastic QED field to zero after FFT into position space --- lib/qcd/action/gauge/Photon.h | 2 ++ 1 file changed, 2 insertions(+) diff --git a/lib/qcd/action/gauge/Photon.h b/lib/qcd/action/gauge/Photon.h index 852ecb3e..ca0a8d40 100644 --- a/lib/qcd/action/gauge/Photon.h +++ b/lib/qcd/action/gauge/Photon.h @@ -172,6 +172,8 @@ namespace QCD{ pokeLorentz(aTilde, r, mu); } fft.FFT_all_dim(out, aTilde, FFT::backward); + + out = 0.5*(out + conjugate(out)); } // template // void Photon::FeynmanGaugeMomentumSpacePropagator_L(GaugeField &out, From a71b69389b6fc7360dbebbf5ed8d4fa3a6952016 Mon Sep 17 00:00:00 2001 From: James Harrison Date: Mon, 14 Nov 2016 18:23:04 +0000 Subject: [PATCH 13/14] QedFVol: calculate square Wilson loops up to 10x10 --- programs/qed-fvol/qed-fvol.cc | 74 +++++++---------------------------- 1 file changed, 14 insertions(+), 60 deletions(-) diff --git a/programs/qed-fvol/qed-fvol.cc b/programs/qed-fvol/qed-fvol.cc index d026057e..31312b1e 100644 --- a/programs/qed-fvol/qed-fvol.cc +++ b/programs/qed-fvol/qed-fvol.cc @@ -64,75 +64,29 @@ int main(int argc, char *argv[]) EmField a(&grid); EmField expA(&grid); - Real avgPlaqAexp, avgWl2x2Aexp, avgWl2x2Aexp_time, avgWl2x2Aexp_space; + Real wlA, logWlA; pRNG.SeedRandomDevice(); photon.StochasticField(a, pRNG); // Exponentiate photon field Complex imag_unit(0, 1); - expA = exp(imag_unit*0.5*(a+conjugate(a))); + expA = exp(imag_unit*a); - // Calculate plaquette from photon field - EmComp plaqA(&grid); - EmComp wlA(&grid); - EmComp tmp(&grid); - std::vector a_comp(4, &grid); - - for (int dir = 0; dir < Nd; dir++) { - a_comp[dir] = PeekIndex(a, dir); + // Calculate Wilson loops + for(int i=1; i<=10; i++){ + LOG(Message) << i << 'x' << i << " Wilson loop" << std::endl; + wlA = NewWilsonLoops::avgWilsonLoop(expA, i, i) * 3; + logWlA = -2*log(wlA); + LOG(Message) << "-2log(W) average: " << logWlA << std::endl; + wlA = NewWilsonLoops::avgTimelikeWilsonLoop(expA, i, i) * 3; + logWlA = -2*log(wlA); + LOG(Message) << "-2log(W) timelike: " << logWlA << std::endl; + wlA = NewWilsonLoops::avgSpatialWilsonLoop(expA, i, i) * 3; + logWlA = -2*log(wlA); + LOG(Message) << "-2log(W) spatial: " << logWlA << std::endl; } - plaqA = zero; - wlA = zero; - - for(int mu = 1; mu < Nd; mu++) { - for(int nu = 0; nu < mu; nu++) { - tmp = a_comp[mu] + Cshift(a_comp[nu], mu, 1) - Cshift(a_comp[mu], nu, 1) - a_comp[nu]; - plaqA = plaqA + cos(tmp); - - tmp = a_comp[mu] + Cshift(a_comp[mu], mu, 1) - + Cshift(a_comp[nu], mu, 2) + Cshift(Cshift(a_comp[nu], mu, 2), nu, 1) - - Cshift(Cshift(a_comp[mu], nu, 2), mu, 1) - Cshift(a_comp[mu], nu, 2) - - Cshift(a_comp[nu], nu, 1) - a_comp[nu]; - wlA = wlA + cos(tmp); - } - } - - Real vol = grid.gSites(); - Real faces = (1.0 * Nd * (Nd - 1)) / 2.0; - - Complex avgPlaqA = sum(trace(plaqA)); - avgPlaqA = avgPlaqA / vol / faces; - - Complex avgWlA = sum(trace(wlA)); - avgWlA = avgWlA / vol / faces; - - TComplex tplaqsite; - LatticeComplex plaqtrace = trace(plaqA); - std::vector site0 = {0,0,0,0}; - peekSite(tplaqsite, plaqtrace, site0); - Complex plaqsite = TensorRemove(tplaqsite); - - // Calculate plaquette from exponentiated photon field - avgPlaqAexp = NewWilsonLoops::avgPlaquette(expA); - avgWl2x2Aexp = NewWilsonLoops::avgWilsonLoop(expA, 2, 2); - avgWl2x2Aexp_time = NewWilsonLoops::avgTimelikeWilsonLoop(expA, 2, 2); - avgWl2x2Aexp_space = NewWilsonLoops::avgSpatialWilsonLoop(expA, 2, 2); - - avgPlaqAexp = avgPlaqAexp*3; - avgWl2x2Aexp = avgWl2x2Aexp*3; - avgWl2x2Aexp_time = avgWl2x2Aexp_time*3; - avgWl2x2Aexp_space = avgWl2x2Aexp_space*3; - - LOG(Message) << "Plaquette average (from A): " << avgPlaqA << std::endl; - LOG(Message) << "Plaquette average (from exp(A)): " << avgPlaqAexp << std::endl; - LOG(Message) << "2x2 Wilson Loop average (from A): " << avgWlA << std::endl; - LOG(Message) << "2x2 Wilson Loop average (from exp(A)): " << avgWl2x2Aexp << std::endl; - LOG(Message) << "2x2 Wilson Loop timelike average (from exp(A)): " << avgWl2x2Aexp_time << std::endl; - LOG(Message) << "2x2 Wilson Loop spatial average (from exp(A)): " << avgWl2x2Aexp_space << std::endl; - LOG(Message) << "Plaquette (one site): " << plaqsite / faces << std::endl; - // epilogue LOG(Message) << "Grid is finalizing now" << std::endl; Grid_finalize(); From 739c2308b5ce9a9464dbbd9057dbe49f6b04cf59 Mon Sep 17 00:00:00 2001 From: James Harrison Date: Tue, 15 Nov 2016 13:07:52 +0000 Subject: [PATCH 14/14] Set imaginary part of stochastic QED field to zero using real() instead of conjugate(). --- lib/qcd/action/gauge/Photon.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/qcd/action/gauge/Photon.h b/lib/qcd/action/gauge/Photon.h index ca0a8d40..b6c1b76f 100644 --- a/lib/qcd/action/gauge/Photon.h +++ b/lib/qcd/action/gauge/Photon.h @@ -173,7 +173,7 @@ namespace QCD{ } fft.FFT_all_dim(out, aTilde, FFT::backward); - out = 0.5*(out + conjugate(out)); + out = real(out); } // template // void Photon::FeynmanGaugeMomentumSpacePropagator_L(GaugeField &out,