From e93cd3dee7ddfb7e688ed300adcfaacc5cbad7b6 Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Sun, 14 Dec 2014 11:41:50 +0100 Subject: change drand48_r to random_r --- packages/base/src/C/vector-aux.c | 21 ++++++++++++++------- 1 file changed, 14 insertions(+), 7 deletions(-) (limited to 'packages') diff --git a/packages/base/src/C/vector-aux.c b/packages/base/src/C/vector-aux.c index f8feb37..54cc413 100644 --- a/packages/base/src/C/vector-aux.c +++ b/packages/base/src/C/vector-aux.c @@ -700,15 +700,15 @@ int saveMatrix(char * file, char * format, KDMAT(a)){ //////////////////////////////////////////////////////////////////////////////// -inline double urandom(struct drand48_data * buffer) { - double res; - drand48_r(buffer,&res); - return res; +inline double urandom(struct random_data * buffer) { + int32_t res; + random_r(buffer,&res); + return (double)res/RAND_MAX; } // http://c-faq.com/lib/gaussian.html -double gaussrand(struct drand48_data *buffer, +double gaussrand(struct random_data *buffer, int *phase, double *pV1, double *pV2, double *pS) { double V1=*pV1, V2=*pV2, S=*pS; @@ -735,8 +735,15 @@ double gaussrand(struct drand48_data *buffer, } int random_vector(unsigned int seed, int code, DVEC(r)) { - struct drand48_data buffer; - srand48_r(seed,&buffer); + struct random_data buffer; + char random_state[128]; + memset(&buffer, 0, sizeof(struct random_data)); + memset(random_state, 0, sizeof(random_state)); + + initstate_r(seed,random_state,sizeof(random_state),&buffer); + // setstate_r(random_state,&buffer); + // srandom_r(seed,&buffer); + int phase = 0; double V1,V2,S; -- cgit v1.2.3 From 453305d525f1938456ef46abf5387725ba026cfe Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Sun, 14 Dec 2014 11:50:10 +0100 Subject: bump version and thanks --- packages/base/THANKS.md | 2 ++ packages/base/hmatrix.cabal | 2 +- 2 files changed, 3 insertions(+), 1 deletion(-) (limited to 'packages') diff --git a/packages/base/THANKS.md b/packages/base/THANKS.md index 90e6199..c220f16 100644 --- a/packages/base/THANKS.md +++ b/packages/base/THANKS.md @@ -177,3 +177,5 @@ module reorganization, monadic mapVectorM, and many other improvements. - "yongqli" reported the bug in randomVector (rand() is not thread-safe). +- "yongqli" and Kiwamu Ishikura reported that drand48_r() is not portable. + diff --git a/packages/base/hmatrix.cabal b/packages/base/hmatrix.cabal index 2900c42..c007509 100644 --- a/packages/base/hmatrix.cabal +++ b/packages/base/hmatrix.cabal @@ -1,5 +1,5 @@ Name: hmatrix -Version: 0.16.1.0 +Version: 0.16.1.1 License: BSD3 License-file: LICENSE Author: Alberto Ruiz -- cgit v1.2.3 From cc38f9a8e6d51b2e96037375d338546cc606640d Mon Sep 17 00:00:00 2001 From: Alberto Ruiz Date: Sun, 21 Dec 2014 12:07:38 +0100 Subject: temporary non thread safe randomVector workaround for OSX --- packages/base/hmatrix.cabal | 3 ++- packages/base/src/C/vector-aux.c | 57 ++++++++++++++++++++++++++++++++++++++++ 2 files changed, 59 insertions(+), 1 deletion(-) (limited to 'packages') diff --git a/packages/base/hmatrix.cabal b/packages/base/hmatrix.cabal index c007509..b3559a9 100644 --- a/packages/base/hmatrix.cabal +++ b/packages/base/hmatrix.cabal @@ -1,5 +1,5 @@ Name: hmatrix -Version: 0.16.1.1 +Version: 0.16.1.2 License: BSD3 License-file: LICENSE Author: Alberto Ruiz @@ -118,6 +118,7 @@ library if arch(i386) cc-options: -arch i386 frameworks: Accelerate + cpp-options: -DOSX if os(freebsd) extra-lib-dirs: /usr/local/lib diff --git a/packages/base/src/C/vector-aux.c b/packages/base/src/C/vector-aux.c index 54cc413..946d069 100644 --- a/packages/base/src/C/vector-aux.c +++ b/packages/base/src/C/vector-aux.c @@ -700,6 +700,61 @@ int saveMatrix(char * file, char * format, KDMAT(a)){ //////////////////////////////////////////////////////////////////////////////// +#ifdef OSX + +#pragma message "randomVector is not thread-safe in OSX" + +double gaussrand(int *phase, double *pV1, double *pV2, double *pS) +{ + double V1=*pV1, V2=*pV2, S=*pS; + double X; + + if(*phase == 0) { + do { + double U1 = (double)random() / (double)RAND_MAX; + double U2 = (double)random() / (double)RAND_MAX; + + V1 = 2 * U1 - 1; + V2 = 2 * U2 - 1; + S = V1 * V1 + V2 * V2; + } while(S >= 1 || S == 0); + + X = V1 * sqrt(-2 * log(S) / S); + } else + X = V2 * sqrt(-2 * log(S) / S); + + *phase = 1 - *phase; + *pV1=V1; *pV2=V2; *pS=S; + + return X; + +} + +int random_vector(unsigned int seed, int code, DVEC(r)) { + int phase = 0; + double V1,V2,S; + + int k; + switch (code) { + case 0: { // uniform + for (k=0; k