Header file for libtpcrand. More...

#include "tpcclibConfig.h"
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <sys/time.h>
#include <stdint.h>
#include <inttypes.h>
#include "tpcextensions.h"

Data Structures
struct	MERTWI

Macros
#define	TPCCLIB_MERTWI_NN 312
#define	TPCCLIB_MERTWI_A UINT64_C(0xB5026F5AA96619E9)

Functions
uint32_t	mertwiSeed32 (void)
	Make uint32_t seed for pseudo-random number generators.
uint64_t	mertwiSeed64 (void)
	Make uint64_t seed for pseudo-random number generators.
void	mertwiInit (MERTWI *mt)
void	mertwiInitWithSeed64 (MERTWI *mt, uint64_t seed)
	Initialize the state vector mt[] inside data structure for Mersenne Twister MT19937 pseudo-random number generator using given seed.
void	mertwiInitByArray64 (MERTWI *mt, uint64_t init_key[], uint64_t key_length)
	Initialize the state vector mt[] inside data structure for Mersenne Twister MT19937 pseudo-random number generator using given array.
uint64_t	mertwiRandomInt64 (MERTWI *mt)
	Generate a random number on [0, 2^64-1]-interval using Mersenne Twister MT19937.
int64_t	mertwiRandomInt63 (MERTWI *mt)
	Generate a random number on [0, 2^63-1]-interval using Mersenne Twister MT19937.
double	mertwiRandomDouble1 (MERTWI *mt)
	Generate a 64-bit double precision floating point pseudo-random number in the range of [0,1] with uniform distribution using Mersenne Twister MT19937.
double	mertwiRandomDouble2 (MERTWI *mt)
	Generate a 64-bit double precision floating point pseudo-random number in the range of [0,1) with uniform distribution using Mersenne Twister MT19937.
double	mertwiRandomDouble3 (MERTWI *mt)
	Generate a 64-bit double precision floating point pseudo-random number in the range of (0,1) with uniform distribution using Mersenne Twister MT19937.
double	mertwiRandomGaussian (MERTWI *mt)
	Generate a 64-bit double precision floating point pseudo-random number with normal (Gaussian) distribution with zero mean 0 and SD 1.
double	mertwiRandomExponential (MERTWI *mt, double mean)
	Generate pseudo-random number with exponential distribution and specified mean.
int	mertwiRandomBetween (MERTWI mt, unsigned int nr, double d, double low, double up, int type)
	Fill the given double array with random numbers with uniform distribution between the specified limits.
unsigned int	drandSeed (short int seed)
	Make and optionally set the seed for rand(), drand, drandRange, and drandGaussian().
double	drand ()
int	drandRange (unsigned int nr, double *d, double low, double up, int type)
double	drandGaussian ()
	Get pseudo-random number with normal (Gaussian) distribution with mean 0 and SD 1.
double	drandExponential (double mean)
	Get pseudo-random number with exponential distribution.
void	randomShuffle (int array, unsigned int n, MERTWI mt)
	Random shuffle.
void	randomShuffleUnsigned (unsigned int array, unsigned int n, MERTWI mt)
	Random shuffle.
void	randomPermutation (int array, unsigned int n, int a, MERTWI mt)
	Random permutation for an array of signed integers.
void	randomPermutationUnsigned (unsigned int array, unsigned int n, unsigned int a, MERTWI mt)
	Random permutation for an array of unsigned integers.
int	haltonElement (const int i, const int dim, double *r)
	Calculation of an element of quasi-random low-discrepancy Halton sequence.
int	haltonPrime (const int N)
	Subfunction for calculation of quasi-random Halton sequence.

Detailed Description

Header file for libtpcrand.

Header file for library libtpcrand.

Author: Vesa Oikonen

Copyright: (c) Turku PET Centre

Definition in file tpcrand.h.

Macro Definition Documentation

◆ TPCCLIB_MERTWI_A

#define TPCCLIB_MERTWI_A UINT64_C(0xB5026F5AA96619E9)

Mersenne Twister required constant

Definition at line 29 of file tpcrand.h.

Referenced by mertwiInit(), and mertwiRandomInt64().

◆ TPCCLIB_MERTWI_NN

#define TPCCLIB_MERTWI_NN 312

Mersenne Twister state vector length

Definition at line 27 of file tpcrand.h.

Referenced by mertwiInit().

Function Documentation

◆ drand()

double drand ( )

extern

Alternative function to rand() which returns a double precision floating point number in the range of [0,1] with uniform distribution.

Precondition: Uses rand(), therefore set seed for a new series of pseudo-random numbers; to produce truly random numbers (not just pseudo-random), do srand(time(0)) before calling this function. If no seed is set, then value 1 is used as default seed by rand().

See also: mertwiRandomDouble1, drandRange, drandSeed

Author: Vesa Oikonen

Returns: Random double value in the range [0,1].

Definition at line 66 of file gaussdev.c.

{
  double d, s;
  s=1.0/(1.0+RAND_MAX);
  do {
    d = ( ( s*rand() + rand() )*s + rand() ) * s;
  } while(d>=1.0);
  return d;
}

Referenced by drandExponential(), drandGaussian(), drandRange(), nloptPowellBrent(), and nloptRandomPoint().

◆ drandExponential()

double drandExponential ( double mean )

extern

Get pseudo-random number with exponential distribution.

Precondition: Set seed for a new series of pseudo-random numbers; to produce truly random numbers (not just pseudo-random), do srand(time(0)) before calling this function. If no seed is set, then value 1 is used as default seed by rand().

See also: mertwiRandomExponential, drandRange, drandSeed

Returns: Returns the random non-negative double value with exponential distribution and given mean.

Parameters

mean	Mean of the exponential distribution.

Definition at line 178 of file gaussdev.c.

  {
  double r;
  do {
    r=drand(); 
  } while(r==0.0);
  return(-mean*log(r)); 
}

