subroutine aptqexv (dx, dy, dz, qc, qx, qy, qz, qxy, qyz, qzx,
     &                    qxx, qyy, qzz, tol,
     &                    nextr, ex, ey, ez, adir, nerr)

ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c                             SUBROUTINE APTQEXV
c
c     call aptqexv (dx, dy, dz, qc, qx, qy, qz, qxy, qyz, qzx,
c    &              qxx, qyy, qzz, tol,
c    &              nextr, ex, ey, ez, adir, nerr)
c
c     Version:  aptqexv  Updated    1994 October 4 16:50.
c               aptqexv  Originated 1994 October 4 16:50.
c
c     Author:   Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c     Purpose:  For the quadric surface specified by the coordinates qc to qzz,
c               to find any extremities e = (ex, ey, ez) in the direction of the
c               vector d = (dx, dy, dz), and the type adir of extremities.
c               The possible types of extrema are minima and maxima (points,
c               lines or planes), saddle points, constant planes, intersection
c               lines, and isolated points.
c 
c               The implicit equation of the quadric curve is:
c
c                 F(x,y,z) = qc + qx  * x     + qy  * y     +
c                                 qxy * x * y + qyz * y * z + qzx * z * x +
c                                 qxx * x**2  + qyy * y**2  + qzz * z**2  = 0  .
c
c               Results may be unpredictable if the coefficients of the quadric
c               surface differ by many orders of magnitude.
c
c     Method:   See APT subroutine aptextr.
c
c     Input:    dx, dy, dz, qc, qx, qy, qz, qxy, qyz, qzx, qxx, qyy, qzz, tol.
c
c     Output:   nextr, ex, ey, ez, adir, nerr.
c
c     Calls: aptmopv, aptqext, aptrois, aptrotv 
c               
c
c     Glossary:
c
c     adir      Output   An ASCII 8-character word, the type of extreme point:
c                          "Minimum", "Maximum", "Plane", "IntPlane",
c                          "Saddle", "Point".
c                          If none, return "unknown".  Size = 2.
c
c     dx,dy,dz  Input    Components x, y and z of the direction in which any
c                          extreme points are to be found.
c
c     ex,ey,ez  Output   The coordinates of any extreme points.  Size = 2.
c
c     nerr      Output   Indicates an input error, if not zero:
c                          1 if all coefficients, except possible qc, are zero,
c                          or the quadric surface type is 17, 19, 20, 23 or 24.
c                          2 if the magnitude of vector "d" is too small.
c
c     nextr     Output   Number of extreme or intersection points found.
c
c     qc        Input    Constant term in the equation of the quadric curve.
c
c     qx,qy,qz  Input    Coefficients of x, y and z, respectively, in the
c                          equation of the quadric curve.
c
c     qxx       Input    Coefficient of x**2 in the equation of the quadric
c                          curve.
c
c     qxy       Input    Coefficient of x * y in the equation of the quadric
c                          curve.
c
c     qyy       Input    Coefficient of y**2 in the equation of the quadric
c                          curve.
c
c     qyz       Input    Coefficient of y * z in the equation of the quadric
c                          curve.
c
c     qzx       Input    Coefficient of z * x in the equation of the quadric
c                          curve.
c
c     qzz       Input    Coefficient of z**2 in the equation of the quadric
c                          curve.
c
c     tol       Input    Numerical tolerance limit.
c                          On Cray computers, recommend 1.e-11.
c                          On 32-bit computers, recommend 1.e-6.
c
c     History:
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.

c.... Dimensioned arguments.

      dimension adir(2)               ! Type of extreme point.
      character*8 adir                ! Type of extreme point.
      dimension ex(2)                 ! Coordinate x of an extreme point.
      dimension ey(2)                 ! Coordinate y of an extreme point.
      dimension ez(2)                 ! Coordinate z of an extreme point.

c.... Local variables.

      common /laptqexv/ abdir(6)      ! Type of extreme point.
      character*8 abdir               ! Type of extreme point.
      common /laptqexv/ rotm(3,3)     ! Rotation matrix, vector "d" to z axis.
      common /laptqexv/ bx(6)         ! Coordinate x of an extreme point.
      common /laptqexv/ by(6)         ! Coordinate y of an extreme point.
      common /laptqexv/ bz(6)         ! Coordinate z of an extreme point.
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptqexv finding extrema of quadric in "d" dir.' /
cbug     &  '  tol=',1pe13.5 /
cbug     &  '  dx,dy,dz=   ',1p3e22.14 )
cbug 9902 format (
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) tol, dx, dy, dz
cbug      write ( 3, 9902) qc, qx, qy, qz, qxy, qyz, qzx, qxx, qyy, qzz
cbugc***DEBUG ends.

c.... Initialize.

      nerr    = 0

      do 110 n = 1, 2
        ex(n)   = -1.e99
        ey(n)   = -1.e99
        ez(n)   = -1.e99
        adir(n) = 'unknown'
  110 continue

c.... Find any extrema in the quadric in the direction (dx, dy, dz).

      nextr  = 0

c.... Find the rotation operator to rotate the axis vector (dx, dy, dz) to the
c....   z axis.

      call aptrotv (dx, dy, dz, 0., 0., 1., tol, rotm, nerra)
