From 077aa728b9e876324082127836ae5434c8d64e60 Mon Sep 17 00:00:00 2001
From: paboyle <paboyle@ph.ed.ac.uk>
Date: Thu, 30 Mar 2017 13:42:09 +0900
Subject: [PATCH] Fix the ZMobius (I think)

---
 lib/qcd/action/fermion/CayleyFermion5D.cc | 31 ++++++++++++++---------
 1 file changed, 19 insertions(+), 12 deletions(-)

diff --git a/lib/qcd/action/fermion/CayleyFermion5D.cc b/lib/qcd/action/fermion/CayleyFermion5D.cc
index 90614295..51677522 100644
--- a/lib/qcd/action/fermion/CayleyFermion5D.cc
+++ b/lib/qcd/action/fermion/CayleyFermion5D.cc
@@ -170,7 +170,6 @@ void CayleyFermion5D<Impl>::Mooee       (const FermionField &psi, FermionField &
   lower[0]   =-mass*lower[0];
   M5D(psi,psi,chi,lower,diag,upper);
 }
-
 template<class Impl>
 void CayleyFermion5D<Impl>::MooeeDag    (const FermionField &psi, FermionField &chi)
 {
@@ -192,7 +191,7 @@ void CayleyFermion5D<Impl>::MooeeDag    (const FermionField &psi, FermionField &
       lower[s]=-cee[s-1];
     }
   }
-  // Conjugate the terms ?
+  // Conjugate the terms 
   for (int s=0;s<Ls;s++){
     diag[s] =conjugate(diag[s]);
     upper[s]=conjugate(upper[s]);
@@ -219,14 +218,22 @@ void CayleyFermion5D<Impl>::MeooeDag5D    (const FermionField &psi, FermionField
   int Ls=this->Ls;
   std::vector<Coeff_t> diag =bs;
   std::vector<Coeff_t> upper=cs;
-  std::vector<Coeff_t> lower=cs;
-  upper[Ls-1]=-mass*upper[Ls-1];
-  lower[0]   =-mass*lower[0];
-  // Conjugate the terms ?
+  std::vector<Coeff_t> lower=cs; 
+
   for (int s=0;s<Ls;s++){
-    diag[s] =conjugate(diag[s]);
-    upper[s]=conjugate(upper[s]);
-    lower[s]=conjugate(lower[s]);
+    if ( s== 0 ) {
+      upper[s] = cs[s+1];
+      lower[s] =-mass*cs[Ls-1];
+    } else if ( s==(Ls-1) ) { 
+      upper[s] =-mass*cs[0];
+      lower[s] = cs[s-1];
+    } else { 
+      upper[s] = cs[s+1];
+      lower[s] = cs[s-1];
+    }
+    upper[s] = conjugate(upper[s]);
+    lower[s] = conjugate(lower[s]);
+    diag[s]  = conjugate(diag[s]);
   }
   M5Ddag(psi,psi,Din,lower,diag,upper);
 }
@@ -313,7 +320,7 @@ void CayleyFermion5D<Impl>::MDeriv  (GaugeField &mat,const FermionField &U,const
     this->DhopDeriv(mat,U,Din,dag);
   } else {
     //      U d/du [D_w D5]^dag V = U D5^dag d/du DW^dag Y // implicit adj on U in call
-    Meooe5D(U,Din);
+    MeooeDag5D(U,Din);
     this->DhopDeriv(mat,Din,V,dag);
   }
 };
@@ -328,7 +335,7 @@ void CayleyFermion5D<Impl>::MoeDeriv(GaugeField &mat,const FermionField &U,const
     this->DhopDerivOE(mat,U,Din,dag);
   } else {
     //      U d/du [D_w D5]^dag V = U D5^dag d/du DW^dag Y // implicit adj on U in call
-      Meooe5D(U,Din);
+      MeooeDag5D(U,Din);
       this->DhopDerivOE(mat,Din,V,dag);
   }
 };
@@ -343,7 +350,7 @@ void CayleyFermion5D<Impl>::MeoDeriv(GaugeField &mat,const FermionField &U,const
     this->DhopDerivEO(mat,U,Din,dag);
   } else {
     //      U d/du [D_w D5]^dag V = U D5^dag d/du DW^dag Y // implicit adj on U in call
-    Meooe5D(U,Din);
+    MeooeDag5D(U,Din);
     this->DhopDerivEO(mat,Din,V,dag);
   }
 };