Referenced by nloptIATGO().

◆ drandGaussian()

double drandGaussian ( )

extern

Get pseudo-random number with normal (Gaussian) distribution with mean 0 and SD 1.

Applies the polar form of Box-Müller transform to produce pseudo-random numbers with Gaussian (normal) distribution which has a zero mean and standard deviation of one.

Box GEP, Muller ME. A note on the generation of random normal deviates. Annals of Mathematical Statistics, Volume 29, Issue 2, 1958, 610-611. Available from JSTOR https://www.jstor.org/

Precondition: Set seed for a new series of pseudo-random numbers; to produce truly random numbers (not just pseudo-random), do srand(time(0)) before calling this function. If no seed is set, then value 1 is used as default seed by rand().

See also: mertwiRandomGaussian, drandRange, drandSeed

Returns: Returns the pseudo-random double value with normal distribution with zero mean and standard deviation 1.

Definition at line 144 of file gaussdev.c.

{
  static int ready=0;
  static double dev;
  double fac, rsq, a, b;
 
  /* If we don't have deviate already, then we'll have to make one */
  if(!ready) {
    do {
      a = 2.*drand() - 1.0; 
      b = 2.*drand() - 1.0; 
      rsq = a*a + b*b; 
    } while (rsq>1.0 || rsq==0.0);
   
    fac = sqrt(-2.0*log(rsq)/rsq);
    dev=a*fac; ready=1; 
    return(b*fac); 
  } else { /* dev is ready so return it */
    ready=0;
    return(dev);
  }
}

Referenced by nloptGaussianPoint().

◆ drandRange()

int drandRange	(	unsigned int	nr,
		double *	d,
		double	low,
		double	up,
		int	type )

extern

Fill the given double array with random numbers with uniform distribution between the specified limits.

With uniform distribution, the SD=(up-low)/sqrt(12), and CV=(up-low)/(sqrt(3)*(low+up)).

Precondition: Set seed for a new series of pseudo-random numbers; to produce truly random numbers (not just pseudo-random), do srand(time(0)) before calling this function. If no seed is set, then value 1 is used as default seed by rand().

See also: mertwiRandomBetween, drand, drandSeed

Author: Vesa Oikonen

Returns: 0 when successful, otherwise <> 0.

Parameters

nr	Nr of values in double array.
d	Pointer to pre-allocated double array.
low	Lower limit for random values.
up	Upper limit for random values.
type	Distribution: 0=even, 1=square-root transformation.

Definition at line 90 of file gaussdev.c.

  {
  unsigned int i;
  double dif, v, stl, stu;
 
  if(nr<1) return 0;
  if(d==NULL || type<0 || type>1) return 1;
 
  dif=up-low; if(dif<0.0) return 2;
  if(dif==0.0) {
    for(i=0; i<nr; i++) d[i]=low;
    return 0;
  }
 
  if(type==0) {
    for(i=0; i<nr; i++) d[i] = drand()*dif + low;
  } else if(type==1) {
    stl=copysign(sqrt(fabs(low)),low); if(!isnormal(stl)) stl=0.0;
    stu=copysign(sqrt(fabs(up)), up); if(!isnormal(stu)) stu=0.0;
    dif=stu-stl;
    for(i=0; i<nr; i++) {v=drand()*dif+stl; d[i]=copysign(v*v, v);}
  }
 
  return 0;
}

◆ drandSeed()

unsigned int drandSeed ( short int seed )

extern

Make and optionally set the seed for rand(), drand, drandRange, and drandGaussian().

Uses microseconds from the computer clock and process ID to reduce the chance of getting the same seed for simultaneously executing program threads and instances.

See also: mertwiSeed64

Returns: Returns the seed for srand().

Parameters

seed	Also sets seed with srand (1) or not (0).

Definition at line 27 of file gaussdev.c.

  {
  unsigned int li;
#if defined HAVE_TIMESPEC_GET
  struct timespec ts;
  timespec_get(&ts, TIME_UTC);
  //li=(unsigned int)ts.tv_sec+(unsigned int)ts.tv_nsec+(unsigned int)getpid();
  li=((ts.tv_sec % 10000)*523 ^ ts.tv_nsec*10) ^ ((getpid() % 1000)*983);
#elif defined HAVE_CLOCK_GETTIME
  struct timespec ts;
  clock_gettime(CLOCK_REALTIME, &ts);
  //li=(unsigned int)ts.tv_sec+(unsigned int)ts.tv_nsec+(unsigned int)getpid();
  li=((ts.tv_sec % 10000)*523 ^ ts.tv_nsec*10) ^ ((getpid() % 1000)*983);
#elif defined HAVE_GETTIMEOFDAY
  struct timeval tv;
  gettimeofday(&tv, 0);
  li=((tv.tv_sec % 10000)*523 ^ tv.tv_usec*13) ^ ((getpid() % 1000)*983);
#else
  li=(unsigned int)time(NULL)+(unsigned int)getpid();
#endif
  li+=(unsigned int)rand();
  if(seed) srand(li);
  return(li);
}

◆ haltonElement()

int haltonElement	(	const int	i,
		const int	dim,
		double *	r )

extern

Calculation of an element of quasi-random low-discrepancy Halton sequence.

Returns: 0 when successful, otherwise <> 0.

Parameters

i	The index of the element to be calculated; 0<=I.
dim	The spatial dimension (size of element; number of parameters).
r	An array for the calculated element, of length DIM, containing values in range [0,1].