cbugc***DEBUG begins.
cbug 9990 format ('  op11,op12,op13  ',1p3e20.12 /
cbug     &        '  op21,op22,op23  ',1p3e20.12 /
cbug     &        '  op31,op32,op33  ',1p3e20.12)
cbug      write ( 3, 9990) ((rotm(i,j), j = 1, 3), i = 1, 3)
cbugc***DEBUG ends.
      if (nerra .ne. 0) then
        nerr = 2
        go to 999
      endif

c.... Rotate the quadric surface.

      ac  = qc
      ax  = qx
      ay  = qy
      az  = qz
      axy = qxy
      ayz = qyz
      azx = qzx
      axx = qxx
      ayy = qyy
      azz = qzz
cbugc***DEBUG begins.
cbug 9960 format ('    ac            ',1pe20.12 )
cbug 9980 format ('    ax,  ay,  az  ',1p3e20.12 )
cbug 9985 format ('    axy, ayz, azx ',1p3e20.12 )
cbug 9991 format ('    axx, ayy, azz ',1p3e20.12 )
cbug      write ( 3, 9960) ac
cbug      write ( 3, 9980) ax, ay, az
cbug      write ( 3, 9985) axy, ayz, azx
cbug      write ( 3, 9991) axx, ayy, azz
cbugc***DEBUG ends.

      call aptrois (rotm, 0, 0., 0., 0., ac, ax, ay, az,
     &              axy, ayz, azx, axx, ayy, azz, tol, nerra)
cbugc***DEBUG begins.
cbug      write ( 3, 9960) ac
cbug      write ( 3, 9980) ax, ay, az
cbug      write ( 3, 9985) axy, ayz, azx
cbug      write ( 3, 9991) axx, ayy, azz
cbugc***DEBUG ends.

c.... Find the extrema in the z direction.

      call aptqext (ac, ax, ay, az, axy, ayz, azx, axx, ayy, azz,
     &              tol, ntyps, cenx, ceny, cenz, nmm, bx, by, bz,
     &              abdir, nerra)
cbugc***DEBUG begins.
cbug 9922 format ('nmm=',i2)
cbug 9923 format ('n=',i2,' abdir="',a8,'"')
cbug      write ( 3, 9922) nmm
cbug      if (nmm .gt. 0) then
cbug        write ( 3, 9923) (n, abdir(n), n = 1, nmm)
cbug      endif
cbugc***DEBUG ends.
      if (nerra .ne. 0) then
        nerr = 1
        go to 999
      endif

c.... Rotate any extrema back to the original direction.

      if (nmm .gt. 0) then

        call aptmopv (rotm, 1, 0., 0., 0., bx, by, bz, nmm, tol, nerra)

      endif

c.... Find any 'Min z' or 'Max z', relabel 'Minimum' or 'Maximum'.
c.... Find any 'Const z' or 'Inter z', relabel 'Plane' or 'IntPlane'.
c.... Find any 'Saddle z' or 'Point z', relabel 'Saddle' or 'Point'.

      do 150 n = 1, nmm               ! Loop over extreme points.

        if     (abdir(n) .eq. 'Min z') then
          nextr       = nextr + 1
          adir(nextr) = 'Minimum'
          ex(nextr)   = bx(n)
          ey(nextr)   = by(n)
          ez(nextr)   = bz(n)
        elseif (abdir(n) .eq. 'Max z') then
          nextr       = nextr + 1
          adir(nextr) = 'Maximum'
          ex(nextr)   = bx(n)
          ey(nextr)   = by(n)
          ez(nextr)   = bz(n)
        elseif (abdir(n) .eq. 'Const z') then
          nextr       = nextr + 1
          adir(nextr) = 'Plane'
          ex(nextr)   = bx(n)
          ey(nextr)   = by(n)
          ez(nextr)   = bz(n)
        elseif (abdir(n) .eq. 'Inter z') then
          nextr       = nextr + 1
          adir(nextr) = 'IntPlane'
          ex(nextr)   = bx(n)
          ey(nextr)   = by(n)
          ez(nextr)   = bz(n)
        elseif (abdir(n) .eq. 'Saddle z') then
          nextr       = nextr + 1
          adir(nextr) = 'Saddle'
          ex(nextr)   = bx(n)
          ey(nextr)   = by(n)
          ez(nextr)   = bz(n)
        elseif (abdir(n) .eq. 'Point z') then
          nextr       = nextr + 1
          adir(nextr) = 'Point'
          ex(nextr)   = bx(n)
          ey(nextr)   = by(n)
          ez(nextr)   = bz(n)
        endif

  150 continue                        ! End of loop over extreme points.

  999 continue
cbugc***DEBUG begins.
cbug 9910 format (/ 'aptqexv results:  nextr=',i2,' nerr=',i2)
cbug 9912 format (a8,'  xyz=',1p3e22.14)
cbug      write ( 3, 9910) nextr, nerr
cbug      do 991 n = 1, 2
cbug        write ( 3, 9912) adir(n), ex(n), ey(n), ez(n)
cbug  991 continue
cbugc***DEBUG ends.

      return

c.... End of subroutine aptqexv.      (+1 line)
      end

UCRL-WEB-209832