mirror of
				https://github.com/paboyle/Grid.git
				synced 2025-11-04 05:54:32 +00:00 
			
		
		
		
	Format, NAMESPACE
This commit is contained in:
		@@ -1,4 +1,4 @@
 | 
			
		||||
    /*************************************************************************************
 | 
			
		||||
/*************************************************************************************
 | 
			
		||||
 | 
			
		||||
    Grid physics library, www.github.com/paboyle/Grid 
 | 
			
		||||
 | 
			
		||||
@@ -23,8 +23,8 @@ Author: Peter Boyle <paboyle@ph.ed.ac.uk>
 | 
			
		||||
    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 */
 | 
			
		||||
*************************************************************************************/
 | 
			
		||||
/*  END LEGAL */
 | 
			
		||||
#ifndef GRID_CARTESIAN_H
 | 
			
		||||
#define GRID_CARTESIAN_H
 | 
			
		||||
 | 
			
		||||
 
 | 
			
		||||
@@ -1,4 +1,4 @@
 | 
			
		||||
    /*************************************************************************************
 | 
			
		||||
/*************************************************************************************
 | 
			
		||||
 | 
			
		||||
    Grid physics library, www.github.com/paboyle/Grid 
 | 
			
		||||
 | 
			
		||||
@@ -25,265 +25,264 @@
 | 
			
		||||
    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 */
 | 
			
		||||
*************************************************************************************/
 | 
			
		||||
/*  END LEGAL */
 | 
			
		||||
#ifndef GRID_CARTESIAN_BASE_H
 | 
			
		||||
#define GRID_CARTESIAN_BASE_H
 | 
			
		||||
 | 
			
		||||
 | 
			
		||||
namespace Grid{
 | 
			
		||||
NAMESPACE_BEGIN(Grid);
 | 
			
		||||
 | 
			
		||||
  //////////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Commicator provides information on the processor grid
 | 
			
		||||
  //////////////////////////////////////////////////////////////////////
 | 
			
		||||
  //    unsigned long _ndimension;
 | 
			
		||||
  //    std::vector<int> _processors; // processor grid
 | 
			
		||||
  //    int              _processor;  // linear processor rank
 | 
			
		||||
  //    std::vector<int> _processor_coor;  // linear processor rank
 | 
			
		||||
  //////////////////////////////////////////////////////////////////////
 | 
			
		||||
  class GridBase : public CartesianCommunicator , public GridThread {
 | 
			
		||||
//////////////////////////////////////////////////////////////////////
 | 
			
		||||
// Commicator provides information on the processor grid
 | 
			
		||||
//////////////////////////////////////////////////////////////////////
 | 
			
		||||
//    unsigned long _ndimension;
 | 
			
		||||
//    std::vector<int> _processors; // processor grid
 | 
			
		||||
//    int              _processor;  // linear processor rank
 | 
			
		||||
//    std::vector<int> _processor_coor;  // linear processor rank
 | 
			
		||||
//////////////////////////////////////////////////////////////////////
 | 
			
		||||
class GridBase : public CartesianCommunicator , public GridThread {
 | 
			
		||||
 | 
			
		||||
public:
 | 
			
		||||
    int dummy;
 | 
			
		||||
    // Give Lattice access
 | 
			
		||||
    template<class object> friend class Lattice;
 | 
			
		||||
  int dummy;
 | 
			
		||||
  // Give Lattice access
 | 
			
		||||
  template<class object> friend class Lattice;
 | 
			
		||||
 | 
			
		||||
    GridBase(const std::vector<int> & processor_grid) : CartesianCommunicator(processor_grid) {};
 | 
			
		||||
    GridBase(const std::vector<int> & processor_grid,
 | 
			
		||||
	     const CartesianCommunicator &parent,
 | 
			
		||||
	     int &split_rank) 
 | 
			
		||||
      : CartesianCommunicator(processor_grid,parent,split_rank) {};
 | 
			
		||||
    GridBase(const std::vector<int> & processor_grid,
 | 
			
		||||
	     const CartesianCommunicator &parent) 
 | 
			
		||||
      : CartesianCommunicator(processor_grid,parent,dummy) {};
 | 
			
		||||
  GridBase(const std::vector<int> & processor_grid) : CartesianCommunicator(processor_grid) {};
 | 
			
		||||
  GridBase(const std::vector<int> & processor_grid,
 | 
			
		||||
	   const CartesianCommunicator &parent,
 | 
			
		||||
	   int &split_rank) 
 | 
			
		||||
    : CartesianCommunicator(processor_grid,parent,split_rank) {};
 | 
			
		||||
  GridBase(const std::vector<int> & processor_grid,
 | 
			
		||||
	   const CartesianCommunicator &parent) 
 | 
			
		||||
    : CartesianCommunicator(processor_grid,parent,dummy) {};
 | 
			
		||||
 | 
			
		||||
    virtual ~GridBase() = default;
 | 
			
		||||
  virtual ~GridBase() = default;
 | 
			
		||||
 | 
			
		||||
    // Physics Grid information.
 | 
			
		||||
    std::vector<int> _simd_layout;// Which dimensions get relayed out over simd lanes.
 | 
			
		||||
    std::vector<int> _fdimensions;// (full) Global dimensions of array prior to cb removal
 | 
			
		||||
    std::vector<int> _gdimensions;// Global dimensions of array after cb removal
 | 
			
		||||
    std::vector<int> _ldimensions;// local dimensions of array with processor images removed
 | 
			
		||||
    std::vector<int> _rdimensions;// Reduced local dimensions with simd lane images and processor images removed 
 | 
			
		||||
    std::vector<int> _ostride;    // Outer stride for each dimension
 | 
			
		||||
    std::vector<int> _istride;    // Inner stride i.e. within simd lane
 | 
			
		||||
    int _osites;                  // _isites*_osites = product(dimensions).
 | 
			
		||||
    int _isites;
 | 
			
		||||
    int _fsites;                  // _isites*_osites = product(dimensions).
 | 
			
		||||
    int _gsites;
 | 
			
		||||
    std::vector<int> _slice_block;// subslice information
 | 
			
		||||
    std::vector<int> _slice_stride;
 | 
			
		||||
    std::vector<int> _slice_nblock;
 | 
			
		||||
  // Physics Grid information.
 | 
			
		||||
  std::vector<int> _simd_layout;// Which dimensions get relayed out over simd lanes.
 | 
			
		||||
  std::vector<int> _fdimensions;// (full) Global dimensions of array prior to cb removal
 | 
			
		||||
  std::vector<int> _gdimensions;// Global dimensions of array after cb removal
 | 
			
		||||
  std::vector<int> _ldimensions;// local dimensions of array with processor images removed
 | 
			
		||||
  std::vector<int> _rdimensions;// Reduced local dimensions with simd lane images and processor images removed 
 | 
			
		||||
  std::vector<int> _ostride;    // Outer stride for each dimension
 | 
			
		||||
  std::vector<int> _istride;    // Inner stride i.e. within simd lane
 | 
			
		||||
  int _osites;                  // _isites*_osites = product(dimensions).
 | 
			
		||||
  int _isites;
 | 
			
		||||
  int _fsites;                  // _isites*_osites = product(dimensions).
 | 
			
		||||
  int _gsites;
 | 
			
		||||
  std::vector<int> _slice_block;// subslice information
 | 
			
		||||
  std::vector<int> _slice_stride;
 | 
			
		||||
  std::vector<int> _slice_nblock;
 | 
			
		||||
 | 
			
		||||
    std::vector<int> _lstart;     // local start of array in gcoors _processor_coor[d]*_ldimensions[d]
 | 
			
		||||
    std::vector<int> _lend  ;     // local end of array in gcoors   _processor_coor[d]*_ldimensions[d]+_ldimensions_[d]-1
 | 
			
		||||
  std::vector<int> _lstart;     // local start of array in gcoors _processor_coor[d]*_ldimensions[d]
 | 
			
		||||
  std::vector<int> _lend  ;     // local end of array in gcoors   _processor_coor[d]*_ldimensions[d]+_ldimensions_[d]-1
 | 
			
		||||
 | 
			
		||||
public:
 | 
			
		||||
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Checkerboarding interface is virtual and overridden by 
 | 
			
		||||
    // GridCartesian / GridRedBlackCartesian
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
    virtual int CheckerBoarded(int dim)=0;
 | 
			
		||||
    virtual int CheckerBoard(const std::vector<int> &site)=0;
 | 
			
		||||
    virtual int CheckerBoardDestination(int source_cb,int shift,int dim)=0;
 | 
			
		||||
    virtual int CheckerBoardShift(int source_cb,int dim,int shift,int osite)=0;
 | 
			
		||||
    virtual int CheckerBoardShiftForCB(int source_cb,int dim,int shift,int cb)=0;
 | 
			
