From cdf407bd7994dbe41a952bf29cacc1a2ca9c722e Mon Sep 17 00:00:00 2001
From: Pjotr Prins
Date: Thu, 6 Sep 2018 11:59:03 +0000
Subject: More debugging and a performance check

---
 src/debug.cpp  | 4 ++--
 src/debug.h    | 3 +++
 src/lapack.cpp | 2 +-
 src/lmm.cpp    | 3 ++-
 src/mvlmm.cpp  | 6 +++++-
 src/param.cpp  | 3 +--
 src/version.h  | 6 +++---
 7 files changed, 17 insertions(+), 10 deletions(-)

(limited to 'src')

diff --git a/src/debug.cpp b/src/debug.cpp
index 2f1f861..9f9f862 100644
--- a/src/debug.cpp
+++ b/src/debug.cpp
@@ -137,7 +137,7 @@ void enable_segfpe() {
   if (!is_fpe_check_mode() || is_legacy_mode()) return;
   #ifdef __GNUC__
     #if defined(__x86_64__)
-       debug_msg("enable segfpe hardware floating point error detection");
+  // debug_msg("enable segfpe hardware floating point error detection");
        signal(SIGFPE, sighandler);
        feenableexcept(FE_INVALID | FE_DIVBYZERO | FE_OVERFLOW | FE_UNDERFLOW);
     #endif
@@ -148,7 +148,7 @@ void disable_segfpe() {
   if (!is_fpe_check_mode() || is_legacy_mode()) return;
   #ifdef __GNUC__
     #if defined(__x86_64__)
-      debug_msg("disable segfpe");
+  // debug_msg("disable segfpe");
       fedisableexcept(FE_INVALID | FE_DIVBYZERO | FE_OVERFLOW | FE_UNDERFLOW);
     #endif
   #endif
diff --git a/src/debug.h b/src/debug.h
index e83c813..a665a37 100644
--- a/src/debug.h
+++ b/src/debug.h
@@ -115,6 +115,9 @@ inline void fail_msg(std::string msg) {
   std::raise(SIGINT); // keep stack trace for gdb
 }
 
+#define info_msg(msg) cerr << "**** INFO: " << msg << "." << endl;
+#define msg(msg) info_msg(msg);
+
 #if defined NDEBUG
   #define __SHOW_FUNC __func__
 
diff --git a/src/lapack.cpp b/src/lapack.cpp
index 125e5a0..79d49fd 100644
--- a/src/lapack.cpp
+++ b/src/lapack.cpp
@@ -319,7 +319,7 @@ void LUDecomp(gsl_matrix *LU, gsl_permutation *p, int *signum) {
 void LUInvert(const gsl_matrix *LU, const gsl_permutation *p, gsl_matrix *ret_inverse) {
   // debug_msg("entering");
   auto det = LULndet(LU);
-  assert(det != 1.0);
+  enforce_msg(det != 1.0,"LU determinant is zero -> LU is not invertable");
 
   enforce_gsl(gsl_linalg_LU_invert(LU, p, ret_inverse));
 }
diff --git a/src/lmm.cpp b/src/lmm.cpp
index d08c90e..80372ee 100644
--- a/src/lmm.cpp
+++ b/src/lmm.cpp
@@ -296,7 +296,8 @@ void CalcPab(const size_t n_cvt, const size_t e_mode, const gsl_vector *Hi_eval,
   write(Hi_eval,"Hi_eval");
   write(Uab,"Uab");
   // write(ab,"ab");
-  assert(!has_nan(Hi_eval));
+  if (is_check_mode())
+    assert(!has_nan(Hi_eval));
   assert(!has_nan(Uab));
   // assert(!has_nan(ab));
 
diff --git a/src/mvlmm.cpp b/src/mvlmm.cpp
index d877302..faea21e 100644
--- a/src/mvlmm.cpp
+++ b/src/mvlmm.cpp
@@ -321,12 +321,15 @@ double CalcQi(const gsl_vector *eval, const gsl_vector *D_l,
 }
 
 // xHiy=\sum_{k=1}^n x_k\otimes ((delta_k*Dl+I)^{-1}Ul^TVe^{-1/2}y.
+//
+// FIXME: mvlmm spends a massive amount of time here
 void CalcXHiY(const gsl_vector *eval, const gsl_vector *D_l,
               const gsl_matrix *X, const gsl_matrix *UltVehiY,
               gsl_vector *xHiy) {
+  debug_msg("enter");
   size_t n_size = eval->size, c_size = X->size1, d_size = D_l->size;
 
-  gsl_vector_set_zero(xHiy);
+  // gsl_vector_set_zero(xHiy);
 
   double x, delta, dl, y, d;
   for (size_t i = 0; i < d_size; i++) {
@@ -342,6 +345,7 @@ void CalcXHiY(const gsl_vector *eval, const gsl_vector *D_l,
       gsl_vector_set(xHiy, j * d_size + i, d);
     }
   }
+  // debug_msg("exit");
 
   return;
 }
diff --git a/src/param.cpp b/src/param.cpp
index 6537379..de0c257 100644
--- a/src/param.cpp
+++ b/src/param.cpp
@@ -1942,8 +1942,7 @@ void PARAM::CheckCvt() {
     indicator_cvt.clear();
     n_cvt = 1;
   } else if (flag_ipt == 0) {
-    cout << "no intercept term is found in the cvt file. "
-         << "a column of 1s is added." << endl;
+    info_msg("no intercept term is found in the cvt file: a column of 1s is added");
     for (vector<int>::size_type i = 0; i < indicator_idv.size(); ++i) {
       if (indicator_idv[i] == 0 || indicator_cvt[i] == 0) {
         continue;
diff --git a/src/version.h b/src/version.h
index 518a9f3..5afd64d 100644
--- a/src/version.h
+++ b/src/version.h
@@ -1,5 +1,5 @@
 // version.h generated by GEMMA scripts/gen_version_info.sh
-#define GEMMA_VERSION "0.98-pre2"
-#define GEMMA_DATE "2018-07-14"
+#define GEMMA_VERSION "0.98-beta1"
+#define GEMMA_DATE "2018-09-06"
 #define GEMMA_YEAR "2018"
-#define GEMMA_PROFILE "/gnu/store/9ahrb1swr06kjm2gr2zg0fsyvps3xqgz-profile"
+#define GEMMA_PROFILE "/gnu/store/8h46wb8hi8n492i2vlbfcz4dm44xwlpc-profile"
-- 
cgit v1.2.3