lloda pushed a commit to branch wip-exception-truncate
in repository guile.
commit f4d7038255ab721142334096bcf3c820e63c38d5
Author: Daniel Llorens <daniel.llor...@bluewin.ch>
Date: Thu Mar 30 14:29:59 2017 +0200
Support general arrays in random:hollow-sphere!
* libguile/random.c (vector_scale_x, vector_sum_squares): Handle general
rank-1 #t or 'f64 arrays.
* test-suite/tests/random.test: Add tests for random:hollow-sphere!.
---
libguile/random.c | 105 ++++++++++++++++++++++++-------------------
test-suite/tests/random.test | 37 ++++++++++++++-
2 files changed, 94 insertions(+), 48 deletions(-)
diff --git a/libguile/random.c b/libguile/random.c
index a8ad075..58791af 100644
--- a/libguile/random.c
+++ b/libguile/random.c
@@ -498,66 +498,77 @@ SCM_DEFINE (scm_random_normal, "random:normal", 0, 1, 0,
}
#undef FUNC_NAME
+/* FIXME see scm_array_handle_ref for handling possible overflow */
static void
vector_scale_x (SCM v, double c)
{
- size_t n;
- if (scm_is_vector (v))
- {
- n = SCM_SIMPLE_VECTOR_LENGTH (v);
- while (n-- > 0)
- SCM_REAL_VALUE (SCM_SIMPLE_VECTOR_REF (v, n)) *= c;
- }
- else
- {
- /* must be a f64vector. */
- scm_t_array_handle handle;
- size_t i, len;
- ssize_t inc;
- double *elts;
-
- elts = scm_f64vector_writable_elements (v, &handle, &len, &inc);
+ scm_t_array_handle handle;
+ scm_t_array_dim const * dims;
+ ssize_t i, inc, ubnd;
- for (i = 0; i < len; i++, elts += inc)
- *elts *= c;
-
- scm_array_handle_release (&handle);
+ scm_array_get_handle (v, &handle);
+ dims = scm_array_handle_dims (&handle);
+ if (1 == scm_array_handle_rank (&handle))
+ {
+ ubnd = dims[0].ubnd;
+ inc = dims[0].inc;
+ if (handle.element_type == SCM_ARRAY_ELEMENT_TYPE_F64)
+ {
+ double *elts = (double *)(handle.writable_elements) + handle.base;
+ for (i = dims[0].lbnd; i <= ubnd; ++i, elts += inc)
+ *elts *= c;
+ return;
+ }
+ else if (handle.element_type == SCM_ARRAY_ELEMENT_TYPE_SCM)
+ {
+ SCM *elts = (SCM *)(handle.writable_elements) + handle.base;
+ for (i = dims[0].lbnd; i <= ubnd; ++i, elts += inc)
+ SCM_REAL_VALUE (*elts) *= c;
+ return;
+ }
}
+ scm_array_handle_release (&handle);
+ scm_misc_error (NULL, "must be a rank-1 array of type #t or 'f64",
scm_list_1 (v));
}
+/* FIXME see scm_array_handle_ref for handling possible overflow */
static double
vector_sum_squares (SCM v)
{
double x, sum = 0.0;
- size_t n;
- if (scm_is_vector (v))
- {
- n = SCM_SIMPLE_VECTOR_LENGTH (v);
- while (n-- > 0)
- {
- x = SCM_REAL_VALUE (SCM_SIMPLE_VECTOR_REF (v, n));
- sum += x * x;
- }
- }
- else
- {
- /* must be a f64vector. */
- scm_t_array_handle handle;
- size_t i, len;
- ssize_t inc;
- const double *elts;
-
- elts = scm_f64vector_elements (v, &handle, &len, &inc);
-
- for (i = 0; i < len; i++, elts += inc)
- {
- x = *elts;
- sum += x * x;
- }
+ scm_t_array_handle handle;
+ scm_t_array_dim const * dims;
+ ssize_t i, inc, ubnd;
- scm_array_handle_release (&handle);
+ scm_array_get_handle (v, &handle);
+ dims = scm_array_handle_dims (&handle);
+ if (1 == scm_array_handle_rank (&handle))
+ {
+ ubnd = dims[0].ubnd;
+ inc = dims[0].inc;
+ if (handle.element_type == SCM_ARRAY_ELEMENT_TYPE_F64)
+ {
+ const double *elts = (const double *)(handle.elements) + handle.base;
+ for (i = dims[0].lbnd; i <= ubnd; ++i, elts += inc)
+ {
+ x = *elts;
+ sum += x * x;
+ }
+ return sum;
+ }
+ else if (handle.element_type == SCM_ARRAY_ELEMENT_TYPE_SCM)
+ {
+ const SCM *elts = (const SCM *)(handle.elements) + handle.base;
+ for (i = dims[0].lbnd; i <= ubnd; ++i, elts += inc)
+ {
+ x = SCM_REAL_VALUE (*elts);
+ sum += x * x;
+ }
+ return sum;
+ }
}
- return sum;
+ scm_array_handle_release (&handle);
+ scm_misc_error (NULL, "must be an array of type #t or 'f64", scm_list_1 (v));
}
/* For the uniform distribution on the solid sphere, note that in
diff --git a/test-suite/tests/random.test b/test-suite/tests/random.test
index ab20b58..1492651 100644
--- a/test-suite/tests/random.test
+++ b/test-suite/tests/random.test
@@ -20,7 +20,8 @@
#:use-module ((system base compile) #:select (compile))
#:use-module (test-suite lib)
#:use-module (srfi srfi-4)
- #:use-module (srfi srfi-4 gnu))
+ #:use-module (srfi srfi-4 gnu)
+ #:use-module ((ice-9 control) #:select (let/ec)))
; see strings.test, arrays.test.
(define exception:wrong-type-arg
@@ -53,3 +54,37 @@
(random:normal-vector! b (random-state-from-platform))
(random:normal-vector! c (random-state-from-platform))
(and (not (equal? a b)) (not (equal? a c)))))))
+
+;;;
+;;; random:hollow-sphere!
+;;;
+
+(with-test-prefix "random:hollow-sphere!"
+
+ (define (sqr a)
+ (* a a))
+ (define (norm a)
+ (sqrt (+ (sqr (array-ref a 0)) (sqr (array-ref a 1)) (sqr (array-ref a
2)))))
+ (define double-eps 1e-15)
+
+ (pass-if "non uniform"
+ (let ((a (transpose-array (make-array 0. 3 10) 1 0)))
+ (let/ec exit
+ (array-slice-for-each 1
+ (lambda (a)
+ (random:hollow-sphere! a)
+ (if (> (magnitude (- 1 (norm a))) double-eps) (exit #f)))
+ a)
+ #t)))
+
+ (pass-if "uniform (f64)"
+ (let ((a (transpose-array (make-array 0. 3 10) 1 0)))
+ (let/ec exit
+ (array-slice-for-each 1
+ (lambda (a)
+ (random:hollow-sphere! a)
+ (if (> (magnitude (- 1 (norm a))) double-eps) (exit #f)))
+ a)
+ #t))))
+
+
