SF_RANDOM

Description

SciFortran module for RNG

Quick access

Subroutines and functions

interface sf_random/mersenne_init(seed)

Generic interface to initialize the Mersenne Twister random number generator with the integer seed (init_genrand), same as mt_init. The 624 words of the state are filled from the seed with the linear recurrence of Knuth (multiplier 1812433253), assuming that the integer overflow does not stop the program. If the generator is not initialized, the first call to mersenne uses the default seed 4357.

Parameters:: seed [integer, in] – seed of the generator

interface sf_random/mt_init(seed)

Generic interface to initialize the Mersenne Twister random number generator with the integer seed (init_genrand), same as mersenne_init. See mersenne_init for the details.

Parameters:: seed [integer, in] – seed of the generator

interface sf_random/mersenne(): Generic interface to the Mersenne Twister random number generator MT19937 (grnd). Every call returns a real(8) random number uniformly distributed between 0 and 1, obtained dividing a 32 bit random integer by \(2^{32}-1\). The 624 words of the state are regenerated every 624 calls. The generator is initialized with mersenne_init, and its state can be saved and read back with mt_save and mt_get. All the mt_* generators use it as the source of uniform numbers.

interface sf_random/mt_random(a)

Generic interface to fill an array of rank 1 to 7 with uniform random numbers between 0 and 1, generated with mersenne. The array is real(8) (d_grnd_1 to d_grnd_7) or complex(8) (c_grnd_1 to c_grnd_7), in which case the real and the imaginary parts are two independent numbers.

Parameters:: a (various shapes) [real/complex] – real(8) array of rank 1, filled with random numbers

interface sf_random/mt_uniform(l, h, a, b)

Generic interface for uniform random numbers generated with mersenne: an integer between l and h, both included (igrnd), or a real(8) number between a and b (dgrnd_uniform).

Parameters:

l [integer, in] – lower bound, included
h [integer, in] – upper bound, included
a [real] – lower bound
b [real] – upper bound

Result:

c [real] – random number in [a,b]

interface sf_random/mt_normal(mean, stdev)

Generic interface for normally distributed random numbers: a standard normal deviate, with zero mean and unit standard deviation (gaussrnd), or a normal deviate with the given mean and standard deviation (normalrnd). The former uses the polar form of the Box-Muller method (Numerical Recipes), generating two deviates at a time and returning the second one at the following call; the latter uses the basic Box-Muller method, one deviate per call.

Parameters:

mean [real] – mean
stdev [real] – standard deviation, positive

Result:

c [real] – normal deviate

interface sf_random/mt_exponential(mean)

Generic interface for random numbers from an exponential distribution with mean mean, positive, obtained by inversion (exponentialrnd). If the mean is not positive a message is printed and the result is undefined.

Parameters:: mean [real] – mean of the distribution, positive
Result:: c [real] – exponential deviate

interface sf_random/mt_gamma(shape, scale)

Generic interface for random numbers from a gamma distribution with parameters shape and scale, both positive, i.e. with mean \(shape\cdot scale\) (gammarnd). It uses the method of Marsaglia and Tsang (ACM TOMS 26, 363, 2000), and for \(shape<1\) the relation with a gamma deviate of shape \(shape+1\). If a parameter is not positive a message is printed and the result is undefined.

Parameters:

shape [real] – shape parameter, positive
scale [real] – scale parameter, positive

Result:

ans [real] – gamma deviate, mean shape*scale

interface sf_random/mt_chi_square(dof)

Generic interface for random numbers from a chi-square distribution with dof degrees of freedom (chi_squarernd), obtained from mt_gamma.

Parameters:: dof [real] – degrees of freedom
Result:: ans [real] – chi-square deviate

interface sf_random/mt_inverse_gamma(shape, scale)

Generic interface for random numbers from an inverse gamma distribution with parameters shape and scale (inverse_gammarnd), obtained as the inverse of a gamma deviate, mt_gamma, with shape shape and scale \(1/scale\).

Parameters:

shape [real] – shape parameter, positive
scale [real] – scale parameter, positive

Result:

