subroutine aptclis (px, py, pz, ux, uy, uz, rad, tol, qc, qx, qy,
& qz, qxy, qyz, qzx, qxx, qyy, qzz, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTCLIS
c
c call aptclis (px, py, pz, ux, uy, uz, rad, tol, qc, qx, qy, qz,
c & qxy, qyz, qzx, qxx, qyy, qzz, nerr)
c
c Version: aptclis Updated 1993 April 29 11:20.
c aptclis Originated 1993 April 29 11:20.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To find the coefficients of the equation for a right circular
c cylinder, given point p = (px, py, pz) on the axis, the axial
c direction vector u = (ux, uy, uz), and the radius rad.
c Flag nerr indicates an input error, if not zero.
c
c The equation for the cylinder is as follows:
c
c f(x, y, z) = qc + qx * x + qy * y + qz * z +
c qxy * x * y + qyz * y * z + qzx * z * x +
c qxx * x**2 + qyy * y**2 + qzz * z**2 = 0
c
c Input: px, py, pz, ux, uy, uz, rad, tol.
c
c Output: qc, qx, qy, qz, qxy, qyz, qzx, qxx, qyy, qzz, nerr.
c
c Calls: aptrois, aptrotv, apttris
c
c
c Glossary:
c
c nerr Output Indicates an input error, if not zero:
c 1 if the magnitude of vector "u" is no more than tol.
c 2 if rad is no more than tol times the length of
c vector "p".
c
c px,py,pz Input The x, y, z coordinates of point "p" on the axis of
c the cylinder.
c
c qc Output Constant term in the equation of the cylinder.
c Zero only if the cylinder passes through the origin.
c
c qx,qy,qz Output Coefficients of x, y, z in the equation of the
c cylinder. All zero only if the axis of the cylinder
c passes through the origin.
c
c q.. Output Coefficients of the second-order terms in the equation
c of the cylinder. qxy, qyz, qzx, qxx, qyy, qzz are
c the coefficients of x*y, y*z, z*x, x**2, y**2, z**2,
c respectively. The coefficients qxy, qyz, qzx are all
c zero only if the axis of the cylinder is aligned with
c a major axis.
c
c tol Input Numerical tolerance limit.
c On Cray computers, recommend 1.e-5 to 1.e-11.
c
c ux,uy,uz Input The x, y, z components of the vector "u" in the
c direction of the axis of the cylinder.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Local variables.
c---- Coordinate transformation matrix.
common /laptclis/ rotm (3,3)
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptclis finding equation of a cylinder.',
cbug & ' tol=',1pe13.5)
cbug 9902 format (/
cbug & ' px,py,pz= ',1p3e22.14 /
cbug & ' us,uy,uz= ',1p3e22.14 /
cbug & ' rad= ',1pe22.14)
cbug write ( 3, 9901) tol
cbug write ( 3, 9902) px, py, pz, ux, uy, uz, rad
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
qc = 0.0
qx = 0.0
qy = 0.0
qz = 0.0
qxy = 0.0
qyz = 0.0
qzx = 0.0
qxx = 0.0
qyy = 0.0
qzz = 0.0
if ((ux**2 + uy**2 + uz**2) .le. tol**2) then
nerr = 1
go to 210
endif
if (rad**2 .le. tol**2 * (px**2 + py**2 + pz**2)) then
nerr = 2
go to 210
endif
c.... Find the rotation matrix to rotate the z axis to direction "u".
call aptrotv (0., 0., 1., ux, uy, uz, tol, rotm, nerr)
c.... Find the coefficients of a right circular cylinder on the z axis, with
c.... radius rad: -rad**2 + x**2 +y**2 = 0.
qc = -rad**2
qxx = 1.0
qyy = 1.0
c.... Rotate the cylinder around the origin to point the axis in direction "u".
call aptrois (rotm, 0, 0., 0., 0., qc, qx, qy, qz, qxy, qyz, qzx,
& qxx, qyy, qzz, tol, nerr)
c.... Translate the origin to point "p".
call apttris (1, px, py, pz, qc, qx, qy, qz, qxy, qyz, qzx,
& qxx, qyy, qzz, tol, nerr)
210 continue
cbugc***DEBUG begins.
cbug 9903 format (/ 'aptclis results: nerr=',i3 /
cbug & ' qc= ',1pe22.14 /
cbug & ' qx,qy,qz= ',1p3e22.14 /
cbug & ' qxy,qyz,qzx=',1p3e22.14 /
cbug & ' qxx,qyy,qzz=',1p3e22.14)
cbug write ( 3, 9903) nerr, qc, qx, qy, qz, qxy, qyz, qzx,
cbug & qxx, qyy, qzz
cbugc***DEBUG ends.
return
c.... End of subroutine aptclis. (+1 line.)
end
UCRL-WEB-209832