From aceaee774c2d7d6ef4ca63899bda53de16ec8060 Mon Sep 17 00:00:00 2001
From: paboyle <paboyle@ph.ed.ac.uk>
Date: Wed, 11 May 2016 14:12:02 +0100
Subject: [PATCH] ExtractSlice / InsertSlice for lower dimensional lattices
 where the lattice is not distributed in the orthogonal direction. Useful for
 fermion 4d/5d etc..

---
 lib/lattice/Lattice_transfer.h | 85 +++++++++++++++++++++++++++++++++-
 1 file changed, 84 insertions(+), 1 deletion(-)

diff --git a/lib/lattice/Lattice_transfer.h b/lib/lattice/Lattice_transfer.h
index f39b4bab..eb6d6515 100644
--- a/lib/lattice/Lattice_transfer.h
+++ b/lib/lattice/Lattice_transfer.h
@@ -44,7 +44,7 @@ inline void subdivides(GridBase *coarse,GridBase *fine)
   }
 }
 
-
+ 
   ////////////////////////////////////////////////////////////////////////////////////////////
   // remove and insert a half checkerboard
   ////////////////////////////////////////////////////////////////////////////////////////////
@@ -325,6 +325,89 @@ inline void blockPromote(const Lattice<iVector<CComplex,nbasis > > &coarseData,
   
 }
 
+template<class vobj>
+void InsertSlice(const Lattice<vobj> &lowDim,Lattice<vobj> & higherDim,int slice, int orthog)
+{
+  typedef typename vobj::scalar_object sobj;
+  sobj s;
+
+  GridBase *lg = lowdim._grid;
+  GridBase *hg = higherdim._grid;
+  int nl = lg->_ndimension;
+  int nh = hg->_ndimension;
+
+  assert(nl+1 == nh);
+  assert(orthog<nh);
+  assert(orthog>=0);
+  assert(hg->_processors[orthog]==0);
+
+  int dl; dl = 0;
+  for(int d=0;d<nh;d++){
+    if ( d != orthog) {
+      assert(lg->_processors[dl]  == hg->_processors[d]);
+      assert(lg->_ldimensions[dl] == hg->_ldimensions[d]);
+      dl++;
+    }
+  }
+
+  // the above should guarantee that the operations are local
+  std::vector<int> lcoor(nl);
+  std::vector<int> hcoor(nh);
+  hcoor[orthog] = slice;
+  for(int idx=0;idx<lg->lSites();idx++){
+    lg->LocalIndexToLocalCoor(idx,lcoor);
+    dl=0;
+    for(int d=0;d<nh;d++){
+      if ( d!orthog ) { 
+	hcoor[d]=lcoor[ld++];
+      }
+    }
+    peekLocalSite(s,higherDim,hcoor);
+    pokeLocalSite(s,lowDim,lcoor);
+  }
+}
+
+template<class vobj>
+void ExtractSlice(Lattice<vobj> &lowDim,const Lattice<vobj> & higherDim,int slice, int dir)
+{
+  typedef typename vobj::scalar_object sobj;
+  sobj s;
+
+  GridBase *lg = lowdim._grid;
+  GridBase *hg = higherdim._grid;
+  int nl = lg->_ndimension;
+  int nh = hg->_ndimension;
+
+  assert(nl+1 == nh);
+  assert(orthog<nh);
+  assert(orthog>=0);
+  assert(hg->_processors[orthog]==0);
+
+  int dl; dl = 0;
+  for(int d=0;d<nh;d++){
+    if ( d != orthog) {
+      assert(lg->_processors[dl]  == hg->_processors[d]);
+      assert(lg->_ldimensions[dl] == hg->_ldimensions[d]);
+      dl++;
+    }
+  }
+  // the above should guarantee that the operations are local
+  std::vector<int> lcoor(nl);
+  std::vector<int> hcoor(nh);
+  hcoor[orthog] = slice;
+  for(int idx=0;idx<lg->lSites();idx++){
+    lg->LocalIndexToLocalCoor(idx,lcoor);
+    dl=0;
+    for(int d=0;d<nh;d++){
+      if ( d!orthog ) { 
+	hcoor[d]=lcoor[ld++];
+      }
+    }
+    peekLocalSite(s,lowDim,lcoor);
+    pokeLocalSite(s,higherDim,hcoor);
+  }
+
+}
 
 template<class vobj>
 void Replicate(Lattice<vobj> &coarse,Lattice<vobj> & fine)