Feature/deflated solvers revisit

CMakeLists.txt

@@ -553,6 +553,7 @@ if(QUDA_ARPACK)
   if(QUDA_DOWNLOAD_ARPACK)
   #TODO: switch to CPM
     include(GNUInstallDirs)
+    include(ExternalProject)
     ExternalProject_Add(
       ARPACK-NG
       GIT_REPOSITORY https://github.com/opencollab/arpack-ng.git

include/blas_quda.h

@@ -1,6 +1,4 @@
-#ifndef _QUDA_BLAS_H
-#define _QUDA_BLAS_H
-
+#pragma once
 #include <quda_internal.h>
 #include <color_spinor_field.h>
 
@@ -35,6 +33,8 @@ namespace quda {
     }
 
     void ax(double a, ColorSpinorField &x);
+    void ax(double a, std::vector<ColorSpinorField *> &x);
+    void ax(double *a, std::vector<ColorSpinorField *> &x); // not a true block-blas routine
 
     void axpbyz(double a, ColorSpinorField &x, double b, ColorSpinorField &y, ColorSpinorField &z);
 
@@ -122,13 +122,14 @@ namespace quda {
 
     double quadrupleCG3InitNorm(double a, ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v);
     double quadrupleCG3UpdateNorm(double a, double b, ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v);
-
+
+    double4 quadrupleEigCGUpdate(double a, double b, ColorSpinorField &x, ColorSpinorField &y, ColorSpinorField &z, ColorSpinorField &w, ColorSpinorField &v);
+
     // multi-blas kernels - defined in multi_blas.cu
-
     /**
        @brief Compute the block "axpy" with over the set of
-              ColorSpinorFields.  E.g., it computes y = x * a + y
-              The dimensions of a can be rectangular, e.g., the width of x and y need not be same.
+       ColorSpinorFields.  E.g., it computes y = x * a + y
+       The dimensions of a can be rectangular, e.g., the width of x and y need not be same.
        @param a[in] Matrix of real coefficients
        @param x[in] vector of input ColorSpinorFields
        @param y[in,out] vector of input/output ColorSpinorFields
@@ -680,5 +681,3 @@ et of
   } // namespace blas
 
 } // namespace quda
-
-#endif // _QUDA_BLAS_H

include/comm_quda.h

@@ -388,9 +388,12 @@ namespace quda
   void comm_wait(MsgHandle *mh);
   int comm_query(MsgHandle *mh);
 
+//  void comm_nonblocking_allreduce_array(MsgHandle *&mh, double *outdata, double *indata, size_t size);
+
   template <typename T> void comm_allreduce_sum(T &v);
   template <typename T> void comm_allreduce_max(T &v);
   template <typename T> void comm_allreduce_min(T &v);
+  template <typename T> void comm_nonblocking_allreduce_sum(MsgHandle *&mh, T &v);//?
 
   void comm_allreduce_int(int &data);
   void comm_allreduce_xor(uint64_t &data);

include/communicator_quda.h

@@ -743,6 +743,8 @@ namespace quda
   }
 
   void comm_allreduce_sum_array(double *data, size_t size);
+  //
+  void comm_nonblocking_allreduce_sum_array(MsgHandle *&mh, double *outdata, double *indata, size_t size);
 
   void comm_allreduce_max_array(double *data, size_t size);
 

include/complex_quda.h

@@ -1054,9 +1054,8 @@ lhs.real()*rhs.imag()+lhs.imag()*rhs.real());
     return polar(::sqrt(abs(z)),arg(z)/ValueType(2));
   }
 
-  template <>
-    __host__ __device__
-    inline complex<float> sqrt(const complex<float>& z){
+  template <> __host__ __device__ inline complex<float> sqrt(const complex<float> &z)
+  {
     return polar(::sqrtf(abs(z)),arg(z)/float(2));
   }
 

include/enum_quda_fortran.h

@@ -102,25 +102,27 @@
 #define QUDA_BICGSTAB_INVERTER 1
 #define QUDA_GCR_INVERTER 2
 #define QUDA_MR_INVERTER 3
-#define QUDA_SD_INVERTER 4
-#define QUDA_PCG_INVERTER 5
-#define QUDA_EIGCG_INVERTER 6
-#define QUDA_INC_EIGCG_INVERTER 7
-#define QUDA_GMRESDR_INVERTER 8
-#define QUDA_GMRESDR_PROJ_INVERTER 9
-#define QUDA_GMRESDR_SH_INVERTER 10
-#define QUDA_FGMRESDR_INVERTER 11
-#define QUDA_MG_INVERTER 12
-#define QUDA_BICGSTABL_INVERTER 13
-#define QUDA_CGNE_INVERTER 14
-#define QUDA_CGNR_INVERTER 15
-#define QUDA_CG3_INVERTER 16
-#define QUDA_CG3NE_INVERTER 17
-#define QUDA_CG3NR_INVERTER 18
-#define QUDA_CA_CG_INVERTER 19
-#define QUDA_CA_CGNE_INVERTER 20
-#define QUDA_CA_CGNR_INVERTER 21
-#define QUDA_CA_GCR_INVERTER 22
+#define QUDA_MPBICGSTAB_INVERTER 4
+#define QUDA_SD_INVERTER 5
+#define QUDA_PCG_INVERTER 7
+#define QUDA_MPCG_INVERTER 8
+#define QUDA_EIGCG_INVERTER 9
+#define QUDA_INC_EIGCG_INVERTER 10
+#define QUDA_GMRESDR_INVERTER 11
+#define QUDA_FGMRESDR_INVERTER 12
+#define QUDA_GCRODR_INVERTER 13
+#define QUDA_FGCRODR_INVERTER 14
+#define QUDA_MG_INVERTER 15
+#define QUDA_BICGSTABL_INVERTER 16
+#define QUDA_CGNE_INVERTER 17 
+#define QUDA_CGNR_INVERTER 18
+#define QUDA_CG3_INVERTER 19
+#define QUDA_CG3NE_INVERTER 20
+#define QUDA_CG3NR_INVERTER 21
+#define QUDA_CA_CG_INVERTER 22
+#define QUDA_CA_CGNE_INVERTER 23
+#define QUDA_CA_CGNR_INVERTER 24
+#define QUDA_CA_GCR_INVERTER 25
 #define QUDA_INVALID_INVERTER QUDA_INVALID_ENUM
 
 #define QudaEigType integer(4)