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