subroutine aptexps (nopt, sigma, np, xran, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTEXPS
c
c call aptexps (nopt, sigma, np, xran, nerr)
c
c Version: aptexps Updated 1994 June 27 16:40.
c aptexps Originated 1990 August 16 15:00.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To randomly sample each of np values xran(n) from an
c exponential probability distribution function (PDF)
c with decay constant sigma(k), where k = 1 (nopt = 0),
c or k = n (nopt = 1):
c
c PDF = sigma(k) * exp (-sigma(k) * x)
c
c where x ranges from zero to plus infinity.
c Flag nerr indicates any input error.
c
c Input: nopt, sigma, np.
c
c Output: xran, nerr.
c
c Glossary:
c
c nerr Output Indicates an input error, if not zero.
c 1 if nopt is not 0 or 1.
c 2 if np is not positive.
c
c nopt Input Indicates the sampling option:
c 0: All np xran values will be randomly sampled from
c one exponential probability distribution function
c with decay constant sigma (a scalar).
c 1: Each value xran(n) will be randomly sampled from
c an exponential probability distribution function
c with decay constant sigma(n) (an array).
c
c np Input Size of arrays xran, and if nopt = 1, sigma.
c
c sigma Input Decay constant. Size 1 (nopt = 0) or np (nopt = 1).
c If x is a collision distance, sigma is a cross
c section (per unit distance) or 1 / mean-free-path.
c If x is a time, sigma is log (2) / half-life.
c The absolute value of sigma is used.
c
c xran Output Randomly sampled values of x. For xran(n), from the
c exponential probability distribution function with
c decay constant sigma, where k = 1 (nopt = 0),
c or k = n (nopt = 1). Size np.
c The expected value of xran is 1.0 / sigma, with
c a standard deviation of 1.0 / sigma.
c
c History:
c
c 1994 June 27 16:40. Corrected bug: output xran values were a factor
c of sigma**2 times the correct values.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Dimensioned arguments.
c---- Decay constant.
dimension sigma (1)
c---- Randomly sampled value of x.
dimension xran (1)
c.... Local variables.
c---- Index in local array.
common /laptexps/ n
c---- First index of subset of data.
common /laptexps/ n1
c---- Last index of subset of data.
common /laptexps/ n2
c---- Index in external array.
common /laptexps/ nn
c---- Size of current subset of data.
common /laptexps/ ns
c---- Random numbers.
common /laptexps/ ranfp (64)
cbugc***DEBUG begins.
cbugc---- Standard deviation from mean.
cbug common /laptexps/ xdev
cbugc---- Mean value.
cbug common /laptexps/ xmean
cbug 9901 format (/ 'aptexps sampling an exponential distribution.',
cbug & ' nopt=',i2,' np=',i5)
cbug 9902 format (' sigma=',1pe22.14)
cbug 9903 format (i5,' sigma=',1pe22.14)
cbug write ( 3, 9901) nopt, np
cbug if (nopt .eq. 0) then
cbug write ( 3, 9902) sigma(1)
cbug elseif (nopt .eq. 1) then
cbug write (3, 9903) (n, sigma(n), n = 1, np)
cbug endif
cbugc***DEBUG ends.
c.... Test for input errors.
nerr = 0
if ((nopt .ne. 0) .and. (nopt .ne. 1)) then
nerr = 1
go to 210
endif
if (np .le. 0) then
nerr = 2
go to 210
endif
c.... Do data sets in groups of 64.
n1 = 1
n2 = min (np, 64)
110 ns = n2 - n1 + 1
c.... Get needed random numbers.
c---- Loop over subset of data.
do 120 n = 1, ns
ranfp(n) = ranf( )
c---- End of loop over subset of data.
120 continue
c.... Randomly sample values of x.
c---- One PDF for all xran.
if (nopt .eq. 0) then
c---- Loop over subset of data.
do 130 n = 1, ns
nn = n + n1 - 1
xran(nn) = -(1.0 / abs (sigma(1))) * alog (ranfp(n))
c---- End of loop over subset of data.
130 continue
c---- One PDF for each xran.
else
c---- Loop over subset of data.
do 140 n = 1, ns
nn = n + n1 - 1
xran(nn) = -(1.0 / abs (sigma(nn))) * alog (ranfp(n))
c---- End of loop over subset of data.
140 continue
c---- Tested nopt.
endif
c.... See if all sets are done.
c---- Do another data set.
if (n2 .lt. np) then
n1 = n2 + 1
n2 = min (np, n1 + 63)
go to 110
endif
cbugc***DEBUG begins.
cbug 9904 format (/ 'aptexps results: xran')
cbug write ( 3, 9904)
cbug call aptmean (xran, np, 1.e-11, xmean, xdev, nerr)
cbugc***DEBUG ends.
210 return
c.... End of subroutine aptexps. (+1 line.)
end
UCRL-WEB-209832