Const correctness patch

by Christoph. Differs from the one by Rudy in BFM since it vectorises the twisted
4D solves in pairs.

.gitignore

@@ -83,6 +83,7 @@ ltmain.sh
 .Trashes
 ehthumbs.db
 Thumbs.db
+.dirstamp
 
 # build directory #
 ###################
@@ -97,11 +98,8 @@ build.sh
 
 # Eigen source #
 ################
-lib/Eigen/*
-
-# FFTW source #
-################
-lib/fftw/*
+Grid/Eigen
+Eigen/*
 
 # libtool macros #
 ##################
@@ -112,21 +110,7 @@ m4/libtool.m4
 ################
 gh-pages/
 
-# Buck files #
-##############
-.buck*
-buck-out
-BUCK
-make-bin-BUCK.sh
-
 # generated sources #
 #####################
-lib/qcd/spin/gamma-gen/*.h
-lib/qcd/spin/gamma-gen/*.cc
-lib/version.h
-
-# vs code editor files #
-########################
-.vscode/
-.vscode/settings.json
-settings.json
+Grid/qcd/spin/gamma-gen/*.h
+Grid/qcd/spin/gamma-gen/*.cc

.travis.yml

@@ -9,6 +9,11 @@ matrix:
     - os:        osx
       osx_image: xcode8.3
       compiler: clang
+      env: PREC=single
+    - os:        osx
+      osx_image: xcode8.3
+      compiler: clang
+      env: PREC=double
       
 before_install:
     - export GRIDDIR=`pwd`
@@ -16,7 +21,7 @@ before_install:
     - if [[ "$TRAVIS_OS_NAME" == "linux" ]] && [[ "$CC" == "clang" ]]; then export PATH="${GRIDDIR}/clang/bin:${PATH}"; fi
     - if [[ "$TRAVIS_OS_NAME" == "linux" ]] && [[ "$CC" == "clang" ]]; then export LD_LIBRARY_PATH="${GRIDDIR}/clang/lib:${LD_LIBRARY_PATH}"; fi
     - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then brew update; fi
-    - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then brew install libmpc; fi
+    - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then brew install libmpc openssl; fi
     
 install:
     - export CWD=`pwd`
@@ -33,6 +38,7 @@ install:
     - which $CXX
     - $CXX --version
     - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then export LDFLAGS='-L/usr/local/lib'; fi
+    - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then export EXTRACONF='--with-openssl=/usr/local/opt/openssl'; fi
     
 script:
     - ./bootstrap.sh
@@ -49,12 +55,7 @@ script:
     - make -j4
     - make install
     - cd $CWD/build
-    - ../configure --enable-precision=single --enable-simd=SSE4 --enable-comms=none --with-lime=$CWD/build/lime/install
+    - ../configure --enable-precision=$PREC --enable-simd=SSE4 --enable-comms=none --with-lime=$CWD/build/lime/install ${EXTRACONF}
     - make -j4 
     - ./benchmarks/Benchmark_dwf --threads 1 --debug-signals
-    - echo make clean
-    - ../configure --enable-precision=double --enable-simd=SSE4 --enable-comms=none --with-lime=$CWD/build/lime/install
-    - make -j4
-    - ./benchmarks/Benchmark_dwf --threads 1 --debug-signals
     - make check
-

Grid/GridCore.h

@@ -48,6 +48,7 @@ Author: paboyle <paboyle@ph.ed.ac.uk>
 #include <Grid/serialisation/Serialisation.h>
 #include <Grid/threads/Threads.h>
 #include <Grid/util/Util.h>
+#include <Grid/util/Sha.h>
 #include <Grid/communicator/Communicator.h> 
 #include <Grid/cartesian/Cartesian.h>    
 #include <Grid/tensors/Tensors.h>      

Grid/Grid_Eigen_Dense.h

@@ -1,10 +1,14 @@
 #pragma once
+// Force Eigen to use MKL if Grid has been configured with --enable-mkl
+#ifdef USE_MKL
 #define EIGEN_USE_MKL_ALL
+#endif
+
 #if defined __GNUC__
 #pragma GCC diagnostic push
 #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
 #endif
-#include <Eigen/Dense>
+#include <Grid/Eigen/Dense>
 #if defined __GNUC__
 #pragma GCC diagnostic pop
 #endif

Grid/Makefile.am

@@ -21,6 +21,32 @@ if BUILD_HDF5
   extra_headers+=serialisation/Hdf5Type.h
 endif
 
+all: version-cache
+
+version-cache:
+	@if [ `git status --porcelain | grep -v '??' | wc -l` -gt 0 ]; then\
+		a="uncommited changes";\
+	else\
+		a="clean";\
+	fi;\
+	echo "`git log -n 1 --format=format:"#define GITHASH \\"%H:%d $$a\\"%n" HEAD`" > vertmp;\
+	if [ -e version-cache ]; then\
+		d=`diff vertmp version-cache`;\
+		if [ "$${d}" != "" ]; then\
+			mv vertmp version-cache;\
+			rm -f Version.h;\
+		fi;\
+	else\
+		mv vertmp version-cache;\
+		rm -f Version.h;\
+	fi;\
+	rm -f vertmp
+
+Version.h:
+	cp version-cache Version.h
+
+.PHONY: version-cache
+
 #
 # Libraries
 #
@@ -30,8 +56,8 @@ include Eigen.inc
 lib_LIBRARIES = libGrid.a
 
 CCFILES += $(extra_sources)
-HFILES  += $(extra_headers)
+HFILES  += $(extra_headers) Config.h Version.h
 
 libGrid_a_SOURCES              = $(CCFILES)
-libGrid_adir                   = $(pkgincludedir)
-nobase_dist_pkginclude_HEADERS = $(HFILES) $(eigen_files) $(eigen_unsupp_files) Config.h
+libGrid_adir                   = $(includedir)/Grid
+nobase_dist_pkginclude_HEADERS = $(HFILES) $(eigen_files) $(eigen_unsupp_files)

Grid/algorithms/iterative/ConjugateGradient.h

@@ -150,13 +150,13 @@ class ConjugateGradient : public OperatorFunction<Field> {
 	std::cout << GridLogMessage << "\tTrue residual " << true_residual<<std::endl;
 	std::cout << GridLogMessage << "\tTarget " << Tolerance << std::endl;
 
-        std::cout << GridLogMessage << "Time breakdown "<<std::endl;
-	std::cout << GridLogMessage << "\tElapsed    " << SolverTimer.Elapsed() <<std::endl;
-	std::cout << GridLogMessage << "\tMatrix     " << MatrixTimer.Elapsed() <<std::endl;
-	std::cout << GridLogMessage << "\tLinalg     " << LinalgTimer.Elapsed() <<std::endl;
-	std::cout << GridLogMessage << "\tInner      " << InnerTimer.Elapsed() <<std::endl;
-	std::cout << GridLogMessage << "\tAxpyNorm   " << AxpyNormTimer.Elapsed() <<std::endl;
-	std::cout << GridLogMessage << "\tLinearComb " << LinearCombTimer.Elapsed() <<std::endl;
+        std::cout << GridLogPerformance << "Time breakdown "<<std::endl;
+	std::cout << GridLogPerformance << "\tElapsed    " << SolverTimer.Elapsed() <<std::endl;
+	std::cout << GridLogPerformance << "\tMatrix     " << MatrixTimer.Elapsed() <<std::endl;
+	std::cout << GridLogPerformance << "\tLinalg     " << LinalgTimer.Elapsed() <<std::endl;
+	std::cout << GridLogPerformance << "\tInner      " << InnerTimer.Elapsed() <<std::endl;
+	std::cout << GridLogPerformance << "\tAxpyNorm   " << AxpyNormTimer.Elapsed() <<std::endl;
+	std::cout << GridLogPerformance << "\tLinearComb " << LinearCombTimer.Elapsed() <<std::endl;
 
         if (ErrorOnNoConverge) assert(true_residual / Tolerance < 10000.0);
 

Grid/algorithms/iterative/Deflation.h

@@ -31,21 +31,13 @@ Author: Peter Boyle <paboyle@ph.ed.ac.uk>
 namespace Grid { 
 
 template<class Field>
-class Guesser {
-public:
-  Guesser(void) = default;
-  virtual ~Guesser(void) = default;
-  virtual void operator()(const Field &src, Field &guess) = 0;
-};
-
-template<class Field>
-class ZeroGuesser: public Guesser<Field> {
+class ZeroGuesser: public LinearFunction<Field> {
 public:
   virtual void operator()(const Field &src, Field &guess) { guess = zero; };
 };
 
 template<class Field>
-class SourceGuesser: public Guesser<Field> {
+class SourceGuesser: public LinearFunction<Field> {
 public:
   virtual void operator()(const Field &src, Field &guess) { guess = src; };
 };
@@ -54,7 +46,7 @@ public:
 // Fine grid deflation
 ////////////////////////////////
 template<class Field>
-class DeflatedGuesser: public Guesser<Field> {
+class DeflatedGuesser: public LinearFunction<Field> {
 private:
   const std::vector<Field> &evec;
   const std::vector<RealD> &eval;
@@ -76,7 +68,7 @@ public:
 };
 
 template<class FineField, class CoarseField>
-class LocalCoherenceDeflatedGuesser: public Guesser<FineField> {
+class LocalCoherenceDeflatedGuesser: public LinearFunction<FineField> {
 private:
   const std::vector<FineField>   &subspace;
   const std::vector<CoarseField> &evec_coarse;

Grid/communicator/SharedMemoryMPI.cc

@@ -132,7 +132,6 @@ int Log2Size(int TwoToPower,int MAXLOG2)
 }
 void GlobalSharedMemory::OptimalCommunicator(const std::vector<int> &processors,Grid_MPI_Comm & optimal_comm)
 {
-#undef HYPERCUBE 
 #ifdef HYPERCUBE
   ////////////////////////////////////////////////////////////////
   // Assert power of two shm_size.
@@ -175,7 +174,7 @@ void GlobalSharedMemory::OptimalCommunicator(const std::vector<int> &processors,
 
   std::string hname(name);
   std::cout << "hostname "<<hname<<std::endl;
-  std::cout << "R " << R << " I " << I << " N "<< N<<
+  std::cout << "R " << R << " I " << I << " N "<< N
             << " hypercoor 0x"<<std::hex<<hypercoor<<std::dec<<std::endl;
 
   //////////////////////////////////////////////////////////////////

Grid/lattice/Lattice_rng.h

@@ -251,7 +251,7 @@ namespace Grid {
 
       dist[0].reset();
       for(int idx=0;idx<words;idx++){
-  fillScalar(buf[idx],dist[0],_generators[0]);
+	fillScalar(buf[idx],dist[0],_generators[0]);
       }
 
       CartesianCommunicator::BroadcastWorld(0,(void *)&l,sizeof(l));
@@ -283,7 +283,7 @@ namespace Grid {
       RealF *pointer=(RealF *)&l;
       dist[0].reset();
       for(int i=0;i<2*vComplexF::Nsimd();i++){
-  fillScalar(pointer[i],dist[0],_generators[0]);
+	fillScalar(pointer[i],dist[0],_generators[0]);
       }
       CartesianCommunicator::BroadcastWorld(0,(void *)&l,sizeof(l));
     }
@@ -291,7 +291,7 @@ namespace Grid {
       RealD *pointer=(RealD *)&l;
       dist[0].reset();
       for(int i=0;i<2*vComplexD::Nsimd();i++){
-  fillScalar(pointer[i],dist[0],_generators[0]);
+	fillScalar(pointer[i],dist[0],_generators[0]);
       }
       CartesianCommunicator::BroadcastWorld(0,(void *)&l,sizeof(l));
     }
@@ -299,7 +299,7 @@ namespace Grid {
       RealF *pointer=(RealF *)&l;
       dist[0].reset();
       for(int i=0;i<vRealF::Nsimd();i++){
-  fillScalar(pointer[i],dist[0],_generators[0]);
+	fillScalar(pointer[i],dist[0],_generators[0]);
       }
       CartesianCommunicator::BroadcastWorld(0,(void *)&l,sizeof(l));
     }
@@ -317,6 +317,19 @@ namespace Grid {
       std::seed_seq src(seeds.begin(),seeds.end());
       Seed(src,0);
     }
+
+    void SeedUniqueString(const std::string &s){
+      std::vector<int> seeds;
+      std::stringstream sha;
+      seeds = GridChecksum::sha256_seeds(s);
+      for(int i=0;i<seeds.size();i++) { 
+        sha << std::hex << seeds[i];
+      }
+      std::cout << GridLogMessage << "Intialising serial RNG with unique string '" 
+                << s << "'" << std::endl;
+      std::cout << GridLogMessage << "Seed SHA256: " << sha.str() << std::endl;
+      SeedFixedIntegers(seeds);
+    }
   };
 
   class GridParallelRNG : public GridRNGbase {
@@ -377,6 +390,14 @@ namespace Grid {
       _time_counter += usecond()- inner_time_counter;
     };
 
+    void SeedUniqueString(const std::string &s){
+      std::vector<int> seeds;
+      seeds = GridChecksum::sha256_seeds(s);
+      std::cout << GridLogMessage << "Intialising parallel RNG with unique string '" 
+                << s << "'" << std::endl;
+      std::cout << GridLogMessage << "Seed SHA256: " << GridChecksum::sha256_string(seeds) << std::endl;
+      SeedFixedIntegers(seeds);
+    }
     void SeedFixedIntegers(const std::vector<int> &seeds){
 
       // Everyone generates the same seed_seq based on input seeds

Grid/lattice/Lattice_transfer.h

@@ -485,7 +485,7 @@ void InsertSliceLocal(const Lattice<vobj> &lowDim, Lattice<vobj> & higherDim,int
 
 
 template<class vobj>
-void ExtractSliceLocal(Lattice<vobj> &lowDim, Lattice<vobj> & higherDim,int slice_lo,int slice_hi, int orthog)
+void ExtractSliceLocal(Lattice<vobj> &lowDim, const Lattice<vobj> & higherDim,int slice_lo,int slice_hi, int orthog)
 {
   typedef typename vobj::scalar_object sobj;
 
@@ -520,7 +520,7 @@ void ExtractSliceLocal(Lattice<vobj> &lowDim, Lattice<vobj> & higherDim,int slic
 
 
 template<class vobj>
-void Replicate(Lattice<vobj> &coarse,Lattice<vobj> & fine)
+void Replicate(const Lattice<vobj> &coarse,Lattice<vobj> & fine)
 {
   typedef typename vobj::scalar_object sobj;
 

Grid/parallelIO/IldgIO.h

@@ -250,8 +250,7 @@ class GridLimeReader : public BinaryIO {
   ////////////////////////////////////////////
   // Read a generic serialisable object
   ////////////////////////////////////////////
-  template<class serialisable_object>
-  void readLimeObject(serialisable_object &object,std::string object_name,std::string record_name)
+  void readLimeObject(std::string &xmlstring,std::string record_name)
   {
     // should this be a do while; can we miss a first record??
     while ( limeReaderNextRecord(LimeR) == LIME_SUCCESS ) { 
@@ -266,15 +265,23 @@ class GridLimeReader : public BinaryIO {
 	limeReaderReadData((void *)&xmlc[0], &nbytes, LimeR);    
 	//	std::cout << GridLogMessage<< " readLimeObject matches XML " << &xmlc[0] <<std::endl;
 
-  std::string xmlstring(&xmlc[0]);
-	XmlReader RD(xmlstring, true, "");
-	read(RD,object_name,object);
+   xmlstring = std::string(&xmlc[0]);
 	return;
       }
 
     }  
     assert(0);
   }
+
+  template<class serialisable_object>
+  void readLimeObject(serialisable_object &object,std::string object_name,std::string record_name)
+  {
+    std::string xmlstring;
+
+    readLimeObject(xmlstring, record_name);
+	  XmlReader RD(xmlstring, true, "");
+	  read(RD,object_name,object);
+  }
 };
 
 class GridLimeWriter : public BinaryIO 

Grid/perfmon/Timer.h

@@ -102,6 +102,9 @@ public:
     assert(running == false);
     return (uint64_t) accumulator.count();
   }
+  bool isRunning(void){
+    return running;
+  }
 };
 
 }

Grid/qcd/QCD.h

@@ -90,17 +90,20 @@ namespace QCD {
     // That probably makes for GridRedBlack4dCartesian grid.
 
     // s,sp,c,spc,lc
-    template<typename vtype> using iSinglet                   = iScalar<iScalar<iScalar<vtype> > >;
-    template<typename vtype> using iSpinMatrix                = iScalar<iMatrix<iScalar<vtype>, Ns> >;
-    template<typename vtype> using iColourMatrix              = iScalar<iScalar<iMatrix<vtype, Nc> > > ;
-    template<typename vtype> using iSpinColourMatrix          = iScalar<iMatrix<iMatrix<vtype, Nc>, Ns> >;
-    template<typename vtype> using iLorentzColourMatrix       = iVector<iScalar<iMatrix<vtype, Nc> >, Nd > ;
-    template<typename vtype> using iDoubleStoredColourMatrix  = iVector<iScalar<iMatrix<vtype, Nc> >, Nds > ;
-    template<typename vtype> using iSpinVector                = iScalar<iVector<iScalar<vtype>, Ns> >;
-    template<typename vtype> using iColourVector              = iScalar<iScalar<iVector<vtype, Nc> > >;
-    template<typename vtype> using iSpinColourVector          = iScalar<iVector<iVector<vtype, Nc>, Ns> >;
-    template<typename vtype> using iHalfSpinVector            = iScalar<iVector<iScalar<vtype>, Nhs> >;
-    template<typename vtype> using iHalfSpinColourVector      = iScalar<iVector<iVector<vtype, Nc>, Nhs> >;
+
+    template<typename vtype> using iSinglet                     = iScalar<iScalar<iScalar<vtype> > >;
+    template<typename vtype> using iSpinMatrix                  = iScalar<iMatrix<iScalar<vtype>, Ns> >;
+    template<typename vtype> using iColourMatrix                = iScalar<iScalar<iMatrix<vtype, Nc> > > ;
+    template<typename vtype> using iSpinColourMatrix            = iScalar<iMatrix<iMatrix<vtype, Nc>, Ns> >;
+    template<typename vtype> using iLorentzColourMatrix         = iVector<iScalar<iMatrix<vtype, Nc> >, Nd > ;
+    template<typename vtype> using iDoubleStoredColourMatrix    = iVector<iScalar<iMatrix<vtype, Nc> >, Nds > ;
+    template<typename vtype> using iSpinVector                  = iScalar<iVector<iScalar<vtype>, Ns> >;
+    template<typename vtype> using iColourVector                = iScalar<iScalar<iVector<vtype, Nc> > >;
+    template<typename vtype> using iSpinColourVector            = iScalar<iVector<iVector<vtype, Nc>, Ns> >;
+    template<typename vtype> using iHalfSpinVector              = iScalar<iVector<iScalar<vtype>, Nhs> >;
+    template<typename vtype> using iHalfSpinColourVector        = iScalar<iVector<iVector<vtype, Nc>, Nhs> >;
+    template<typename vtype> using iSpinColourSpinColourMatrix  = iScalar<iMatrix<iMatrix<iMatrix<iMatrix<vtype, Nc>, Ns>, Nc>, Ns> >;
+
 
     template<typename vtype> using iGparitySpinColourVector       = iVector<iVector<iVector<vtype, Nc>, Ns>, Ngp >;
     template<typename vtype> using iGparityHalfSpinColourVector   = iVector<iVector<iVector<vtype, Nc>, Nhs>, Ngp >;
@@ -127,10 +130,28 @@ namespace QCD {
     typedef iSpinColourMatrix<Complex  >    SpinColourMatrix;
     typedef iSpinColourMatrix<ComplexF >    SpinColourMatrixF;
     typedef iSpinColourMatrix<ComplexD >    SpinColourMatrixD;
-
+    
     typedef iSpinColourMatrix<vComplex >    vSpinColourMatrix;
     typedef iSpinColourMatrix<vComplexF>    vSpinColourMatrixF;
     typedef iSpinColourMatrix<vComplexD>    vSpinColourMatrixD;
+    
+    // SpinColourSpinColour matrix
+    typedef iSpinColourSpinColourMatrix<Complex  >    SpinColourSpinColourMatrix;
+    typedef iSpinColourSpinColourMatrix<ComplexF >    SpinColourSpinColourMatrixF;
+    typedef iSpinColourSpinColourMatrix<ComplexD >    SpinColourSpinColourMatrixD;
+
+    typedef iSpinColourSpinColourMatrix<vComplex >    vSpinColourSpinColourMatrix;
+    typedef iSpinColourSpinColourMatrix<vComplexF>    vSpinColourSpinColourMatrixF;
+    typedef iSpinColourSpinColourMatrix<vComplexD>    vSpinColourSpinColourMatrixD;
+
+    // SpinColourSpinColour matrix
+    typedef iSpinColourSpinColourMatrix<Complex  >    SpinColourSpinColourMatrix;
+    typedef iSpinColourSpinColourMatrix<ComplexF >    SpinColourSpinColourMatrixF;
+    typedef iSpinColourSpinColourMatrix<ComplexD >    SpinColourSpinColourMatrixD;
+
+    typedef iSpinColourSpinColourMatrix<vComplex >    vSpinColourSpinColourMatrix;
+    typedef iSpinColourSpinColourMatrix<vComplexF>    vSpinColourSpinColourMatrixF;
+    typedef iSpinColourSpinColourMatrix<vComplexD>    vSpinColourSpinColourMatrixD;
 
     // LorentzColour
     typedef iLorentzColourMatrix<Complex  > LorentzColourMatrix;
@@ -229,6 +250,9 @@ namespace QCD {
     typedef Lattice<vSpinColourMatrixF>     LatticeSpinColourMatrixF;
     typedef Lattice<vSpinColourMatrixD>     LatticeSpinColourMatrixD;
 
+    typedef Lattice<vSpinColourSpinColourMatrix>      LatticeSpinColourSpinColourMatrix;
+    typedef Lattice<vSpinColourSpinColourMatrixF>     LatticeSpinColourSpinColourMatrixF;
+    typedef Lattice<vSpinColourSpinColourMatrixD>     LatticeSpinColourSpinColourMatrixD;
 
     typedef Lattice<vLorentzColourMatrix>  LatticeLorentzColourMatrix;
     typedef Lattice<vLorentzColourMatrixF> LatticeLorentzColourMatrixF;