subroutine aptrext (a, b, c, tol,
& next, fa, aa, fb, bb, fc, cc, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTREXT
c
c call aptrext (a, b, c, tol, next, fa, aa, fb, bb, fc, cc, nerr)
c
c Version: aptrext Updated 1999 September 13 14:00.
c aptrext Originated 1999 August 30 13:00.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: For the triangle with edges a, b and c, to find the number next
c of extended triangles formed by extending sides a and b
c (to lengths aa and bb), and connecting their new endpoints
c (with new edge cc), that add a new section with the same
c perimeter and area as the original triangle, and to find the
c ratios of side lengths fa = a / aa, fb = b / bb and fc = c / cc.
c Flag nerr indicates any input error.
c
c Input: a, b, c, tol.
c
c Output: next, fa, aa, fb, bb, fc, cc, nerr.
c
c Calls: aptqrts
c
c Glossary:
c
c a Input The length of edge a of the original triangle.
c
c b Input The length of edge b of the original triangle.
c
c c Input The length of edge c of the original triangle.
c
c aa Output The length of the extension of side a. Size 2.
c
c bb Output The length of the extension of side b. Size 2.
c
c cc Output The length of the new edge joining the ends of
c edges aa and bb.
c
c fa Output The ratio a / aa. Size 2.
c
c fb Output The ratio b / bb. Size 2.
c
c fc Output The ratio c / cc. Size 2.
c
c next Output The number of extended triangles (0, 1 or 2).
c
c nerr Output Indicates an input error, if not 0.
c 1 if an edge length is too small or negative.
c 2 if the initial triangle is impossible (an edge
c length is too large (not less than the sum of the
c other two edge lengths).
c
c tol Input Numerical tolerance limit.
c On Cray computers, recommend 1.e-5 to 1.e-11.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
implicit none
c.... Arguments.
real a ! Length of edge a of original triangle.
real b ! Length of edge b of original triangle.
real c ! Length of edge c of original triangle.
real aa(2) ! Length of extended edge aa = a / fa.
real bb(2) ! Length of extended edge bb = b / fb.
real cc(2) ! Length of new edge cc = c / fc.
real fa(2) ! Ratio a / aa.
real fb(2) ! Ratio b / bb.
real fc(2) ! Ratio c / cc.
integer next ! Number of extended triangles.
integer nerr ! Error indicator.
real tol ! Numerical tolerence limit.
c.... Local variables.
real aax ! New edge length aa.
real bbx ! New edge length bb.
real ccx ! New edge length cc.
real c2 ! Coefficient of t**2 in quadratic eq.
real c1 ! Coefficient of t in quadratic eq.
real c0 ! Coefficient of 1 in quadratic eq.
real fax ! Ratio a / aa.
real fbx ! Ratio b / bb.
real fcx ! Ratio c / cc.
integer itrun ! Truncation flag from aptqrts.
integer n ! Index of cut line.
integer nerra ! Error flag from apt subroutine.
integer nroots ! # of real roots of quadratic equation.
real per ! Perimeter of extended triangle.
real qq ! c1**2 - 4*c2*c1.
real root1 ! First real root of quadatic equation.
real root2 ! Second real root of quadratic equation.
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptrext finding triangle extensions',
cbug & ' tol=',1pe22.14 /
cbug & ' a ,b ,c = ',1p3e22.14 )
cbug write ( 3, 9901) tol, a, b, c
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
next = 0
do 110 n = 1, 2
fa(n) = -1.e99
fb(n) = -1.e99
fc(n) = -1.e99
aa(n) = -1.e99
bb(n) = -1.e99
cc(n) = -1.e99
110 continue
c=======================================================================********
c.... Test for input errors.
call aptrich (a, b, c, tol, nerr)
if (nerr .ne. 0) go to 410
c.... Find any extended triangles.
per = 2.0 * (a + b) ! Perimeter of extended triangle.
ccx = c**2 / per + 0.25 * per
fcx = c / ccx
c2 = 1.0
c1 = ccx - per
c0 = 2.0 * a * b
cbugcbugc***DEBUG begins.
cbugcbug 9804 format ('DEBUG 4 a,b,c=',1p3e22.14 /
cbugcbug & 'per,fcx,ccx =',1p3e22.14 /
cbugcbug & 'c0, c1, c2 =',1p3e22.14 )
cbugcbug write ( 3, 9804) a, b, c, per, fcx, ccx, c0, c1, c2
cbugcbugc***DEBUG ends.
call aptqrts (0, c2, c1, c0, qq, tol,
& nroots, root1, root2, itrun)
if (nroots .le. 0) then
go to 410
elseif (nroots .eq. 1) then
aax = root1
bbx = root1
fax = a / aax
fbx = b / bbx
fcx = c / ccx
cbugcbugc***DEBUG begins.
cbugcbug 9701 format ('DEBUG 1 a,b,c=',1p3e22.14 /
cbugcbug & 'fax,fbx,fcx= ',1p3e22.14 /
cbugcbug & 'aax,bbx,ccx= ',1p3e22.14 )
cbugcbug write ( 3, 9701) a, b, c, fax, fbx, fcx, aax, bbx, ccx
cbugcbugc***DEBUG ends.
if ((fax .lt. (0.5 - tol)) .or.
& (fax .gt. (1.0 + tol))) go to 410
if ((fbx .lt. (0.5 - tol)) .or.
& (fbx .gt. (1.0 + tol))) go to 410
if ((aax .ge. (bbx + ccx) * (1.0 - tol)) .or.
& (bbx .ge. (ccx + aax) * (1.0 - tol)) .or.
& (ccx .ge. (aax + bbx) * (1.0 - tol))) then
go to 410
endif
next = next + 1
fa(next) = fax
fb(next) = fbx
fc(next) = fcx
aa(next) = aax
bb(next) = bbx
cc(next) = ccx
elseif (nroots .eq. 2) then
aax = root1
bbx = root2
fax = a / aax
fbx = b / bbx
fcx = c / ccx
cbugcbugc***DEBUG begins.
cbugcbug 9702 format ('DEBUG 2 a,b,c=',1p3e22.14 /
cbugcbug & 'fax,fbx,fcx= ',1p3e22.14 /
cbugcbug & 'aax,bbx,ccx= ',1p3e22.14 )
cbugcbug write ( 3, 9702) a, b, c, fax, fbx, fcx, aax, bbx, ccx
cbugcbugc***DEBUG ends.
if ((fax .lt. (0.5 - tol)) .or.
& (fax .gt. (1.0 + tol))) go to 350
if ((fbx .lt. (0.5 - tol)) .or.
& (fbx .gt. (1.0 + tol))) go to 350
if ((aax .ge. (bbx + ccx) * (1.0 - tol)) .or.
& (bbx .ge. (ccx + aax) * (1.0 - tol)) .or.
& (ccx .ge. (aax + bbx) * (1.0 - tol))) then
go to 350
endif
next = next + 1
fa(next) = fax
fb(next) = fbx
fc(next) = fcx
aa(next) = aax
bb(next) = bbx
cc(next) = ccx
350 aax = root2
bbx = root1
fax = a / aax
fbx = b / bbx
fcx = c / ccx
cbugcbugc***DEBUG begins.
cbugcbug 9703 format ('DEBUG 3 a,b,c=',1p3e22.14 /
cbugcbug & 'fax,fbx,fcx= ',1p3e22.14 /
cbugcbug & 'aax,bbx,ccx= ',1p3e22.14 )
cbugcbug write ( 3, 9703) a, b, c, fax, fbx, fcx, aax, bbx, ccx
cbugcbugc***DEBUG ends.
if ((fax .lt. (0.5 - tol)) .or.
& (fax .gt. (1.0 + tol))) go to 410
if ((fbx .lt. (0.5 - tol)) .or.
& (fbx .gt. (1.0 + tol))) go to 410
if ((aax .ge. (bbx + ccx) * (1.0 - tol)) .or.
& (bbx .ge. (ccx + aax) * (1.0 - tol)) .or.
& (ccx .ge. (aax + bbx) * (1.0 - tol))) then
go to 410
endif
next = next + 1
fa(next) = fax
fb(next) = fbx
fc(next) = fcx
aa(next) = aax
bb(next) = bbx
cc(next) = ccx
endif ! Cut endpoints are real.
c=======================================================================********
410 continue
cbugc***DEBUG begins.
cbug 9916 format (/ 'aptrext results. nerr = ',i2,' next = ',i1 )
cbug 9918 format ('next =',i2 /
cbug & 'fa, fb, fc = ',1p3e22.12 /
cbug & 'aa, bb, cc = ',1p3e22.12 )
cbug write ( 3, 9916) nerr, next
cbug
cbug if (next .gt. 0) then
cbug do 420 n = 1, next
cbug write ( 3, 9918) n, fa(n), fb(n), fc(n), aa(n), bb(n), cc(n)
cbug 420 continue
cbug endif
cbugc***DEBUG ends.
return
c.... End of subroutine aptrext. (+1 line.)
end
UCRL-WEB-209832