From db86cdd7bd9ac6e6d355fc6239daf5939ed3cda1 Mon Sep 17 00:00:00 2001 From: fionnoh Date: Tue, 10 Jul 2018 13:30:45 +0100 Subject: [PATCH] Possible trash commit --- extras/Hadrons/AllToAllVectors.hpp | 1 - .../Modules/MContraction/MesonFieldGamma.hpp | 5 +- lib/lattice/Lattice_reduction.h | 121 +++++++++++++++++- lib/tensors/Tensor_inner.h | 26 ++++ 4 files changed, 148 insertions(+), 5 deletions(-) diff --git a/extras/Hadrons/AllToAllVectors.hpp b/extras/Hadrons/AllToAllVectors.hpp index d9f3c281..a514beac 100644 --- a/extras/Hadrons/AllToAllVectors.hpp +++ b/extras/Hadrons/AllToAllVectors.hpp @@ -102,7 +102,6 @@ class A2AModesSchurDiagTwo void low_mode_v(Matrix &action, const Field &evec, const RealD &eval, Field &vout_5d, Field &vout_4d) { - GridBase *grid = action.RedBlackGrid(); Field src_o(grid); Field sol_e(grid); diff --git a/extras/Hadrons/Modules/MContraction/MesonFieldGamma.hpp b/extras/Hadrons/Modules/MContraction/MesonFieldGamma.hpp index 59c90c08..0b81f607 100644 --- a/extras/Hadrons/Modules/MContraction/MesonFieldGamma.hpp +++ b/extras/Hadrons/Modules/MContraction/MesonFieldGamma.hpp @@ -39,7 +39,8 @@ class TMesonFieldGamma : public Module public: GRID_SERIALIZABLE_CLASS_MEMBERS(Result, Gamma::Algebra, gamma, - std::vector>>, MesonField); + std::vector>>, MesonField, + ComplexD, last); }; public: @@ -178,7 +179,7 @@ void TMesonFieldGamma::execute(void) LOG(Message) << "MF for i = " << i << " of " << N << std::endl; } } - + result[0].last = MesonField_ij[7]; saveResult(par().output, "meson", result); } diff --git a/lib/lattice/Lattice_reduction.h b/lib/lattice/Lattice_reduction.h index 3be4b6cb..12bfebed 100644 --- a/lib/lattice/Lattice_reduction.h +++ b/lib/lattice/Lattice_reduction.h @@ -39,8 +39,9 @@ template inline RealD norm2(const Lattice &arg){ // Double inner product template -inline ComplexD innerProduct(const Lattice &left,const Lattice &right) +inline ComplexD innerProduct(const Lattice &left,const Lattice &right) { + std::cout << GridLogMessage << "Start alloc innerProduct" << std::endl; typedef typename vobj::scalar_type scalar_type; typedef typename vobj::vector_typeD vector_type; GridBase *grid = left._grid; @@ -49,6 +50,8 @@ inline ComplexD innerProduct(const Lattice &left,const Lattice &righ ComplexD inner; Vector sumarray(grid->SumArraySize()*pad); + std::cout << GridLogMessage << "End alloc innerProduct" << std::endl; + std::cout << GridLogMessage << "Start parallel for innerProduct" << std::endl; parallel_for(int thr=0;thrSumArraySize();thr++){ int nwork, mywork, myoff; GridThread::GetWork(left._grid->oSites(),thr,mywork,myoff); @@ -62,13 +65,18 @@ inline ComplexD innerProduct(const Lattice &left,const Lattice &righ ComplexD tmp = Reduce(TensorRemove(vinner)) ; vstream(sumarray[thr*pad],tmp); } - + std::cout << GridLogMessage << "End parallel for innerProduct" << std::endl; + + std::cout << GridLogMessage << "Start inner sum innerProduct" << std::endl; inner=0.0; for(int i=0;iSumArraySize();i++){ inner = inner+sumarray[i*pad]; } right._grid->GlobalSum(inner); return inner; + std::cout << GridLogMessage << "End inner sum innerProduct" << std::endl; + + std::cout << GridLogMessage << "End innerProduct" << std::endl; } ///////////////////////// @@ -274,6 +282,115 @@ template inline void sliceSum(const Lattice &Data,std::vector< } } +template +static void mySliceInnerProductVector( std::vector & result, const Lattice &lhs,const Lattice &rhs,int orthogdim) +{ + std::cout << GridLogMessage << "Start mySsliceInnerProductVector" << std::endl; + + typedef typename vobj::scalar_type scalar_type; + std::vector lsSum; + localSliceInnerProductVector(result, lhs, rhs, lsSum, orthogdim); + globalSliceInnerProductVector(result, lhs, lsSum, orthogdim); + std::cout << GridLogMessage << "End mySliceInnerProductVector" << std::endl; +} + +template +static void localSliceInnerProductVector(std::vector &result, const Lattice &lhs, const Lattice &rhs, std::vector &lsSum, int orthogdim) +{ + std::cout << GridLogMessage << "Start prep" << std::endl; + typedef typename vobj::vector_type vector_type; + typedef typename vobj::scalar_type scalar_type; + GridBase *grid = lhs._grid; + assert(grid!=NULL); + conformable(grid,rhs._grid); + + const int Nd = grid->_ndimension; + const int Nsimd = grid->Nsimd(); + + assert(orthogdim >= 0); + assert(orthogdim < Nd); + + int fd=grid->_fdimensions[orthogdim]; + int ld=grid->_ldimensions[orthogdim]; + int rd=grid->_rdimensions[orthogdim]; + std::cout << GridLogMessage << "Start alloc" << std::endl; + + std::vector > lvSum(rd); // will locally sum vectors first + lsSum.resize(ld,scalar_type(0.0)); // sum across these down to scalars + std::vector> extracted(Nsimd); // splitting the SIMD + std::cout << GridLogMessage << "End alloc" << std::endl; + + result.resize(fd); // And then global sum to return the same vector to every node for IO to file + for(int r=0;r_slice_nblock[orthogdim]; + int e2= grid->_slice_block [orthogdim]; + int stride=grid->_slice_stride[orthogdim]; + std::cout << GridLogMessage << "End prep" << std::endl; + std::cout << GridLogMessage << "Start parallel inner product, _rd = " << rd << std::endl; + vector_type vv; + parallel_for(int r=0;r_ostride[orthogdim]; // base offset for start of plane + + for(int n=0;n icoor(Nd); + for(int rt=0;rt temp; + temp._internal = lvSum[rt]; + extract(temp,extracted); + + for(int idx=0;idxiCoorFromIindex(icoor,idx); + + int ldx =rt+icoor[orthogdim]*rd; + + lsSum[ldx]=lsSum[ldx]+extracted[idx]._internal; + + } + } + std::cout << GridLogMessage << "End sum over simd lanes" << std::endl; +} +template +static void globalSliceInnerProductVector(std::vector &result, const Lattice &lhs, std::vector &lsSum, int orthogdim) +{ + typedef typename vobj::scalar_type scalar_type; + GridBase *grid = lhs._grid; + int fd = result.size(); + int ld = lsSum.size(); + // sum over nodes. + std::vector gsum; + gsum.resize(fd, scalar_type(0.0)); + std::cout << GridLogMessage << "Start of gsum[t] creation:" << std::endl; + for(int t=0;t_processor_coor[orthogdim] ) { + gsum[t]=lsSum[lt]; + } + } + std::cout << GridLogMessage << "End of gsum[t] creation:" << std::endl; + std::cout << GridLogMessage << "Start of GlobalSumVector:" << std::endl; + grid->GlobalSumVector(&gsum[0], fd); + std::cout << GridLogMessage << "End of GlobalSumVector:" << std::endl; + + result = gsum; +} template static void sliceInnerProductVector( std::vector & result, const Lattice &lhs,const Lattice &rhs,int orthogdim) { diff --git a/lib/tensors/Tensor_inner.h b/lib/tensors/Tensor_inner.h index 46185652..b7a8417d 100644 --- a/lib/tensors/Tensor_inner.h +++ b/lib/tensors/Tensor_inner.h @@ -106,6 +106,7 @@ inline vRealD innerProductD(const vRealF &l,const vRealF &r){ typedef decltype(innerProduct(lhs._internal[0],rhs._internal[0])) ret_t; iScalar ret; ret=zero; + // std::cout << GridLogMessage << "innerProduct iVector" << std::endl; for(int c1=0;c1 ret; + // std::cout << GridLogMessage << "innerProduct iScalar" << std::endl; + ret._internal = innerProduct(lhs._internal,rhs._internal); return ret; } + template inline + auto myInnerProduct (const iVector& lhs,const iVector& rhs) -> iScalar + { + typedef decltype(innerProduct(lhs._internal[0],rhs._internal[0])) ret_t; + iScalar ret; + ret=zero; + std::cout << GridLogMessage << "myInnerProduct iVector, N = " << N << std::endl; + for(int c1=0;c1 inline + auto myInnerProduct (const iScalar& lhs,const iScalar& rhs) -> iScalar + { + typedef decltype(innerProduct(lhs._internal,rhs._internal)) ret_t; + iScalar ret; + std::cout << GridLogMessage << "myInnerProduct iScalar" << std::endl; + + ret._internal = myInnerProduct(lhs._internal,rhs._internal); + return ret; + } } #endif