subroutine aptafac (inum, nfactm, ifact, nfact, inumr, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTAFAC
c
c call aptafac (inum, nfactm, ifact, nfact, inumr, nerr)
c
c Version: aptafac Updated 2005 September 23 15:45.
c aptafac Originated 2005 September 23 15:45.
c
c Authors: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: For the integer inum, find all factors ifact(n) > 1,
c n = 1, nfact, with nfact limited to nfactm, and with each
c factor used only once, except possibly the last factor.
c If not all factors of inum are returned, because nfactm is
c insufficient, then inumr, the remaining unfactored part of inum,
c is returned with a value other than 1, with nerr = 3.
c
c The method is slow, but avoids storing a large number of
c numbers. Each successive integer, 2, 3, 4, 5, 6, is tested as a
c factor one time, and the result of dividing inum by each factor
c used to limit the search to integers no larger than the square
c root of that result.
c
c Flag nerr indicates any input error or limited results.
c
c Input: inum, nfactm.
c
c Output: ifact, nfact, inumr, nerr.
c
c Calls: (none)
c
c Glossary:
c
c ifact Output Factor ifact(n) is the n'th smallest factor of inum,
c with each factor only used once, except possibly
c the last.
c If inum = 1, then ifact(1) = 1.
c Size at least nfactm. Initially set to -999999999.
c
c inum Input Any positive integer, for which all factors are to be
c found, up to a maximum of nfactm factors.
c
c inumr Output Any remaining unfactored part of inum, when nfact
c is prevented from exceeding nfactm. Normally 1.
c
c nerr Output If not zero, indicates an error.
c 1 if inum is not 1 or more.
c 2 if nfactm is not positive.
c 3 if more than nfactm factors exist.
c
c nfact Output The number of factors found for integer inum.
c None may be used more than once.
c If nfactm are found, and more exist, nerr = 3.
c
c nfactm Input The maximum allowed size of arrays ifact.
c On a machine with 64-bit integers, 18 is enough.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Declarations for arguments.
integer ifact(1) ! A factor of inum.
integer inum ! Any positive integer.
integer inumr ! Result of inum divided by all factors.
integer nerr ! Indicates an error, if not zero.
integer nfact ! How many factors of inum found.
integer nfactm ! The maximum allowed value of nfact.
c.... Declarations for local variables.
integer itest ! Integer to test as factor.
integer itestm ! Maximum integer to test as factor.
integer n ! Index in arrays ifact.
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptafac finding factors of inum =',i20 /
cbug & ' but no more than',i6,' factors.' )
cbug write ( 3, 9901) inum, nfactm
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
nfact = 0
inumr = inum
if (nfactm .gt. 0) then
do 100 n = 1, nfactm
ifact(n) = -999999999
100 continue
endif
c.... Test for input errors.
if (inum .lt. 1) then
nerr = 1
go to 210
endif
if (nfactm .le. 0) then
nerr = 2
go to 210
endif
if (inum .eq. 1) then
nfact = 1
ifact(1) = 1
go to 210
endif
c.... Find maximum possible factor.
itestm = sqrt (float (inumr))
cbugcbugc***DEBUG begins.
cbugcbug 9700 format ('DEBUG 0. nfact, itestm =',i3,21x, i18)
cbugcbug write ( 3, 9700) nfact, itestm
cbugcbugc***DEBUG ends.
c.... Test 2.
itest = 2
if (itest .gt. itestm) then ! 2 is too high.
nfact = nfact + 1
cbugcbugc***DEBUG begins.
cbugcbug 9801 format ('DEBUG A. nfact, itestm =',i3,21x, i18)
cbugcbug write ( 3, 9801) nfact, itestm
cbugcbugc***DEBUG ends.
ifact(nfact) = inumr
inumr = 1
go to 210
endif ! 2 is too high.
if (mod (inumr, itest) .eq. 0) then ! Residual inumr is even.
if (nfact .ge. nfactm) then
nerr = 3
go to 210
endif
nfact = nfact + 1
cbugcbugc***DEBUG begins.
cbugcbug 9802 format ('DEBUG B. nfact, itestm =',i3,21x, i18)
cbugcbug write ( 3, 9802) nfact, itestm
cbugcbugc***DEBUG ends.
ifact(nfact) = itest
inumr = inumr / itest
if (inumr .eq. 1) go to 210
itestm = sqrt (float (inumr))
cbugcbugc***DEBUG begins.
cbugcbug 9701 format ('DEBUG 1. nfact,ifact,itestm =',i3,i18,i18)
cbugcbug write ( 3, 9701) nfact, ifact(nfact), itestm
cbugcbugc***DEBUG ends.
endif
if (nfact .ge. nfactm) then
nerr = 3
go to 210
endif
c.... Test each higher integer.
itest = 3
130 if (itest .gt. itestm) then
nfact = nfact + 1
cbugcbugc***DEBUG begins.
cbugcbug 9803 format ('DEBUG C. nfact, itestm =',i3,21x, i18)
cbugcbug write ( 3, 9803) nfact, itestm
cbugcbugc***DEBUG ends.
ifact(nfact) = inumr
inumr = 1
go to 210
endif
if (mod (inumr, itest) .eq. 0) then ! Test next integer.
nfact = nfact + 1
cbugcbugc***DEBUG begins.
cbugcbug 9804 format ('DEBUG C. nfact, itestm =',i3,21x, i18)
cbugcbug write ( 3, 9804) nfact, itestm
cbugcbugc***DEBUG ends.
ifact(nfact) = itest
inumr = inumr / itest
if (inumr .eq. 1) go to 210
itestm = sqrt (float (inumr))
cbugcbugc***DEBUG begins.
cbugcbug 9703 format ('DEBUG 3. nfact,ifact,itestm =',i3,i18,i18)
cbugcbug write ( 3, 9703) nfact, ifact(nfact), itestm
cbugcbugc***DEBUG ends.
endif ! Tested next odd integer.
if (nfact .ge. nfactm) then
nerr = 3
go to 210
endif
itest = itest + 1
go to 130
210 continue
cbugc***DEBUG begins.
cbug 9902 format (/ 'aptafac found',i12,' factors. nerr =',i3,'.' )
cbug 9903 format (' n =',i3,' ifact(n) =',i22)
cbug 9905 format (' Unfactored residue =',i20)
cbug write ( 3, 9902) nfact, nerr
cbug if (nfact .gt. 0) then
cbug write ( 3, 9903) (n, ifact(n), n = 1, nfact)
cbug endif
cbug if (nerr .eq. 3) then
cbug write ( 3, 9905) inumr
cbug endif
cbugc***DEBUG ends.
return
c.... End of subroutine aptafac. (+1 line.)
end
UCRL-WEB-209832