public inbox for gsl-discuss@sourceware.org
 help / color / mirror / Atom feed
* RN generators Seed
@ 2003-08-26  7:50 Maura Edeweiss Monville
  2003-08-26 13:28 ` Jeff Spirko
  2003-08-26 13:44 ` Robert G. Brown
  0 siblings, 2 replies; 6+ messages in thread
From: Maura Edeweiss Monville @ 2003-08-26  7:50 UTC (permalink / raw)
  To: gsl-discuss

I'm currently using a RNG from the GNU Scietific Library in my simulation
I need to rerun my simulation many times with different seeds.
A friend of mine states the following:
"You can use the final random seed of a run as the starting random seed of
the next run (as you probably know, the random generators usually start
with a seed and then modify it in the process of generating random
numbers, so essentially the seed becomes a random number itself)"
I wonder how the final seed value can be obtained when using GSL Random  
Number Generators .... I only found the possibility to read and save the 
RNG "status". I do not know the relation between the RNG "status" and the
RNG "seed" ...     Can anyone help in this regards ?

Thank you
Maura

^ permalink raw reply	[flat|nested] 6+ messages in thread

* Re: RN generators Seed
  2003-08-26  7:50 RN generators Seed Maura Edeweiss Monville
@ 2003-08-26 13:28 ` Jeff Spirko
  2003-08-26 14:14   ` Nicolas Bock
  2003-08-26 13:44 ` Robert G. Brown
  1 sibling, 1 reply; 6+ messages in thread
From: Jeff Spirko @ 2003-08-26 13:28 UTC (permalink / raw)
  To: Maura Edeweiss Monville; +Cc: gsl-discuss

On Tue, Aug 26, 2003 at 12:15:27AM -0500, Maura Edeweiss Monville wrote:
> I'm currently using a RNG from the GNU Scietific Library in my simulation
> I need to rerun my simulation many times with different seeds.
> A friend of mine states the following:
> "You can use the final random seed of a run as the starting random seed of
> the next run (as you probably know, the random generators usually start
> with a seed and then modify it in the process of generating random
> numbers, so essentially the seed becomes a random number itself)"
> I wonder how the final seed value can be obtained when using GSL Random  
> Number Generators .... I only found the possibility to read and save the 
> RNG "status". I do not know the relation between the RNG "status" and the
> RNG "seed" ...     Can anyone help in this regards ?

I don't know the exact answer to your question, but I do have some
suggestions.  I often use the current time in seconds to seed an
RNG.  There is a slight problem with this.  If one simulation (or in
my case, setting up a simulation) takes less than a second, two runs
may use the same seed.  I ended up switching to using the bash
shell's $RANDOM variable.  It's limited to 16-bits, but it could be
combined with the time to get a better starting seed.

-- 
Jeff Spirko   spirko@lehigh.edu   spirko@yahoo.com   WD3V   |=>

The study of non-linear physics is like the study of non-elephant biology.

All theoretical chemistry is really physics;
and all theoretical chemists know it. -- Richard P. Feynman 

^ permalink raw reply	[flat|nested] 6+ messages in thread

* Re: RN generators Seed
  2003-08-26  7:50 RN generators Seed Maura Edeweiss Monville
  2003-08-26 13:28 ` Jeff Spirko