Definition at line 38 of file halton.c.

  {
  if(i<0 || dim<1 || r==NULL) return(1);
 
  int t[dim];
  for(int j=0; j<dim; j++) t[j]=1+i; // Plus one by VO, because otherwise results are set to zero for i=0
 
  double prime_inv[dim];
  for(int j=0; j<dim; j++) prime_inv[j] = 1.0 / (double)haltonPrime(1+j);
 
  for(int j=0; j<dim; j++) r[j]=0.0;
  while(1) {
    for(int j=0; j<dim; j++) {
      int d = (t[j] % haltonPrime(1+j));
      r[j] += (double)d * prime_inv[j];
      prime_inv[j] /= (double)haltonPrime(1+j);
      t[j] /= haltonPrime(1+j);
    }
    int tmax=t[0];
    for(int j=1; j<dim; j++) if(t[j]>tmax) tmax=t[j];
    if(tmax<1) break;
  }
  return(0);
}

Referenced by nloptSimplexMS().

◆ haltonPrime()

int haltonPrime ( const int N )

extern

Subfunction for calculation of quasi-random Halton sequence.

See also: haltonElement

Returns: Returns the required prime number (>0), or 0 in case of an error.

Parameters

N	The index of the desired prime number; 0 <= N <= iMAX. N<0 returns iMAX, the index of the largest prime available. N=0 is legal, returning prime 1. N>iMAX returns zero, which should be interpreted as an error.