ans [real] – inverse gamma deviate If X is gamma(shape, scale) then 1/Y is inverse gamma(shape, 1/scale)

interface sf_random/mt_weibull(shape, scale)

Generic interface for random numbers from a Weibull distribution with parameters shape and scale, both positive (weibullrnd), obtained by inversion, \(scale\,(-\ln u)^{1/shape}\). If a parameter is not positive a message is printed and the result is undefined.

Parameters:

shape [real] – shape parameter, positive
scale [real] – scale parameter, positive

Result:

ans [real] – Weibull deviate

interface sf_random/mt_cauchy(median, scale)

Generic interface for random numbers from a Cauchy distribution with median median and scale scale, positive (cauchyrnd), obtained by inversion. If the scale is not positive a message is printed and the result is undefined.

Parameters:

median [real] – median of the distribution
scale [real] – scale parameter, positive

Result:

ans [real] – Cauchy deviate

interface sf_random/mt_student_t(dof)

Generic interface for random numbers from a Student t distribution with dof degrees of freedom, positive (student_trnd), obtained as the ratio of a standard normal deviate, mt_normal, and the square root of a chi-square deviate, mt_chi_square, divided by the degrees of freedom. If dof is not positive a message is printed and the result is undefined.

Parameters:: dof [real] – degrees of freedom, positive
Result:: ans [real] – Student t deviate

interface sf_random/mt_laplace(mean, scale)

Generic interface for random numbers from a Laplace (double exponential) distribution with mean mean and scale scale, positive (laplacernd), obtained by inversion. If the scale is not positive a message is printed and the result is undefined.

Parameters:

mean [real] – mean of the distribution
scale [real] – scale parameter, positive

Result:

ans [real] – Laplace deviate

interface sf_random/mt_log_normal(mu, sigma)

Generic interface for random numbers from a log-normal distribution (log_normalrnd), i.e. the exponential of a normal deviate, mt_normal, with mean mu and standard deviation sigma. Note that mu and sigma are the mean and standard deviation of the logarithm of the variable, not of the variable.

Parameters:

mu [real] – mean of the logarithm of the variable
sigma [real] – standard deviation of the logarithm

Result:

ans [real] – log-normal deviate

interface sf_random/mt_beta(a, b)

Generic interface for random numbers from a beta distribution with positive shape parameters a and b (betarnd), obtained as \(u/(u+v)\) with \(u\) and \(v\) gamma deviates, mt_gamma, of shapes a and b. If a parameter is not positive a message is printed and the result is undefined.

Parameters:

a [real] – first shape parameter, positive
b [real] – second shape parameter, positive

Result:

ans [real] – beta deviate, in [0,1]

interface sf_random/mt_save(fname, forma, unum)

Generic interface to save the state of the Mersenne Twister generator, i.e. the index and the 624 words of the state, to a file (mtsavef) or to an open unit (mtsaveu). With forma equal to 'u' or 'U' the state is written unformatted, otherwise formatted. NOTE: a file is not overwritten, the state is APPENDED to it, and the unit 10 is used. The state can be read back with mt_get.

Parameters:

fname [character(len=*), in] – name of the file, the state is appended
forma [character, in] – ‘u’ or ‘U’: unformatted, otherwise formatted
unum [integer, in] – open unit where the state is written

interface sf_random/mt_get(fname, forma, unum)

Generic interface to read the state of the Mersenne Twister generator, saved with mt_save, from a file (mtgetf) or from an open unit (mtgetu). With forma equal to 'u' or 'U' the state is read unformatted, otherwise formatted. The file must exist and the first saved state is read; the unit 10 is used for the file.

Parameters:

fname [character(len=*), in] – name of the file, must exist
forma [character, in] – ‘u’ or ‘U’: unformatted, otherwise formatted
unum [integer, in] – open unit from which the state is read

subroutine sf_random/random_number_init([shift])

This subroutine initializes the random number generator of the compiler, the one used by the intrinsic random_number, with a seed array built from the system clock, seed(i) = clock + 37*(i-1). The optional argument shift is added to every element of the seed, and it can be used to make different the seeds of processes started at the same time.

Options:: shift [integer] – optional: integer added to every element of the seed

