subroutine aptcmad (int, nint, id, im, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTCMAD
c
c call aptcmad (int, nint, id, im, nerr)
c
c Version: aptcmad Updated 1903 September 12 13:50.
c aptcmad Originated 1903 September 10 18:00.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: For the set of integers int(n), n = 1, nint, find id, the
c greatest common divisor (GCD) and im, the least common multiple
c (LCM). To find the factors of any of these, use subroutine
c aptfact.
c
c Flag nerr indicates any input error, or lack of result.
c
c Input: int, nint.
c
c Output: id, im, nerr.
c
c Glossary:
c
c id Output The GCD of the integers in set int.
c It contains only those prime factors, including 1,
c found in every integer in the set int, with the
c lowest power found for each such factor.
c
c im Output The LCM of the integers in set int.
c It contains every prime factor found in the set int,
c with the highest power found for each such factor.
c If im is larger than 1.e18, nerr = 4, and im is
c replaced with -999999999.
c
c int Input A set of integers int(n), n = 1, nint, for which
c id, the GCD, and im, the LCM are to be found.
c Must all be nonzero. If not, nerr = 2.
c If any is not factorable, nerr = 3.
c Size at least nint.
c
c nerr Output Indicates an input or result error, if not zero.
c 1 if nint < 1.
c 2 if any number in set int is not positive.
c 3 if any number in set int can not be factored.
c 4 if im exceeds the largest machine integer.
c
c nint Output The number of integers in the set int. Must be
c positive. If not, nerr = 1.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
implicit none
c.... Arguments.
integer int(1) ! Nonzero integers int(n), n = 1, nint.
integer id ! GCD of set int.
integer im ! LCM of set int.
integer nerr ! Error flag, none if zero.
integer nint ! Size of set of integers int.
c.... Local variables.
integer ifd(20) ! A prime factor of id(nd), nd = 1, nfd.
integer ifi(20) ! A prime factor of an integer int(n).
integer ifm(20) ! A prime factor of im(nm), nm = 1, nfm.
integer imsave ! Value of im before testing size.
integer inta ! Positive integers abs (int(n)).
integer ipd(20) ! Power of ifd(nd), nd = 1, nfd.
integer ipi(20) ! Power of ifi(ni), ni = 1, nfi.
integer ipm(20) ! Power of ifm(nm), nm = 1, nfm.
integer ir ! Residue of imcomplete factorization.
integer it ! Totient function. Not used.
integer n ! Index, 1 to nint.
integer nd ! Index, 1 to nfd.
integer ni ! Index, 1 to nfi.
integer nm ! Index, 1 to nfm.
integer nerra ! Error flag from aptpfac, none if zero.
integer nfd ! How many prime factors of id found.
integer nfi ! How many prime factors of int found.
integer nfm ! How many prime factors of im found.
real zim ! Floating point value of im.
cbugc***DEBUG begins.
cbug 9001 format (/ 'aptcmad finding GCD and LCM of ',i6,' values:' /
cbug & (i20))
cbug write ( 3, 9001) nint, (int(n), n = 1, nint)
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
id = -999999999
im = -999999999
c.... Test for errors, limits.
if (nint .lt. 1) then
nerr = 1
go to 210
endif
do 110 n = 1, nint
if (int(n) .eq. 0) then
nerr = 2
go to 210
endif
110 continue
c.... Initialize.
nfi = 0
nfd = 0
nfm = 0
c.... Find prime factors of each integer in set int.
do 170 n = 1, nint ! Loop over set int.
inta = abs (int(n))
call aptpfac (inta, 20, ifi, ipi, nfi, it, ir, nerra)
if (nerra .ne. 0) then
nerr = 3
go to 210
endif
cbugcbugc***DEBUG begins.
cbugcbug 9500 format ('II* n =',i3,' INT factor =',i20,'^',i2)
cbugcbug write ( 3, 9010)
cbugcbug write ( 3, 9500) (n, ifi(ni), ipi(ni), ni = 1, nfi)
cbugcbugc***DEBUG ends.
c.... Add any new factors to im.
if (n .eq. 1) then ! Only when n = 1.
c.... Initialize factors of GCD and LCM to those of first integer.
nfd = nfi
nfm = nfi
do 120 ni = 1, nfi ! Loop over factors of int.
ifd(ni) = ifi(ni)
ipd(ni) = ipi(ni)
ifm(ni) = ifi(ni)
ipm(ni) = ipi(ni)
120 continue ! End of loop over factors of int.
cbugcbugc***DEBUG begins.
cbugcbug 9501 format ('L1* n =',i3,' LCM factor =',i20,'^',i2)
cbugcbug write ( 3, 9010)
cbugcbug write ( 3, 9501) (n, ifm(nm), ipm(nm), nm = 1, nfm)
cbugcbugc***DEBUG ends.
cbugcbugc***DEBUG begins.
cbugcbug 9502 format ('B1* n =',i3,' GCD factor =',i20,'^',i2)
cbugcbug write ( 3, 9010)
cbugcbug write ( 3, 9502) (n, ifd(nd), ipd(nd), nd = 1, nfd)
cbugcbugc***DEBUG ends.
else ! When n is not 1.
c.... Find any unshared factors of GCD, use lower power on shared factors.
do 140 nd = 1, nfd ! Loop over factors of GCD.
do 130 ni = 1, nfi ! Loop over factors of int.
if (ifi(ni) .eq. ifd(nd)) then ! Shared factor.
if (ipi(ni) .lt. ipd(nd)) then ! Lower power.
ipd(nd) = ipi(ni)
endif
go to 140
endif ! End of testing factor.
130 continue ! End of loop over factors of int.
c.... This factor not in this integer. Remove from GCD.
ifd(nd) = 1
ipd(nd) = 1
140 continue ! End of loop over factors of GCD.
cbugcbugc***DEBUG begins.
cbugcbug 9504 format ('D2* n =',i3,' GCD factor =',i20,'^',i2)
cbugcbug write ( 3, 9010)
cbugcbug write ( 3, 9504) (n, ifd(nd), ipd(nd), nd = 1, nfd)
cbugcbugc***DEBUG ends.
c.... Find any new factors of LCM, use higher power on shared factors.
do 160 ni = 1, nfi ! Loop over factors of int.
do 150 nm = 1, nfm ! Loop over factors of LCM.
if (ifi(ni) .eq. ifm(nm)) then ! Shared factor.
if (ipi(ni) .gt. ipm(nm)) then ! Higher power.
ipm(nm) = ipi(ni)
endif
go to 160
endif ! End of testing factor.
150 continue ! End of loop over factors of LCM.
c.... This factor not in LCM. Add to LCM.
nfm = nfm + 1
ifm(nfm) = ifi(ni)
ipm(nfm) = ipi(ni)
160 continue ! End of loop over factors of int.
cbugcbugc***DEBUG begins.
cbugcbug 9503 format ('L2* n =',i3,' LCM factor =',i20,'^',i2)
cbugcbug write ( 3, 9010)
cbugcbug write ( 3, 9503) (n, ifm(nm), ipm(nm), nm = 1, nfm)
cbugcbugc***DEBUG ends.
endif ! Tested n.
170 continue ! End of loop over set int.
c.... Construct GCD.
id = 1
do 180 nd = 1, nfd
id = id * ifd(nd)**ipd(nd)
180 continue
c.... Construct LCM.
im = 1
zim = 1
do 190 nm = 1, nfm
im = im * ifm(nm)**ipm(nm)
zim = zim * ifm(nm)**ipm(nm)
190 continue
imsave = im
if (zim .gt. 1.e18) then
im = -999999999
nerr = 4
endif
210 continue
cbugc***DEBUG begins.
cbug 9002 format (/ 'aptcmad results: nerr =',i3 /
cbug & ' GCD =',i20 /
cbug & ' LCM =',i20 )
cbug 9003 format (' INT factor =',i20,'^',i2)
cbug 9004 format (' GCD factor =',i20,'^',i2)
cbug 9005 format (' LCM factor =',i20,'^',i2)
cbug 9006 format (' LCM (bad?) =',i20 )
cbug 9010 format (/)
cbug write ( 3, 9002) nerr, id, im
cbug write ( 3, 9010)
cbug write ( 3, 9004) (ifd(nd), ipd(nd), nd = 1, nfd)
cbug write ( 3, 9010)
cbug write ( 3, 9005) (ifm(nm), ipm(nm), nm = 1, nfm)
cbug if (nerr .eq. 4) then
cbug write ( 3, 9010)
cbug write ( 3, 9006) imsave
cbug endif
cbugc***DEBUG ends.
return
c.... End of subroutine aptcmad. (+1 line.)
end
UCRL-WEB-209832