subroutine aptvpap (ax, ay, az, bx, by, bz, np, tol,
& cx, cy, cz, clen, dx, dy, dz, dlen, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTVPAP
c
c call aptvpap (ax, ay, az, bx, by, bz, np, tol,
c cx, cy, cz, clen, dx, dy, dz, dlen, nerr)
c
c Version: aptvpap Updated 1997 April 18 12:00.
c aptvpap Originated 1997 April 18 12:00.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To find, for each of np pairs of vectors a = (ax, ay, az) and
c b = (bx, by, bz), the components of vector "a" parallel and
c perpendicular to vector "b", returned as vectors
c c = (cx, cy, cz) and d = (dx, dy, dz), with magnitudes clen
c and dlen.
c If the magnitude of vector "b" does not exceed tol, its
c orientation is indeterminate, and the returned values of vectors
c "c" and "d" will be zero, and clen and dlen will be -1.e+99.
c The components of vectors "c" or "d" may be truncated to zero,
c if less than the estimated error in their calculation, based on
c tol.
c Flag nerr indicates any input error.
c
c Input: ax, ay, az, bx, by, bz, np, tol.
c
c Output: cx, cy, cz, clen, dx, dy, dz, dlen, nerr.
c
c Calls: aptvdis, aptvdot, aptvunb
c
c
c Glossary:
c
c ax,ay,az Input The x, y and z components of vector "a". Size np.
c
c bx,by,bz Input The x, y and z components of vector "b". Size np.
c If the length of vector "b" does not exceed tol,
c clen and dlen will be returned as -1.e+99.
c
c clen Output The magnitude of vector "c". Size np.
c Zero if vector "a" is perpendicular to vector "b".
c Will be -1.e+99 if the length of vector "b" does
c not exceed tol.
c
c cx,cy,cz Output The x, y and z components of vector "c", the component
c of vector "a" parallel to vector "b". Size np.
c Each component may be truncated to zero, if less than
c the estimated error in its calculation, based on tol.
c Will be zero if the length of vector "b" does
c not exceed tol, but clen will be -1.e99l.
c
c dlen Output The magnitude of vector "d". Size np.
c Zero if vector "a" is parallel to vector "b".
c Will be -1.e+99 if the length of vector "b" does
c not exceed tol.
c
c dx,dy,dz Output The x, y and z components of vector "d", the component
c of vector "a" perpendicular to vector "b". Size np.
c Each component may be truncated to zero, if less than
c the estimated error in its calculation, based on tol.
c Will be zero if the length of vector "b" does
c not exceed tol, but dlen will be -1.e99.
c
c nerr Output Indicates an input error, if not 0. See clen, dlen.
c 1 if np is not positive.
c
c np Input Size of arrays.
c
c tol Input Numerical tolerance limit.
c On Cray computers, recommend 1.e-5 to 1.e-11.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Dimensioned arguments.
dimension ax (1) ! Component x of vector "a".
dimension ay (1) ! Component y of vector "a"
dimension az (1) ! Component z of vector "a"
dimension bx (1) ! Component x of vector "b".
dimension by (1) ! Component y of vector "b"
dimension bz (1) ! Component z of vector "b"
dimension cx (1) ! Component x of vector "c".
dimension cy (1) ! Component y of vector "c"
dimension cz (1) ! Component z of vector "c"
dimension clen (1) ! Length of vector "c".
dimension dlen (1) ! Length of vector "d".
dimension dx (1) ! Component x of vector "d".
dimension dy (1) ! Component y of vector "d"
dimension dz (1) ! Component z of vector "d"
c.... Local variables.
common /laptvpap/ n1 ! First index of subset of data.
common /laptvpap/ n2 ! Last index of subset of data.
common /laptvpap/ n ! Index in local array.
common /laptvpap/ nn ! Index in external array.
common /laptvpap/ ns ! Size of current subset of data.
common /laptvpap/ blen (64) ! Magnitude of vector "b".
common /laptvpap/ ubx (64) ! Component x of unit vector "ub".
common /laptvpap/ uby (64) ! Component y of unit vector "ub".
common /laptvpap/ ubz (64) ! Component z of unit vector "ub".
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptvpap finding components of A w/o B:' /
cbug & (i3,' ax,ay,az= ',1p3e22.14 /
cbug & ' bx,by,bz= ',1p3e22.14))
cbug write ( 3, 9901) (n, ax(n), ay(n), az(n),
cbug & bx(n), by(n), bz(n), n = 1, np)
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
c.... Test for input errors.
if (np .le. 0) then
nerr = 1
go to 999
endif
c.... Set up the indices of the first subset of data.
n1 = 1
n2 = min (np, 64)
110 ns = n2 - n1 + 1
c.... Find the unit vector "ub" parallel to vector "b", and
c.... blen, the magnitude of vector "b".
call aptvunb (bx(n1), by(n1), bz(n1), ns, 0.,
& ubx, uby, ubz, blen, nerr)
c.... Find the magnitude of vector "c".
call aptvdot (ax(n1), ay(n1), az(n1), ubx, uby, ubz, ns, tol,
& clen(n1), nerr)
c.... Find the components of vector "a" parallel to vector "b".
do 120 n = 1, ns ! Loop over subset of data.
nn = n + n1 - 1
cx(nn) = clen(nn) * ubx(n)
cy(nn) = clen(nn) * uby(n)
cz(nn) = clen(nn) * ubz(n)
if (blen(n) .le. tol) then
clen(nn) = -1.e99
endif
120 continue ! End of loop over subset of data.
c.... Find the components of vector "a" perpendicular to vector "b".
call aptvdis (cx(n1), cy(n1), cz(n1), ax(n1), ay(n1), az(n1),
& ns, tol, dx(n1), dy(n1), dz(n1), dlen(n1), nerr)
do 130 n = 1, ns ! Loop over subset of data.
nn = n + n1 - 1
if (blen(n) .le. tol) then
dx(nn) = 0.0
dy(nn) = 0.0
dz(nn) = 0.0
dlen(nn) = -1.e99
endif
130 continue ! End of loop over subset of data.
c.... See if all data subsets are done.
if (n2 .lt. np) then ! Do another subset of data.
n1 = n2 + 1
n2 = min (np, n1 + 63)
go to 110
endif
999 continue
cbugc***DEBUG begins.
cbug 9902 format (/ 'aptvpap results: nerr=',i3 /
cbug & (i3,' clen= ',1pe22.14 /
cbug & ' cx,cy,cz=',1p3e22.14 /
cbug & ' dlen= ',1pe22.14 /
cbug & ' dx,dy,dz=',1p3e22.14 ))
cbug write ( 3, 9902) nerr, (n, clen(n), cx(n), cy(n), cz(n),
cbug & dlen(n), dx(n), dy(n), dz(n), n = 1, np)
cbugc***DEBUG ends.
return
c.... End of subroutine aptvpap. (+1 line.)
end
UCRL-WEB-209832