3 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 James Theiler, Brian Gough
5 * This program is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 3 of the License, or (at
8 * your option) any later version.
10 * This program is distributed in the hope that it will be useful, but
11 * WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 * General Public License for more details.
15 * You should have received a copy of the GNU General Public License
16 * along with this program; if not, write to the Free Software
17 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
24 #include <gsl/gsl_errno.h>
25 #include <gsl/gsl_rng.h>
28 gsl_rng_alloc (const gsl_rng_type * T)
31 gsl_rng *r = (gsl_rng *) malloc (sizeof (gsl_rng));
35 GSL_ERROR_VAL ("failed to allocate space for rng struct",
39 r->state = malloc (T->size);
43 free (r); /* exception in constructor, avoid memory leak */
45 GSL_ERROR_VAL ("failed to allocate space for rng state",
51 gsl_rng_set (r, gsl_rng_default_seed); /* seed the generator */
57 gsl_rng_memcpy (gsl_rng * dest, const gsl_rng * src)
59 if (dest->type != src->type)
61 GSL_ERROR ("generators must be of the same type", GSL_EINVAL);
64 memcpy (dest->state, src->state, src->type->size);
70 gsl_rng_clone (const gsl_rng * q)
72 gsl_rng *r = (gsl_rng *) malloc (sizeof (gsl_rng));
76 GSL_ERROR_VAL ("failed to allocate space for rng struct",
80 r->state = malloc (q->type->size);
84 free (r); /* exception in constructor, avoid memory leak */
86 GSL_ERROR_VAL ("failed to allocate space for rng state",
92 memcpy (r->state, q->state, q->type->size);
98 gsl_rng_set (const gsl_rng * r, unsigned long int seed)
100 (r->type->set) (r->state, seed);
103 #ifndef HIDE_INLINE_STATIC
105 gsl_rng_get (const gsl_rng * r)
107 return (r->type->get) (r->state);
111 gsl_rng_uniform (const gsl_rng * r)
113 return (r->type->get_double) (r->state);
117 gsl_rng_uniform_pos (const gsl_rng * r)
122 x = (r->type->get_double) (r->state) ;
129 /* Note: to avoid integer overflow in (range+1) we work with scale =
130 range/n = (max-min)/n rather than scale=(max-min+1)/n, this reduces
131 efficiency slightly but avoids having to check for the out of range
132 value. Note that range is typically O(2^32) so the addition of 1
133 is negligible in most usage. */
136 gsl_rng_uniform_int (const gsl_rng * r, unsigned long int n)
138 unsigned long int offset = r->type->min;
139 unsigned long int range = r->type->max - offset;
140 unsigned long int scale;
143 if (n > range || n == 0)
145 GSL_ERROR_VAL ("invalid n, either 0 or exceeds maximum value of generator",
153 k = (((r->type->get) (r->state)) - offset) / scale;
162 gsl_rng_max (const gsl_rng * r)
168 gsl_rng_min (const gsl_rng * r)
174 gsl_rng_name (const gsl_rng * r)
176 return r->type->name;
180 gsl_rng_size (const gsl_rng * r)
182 return r->type->size;
186 gsl_rng_state (const gsl_rng * r)
192 gsl_rng_print_state (const gsl_rng * r)
195 unsigned char *p = (unsigned char *) (r->state);
196 const size_t n = r->type->size;
198 for (i = 0; i < n; i++)
200 /* FIXME: we're assuming that a char is 8 bits */
201 printf ("%.2x", *(p + i));
207 gsl_rng_free (gsl_rng * r)