function sf_random/random_number_seed([info, file])

This function returns a positive integer to be used as the seed of a random number generator, read from the device /dev/urandom (default) or /dev/random, with stream i/o. The most significant bit is cleared to get a positive number. If the device cannot be opened or read, the seed is computed from the current date and time. The simplest way to call it is random_number_seed(). Note that /dev/random may be slow, since it can wait for enough entropy.

Options:

info [integer, in] – optional: if /= 0 print the origin of the seed to the standard output
file [integer, in] – optional: 0 use /dev/urandom (default), /= 0 use /dev/random

Result:

random_number_seed [integer]

function sf_random/random_normal()

This function returns a single precision normally distributed random number with zero mean and unit variance, generated with the ratio of uniforms method of Kinderman and Monahan, augmented with quadratic bounding curves (ACM Algorithm 712). The uniform numbers are generated with the intrinsic random_number.

Result:: fn_val [real] – normal deviate, mean 0 and variance 1

function sf_random/random_gamma(s, first)

This function returns a single precision random number from a gamma distribution with shape parameter s and unit scale, calling random_gamma1 if \(s>1\), random_exponential if \(s=1\) and random_gamma2 if \(s<1\). The program stops if s is not positive.

Parameters:

s [real, in] – shape parameter, positive
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – gamma deviate

function sf_random/random_gamma1(s, first)

This function returns a single precision random number from a gamma distribution with shape parameter s not smaller than 1 and unit scale, with the method of Marsaglia and Tsang (ACM TOMS 26, 363, 2000). It is called by random_gamma.

Parameters:

s [real, in] – shape parameter, >= 1
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – gamma deviate

function sf_random/random_gamma2(s, first)

This function returns a single precision random number from a gamma distribution with shape parameter s smaller than 1 and unit scale, with a switching method (Dagpunar, 1988). It is called by random_gamma. The program stops if s is not in the interval (0,1).

Parameters:

s [real, in] – shape parameter, in (0,1)
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – gamma deviate

function sf_random/random_chisq(ndf, first)

This function returns a single precision random number from a chi-square distribution with ndf degrees of freedom, computed as twice a gamma deviate of shape \(ndf/2\), random_gamma.

Generates a random variate from the chi-squared distribution with ndf degrees of freedom

Parameters:

ndf [integer, in] – degrees of freedom
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – chi-square deviate

function sf_random/random_exponential()

This function returns a single precision random number from an exponential distribution with unit mean, density \(e^{-x}\), obtained by inversion.

Result:: fn_val [real] – exponential deviate, unit mean

function sf_random/random_weibull(a)

This function returns a single precision random number from a Weibull distribution with shape parameter a and unit scale, density \(a x^{a-1} e^{-x^a}\). For speed, there is no check that a is not zero or very small.

Parameters:: a [real, in] – shape parameter, positive
Result:: fn_val [real] – Weibull deviate

function sf_random/random_beta(aa, bb, first)

This function returns a single precision random number in [0,1] from a beta distribution with density proportional to \(x^{aa-1}(1-x)^{bb-1}\), with Cheng’s log-logistic method (Dagpunar, 1988). The program stops if a shape parameter is not positive.

Parameters:

aa [real, in] – first shape parameter, positive
bb [real, in] – second shape parameter, positive
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – beta deviate, in [0,1]

function sf_random/random_inv_gauss(h, b, first)

This function returns a single precision random number from a reparameterised generalised inverse Gaussian distribution, with density proportional to \(x^{h-1}\exp[-b(x+1/x)/2]\), with a ratio method (Dagpunar, 1988). The program stops if \(h<0\), if \(b\le 0\), or if the ratio \(h:b\) is too small.

Parameters:

h [real, in] – parameter of the distribution, >= 0
b [real, in] – parameter of the distribution, > 0
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – deviate, in [0,inf)

function sf_random/random_poisson(mu, first)

This function returns an integer random number from a Poisson distribution with mean mu, with the method of Ahrens and Dieter (ACM TOMS 8, 163, 1982), based on modified normal distributions.

Parameters:

mu [real, in] – mean of the distribution
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

ival [integer] – Poisson deviate

