From 6c9c68c93dfdb143422f26868abab97ffd9d727a Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Fran=C3=A7ois=20Coppens?= <mail@fmgjcoppens.nl>
Date: Thu, 20 May 2021 19:21:59 +0200
Subject: [PATCH] Added  LAPACK option to test_h5.cpp to compare residuals and
 number of operations.

---
 Makefile               | 14 ++++++++------
 include/SM_Helpers.hpp |  3 +++
 src/SM_Helpers.cpp     | 18 ++++++++++++++++++
 tests/test_h5.cpp      |  3 +++
 4 files changed, 32 insertions(+), 6 deletions(-)

diff --git a/Makefile b/Makefile
index f144ad7..9cf057d 100644
--- a/Makefile
+++ b/Makefile
@@ -21,11 +21,13 @@ else
     $(error No valid compiler environment set in $$ENV. \
 	First run: $$ source smvars.sh {intel | llvm | gnu})
 endif
-CXXFLAGS = $(OPT) $(ARCH) $(DEBUG) -fPIC $(THRESHOLD)
+CXXFLAGS = $(OPT) $(ARCH) $(DEBUG) $(THRESHOLD) -fPIC
 FFLAGS = $(CXXFLAGS)
 H5CXX = h5c++
 H5CXXFLAGS = $(CXXFLAGS)
 FLIBS = -lstdc++
+LFLAGS = -mkl=sequential
+H5LFLAGS = -L${MKLROOT}/lib/intel64 -Wl,--no-as-needed -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl
 
 ## Includes and dependencies
 INCLUDE = -I $(INC_DIR)/
@@ -93,16 +95,16 @@ $(OBJ_DIR)/%.o: $(TST_DIR)/%.f90 | $(OBJ_DIR)
 
 #### LINKING
 $(BIN_DIR)/cMaponiA3_test_3x3_3: $(OBJ_DIR)/cMaponiA3_test_3x3_3.o $(DEPS_CXX) | $(BIN_DIR)
-	$(CXX) -o $@ $^
+	$(CXX) $(LFLAGS) -o $@ $^
 
 $(BIN_DIR)/test_h5: $(OBJ_DIR)/test_h5.o $(DEPS_CXX) | $(BIN_DIR)
-	$(H5CXX) -o $@ $^
+	$(H5CXX) $(H5LFLAGS) -o $@ $^
 
 $(BIN_DIR)/fMaponiA3_test_3x3_3: $(DEPS_F) $(OBJ_DIR)/fMaponiA3_test_3x3_3.o  | $(BIN_DIR)
-	$(FC) $(FLIBS) -o $@ $^
+	$(FC) $(LFLAGS) $(FLIBS) -o $@ $^
 
 $(BIN_DIR)/fMaponiA3_test_4x4_2: $(DEPS_F) $(OBJ_DIR)/fMaponiA3_test_4x4_2.o  | $(BIN_DIR)
-	$(FC) $(FLIBS) -o $@ $^
+	$(FC) $(LFLAGS) $(FLIBS) -o $@ $^
 
 $(BIN_DIR)/QMCChem_dataset_test: $(DEPS_F) $(OBJ_DIR)/QMCChem_dataset_test.o  | $(BIN_DIR)
-	$(FC) $(FLIBS) -o $@ $^
+	$(FC) $(LFLAGS) $(FLIBS) -o $@ $^
diff --git a/include/SM_Helpers.hpp b/include/SM_Helpers.hpp
index 0840dcc..334af69 100644
--- a/include/SM_Helpers.hpp
+++ b/include/SM_Helpers.hpp
@@ -4,6 +4,7 @@
 #include <cstring>
 #include <iostream>
 #include <string>
+#include <mkl_lapacke.h>
 
 // #define DEBUG
 
@@ -18,6 +19,8 @@ void selectLargestDenominator(unsigned int l, unsigned int N_updates,
                               unsigned int *Updates_index, unsigned int *p,
                               double ***ylk);
 
+lapack_int inverse(double *A, unsigned n);
+
 template <typename T> void showScalar(T scalar, std::string name) {
   std::cout << name << " = " << scalar << std::endl << std::endl;
 }
diff --git a/src/SM_Helpers.cpp b/src/SM_Helpers.cpp
index 157bac4..5bc9e1c 100644
--- a/src/SM_Helpers.cpp
+++ b/src/SM_Helpers.cpp
@@ -39,3 +39,21 @@ void selectLargestDenominator(unsigned int l, unsigned int N_updates,
   }
   Switch(p, l, lbar);
 }
+
+// Inplace inverse n x n matrix A.
+// returns:
+//   ret = 0 on success
+//   ret < 0 illegal argument value
+//   ret > 0 singular matrix
+lapack_int inverse(double *A, unsigned n) {
+  int ipiv[n + 1];
+  lapack_int ret;
+
+  ret = LAPACKE_dgetrf(LAPACK_COL_MAJOR, n, n, A, n, ipiv);
+
+  if (ret != 0)
+    return ret;
+
+  ret = LAPACKE_dgetri(LAPACK_COL_MAJOR, n, A, n, ipiv);
+  return ret;
+}
diff --git a/tests/test_h5.cpp b/tests/test_h5.cpp
index 2724512..7cf1ec5 100644
--- a/tests/test_h5.cpp
+++ b/tests/test_h5.cpp
@@ -81,6 +81,9 @@ int test_cycle(H5File file, int cycle, std::string version, double tolerance) {
     SM3(slater_inverse, dim, nupdates, u, col_update_index);
   } else if (version == "sm4") {
     SM4(slater_inverse, dim, nupdates, u, col_update_index);
+  } else if (version == "lapack") {
+    memcpy(slater_inverse, slater_matrix, dim * dim * sizeof(double));
+    inverse(slater_inverse, dim);
   } else {
     std::cerr << "Unknown version " << version << std::endl;
     exit(1);