Definition at line 75 of file halton.c.

  {
#define iMAX 1600
 
  if(N<0) return(iMAX);
  if(N==0) return(1);
  if(N>iMAX) return(0);
  const int npvec[iMAX] = {
        2,    3,    5,    7,   11,   13,   17,   19,   23,   29,
       31,   37,   41,   43,   47,   53,   59,   61,   67,   71,
       73,   79,   83,   89,   97,  101,  103,  107,  109,  113,
      127,  131,  137,  139,  149,  151,  157,  163,  167,  173,
      179,  181,  191,  193,  197,  199,  211,  223,  227,  229,
      233,  239,  241,  251,  257,  263,  269,  271,  277,  281,
      283,  293,  307,  311,  313,  317,  331,  337,  347,  349,
      353,  359,  367,  373,  379,  383,  389,  397,  401,  409,
      419,  421,  431,  433,  439,  443,  449,  457,  461,  463,
      467,  479,  487,  491,  499,  503,  509,  521,  523,  541,
      547,  557,  563,  569,  571,  577,  587,  593,  599,  601,
      607,  613,  617,  619,  631,  641,  643,  647,  653,  659,
      661,  673,  677,  683,  691,  701,  709,  719,  727,  733,
      739,  743,  751,  757,  761,  769,  773,  787,  797,  809,
      811,  821,  823,  827,  829,  839,  853,  857,  859,  863,
      877,  881,  883,  887,  907,  911,  919,  929,  937,  941,
      947,  953,  967,  971,  977,  983,  991,  997, 1009, 1013,
     1019, 1021, 1031, 1033, 1039, 1049, 1051, 1061, 1063, 1069,
     1087, 1091, 1093, 1097, 1103, 1109, 1117, 1123, 1129, 1151,
     1153, 1163, 1171, 1181, 1187, 1193, 1201, 1213, 1217, 1223,
     1229, 1231, 1237, 1249, 1259, 1277, 1279, 1283, 1289, 1291,
     1297, 1301, 1303, 1307, 1319, 1321, 1327, 1361, 1367, 1373,
     1381, 1399, 1409, 1423, 1427, 1429, 1433, 1439, 1447, 1451,
     1453, 1459, 1471, 1481, 1483, 1487, 1489, 1493, 1499, 1511,
     1523, 1531, 1543, 1549, 1553, 1559, 1567, 1571, 1579, 1583,
     1597, 1601, 1607, 1609, 1613, 1619, 1621, 1627, 1637, 1657,
     1663, 1667, 1669, 1693, 1697, 1699, 1709, 1721, 1723, 1733,
     1741, 1747, 1753, 1759, 1777, 1783, 1787, 1789, 1801, 1811,
     1823, 1831, 1847, 1861, 1867, 1871, 1873, 1877, 1879, 1889,
     1901, 1907, 1913, 1931, 1933, 1949, 1951, 1973, 1979, 1987,
     1993, 1997, 1999, 2003, 2011, 2017, 2027, 2029, 2039, 2053,
     2063, 2069, 2081, 2083, 2087, 2089, 2099, 2111, 2113, 2129,
     2131, 2137, 2141, 2143, 2153, 2161, 2179, 2203, 2207, 2213,
     2221, 2237, 2239, 2243, 2251, 2267, 2269, 2273, 2281, 2287,
     2293, 2297, 2309, 2311, 2333, 2339, 2341, 2347, 2351, 2357,
     2371, 2377, 2381, 2383, 2389, 2393, 2399, 2411, 2417, 2423,
     2437, 2441, 2447, 2459, 2467, 2473, 2477, 2503, 2521, 2531,
     2539, 2543, 2549, 2551, 2557, 2579, 2591, 2593, 2609, 2617,
     2621, 2633, 2647, 2657, 2659, 2663, 2671, 2677, 2683, 2687,
     2689, 2693, 2699, 2707, 2711, 2713, 2719, 2729, 2731, 2741,
     2749, 2753, 2767, 2777, 2789, 2791, 2797, 2801, 2803, 2819,
     2833, 2837, 2843, 2851, 2857, 2861, 2879, 2887, 2897, 2903,
     2909, 2917, 2927, 2939, 2953, 2957, 2963, 2969, 2971, 2999,
     3001, 3011, 3019, 3023, 3037, 3041, 3049, 3061, 3067, 3079,
     3083, 3089, 3109, 3119, 3121, 3137, 3163, 3167, 3169, 3181,
     3187, 3191, 3203, 3209, 3217, 3221, 3229, 3251, 3253, 3257,
     3259, 3271, 3299, 3301, 3307, 3313, 3319, 3323, 3329, 3331,
     3343, 3347, 3359, 3361, 3371, 3373, 3389, 3391, 3407, 3413,
     3433, 3449, 3457, 3461, 3463, 3467, 3469, 3491, 3499, 3511,
     3517, 3527, 3529, 3533, 3539, 3541, 3547, 3557, 3559, 3571,
     3581, 3583, 3593, 3607, 3613, 3617, 3623, 3631, 3637, 3643,
     3659, 3671, 3673, 3677, 3691, 3697, 3701, 3709, 3719, 3727,
     3733, 3739, 3761, 3767, 3769, 3779, 3793, 3797, 3803, 3821,
     3823, 3833, 3847, 3851, 3853, 3863, 3877, 3881, 3889, 3907,
     3911, 3917, 3919, 3923, 3929, 3931, 3943, 3947, 3967, 3989,
     4001, 4003, 4007, 4013, 4019, 4021, 4027, 4049, 4051, 4057,
     4073, 4079, 4091, 4093, 4099, 4111, 4127, 4129, 4133, 4139,
     4153, 4157, 4159, 4177, 4201, 4211, 4217, 4219, 4229, 4231,
     4241, 4243, 4253, 4259, 4261, 4271, 4273, 4283, 4289, 4297,
     4327, 4337, 4339, 4349, 4357, 4363, 4373, 4391, 4397, 4409,
     4421, 4423, 4441, 4447, 4451, 4457, 4463, 4481, 4483, 4493,
     4507, 4513, 4517, 4519, 4523, 4547, 4549, 4561, 4567, 4583,
     4591, 4597, 4603, 4621, 4637, 4639, 4643, 4649, 4651, 4657,
     4663, 4673, 4679, 4691, 4703, 4721, 4723, 4729, 4733, 4751,
     4759, 4783, 4787, 4789, 4793, 4799, 4801, 4813, 4817, 4831,
     4861, 4871, 4877, 4889, 4903, 4909, 4919, 4931, 4933, 4937,
     4943, 4951, 4957, 4967, 4969, 4973, 4987, 4993, 4999, 5003,
     5009, 5011, 5021, 5023, 5039, 5051, 5059, 5077, 5081, 5087,
     5099, 5101, 5107, 5113, 5119, 5147, 5153, 5167, 5171, 5179,
     5189, 5197, 5209, 5227, 5231, 5233, 5237, 5261, 5273, 5279,
     5281, 5297, 5303, 5309, 5323, 5333, 5347, 5351, 5381, 5387,
     5393, 5399, 5407, 5413, 5417, 5419, 5431, 5437, 5441, 5443,
     5449, 5471, 5477, 5479, 5483, 5501, 5503, 5507, 5519, 5521,
     5527, 5531, 5557, 5563, 5569, 5573, 5581, 5591, 5623, 5639,
     5641, 5647, 5651, 5653, 5657, 5659, 5669, 5683, 5689, 5693,
     5701, 5711, 5717, 5737, 5741, 5743, 5749, 5779, 5783, 5791,
     5801, 5807, 5813, 5821, 5827, 5839, 5843, 5849, 5851, 5857,
     5861, 5867, 5869, 5879, 5881, 5897, 5903, 5923, 5927, 5939,
     5953, 5981, 5987, 6007, 6011, 6029, 6037, 6043, 6047, 6053,
     6067, 6073, 6079, 6089, 6091, 6101, 6113, 6121, 6131, 6133,
     6143, 6151, 6163, 6173, 6197, 6199, 6203, 6211, 6217, 6221,
     6229, 6247, 6257, 6263, 6269, 6271, 6277, 6287, 6299, 6301,
     6311, 6317, 6323, 6329, 6337, 6343, 6353, 6359, 6361, 6367,
     6373, 6379, 6389, 6397, 6421, 6427, 6449, 6451, 6469, 6473,
     6481, 6491, 6521, 6529, 6547, 6551, 6553, 6563, 6569, 6571,
     6577, 6581, 6599, 6607, 6619, 6637, 6653, 6659, 6661, 6673,
     6679, 6689, 6691, 6701, 6703, 6709, 6719, 6733, 6737, 6761,
     6763, 6779, 6781, 6791, 6793, 6803, 6823, 6827, 6829, 6833,
     6841, 6857, 6863, 6869, 6871, 6883, 6899, 6907, 6911, 6917,
     6947, 6949, 6959, 6961, 6967, 6971, 6977, 6983, 6991, 6997,
     7001, 7013, 7019, 7027, 7039, 7043, 7057, 7069, 7079, 7103,
     7109, 7121, 7127, 7129, 7151, 7159, 7177, 7187, 7193, 7207,
     7211, 7213, 7219, 7229, 7237, 7243, 7247, 7253, 7283, 7297,
     7307, 7309, 7321, 7331, 7333, 7349, 7351, 7369, 7393, 7411,
     7417, 7433, 7451, 7457, 7459, 7477, 7481, 7487, 7489, 7499,
     7507, 7517, 7523, 7529, 7537, 7541, 7547, 7549, 7559, 7561,
     7573, 7577, 7583, 7589, 7591, 7603, 7607, 7621, 7639, 7643,
     7649, 7669, 7673, 7681, 7687, 7691, 7699, 7703, 7717, 7723,
     7727, 7741, 7753, 7757, 7759, 7789, 7793, 7817, 7823, 7829,
     7841, 7853, 7867, 7873, 7877, 7879, 7883, 7901, 7907, 7919,
     7927, 7933, 7937, 7949, 7951, 7963, 7993, 8009, 8011, 8017,
     8039, 8053, 8059, 8069, 8081, 8087, 8089, 8093, 8101, 8111,
     8117, 8123, 8147, 8161, 8167, 8171, 8179, 8191, 8209, 8219,
     8221, 8231, 8233, 8237, 8243, 8263, 8269, 8273, 8287, 8291,
     8293, 8297, 8311, 8317, 8329, 8353, 8363, 8369, 8377, 8387,
     8389, 8419, 8423, 8429, 8431, 8443, 8447, 8461, 8467, 8501,
     8513, 8521, 8527, 8537, 8539, 8543, 8563, 8573, 8581, 8597,
     8599, 8609, 8623, 8627, 8629, 8641, 8647, 8663, 8669, 8677,
     8681, 8689, 8693, 8699, 8707, 8713, 8719, 8731, 8737, 8741,
     8747, 8753, 8761, 8779, 8783, 8803, 8807, 8819, 8821, 8831,
     8837, 8839, 8849, 8861, 8863, 8867, 8887, 8893, 8923, 8929,
     8933, 8941, 8951, 8963, 8969, 8971, 8999, 9001, 9007, 9011,
     9013, 9029, 9041, 9043, 9049, 9059, 9067, 9091, 9103, 9109,
     9127, 9133, 9137, 9151, 9157, 9161, 9173, 9181, 9187, 9199,
     9203, 9209, 9221, 9227, 9239, 9241, 9257, 9277, 9281, 9283,
     9293, 9311, 9319, 9323, 9337, 9341, 9343, 9349, 9371, 9377,
     9391, 9397, 9403, 9413, 9419, 9421, 9431, 9433, 9437, 9439,
     9461, 9463, 9467, 9473, 9479, 9491, 9497, 9511, 9521, 9533,
     9539, 9547, 9551, 9587, 9601, 9613, 9619, 9623, 9629, 9631,
     9643, 9649, 9661, 9677, 9679, 9689, 9697, 9719, 9721, 9733,
     9739, 9743, 9749, 9767, 9769, 9781, 9787, 9791, 9803, 9811,
     9817, 9829, 9833, 9839, 9851, 9857, 9859, 9871, 9883, 9887,
     9901, 9907, 9923, 9929, 9931, 9941, 9949, 9967, 9973,10007,
    10009,10037,10039,10061,10067,10069,10079,10091,10093,10099,
    10103,10111,10133,10139,10141,10151,10159,10163,10169,10177,
    10181,10193,10211,10223,10243,10247,10253,10259,10267,10271,
    10273,10289,10301,10303,10313,10321,10331,10333,10337,10343,
    10357,10369,10391,10399,10427,10429,10433,10453,10457,10459,
    10463,10477,10487,10499,10501,10513,10529,10531,10559,10567,
    10589,10597,10601,10607,10613,10627,10631,10639,10651,10657,
    10663,10667,10687,10691,10709,10711,10723,10729,10733,10739,
    10753,10771,10781,10789,10799,10831,10837,10847,10853,10859,
    10861,10867,10883,10889,10891,10903,10909,10937,10939,10949,
    10957,10973,10979,10987,10993,11003,11027,11047,11057,11059,
    11069,11071,11083,11087,11093,11113,11117,11119,11131,11149,
    11159,11161,11171,11173,11177,11197,11213,11239,11243,11251,
    11257,11261,11273,11279,11287,11299,11311,11317,11321,11329,
    11351,11353,11369,11383,11393,11399,11411,11423,11437,11443,
    11447,11467,11471,11483,11489,11491,11497,11503,11519,11527,
    11549,11551,11579,11587,11593,11597,11617,11621,11633,11657,
    11677,11681,11689,11699,11701,11717,11719,11731,11743,11777,
    11779,11783,11789,11801,11807,11813,11821,11827,11831,11833,
    11839,11863,11867,11887,11897,11903,11909,11923,11927,11933,
    11939,11941,11953,11959,11969,11971,11981,11987,12007,12011,
    12037,12041,12043,12049,12071,12073,12097,12101,12107,12109,
    12113,12119,12143,12149,12157,12161,12163,12197,12203,12211,
    12227,12239,12241,12251,12253,12263,12269,12277,12281,12289,
    12301,12323,12329,12343,12347,12373,12377,12379,12391,12401,
    12409,12413,12421,12433,12437,12451,12457,12473,12479,12487,
    12491,12497,12503,12511,12517,12527,12539,12541,12547,12553,
    12569,12577,12583,12589,12601,12611,12613,12619,12637,12641,
    12647,12653,12659,12671,12689,12697,12703,12713,12721,12739,
    12743,12757,12763,12781,12791,12799,12809,12821,12823,12829,
    12841,12853,12889,12893,12899,12907,12911,12917,12919,12923,
    12941,12953,12959,12967,12973,12979,12983,13001,13003,13007,
    13009,13033,13037,13043,13049,13063,13093,13099,13103,13109,
    13121,13127,13147,13151,13159,13163,13171,13177,13183,13187,
    13217,13219,13229,13241,13249,13259,13267,13291,13297,13309,
    13313,13327,13331,13337,13339,13367,13381,13397,13399,13411,
    13417,13421,13441,13451,13457,13463,13469,13477,13487,13499
  };
  return(npvec[N-1]);
# undef iMAX
}

