Generates pseudorandom numbers from a general discrete distribution using an alias method.
IOPT — Indicator of whether the alias vectors are to be initialized. (Input)
IOPT Action
0 The alias vectors are to be initialized using the probabilities in PROBS. IOPT is set to 0 on the first call to RNGDA.
1 The alias vectors IWK and WK are used but PROBS is not used.
IMIN — Smallest
value the random deviate can assume. (Input)
This is the value
corresponding to the probability in PROBS(1).
PROBS — Vector of
length NMASS
containing probabilities associated with the individual mass points.
(Input)
The elements of PROBS must be
nonnegative and must sum to 1.0.
IWK — Index vector of length NMASS. (Input, if IOPT = 1; output, if IOPT = 0) IWK is a work vector.
WK — Index vector of length NMASS. (Input, if IOPT = 1; output, if IOPT = 0) WK is a work vector.
IR — Vector of length NR containing the random discrete deviates. (Output)
NR — Number of
random numbers to generate. (Input)
Default: NR = size (IR,1).
NMASS — Number of
mass points in the discrete distribution. (Input)
Default: NMASS = size (PROBS,1).
Generic: CALL RNGDA (IOPT, IMIN, PROBS, IWK, WK, IR [,…])
Specific: The specific interface names are S_RNGDA and D_RNGDA.
Single: CALL RNGDA (NR, IOPT, IMIN, NMASS, PROBS, IWK, WK, IR)
Double: The double precision name is DRNGDA.
Routine RNGDA generates pseudorandom numbers from a discrete distribution with probability function given in the vector PROBS; that is
Pr(X = i) = pj
for
The algorithm is the alias method, due to Walker (1974), with modifications suggested by Kronmal and Peterson (1979). The method involves a setup phase, in which the vectors IWK and WK are filled. After the vectors are filled, the generation phase is very fast.
1.
In the interest of efficiency, this routine does only limited error checking
when
IOPT =
1.
2. The routine RNSET can be used to initialize the seed of the random number generator. The routine RNOPT can be used to select the form of the generator.
In this example, RNGDA is used to generate five pseudorandom variates from the discrete distribution:
Pr(X = 1) = .05
Pr(X = 2) = .45
Pr(X = 3) = .31
Pr(X = 4) = .04
Pr(X = 5) = .15
When RNGDA is called the first time, IOPT is input as 0. This causes the work arrays to be initialized. In the next call, IOPT is 1, so the setup phase is bypassed.
USE RNGDA_INT
USE UMACH_INT
USE RNSET_INT
IMPLICIT NONE
INTEGER NMASS, NR
PARAMETER (NMASS=5, NR=5)
!
INTEGER IMIN, IOPT, IR(NR), ISEED, IWK(NMASS), NOUT
REAL PROBS(NMASS), WK(NMASS)
!
CALL UMACH (2, NOUT)
IMIN = 1
PROBS(1) = 0.05
PROBS(2) = 0.45
PROBS(3) = 0.31
PROBS(4) = 0.04
PROBS(5) = 0.15
IOPT = 0
ISEED = 123457
CALL RNSET (ISEED)
CALL RNGDA (IOPT, IMIN, PROBS, IWK, WK, IR)
WRITE (NOUT,99998) IR
99998 FORMAT (' Random deviates: ', 5I4)
IOPT = 1
CALL RNGDA (IOPT, IMIN, PROBS, IWK, WK, IR)
WRITE (NOUT,99999) IR
99999 FORMAT (' ', 5I4)
END
Random deviates: 3 2 2 3 5
1 3 4 5 3
In this example, RNGDA is used to generate five pseudorandom binomial variates with parameters 20 and 0.5.
USE UMACH_INT
USE RNSET_INT
USE RNGDA_INT
USE BINPR_INT
IMPLICIT NONE
INTEGER NMASS, NR
PARAMETER (NMASS=21, NR=5)
!
INTEGER IMIN, IOPT, IR(NR), ISEED, IWK(NMASS), K, N, NOUT
REAL P, PROBS(NMASS), WK(NMASS)
!
CALL UMACH (2, NOUT)
N = 20
P = 0.5
IMIN = 0
DO 10 K=1, NMASS
PROBS(K) = BINPR(K-1,N,P)
10 CONTINUE
IOPT = 0
ISEED = 123457
CALL RNSET (ISEED)
CALL RNGDA (IOPT, IMIN, PROBS, IWK, WK, IR)
WRITE (NOUT,99999) IR
99999 FORMAT (' Binomial (20, .5) deviates: ', 5I4)
END
Binomial (20, .5) deviates: 12 10 16 12 11
PHONE: 713.784.3131 FAX:713.781.9260 |