17 Pseudo-Random Numbers

 17.1 Setting the seed (replacements for G05CBF and G05CCF)
 17.2 Data type of the random numbers
 17.3 Replacements for G05CAF and G05DAF
 17.4 Replacement for G05DBF
 17.5 Replacements for G05DDF, G05DRF, and G05FFF

The routines for creating pseudo-random numbers in this library all have a period of 226 and 6–7 digits accuracy. They are based upon code by Ahrens, Dieter, & Grube. They use a multiplicative congruential generator which is certainly not the state of the art and may not be suitable for critical or sophisticated use.

The routines are:

17.1 Setting the seed (replacements for G05CBF and G05CCF)

Before any random numbers can be selected, a seed must be set using PDA_RNSED. The integer seed should satisfy the relationship

seed = 4 ∗k + 1

where k is a non-negative integer. A fixed seed gives rise to a reproducible sequence of pseudo-random numbers.

For a non-repeatable sequence, there is no equivalent to NAG routine G05CCF because the system clock used to create the seed is not accessible portably in Fortran, and PDA is independent of other libraries. However, the following code has the desired effect.

        INTEGER SEED, STATUS, TICKS, PID
        INCLUDE ’PRM_PAR’
  
        CALL PSX_TIME( TICKS, STATUS )
        CALL PSX_GETPID( PID, STATUS )
        SEED = TICKS + PID
        SEED = MOD( SEED, VAL__MAXI / 4 ) * 4 + 1
        SEED = MOD( SEED, 2**28 )
  
        CALL PDA_RNSED( SEED )

PSX_TIME returns the time in units of clock ticks since some arbitrary time. See SUN/121 for more details and linking instructions. The above code also permits storage of the chosen seed.

17.2 Data type of the random numbers

There is a major difference between the PDA random-number routines and those provided in the standard NAG library: in general the former are single-precision functions, whereas the latter are double precision. However, PDA_RNPOI and the corresponding G05DRF are both integer functions.

17.3 Replacements for G05CAF and G05DAF

Like G05CAF, PDA_RAND has a dummy argument demanded by the Fortran standard. It is convenient to set it to zero. Here is an example where two random numbers are drawn from a uniform distribution between 0 and 1. In this example a fixed seed is used, but you could use the computer’s clock to create a random seed (see Section 17.1).

        INTEGER SEED
        EXTERNAL PDA_RAND
        REAL PDA_RAND, VALUES( 2 )
  
  *  Use a fixed seed of 1.
        SEED = 1
        CALL PDA_RNSED( SEED )
  
  *  Obtain two random numbers from a uniform distribution between 0
  *  and 1.
        VALUE( 1 ) = PDA_RAND( 0.0 )
        VALUE( 2 ) = PDA_RAND( 0.0 )

The EXTERNAL statement is recommended, although in many cases it will be unnecessary. To obtain in the range [a,b] as provided by G05DAF, merely apply the following relationship.

random value = (b−a) ∗ PDA_RAND(0.0) + a

17.4 Replacement for G05DBF

PDA_RNEXP is only a partial replacement for G05DBF in that it computes pseudo-random numbers from e−x, whereas G05DBF uses the function 1 ae−x/a. Thus its argument is also a dummy mandated by the Fortran standard.

17.5 Replacements for G05DDF, G05DRF, and G05FFF

The following code shows the remaining three routines in action.

        INTEGER SEED
        EXTERNAL PDA_RNGAM, PDA_RNNOR, PDA_RNPOI
        REAL PDA_RNGAM, PDA_RNNOR, PDA_RNPOI, VALUES( 3 )
  
  *  Use a fixed seed of 1001.
        SEED = 1001
        CALL PDA_RNSED( SEED )
  
  *  Obtain a random number from a Normal distribution of mean 4.2 and
  *  standard deviation 0.15
        VALUE( 1 ) = PDA_RNNOR( 4.2, 0.15 )
  
  *  Obtain a random number from a Poisson distribution of mean 3.4.
        VALUE( 2 ) = PDA_RNPOI( 3.4 )
  
  *  Obtain a random number from a Gamma-function distribution of mean
  *  1.2.
        VALUE( 3 ) = PDA_RNGAM( 1.2 )

Apart from the change of data type, calls to G05DDF can be replaced with PDA_RNNOR using the same arguments. PDA_RNPOI is in effect a renamed G05DRF.

PDA_RNGAM only has one argument—the mean—of the Gamma function, whereas G05FFF has a second scaling parameter similar in role to the a argument of G05DBF. G05FFF also generates a vector of pseudo-random numbers.