Referenced by haltonElement().

◆ mertwiInit()

void mertwiInit ( MERTWI * mt )

extern

Prepare the data structure of Mersenne Twister MT19937 for usage. Do not call any mertwi* functions before calling this function!

See also: mertwiSeed64, mertwiInitWithSeed64, mertwiInitByArray64

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 28 of file mertwi.c.

  {
  mt->n=TPCCLIB_MERTWI_NN; // Length of mt array
  mt->m=mt->n/2; // N/2
  mt->a=TPCCLIB_MERTWI_A;
  mt->um=UINT64_C(0xFFFFFFFF80000000);
  mt->lm=UINT64_C(0x7FFFFFFF);
  mt->mti=mt->n+1; // means that mt[] is not initialized
}

Referenced by nloptMPSO().

◆ mertwiInitByArray64()

void mertwiInitByArray64	(	MERTWI *	mt,
		uint64_t	init_key[],
		uint64_t	key_length )

extern

Initialize the state vector mt[] inside data structure for Mersenne Twister MT19937 pseudo-random number generator using given array.

Call either this or mertwiInitWithSeed64 before generating random numbers with MT19937 using functions.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before this.

See also: mertwiSeed64, mertwiInitWithSeed64, mertwiInit

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.
init_key	The array for initializing keys.
key_length	Length of initialization array init_key[].

