subroutine aptlnic (au, av, bu, bv, rc, cu, cv, np, tol,
& du, dv, eu, ev, nint, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTLNIC
c
c call aptlnic (au, av, bu, bv, rc, cu, cv, np, tol,
c & du, dv, eu, ev, nint, nerr)
c
c Version: aptlnic Updated 1993 December 6 15:20..
c aptlnic Originated 1990 March 21 13:40.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To find, for each of np sets of input data, any intersection
c points d = (du, dv) and e = (eu, ev) of the line through
c points a = (au, av) and b = (bu, bv), and the circle centered
c at point c = (cu, cv) with radius rc, all in the uv plane.
c Flag nint indicates the number or type of intersection points.
c Flag nerr indicates any input error.
c
c Input: au, av, bu, bv, rc, cu, cv, np, tol.
c
c Output: du, dv, eu, ev, nint, nerr.
c
c Calls: aptqrtv, aptvadc, aptvdic, aptvdoc
c
c
c Glossary:
c
c au, av Input The u and v coordinates of point "a" on the line "ab"
c in the uv plane. Size np.
c
c bu, bv Input The u and v coordinates of point "b" on the line "ab"
c in the uv plane. Size np.
c
c cu, cv Input The u and v coordinates of point "c" at the center of
c the circle in the uv plane with radius rc. Size np.
c
c du, dv Output The u and v coordinates of point "d" at an intersection
c of the line "ab" and the circle of radius rc centered
c at point "c", if an intersection occurs (nint = 1 or
c nint = 2). Meaningless if nint = -1 or 0. Size np.
c
c eu, ev Output The u and v coordinates of point "e" at an intersection
c of the line "ab" and the circle or radius rc centered
c at point "c", if an intersection occurs (nint = 2).
c Meaningless if nint = -1, 0 or 1. Size np.
c
c nerr Output Indicates an input error, if not 0.
c 1 if np is not positive.
c
c nint Output Indicates the type and number of intersection points:
c -1 if points "a" and "b" coincide, based on tol.
c Points "d" and "e" are meaningless.
c 0 if no intersection occurs.
c Points "d" and "e" are meaningless.
c 1 if the line and circle are tangent at the single
c point "d". Point "e" is meaningless.
c 2 if the line and circle intersect at the two
c distinct points "d" and "e".
c Size np.
c
c np Input Size of arrays.
c
c rc Input The radius of the circle in the uv plane, centered at
c point "c". The absolute value is used. Size np.
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---- Coordinate u of point "a".
dimension au (1)
c---- Coordinate v of point "a".
dimension av (1)
c---- Coordinate u of point "b".
dimension bu (1)
c---- Coordinate v of point "b".
dimension bv (1)
c---- Coordinate u of point "c".
dimension cu (1)
c---- Coordinate v of point "c".
dimension cv (1)
c---- Coordinate u of point "d".
dimension du (1)
c---- Coordinate v of point "d".
dimension dv (1)
c---- Coordinate u of point "e".
dimension eu (1)
c---- Coordinate v of point "e".
dimension ev (1)
c---- Type and number of intersections.
dimension nint (1)
c---- Radius of circle at "c".
dimension rc (1)
c.... Local variables.
c---- Coefficient of fad**2 and fae**2.
common /laptlnic/ aa (64)
c---- Component u of vector "ab".
common /laptlnic/ abu (64)
c---- Component v of vector "ab".
common /laptlnic/ abv (64)
c---- Value of dac - rc.
common /laptlnic/ acmr
c---- Value of dac + rc.
common /laptlnic/ acpr
c---- Component u of vector "ac".
common /laptlnic/ acu (64)
c---- Component v of vector "ac".
common /laptlnic/ acv (64)
c---- Coefficient of fad and fae.
common /laptlnic/ bb (64)
c---- Constant term in quadratic.
common /laptlnic/ cc (64)
c---- Distance from "a" to "b".
common /laptlnic/ dab (64)
c---- Distance from "a" to "c".
common /laptlnic/ dac (64)
c---- Fractional distance along "ab" of "d".
common /laptlnic/ fad (64)
c---- Fractional distance along "ab" of "e".
common /laptlnic/ fae (64)
c---- 1 if qq truncated to zero.
common /laptlnic/ itrun (64)
c---- Index in arrays.
common /laptlnic/ n
c---- First index of subset of data.
common /laptlnic/ n1
c---- Last index of subset of data.
common /laptlnic/ n2
c---- Index in external array.
common /laptlnic/ nn
c---- Size of current subset of data.
common /laptlnic/ ns
c---- Value of bb**2 - 4.0 * aa * cc.
common /laptlnic/ qq (64)
c---- Dot product of vectors "ab" and "ac".
common /laptlnic/ vdot (64)
c---- Magnitude of a vector.
common /laptlnic/ vlen (64)
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptlnic finding intersection of line ab and circle',
cbug & ' c, r.' / ' tol= ',1pe22.14)
cbug 9902 format (i3,' au,av= ',1p2e22.14 /
cbug & ' bu,bv= ',1p2e22.14 /
cbug & ' cu,cv,rc=',1p3e22.14)
cbug write ( 3, 9901) tol
cbug write ( 3, 9902) (n, au(n), av(n), bu(n), bv(n),
cbug & cu(n), cv(n), rc(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 210
endif
c.... Set up the indices of the first subset of data.
n1 = 1
n2 = min (np, 64)
c.... Loop over subsets of data.
110 ns = n2 - n1 + 1
c.... Find the vector "ab" from point "a" to point "b".
c.... The magnitude dab may be truncated to zero.
call aptvdic (au(n1), av(n1), bu(n1), bv(n1), ns, tol,
& abu, abv, dab, nerr)
c.... Find the vector "ac" from point "a" to point "c".
c.... The magnitude dab may be truncated to zero.
call aptvdic (au(n1), av(n1), cu(n1), cv(n1), ns, tol,
& acu, acv, dac, nerr)
c.... Find the dot product of vectors "ab" and "ac".
call aptvdoc (abu, abv, acu, acv, ns, tol,
& vdot, nerr)
c.... Find the coefficients of the quadratic equation for fad, fae.
c---- Loop over subset of data.
do 120 n = 1, ns
nn = n + n1 - 1
aa(n) = dab(n)**2
bb(n) = -2.0 * vdot(n)
acpr = dac(n) + abs (rc(nn))
acmr = dac(n) - abs (rc(nn))
if (abs (acmr) .lt. tol * acpr) then
acmr = 0.0
endif
cc(n) = acpr * acmr
c---- End of loop over subset of data.
120 continue
c.... Find the fractional distances fad and fae.
call aptqrtv (0, aa, bb, cc, qq, ns, tol,
& nint, fad, fae, itrun, nerr)
c.... Find the intersection points "d" and "e".
call aptvadc (au(n1), av(n1), 1.0, fad, abu, abv, ns, tol,
& du(n1), dv(n1), vlen, nerr)
call aptvadc (au(n1), av(n1), 1.0, fae, abu, abv, ns, tol,
& eu(n1), ev(n1), vlen, nerr)
c.... See if all data subsets are done.
c---- Do another subset of data.
if (n2 .lt. np) then
n1 = n2 + 1
n2 = min (np, n1 + 63)
go to 110
endif
cbugc***DEBUG begins.
cbug 9903 format (/ 'aptlnic results:' /
cbug & (i3,' nint=',i2 /
cbug & ' du,dv=',1p2e22.14 /
cbug & ' eu,ev=',1p2e22.14))
cbug write ( 3, 9903) (n, nint(n),
cbug & du(n), dv(n), eu(n), ev(n), n = 1, np)
cbugc***DEBUG ends.
210 return
c.... End of subroutine aptlnic. (+1 line.)
end
UCRL-WEB-209832