From 1ed69816b9ecb9838c0a9d613e08fc2fc042afa2 Mon Sep 17 00:00:00 2001
From: Guido Cossu <guido.cossu@ed.ac.uk>
Date: Thu, 30 Mar 2017 11:14:27 +0100
Subject: [PATCH] First steps for the force term

---
 lib/qcd/action/fermion/WilsonCloverFermion.cc | 13 ++++++++++++-
 lib/qcd/utils/WilsonLoops.h                   |  2 +-
 2 files changed, 13 insertions(+), 2 deletions(-)

diff --git a/lib/qcd/action/fermion/WilsonCloverFermion.cc b/lib/qcd/action/fermion/WilsonCloverFermion.cc
index b94c72c0..7e51dcfe 100644
--- a/lib/qcd/action/fermion/WilsonCloverFermion.cc
+++ b/lib/qcd/action/fermion/WilsonCloverFermion.cc
@@ -82,8 +82,10 @@ namespace QCD {
     WilsonLoops<Impl>::FieldStrength(Ex, _Umu, Tdir, Xdir);
     WilsonLoops<Impl>::FieldStrength(Ey, _Umu, Tdir, Ydir);
     WilsonLoops<Impl>::FieldStrength(Ez, _Umu, Tdir, Zdir);
+    // Save the contracted term with sigma 
+    // into a dense matrix site by site
 
-    // Invert the Moo, Mee terms (?)
+    // Invert the Moo, Mee terms (using Eigen)
   }
 
 
@@ -110,6 +112,14 @@ namespace QCD {
   template<class Impl>
   void WilsonCloverFermion<Impl>::MDeriv(GaugeField&mat, const FermionField&U, const FermionField&V, int dag){
     GaugeField tmp(mat._grid);
+
+    conformable(U._grid, _grid);
+    conformable(U._grid, V._grid);
+    conformable(U._grid, mat._grid);
+
+    mat.checkerboard = U.checkerboard;
+    tmp.checkerboard = U.checkerboard;
+
     this->DhopDeriv(mat, U, V, dag);
     MooDeriv(tmp, U, V, dag);
     mat += tmp;
@@ -118,6 +128,7 @@ namespace QCD {
  // Derivative parts
   template<class Impl>
   void WilsonCloverFermion<Impl>::MooDeriv(GaugeField&mat, const FermionField&U, const FermionField&V, int dag){
+    
     assert(0); // not implemented yet
   }
 
diff --git a/lib/qcd/utils/WilsonLoops.h b/lib/qcd/utils/WilsonLoops.h
index ca2b2b8d..1984d2b8 100644
--- a/lib/qcd/utils/WilsonLoops.h
+++ b/lib/qcd/utils/WilsonLoops.h
@@ -301,7 +301,7 @@ public:
       GaugeMat v = adj(Vup) - adj(Vdn);
       GaugeMat u = PeekIndex<LorentzIndex>(Umu, mu);  // some redundant copies
       GaugeMat vu = v*u;
-      FS = 0.25*Ta(u*v - Cshift(vu, mu, +1));
+      FS = 0.25*Ta(u*v + Cshift(vu, mu, +1));// need jsut the antihermite part
   }