Definition at line 113 of file mertwi.c.

  {
  unsigned int i, j;
  uint64_t k;
  mertwiInitWithSeed64(mt, UINT64_C(19650218));
  i=1; j=0; if(mt->n>key_length) k=mt->n; else k=key_length;
  for(; k; k--) {
    mt->mt[i] = 
      (mt->mt[i] ^ ((mt->mt[i-1] ^ (mt->mt[i-1]>>62)) * UINT64_C(3935559000370003845)))
      + init_key[j] + j;
    i++; j++;
    if(i>=mt->n) {mt->mt[0]=mt->mt[mt->n-1]; i=1;}
    if(j>=key_length) j=0;
  }
  for(k=mt->n-1; k; k--) {
    mt->mt[i]=
      (mt->mt[i] ^ ((mt->mt[i-1] ^ (mt->mt[i-1] >> 62)) * UINT64_C(2862933555777941757)))
       - i;
    i++;
    if(i>=mt->n) {mt->mt[0]=mt->mt[mt->n-1]; i=1;}
  }
  mt->mt[0] = UINT64_C(1) << 63; /* MSB is 1; assuring non-zero initial array */ 
}

◆ mertwiInitWithSeed64()

void mertwiInitWithSeed64	(	MERTWI *	mt,
		uint64_t	seed )

extern

Initialize the state vector mt[] inside data structure for Mersenne Twister MT19937 pseudo-random number generator using given seed.

Call either this or mertwiInitByArray64 before generating random numbers with MT19937 using functions.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before this.

See also: mertwiSeed64, mertwiInitByArray64, mertwiInit

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.
seed	Seed, for example from mertwiSeed64().

Definition at line 94 of file mertwi.c.

  {
  mt->mt[0]=seed;
  for(mt->mti=1; mt->mti<mt->n; mt->mti++) 
    mt->mt[mt->mti]=
      (UINT64_C(6364136223846793005)*(mt->mt[mt->mti-1] ^ (mt->mt[mt->mti-1]>>62)) + mt->mti);
}

Referenced by mertwiInitByArray64(), mertwiRandomInt64(), and nloptMPSO().

◆ mertwiRandomBetween()

int mertwiRandomBetween	(	MERTWI *	mt,
		unsigned int	nr,
		double *	d,
		double	low,
		double	up,
		int	type )

extern

Fill the given double array with random numbers with uniform distribution between the specified limits.

Applies Mersenne Twister MT19937 pseudo-random number generator. With uniform distribution, the SD=(up-low)/sqrt(12), and CV=(up-low)/(sqrt(3)*(low+up)).

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64, mertwiRandomDouble1

Author: Vesa Oikonen

Returns: 0 when successful, otherwise <> 0.

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.
nr	Number of values in double array d[].
d	Pointer to pre-allocated double array, to be filled here.
low	Lower limit for random values.
up	Upper limit for random values.
type	Distribution: 0=even, 1=square-root transformation.

Definition at line 270 of file mertwi.c.

  {
  unsigned int i;
  double dif, v, stl, stu;
 
  if(nr<1) return(0);
  if(mt==NULL || d==NULL || type<0 || type>1) return(1);
 
  dif=up-low; if(dif<0.0) return(2);
  if(dif==0.0) {
    for(i=0; i<nr; i++) d[i]=low;
    return(0);
  }
 
  if(type==0) {
    for(i=0; i<nr; i++) d[i] = mertwiRandomDouble1(mt)*dif + low;
  } else if(type==1) {
    stl=copysign(sqrt(fabs(low)),low); if(!isnormal(stl)) stl=0.0;
    stu=copysign(sqrt(fabs(up)), up); if(!isnormal(stu)) stu=0.0;
    dif=stu-stl;
    for(i=0; i<nr; i++) {
      v=mertwiRandomDouble1(mt)*dif+stl; d[i]=copysign(v*v, v);
    }
  }
 
  return(0);
}

◆ mertwiRandomDouble1()

double mertwiRandomDouble1 ( MERTWI * mt )

extern

Generate a 64-bit double precision floating point pseudo-random number in the range of [0,1] with uniform distribution using Mersenne Twister MT19937.

With uniform distribution, the SD=(up-low)/sqrt(12), and CV=(up-low)/(sqrt(3)*(low+up)).

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64

Returns: Returns the random double value in the range [0,1].

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 216 of file mertwi.c.

  {
  return(mertwiRandomInt64(mt) >> 11) * (1.0/9007199254740991.0);
}

Referenced by mertwiRandomBetween(), mertwiRandomExponential(), mertwiRandomGaussian(), nloptMPSO(), and nloptRandomPoint().

◆ mertwiRandomDouble2()

double mertwiRandomDouble2 ( MERTWI * mt )

extern

Generate a 64-bit double precision floating point pseudo-random number in the range of [0,1) with uniform distribution using Mersenne Twister MT19937.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64

Returns: Returns the random double value in the range [0,1).

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 232 of file mertwi.c.

  {
  return(mertwiRandomInt64(mt) >> 11) * (1.0/9007199254740992.0);
}

◆ mertwiRandomDouble3()

double mertwiRandomDouble3 ( MERTWI * mt )

extern

Generate a 64-bit double precision floating point pseudo-random number in the range of (0,1) with uniform distribution using Mersenne Twister MT19937.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64

Returns: Returns the random double value in the range (0,1).

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 248 of file mertwi.c.

  {
  return((mertwiRandomInt64(mt) >> 12) + 0.5) * (1.0/4503599627370496.0);
}

◆ mertwiRandomExponential()

