subroutine aptcomf (n, m, tol, ncomb, comb, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTCOMF
c
c call aptcomf (n, m, tol, ncomb, comb, nerr)
c
c Version: aptcomf Updated 2001 November 28 16:00.
c aptcomf Originated 2001 November 28 16:00.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: to find C(n,m), the number of combinations of n things taken
c m at a time, and return the integer and floating point values
c of C(n,m), ncomb and comb, respectively.
c
c Note: The formula is C(n,m) = ncomb = comb = n! / ((n - m)! * m!),
c or C(n,m) = (n - m + 1) * (n - m + 2) * ... * (n - 1) * n /
c 1 * 2 * ... * (m - 1) * m.
c
c C(n,m), (m=0,n) are the coefficients of the binomial
c expansion (x + y)**n = sum(m=0,n) {C(n,m) * x**(n-m) * y**m}.
c
c Input: n, m, tol.
c
c Output: ncomb, comb, nerr.
c
c Glossary:
c
c comb Output The floating point value of n! / (n - m)!*m!).
c If n < 0, returned as -1.e99, with nerr = 1.
c WARNING: may overflow largest machine floating point
c number, if n is too large.
c
c n Input Total number of things.
c
c nerr Output Indicates an input or result error, if not zero.
c -1 if ncomb exceeds the largest machine integer.
c 1 if n < 0.
c 2 if m < 0, or m > n.
c
c ncomb Output The integer value of n! / (n - m)!.
c If n < 0, returned as -999999, with nerr = 1.
c If larger then the largest machine integer,
c returned as -999999, with nerr = -1, but comb will
c still be accurate, unless comb overflows the largest
c machine floating point number.
c
c m Input The number of things selected out of n, for each
c combination.
c
c tol Input The precision with which ncomb and comb must agree,
c to insure that ncomb does not exceed the largest
c machine integer.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
implicit none
c.... Arguments.
integer n ! Total number of things.
integer nerr ! Error indicator.
integer ncomb ! Integer value of n! / (n - m)!.
integer m ! Items taken out of n, each combination.
real comb ! Real value of n! / (n - m)!.
real tol ! Reqired precision of comb, ncomb match.
c.... Local variables.
integer nbeg ! n - m + 1.
integer nn ! Index, nbeg to n.
real fden ! Real value of m! or (n - m)!..
real fn ! Real value of nn - nbeg + 1.
real fnn ! Real value of nn.
real fncomb ! Real value of ncomb.
cbugc***DEBUG begins.
cbug 9001 format (/ 'aptcomf finding combinations of (',i20,',',i20,').' /
cbug & 'tol = ',1pe22.14)
cbug write ( 3, 9001) n, m, tol
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
ncomb = -999999
fncomb = -999999.0
comb = -1.e99
c.... Test for errors, limits.
if (n .lt. 0) then
nerr = 1
go to 210
endif
if ((m .lt. 0) .or. (m .gt. n)) then
nerr = 2
go to 210
endif
c.... Find number of combinations of n things taken m at a time.
ncomb = 1.0
fncomb = 1.0
comb = 1.0
fden = 1.0
if (m .eq. n) go to 210
if (n .gt. 1) then
nbeg = max (m, n - m) + 1
fn = 0.0
do 110 nn = nbeg, n
fn = fn + 1.0
fden = fden * fn
fnn = nn
comb = comb * fnn
110 continue
comb = comb / fden
ncomb = nint (comb)
fncomb = ncomb
if (abs (fncomb - comb) .gt. tol * comb) then
nerr = -1
ncomb = -999999
endif
endif ! Tested n.
210 continue
cbugc***DEBUG begins.
cbug 9002 format (/ 'aptcomf found # of combinations. nerr = ',i3 /
cbug & 'nc, fnc, chk =',i20,2x,1p2e22.14)
cbug write ( 3, 9002) nerr, ncomb, comb, fncomb
cbugc***DEBUG ends.
return
end
UCRL-WEB-209832