@ 2003-08-26 13:44 ` Robert G. Brown
  1 sibling, 0 replies; 6+ messages in thread
From: Robert G. Brown @ 2003-08-26 13:44 UTC (permalink / raw)
  To: Maura Edeweiss Monville; +Cc: gsl-discuss

[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #1: Type: TEXT/PLAIN; charset=X-UNKNOWN, Size: 2976 bytes --]

On Tue, 26 Aug 2003, Maura Edeweiss Monville wrote:

> I'm currently using a RNG from the GNU Scietific Library in my simulation
> I need to rerun my simulation many times with different seeds.
> A friend of mine states the following:
> "You can use the final random seed of a run as the starting random seed of
> the next run (as you probably know, the random generators usually start
> with a seed and then modify it in the process of generating random
> numbers, so essentially the seed becomes a random number itself)"
> I wonder how the final seed value can be obtained when using GSL Random  
> Number Generators .... I only found the possibility to read and save the 
> RNG "status". I do not know the relation between the RNG "status" and the
> RNG "seed" ...     Can anyone help in this regards ?

This still leaves you with several problems -- generating the FIRST seed
to start your computation(s) and the possibility of preserving any
residual correlations (however weak) between one run and the next.
Similar objections are often raised against the use of the clock as a
starting seed -- in PRINCIPLE a rng started with two distinct seeds
produces random sequences that are not correlated, but in practice they
may well do so.  A better way (in my opinion) to seed your rng is to use
the built in entropy based rng, /dev/random.

I actually have built a gsl-compatible call wrapper for /dev/random just
so I could test it in the context of testing other gsl random number
routines and use it exclusively to seed rng's.  Given the fairly trivial
wrapper, /dev/random could be added directly to the existing gsl library
easily enough, perhaps with a test to see if /dev/random exists before
actually registering the method (it does in all modern linuces but I'm
not sure about all OS's on which gsl runs:-).

Alternatively the direct utilization code is trivial -- open the "file"
/dev/random, read out a seed's worth of data, close it, you're done.  I
also have code for that if you like.

Note that /dev/random is slow as molasses (relying as it does on the
gradual accumulation of entropy from many "random" hardware sources) and
so is not a good general purpose rng for most purposes.  It is excellent
for generating seeds though, as the uniform deviates are of high enough
quality themselves to pass all the diehard tests I've subjected them to
so far (which is still not all of them, but they aren't egregiously
bit-correlated at any rate).

If there is any interest in the gsl-compatible call wrapper for
/dev/random, the code to just open and read it, or putting /dev/random
into the gsl rng interface directly on systems where it exists, let me
know.

    rgb

> 
> Thank you
> Maura
> 

Robert G. Brown	                       http://www.phy.duke.edu/~rgb/
Duke University Dept. of Physics, Box 90305
Durham, N.C. 27708-0305
Phone: 1-919-660-2567  Fax: 919-660-2525     email:rgb@phy.duke.edu



^ permalink raw reply	[flat|nested] 6+ messages in thread

* Re: RN generators Seed
  2003-08-26 13:28 ` Jeff Spirko
@ 2003-08-26 14:14   ` Nicolas Bock
  0 siblings, 0 replies; 6+ messages in thread
From: Nicolas Bock @ 2003-08-26 14:14 UTC (permalink / raw)
  To: gsl-discuss

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

On Tuesday 26 August 2003 09:28, Jeff Spirko wrote:

> I don't know the exact answer to your question, but I do have some
> suggestions.  I often use the current time in seconds to seed an
> RNG.  There is a slight problem with this.  If one simulation (or in
> my case, setting up a simulation) takes less than a second, two runs
> may use the same seed.  I ended up switching to using the bash
> shell's $RANDOM variable.  It's limited to 16-bits, but it could be
> combined with the time to get a better starting seed.

why don't you use /dev/random for the seed? This should give you completely 
random runs, no matter how close they are in time.

		nick
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.2.1 (GNU/Linux)

iD8DBQE/S2M0f1nNyq3gDhkRArkoAJ0TDY3m0cAKf6VGvFDi601mX0X98wCfbW1q
7udrcTxcuvO2LH3JXOjZ5JY=
=habO
-----END PGP SIGNATURE-----

^ permalink raw reply	[flat|nested] 6+ messages in thread

* Re: RN generators Seed
  2003-08-27 19:04 ` Robert G. Brown
