subroutine aptqnor (ax, ay, az, qc, qx, qy, qz, qxy, qyz, qzx,
& qxx, qyy, qzz, tol, sc, sx, sy, sz, slen,
& nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTQNOR
c
c call aptqnor (ax, ay, az, qc, qx, qy, qz, qxy, qyz, qzx,
c qxx, qyy, qzz, tol, sc, sx, sy, sz, slen, nerr)
c
c Version: aptqnor Updated 1997 October 13 10:50.
c aptqnor Originated 1997 October 13 10:50.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: For point "a" = (ax, ay, az), and the family of quadric
c surfaces represented by f(x,y,z) = c (for arbitrary c):
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,
c
c to find the components sx, sy and sz and magnitude slen
c of the normal vector "s" = (sx,sy,sz), where s = grad f:
c
c sx = qx + 2.0 * qxx * x + qxy * y + qzx * z
c sy = qy + qxy * x + 2.0 * qyy * y + qyz * z
c sz = qz + qzx * x + qyz * y + 2.0 * qzz * z,
c
c and the constant term sc in the equation of the tangent plane
c at point "a",
c
c sc + sx * x + sy * y + sz * z = 0
c
c where sc = -a dot s.
c
c Any component of "s" less than the estimated error in its
c calculation, based on tol, will be truncated to zero.
c
c Input: ax, ay, az, qc, qx, qy, qz, qxy, qyz, qzx,
c qxx, qyy, qzz, tol.
c
c Output: sx, sy, sz, slen, sc.
c
c Calls: aptvunz, aptvdot
c
c Glossary:
c
c ax,ay,az Input The x, y, z coordinates of point "a".
c
c sc Output Constant term in implicit equation of the plane
c tangent to the quadric surface at point "a".
c If normal vector is zero, sc is returned as -1.e99.
c
c sx,sy,sz Output The x, y, z components of the normal vector "b", the
c gradient of function f(x,y,z), at point "a".
c The coefficients of x, y and z, respectively,
c in the implicit equation of the tangent plane.
c
c nerr Output Flag nerr indicates any input or result error:
c 1 if vector "b" is zero, so no tangent plane exists.
c
c tol Input Numerical tolerance limit.
c On computers with 64-bit floating point words,
c recommend 1.e-5 to 1.e-15.
c
c History:
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptqnor finding the normal vector at point "a"' /
cbug & ' of the family of quadric surfaces with the',
cbug & ' coefficients qc-qzz:' /
cbug & ' ax, ay, az =',1p3e22.14 /
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, 9901) ax, ay, az,
cbug & qc, qx, qy, qz, qxy, qyz, qzx, qxx, qyy, qzz
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
sc = -1.e99
sx = -1.e99
sy = -1.e99
sz = -1.e99
c.... Find the normal vector "b" at point "a".
sx = qx + 2.0 * qxx * ax + qxy * ay + qzx * az
sy = qy + qxy * ax + 2.0 * qyy * ay + qyz * az
sz = qz + qzx * ax + qyz * ay + 2.0 * qzz * az
c.... Truncate the components of "b" to zero if less than the estimated error
c.... in their calculation, based on tol.
sxerr = tol * (abs (qx) + 4.0 * abs (qxx * ax) +
& 2.0 * (abs (qxy * ay) + abs (qyz * az)))
syerr = tol * (abs (qy) + 4.0 * abs (qyy * ay) +
& 2.0 * (abs (qyz * az) + abs (qzx * ax)))
szerr = tol * (abs (qz) + 4.0 * abs (qzz * az) +
& 2.0 * (abs (qzx * ax) + abs (qxy * ay)))
if (abs (sx) .le. sxerr) sx = 0.0
if (abs (sy) .le. syerr) sy = 0.0
if (abs (sz) .le. szerr) sz = 0.0
c.... Find the magnitude of normal vector "b", after truncating any
c.... components to zero if very small relative to the magnitude.
call aptvunz (sx, sy, sz, tol, slen)
if (slen .eq. 0.0) then
nerr = 1
go to 999
endif
sx = slen * sx
sy = slen * sy
sz = slen * sz
c.... Find the equation of the tangent plane.
call aptvdot (ax, ay, az, sx, sy, sz, 1, tol, adotb, nerr)
sc = -adotb
999 continue
cbugc***DEBUG begins.
cbug 9905 format (/ 'aptqnor results: nerr=',i3 /
cbug & ' sc =',1pe22.14 /
cbug & ' sx, sy, sz =',1p3e22.14 )
cbug write ( 3, 9905) nerr, sc, sx, sy, sz
cbugc***DEBUG ends.
return
c.... End of subroutine aptqnor. (+1 line.)
end
UCRL-WEB-209832