double mertwiRandomExponential	(	MERTWI *	mt,
		double	mean )

extern

Generate pseudo-random number with exponential distribution and specified mean.

Applies Mersenne Twister MT19937 pseudo-random number generator.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64, mertwiRandomDouble1

Returns: Returns the random non-negative double value with exponential distribution and given mean.

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.
mean	Mean of the exponential distribution.

Definition at line 322 of file mertwi.c.

  {
  double r;
  do {r=mertwiRandomDouble1(mt);} while(r==0.0);
  return(-mean*log(r)); 
}

Referenced by nloptMPSO().

◆ mertwiRandomGaussian()

double mertwiRandomGaussian ( MERTWI * mt )

extern

Generate a 64-bit double precision floating point pseudo-random number with normal (Gaussian) distribution with zero mean 0 and SD 1.

Applies Mersenne Twister MT19937 pseudo-random number generator, and the polar form of Box-Müller transform to produce numbers with Gaussian (normal) distribution which has a zero mean and standard deviation of one.

Box GEP, Muller ME. A note on the generation of random normal deviates. Annals of Mathematical Statistics, Volume 29, Issue 2, 1958, 610-611. Available from JSTOR https://www.jstor.org/

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64, mertwiRandomDouble1, mertwiRandomBetween

Returns: The random double value with normal distribution with zero mean and standard deviation 1.

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 354 of file mertwi.c.

  {
  static int ready=0;
  static double dev;
  double fac, rsq, a, b;
 
  /* If we don't have deviate already, then we'll have to make one */
  if(!ready) {
    do {
      a = 2.*mertwiRandomDouble1(mt) - 1.0; 
      b = 2.*mertwiRandomDouble1(mt) - 1.0; 
      rsq = a*a + b*b; 
    } while (rsq>1.0 || rsq==0.0);
   
    fac = sqrt(-2.0*log(rsq)/rsq);
    dev=a*fac; ready=1; 
    return(b*fac); 
  } else { /* dev is ready so return it */
    ready=0;
    return(dev);
  }
}

Referenced by nloptGaussianPoint().

◆ mertwiRandomInt63()

int64_t mertwiRandomInt63 ( MERTWI * mt )

extern

Generate a random number on [0, 2^63-1]-interval using Mersenne Twister MT19937.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64

Returns: Returns the random number with range from 0 to INT64_MAX.

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 197 of file mertwi.c.

  {
  return(int64_t)(mertwiRandomInt64(mt) >> 1);
}

Referenced by nloptMPSO().

◆ mertwiRandomInt64()

uint64_t mertwiRandomInt64 ( MERTWI * mt )

extern

Generate a random number on [0, 2^64-1]-interval using Mersenne Twister MT19937.

Precondition: Initialize MERTWI data structure by calling mertwiInit() before using this function. Preferably also initialize state vector by calling either mertwiInitWithSeed64() or mertwiInitByArray64; if you do not do that, mertwiInitWithSeed64() is called automatically with predetermined seed.

See also: mertwiInit, mertwiInitWithSeed64, mertwiInitByArray64

Returns: Returns the random number with range from 0 to UINT64_MAX.

Parameters

mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 153 of file mertwi.c.

  {
  unsigned int i;
  uint64_t x;
  static uint64_t mag01[2]={UINT64_C(0), TPCCLIB_MERTWI_A};
 
  if(mt->mti>=mt->n) { /* generate NN words at one time */
 
    /* if init_genrand64() has not been called, a default initial seed is used */
    if(mt->mti==mt->n+1) mertwiInitWithSeed64(mt, UINT64_C(5489)); 
 
    for(i=0; i<mt->n-mt->m; i++) {
      x=(mt->mt[i]&mt->um)|(mt->mt[i+1]&mt->lm);
      mt->mt[i] = mt->mt[i+mt->m] ^ (x>>1) ^ mag01[(int)(x&UINT64_C(1))];
    }
    for(; i<mt->n-1; i++) {
      x = (mt->mt[i]&mt->um)|(mt->mt[i+1]&mt->lm);
      mt->mt[i] = mt->mt[i+(mt->m-mt->n)] ^ (x>>1) ^ mag01[(int)(x&UINT64_C(1))];
    }
    x = (mt->mt[mt->n-1]&mt->um)|(mt->mt[0]&mt->lm);
    mt->mt[mt->n-1] = mt->mt[mt->m-1] ^ (x>>1) ^ mag01[(int)(x&UINT64_C(1))];
    mt->mti=0;
  }
  x = mt->mt[mt->mti];
 
  x ^= (x>>29) & UINT64_C(0x5555555555555555);
  x ^= (x<<17) & UINT64_C(0x71D67FFFEDA60000);
  x ^= (x<<37) & UINT64_C(0xFFF7EEE000000000);
  x ^= (x>>43);
  mt->mti++;
 
  return(x);
}

Referenced by mertwiRandomDouble1(), mertwiRandomDouble2(), mertwiRandomDouble3(), mertwiRandomInt63(), randomShuffle(), and randomShuffleUnsigned().

◆ mertwiSeed32()

uint32_t mertwiSeed32 ( void )

extern

Make uint32_t seed for pseudo-random number generators.

Uses microseconds from the computer clock and process ID to reduce the chance of getting the same seed for simultaneously executing program threads and instances.

Returns: Returns the seed.

See also: mertwiSeed64

Definition at line 48 of file mertwi.c.