@ 2003-08-27 19:11   ` Robert G. Brown
  0 siblings, 0 replies; 6+ messages in thread
From: Robert G. Brown @ 2003-08-27 19:11 UTC (permalink / raw)
  To: Maura Edelweiss Monville; +Cc: GSL Discussion list

On Wed, 27 Aug 2003, Robert G. Brown wrote:

>  /*
>   * and add the new one(s)...
>   */
>  if(i < N){
>    i = devrandom; /* save the index/type of the new generator */
>    types[i++] = (gsl_rng_dev_random);
>  }

Sheesh, adding bugs to my own code when I add the "illustration" -- this
should be:

>    devrandom = i; /* save the index/type of the new generator */

> 
> %< Snip snip snip ==================================================================
> 
> Thereafter you should be able to use /dev/random just like any other gsl rng,
> e.g.
>        gsl_rng *random;
>        random = gsl_rng_alloc (types[randnum]);

and this should be

>        random = gsl_rng_alloc (types[devrandom]);

Hopefully there aren't too many other typos...

   rgb

-- 
Robert G. Brown	                       http://www.phy.duke.edu/~rgb/
Duke University Dept. of Physics, Box 90305
Durham, N.C. 27708-0305
Phone: 1-919-660-2567  Fax: 919-660-2525     email:rgb@phy.duke.edu



^ permalink raw reply	[flat|nested] 6+ messages in thread

* Re: RN generators Seed
       [not found] <3F4B7547.9D138C75@artsci.wustl.edu>
@ 2003-08-27 19:04 ` Robert G. Brown
  2003-08-27 19:11   ` Robert G. Brown
  0 siblings, 1 reply; 6+ messages in thread
From: Robert G. Brown @ 2003-08-27 19:04 UTC (permalink / raw)
  To: Maura Edelweiss Monville; +Cc: GSL Discussion list

On Tue, 26 Aug 2003, Maura Edelweiss Monville wrote:

> Thank you very much. 
> Since I 'm using a RN generator from GSL I'm definitely interested in
> the gsl-compatible call wrapper for /dev/random. 
> I guess it is not in the same web site as the othe GSL goodies ... 
> I would appreciate being given the possibility to download it. 
> I develop on SuSE 8.2.

I'm including BOTH the code to generate a good random number seed from
/dev/random as a direct subroutine call (at the very end of this message) 
AND the wrapper code below.  The gsl rng god (if listening) is welcome to 
add /dev/random  to the gsl rng list of supported generators if they so desire.

[If you do you'll probably want to make error handling gsl-consistent (I'm
sloppy) and perhaps insert a stat/test for /dev/random's existence at the
point where the generator is actually added to the list of generators, as it
may not exist on some architectures where gsl runs.  You also MAY want to add
a gsl_rng_set_dev_random(random) function that just sets the seed from
/dev/random -- or not.]



To use the wrapper in your gsl routine, save dev_random.c below, add it to your 
makefile so it will compile and link.  Add the following into your header
file(s) (ordinarily they'd be in <gsl/gsl_rng.h>).

 const gsl_rng_type **types;    /* where all the rng types go */
 gsl_rng *random;               /* global gsl random number generator */

You'll then need a fragment like the following in your startup code to add the new type to the gsl
list of rng's BEFORE trying to use the new generator in any way:

%< Snip snip snip ==================================================================
 /*
  * Now add new RNG's onto the GSL types list.  Apparently this is all that is
  * needed -- they subsequently just "work" in the gsl call format.  This is code
  * taken from gsl types.c where N is apparently a hard-coded value.  This means
  * that we cannot CURRENTLY go over 100 generators total without e.g.
  * increasing N in the gsl sources.
  */
#define N 100
#define GSL_VAR

/* List new rng types to be added. */
 GSL_VAR const gsl_rng_type *gsl_rng_dev_random;

 int i,devrandom;

 /*
  * get to the end of the list of built in gsl generators
  */
 types = gsl_rng_types_setup ();
 i = 0;
 while(types[i] != NULL){
   i++;
   /* insert list code if desired */
 }
 /*
  * and add the new one(s)...
  */
 if(i < N){
   i = devrandom; /* save the index/type of the new generator */
   types[i++] = (gsl_rng_dev_random);
 }

%< Snip snip snip ==================================================================

Thereafter you should be able to use /dev/random just like any other gsl rng,
e.g.
       gsl_rng *random;
       random = gsl_rng_alloc (types[randnum]);
       random_max = gsl_rng_max(random);
       for(i = 0;i<100;i++){
         printf("%s returns a random int = %d\n",
           gsl_rng_name(random),    /* should return string "/dev/random", */
           gsl_rng_get(random));     /* should return a random int from /dev/random */
       }

(or the like, I'm not testing this code for typos:-) should return 100 random
ints from /dev/random.  If you run long enough loops of this, you can actually
see /dev/random "hang" from time to time waiting for entropy...

Really this is a pretty cool interface and it's remarkably easy to add your own rng's 
to it.

All code in this note is GPL v2, BTW and can be freely added to the GSL or
anybody's personal code, just to complete the formalities there.

   rgb

-- 
Robert G. Brown	                       http://www.phy.duke.edu/~rgb/
Duke University Dept. of Physics, Box 90305
Durham, N.C. 27708-0305
Phone: 1-919-660-2567  Fax: 919-660-2525     email:rgb@phy.duke.edu


%< Snip snip snip =============== dev_random.c ================================
/* dev_random
 * 
 * Copyright (C) 2003 Robert G. Brown
 * 
 * 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 of the License, or (at
 * your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 */

#include <gsl/gsl_rng.h>

/*
 * This is a wrapping of the /dev/random hardware rng
 */

static unsigned long int dev_random_get (void *vstate);
static double dev_random_get_double (void *vstate);
static void dev_random_set (void *vstate, unsigned long int s);

typedef struct
  {
    FILE *fp;
  }
dev_random_state_t;

static unsigned long int
dev_random_get (void *vstate)
{
  dev_random_state_t *state = (dev_random_state_t *) vstate;
  unsigned long int j;

  if(state->fp != NULL) {
    fread(&j,sizeof(j),1,state->fp);
    return j;
  } else {
    fprintf(stderr,"Error: /dev/random not open.  Exiting.\n");
    exit(0);
  }

}

static double
dev_random_get_double (void *vstate)
{
  return dev_random_get (vstate) / (double) UINT_MAX;
}

static void
dev_random_set (void *vstate, unsigned long int s)
{
  dev_random_state_t *state = (dev_random_state_t *) vstate;

 if ((state->fp = fopen("/dev/random","r")) == NULL) {
   fprintf(stderr,"Error: Cannot open /dev/random, exiting.\n");
   exit(0);
 }

 return;

}

static const gsl_rng_type dev_random_type =
{"/dev/random",			/* name */
 UINT_MAX,			/* RAND_MAX */
 0,				/* RAND_MIN */
 sizeof (dev_random_state_t),
 &dev_random_set,
 &dev_random_get,
 &dev_random_get_double};

const gsl_rng_type *gsl_rng_dev_random = &dev_random_type;

%< Snip snip snip =============== end of dev_random.c   ============================

%< Snip snip snip =============== random_seed.c  ===================================
/*
 *========================================================================
 * $Id: random_seed.c,v 1.2 2003/06/10 15:21:15 rgb Exp $
 *
 * See copyright in copyright.h and the accompanying file COPYING
 *========================================================================
 */

/*
 * *========================================================================
 * NOTE WELL:  This routine automagically falls back on setting from clock if
 * /dev/random is not available.  It also is missing some includes -- presuming
 * it will be put somewhere with either a shared header (where e.g. verbose is
 * defined) that has them or you'll add them by hand so it can standalone...
 * *========================================================================
 */

unsigned long int random_seed()
{

 unsigned long int seed;
 struct timeval tv;
 FILE *devrandom;

 if ((devrandom = fopen("/dev/random","r")) == NULL) {
   gettimeofday(&tv,0);
   seed = tv.tv_sec + tv.tv_usec;
   if(verbose) printf("Got seed %u from gettimeofday()\n",seed);
 } else {
   fread(&seed,sizeof(seed),1,devrandom);
   if(verbose) printf("Got seed %u from /dev/random\n",seed);
   fclose(devrandom);
 }

 return(seed);

}

%< Snip snip snip ==================================================================


^ permalink raw reply	[flat|nested] 6+ messages in thread

end of thread, other threads:[~2003-08-27 19:11 UTC | newest]

Thread overview: 6+ messages (download: mbox.gz / follow: Atom feed)
-- links below jump to the message on this page --
2003-08-26  7:50 RN generators Seed Maura Edeweiss Monville
2003-08-26 13:28 ` Jeff Spirko
2003-08-26 14:14   ` Nicolas Bock
2003-08-26 13:44 ` Robert G. Brown
     [not found] <3F4B7547.9D138C75@artsci.wustl.edu>
2003-08-27 19:04 ` Robert G. Brown
2003-08-27 19:11   ` Robert G. Brown

This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).