From 4a56c11c95e9f13670c906b353fe9360344eb913 Mon Sep 17 00:00:00 2001
From: Pjotr Prins
Date: Sat, 28 Nov 2020 14:39:05 +0000
Subject: Sane random generator handling

---
 src/param.cpp | 42 ++++++++++++++++++++++++++----------------
 1 file changed, 26 insertions(+), 16 deletions(-)

(limited to 'src/param.cpp')

diff --git a/src/param.cpp b/src/param.cpp
index e4781c4..db6d8d5 100644
--- a/src/param.cpp
+++ b/src/param.cpp
@@ -106,6 +106,10 @@ PARAM::PARAM(void)
       time_total(0.0), time_G(0.0), time_eigen(0.0), time_UtX(0.0),
       time_UtZ(0.0), time_opt(0.0), time_Omega(0.0) {}
 
+PARAM::~PARAM() {
+  gsl_rng_free(gsl_r);
+}
+
 // Read files: obtain ns_total, ng_total, ns_test, ni_test.
 void PARAM::ReadFiles(void) {
   string file_str;
@@ -817,6 +821,28 @@ void PARAM::CheckParam(void) {
     flag++;
   }
 
+  // Always set up random environment.
+  gsl_rng_env_setup(); // sets gsl_rng_default_seed
+  const gsl_rng_type *T = gsl_rng_default; // pick up environment GSL_RNG_SEED
+
+  if (randseed >= 0)
+    gsl_rng_default_seed = randseed; // CLI option used
+  else if (gsl_rng_default_seed == 0) { // by default we will randomize the seed
+    time_t rawtime;
+    time(&rawtime);
+    tm *ptm = gmtime(&rawtime);
+
+    gsl_rng_default_seed =
+      (unsigned)(ptm->tm_hour % 24 * 3600 + ptm->tm_min * 60 + ptm->tm_sec);
+  }
+  gsl_r = gsl_rng_alloc(T);
+
+  if (is_debug_mode()) {
+    printf ("GSL random generator type: %s; ", gsl_rng_name (gsl_r));
+    printf ("seed = %lu (option %li); ", gsl_rng_default_seed, randseed);
+    printf ("first value = %lu\n", gsl_rng_get (gsl_r));
+  }
+
   if (flag != 1 && a_mode != 15 && a_mode != 27 && a_mode != 28 &&
       a_mode != 43 && a_mode != 5 && a_mode != 61 && a_mode != 62 &&
       a_mode != 63 && a_mode != 66 && a_mode != 67) {
@@ -2015,22 +2041,6 @@ void PARAM::ProcessCvtPhen() {
            << "analyzed individuals. " << endl;
     } else {
 
-      // Set up random environment.
-      gsl_rng_env_setup();
-      gsl_rng *gsl_r;
-      const gsl_rng_type *gslType;
-      gslType = gsl_rng_default;
-      if (randseed < 0) {
-        time_t rawtime;
-        time(&rawtime);
-        tm *ptm = gmtime(&rawtime);
-
-        randseed = (unsigned)(ptm->tm_hour % 24 * 3600 + ptm->tm_min * 60 +
-                              ptm->tm_sec);
-      }
-      gsl_r = gsl_rng_alloc(gslType);
-      gsl_rng_set(gsl_r, randseed);
-
       // From ni_test, sub-sample ni_subsample.
       vector<size_t> a, b;
       for (size_t i = 0; i < ni_subsample; i++) {
-- 
cgit v1.2.3