subroutine aptrint (np, imin, imax, frstr, iran, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTRINT
c
c call aptrint (np, imin, imax, frstr, iran, nerr)
c
c Version: aptrint Updated 1990 December 11 15:20.
c aptrint Originated 1990 December 11 15:20.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To randomly sample np integers from imin to imax, where each
c integer has an equal probability.
c Any fraction frstr of sampling may be striated, from 0.0 to 1.0.
c The results are stored in iran.
c Flag nerr indicates any input error.
c
c WARNING: STRIATED SAMPLING INTRODUCES CORRELATIONS BETWEEN SAMPLE INDICES
c AND INTEGER VALUES. USE ONLY WHEN SUCH CORRELATIONS ACCEPTABLE.
c HOWEVER, THE RESULTS MAY BE RANDOMLY REORDERED WITH APTRAND.
c
c Input: np, imin, imax, frstr.
c
c Output: iran, nerr.
c
c Glossary:
c
c frstr Input The fraction (from 0.0 to 1.0) of samples to striate
c over the integers. Remaining samples will be
c randomly sampled over the integers.
c
c imin Input The minimum integer to be sampled.
c
c imax Input The maximum integer to be sampled.
c
c iran Output The randomly sampled integers. Size np.
c
c nerr Output Indicates an input error, if not 0.
c 1 if np is not positive.
c 2 if frstr is less than 0.0 or more than 1.0.
c
c np Input Size of array iran. Number of samples.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Dimensioned arguments.
c---- Randomly sampled integer.
dimension iran (1)
c.... Local variables.
c---- Floating point result.
common /laptrint/ fran
c---- Index in internal array.
common /laptrint/ n
c---- First index of subset of data.
common /laptrint/ n1
c---- Last index of subset of data.
common /laptrint/ n2
c---- Number of integers in set.
common /laptrint/ nint
c---- Index in external array.
common /laptrint/ nn
c---- Size of current subset of data.
common /laptrint/ ns
c---- Number of striated samples.
common /laptrint/ nstrd
c---- Number of unstriated samples.
common /laptrint/ nunst
c---- Random numbers in range 0.0 to 1.0.
common /laptrint/ ranfp (64)
c---- Striated random number, 0.0 to 1.0.
common /laptrint/ ranst
cbugc***DEBUG begins.
cbugc---- Fractional sampling of integer.
cbug common /laptrint/ fb
cbugc---- Integer value.
cbug common /laptrint/ i
cbug 9901 format (/ 'aptrint sampling from integers.' /
cbug & ' np=',i5,' imin=',i6,' imax=',i6,' frstr=',1pe22.14)
cbug write ( 3, 9901) np, imin, imax, frstr
cbugc***DEBUG ends.
c.... initialize.
nerr = 0
if (imin .le. imax) then
iminn = imin
imaxx = imax
else
iminn = imax
imaxx = imin
endif
nint = imaxx - iminn + 1
c.... Test for input errors.
if (np .le. 0) then
nerr = 1
go to 210
endif
if ((frstr .lt. 0.0) .or. (frstr .gt. 1.0)) then
nerr = 2
go to 210
endif
c.... Find the number of striated and unstriated samples.
nstrd = frstr * np + 0.5
nunst = np - nstrd
c.... See if any unstriated samples are needed.
c---- Need unstriated samples.
if (nunst .gt. 0) then
c.... Set up the indices of the first subset of samples.
n1 = 1
n2 = min (nunst, 64)
c.... Loop over the subset of samples.
c---- Loop over subset of samples.
110 ns = n2 - n1 + 1
c.... Generate the needed random numbers.
c---- Loop over samples.
do 120 n = 1, ns
ranfp(n) = ranf( )
c---- End of loop over samples.
120 continue
c.... Randomly sample integers.
c---- Loop over samples.
do 130 n = 1, ns
nn = n + n1 - 1
fran = iminn + ranfp(n) * nint
fran = fran + min (0.0, sign (1.0, fran))
iran(nn) = fran
c---- End of loop over samples.
130 continue
c.... See if all subsets of samples are done.
c---- Do another subset of data.
if (n2 .lt. nunst) then
n1 = n2 + 1
n2 = min (nunst, n1 + 63)
c---- End of loop over subset of samples.
go to 110
endif
c---- Tested nunst.
endif
c.... See if any striated samples are needed.
c---- Need striated samples.
if (nstrd .gt. 0) then
c.... Set up the indices of the first subset of samples.
n1 = 1 + nunst
n2 = min (np, n1 + 63)
c.... Loop over the subset of samples.
c---- Loop over subset of samples.
140 ns = n2 - n1 + 1
c.... Generate the needed random numbers.
c---- Loop over samples.
do 150 n = 1, ns
ranfp(n) = ranf( )
c---- End of loop over samples.
150 continue
c.... Randomly sample integers, using striated sampling.
c---- Loop over samples.
do 160 n = 1, ns
nn = n + n1 - 1
ranst = (nn - nunst - ranfp(n)) / nstrd
fran = iminn + ranst * nint
fran = fran + min (0.0, sign (1.0, fran))
iran(nn) = fran
c---- End of loop over samples.
160 continue
c.... See if all subsets of samples are done.
c---- Do another subset of data.
if (n2 .lt. np) then
n1 = n2 + 1
n2 = min (np, n1 + 63)
c---- End of loop over subset of samples.
go to 140
endif
c---- Tested nstrd.
endif
210 continue
cbugc***DEBUG begins.
cbug 9904 format (/ 'aptrint results: nerr=',i3,
cbug & ' nunst=',i5,' nstrd=',i5)
cbug 9905 format (4(' n=',i5,' iran=',i5))
cbug write ( 3, 9904) nerr, nunst, nstrd
cbug if (nerr .ne. 0) return
cbug write ( 3, 9905) (n, iran(n), n = 1, np)
cbug
cbug 9906 format (/ ' summary by integers:')
cbug 9907 format (i5,' samples=',i5,' fraction=',1pe22.14)
cbug write ( 3, 9906)
cbugc---- Loop over integers.
cbug do 190 i = iminn, imaxx
cbug ns = 0
cbugc---- Loop over samples.
cbug do 185 n = 1, np
cbug if (iran(n) .eq. i) then
cbug ns = ns + 1
cbug endif
cbugc---- End of loop over samples.
cbug 185 continue
cbug fb = ns / (np + 1.e-99)
cbug if (ns .ge. 1) then
cbug write ( 3, 9907) i, ns, fb
cbug endif
cbugc---- End of loop over integers.
cbug 190 continue
cbugc***DEBUG ends.
return
c.... End of subroutine aptrint. (+1 line.)
end
UCRL-WEB-209832