		||||
    virtual int CheckerBoardFromOindex (int Oindex)=0;
 | 
			
		||||
    virtual int CheckerBoardFromOindexTable (int Oindex)=0;
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Checkerboarding interface is virtual and overridden by 
 | 
			
		||||
  // GridCartesian / GridRedBlackCartesian
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
  virtual int CheckerBoarded(int dim)=0;
 | 
			
		||||
  virtual int CheckerBoard(const std::vector<int> &site)=0;
 | 
			
		||||
  virtual int CheckerBoardDestination(int source_cb,int shift,int dim)=0;
 | 
			
		||||
  virtual int CheckerBoardShift(int source_cb,int dim,int shift,int osite)=0;
 | 
			
		||||
  virtual int CheckerBoardShiftForCB(int source_cb,int dim,int shift,int cb)=0;
 | 
			
		||||
  virtual int CheckerBoardFromOindex (int Oindex)=0;
 | 
			
		||||
  virtual int CheckerBoardFromOindexTable (int Oindex)=0;
 | 
			
		||||
 | 
			
		||||
    //////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Local layout calculations
 | 
			
		||||
    //////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    // These routines are key. Subdivide the linearised cartesian index into
 | 
			
		||||
    //      "inner" index identifying which simd lane of object<vFcomplex> is associated with coord
 | 
			
		||||
    //      "outer" index identifying which element of _odata in class "Lattice" is associated with coord.
 | 
			
		||||
    //
 | 
			
		||||
    // Compared to, say, Blitz++ we simply need to store BOTH an inner stride and an outer
 | 
			
		||||
    // stride per dimension. The cost of evaluating the indexing information is doubled for an n-dimensional
 | 
			
		||||
    // coordinate. Note, however, for data parallel operations the "inner" indexing cost is not paid and all
 | 
			
		||||
    // lanes are operated upon simultaneously.
 | 
			
		||||
  //////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Local layout calculations
 | 
			
		||||
  //////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
  // These routines are key. Subdivide the linearised cartesian index into
 | 
			
		||||
  //      "inner" index identifying which simd lane of object<vFcomplex> is associated with coord
 | 
			
		||||
  //      "outer" index identifying which element of _odata in class "Lattice" is associated with coord.
 | 
			
		||||
  //
 | 
			
		||||
  // Compared to, say, Blitz++ we simply need to store BOTH an inner stride and an outer
 | 
			
		||||
  // stride per dimension. The cost of evaluating the indexing information is doubled for an n-dimensional
 | 
			
		||||
  // coordinate. Note, however, for data parallel operations the "inner" indexing cost is not paid and all
 | 
			
		||||
  // lanes are operated upon simultaneously.
 | 
			
		||||
  
 | 
			
		||||
    virtual int oIndex(std::vector<int> &coor)
 | 
			
