subroutine aptbpow (npow, nbase, idiga, ndiga,
& idigbp, idigbw, ndigm,
& idigb, ndigb, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTBPOW
c
c call aptbpow (npow, nbase, idiga, ndiga, idigbp, idigbw,
c & ndigm, idigb, ndigb, nerr)
c
c Version: aptbpow Updated 2006 July 6 15:30.
c aptbpow Originated 2005 May 23 14:50.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To find the npow'th power of the big integer "a", stored as
c the array of non-negative base nbase digits idiga, of length
c ndiga, and store the result as the array of non-negative base
c nbase digits idigb, of length ndigb. Integer arrays idigbp and
c idigbw, of length ndigm, are needed for intermediate results.
c
c The digits of idiga and idigb are the base nbase
c digits representing the non-negative decimal integers ideca
c and idecb, in order from most significant to least
c significant, using the equations:
c ideca = sum (idiga(n) * nbase^(N-n), n = 1, N = ndiga),
c idecb = sum (idigb(n) * nbase^(N-n), n = 1, N = ndigb),
c idecb = ideca^npow.
c See aptbrev to reverse the order of the digits.
c
c Flag nerr indicates any input error.
c
c See aptbtod, aptdtob, aptbadd, aptbsub, aptbmul, aptbdiv,
c aptbrev, aptbrtn, aptbfac.
c
c Input: nbase, idiga, ndiga, idigbp, idigbw, ndigm.
c
c Output: idigb, ndigb, nerr.
c
c Calls: aptbmul, aptbrev
c
c Glossary:
c
c idiga Input A big integer "a", stored as an array of ndiga
c non-negative base nbase digits, each representing a
c single "digit" in the base nbase representation of
c "a", from most to least significant. If nbase > 10,
c each "digit" may require 2 or more integers.
c For example, for decimal integer 4821 and
c nbase = 16 (hexadecimal), idiga(n) = (5, 13, 2, 1),
c with ndiga = 4, or
c 4821 (dec) = 5 * 1 + 13 * 16 + 2 * 256 + 1 * 4096
c
c idigb Input See idiga. Array idigb must not be idiga.
c Big integer "b" is the npow'th power of "a".
c
c idigbp Input See idiga. Temporary working space for idigb in
c subroutine aptbpow. Array idigbp may not be idiga,
c idigbw or idigb, but must have its own memory space,
c which must be at least npow * ndiga.
c
c idigbw Input See idiga. Temporary working space for idigb in
c subroutine aptbmul. Array idigbw must not be
c idiga or idigb, but must have its own memory space,
c which must be at least npow * ndiga.
c
c idigw Input Temporary for subroutine aptbmul. Must have its own
c memory space, at least ndigm.
c
c nbase Input The number base for which the digits of integer arrays
c idiga and idigb are the coefficients of
c the powers of nbase, in order from least to most
c significant.
c
c ndiga Input The length of the integer array idiga.
c
c ndigb Output The length of the integer array idigb.
c Memory space must be at least npow * ndiga.
c
c ndigm Input The maximum number of words allowed in integer arrays
c idigbp and idigpw. Memory space for idigb must be
c at least ndigm, and ndigm must be at least as big as
c npow * ndiga.
c
c nerr Output Indicates an input or a result error, if not zero.
c 1 if npow is negative.
c 2 if nbase is less than 2.
c 3 if ndiga is not positive.
c 4 if any digits idiga are negative.
c 5 if idiga and npow are both zero.
c 6 if ndiga exceeds ndigm.
c 7 if ndigb exceeds ndigm.
c
c npow Input The power to which big integer "a" is to be raised.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Dimensioned arguments.
integer idiga(1)
integer idigb(1)
integer idigbp(1)
integer idigbw(1)
c.... Local variables.
cbugc***DEBUG begins.
cbug 9900 format (/)
cbug 9901 format (/ 'aptbpow finding power',i7,' of the base nbase =',i7 /
cbug & ' digit array idiga of length',i7 /
cbug & ' to get base nbase digit array idigb.' /
cbug & ' ndigm =',i7,'.')
cbug 9902 format (' a n =',i7,', idiga =',i7,'.')
cbug 9904 format (' b n =',i7,', idigb =',i7,'.')
cbug 9906 format (' bp n =',i7,', idigbp =',i7,'.')
cbug write ( 3, 9901) npow, nbase, ndiga, ndigm
cbug write ( 3, 9900)
cbug write ( 3, 9902) (n, idiga(n), n = 1, ndiga)
cbugc***DEBUG ends.
c.... Initialize.
c.... Change "a" to least significant digit first.
call aptbrev (idiga, ndiga, nerrb)
cbugcbugc***DEBUG begins.
cbugcbug 9913 format (/,'aptbpow nbug=',i3,' nerr=',i3)
cbugcbug nbug = 1
cbugcbug write ( 3, 9913) nbug, nerrb
cbugcbugc***DEBUG ends.
nerr = 0
c.... Test for input errors.
if (npow .lt. 0) then
nerr = 1
go to 210
endif
if (nbase .lt. 2) then
nerr = 2
go to 210
endif
if (ndiga .lt. 0) then
nerr = 3
go to 210
endif
if (ndiga .ge. 1) then
do n = 1, ndiga
if (idiga(n) .lt. 0) then
nerr = 4
go to 210
endif
enddo
endif
if ((npow .eq. 0) .and. (ndiga .le. 1) .and.
& (idiga(1) .eq. 0)) then
nerr = 5
go to 210
endif
if (ndiga .gt. ndigm) then
nerr = 6
goto 210
endif
c.... Initialize "bp" and "bw" to zero.
do n = 1, ndigm
idigbp(n) = 0
idigbw(n) = 0
enddo
c.... Find the specified power of "a".
if (npow .eq. 0) then
ndigbp = 1
idigbp(1) = 1
go to 210
endif
if (npow .eq. 1) then
ndigbp = ndiga
do n = 1, ndigbp
idigbp(n) = idiga(n)
enddo
go to 210
endif
c.... Multiply "a" by itself to get "bp" = "a" squared.
c.... Change "a" to most significant digit first.
call aptbrev (idiga, ndiga, nerrb)
cbugcbugc***DEBUG begins.
cbugcbug nbug = 22
cbugcbug write ( 3, 9913) nbug, nerrb
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9902) (n, idiga(n), n = 1, ndiga)
cbugcbugc***DEBUG ends.
call aptbmul (nbase, idiga, ndiga, idiga, ndiga,
& idigbw, ndigm, idigbp, ndigbp, nerra)
cbugcbugc***DEBUG begins.
cbugcbug nbug = 4
cbugcbug write ( 3, 9913) nbug, nerra
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9902) (n, idiga(n), n = 1, ndiga)
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9906) (n, idigbp(n), n = 1, ndigbp)
cbugcbugc***DEBUG ends.
c.... Change "a" and "bp" to least significant digit first.
call aptbrev (idiga, ndiga, nerrb)
call aptbrev (idigbp, ndigbp, nerrb)
if (nerra .ne. 0) then
nerr = 7
go to 210
endif
if (npow .ge. 3) then ! Need cube or higher power.
do np = 3, npow ! Loop over powers.
c.... Multiple last power of "a" by "a" to get next power.
c.... Change "a" and "bp" to most significant digit first.
call aptbrev (idiga, ndiga, nerrb)
call aptbrev (idigbp, ndigbp, nerrb)
cbugcbugc***DEBUG begins.
cbugcbug nbug = 6
cbugcbug write ( 3, 9913) nbug, nerrb
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9902) (n, idiga(n), n = 1, ndiga)
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9906) (n, idigbp(n), n = 1, ndigbp)
cbugcbugc***DEBUG ends.
call aptbmul (nbase, idiga, ndiga, idigbp, ndigbp,
& idigbw, ndigm, idigbp, ndigbp, nerra)
cbugcbugc***DEBUG begins.
cbugcbug nbug = 61
cbugcbug write ( 3, 9913) nbug, nerra
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9902) (n, idiga(n), n = 1, ndiga)
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9906) (n, idigbp(n), n = 1, ndigbp)
cbugcbugc***DEBUG ends.
c.... Change "a" and "bp" to least significant digit first.
call aptbrev (idiga, ndiga, nerrb)
call aptbrev (idigbp, ndigbp, nerrb)
if (nerra .ne. 0) then
nerr = 7
go to 210
endif
enddo ! End of loop over powers.
endif ! Needed cube or higher power.
210 continue
c.... Change "a", "bp" and "b" to most significant digit first.
call aptbrev (idiga, ndiga, nerrb)
call aptbrev (idigbp, ndigbp, nerrb)
c.... Copy "bp" into "b".
if (nerr .eq. 0) then
ndigb = ndigbp
do n = 1, ndigb
idigb(n) = idigbp(n)
enddo
endif
cbugcbugc***DEBUG begins.
cbugcbug nbug = 8
cbugcbug write ( 3, 9913) nbug, nerrb
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9902) (n, idiga(n), n = 1, ndiga)
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9906) (n, idigbp(n), n = 1, ndigbp)
cbugcbug write ( 3, 9900)
cbugcbug write ( 3, 9904) (n, idigb(n), n = 1, ndigb)
cbugcbugc***DEBUG ends.
if (nerr .ne. 0) then
ndigb = 0
do n = 1, ndigm
idigb(n) = -999999
enddo
endif
cbugc***DEBUG begins.
cbug 9911 format (/ 'aptbpow results: nerr=',i2,', ndigb =',i7,'.')
cbug 9914 format (' a^npow n =',i7,', idigb =',i7,'.')
cbug write ( 3, 9911) nerr, ndigb
cbug write ( 3, 9900)
cbug write ( 3, 9914) (n, idigb(n), n = 1, ndigb)
cbugc***DEBUG ends.
return
c.... End of subroutine aptbpow. (+1 line.)
end
UCRL-WEB-209832