subroutine aptrotc (au, av, bu, bv, cu, cv, np, tol, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTROTC
c
c call aptrotc (au, av, bu, bv, cu, cv, np, tol, nerr)
c
c Version: aptrotc Updated 1990 December 3 14:20.
c aptrotc Originated 1989 December 29 10:00.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To rotate the np 2-D vectors c = (cu, cv) around the w axis by
c the angle needed to make 2-D vector a = (au, av) parallel to
c 2-D vector b = (bu, bv). All are in the uv plane. Directions
c u, v and w are orthogonal. Any new components of c = (cu, cv)
c that are smaller than their estimated error, based on tol, will
c be truncated to zero. Flag nerr indicates any input error.
c
c History: 1990 March 14. Changed tol to 0.0 in call to unit vector
c subroutine. Allows small magnitudes.
c
c Input: au, av, bu, bv, cu, cv, tol.
c
c Output: cu, cv, nerr.
c
c Calls: aptvaxc, aptvdoc, aptvubc
c
c
c Glossary:
c
c au, av Input The u and v components of the first vector.
c The w components are zero.
c
c bu, bv Input The u and v components of the second vector.
c The w components are zero.
c
c cu, cv Input The u and v components of vector "c". Size np.
c The w components are zero.
c
c cu, cv Output The u and v components of vector "c", after
c rotation. The w components are zero. Will be
c truncated to zero, if smaller than the estimated
c error in their calculation, based on tol. Size np.
c
c nerr Output Indicates an input error, if not 0.
c 1 if np is not positive.
c 2 if vector "a" is too short.
c 3 if vector "b" is too short.
c
c np Input The size of arrays cu, cv.
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.
c---- Component u of vector "c".
dimension cu (1)
c---- Component v of vector "c".
dimension cv (1)
c.... Local variables.
c---- Cosine of angle from "a" to "b".
common /laptrotc/ costh
c---- Estimated error in cu.
common /laptrotc/ cuerr
c---- Initial value of cu.
common /laptrotc/ cus
c---- Estimated error in cv.
common /laptrotc/ cverr
c---- Initial value of cv.
common /laptrotc/ cvs
c---- Component u of first unit vector.
common /laptrotc/ du
c---- Component v of first unit vector.
common /laptrotc/ dv
c---- Component u of second unit vector.
common /laptrotc/ eu
c---- Component v of second unit vector.
common /laptrotc/ ev
c---- Index in local array.
common /laptrotc/ n
c---- Sine of angle from "a" to "b".
common /laptrotc/ sinth
c---- Magnitude of vector "a".
common /laptrotc/ vlena
c---- Magnitude of vector "b".
common /laptrotc/ vlenb
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptrotc. Rotate (au, av) parallel to (bu, bv):' /
cbug & ' au,av=',1p2e22.14 /
cbug & ' bu,bv=',1p2e22.14)
cbug 9902 format (i3,' cu,cv=',1p2e22.14)
cbug write ( 3, 9901) au, av, bu, bv
cbug write ( 3, 9902) (n, cu(n), cv(n), n = 1, np)
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
c.... Test for input errors.
c---- No input data.
if (np .le. 0) then
nerr = 1
go to 210
endif
c.... Find the unit vector "d" parallel to 2-D vector "a".
call aptvubc (au, av, 1, 0., du, dv, vlena, nerr)
c---- Vector "a" is too short.
if (vlena .le. tol) then
nerr = 2
go to 210
endif
c.... Find the unit vector "e" parallel to 2-D vector "b".
call aptvubc (bu, bv, 1, 0., eu, ev, vlenb, nerr)
c---- Vector "b" is too short.
if (vlenb .le. tol) then
nerr = 3
go to 210
endif
c.... Rotate 2-D vector "c".
c---- No truncation test.
if (tol .eq. 0.0) then
costh = du * eu + dv * ev
sinth = du * ev - dv * eu
c---- Loop over vectors.
do 110 n = 1, np
cus = cu(n)
cvs = cv(n)
cu(n) = costh * cus - sinth * cvs
cv(n) = sinth * cus + costh * cvs
c---- End of loop over vectors.
110 continue
c---- Truncate small components to zero.
else
c.... Find the cosine of the angle of rotation.
c.... costh = d .dot. e (scalar product).
call aptvdoc (du, dv, eu, ev, 1, tol, costh, nerr)
c.... Find the sine of the angle of rotation.
c.... sinth = w component of d x e (cross product).
call aptvaxc (du, dv, eu, ev, 1, tol, sinth, nerr)
c---- Loop over vectors.
do 120 n = 1, np
cus = cu(n)
cvs = cv(n)
cu(n) = costh * cus - sinth * cvs
cv(n) = sinth * cus + costh * cvs
cuerr = tol * (abs (costh * cus) + abs (sinth * cvs))
cverr = tol * (abs (sinth * cus) + abs (costh * cvs))
if (abs (cu(n)) .lt. cuerr) then
cu(n) = 0.0
endif
if (abs (cv(n)) .lt. cverr) then
cv(n) = 0.0
endif
c---- End of loop over vectors.
120 continue
c---- Tested tol.
endif
cbugc***DEBUG begins.
cbug 9903 format (/ 'aptrotc results:')
cbug write ( 3, 9903)
cbug write ( 3, 9902) (n, cu(n), cv(n), n = 1, np)
cbugc***DEBUG ends.
210 return
c.... End of subroutine aptrotc. (+1 line.)
end
UCRL-WEB-209832