{
  uint32_t li;
#if defined HAVE_TIMESPEC_GET
  struct timespec ts;
  timespec_get(&ts, TIME_UTC);
  li=((ts.tv_sec % 10000)*523 ^ ts.tv_nsec*10) ^ ((getpid() % 1000)*983);
#elif defined HAVE_CLOCK_GETTIME
  struct timespec ts;
  clock_gettime(CLOCK_REALTIME, &ts);
  li=((ts.tv_sec % 10000)*523 ^ ts.tv_nsec*10) ^ ((getpid() % 1000)*983);
#elif defined HAVE_GETTIMEOFDAY
  struct timeval tv;
  gettimeofday(&tv, 0);
  li=((tv.tv_sec % 10000)*523 ^ tv.tv_usec*13) ^ ((getpid() % 1000)*983);
#else
  li=(unsigned int)time(NULL)+(unsigned int)getpid();
#endif
  li+=(uint32_t)rand();
  return(li);
}

Referenced by mertwiSeed64().

◆ mertwiSeed64()

uint64_t mertwiSeed64 ( void )

extern

Make uint64_t seed for pseudo-random number generators.

Uses microseconds from the computer clock and process ID to reduce the chance of getting the same seed for simultaneously executing program threads and instances.

Returns: Returns the seed.

See also: mertwiSeed32

Definition at line 76 of file mertwi.c.

{
  uint32_t li;
  uint64_t lli;
  lli=li=mertwiSeed32();
  lli=li; lli<<=32; lli+=li;
  return(lli);
}

Referenced by nloptMPSO().

◆ randomPermutation()

void randomPermutation	(	int *	array,
		unsigned int	n,
		int	a,
		MERTWI *	mt )

extern

Random permutation for an array of signed integers.

Given (allocated) array of length n is filled with random permutation of numbers (signed integers) in the range [a:n+a-1]; that is, each number once in random order.

Precondition: Prepare the Mersenne Twister, for example using code MERTWI mt; mertwiInit(&mt); mertwiInitWithSeed64(&mt, mertwiSeed64());

See also: randomPermutationUnsigned, randomShuffle

Parameters

array	Allocated array to be filled.
n	Length of array.
a	Smallest integer value.
mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 80 of file shuffle.c.

  {
  if(n<1 || array==NULL) return;
  int i;
  for(i=0; i<(int)n; i++) array[i]=i+a;
  randomShuffle(array, n, mt);
}

◆ randomPermutationUnsigned()

void randomPermutationUnsigned	(	unsigned int *	array,
		unsigned int	n,
		unsigned int	a,
		MERTWI *	mt )

extern

Random permutation for an array of unsigned integers.

Given (allocated) array of length n is filled with random permutation of numbers (unsigned integers) in the range [a:n+a-1]; that is, each number once in random order.

Precondition: Prepare the Mersenne Twister, for example using code MERTWI mt; mertwiInit(&mt); mertwiInitWithSeed64(&mt, mertwiSeed64());

See also: randomPermutation, randomShuffleUnsigned

Parameters

array	Allocated array to be filled.
n	Length of array.
a	Smallest integer value.
mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 106 of file shuffle.c.

  {
  if(n<1 || array==NULL) return;
  unsigned int i;
  for(i=0; i<n; i++) array[i]=i+a;
  randomShuffleUnsigned(array, n, mt);
}

◆ randomShuffle()

void randomShuffle	(	int *	array,
		unsigned int	n,
		MERTWI *	mt )

extern

Random shuffle.

Arrange the n elements of integer array in random order. Only effective if N is much smaller than UINT64_MAX.

Precondition: Prepare the Mersenne Twister, for example using code MERTWI mt; mertwiInit(&mt); mertwiInitWithSeed64(&mt, mertwiSeed64());

See also: randomShuffleUnsigned

Parameters

array	Integer array to be shuffled.
n	Array length.
mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 27 of file shuffle.c.

  {
  if(n<=1 || array==NULL) return;
  unsigned int i, j;
  int tmp;
  for(i=0; i<n-1; i++) {
    j=i+mertwiRandomInt64(mt)/(UINT64_MAX/(n-i)+1);
    tmp=array[j]; array[j]=array[i]; array[i]=tmp;
  }
}

Referenced by randomPermutation().

◆ randomShuffleUnsigned()

void randomShuffleUnsigned	(	unsigned int *	array,
		unsigned int	n,
		MERTWI *	mt )

extern

Random shuffle.

Arrange the n elements of unsigned integer array in random order. Only effective if N is much smaller than UINT64_MAX.

Precondition: Prepare the Mersenne Twister, for example using code MERTWI mt; mertwiInit(&mt); mertwiInitWithSeed64(&mt, mertwiSeed64());

See also: randomShuffle, randomPermutation

Parameters

array	Integer array to be shuffled.
n	Array length.
mt	Pointer to data structure for Mersenne Twister MT19937 pseudo-random number generator.

Definition at line 54 of file shuffle.c.

  {
  if(n<=1 || array==NULL) return;
  unsigned int i, j, tmp;
  for(i=0; i<n-1; i++) {
    j=i+mertwiRandomInt64(mt)/(UINT64_MAX/(n-i)+1);
    tmp=array[j]; array[j]=array[i]; array[i]=tmp;
  }
}

Referenced by randomPermutationUnsigned().

Data Structures

Macros

Functions

Detailed Description

Macro Definition Documentation

◆ TPCCLIB_MERTWI_A

◆ TPCCLIB_MERTWI_NN

Function Documentation

◆ drand()

◆ drandExponential()

◆ drandGaussian()

◆ drandRange()

◆ drandSeed()

◆ haltonElement()

◆ haltonPrime()

◆ mertwiInit()

◆ mertwiInitByArray64()

◆ mertwiInitWithSeed64()

◆ mertwiRandomBetween()

◆ mertwiRandomDouble1()

◆ mertwiRandomDouble2()

◆ mertwiRandomDouble3()

◆ mertwiRandomExponential()

◆ mertwiRandomGaussian()

◆ mertwiRandomInt63()

◆ mertwiRandomInt64()

◆ mertwiSeed32()

◆ mertwiSeed64()

◆ randomPermutation()

◆ randomPermutationUnsigned()

◆ randomShuffle()

◆ randomShuffleUnsigned()