include/gmx_random.h

   1
   2 #ifndef _GMX_RANDOM_H_
   3 #define _GMX_RANDOM_H_
   4
   5 #include <stdio.h>
   6 #include <types/simple.h>
   7
   8 #ifdef HAVE_CONFIG_H
   9 #include <config.h>
  10 #define gmx_inline inline
  11 #else
  12 #ifdef __GNUC__
  13 #define gmx_inline __inline
  14 #else
  15 #define inline
  16 #endif
  17 #endif
  18
  19
  20 /*! \brief Abstract datatype for a random number generator
  21  *
  22  * This is a handle to the full state of a random number generator.
  23  * You can not access anything inside the gmx_rng structure outside this
  24  * file.
  25  */
  26 typedef struct gmx_rng *
  27 gmx_rng_t;
  28
  29
  30 /*! \brief Create a new RNG, seeded from a single integer.
  31  *
  32  * If you dont want to pick a seed, just call it as
  33  * rng=gmx_rng_init(gmx_rng_make_seed()) to seed it from
  34  * the system time or a random device.
  35  *
  36  * \param seed Random seed, unsigned 32-bit integer.
  37  *
  38  * \return Reference to a random number generator, or NULL if there was an
  39  *         error.
  40  *
  41  * \threadsafe Yes.
  42  */
  43 gmx_rng_t
  44 gmx_rng_init(unsigned int seed);
  45
  46
  47 /*! \brief Generate a 'random' RNG seed.
  48  *
  49  * This routine tries to get a seed from /dev/random if present,
  50  * and if not it uses time-of-day and process id to generate one.
  51  *
  52  * \return 32-bit unsigned integer random seed.
  53  *
  54  * Tip: If you use this in your code, it is a good idea to write the
  55  * returned random seed to a logfile, so you can recreate the exact sequence
  56  * of random number if you need to reproduce your run later for one reason
  57  * or another.
  58  *
  59  * \threadsafe Yes.
  60  */
  61 unsigned int
  62 gmx_rng_make_seed(void);
  63
  64
  65 /*! \brief Initialize a RNG with 624 integers (>32 bits of entropy).
  66  *
  67  *  The Mersenne twister RNG used in Gromacs has an extremely long period,
  68  *  but when you only initialize it with a 32-bit integer there are only
  69  *  2^32 different possible sequences of number - much less than the generator
  70  *  is capable of.
  71  *
  72  *  If you really need the full entropy, this routine makes it possible to
  73  *  initialize the RNG with up to 624 32-bit integers, which will give you
  74  *  up to 2^19968 bits of entropy.
  75  *
  76  *  \param seed Array of unsigned integers to form a seed
  77  *  \param seed_length Number of integers in the array, up to 624 are used.
  78  *
  79  * \return Reference to a random number generator, or NULL if there was an
  80  *         error.
  81  *
  82  * \threadsafe Yes.
  83  */
  84 gmx_rng_t
  85 gmx_rng_init_array(unsigned int    seed[],
  86                    int             seed_length);
  87
  88
  89 /*! \brief Release resources of a RNG
  90  *
  91  *  This routine destroys a random number generator and releases all
  92  *  resources allocated by it.
  93  *
  94  *  \param rng Handle to random number generator previously returned by
  95  *                     gmx_rng_init() or gmx_rng_init_array().
  96  *
  97  * \threadsafe Function itself is threadsafe, but you should only destroy a
  98  *             certain RNG once (i.e. from one thread).
  99  */
 100 void
 101 gmx_rng_destroy(gmx_rng_t rng);
 102
 103
 104 /*! \brief Random 32-bit integer from a uniform distribution
 105  *
 106  *  This routine returns a random integer from the random number generator
 107  *  provided, and updates the state of that RNG.
 108  *
 109  *  \param rng Handle to random number generator previously returned by
 110  *                     gmx_rng_init() or gmx_rng_init_array().
 111  *
 112  *  \return 32-bit unsigned integer from a uniform distribution.
 113  *
 114  *  \threadsafe Function yes, input data no. You should not call this function
 115  *              from two different threads using the same RNG handle at the
 116  *              same time. For performance reasons we cannot lock the handle
 117  *              with a mutex every time we need a random number - that would
 118  *              slow the routine down a factor 2-5. There are two simple
 119  *              solutions: either use a mutex and lock it before calling
 120  *              the function, or use a separate RNG handle for each thread.
 121  */
 122 gmx_inline unsigned int
 123 gmx_rng_uniform_uint32(gmx_rng_t rng);
 124
 125
 126 /*! \brief Random gmx_real_t 0<=x<1 from a uniform distribution
 127  *
 128  *  This routine returns a random floating-point number from the
 129  *  random number generator provided, and updates the state of that RNG.
 130  *
 131  *  \param rng Handle to random number generator previously returned by
 132  *                     gmx_rng_init() or gmx_rng_init_array().
 133  *
 134  *  \return floating-point number 0<=x<1 from a uniform distribution.
 135  *
 136  *  \threadsafe Function yes, input data no. You should not call this function
 137  *              from two different threads using the same RNG handle at the
 138  *              same time. For performance reasons we cannot lock the handle
 139  *              with a mutex every time we need a random number - that would
 140  *              slow the routine down a factor 2-5. There are two simple
 141  *              solutions: either use a mutex and lock it before calling
 142  *              the function, or use a separate RNG handle for each thread.
 143  */
 144 gmx_inline real
 145 gmx_rng_uniform_real(gmx_rng_t rng);
 146
 147
 148 /*! \brief Random gmx_real_t from a gaussian distribution
 149  *
 150  *  This routine returns a random floating-point number from the
 151  *  random number generator provided, and updates the state of that RNG.
 152  *
 153  *  The Box-Muller algorithm is used to provide gaussian random numbers. This
 154  *  is not the fastest known algorithm for gaussian numbers, but in contrast
 155  *  to the alternatives it is very well studied and you can trust the returned
 156  *  random numbers to have good properties and no correlations.
 157  *
 158  *  \param rng Handle to random number generator previously returned by
 159  *                        gmx_rng_init() or gmx_rng_init_array().
 160  *
 161  *  \return Gaussian random floating-point number with average 0.0 and
 162  *          standard deviation 1.0. You can get any average/mean you want
 163  *          by first multiplying with the desired average and then adding
 164  *          the average you want.
 165  *
 166  *  \threadsafe Function yes, input data no. You should not call this function
 167  *              from two different threads using the same RNG handle at the
 168  *              same time. For performance reasons we cannot lock the handle
 169  *              with a mutex every time we need a random number - that would
 170  *              slow the routine down a factor 2-5. There are two simple
 171  *              solutions: either use a mutex and lock it before calling
 172  *              the function, or use a separate RNG handle for each thread.
 173  *
 174  *  It works perfectly to mix calls to get uniform and gaussian random numbers
 175  *  from the same generator, but since it will affect the sequence of returned
 176  *  numbers it is probably better to use separate random number generator
 177  *  structures.
 178  */
 179 real
 180 gmx_rng_gaussian_real(gmx_rng_t rng);
 181
 182
 183
 184 /* Return a new gaussian random number with expectation value
 185  * 0.0 and standard deviation 1.0. This routine uses a table
 186  * lookup for maximum speed.
 187  *
 188  * WARNING: The lookup table is 16k by default, which means
 189  *          the granularity of the random numbers is coarser
 190  *          than what you get from gmx_rng_gauss_real().
 191  *          In most cases this is no problem whatsoever,
 192  *          and it is particularly true for BD/SD integration.
 193  *          Note that you will NEVER get any really extreme
 194  *          numbers: the maximum absolute value returned is
 195  *          4.0255485.
 196  *
 197  * threadsafe: yes
 198  */
 199 real
 200 gmx_rng_gaussian_table(gmx_rng_t rng);
 201
 202
 203 #endif /* _GMX_RANDOM_H_ */
 204