-/* Copyright (C) 1999 Free Software Foundation, Inc.
+/* Copyright (C) 1999, 2000 Free Software Foundation, Inc.
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2, or (at your option)
* whether to permit this exception to apply to your modifications.
* If you do not wish that, delete this exception notice. */
+/* Software engineering face-lift by Greg J. Badros, 11-Dec-1999,
+ gjb@cs.washington.edu, http://www.cs.washington.edu/homes/gjb */
+
+
/* Author: Mikael Djurfeldt <djurfeldt@nada.kth.se> */
#include "_scm.h"
#include <stdio.h>
#include <math.h>
-#include "genio.h"
#include "smob.h"
#include "numbers.h"
#include "feature.h"
+#include "validate.h"
#include "random.h"
\f
*/
scm_rstate *
-scm_i_make_rstate (char *seed, int n)
+scm_c_make_rstate (char *seed, int n)
{
scm_rstate *state = malloc (scm_the_rng.rstate_size);
if (state == 0)
return state;
}
+scm_rstate *
+scm_c_default_rstate ()
+{
+ SCM state = SCM_CDR (scm_var_random_state);
+ SCM_ASSERT (SCM_RSTATEP (state),
+ state, "*random-state* contains bogus random state", 0);
+ return SCM_RSTATE (state);
+}
+
inline double
-scm_i_uniform01 (scm_rstate *state)
+scm_c_uniform01 (scm_rstate *state)
{
double x = (double) scm_the_rng.random_bits (state) / (double) 0xffffffffUL;
return ((x + (double) scm_the_rng.random_bits (state))
}
double
-scm_i_normal01 (scm_rstate *state)
+scm_c_normal01 (scm_rstate *state)
{
if (state->reserved0)
{
{
double r, a, n;
- r = sqrt (-2.0 * log (scm_i_uniform01 (state)));
- a = 2.0 * M_PI * scm_i_uniform01 (state);
+ r = sqrt (-2.0 * log (scm_c_uniform01 (state)));
+ a = 2.0 * M_PI * scm_c_uniform01 (state);
n = r * sin (a);
state->reserved1 = r * cos (a);
}
double
-scm_i_exp1 (scm_rstate *state)
+scm_c_exp1 (scm_rstate *state)
{
- return - log (scm_i_uniform01 (state));
+ return - log (scm_c_uniform01 (state));
}
unsigned char scm_masktab[256];
unsigned long
-scm_i_random (unsigned long m, scm_rstate *state)
+scm_c_random (scm_rstate *state, unsigned long m)
{
unsigned int r, mask;
mask = (m < 0x100
}
SCM
-scm_i_random_bignum (SCM m, scm_rstate *state)
+scm_c_random_bignum (scm_rstate *state, SCM m)
{
SCM b;
int i, nd;
static SCM
make_rstate (scm_rstate *state)
{
- SCM cell;
- SCM_NEWCELL (cell);
- SCM_ENTER_A_SECTION;
- SCM_SETCDR (cell, state);
- SCM_SETCAR (cell, scm_tc16_rstate);
- SCM_EXIT_A_SECTION;
- return cell;
-}
-
-static int
-print_rstate (SCM rstate, SCM port, scm_print_state *pstate)
-{
- scm_puts ("#<random-state ", port);
- scm_intprint ((long) SCM_RSTATE (rstate), 16, port);
- scm_putc ('>', port);
- return 1;
+ SCM_RETURN_NEWSMOB (scm_tc16_rstate, state);
}
static scm_sizet
return scm_the_rng.rstate_size;
}
-static scm_smobfuns rstate_smob = { 0, free_rstate, print_rstate, 0};
-
/*
* Scheme level interface.
*/
SCM_GLOBAL_VCELL_INIT (scm_var_random_state, "*random-state*", scm_seed_to_random_state (scm_makfrom0str ("URL:http://stat.fsu.edu/~geo/diehard.html")));
-SCM_PROC (s_random, "random", 1, 1, 0, scm_random);
-
-SCM
-scm_random (SCM n, SCM state)
+SCM_DEFINE (scm_random, "random", 1, 1, 0,
+ (SCM n, SCM state),
+ "Return a number in [0,N).\n"
+ "\n"
+ "Accepts a positive integer or real n and returns a \n"
+ "number of the same type between zero (inclusive) and \n"
+ "N (exclusive). The values returned have a uniform \n"
+ "distribution.\n"
+ "\n"
+ "The optional argument STATE must be of the type produced by\n"
+ "`seed->andom-state'. It defaults to the value of the variable\n"
+ "*random-state*. This object is used to maintain the state of\n"
+ "the pseudo-random-number generator and is altered as a side\n"
+ "effect of the random operation.\n"
+ "")
+#define FUNC_NAME s_scm_random
{
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state, SCM_ARG2, s_random);
+ SCM_VALIDATE_RSTATE (2,state);
if (SCM_INUMP (n))
{
unsigned long m = SCM_INUM (n);
- SCM_ASSERT (m > 0, n, SCM_ARG1, s_random);
- return SCM_MAKINUM (scm_i_random (m, SCM_RSTATE (state)));
+ SCM_ASSERT_RANGE (1,n,m > 0);
+ return SCM_MAKINUM (scm_c_random (SCM_RSTATE (state), m));
}
- SCM_ASSERT (SCM_NIMP (n), n, SCM_ARG1, s_random);
+ SCM_VALIDATE_NIM (1,n);
if (SCM_REALP (n))
- return scm_makdbl (SCM_REALPART (n) * scm_i_uniform01 (SCM_RSTATE (state)),
+ return scm_makdbl (SCM_REALPART (n) * scm_c_uniform01 (SCM_RSTATE (state)),
0.0);
- SCM_ASSERT (SCM_TYP16 (n) == scm_tc16_bigpos, n, SCM_ARG1, s_random);
- return scm_i_random_bignum (n, SCM_RSTATE (state));
+ SCM_VALIDATE_SMOB (1, n, big);
+ return scm_c_random_bignum (SCM_RSTATE (state), n);
}
+#undef FUNC_NAME
-SCM_PROC (s_copy_random_state, "copy-random-state", 0, 1, 0, scm_copy_random_state);
-
-SCM
-scm_copy_random_state (SCM state)
+SCM_DEFINE (scm_copy_random_state, "copy-random-state", 0, 1, 0,
+ (SCM state),
+ "Return a copy of the random state STATE.")
+#define FUNC_NAME s_scm_copy_random_state
{
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG1,
- s_copy_random_state);
+ SCM_VALIDATE_RSTATE (1,state);
return make_rstate (scm_the_rng.copy_rstate (SCM_RSTATE (state)));
}
+#undef FUNC_NAME
-SCM_PROC (s_seed_to_random_state, "seed->random-state", 1, 0, 0, scm_seed_to_random_state);
-
-SCM
-scm_seed_to_random_state (SCM seed)
+SCM_DEFINE (scm_seed_to_random_state, "seed->random-state", 1, 0, 0,
+ (SCM seed),
+ "Return a new random state using SEED.")
+#define FUNC_NAME s_scm_seed_to_random_state
{
if (SCM_NUMBERP (seed))
seed = scm_number_to_string (seed, SCM_UNDEFINED);
- SCM_ASSERT (SCM_NIMP (seed) && SCM_STRINGP (seed),
- seed,
- SCM_ARG1,
- s_seed_to_random_state);
- return make_rstate (scm_i_make_rstate (SCM_ROCHARS (seed),
+ SCM_VALIDATE_STRING (1,seed);
+ return make_rstate (scm_c_make_rstate (SCM_ROCHARS (seed),
SCM_LENGTH (seed)));
}
+#undef FUNC_NAME
-SCM_PROC (s_random_uniform, "random:uniform", 0, 1, 0, scm_random_uniform);
-
-SCM
-scm_random_uniform (SCM state)
+SCM_DEFINE (scm_random_uniform, "random:uniform", 0, 1, 0,
+ (SCM state),
+ "Returns a uniformly distributed inexact real random number in [0,1).")
+#define FUNC_NAME s_scm_random_uniform
{
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG1,
- s_random_uniform);
- return scm_makdbl (scm_i_uniform01 (SCM_RSTATE (state)), 0.0);
+ SCM_VALIDATE_RSTATE (1,state);
+ return scm_makdbl (scm_c_uniform01 (SCM_RSTATE (state)), 0.0);
}
+#undef FUNC_NAME
+
+SCM_DEFINE (scm_random_normal, "random:normal", 0, 1, 0,
+ (SCM state),
+ "Returns an inexact real in a normal distribution.\n"
+ "The distribution used has mean 0 and standard deviation 1.\n"
+ "For a normal distribution with mean m and standard deviation\n"
+ "d use @code{(+ m (* d (random:normal)))}.\n"
+ "")
+#define FUNC_NAME s_scm_random_normal
+{
+ if (SCM_UNBNDP (state))
+ state = SCM_CDR (scm_var_random_state);
+ SCM_VALIDATE_RSTATE (1,state);
+ return scm_makdbl (scm_c_normal01 (SCM_RSTATE (state)), 0.0);
+}
+#undef FUNC_NAME
+
+#ifdef HAVE_ARRAYS
static void
vector_scale (SCM v, double c)
return sum;
}
-SCM_PROC (s_random_solid_sphere_x, "random:solid-sphere!", 1, 1, 0, scm_random_solid_sphere_x);
-
/* For the uniform distribution on the solid sphere, note that in
* this distribution the length r of the vector has cumulative
* distribution r^n; i.e., u=r^n is uniform [0,1], so r can be
* generated as r=u^(1/n).
*/
-SCM
-scm_random_solid_sphere_x (SCM v, SCM state)
+SCM_DEFINE (scm_random_solid_sphere_x, "random:solid-sphere!", 1, 1, 0,
+ (SCM v, SCM state),
+ "Fills vect with inexact real random numbers\n"
+ "the sum of whose squares is less than 1.0.\n"
+ "Thinking of vect as coordinates in space of \n"
+ "dimension n = (vector-length vect), the coordinates \n"
+ "are uniformly distributed within the unit n-shere.\n"
+ "The sum of the squares of the numbers is returned.\n"
+ "")
+#define FUNC_NAME s_scm_random_solid_sphere_x
{
- SCM_ASSERT (SCM_NIMP (v)
- && (SCM_VECTORP (v) || SCM_TYP7 (v) == scm_tc7_dvect),
- v, SCM_ARG1, s_random_solid_sphere_x);
+ SCM_VALIDATE_VECTOR_OR_DVECTOR (1,v);
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG2,
- s_random_solid_sphere_x);
+ SCM_VALIDATE_RSTATE (2,state);
scm_random_normal_vector_x (v, state);
vector_scale (v,
- pow (scm_i_uniform01 (SCM_RSTATE (state)),
+ pow (scm_c_uniform01 (SCM_RSTATE (state)),
1.0 / SCM_LENGTH (v))
/ sqrt (vector_sum_squares (v)));
return SCM_UNSPECIFIED;
}
-
-SCM_PROC (s_random_hollow_sphere_x, "random:hollow-sphere!", 1, 1, 0, scm_random_hollow_sphere_x);
-
-SCM
-scm_random_hollow_sphere_x (SCM v, SCM state)
+#undef FUNC_NAME
+
+SCM_DEFINE (scm_random_hollow_sphere_x, "random:hollow-sphere!", 1, 1, 0,
+ (SCM v, SCM state),
+ "Fills vect with inexact real random numbers\n"
+ "the sum of whose squares is equal to 1.0.\n"
+ "Thinking of vect as coordinates in space of \n"
+ "dimension n = (vector-length vect), the coordinates\n"
+ "are uniformly distributed over the surface of the \n"
+ "unit n-shere.\n"
+ "")
+#define FUNC_NAME s_scm_random_hollow_sphere_x
{
- SCM_ASSERT (SCM_NIMP (v)
- && (SCM_VECTORP (v) || SCM_TYP7 (v) == scm_tc7_dvect),
- v, SCM_ARG1, s_random_solid_sphere_x);
+ SCM_VALIDATE_VECTOR_OR_DVECTOR (1,v);
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG2,
- s_random_hollow_sphere_x);
+ SCM_VALIDATE_RSTATE (2,state);
scm_random_normal_vector_x (v, state);
vector_scale (v, 1 / sqrt (vector_sum_squares (v)));
return SCM_UNSPECIFIED;
}
+#undef FUNC_NAME
-SCM_PROC (s_random_normal, "random:normal", 0, 1, 0, scm_random_normal);
-SCM
-scm_random_normal (SCM state)
-{
- if (SCM_UNBNDP (state))
- state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG1,
- s_random_normal);
- return scm_makdbl (scm_i_normal01 (SCM_RSTATE (state)), 0.0);
-}
-
-SCM_PROC (s_random_normal_vector_x, "random:normal-vector!", 1, 1, 0, scm_random_normal_vector_x);
-
-SCM
-scm_random_normal_vector_x (SCM v, SCM state)
+SCM_DEFINE (scm_random_normal_vector_x, "random:normal-vector!", 1, 1, 0,
+ (SCM v, SCM state),
+ "Fills vect with inexact real random numbers that are\n"
+ "independent and standard normally distributed\n"
+ "(i.e., with mean 0 and variance 1).\n"
+ "")
+#define FUNC_NAME s_scm_random_normal_vector_x
{
int n;
- SCM_ASSERT (SCM_NIMP (v)
- && (SCM_VECTORP (v) || SCM_TYP7 (v) == scm_tc7_dvect),
- v, SCM_ARG1, s_random_solid_sphere_x);
+ SCM_VALIDATE_VECTOR_OR_DVECTOR (1,v);
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG2,
- s_random_normal_vector_x);
+ SCM_VALIDATE_RSTATE (2,state);
n = SCM_LENGTH (v);
if (SCM_VECTORP (v))
while (--n >= 0)
- SCM_VELTS (v)[n] = scm_makdbl (scm_i_normal01 (SCM_RSTATE (state)), 0.0);
+ SCM_VELTS (v)[n] = scm_makdbl (scm_c_normal01 (SCM_RSTATE (state)), 0.0);
else
while (--n >= 0)
- ((double *) SCM_VELTS (v))[n] = scm_i_normal01 (SCM_RSTATE (state));
+ ((double *) SCM_VELTS (v))[n] = scm_c_normal01 (SCM_RSTATE (state));
return SCM_UNSPECIFIED;
}
+#undef FUNC_NAME
-SCM_PROC (s_random_exp, "random:exp", 0, 1, 0, scm_random_exp);
+#endif /* HAVE_ARRAYS */
-SCM
-scm_random_exp (SCM state)
+SCM_DEFINE (scm_random_exp, "random:exp", 0, 1, 0,
+ (SCM state),
+ "Returns an inexact real in an exponential distribution with mean 1.\n"
+ "For an exponential distribution with mean u use (* u (random:exp)).\n"
+ "")
+#define FUNC_NAME s_scm_random_exp
{
if (SCM_UNBNDP (state))
state = SCM_CDR (scm_var_random_state);
- SCM_ASSERT (SCM_NIMP (state) && SCM_RSTATEP (state),
- state,
- SCM_ARG1,
- s_random_exp);
- return scm_makdbl (scm_i_exp1 (SCM_RSTATE (state)), 0.0);
+ SCM_VALIDATE_RSTATE (1,state);
+ return scm_makdbl (scm_c_exp1 (SCM_RSTATE (state)), 0.0);
}
+#undef FUNC_NAME
void
scm_init_random ()
};
scm_the_rng = rng;
- scm_tc16_rstate = scm_newsmob (&rstate_smob);
+ scm_tc16_rstate = scm_make_smob_type_mfpe ("random-state", 0,
+ NULL, free_rstate, NULL, NULL);
for (m = 1; m <= 0x100; m <<= 1)
for (i = m >> 1; i < m; ++i)