function sf_random/random_binomial1(n, p, first)

This function returns an integer random number from a binomial distribution of n trials with success probability p, with the modal method of Kemp (1986). It is suitable when many variates with the SAME parameters are required; if the parameters change often use random_binomial2.

Parameters:

n [integer, in] – number of trials, >= 1
p [real, in] – probability of success of each trial, in [0,1]
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

ival [integer] – binomial deviate, in [0,n]

function sf_random/bin_prob(n, p, r)

This function returns the binomial probability of \(r\) successes in \(n\) trials with success probability \(p\), \(\binom{n}{r}p^r(1-p)^{n-r}\), computed with the logarithm of the gamma function, lngamma.

Calculate a binomial probability

Parameters:

n [integer, in] – number of trials
p [real, in] – probability of success of each trial
r [integer, in] – number of successes

Result:

fn_val [real] – probability of r successes

function sf_random/lngamma(x)

This function returns the natural logarithm of the gamma function \(\ln\Gamma(x)\), accurate to about \(10^{-14}\). It is not defined for \(x\) equal to zero or to a negative integer. Programmer: Alan Miller.

Logarithm to base e of the gamma function.

Accurate to about 1.e-14. Programmer: Alan Miller Latest revision of Fortran 77 version - 28 February 1988

Parameters:: x [real, in] – argument, not zero or a negative integer
Result:: fn_val [real] – natural logarithm of the gamma function at x

function sf_random/random_binomial2(n, pp, first)

This function returns an integer random number from a binomial distribution of n trials with success probability pp, with the BTPE algorithm of Kachitvichyanukul and Schmeiser (Communications of the ACM 31, 216, 1988). It is suitable when the parameters change often; if many variates with the SAME parameters are required use random_binomial1.

Parameters:

n [integer, in] – number of trials
pp [real, in] – probability of success of each trial
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

ival [integer] – binomial deviate, in [0,n]

function sf_random/random_neg_binomial(sk, p)

This function returns an integer random number from a negative binomial distribution, with unstored inversion and/or the reproductive property (Dagpunar, 1988). The program stops if \(sk\le 0\), or if p is not in the open interval (0,1).

Parameters:

sk [real, in] – number of failures required (power parameter), > 0
p [real, in] – probability of success, in (0,1)

Result:

ival [integer] – negative binomial deviate

function sf_random/random_von_mises(k, first)

This function returns a single precision random number in \((-\pi,\pi]\) from a von Mises distribution with concentration parameter k, with algorithm VMD of Dagpunar (1990). The parameter \(k\) must be non-negative and not larger than 9.5: otherwise an error message is printed and the result is undefined.

Parameters:

k [real, in] – concentration parameter, 0 <= k <= 9.5
first [logical, in] – set .true. on the first call with new parameter values, .false. if unchanged

Result:

fn_val [real] – von Mises deviate, in (-pi,pi]

function sf_random/random_cauchy()

This function returns a single precision random number from the standard Cauchy distribution, with median 0 and scale 1, as the ratio of two uniform deviates inside the unit circle.

Generate a random deviate from the standard Cauchy distribution

Result:: fn_val [real] – Cauchy deviate

function sf_random/nrand([dseed_])

This function returns a uniform random number in the open interval (0,1), generated with the algorithm ran2 of Numerical Recipes, i.e. the combination of two congruential generators of L’Ecuyer with a Bays-Durham shuffle. The optional seed dseed_ initializes the generator if it is not positive, with the seed -dseed_ (or 1). Only the second generator and the shuffle table are saved between calls, the state of the first one is not.

Options:: dseed [integer] – optional: a value <= 0 initializes the generator with seed -dseed_
Result:: nrand [real]

subroutine sf_random/random_order(order, n)

This subroutine generates a random permutation of the integers 1,…,:f:var:n, with the Fisher-Yates shuffle: starting from the end, the last element is swapped with one randomly chosen among the preceding ones. The random numbers are generated with the intrinsic random_number.

generate a random ordering of the integers 1 … n.

Parameters:

order (n) [integer, out] – random permutation of 1,…,n, size n local variables
n [integer, in] – number of integers