		||||
    {
 | 
			
		||||
        int idx=0;
 | 
			
		||||
        // Works with either global or local coordinates
 | 
			
		||||
        for(int d=0;d<_ndimension;d++) idx+=_ostride[d]*(coor[d]%_rdimensions[d]);
 | 
			
		||||
        return idx;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int iIndex(std::vector<int> &lcoor)
 | 
			
		||||
    {
 | 
			
		||||
        int idx=0;
 | 
			
		||||
        for(int d=0;d<_ndimension;d++) idx+=_istride[d]*(lcoor[d]/_rdimensions[d]);
 | 
			
		||||
        return idx;
 | 
			
		||||
    }
 | 
			
		||||
    inline int oIndexReduced(std::vector<int> &ocoor)
 | 
			
		||||
    {
 | 
			
		||||
      int idx=0; 
 | 
			
		||||
      // ocoor is already reduced so can eliminate the modulo operation
 | 
			
		||||
      // for fast indexing and inline the routine
 | 
			
		||||
      for(int d=0;d<_ndimension;d++) idx+=_ostride[d]*ocoor[d];
 | 
			
		||||
      return idx;
 | 
			
		||||
    }
 | 
			
		||||
    inline void oCoorFromOindex (std::vector<int>& coor,int Oindex){
 | 
			
		||||
      Lexicographic::CoorFromIndex(coor,Oindex,_rdimensions);
 | 
			
		||||
    }
 | 
			
		||||
  virtual int oIndex(std::vector<int> &coor)
 | 
			
		||||
  {
 | 
			
		||||
    int idx=0;
 | 
			
		||||
    // Works with either global or local coordinates
 | 
			
		||||
    for(int d=0;d<_ndimension;d++) idx+=_ostride[d]*(coor[d]%_rdimensions[d]);
 | 
			
		||||
    return idx;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int iIndex(std::vector<int> &lcoor)
 | 
			
		||||
  {
 | 
			
		||||
    int idx=0;
 | 
			
		||||
    for(int d=0;d<_ndimension;d++) idx+=_istride[d]*(lcoor[d]/_rdimensions[d]);
 | 
			
		||||
    return idx;
 | 
			
		||||
  }
 | 
			
		||||
  inline int oIndexReduced(std::vector<int> &ocoor)
 | 
			
		||||
  {
 | 
			
		||||
    int idx=0; 
 | 
			
		||||
    // ocoor is already reduced so can eliminate the modulo operation
 | 
			
		||||
    // for fast indexing and inline the routine
 | 
			
		||||
    for(int d=0;d<_ndimension;d++) idx+=_ostride[d]*ocoor[d];
 | 
			
		||||
    return idx;
 | 
			
		||||
  }
 | 
			
		||||
  inline void oCoorFromOindex (std::vector<int>& coor,int Oindex){
 | 
			
		||||
    Lexicographic::CoorFromIndex(coor,Oindex,_rdimensions);
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
    inline void InOutCoorToLocalCoor (std::vector<int> &ocoor, std::vector<int> &icoor, std::vector<int> &lcoor) {
 | 
			
		||||
      lcoor.resize(_ndimension);
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
        lcoor[d] = ocoor[d] + _rdimensions[d] * icoor[d];
 | 
			
		||||
    }
 | 
			
		||||
  inline void InOutCoorToLocalCoor (std::vector<int> &ocoor, std::vector<int> &icoor, std::vector<int> &lcoor) {
 | 
			
		||||
    lcoor.resize(_ndimension);
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      lcoor[d] = ocoor[d] + _rdimensions[d] * icoor[d];
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
    //////////////////////////////////////////////////////////
 | 
			
		||||
    // SIMD lane addressing
 | 
			
		||||
    //////////////////////////////////////////////////////////
 | 
			
		||||
    inline void iCoorFromIindex(std::vector<int> &coor,int lane)
 | 
			
		||||
    {
 | 
			
		||||
      Lexicographic::CoorFromIndex(coor,lane,_simd_layout);
 | 
			
		||||
    }
 | 
			
		||||
  //////////////////////////////////////////////////////////
 | 
			
		||||
  // SIMD lane addressing
 | 
			
		||||
  //////////////////////////////////////////////////////////
 | 
			
		||||
  inline void iCoorFromIindex(std::vector<int> &coor,int lane)
 | 
			
		||||
  {
 | 
			
		||||
    Lexicographic::CoorFromIndex(coor,lane,_simd_layout);
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
    inline int PermuteDim(int dimension){
 | 
			
		||||
      return _simd_layout[dimension]>1;
 | 
			
		||||
    }
 | 
			
		||||
    inline int PermuteType(int dimension){
 | 
			
		||||
      int permute_type=0;
 | 
			
		||||
      //
 | 
			
		||||
      // FIXME:
 | 
			
		||||
      //
 | 
			
		||||
      // Best way to encode this would be to present a mask 
 | 
			
		||||
      // for which simd dimensions are rotated, and the rotation
 | 
			
		||||
      // size. If there is only one simd dimension rotated, this is just 
 | 
			
		||||
      // a permute. 
 | 
			
		||||
      //
 | 
			
		||||
      // Cases: PermuteType == 1,2,4,8
 | 
			
		||||
      // Distance should be either 0,1,2..
 | 
			
		||||
      //
 | 
			
		||||
      if ( _simd_layout[dimension] > 2 ) { 
 | 
			
		||||
        for(int d=0;d<_ndimension;d++){
 | 
			
		||||
          if ( d != dimension ) assert ( (_simd_layout[d]==1)  );
 | 
			
		||||
        }
 | 
			
		||||
        permute_type = RotateBit; // How to specify distance; this is not just direction.
 | 
			
		||||
        return permute_type;
 | 
			
		||||
      }
 | 
			
		||||
 | 
			
		||||
      for(int d=_ndimension-1;d>dimension;d--){
 | 
			
		||||
        if (_simd_layout[d]>1 ) permute_type++;
 | 
			
		||||
  inline int PermuteDim(int dimension){
 | 
			
		||||
    return _simd_layout[dimension]>1;
 | 
			
		||||
  }
 | 
			
		||||
  inline int PermuteType(int dimension){
 | 
			
		||||
    int permute_type=0;
 | 
			
		||||
    //
 | 
			
		||||
    // FIXME:
 | 
			
		||||
    //
 | 
			
		||||
    // Best way to encode this would be to present a mask 
 | 
			
		||||
    // for which simd dimensions are rotated, and the rotation
 | 
			
		||||
    // size. If there is only one simd dimension rotated, this is just 
 | 
			
		||||
    // a permute. 
 | 
			
		||||
    //
 | 
			
		||||
    // Cases: PermuteType == 1,2,4,8
 | 
			
		||||
    // Distance should be either 0,1,2..
 | 
			
		||||
    //
 | 
			
		||||
    if ( _simd_layout[dimension] > 2 ) { 
 | 
			
		||||
      for(int d=0;d<_ndimension;d++){
 | 
			
		||||
	if ( d != dimension ) assert ( (_simd_layout[d]==1)  );
 | 
			
		||||
      }
 | 
			
		||||
      permute_type = RotateBit; // How to specify distance; this is not just direction.
 | 
			
		||||
      return permute_type;
 | 
			
		||||
    }
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Array sizing queries
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
 | 
			
		||||
    inline int iSites(void) const { return _isites; };
 | 
			
		||||
    inline int Nsimd(void)  const { return _isites; };// Synonymous with iSites
 | 
			
		||||
    inline int oSites(void) const { return _osites; };
 | 
			
		||||
    inline int lSites(void) const { return _isites*_osites; }; 
 | 
			
		||||
    inline int gSites(void) const { return _isites*_osites*_Nprocessors; }; 
 | 
			
		||||
    inline int Nd    (void) const { return _ndimension;};
 | 
			
		||||
 | 
			
		||||
    inline const std::vector<int> LocalStarts(void)             { return _lstart;    };
 | 
			
		||||
    inline const std::vector<int> &FullDimensions(void)         { return _fdimensions;};
 | 
			
		||||
    inline const std::vector<int> &GlobalDimensions(void)       { return _gdimensions;};
 | 
			
		||||
    inline const std::vector<int> &LocalDimensions(void)        { return _ldimensions;};
 | 
			
		||||
    inline const std::vector<int> &VirtualLocalDimensions(void) { return _ldimensions;};
 | 
			
		||||
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Utility to print the full decomposition details 
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
 | 
			
		||||
    void show_decomposition(){
 | 
			
		||||
      std::cout << GridLogMessage << "\tFull Dimensions    : " << _fdimensions << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tSIMD layout        : " << _simd_layout << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tGlobal Dimensions  : " << _gdimensions << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tLocal Dimensions   : " << _ldimensions << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tReduced Dimensions : " << _rdimensions << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tOuter strides      : " << _ostride << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tInner strides      : " << _istride << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tiSites             : " << _isites << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\toSites             : " << _osites << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tlSites             : " << lSites() << std::endl;        
 | 
			
		||||
      std::cout << GridLogMessage << "\tgSites             : " << gSites() << std::endl;
 | 
			
		||||
      std::cout << GridLogMessage << "\tNd                 : " << _ndimension << std::endl;             
 | 
			
		||||
    } 
 | 
			
		||||
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Global addressing
 | 
			
		||||
    ////////////////////////////////////////////////////////////////
 | 
			
		||||
    void GlobalIndexToGlobalCoor(int gidx,std::vector<int> &gcoor){
 | 
			
		||||
      assert(gidx< gSites());
 | 
			
		||||
      Lexicographic::CoorFromIndex(gcoor,gidx,_gdimensions);
 | 
			
		||||
    for(int d=_ndimension-1;d>dimension;d--){
 | 
			
		||||
      if (_simd_layout[d]>1 ) permute_type++;
 | 
			
		||||
    }
 | 
			
		||||
    void LocalIndexToLocalCoor(int lidx,std::vector<int> &lcoor){
 | 
			
		||||
      assert(lidx<lSites());
 | 
			
		||||
      Lexicographic::CoorFromIndex(lcoor,lidx,_ldimensions);
 | 
			
		||||
    return permute_type;
 | 
			
		||||
  }
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Array sizing queries
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
 | 
			
		||||
  inline int iSites(void) const { return _isites; };
 | 
			
		||||
  inline int Nsimd(void)  const { return _isites; };// Synonymous with iSites
 | 
			
		||||
  inline int oSites(void) const { return _osites; };
 | 
			
		||||
  inline int lSites(void) const { return _isites*_osites; }; 
 | 
			
		||||
  inline int gSites(void) const { return _isites*_osites*_Nprocessors; }; 
 | 
			
		||||
  inline int Nd    (void) const { return _ndimension;};
 | 
			
		||||
 | 
			
		||||
  inline const std::vector<int> LocalStarts(void)             { return _lstart;    };
 | 
			
		||||
  inline const std::vector<int> &FullDimensions(void)         { return _fdimensions;};
 | 
			
		||||
  inline const std::vector<int> &GlobalDimensions(void)       { return _gdimensions;};
 | 
			
		||||
  inline const std::vector<int> &LocalDimensions(void)        { return _ldimensions;};
 | 
			
		||||
  inline const std::vector<int> &VirtualLocalDimensions(void) { return _ldimensions;};
 | 
			
		||||
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Utility to print the full decomposition details 
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
 | 
			
		||||
  void show_decomposition(){
 | 
			
		||||
    std::cout << GridLogMessage << "\tFull Dimensions    : " << _fdimensions << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tSIMD layout        : " << _simd_layout << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tGlobal Dimensions  : " << _gdimensions << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tLocal Dimensions   : " << _ldimensions << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tReduced Dimensions : " << _rdimensions << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tOuter strides      : " << _ostride << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tInner strides      : " << _istride << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tiSites             : " << _isites << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\toSites             : " << _osites << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tlSites             : " << lSites() << std::endl;        
 | 
			
		||||
    std::cout << GridLogMessage << "\tgSites             : " << gSites() << std::endl;
 | 
			
		||||
    std::cout << GridLogMessage << "\tNd                 : " << _ndimension << std::endl;             
 | 
			
		||||
  } 
 | 
			
		||||
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Global addressing
 | 
			
		||||
  ////////////////////////////////////////////////////////////////
 | 
			
		||||
  void GlobalIndexToGlobalCoor(int gidx,std::vector<int> &gcoor){
 | 
			
		||||
    assert(gidx< gSites());
 | 
			
		||||
    Lexicographic::CoorFromIndex(gcoor,gidx,_gdimensions);
 | 
			
		||||
  }
 | 
			
		||||
  void LocalIndexToLocalCoor(int lidx,std::vector<int> &lcoor){
 | 
			
		||||
    assert(lidx<lSites());
 | 
			
		||||
    Lexicographic::CoorFromIndex(lcoor,lidx,_ldimensions);
 | 
			
		||||
  }
 | 
			
		||||
  void GlobalCoorToGlobalIndex(const std::vector<int> & gcoor,int & gidx){
 | 
			
		||||
    gidx=0;
 | 
			
		||||
    int mult=1;
 | 
			
		||||
    for(int mu=0;mu<_ndimension;mu++) {
 | 
			
		||||
      gidx+=mult*gcoor[mu];
 | 
			
		||||
      mult*=_gdimensions[mu];
 | 
			
		||||
    }
 | 
			
		||||
    void GlobalCoorToGlobalIndex(const std::vector<int> & gcoor,int & gidx){
 | 
			
		||||
      gidx=0;
 | 
			
		||||
      int mult=1;
 | 
			
		||||
      for(int mu=0;mu<_ndimension;mu++) {
 | 
			
		||||
        gidx+=mult*gcoor[mu];
 | 
			
		||||
        mult*=_gdimensions[mu];
 | 
			
		||||
      }
 | 
			
		||||
  }
 | 
			
		||||
  void GlobalCoorToProcessorCoorLocalCoor(std::vector<int> &pcoor,std::vector<int> &lcoor,const std::vector<int> &gcoor)
 | 
			
		||||
  {
 | 
			
		||||
    pcoor.resize(_ndimension);
 | 
			
		||||
    lcoor.resize(_ndimension);
 | 
			
		||||
    for(int mu=0;mu<_ndimension;mu++){
 | 
			
		||||
      int _fld  = _fdimensions[mu]/_processors[mu];
 | 
			
		||||
      pcoor[mu] = gcoor[mu]/_fld;
 | 
			
		||||
      lcoor[mu] = gcoor[mu]%_fld;
 | 
			
		||||
    }
 | 
			
		||||
    void GlobalCoorToProcessorCoorLocalCoor(std::vector<int> &pcoor,std::vector<int> &lcoor,const std::vector<int> &gcoor)
 | 
			
		||||
    {
 | 
			
		||||
      pcoor.resize(_ndimension);
 | 
			
		||||
      lcoor.resize(_ndimension);
 | 
			
		||||
      for(int mu=0;mu<_ndimension;mu++){
 | 
			
		||||
        int _fld  = _fdimensions[mu]/_processors[mu];
 | 
			
		||||
        pcoor[mu] = gcoor[mu]/_fld;
 | 
			
		||||
        lcoor[mu] = gcoor[mu]%_fld;
 | 
			
		||||
      }
 | 
			
		||||
    }
 | 
			
		||||
    void GlobalCoorToRankIndex(int &rank, int &o_idx, int &i_idx ,const std::vector<int> &gcoor)
 | 
			
		||||
    {
 | 
			
		||||
      std::vector<int> pcoor;
 | 
			
		||||
      std::vector<int> lcoor;
 | 
			
		||||
      GlobalCoorToProcessorCoorLocalCoor(pcoor,lcoor,gcoor);
 | 
			
		||||
      rank = RankFromProcessorCoor(pcoor);
 | 
			
		||||
      /*
 | 
			
		||||
  }
 | 
			
		||||
  void GlobalCoorToRankIndex(int &rank, int &o_idx, int &i_idx ,const std::vector<int> &gcoor)
 | 
			
		||||
  {
 | 
			
		||||
    std::vector<int> pcoor;
 | 
			
		||||
    std::vector<int> lcoor;
 | 
			
		||||
    GlobalCoorToProcessorCoorLocalCoor(pcoor,lcoor,gcoor);
 | 
			
		||||
    rank = RankFromProcessorCoor(pcoor);
 | 
			
		||||
    /*
 | 
			
		||||
      std::vector<int> cblcoor(lcoor);
 | 
			
		||||
      for(int d=0;d<cblcoor.size();d++){
 | 
			
		||||
        if( this->CheckerBoarded(d) ) {
 | 
			
		||||
          cblcoor[d] = lcoor[d]/2;
 | 
			
		||||
        }
 | 
			
		||||
      if( this->CheckerBoarded(d) ) {
 | 
			
		||||
      cblcoor[d] = lcoor[d]/2;
 | 
			
		||||
      }
 | 
			
		||||
      */
 | 
			
		||||
      i_idx= iIndex(lcoor);
 | 
			
		||||
      o_idx= oIndex(lcoor);
 | 
			
		||||
    }
 | 
			
		||||
      }
 | 
			
		||||
    */
 | 
			
		||||
    i_idx= iIndex(lcoor);
 | 
			
		||||
    o_idx= oIndex(lcoor);
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
    void RankIndexToGlobalCoor(int rank, int o_idx, int i_idx , std::vector<int> &gcoor)
 | 
			
		||||
    {
 | 
			
		||||
      gcoor.resize(_ndimension);
 | 
			
		||||
      std::vector<int> coor(_ndimension);
 | 
			
		||||
  void RankIndexToGlobalCoor(int rank, int o_idx, int i_idx , std::vector<int> &gcoor)
 | 
			
		||||
  {
 | 
			
		||||
    gcoor.resize(_ndimension);
 | 
			
		||||
    std::vector<int> coor(_ndimension);
 | 
			
		||||
 | 
			
		||||
      ProcessorCoorFromRank(rank,coor);
 | 
			
		||||
      for(int mu=0;mu<_ndimension;mu++) gcoor[mu] = _ldimensions[mu]*coor[mu];
 | 
			
		||||
    ProcessorCoorFromRank(rank,coor);
 | 
			
		||||
    for(int mu=0;mu<_ndimension;mu++) gcoor[mu] = _ldimensions[mu]*coor[mu];
 | 
			
		||||
 | 
			
		||||
      iCoorFromIindex(coor,i_idx);
 | 
			
		||||
      for(int mu=0;mu<_ndimension;mu++) gcoor[mu] += _rdimensions[mu]*coor[mu];
 | 
			
		||||
    iCoorFromIindex(coor,i_idx);
 | 
			
		||||
    for(int mu=0;mu<_ndimension;mu++) gcoor[mu] += _rdimensions[mu]*coor[mu];
 | 
			
		||||
 | 
			
		||||
      oCoorFromOindex (coor,o_idx);
 | 
			
		||||
      for(int mu=0;mu<_ndimension;mu++) gcoor[mu] += coor[mu];
 | 
			
		||||
    oCoorFromOindex (coor,o_idx);
 | 
			
		||||
    for(int mu=0;mu<_ndimension;mu++) gcoor[mu] += coor[mu];
 | 
			
		||||
      
 | 
			
		||||
  }
 | 
			
		||||
  void RankIndexCbToFullGlobalCoor(int rank, int o_idx, int i_idx, int cb,std::vector<int> &fcoor)
 | 
			
		||||
  {
 | 
			
		||||
    RankIndexToGlobalCoor(rank,o_idx,i_idx ,fcoor);
 | 
			
		||||
    if(CheckerBoarded(0)){
 | 
			
		||||
      fcoor[0] = fcoor[0]*2+cb;
 | 
			
		||||
    }
 | 
			
		||||
    void RankIndexCbToFullGlobalCoor(int rank, int o_idx, int i_idx, int cb,std::vector<int> &fcoor)
 | 
			
		||||
    {
 | 
			
		||||
      RankIndexToGlobalCoor(rank,o_idx,i_idx ,fcoor);
 | 
			
		||||
      if(CheckerBoarded(0)){
 | 
			
		||||
        fcoor[0] = fcoor[0]*2+cb;
 | 
			
		||||
      }
 | 
			
		||||
    }
 | 
			
		||||
    void ProcessorCoorLocalCoorToGlobalCoor(std::vector<int> &Pcoor,std::vector<int> &Lcoor,std::vector<int> &gcoor)
 | 
			
		||||
    {
 | 
			
		||||
      gcoor.resize(_ndimension);
 | 
			
		||||
      for(int mu=0;mu<_ndimension;mu++) gcoor[mu] = Pcoor[mu]*_ldimensions[mu]+Lcoor[mu];
 | 
			
		||||
    }
 | 
			
		||||
  }
 | 
			
		||||
  void ProcessorCoorLocalCoorToGlobalCoor(std::vector<int> &Pcoor,std::vector<int> &Lcoor,std::vector<int> &gcoor)
 | 
			
		||||
  {
 | 
			
		||||
    gcoor.resize(_ndimension);
 | 
			
		||||
    for(int mu=0;mu<_ndimension;mu++) gcoor[mu] = Pcoor[mu]*_ldimensions[mu]+Lcoor[mu];
 | 
			
		||||
  }
 | 
			
		||||
};
 | 
			
		||||
 | 
			
		||||
 | 
			
		||||
}
 | 
			
		||||
NAMESPACE_END(Grid);
 | 
			
		||||
#endif
 | 
			
		||||
 
 | 
			
		||||
@@ -1,4 +1,4 @@
 | 
			
		||||
    /*************************************************************************************
 | 
			
		||||
/*************************************************************************************
 | 
			
		||||
 | 
			
		||||
    Grid physics library, www.github.com/paboyle/Grid 
 | 
			
		||||
 | 
			
		||||
@@ -23,12 +23,12 @@ Author: Peter Boyle <paboyle@ph.ed.ac.uk>
 | 
			
		||||
    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 */
 | 
			
		||||
*************************************************************************************/
 | 
			
		||||
/*  END LEGAL */
 | 
			
		||||
#ifndef GRID_CARTESIAN_FULL_H
 | 
			
		||||
#define GRID_CARTESIAN_FULL_H
 | 
			
		||||
 | 
			
		||||
namespace Grid{
 | 
			
		||||
NAMESPACE_BEGIN(Grid);
 | 
			
		||||
    
 | 
			
		||||
/////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
// Grid Support.
 | 
			
		||||
@@ -38,81 +38,81 @@ namespace Grid{
 | 
			
		||||
class GridCartesian: public GridBase {
 | 
			
		||||
 | 
			
		||||
public:
 | 
			
		||||
    int dummy;
 | 
			
		||||
    virtual int  CheckerBoardFromOindexTable (int Oindex) {
 | 
			
		||||
      return 0;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int  CheckerBoardFromOindex (int Oindex)
 | 
			
		||||
    {
 | 
			
		||||
      return 0;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoarded(int dim){
 | 
			
		||||
      return 0;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoard(const std::vector<int> &site){
 | 
			
		||||
        return 0;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoardDestination(int cb,int shift,int dim){
 | 
			
		||||
        return 0;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoardShiftForCB(int source_cb,int dim,int shift, int ocb){
 | 
			
		||||
      return shift;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoardShift(int source_cb,int dim,int shift, int osite){
 | 
			
		||||
      return shift;
 | 
			
		||||
    }
 | 
			
		||||
    /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Constructor takes a parent grid and possibly subdivides communicator.
 | 
			
		||||
    /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    GridCartesian(const std::vector<int> &dimensions,
 | 
			
		||||
		  const std::vector<int> &simd_layout,
 | 
			
		||||
		  const std::vector<int> &processor_grid,
 | 
			
		||||
		  const GridCartesian &parent) : GridBase(processor_grid,parent,dummy)
 | 
			
		||||
    {
 | 
			
		||||
      Init(dimensions,simd_layout,processor_grid);
 | 
			
		||||
    }
 | 
			
		||||
    GridCartesian(const std::vector<int> &dimensions,
 | 
			
		||||
		  const std::vector<int> &simd_layout,
 | 
			
		||||
		  const std::vector<int> &processor_grid,
 | 
			
		||||
		  const GridCartesian &parent,int &split_rank) : GridBase(processor_grid,parent,split_rank)
 | 
			
		||||
    {
 | 
			
		||||
      Init(dimensions,simd_layout,processor_grid);
 | 
			
		||||
    }
 | 
			
		||||
    /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    // Construct from comm world
 | 
			
		||||
    /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    GridCartesian(const std::vector<int> &dimensions,
 | 
			
		||||
		  const std::vector<int> &simd_layout,
 | 
			
		||||
		  const std::vector<int> &processor_grid) : GridBase(processor_grid)
 | 
			
		||||
    {
 | 
			
		||||
      Init(dimensions,simd_layout,processor_grid);
 | 
			
		||||
    }
 | 
			
		||||
  int dummy;
 | 
			
		||||
  virtual int  CheckerBoardFromOindexTable (int Oindex) {
 | 
			
		||||
    return 0;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int  CheckerBoardFromOindex (int Oindex)
 | 
			
		||||
  {
 | 
			
		||||
    return 0;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoarded(int dim){
 | 
			
		||||
    return 0;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoard(const std::vector<int> &site){
 | 
			
		||||
    return 0;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoardDestination(int cb,int shift,int dim){
 | 
			
		||||
    return 0;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoardShiftForCB(int source_cb,int dim,int shift, int ocb){
 | 
			
		||||
    return shift;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoardShift(int source_cb,int dim,int shift, int osite){
 | 
			
		||||
    return shift;
 | 
			
		||||
  }
 | 
			
		||||
  /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Constructor takes a parent grid and possibly subdivides communicator.
 | 
			
		||||
  /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
  GridCartesian(const std::vector<int> &dimensions,
 | 
			
		||||
		const std::vector<int> &simd_layout,
 | 
			
		||||
		const std::vector<int> &processor_grid,
 | 
			
		||||
		const GridCartesian &parent) : GridBase(processor_grid,parent,dummy)
 | 
			
		||||
  {
 | 
			
		||||
    Init(dimensions,simd_layout,processor_grid);
 | 
			
		||||
  }
 | 
			
		||||
  GridCartesian(const std::vector<int> &dimensions,
 | 
			
		||||
		const std::vector<int> &simd_layout,
 | 
			
		||||
		const std::vector<int> &processor_grid,
 | 
			
		||||
		const GridCartesian &parent,int &split_rank) : GridBase(processor_grid,parent,split_rank)
 | 
			
		||||
  {
 | 
			
		||||
    Init(dimensions,simd_layout,processor_grid);
 | 
			
		||||
  }
 | 
			
		||||
  /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
  // Construct from comm world
 | 
			
		||||
  /////////////////////////////////////////////////////////////////////////
 | 
			
		||||
  GridCartesian(const std::vector<int> &dimensions,
 | 
			
		||||
		const std::vector<int> &simd_layout,
 | 
			
		||||
		const std::vector<int> &processor_grid) : GridBase(processor_grid)
 | 
			
		||||
  {
 | 
			
		||||
    Init(dimensions,simd_layout,processor_grid);
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
    virtual ~GridCartesian() = default;
 | 
			
		||||
  virtual ~GridCartesian() = default;
 | 
			
		||||
 | 
			
		||||
    void Init(const std::vector<int> &dimensions,
 | 
			
		||||
	      const std::vector<int> &simd_layout,
 | 
			
		||||
	      const std::vector<int> &processor_grid)
 | 
			
		||||
    {
 | 
			
		||||
      ///////////////////////
 | 
			
		||||
      // Grid information
 | 
			
		||||
      ///////////////////////
 | 
			
		||||
      _ndimension = dimensions.size();
 | 
			
		||||
  void Init(const std::vector<int> &dimensions,
 | 
			
		||||
	    const std::vector<int> &simd_layout,
 | 
			
		||||
	    const std::vector<int> &processor_grid)
 | 
			
		||||
  {
 | 
			
		||||
    ///////////////////////
 | 
			
		||||
    // Grid information
 | 
			
		||||
    ///////////////////////
 | 
			
		||||
    _ndimension = dimensions.size();
 | 
			
		||||
 | 
			
		||||
      _fdimensions.resize(_ndimension);
 | 
			
		||||
      _gdimensions.resize(_ndimension);
 | 
			
		||||
      _ldimensions.resize(_ndimension);
 | 
			
		||||
      _rdimensions.resize(_ndimension);
 | 
			
		||||
      _simd_layout.resize(_ndimension);
 | 
			
		||||
      _lstart.resize(_ndimension);
 | 
			
		||||
      _lend.resize(_ndimension);
 | 
			
		||||
    _fdimensions.resize(_ndimension);
 | 
			
		||||
    _gdimensions.resize(_ndimension);
 | 
			
		||||
    _ldimensions.resize(_ndimension);
 | 
			
		||||
    _rdimensions.resize(_ndimension);
 | 
			
		||||
    _simd_layout.resize(_ndimension);
 | 
			
		||||
    _lstart.resize(_ndimension);
 | 
			
		||||
    _lend.resize(_ndimension);
 | 
			
		||||
 | 
			
		||||
      _ostride.resize(_ndimension);
 | 
			
		||||
      _istride.resize(_ndimension);
 | 
			
		||||
    _ostride.resize(_ndimension);
 | 
			
		||||
    _istride.resize(_ndimension);
 | 
			
		||||
 | 
			
		||||
      _fsites = _gsites = _osites = _isites = 1;
 | 
			
		||||
    _fsites = _gsites = _osites = _isites = 1;
 | 
			
		||||
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        _fdimensions[d] = dimensions[d];   // Global dimensions
 | 
			
		||||
        _gdimensions[d] = _fdimensions[d]; // Global dimensions
 | 
			
		||||
@@ -134,30 +134,30 @@ public:
 | 
			
		||||
 | 
			
		||||
        // Addressing support
 | 
			
		||||
        if (d == 0)
 | 
			
		||||
        {
 | 
			
		||||
          _ostride[d] = 1;
 | 
			
		||||
          _istride[d] = 1;
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    _ostride[d] = 1;
 | 
			
		||||
	    _istride[d] = 1;
 | 
			
		||||
	  }
 | 
			
		||||
        else
 | 
			
		||||
        {
 | 
			
		||||
          _ostride[d] = _ostride[d - 1] * _rdimensions[d - 1];
 | 
			
		||||
          _istride[d] = _istride[d - 1] * _simd_layout[d - 1];
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    _ostride[d] = _ostride[d - 1] * _rdimensions[d - 1];
 | 
			
		||||
	    _istride[d] = _istride[d - 1] * _simd_layout[d - 1];
 | 
			
		||||
	  }
 | 
			
		||||
      }
 | 
			
		||||
 | 
			
		||||
      ///////////////////////
 | 
			
		||||
      // subplane information
 | 
			
		||||
      ///////////////////////
 | 
			
		||||
      _slice_block.resize(_ndimension);
 | 
			
		||||
      _slice_stride.resize(_ndimension);
 | 
			
		||||
      _slice_nblock.resize(_ndimension);
 | 
			
		||||
    ///////////////////////
 | 
			
		||||
    // subplane information
 | 
			
		||||
    ///////////////////////
 | 
			
		||||
    _slice_block.resize(_ndimension);
 | 
			
		||||
    _slice_stride.resize(_ndimension);
 | 
			
		||||
    _slice_nblock.resize(_ndimension);
 | 
			
		||||
 | 
			
		||||
      int block = 1;
 | 
			
		||||
      int nblock = 1;
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
        nblock *= _rdimensions[d];
 | 
			
		||||
    int block = 1;
 | 
			
		||||
    int nblock = 1;
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      nblock *= _rdimensions[d];
 | 
			
		||||
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        nblock /= _rdimensions[d];
 | 
			
		||||
        _slice_block[d] = block;
 | 
			
		||||
@@ -165,7 +165,8 @@ public:
 | 
			
		||||
        _slice_nblock[d] = nblock;
 | 
			
		||||
        block = block * _rdimensions[d];
 | 
			
		||||
      }
 | 
			
		||||
    };
 | 
			
		||||
  };
 | 
			
		||||
};
 | 
			
		||||
}
 | 
			
		||||
 | 
			
		||||
NAMESPACE_END(Grid);
 | 
			
		||||
#endif
 | 
			
		||||
 
 | 
			
		||||
@@ -1,4 +1,4 @@
 | 
			
		||||
    /*************************************************************************************
 | 
			
		||||
/*************************************************************************************
 | 
			
		||||
 | 
			
		||||
    Grid physics library, www.github.com/paboyle/Grid 
 | 
			
		||||
 | 
			
		||||
@@ -24,179 +24,179 @@ Author: Peter Boyle <paboyle@ph.ed.ac.uk>
 | 
			
		||||
    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 */
 | 
			
		||||
*************************************************************************************/
 | 
			
		||||
/*  END LEGAL */
 | 
			
		||||
#ifndef GRID_CARTESIAN_RED_BLACK_H
 | 
			
		||||
#define GRID_CARTESIAN_RED_BLACK_H
 | 
			
		||||
 | 
			
		||||
 | 
			
		||||
namespace Grid {
 | 
			
		||||
NAMESPACE_BEGIN(Grid);
 | 
			
		||||
 | 
			
		||||
  static const int CbRed  =0;
 | 
			
		||||
  static const int CbBlack=1;
 | 
			
		||||
  static const int Even   =CbRed;
 | 
			
		||||
  static const int Odd    =CbBlack;
 | 
			
		||||
static const int CbRed  =0;
 | 
			
		||||
static const int CbBlack=1;
 | 
			
		||||
static const int Even   =CbRed;
 | 
			
		||||
static const int Odd    =CbBlack;
 | 
			
		||||
    
 | 
			
		||||
// Specialise this for red black grids storing half the data like a chess board.
 | 
			
		||||
class GridRedBlackCartesian : public GridBase
 | 
			
		||||
{
 | 
			
		||||
public:
 | 
			
		||||
    std::vector<int> _checker_dim_mask;
 | 
			
		||||
    int              _checker_dim;
 | 
			
		||||
    std::vector<int> _checker_board;
 | 
			
		||||
  std::vector<int> _checker_dim_mask;
 | 
			
		||||
  int              _checker_dim;
 | 
			
		||||
  std::vector<int> _checker_board;
 | 
			
		||||
 | 
			
		||||
    virtual int CheckerBoarded(int dim){
 | 
			
		||||
      if( dim==_checker_dim) return 1;
 | 
			
		||||
      else return 0;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoard(const std::vector<int> &site){
 | 
			
		||||
      int linear=0;
 | 
			
		||||
      assert(site.size()==_ndimension);
 | 
			
		||||
      for(int d=0;d<_ndimension;d++){ 
 | 
			
		||||
	if(_checker_dim_mask[d])
 | 
			
		||||
	  linear=linear+site[d];
 | 
			
		||||
      }
 | 
			
		||||
      return (linear&0x1);
 | 
			
		||||
  virtual int CheckerBoarded(int dim){
 | 
			
		||||
    if( dim==_checker_dim) return 1;
 | 
			
		||||
    else return 0;
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoard(const std::vector<int> &site){
 | 
			
		||||
    int linear=0;
 | 
			
		||||
    assert(site.size()==_ndimension);
 | 
			
		||||
    for(int d=0;d<_ndimension;d++){ 
 | 
			
		||||
      if(_checker_dim_mask[d])
 | 
			
		||||
	linear=linear+site[d];
 | 
			
		||||
    }
 | 
			
		||||
    return (linear&0x1);
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
 | 
			
		||||
    // Depending on the cb of site, we toggle source cb.
 | 
			
		||||
    // for block #b, element #e = (b, e)
 | 
			
		||||
    // we need 
 | 
			
		||||
    virtual int CheckerBoardShiftForCB(int source_cb,int dim,int shift,int ocb){
 | 
			
		||||
      if(dim != _checker_dim) return shift;
 | 
			
		||||
  // Depending on the cb of site, we toggle source cb.
 | 
			
		||||
  // for block #b, element #e = (b, e)
 | 
			
		||||
  // we need 
 | 
			
		||||
  virtual int CheckerBoardShiftForCB(int source_cb,int dim,int shift,int ocb){
 | 
			
		||||
    if(dim != _checker_dim) return shift;
 | 
			
		||||
 | 
			
		||||
      int fulldim =_fdimensions[dim];
 | 
			
		||||
      shift = (shift+fulldim)%fulldim;
 | 
			
		||||
    int fulldim =_fdimensions[dim];
 | 
			
		||||
    shift = (shift+fulldim)%fulldim;
 | 
			
		||||
 | 
			
		||||
      // Probably faster with table lookup;
 | 
			
		||||
      // or by looping over x,y,z and multiply rather than computing checkerboard.
 | 
			
		||||
    // Probably faster with table lookup;
 | 
			
		||||
    // or by looping over x,y,z and multiply rather than computing checkerboard.
 | 
			
		||||
	  
 | 
			
		||||
      if ( (source_cb+ocb)&1 ) {
 | 
			
		||||
	return (shift)/2;
 | 
			
		||||
      } else {
 | 
			
		||||
	return (shift+1)/2;
 | 
			
		||||
      }
 | 
			
		||||
    if ( (source_cb+ocb)&1 ) {
 | 
			
		||||
      return (shift)/2;
 | 
			
		||||
    } else {
 | 
			
		||||
      return (shift+1)/2;
 | 
			
		||||
    }
 | 
			
		||||
    virtual int  CheckerBoardFromOindexTable (int Oindex) {
 | 
			
		||||
      return _checker_board[Oindex];
 | 
			
		||||
    }
 | 
			
		||||
    virtual int  CheckerBoardFromOindex (int Oindex)
 | 
			
		||||
    {
 | 
			
		||||
      std::vector<int> ocoor;
 | 
			
		||||
      oCoorFromOindex(ocoor,Oindex);
 | 
			
		||||
      return CheckerBoard(ocoor);
 | 
			
		||||
    }
 | 
			
		||||
    virtual int CheckerBoardShift(int source_cb,int dim,int shift,int osite){
 | 
			
		||||
  }
 | 
			
		||||
  virtual int  CheckerBoardFromOindexTable (int Oindex) {
 | 
			
		||||
    return _checker_board[Oindex];
 | 
			
		||||
  }
 | 
			
		||||
  virtual int  CheckerBoardFromOindex (int Oindex)
 | 
			
		||||
  {
 | 
			
		||||
    std::vector<int> ocoor;
 | 
			
		||||
    oCoorFromOindex(ocoor,Oindex);
 | 
			
		||||
    return CheckerBoard(ocoor);
 | 
			
		||||
  }
 | 
			
		||||
  virtual int CheckerBoardShift(int source_cb,int dim,int shift,int osite){
 | 
			
		||||
 | 
			
		||||
      if(dim != _checker_dim) return shift;
 | 
			
		||||
    if(dim != _checker_dim) return shift;
 | 
			
		||||
 | 
			
		||||
      int ocb=CheckerBoardFromOindex(osite);
 | 
			
		||||
    int ocb=CheckerBoardFromOindex(osite);
 | 
			
		||||
      
 | 
			
		||||
      return CheckerBoardShiftForCB(source_cb,dim,shift,ocb);
 | 
			
		||||
    }
 | 
			
		||||
    return CheckerBoardShiftForCB(source_cb,dim,shift,ocb);
 | 
			
		||||
  }
 | 
			
		||||
    
 | 
			
		||||
    virtual int CheckerBoardDestination(int source_cb,int shift,int dim){
 | 
			
		||||
      if ( _checker_dim_mask[dim]  ) {
 | 
			
		||||
	// If _fdimensions[checker_dim] is odd, then shifting by 1 in other dims
 | 
			
		||||
	// does NOT cause a parity hop.
 | 
			
		||||
	int add=(dim==_checker_dim) ? 0 : _fdimensions[_checker_dim];
 | 
			
		||||
        if ( (shift+add) &0x1) {
 | 
			
		||||
            return 1-source_cb;
 | 
			
		||||
        } else {
 | 
			
		||||
            return source_cb;
 | 
			
		||||
        }
 | 
			
		||||
  virtual int CheckerBoardDestination(int source_cb,int shift,int dim){
 | 
			
		||||
    if ( _checker_dim_mask[dim]  ) {
 | 
			
		||||
      // If _fdimensions[checker_dim] is odd, then shifting by 1 in other dims
 | 
			
		||||
      // does NOT cause a parity hop.
 | 
			
		||||
      int add=(dim==_checker_dim) ? 0 : _fdimensions[_checker_dim];
 | 
			
		||||
      if ( (shift+add) &0x1) {
 | 
			
		||||
	return 1-source_cb;
 | 
			
		||||
      } else {
 | 
			
		||||
	return source_cb;
 | 
			
		||||
 | 
			
		||||
      }
 | 
			
		||||
    };
 | 
			
		||||
    } else {
 | 
			
		||||
      return source_cb;
 | 
			
		||||
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    // Create Redblack from original grid; require full grid pointer ?
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    GridRedBlackCartesian(const GridBase *base) : GridBase(base->_processors,*base)
 | 
			
		||||
    {
 | 
			
		||||
      int dims = base->_ndimension;
 | 
			
		||||
      std::vector<int> checker_dim_mask(dims,1);
 | 
			
		||||
      int checker_dim = 0;
 | 
			
		||||
      Init(base->_fdimensions,base->_simd_layout,base->_processors,checker_dim_mask,checker_dim);
 | 
			
		||||
    };
 | 
			
		||||
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    // Create redblack from original grid, with non-trivial checker dim mask
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    GridRedBlackCartesian(const GridBase *base,
 | 
			
		||||
			  const std::vector<int> &checker_dim_mask,
 | 
			
		||||
			  int checker_dim
 | 
			
		||||
			  ) :  GridBase(base->_processors,*base) 
 | 
			
		||||
    {
 | 
			
		||||
      Init(base->_fdimensions,base->_simd_layout,base->_processors,checker_dim_mask,checker_dim)  ;
 | 
			
		||||
    }
 | 
			
		||||
  };
 | 
			
		||||
 | 
			
		||||
    virtual ~GridRedBlackCartesian() = default;
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  // Create Redblack from original grid; require full grid pointer ?
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  GridRedBlackCartesian(const GridBase *base) : GridBase(base->_processors,*base)
 | 
			
		||||
  {
 | 
			
		||||
    int dims = base->_ndimension;
 | 
			
		||||
    std::vector<int> checker_dim_mask(dims,1);
 | 
			
		||||
    int checker_dim = 0;
 | 
			
		||||
    Init(base->_fdimensions,base->_simd_layout,base->_processors,checker_dim_mask,checker_dim);
 | 
			
		||||
  };
 | 
			
		||||
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  // Create redblack from original grid, with non-trivial checker dim mask
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  GridRedBlackCartesian(const GridBase *base,
 | 
			
		||||
			const std::vector<int> &checker_dim_mask,
 | 
			
		||||
			int checker_dim
 | 
			
		||||
			) :  GridBase(base->_processors,*base) 
 | 
			
		||||
  {
 | 
			
		||||
    Init(base->_fdimensions,base->_simd_layout,base->_processors,checker_dim_mask,checker_dim)  ;
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
  virtual ~GridRedBlackCartesian() = default;
 | 
			
		||||
#if 0
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    // Create redblack grid ;; deprecate these. Should not
 | 
			
		||||
    // need direct creation of redblack without a full grid to base on
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    GridRedBlackCartesian(const GridBase *base,
 | 
			
		||||
			  const std::vector<int> &dimensions,
 | 
			
		||||
			  const std::vector<int> &simd_layout,
 | 
			
		||||
			  const std::vector<int> &processor_grid,
 | 
			
		||||
			  const std::vector<int> &checker_dim_mask,
 | 
			
		||||
			  int checker_dim
 | 
			
		||||
			  ) :  GridBase(processor_grid,*base) 
 | 
			
		||||
    {
 | 
			
		||||
      Init(dimensions,simd_layout,processor_grid,checker_dim_mask,checker_dim);
 | 
			
		||||
    }
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  // Create redblack grid ;; deprecate these. Should not
 | 
			
		||||
  // need direct creation of redblack without a full grid to base on
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  GridRedBlackCartesian(const GridBase *base,
 | 
			
		||||
			const std::vector<int> &dimensions,
 | 
			
		||||
			const std::vector<int> &simd_layout,
 | 
			
		||||
			const std::vector<int> &processor_grid,
 | 
			
		||||
			const std::vector<int> &checker_dim_mask,
 | 
			
		||||
			int checker_dim
 | 
			
		||||
			) :  GridBase(processor_grid,*base) 
 | 
			
		||||
  {
 | 
			
		||||
    Init(dimensions,simd_layout,processor_grid,checker_dim_mask,checker_dim);
 | 
			
		||||
  }
 | 
			
		||||
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    // Create redblack grid
 | 
			
		||||
    ////////////////////////////////////////////////////////////
 | 
			
		||||
    GridRedBlackCartesian(const GridBase *base,
 | 
			
		||||
			  const std::vector<int> &dimensions,
 | 
			
		||||
			  const std::vector<int> &simd_layout,
 | 
			
		||||
			  const std::vector<int> &processor_grid) : GridBase(processor_grid,*base) 
 | 
			
		||||
    {
 | 
			
		||||
      std::vector<int> checker_dim_mask(dimensions.size(),1);
 | 
			
		||||
      int checker_dim = 0;
 | 
			
		||||
      Init(dimensions,simd_layout,processor_grid,checker_dim_mask,checker_dim);
 | 
			
		||||
    }
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  // Create redblack grid
 | 
			
		||||
  ////////////////////////////////////////////////////////////
 | 
			
		||||
  GridRedBlackCartesian(const GridBase *base,
 | 
			
		||||
			const std::vector<int> &dimensions,
 | 
			
		||||
			const std::vector<int> &simd_layout,
 | 
			
		||||
			const std::vector<int> &processor_grid) : GridBase(processor_grid,*base) 
 | 
			
		||||
  {
 | 
			
		||||
    std::vector<int> checker_dim_mask(dimensions.size(),1);
 | 
			
		||||
    int checker_dim = 0;
 | 
			
		||||
    Init(dimensions,simd_layout,processor_grid,checker_dim_mask,checker_dim);
 | 
			
		||||
  }
 | 
			
		||||
#endif
 | 
			
		||||
 | 
			
		||||
    void Init(const std::vector<int> &dimensions,
 | 
			
		||||
              const std::vector<int> &simd_layout,
 | 
			
		||||
              const std::vector<int> &processor_grid,
 | 
			
		||||
              const std::vector<int> &checker_dim_mask,
 | 
			
		||||
              int checker_dim)
 | 
			
		||||
    {
 | 
			
		||||
      ///////////////////////
 | 
			
		||||
      // Grid information
 | 
			
		||||
      ///////////////////////
 | 
			
		||||
      _checker_dim = checker_dim;
 | 
			
		||||
      assert(checker_dim_mask[checker_dim] == 1);
 | 
			
		||||
      _ndimension = dimensions.size();
 | 
			
		||||
      assert(checker_dim_mask.size() == _ndimension);
 | 
			
		||||
      assert(processor_grid.size() == _ndimension);
 | 
			
		||||
      assert(simd_layout.size() == _ndimension);
 | 
			
		||||
  void Init(const std::vector<int> &dimensions,
 | 
			
		||||
	    const std::vector<int> &simd_layout,
 | 
			
		||||
	    const std::vector<int> &processor_grid,
 | 
			
		||||
	    const std::vector<int> &checker_dim_mask,
 | 
			
		||||
	    int checker_dim)
 | 
			
		||||
  {
 | 
			
		||||
    ///////////////////////
 | 
			
		||||
    // Grid information
 | 
			
		||||
    ///////////////////////
 | 
			
		||||
    _checker_dim = checker_dim;
 | 
			
		||||
    assert(checker_dim_mask[checker_dim] == 1);
 | 
			
		||||
    _ndimension = dimensions.size();
 | 
			
		||||
    assert(checker_dim_mask.size() == _ndimension);
 | 
			
		||||
    assert(processor_grid.size() == _ndimension);
 | 
			
		||||
    assert(simd_layout.size() == _ndimension);
 | 
			
		||||
 | 
			
		||||
      _fdimensions.resize(_ndimension);
 | 
			
		||||
      _gdimensions.resize(_ndimension);
 | 
			
		||||
      _ldimensions.resize(_ndimension);
 | 
			
		||||
      _rdimensions.resize(_ndimension);
 | 
			
		||||
      _simd_layout.resize(_ndimension);
 | 
			
		||||
      _lstart.resize(_ndimension);
 | 
			
		||||
      _lend.resize(_ndimension);
 | 
			
		||||
    _fdimensions.resize(_ndimension);
 | 
			
		||||
    _gdimensions.resize(_ndimension);
 | 
			
		||||
    _ldimensions.resize(_ndimension);
 | 
			
		||||
    _rdimensions.resize(_ndimension);
 | 
			
		||||
    _simd_layout.resize(_ndimension);
 | 
			
		||||
    _lstart.resize(_ndimension);
 | 
			
		||||
    _lend.resize(_ndimension);
 | 
			
		||||
 | 
			
		||||
      _ostride.resize(_ndimension);
 | 
			
		||||
      _istride.resize(_ndimension);
 | 
			
		||||
    _ostride.resize(_ndimension);
 | 
			
		||||
    _istride.resize(_ndimension);
 | 
			
		||||
 | 
			
		||||
      _fsites = _gsites = _osites = _isites = 1;
 | 
			
		||||
    _fsites = _gsites = _osites = _isites = 1;
 | 
			
		||||
 | 
			
		||||
      _checker_dim_mask = checker_dim_mask;
 | 
			
		||||
    _checker_dim_mask = checker_dim_mask;
 | 
			
		||||
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        _fdimensions[d] = dimensions[d];
 | 
			
		||||
        _gdimensions[d] = _fdimensions[d];
 | 
			
		||||
@@ -204,11 +204,11 @@ public:
 | 
			
		||||
        _gsites = _gsites * _gdimensions[d];
 | 
			
		||||
 | 
			
		||||
        if (d == _checker_dim)
 | 
			
		||||
        {
 | 
			
		||||
          assert((_gdimensions[d] & 0x1) == 0);
 | 
			
		||||
          _gdimensions[d] = _gdimensions[d] / 2; // Remove a checkerboard
 | 
			
		||||
	  _gsites /= 2;
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    assert((_gdimensions[d] & 0x1) == 0);
 | 
			
		||||
	    _gdimensions[d] = _gdimensions[d] / 2; // Remove a checkerboard
 | 
			
		||||
	    _gsites /= 2;
 | 
			
		||||
	  }
 | 
			
		||||
        _ldimensions[d] = _gdimensions[d] / _processors[d];
 | 
			
		||||
        assert(_ldimensions[d] * _processors[d] == _gdimensions[d]);
 | 
			
		||||
        _lstart[d] = _processor_coor[d] * _ldimensions[d];
 | 
			
		||||
@@ -223,42 +223,42 @@ public:
 | 
			
		||||
        // all elements of a simd vector must have same checkerboard.
 | 
			
		||||
        // If Ls vectorised, this must still be the case; e.g. dwf rb5d
 | 
			
		||||
        if (_simd_layout[d] > 1)
 | 
			
		||||
        {
 | 
			
		||||
          if (checker_dim_mask[d])
 | 
			
		||||
          {
 | 
			
		||||
            assert((_rdimensions[d] & 0x1) == 0);
 | 
			
		||||
          }
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    if (checker_dim_mask[d])
 | 
			
		||||
	      {
 | 
			
		||||
		assert((_rdimensions[d] & 0x1) == 0);
 | 
			
		||||
	      }
 | 
			
		||||
	  }
 | 
			
		||||
 | 
			
		||||
        _osites *= _rdimensions[d];
 | 
			
		||||
        _isites *= _simd_layout[d];
 | 
			
		||||
 | 
			
		||||
        // Addressing support
 | 
			
		||||
        if (d == 0)
 | 
			
		||||
        {
 | 
			
		||||
          _ostride[d] = 1;
 | 
			
		||||
          _istride[d] = 1;
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    _ostride[d] = 1;
 | 
			
		||||
	    _istride[d] = 1;
 | 
			
		||||
	  }
 | 
			
		||||
        else
 | 
			
		||||
        {
 | 
			
		||||
          _ostride[d] = _ostride[d - 1] * _rdimensions[d - 1];
 | 
			
		||||
          _istride[d] = _istride[d - 1] * _simd_layout[d - 1];
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    _ostride[d] = _ostride[d - 1] * _rdimensions[d - 1];
 | 
			
		||||
	    _istride[d] = _istride[d - 1] * _simd_layout[d - 1];
 | 
			
		||||
	  }
 | 
			
		||||
      }
 | 
			
		||||
 | 
			
		||||
      ////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
      // subplane information
 | 
			
		||||
      ////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
      _slice_block.resize(_ndimension);
 | 
			
		||||
      _slice_stride.resize(_ndimension);
 | 
			
		||||
      _slice_nblock.resize(_ndimension);
 | 
			
		||||
    ////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    // subplane information
 | 
			
		||||
    ////////////////////////////////////////////////////////////////////////////////////////////
 | 
			
		||||
    _slice_block.resize(_ndimension);
 | 
			
		||||
    _slice_stride.resize(_ndimension);
 | 
			
		||||
    _slice_nblock.resize(_ndimension);
 | 
			
		||||
 | 
			
		||||
      int block = 1;
 | 
			
		||||
      int nblock = 1;
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
        nblock *= _rdimensions[d];
 | 
			
		||||
    int block = 1;
 | 
			
		||||
    int nblock = 1;
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      nblock *= _rdimensions[d];
 | 
			
		||||
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        nblock /= _rdimensions[d];
 | 
			
		||||
        _slice_block[d] = block;
 | 
			
		||||
@@ -267,55 +267,55 @@ public:
 | 
			
		||||
        block = block * _rdimensions[d];
 | 
			
		||||
      }
 | 
			
		||||
 | 
			
		||||
      ////////////////////////////////////////////////
 | 
			
		||||
      // Create a checkerboard lookup table
 | 
			
		||||
      ////////////////////////////////////////////////
 | 
			
		||||
      int rvol = 1;
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
    ////////////////////////////////////////////////
 | 
			
		||||
    // Create a checkerboard lookup table
 | 
			
		||||
    ////////////////////////////////////////////////
 | 
			
		||||
    int rvol = 1;
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        rvol = rvol * _rdimensions[d];
 | 
			
		||||
      }
 | 
			
		||||
      _checker_board.resize(rvol);
 | 
			
		||||
      for (int osite = 0; osite < _osites; osite++)
 | 
			
		||||
    _checker_board.resize(rvol);
 | 
			
		||||
    for (int osite = 0; osite < _osites; osite++)
 | 
			
		||||
      {
 | 
			
		||||
        _checker_board[osite] = CheckerBoardFromOindex(osite);
 | 
			
		||||
      }
 | 
			
		||||
    };
 | 
			
		||||
  };
 | 
			
		||||
 | 
			
		||||
  protected:
 | 
			
		||||
    virtual int oIndex(std::vector<int> &coor)
 | 
			
		||||
    {
 | 
			
		||||
      int idx = 0;
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
protected:
 | 
			
		||||
  virtual int oIndex(std::vector<int> &coor)
 | 
			
		||||
  {
 | 
			
		||||
    int idx = 0;
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        if (d == _checker_dim)
 | 
			
		||||
        {
 | 
			
		||||
          idx += _ostride[d] * ((coor[d] / 2) % _rdimensions[d]);
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    idx += _ostride[d] * ((coor[d] / 2) % _rdimensions[d]);
 | 
			
		||||
	  }
 | 
			
		||||
        else
 | 
			
		||||
        {
 | 
			
		||||
          idx += _ostride[d] * (coor[d] % _rdimensions[d]);
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    idx += _ostride[d] * (coor[d] % _rdimensions[d]);
 | 
			
		||||
	  }
 | 
			
		||||
      }
 | 
			
		||||
      return idx;
 | 
			
		||||
    };
 | 
			
		||||
    return idx;
 | 
			
		||||
  };
 | 
			
		||||
 | 
			
		||||
    virtual int iIndex(std::vector<int> &lcoor)
 | 
			
		||||
    {
 | 
			
		||||
      int idx = 0;
 | 
			
		||||
      for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
  virtual int iIndex(std::vector<int> &lcoor)
 | 
			
		||||
  {
 | 
			
		||||
    int idx = 0;
 | 
			
		||||
    for (int d = 0; d < _ndimension; d++)
 | 
			
		||||
      {
 | 
			
		||||
        if (d == _checker_dim)
 | 
			
		||||
        {
 | 
			
		||||
          idx += _istride[d] * (lcoor[d] / (2 * _rdimensions[d]));
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    idx += _istride[d] * (lcoor[d] / (2 * _rdimensions[d]));
 | 
			
		||||
	  }
 | 
			
		||||
        else
 | 
			
		||||
        {
 | 
			
		||||
          idx += _istride[d] * (lcoor[d] / _rdimensions[d]);
 | 
			
		||||
        }
 | 
			
		||||
	  {
 | 
			
		||||
	    idx += _istride[d] * (lcoor[d] / _rdimensions[d]);
 | 
			
		||||
	  }
 | 
			
		||||
      }
 | 
			
		||||
      return idx;
 | 
			
		||||
    }
 | 
			
		||||
    return idx;
 | 
			
		||||
  }
 | 
			
		||||
};
 | 
			
		||||
}
 | 
			
		||||
NAMESPACE_END(Grid);
 | 
			
		||||
#endif
 | 
			
		||||
 
 | 
			
		||||
		Reference in New Issue
	
	Block a user