      SUBROUTINE  AT (N)
C     * * * * * * * * * * * * * * *
C      A DEBUGGING AID
C     * * * * * * * * * * * * * * *
      WRITE (6,10) N
      WRITE (12,10) N
 10   FORMAT('----> HERE AT ',I8)
      RETURN
      END
      SUBROUTINE  BANWTV (NE, N, NTOLPT, NODES, LNODE,
     1                    LBAND, IOPT, NNOW)
C      * * * * * * * * * * * * * * * * * * * * * * * * *
C      DETERMINE SYSTEM BANDWIDTH WITH NODES IN VECTOR
C          (USE NEW NODE NUMBERS IF IOPT .GT. 1)
C      * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  NTOLPT, NODES, LNODE, NNOW, LBAND
      DIMENSION  NODES(1), LNODE(N), NTOLPT(NE+1),
     1           LBAND(NE), NNOW(IOPT)
C     LBAND     = BANDWIDTH OF ELEMENT
C     LNODE     = INCIDENCES ARRAY OF A SINGLE ELEMENT
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = NUMBER OF ELEMENTS 
C     NNOW(I)   = NEW NUMBER OF OLD NODE I
C     NODES     = SYSTEM VECTOR OF ELEMENT INCIDENCES
C     NTOLPT(I) = PT WHERE INCIDENCES FOR ELEM I BEGIN
C      LOOP OVER ELEMENTS
      DO 30 IE = 1,NE
C        EXTRACT INCIDENCES LIST
        L1 = NTOLPT(IE)
        LN = NTOLPT(IE + 1) - L1
        CALL  LNODEV (LN, NODES(L1), LNODE)
C        CONVERT TO NEW NODE NUMBERS?
        IF ( IOPT .NE. 1 )  THEN
          DO 10 IN = 1,LN
            NOLD = LNODE(IN)
            IF ( NOLD .LT. 0 ) STOP 'ERROR IN BANWTV'
   10     LNODE(IN) = NNOW(NOLD)
        ENDIF
C        FIND ELEMENT BANDWIDTH
   30 CALL  ELBAND (LN, 1, LBAND(IE), LNODE)
      RETURN
      END
      SUBROUTINE  CONTRL (M, NE, N, NEED, NUSER, NTRY, LNOPT, 
     1                    MESH, LISTLB, LISTLW, LISTLL, LNEWN, 
     2                    LNEWL, ISAVE, IPICK, INHEAD, L2L)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * *
C                READ TITLE AND CONTROL DATA
C     * * * * * * * * * * * * * * * * * * * * * * * * * * *
      DIMENSION  TITLE(15)
C     INHEAD    = NUMBER OF HEADERS BEFORE CONNECTIVITY
C     IPICK     = 0 RMS VALUE GOVERNS, = 1 MAX GOVERNS
C     ISAVE     = RESTART CONTROL
C     LNOPT     = 0 OUTPUT ELEMENTS, =1 OUTPUT NODES
C     M         = NUMBER OF NODES
C     MESH      = MESH ECHO FLAG, 0 = NONE
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = NUMBER OF ELEMENTS
C     NEED      = NUMBER OF COMMON NODES TO BE NEIGHBOR
C     NTRY    = NUMBER OF PROG SELECTED STARTING ELEMENTS
C     NUSER   = NUMBER OF USER SELECTED STARTING ELEMENTS
      READ  (5,5010) TITLE
 5010 FORMAT ( 15A4 )
      WRITE (6,5010) TITLE
C
C....     ** READ AND ECHO CONTROLS **
C
      READ  (5,5030)  M, NE, N, NEED, NUSER, NTRY, LNOPT, MESH, 
     1                ISAVE, IPICK, INHEAD
      READ  (5,5030)  LISTLB, LISTLW, LISTLL, LNEWN, LNEWL, L2L
 5030 FORMAT ( 16I5 )
      IF ( N      .LT. 1 )  N      = 36
      IF ( NEED   .LT. 1 )  NEED   = 1
      IF ( NEED   .GT. N )  NEED   = 1
      IF ( NTRY   .EQ. 0 )  NTRY   = 5 + NE/1000
      IF ( NTRY   .LT. 0 )  NTRY   = 0
      IF ( ISAVE  .EQ. 2 )  NTRY   = 0
      WRITE (6,5040)  M, NE, N, NEED, NUSER, NTRY, LNOPT, MESH, 
     1                ISAVE, IPICK, INHEAD
 5040 FORMAT ( /,'CONTROL PARAMETERS:            (DEFAULT)', /,
     1 'MAXIMUM NODE NUMBER .....................',I7,/,
     2 'NUMBER OF ELEMENTS ......................',I7,/,
     3 'MAX NUMBER OF NODES PER ELEMENT...(36)...',I7,/,
     4 'COMMON NODES TO BE NEIGHBOR .......(1)...',I7,/,
     5 'NUMBER OF USER STARTING ELEMENTS ..(0)...',I7,/,
     6 'NUMBER OF SYSTEM SELECTED STARTS ..(5)...',I7,/,
     7 'REORDER ELEMENTS=0, OR NODES=1 ....(0)...',I7,/,
     8 'MESH ECHO FLAG, 0=NONE ............(0)...',I7,/,
     9 'RESTART, 0=NO, 1=SAVE, 2=READ......(0)...',I7,/,
     1 'SELECTION METHOD, 0=RMS, 1=MAX.....(0)...',I7,/,
     2 'NUMBER OF HEADERS BEFORE NODES.....(0)...',I7,/)
C          * CITE DATA TO BE PRINTED *
      WRITE (6,5045)  LISTLB, LISTLW, LISTLL, LNEWN, LNEWL
 5045 FORMAT ( 'PRINT FLAGS, 0=OFF 1=ON',/,
     1 'ELEMENT BANDWIDTHS ................(0)...',I7,/,
     2 'ELEMENT WAVEFRONTS ................(0)...',I7,/,
     3 'ELEMENT ADJACENT TO ELEMENTS ......(0)...',I7,/,
     4 'NEW NODE ORDER LIST ...............(0)...',I7,/,
     5 'NEW ELEMENT ORDER LIST ............(0)...',I7,/)
      WRITE (6,5050) L2L
 5050 FORMAT (
     1 'STORE NEIGHBORS -1=OLD 0=NO 1=NEW .(0)...',I7,/)
      RETURN
      END
      SUBROUTINE  DRIVER (IFIXED, ID)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * *
C       ELEMENT AND NODAL RESEQUENCING FOR FRONTAL
C               AND BANDWIDTH SOLUTIONS
C         Copyright J. E. Akin, 1981, Revised 1991
C     * * * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  ID
      DIMENSION  ID(IFIXED)
C      ----- NOTATION ------
C     IFIXED    = AVAILABLE STORAGE
C     INHEAD    = NUMBER OF ELEMENT HEADERS IN LHEAD
C     ISAVE     = RESTART FLAG
C               0 - DIRECT INPUT ONLY
C               1 - BINARY SAVE FOR RESTART
C               2 - RESTART, NTRY=0, READ NUSER STARTS
C     L1OFN(I)  = ELEM OF 1-ST APPEARANCE OF NODE I
C     LADJL     = PACKED LIST OF ELEMENTS ADJACENT TO ELEMS
C     LBAND(I)  = BANDWIDTH OF ELEMENT I
C     LFIRST(I) = ELEMENT OF FIRST APPEARENCE OF NODE I
C     LHEAD     = APPLICATION DEPENDENT INTEGERS BEFORE NODES
C     LHERE     = STAGE OF LAST BINARY RESTART SAVE
C     LLAST(I)  = ELEMENT OF LAST APPEARENCE OF NODE I
C     LLMAX     = DIMENSION OF ARRAY LADJL
C     LNODE     = INCIDENCES ARRAY OF A SINGLE ELEMENT
C     LNOW(I)   = NEW ELEMENT NUMBER OF OLD ELEMENT I
C     LSTART    = SPECIFIED STARTING ELEMENT.
C     LTOLPT(I) = PT WHERE ELEM NEIGHBORS OF ELEM I BEGIN
C     LUN2L(I)  = NO. OF UNNUMBERED ELEM ADJ TO ELEMENT I
C     LUN2N(I)  = NO. OF UNNUMBERED ELEM ADJ TO NODE I
C     LWAS(I)   = OLD NUMBER OF NEW ELEMEMT I
C     LWAVE(I)  = WAVEFRONT AT ELEMENT I
C     M         = NUMBER OF NODES
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = NUMBER OF ELEMENTS
C     NNMAX     = DIMENSION OF ARRAY NODES
C     NNOW(I)   = NEW NUMBER OF OLD NODE I
C     NOADJL(I) = NUMBER OF ELEMENTS ADJACENT TO NODE I
C     NODES     = SYSTEM VECTOR OF ELEMENT INCIDENCES
C     NTOLPT(I) = PT WHERE INCIDENCES FOR ELEM I BEGIN
C     NWAS(I)   = OLD NODE NUMBER OF NEW NODE I
C      ---  I/O UNITS USED ----
      OPEN (UNIT=5,  FILE='CTL',      STATUS='OLD')
      OPEN (UNIT=6,  FILE='OUT',      STATUS='UNKNOWN')
      OPEN (UNIT=7,  FILE='TOPO',     STATUS='OLD')
      OPEN (UNIT=8,  FILE='ECHO',     STATUS='UNKNOWN')
      OPEN (UNIT=9,  FILE='ABSTRACT', STATUS='UNKNOWN')
      OPEN (UNIT=10, FILE='ORDER',    STATUS='UNKNOWN',
     1               FORM='UNFORMATTED')
      OPEN (UNIT=11, FILE='NEW',      STATUS='UNKNOWN')
      OPEN (UNIT=12, FILE='DEBUG',    STATUS='UNKNOWN')
      OPEN (UNIT=13, FILE='RESTART',  STATUS='UNKNOWN',
     1               FORM='UNFORMATTED')
      OPEN (UNIT=14, FILE='L2L',      STATUS='UNKNOWN')
C      GET CONTROL DATA
      CALL  CONTRL (M, NE, N, NEED, NUSER, NTRY, LNOPT, 
     1              MESH, LISTLB, LISTLW, LISTLL, LNEWN, 
     2              LNEWL, ISAVE, IPICK, INHEAD, L2L)
C..       RESTART CONTROLS
      ID(1) = IFIXED
      ID(2) = 0
      IF ( ISAVE .EQ. 2 )  THEN
        READ (13) LFINAL, LHERE
        IF ( LFINAL .GT. IFIXED ) THEN
          STOP 'IFIXED TOO SMALL FOR RESTART'
        ELSE
          ID(1) = LFINAL
          ID(2) = LHERE
        ENDIF
        READ (13) ( ID(K), K = 3,LFINAL )
        MAXTRY = ID(11)
      ELSE
        MAXTRY = NUSER + NTRY
      ENDIF
C
C....         ** CALCULATE ARRAY POINTERS **
C
C      SAVE 20 LOCATIONS FOR RESTART INFO
C      1-LFINAL 2-LHERE 3-NNMAX 4-ISTART 5-LKEEP
C      6-MKEEP  7-RKEEP 8-LLMAX 9-AVELL 10-L18
C      11-MAXTRY
C
C     ORDER OF ARRAY STORAGE IN ARRAY ID, L(I):
C      FIXED NODAL ARRAYS; i=1-LFIRST(M), 2-LLAST(M),
C       3-L1OFN(M)=LFIRST, 4-LUN2N(M)=LLAST,
C       5-NOADJL(M), 6-NNOW(M), 7-NWAS(M)=NOADJL
      L1 = 21
      L2 = L1 + M
      L3 = L1
      L4 = L2
      L5 = L2 + M
      L6 = L5 + M
      L7 = L5
      L8 = L6 + M
C      FIXED ELEMENT ARRAYS; 8-LNODE(N),
C       9-LWAVE(NE), 10-LBAND(NE)=LWAVE,
C       11-LNOW(NE)=LWAVE, 12-LWAS(NE),
c       11-LNOW(NE), 12-LWAS(NE),
C       13-LUN2L(NE)=LNOW, 14-LTOLPT(NE+1),
C       15-NTOLPT(NE+1)
c       15a-LHEAD(NE)=LWAVE, 15b-LTAIL(NE)=LNOW 
      L9  = L8 +  N
      L10 = L9
      L11 = L9
c     L11 = L9 +  NE
      L12 = L9  + NE
c     L12 = L11 + NE
      L13 = L12 + NE
c     L13 = L11
c     L14 = L12 + NE
      L14 = L13 + NE
      L15 = L14 + NE + 1
      L15a= L9
      L15b= L11
      L16 = L15 + NE + 1
C      FIXED MISCELLANEOUS ARRAYS
C       16-LSTART(NUSER+NTRY), 17-LHEAD
c     L17 = L16 + NTRY + NUSER
      L17 = L16 + MAXTRY
      IF ( INHEAD .GT. 0 )  THEN
        L18 = L17 + NE*INHEAD
      ELSE
        L18 = L17
      ENDIF
C      VARIABLE SIZES TO BE FOUND LATER; 18-NODES(NNMAX),
C       19-LADJL(LLMAX).
      IFREE = IFIXED - L18
      IF ( IFREE .LT. NE ) STOP 'NEED MORE SPACE IN DRIVER'
C....     CALL THE REAL MAIN PROGRAM
C     CALL  REDUCE (M, NE, N, NEED, NUSER, NTRY, LNOPT,
C    1              LISTLB, LISTLW, LISTLL, LNEWN, LNEWL,
C    2              LFIRST, LLAST, L1OFN, LUN2N, NOADJL,
C    3              NNOW, NWAS, LNODE, LWAVE, LBAND,
C    4              LNOW, LWAS, LUN2L, LTOLPT, NTOLPT,
C    5              LSTART, LHEAD, NODES, ID, L18, IFREE, 
C    6              IFIXED, MESH, ISAVE, LHERE, IPICK, INHEAD, L2L)
      CALL  REDUCE (M, NE, N, NEED, NUSER, NTRY, LNOPT,
     1              LISTLB, LISTLW, LISTLL, LNEWN, LNEWL,
     2              ID(L1), ID(L2), ID(L3), ID(L4), ID(L5),
     3              ID(L6), ID(L7), ID(L8), ID(L9), ID(L10),
     4              ID(L11), ID(L12), ID(L13), ID(L14), ID(L15),
     5              ID(L16), ID(L17), ID(L18), ID, L18, IFREE, 
     6              IFIXED, MESH, ISAVE, LHERE, IPICK, INHEAD, L2L)
      WRITE (6,*) 'NORMAL EXIT OF REDUCE'
      CLOSE (UNIT=5)
      CLOSE (UNIT=6)
      CLOSE (UNIT=7)
      CLOSE (UNIT=8)
      CLOSE (UNIT=9)
      CLOSE (UNIT=10)
      CLOSE (UNIT=11)
      CLOSE (UNIT=12)
      CLOSE (UNIT=13)
      CLOSE (UNIT=14)
      RETURN
      END
      SUBROUTINE  ELBAND (N, NG, IBW, LNODE)
C     * * * * * * * * * * * * * * * * * * * * * * * * * *
C     THIS SUBROUTINE DETERMINES THE BANDWITH
C     ASSOCIATED WITH A SINGLE ELEMENT
C     * * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  LNODE
      DIMENSION  LNODE(N)
C     LNODE = ELEMENT INCIDENCES
C     N     = NUMBER OF NODES PER ELEMENT
C     NG    = NUMBER OF PARAMETERS PER NODE
C     IBW   = HALF BANDWIDTH INCLUDING THE DIAGONAL
      IBW = 1
      NLESS = N - 1
      DO 20 I = 1,NLESS
        II = I + 1
        LNI = LNODE(I)
        DO 10 J = II,N
          IDIF = LNODE(J) - LNI
          NEW = NG*(IABS(IDIF) + 1)
          IF (NEW .GT. IBW)  IBW = NEW
   10   CONTINUE
   20 CONTINUE
      RETURN
      END
      SUBROUTINE  ELWAVE (M, NE, LFIRST, LLAST, LWAVE)
C     * * * * * * * * * * * * * * * * * * * * * * * * * *
C     BUILD WAVEFRONT LIST FOR EACH ELEMENT IN THE SYSTEM
C               ASSUMING 1 DOF PER NODE
C     * * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  LFIRST, LLAST, LWAVE
      DIMENSION  LFIRST(M), LLAST(M), LWAVE(NE)
C     M         = NUMBER OF NODES 
C     NE        = NUMBER OF ELEMENTS
C     LLAST(I)  = ELEMENT OF LAST APPEARENCE OF NODE I
C     LFIRST(I) = ELEMENT OF FIRST APPEARENCE OF NODE I
C     LWAVE(I)  = WAVEFRONT AT ELEMENT I
      LSTWAV = 0
      DO 10  I = 1,NE
   10 LWAVE(I) = 0
C      STORE NET FRONT INCREASE IN LWAVE
      DO 30  I = 1,M
        L = LFIRST(I)
C        IS THE NODE USED?
        IF ( L .EQ. 0 )  GO TO 30
C          INCREASE NUMBER OF ARRIVALS
          LWAVE(L) = LWAVE(L) + 1
C          CHECK DEPARTING NODES
          LP1 = LLAST(I) + 1
          IF ( LP1 .GT. NE )  GO TO 20
C            INCREASE NUMBER OF DEPARTURES
            LWAVE(LP1) = LWAVE(LP1) - 1
            GO TO 30
   20     LSTWAV = LSTWAV + 1
   30 CONTINUE
C      CALCULATE ELEMENT WAVEFRONTS
      DO 40  I = 2,NE
   40 LWAVE(I) = LWAVE(I) + LWAVE(I-1)
C      CHECK FOR FATAL INPUT ERRORS. IS EXITING WAVE ZERO?
      IF ( (LWAVE(NE) - LSTWAV) .EQ. 0 )  RETURN
      WRITE (12,*) 'LLAST OR LFIRST INCORRECT IN ELWAVE'
      STOP         'LLAST OR LFIRST INCORRECT IN ELWAVE'
      END
      SUBROUTINE  GETL2L (NE, LLMAX, LTOLPT, LADJL, LLUNIT)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C       INPUT THE ELEMENT NEIGHBOR LISTS USING THE NEW
C               ELEMENT ORDER NUMBERS
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
      DIMENSION  LTOLPT(NE+1), LADJL(1)
C     DIMENSION  LTOLPT(NE+1), LADJL(LLMAX)
C     LADJL     = PACKED LIST OF ELEMENTS ADJACENT TO ELEMENTS
C     LLMAX     = DIMENSION OF ARRAY LADJL = LTOLPT(NE+1)-1
C     LLUNIT    = UNIT NUMBER CONTAINING STORED DATA
C     LTOLPT(I) = POINTER WHERE NEIGHBORS OF I BEGIN IN LADJL
C     NE        = NUMBER OF ELEMENTS
C
C      WRITE ARRAY SIZES
      REWIND (LLUNIT)
      READ   (LLUNIT,*) NE, LLMAX
C     WRITE  (6,*)      NE, LLMAX
C      GET THE POINTER DATA
      READ  (LLUNIT,*) LTOLPT
      IF ( LTOLPT(NE+1) .NE. (LLMAX+1) ) STOP 'FATAL ERROR GETL2L'
C      INPUT NEIGHBORS WITH NEW ELEMENT NAMES
      KOUNT = 0
      DO 20  IE = 1, NE
        JSTART = LTOLPT(IE) - 1
        NEIGHB = LTOLPT(IE+1) - LTOLPT(IE)
        IF ( NEIGHB .GT. 0 )  THEN
          KOUNT = KOUNT + NEIGHB
C          GET NEIGHBORS 
          READ (LLUNIT,*) ( LADJL(JSTART+J), J = 1, NEIGHB )
        ELSE
          WRITE (6,*) 'WARNING: NO NEIGHBORS AT ELEMENT ', IE
        ENDIF
  20  CONTINUE
      IF ( KOUNT .NE. LLMAX ) STOP 'FATAL ERROR2 GETL2L'
      WRITE (6,*) 'ELEMENT NEIGHBORS READ'
      RETURN
      END
      SUBROUTINE  ICOMB (LIST, LONG)
C     ------------------------------------------
C     COMB SORT OF INTEGER LIST, SEE BYTE APR 91
C     ------------------------------------------
      DIMENSION  LIST(LONG)
      IGAP = LONG*10
  10  IGAP = MAX0 ( IGAP/13, 1 )
      ISWAP = 0
      DO 20  I = 1, (LONG-IGAP)
        J = I + IGAP
        IF ( LIST(I) .GT. LIST(J) )  THEN
C         SWAP THE ELEMENTS
          IHOLD   = LIST(I)
          LIST(I) = LIST(J)
          LIST(J) = IHOLD
          ISWAP   = ISWAP + 1
        ENDIF
  20  CONTINUE
      IF ( ISWAP .EQ. 0 .AND. IGAP .EQ. 1 )  THEN
        RETURN
      ELSE
        GO TO 10
      ENDIF
      END
      SUBROUTINE  INPUTV (NE, N, NNMAX, NTOLPT, NODES,
     1                    IFREE, MESH, INHEAD, LHEAD )
C     * * * * * * * * * * * * * * * * * * * * * * * * *
C             READ  INPUT DATA FOR MESH TOPOLOGY
C                PACK NODES INTO VECTOR ARRAY
C     * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  ITEMP, NTOLPT, NODES
      DIMENSION  ITEMP(36), NTOLPT(NE+1), NODES(1),
     1           LHEAD(NE,INHEAD)
      DATA       ITEMP / 36*0 /
C     IFREE     = SPACE AVAILABLE TO HOLD NODES ARRAY
C     INHEAD    = NUMBER OF INTEGERS BEFORE CONNECTIVITY
C     MESH      = ECHO TOPOLOGY INPUT IF > 0
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = TOTAL NUMBER OF ELEMENTS
C     NNMAX     = DIMENSION OF ARRAY NODES
C     NODES     = SYSTEM VECTOR OF ELEMENT INCIDENCES
C     NTOLPT(I) = PT WHERE INCIDENCES FOR ELEM I BEGIN
      IF ( MESH .GT. 0 )  THEN
        WRITE (8,*) 'ELEMENT   CONNECTIVITY DATA' 
      ENDIF
      NTOLPT(1) = 1
      NUMELM    = 0
      MAXM      = 0
      KOUNT     = 0
C-->   READ ELEMENT INCIDENCES DATA
 10   CONTINUE
      NUMELM = NUMELM + 1
      IF ( INHEAD .GT. 0 )  THEN
        READ (7,5000,END=99)  (LHEAD(NUMELM,J), J=1,INHEAD),
     1                        (ITEMP(I), I=1,12)
 5000   FORMAT ( 16I5 )
      ELSE
        READ (7,5000,END=99)  (ITEMP(I), I=1,12)
      ENDIF
      IF ( ITEMP(1) .EQ. 0 )  THEN     
C       A BLANK CARD ALSO TERMINATES INPUT
       ISUM = ITEMP(1)+ITEMP(2)+ITEMP(3)+ITEMP(4)+ITEMP(5)
     1      + ITEMP(6)+ITEMP(7)+ITEMP(8)+ITEMP(9)+ITEMP(10)
     2      + ITEMP(11)+ITEMP(12)
        IF (  ISUM .EQ. 0 ) GO TO 99
      ENDIF
C      MORE THAN 12 NODES
      IF ( N .GT. 12 )  THEN     
        READ  (7,5000,END=99)  (ITEMP(I), I=13,24)
      ENDIF
C      MORE THAN 24 NODES
      IF ( N .GT. 24 )  THEN     
        READ  (7,5000,END=99)  (ITEMP(I), I=25,36)
      ENDIF
C      FEWER NODES THAN N
      NUMNOD = N
   60 IF ( ITEMP(NUMNOD) .EQ. 0 )  THEN
        NUMNOD = NUMNOD - 1
        GO TO 60
      ENDIF   
      IF ( MESH .GT. 0 ) WRITE (8,5010) 
     1                   (NUMELM, ITEMP(I), I = 1, NUMNOD)
 5010 FORMAT ( 13I5, (/,12I5) )
      IF ( NUMNOD .EQ. 4 .AND. 
     1    (ITEMP(4) .EQ. ITEMP(3)) ) NUMNOD = 3
C      UPDATE POINTER, ADD ELEM NODES TO NODES VECTOR
      NTOLPT(NUMELM+1) = NTOLPT(NUMELM) + NUMNOD
      NPTL1 = NTOLPT(NUMELM) - 1
      DO 80  I = 1,NUMNOD
        KOUNT = KOUNT + 1
        IF ( KOUNT .GT. IFREE ) THEN
          WRITE (12,*) 'STORAGE EXCEEDED AT ELEMENT', NUMELM
          STOP         'STORAGE EXCEEDED'
        ENDIF
        IF ( ITEMP(I) .GT. MAXM )  MAXM = ITEMP(I)
   80 NODES(NPTL1+I) = ITEMP(I)
      GO TO 10
C<--   GET NEXT ELEM, OR LIST MAXIMUMS
   99 MAXELM = NUMELM - 1
      NNMAX  = NTOLPT(MAXELM + 1) - 1
      WRITE (6,*)  ' '
      WRITE (6,*)  'MAX ELEMENT NUMBER FOUND =', MAXELM
      WRITE (9,*)  'MAX ELEMENT NUMBER FOUND =', MAXELM
      WRITE (6,*)  'MAX NODE NUMBER FOUND    =', MAXM
      WRITE (9,*)  'MAX NODE NUMBER FOUND    =', MAXM
      IF ( MAXELM .LT. NE )  NE = MAXELM
      RETURN
      END
      SUBROUTINE  L2NOW (NE, LWAS, LNOW)
C     * * * * * * * * * * * * * * * * * * * * * * * *
C     CONVERT THE BEST LWAS TO LNOW 
C     * * * * * * * * * * * * * * * * * * * * * * * *
C     LNOW(I)   = NEW ELEMENT NUMBER OF OLD ELEMENT I
C     LWAS(I)   = OLD ELEMENT NUMBER OF NEW ELEMEMT I
C     NE        = NUMBER OF ELEMENTS 
      DIMENSION  LWAS(NE), LNOW(NE)
      DO 10  K = 1,NE
        KK = LWAS(K)
        IF ( KK .EQ. 0 )  GO TO 10
          LNOW(KK) = K
   10 CONTINUE
      RETURN
      END
      SUBROUTINE  LASTLV (NE, LFIRST, LLAST, M, N, NODES,
     1                    LNODE, NTOLPT, NNMAX, NOADJL, MAX)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C     TABULATE ELEMENT OF FIRST & LAST APPEARENCE OF EACH
C          NODE USING VECTOR STORAGE OF NODES
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  NODES, LFIRST, LLAST, LNODE, NTOLPT, NOADJL
      DIMENSION  NODES(NNMAX), LFIRST(M), LLAST(M), LNODE(N),
     1           NTOLPT(NE+1),  NOADJL(M)
C     LNODE     = INCIDENCES ARRAY OF A SINGLE ELEMENT
C     LFIRST(I) = ELEMENT OF FIRST APPEARENCE OF NODE I
C     LLAST(I)  = ELEMENT OF LAST APPEARENCE OF NODE I
C     M         = NUMBER NODES
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = NUMBER OF ELEMENTS 
C     NODES     = SYSTEM VECTOR OF ELEMENT INCIDENCES
C     NOADJL(I) = NUMBER OF ELEMENTS ADJACENT TO NODE I
C     NTOLPT    = POINTER ARRAY FOR ARRAY NODES
      DO 10  I = 1,M
        NOADJL(I) = 0
        LFIRST(I) = 0
   10 LLAST(I)  = 0
C-->   LOOP OVER ELEMENTS
      DO 30  J = 1,NE
        JJ = J
        L1 = NTOLPT(JJ)
        LN = NTOLPT(JJ+1) - L1
        CALL  LNODEV (LN, NODES(L1), LNODE)
        DO 20  I = 1,LN
          L = LNODE(I)
          IF ( L .LT. 1 )  GO TO 20
            IF ( LFIRST(L) .EQ. 0 )  LFIRST(L) = J
            LLAST(L)  = J
            NOADJL(L) = NOADJL(L) + 1
   20   CONTINUE
   30 CONTINUE
C<--   WARN OF NODE WITH NO ELEMENT CONNECTIONS
      KOUNT = 0
      MIN   = 0
      AVE   = 0.
      DO 40  I = 1,M
        AVE = AVE + NOADJL(I)
        IF ( NOADJL(I) .LT. MIN .AND.
     1       NOADJL(I) .GT. 1 )  MIN = NOADJL(I)
        IF ( LFIRST(I) .GT. 0 )  GO TO 40
        KOUNT = KOUNT + 1
        WRITE (12,*) 'NODE ', I,' NOT ON ANY ELEMENT'
   40 CONTINUE
      IF ( KOUNT .GT. 0 )  WRITE (12,*)  
     1    'TOTAL NUMBER OF OMITTED NODES =', KOUNT
C      ESTIMATE SIZE OF LADJL
      AVE = ( (AVE - 1)/FLOAT(M) + MIN )/2
      MAX = IFIX(AVE)*N*NE
      WRITE (6,*) 'LADJL SIZE ESTIMATE =', MAX
      RETURN
      END
      SUBROUTINE  LNODEV (N, NODES, LNODE)
C     * * * * * * * * * * * * * * * * * * * * * * * *
C     EXTRACT NODES FOR A SINGLE ELEMENT FROM 
C                 VECTOR MODE OPTION
C     * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  NODES, LNODE
      DIMENSION  LNODE(N), NODES(N)
C     NODES = COMPACTED SYSTEM ARRAY IN VECTOR FORM
C     LNODE = INCIDENCES LIST FOR SINGLE ELEMENT
C     N     = NUMBER OF NODES ON ELEMENT
      DO 10  I = 1,N
   10 LNODE(I) = NODES(I)
      RETURN
      END
      SUBROUTINE  LRESEQ (M, NL, NNMAX, LLMAX, LSTART, NNOW,
     1                    NOADJL, LWAS, LWAVE, NTOLPT, NODES, 
     2                    LTOLPT, LADJL, L1OFN, LUN2N, LUN2L )
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C            ELEMENT AND NODAL RESEQUENCING FOR FRONTAL
C                    AND BANDWIDTH SOLUTIONS
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  NOADJL, LWAS, LWAVE, NTOLPT, 
C I2 1           NODES, LTOLPT, LADJL, NNOW
      DIMENSION  NOADJL(M), LWAS(NL), LWAVE(NL), NTOLPT(NL+1), 
     1           NODES(NNMAX), LTOLPT(NL+1), LADJL(LLMAX), NNOW(M), 
     2           L1OFN(M), LUN2N(M), LUN2L(NL)
C...       --- INPUT DATA ---
C     LADJL     = LIST OF ELEMENTS ORDER ADJACENT TO EACH ELEM.
C     LLMAX     = DIMENSION OF ARRAY LADJL
C     LSTART    = SPECIFIED STARTING ELEMENT. 
C     LTOLPT(I) = PT IN LADJL WHERE ELEM NEIGHBORS OF ELEM I BEGIN
C     M         = TOTAL NUMBER OF NODES
C     NL        = TOTAL NUMBER OF ELEMENTS
C     NNMAX     = DIMENSION OF ARRAY NODES
C     NOADJL(I) = NO. OF ELEMENTS ADJACENT TO NODE I
C     NODES     = NODAL INCIDENCES ARRAY FOR ALL ELEMENTS
C     NTOLPT(I) = POINT IN NODES WHERE INCIDENCES FOR ELEM I BEGIN
C...      --- GENERATED ARRAYS ---
C     LWAS(I)  = OLD NUMBER OF NEW ELEMEMT I
C     LWAVE(I) = WAVEFRONT OF NEW ELEMENT I
C     NNOW(I)  = NEW NUMBER OF OLD NODE I
C...        --- TEMPORARY SCRATCH ARRAYS ---
C     LUN2L(I) = NO. OF UNNUMBERED ELEM ADJ TO ELEMENT I
C     LUN2N(I) = NO. OF UNNUMBERED ELEM ADJ TO NODE I
C     L1OFN(I) = ELEM OF 1-ST APPEARANCE OF NODE I
C        *** INITIALIZATION ***
C       LOOP COUNTERS
      LOPOLD = 1
      LOOP   = 0
      NLOOP  = 0
C      WAVE & LEAVING NODES @ I - 1
      IM1W   = 0
      IM1LN  = 0
C       SCRATCH ARRAYS 
      DO 10  I = 1,M
        LUN2N(I) = NOADJL(I)
        NNOW(I)  = 0
   10 L1OFN(I) = 0
      DO 20  I = 1,NL
        LWAS(I)  = 0
   20 LUN2L(I) = LTOLPT(I+1) - LTOLPT(I)
C       STARTING ITEMS
      LOLD = LSTART
      IT   = LSTART
      IAN  = NTOLPT(LSTART+1) - NTOLPT(LSTART)
C==>>     *** BEGIN ELEMENT RESEQUENCING LOOP ***
   60 LOOP = LOOP + 1
C      CALCULATE CURRENT FRONT SIZE
      IM1W  = IM1W + IAN - IM1LN
      IM1LN = 0
      LWAS(LOOP)  = IT
      LWAVE(LOOP) = IM1W
C      EXTRACT ITS' ELEMENT INCIDENCES LIST
      L1   = NTOLPT(IT)
      L1M1 = L1 - 1
      LN   = NTOLPT(IT+1) - L1
C-->    UPDATE 1ST APPEARANCE & ADJ EL LIST FOR IT'S NODES
      DO 70  II = 1,LN
C        OLD NODE NUMBER
        I = NODES(L1M1+II)
C        FLAG FIRST APPEARANCE ON NEW ELEMENT
        IF ( L1OFN(I) .EQ. 0 )  L1OFN(I) = LOOP
C        REDUCE NUMBER OF UN-NUMBERED ADJACENT ELEMENTS
        LUN2N(I) = LUN2N(I) - 1
        IF ( LUN2N(I) .EQ. 0 )  THEN    
C          NUMBER EXITING NODES AND STORE FOR NNOW
          IM1LN   = IM1LN + 1
          NLOOP   = NLOOP + 1
          NNOW(I) = NLOOP
        ENDIF
   70 CONTINUE
C<--    EXTRACT ELEMENT NEIGHBORS OF IT
      IF ( LOOP .EQ. NL )  GO TO 250
      K1 = LTOLPT(IT)
      K2 = LTOLPT(IT+1) - 1
C-->   LOOP OVER ELEM NEIGHBORS OF IT
      DO 80  K = K1,K2
C        CANDIDATE ELEMENT NUMBER
        KLEM = IABS( LADJL(K) )
C        UPDATE NO. OF UN-NUMBERED ELEM. NEIGHBORS
        LUN2L(KLEM) = LUN2L(KLEM) - 1
C        FIND ITS LOCATION IN NEIGHBOR LIST
        L1 = LTOLPT(KLEM)
        L2 = LTOLPT(KLEM+1) - 1
        LSIZE = L2 - L1 + 1
        CALL  WHERE (LSIZE, IT, LOC, LADJL(L1))
        LOC = L1 - LOC - 1
C        SET IT NEGATIVE IN LADJL TO FLAG AS NEW NUMBER
   80 LADJL(LOC) = -IT
C<--    ARE LOLD'S NEIGHBORS GONE
      IF ( LUN2L(LOLD) .EQ. 0 )  THEN
C-->     SELECT NEXT LOLD THAT IS STILL ACTIVE
        DO 90  KL = LOPOLD,LOOP
          KK  = KL
C          ORIGINAL ELEMENT NUMBER
          LKL = LWAS(KL)
C          USE FIRST WITH REMAINING NEIGHBORS
          IF ( LUN2L(LKL) .GT. 0 )  GO TO 100
   90   CONTINUE
C<--   UPDATE THE OLDEST VALUES
  100   LOPOLD = KK
        LOLD   = LWAS(KK)
      ENDIF
C        *** SELECT NEXT ELEMENT FROM NEIGHBORS OF LOLD ***
C      INITIALIZE COUNTER TESTS
      NEIGH  = NL+1
      MAXOUT = 0
      MAXOLD = 0
      MININ  = M+1
      L1OUT  = 0
      L2OUT  = 0
      L3OUT  = 0
C       FIND ELEMENT NEIGHBORS OF LOLD
      K1 = LTOLPT(LOLD)
      K2 = LTOLPT(LOLD+1) - 1
C-->>  LOOP OVER ELEMENT NEIGHBORS
      DO 240  K = K1,K2
        KLEM = LADJL(K)
C        HAS THIS CANDIDATE BEEN RENUMBERED?
        IF ( KLEM .LT. 0 )  GO TO 240
C        IF NOT, FIND INCIDENCES LIST FOR KLEM
        L1 = NTOLPT(KLEM)
        L1M1 = L1 - 1
        LN = NTOLPT(KLEM+1) - L1
C     ** COMPUTE COUNTERS FOR CANDIDATE KLEM **
        KOUNT1 = 0
        KOUNT2 = 0
        LOOK1  = 0
        LOOK2  = 0
        LOOK3  = 0
C->        LOOP OVER NODES OF CANDIDATE ELEMENT
        DO 120  II = 1,LN
          I = NODES(L1M1 + II)
C          FIRST APPEARANCE OF NODE
          IF ( L1OFN(I) .EQ. 0 )  KOUNT1 = KOUNT1 + 1
C          LAST APPEARANCE OF NODE
          IF ( LUN2N(I) .EQ. 1 )  KOUNT2 = KOUNT2 + 1
C          LOOK AHEAD FOR LAST APPEARANCE OF NODES
          IF ( LUN2N(I) .EQ. 2 )  LOOK1  = LOOK1 + 1
          IF ( LUN2N(I) .EQ. 3 )  LOOK2  = LOOK2 + 1
          IF ( LUN2N(I) .EQ. 4 )  LOOK3  = LOOK3 + 1
  120   CONTINUE
C<-   NOTE: WE COULD COMPUTE A WEIGHTED AVERAGE OF THE
C           ABOVE KOUNTS AND LOOKS AT THIS POINT.
C
C...     *** CURRENT SELECTION FOR NEXT ELEMENT ***
C
C        TRY TIE BREAKERS (EXPECT MANY TIES):
C         ARRIVING NEW NODES
        IF ( KOUNT1 - MININ  )  170,130,230
C         LEAVING NODES
  130   IF ( KOUNT2 - MAXOUT )  230,140,180
C         UN-NUMBERED ELEMENT NEIGHBORS
  140   KOUNT3 = LUN2L(KLEM)
        IF ( KOUNT3 - NEIGH  )  190,150,230
C          AGE OF ACTIVE NODES
  150   KOUNT4 = LN - KOUNT1
        IF ( KOUNT4 - MAXOLD )  230,160,200
C          NODES POSSIBLY LEAVING NEXT TIME
  160   IF ( LOOK1  - L1OUT  )  230,165,210
C          NODES POSSIBLY LEAVING TIME AFTER NEXT
  165   IF ( LOOK2  - L2OUT  )  230,166,215
C          LEAVING THIRD TIME
  166   IF ( LOOK3  - L3OUT  )  230,220,216
C        *** UPDATE TIE BREAKER VALUES ***
  170   MININ  = KOUNT1
  180   MAXOUT = KOUNT2
  190   NEIGH  = KOUNT3
  200   MAXOLD = KOUNT4
  210   L1OUT  = LOOK1
  215   L2OUT  = LOOK2
  216   L3OUT  = LOOK3
C        CANDIDATE SELECTED, ASSIGN TO IT
  220   IT     = KLEM
C        CANDIDATE REJECTED
  230   CONTINUE
C<<--   END OF ELEM NEIGHBOR LOOP
  240 CONTINUE
C      RECORD NUMBER OF NODES ENTERING AND LEAVING
      IAN = MININ
C     ILN = MAXOUT
C<<==   RETURN FOR THE NEXT REMAINING ELEMENT
      GO TO 60
C
C      *** ELEMENT RESEQUENCING COMPLETED ***
C
C      CHECK FOR ZERO FRONT UPON EXIT
  250 IF ( LWAVE(NL) .NE. IM1LN )  STOP
     1    'FRONT DOES NOT VANISH WITH LAST ELEMENT!'
C      RESET ARRAY LADJL FOR NEXT STARTING ELEMENT
      DO 280 K = 1,LLMAX
  280 LADJL(K) = IABS( LADJL(K) )
      RETURN
      END
      SUBROUTINE  LRMS (NE, LARRAY, MAX, RMS)
C     * * * * * * * * * * * * * * * * * * * * * * *
C     FIND MAXIMUM AND ROOT MEAN SQUARE VALUES
C     OF ELEMENT BANDWIDTH OR WAVEFRONT
C     * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  LARRAY
      DIMENSION  LARRAY(NE)
C     LARRAY  = BANDWIDTH OR WAVEFRONT OF EACH ELEM
C     MAX     = MAXIMUM VALUE
C     NE      = NUMBER OF ELEMENTS
C     RMS     = ROOT MEAN SQUARE VALUE
      MAX = LARRAY(1)
      RMS = MAX*MAX
      DO 10  I = 2,NE
        LTEST = LARRAY(I)
        IF ( LTEST .GT. MAX )  MAX = LTEST
   10 RMS = RMS + LTEST*LTEST
      RMS = SQRT( RMS/NE )
      RETURN
      END
      SUBROUTINE  LTOLV (NE, N, M, LFIRST, LLAST, NODES, LNODE, 
     1                   LADJL, LTOLPT, LLMAX, AVELL, NTOLPT,
     2                   NEED )
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C         BUILD LIST OF ELEMENTS ADJACENT TO OTHER ELEMENTS
C                 WITH NODES STORED AS VECTOR
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C                   SEE ORIGINAL ELTOEL
C I2  INTEGER*2  LFIRST, LLAST, NODES, LNODE, LTOLPT, LADJL, 
C I2 1           NTOLPT
      DIMENSION  LFIRST(M), LLAST(M), NODES(1), LNODE(N),
     1           NTOLPT(NE+1), LTOLPT(NE+1), LADJL(1)
C     LADJL     = PACKED LIST OF ELEMENTS ADJACENT TO ELEMENTS
C                NEIGHBORS OF ELEM L RUN FROM LADJL(LTOLPT(L))
C                TO LADJL(LTOLPT(L+1)-1)
C     LFIRST(I) = ELEMENT WHERE NODE I FIRST APPEARED  
C     LLAST(I)  = ELEMENT WHERE NODE I LAST APPEARED  
C     LLMAX     = SIZE OF LADJL (RETURNED), ESTIMATED (INPUT)
C     LNODE     = INCIDENCES ARRAY FOR A SINGLE ELEMENT
C     LTOLPT    = POINTER ARRAY FOR COMPRESSED ARRAY LADJL
C     M         = TOTAL NUMBER OF NODES
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = TOTAL NUMBER OF ELEMENTS
C     NEED      = COMMON NODES NEEDED TO BE A NEIGHBOR
C             1 - ALL (CORNER, EDGE, FACE)
C          >= 2 - ONLY EDGES AND FACES
C          OTHER- FACES ONLY
C     NODES     = SYSTEM VECTOR OF INCIDENCES OF ALL ELEMENTS
C     NTOLPT    = POINTER ARRAY FOR COMPRESSED ARRAY NODES
      IF ( NEED .LT. 1 )  NEED = 1
      KOUNT = 1
C==>  LOOP OVER ALL ELEMENTS
      DO 60  IE = 1,NE
C        SET BEGINNING POINTER FOR ELEMENT IE
        LTOLPT(IE) = KOUNT
C        EXTRACT INCIDENCES LIST FOR ELEMENT IE
        L1     = NTOLPT(IE)
        LN     = NTOLPT(IE+1)-L1
        CALL  LNODEV (LN, NODES(L1), LNODE)
        LIN    = LNODE(1)
        LSTART = LFIRST(LIN)
        LSTOP  = LLAST(LIN)
C        ESTABLISH RANGE OF POSSIBLE ELEMENT NEIGHBORS
        DO 10  IN = 2,LN
          LIN = LNODE(IN)
          IF ( LSTART .GT. LFIRST(LIN) )  LSTART = LFIRST(LIN)
          IF ( LSTOP  .LT. LLAST(LIN)  )  LSTOP  = LLAST(LIN)
   10   CONTINUE
C-->     LOOP OVER POSSIBLE ELEMENT NEIGHBORS
        DO 50  L = LSTART,LSTOP
          IF ( L .EQ. IE )  GO TO 50
C          LOOP OVER INCIDENCES OF POSSIBLE ELEMENT NEIGHBOR
          L1     = NTOLPT(L) - 1
          NPN    = NTOLPT(L+1) - L1 - 1
          NFOUND = 0
          DO 30  IN = 1,NPN
            NTEST = NODES(L1+IN)
            IF ( LLAST(NTEST)  .LT. IE )  GO TO 30
            IF ( LFIRST(NTEST) .GT. IE )  GO TO 30
C            COMPARE WITH INCIDENCES OF ELEMENT IE
            DO 20  JN = 1,LN
c             IF ( LNODE(JN) .EQ. NTEST )  GO TO 40
              IF ( LNODE(JN) .EQ. NTEST )  THEN
C                FOUND ENOUGH NODES?
                NFOUND = NFOUND + 1
                IF ( NFOUND .EQ. NEED )  GO TO 40
              ENDIF
   20       CONTINUE
   30     CONTINUE
C          NOT A NEIGHBOR
          GO TO 50
C          FOUND ONE, INSERT INTO LADJL, UPDATE POINTER
   40     LADJL(KOUNT) = L
          KOUNT = KOUNT + 1
C          CHECK DIMENSION LIMITS
          IF ( KOUNT .GT. LLMAX )  THEN
            WRITE (12,*) 'DIMENSION EXCEEDED IN LTOLV',
     1                   'FOR ELEMENT', IE
            WRITE (12,*)  KOUNT, LLMAX
            STOP         'DIMENSION EXCEEDED IN LTOLV'
          ENDIF
   50   CONTINUE
C<--
   60 CONTINUE
      IF ( KOUNT .EQ. 1 ) STOP 'NO ELEMENT NEIGHBORS'
C<==   LOCATE END OF ARRAY
      LLMAX = KOUNT - 1
      LTOLPT(NE+1) = KOUNT
      DO 70  I = 1,NE
      IF ( (LTOLPT(I+1) - LTOLPT(I)) .LE. 0 )  WRITE (12,*)
     1   'WARNING: ELEMENT', I,' HAS NO ELEMENT NEIGHBORS'
   70 CONTINUE
      AVELL = FLOAT( LLMAX )/NE
      RETURN
      END
      PROGRAM  MAIN
C     DUMMY MAIN TO SET UP STORAGE FOR REDUCE
CIBM  INTEGER*2 ID
      PARAMETER  ( IFIXED = 75000 )
      DIMENSION  ID(IFIXED)
      CALL  DRIVER (IFIXED, ID)
      STOP 'NORMAL EXIT OF REDUCE'
      END
      SUBROUTINE  N2WAS (M, NNOW, NWAS)
C     * * * * * * * * * * * * * * * * * * * *
C     CONVERT THE BEST NNOW TO NWAS 
C     * * * * * * * * * * * * * * * * * * * *
      DIMENSION  NNOW(M), NWAS(M)
C     M         = NUMBER OF NODES 
C     NNOW(I)   = NEW NUMBER OF OLD NODE I
C     NWAS(I) = OLD NODE NUMBER OF NEW NODE I
      DO 10  K = 1,M
        KK = NNOW(K)
        IF ( KK .EQ. 0 )  GO TO 10
          NWAS(KK) = K
   10 CONTINUE
      RETURN
      END
      SUBROUTINE NEWDAT (M, NE, N, NNMAX, IOPT, NTOLPT, NODES, 
     1                   NNOW, LWAS, LNODE, INHEAD, LHEAD)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C     OUTPUT NEW ELEMENT (AND HEADER) AND/OR NODE LISTS
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * *
      DIMENSION NTOLPT(NE+1), NODES(NNMAX), NNOW(M), LWAS(NE), 
     1          LNODE(N), LHEAD(NE,INHEAD)
C     IOPT = OUTPUT OPTIONS
C          0 - NEW ELEMENT ORDER WITH NEW NODE NUMBERS
C          1 - NEW ELEMENT ORDER WITH OLD NODE NUMBERS
C          2 - OLD ELEMENT ORDER WITH NEW NODE NUMBERS
C
C-->   LOOP OVER ELEMENTS
      DO 10 L = 1,NE
        IF ( IOPT .LT. 2 ) THEN
           LOUT = LWAS(L)
        ELSE
           LOUT = L
        ENDIF
C        LOCATE ORIGINAL NODES ON ELEMENT LOUT
        K1    = NTOLPT(LOUT)
c       K1    = NTOLPT(LOUT) - 1
        KSIZE = NTOLPT(LOUT+1) - NTOLPT(LOUT)
        CALL  LNODEV (KSIZE, NODES(K1), LNODE)
C->      LOOP OVER ORIGINAL NODES OF ELEMENT
        DO 20 K = 1,KSIZE
          NOUT = LNODE(K)
c         NOUT = NODES(K1+K)
          IF ( IOPT .NE. 1 .AND. NOUT .GT. 0 )
     1      NOUT = NNOW(NOUT)
          LNODE(K)= NOUT
   20   CONTINUE
C<-      OUTPUT THE ELEMENT HEADERS AND ITS TOPOLOGY
        IF ( INHEAD .GT. 0 )  THEN
          WRITE (11,5000)  ( LHEAD(LOUT,J), J = 1,INHEAD),
     1                     ( LNODE(K), K = 1,KSIZE )
        ELSE
          WRITE (11,5000)  ( LNODE(K), K = 1,KSIZE )
 5000     FORMAT ( 16I5, (/, 16I5) )
        ENDIF
   10   CONTINUE
C<--
      RETURN
      END
      SUBROUTINE  NEWL2L (NE, LLMAX, LTOLPT, LADJL, LNOW, 
     1                    LWAS, LTEMP)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C      OUTPUT THE ELEMENT NEIGHBOR LISTS USING THE NEW
C               ELEMENT ORDER NUMBERS
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
      DIMENSION  LNOW(NE), LWAS(NE), LTOLPT(NE+1), 
     1           LADJL(LLMAX), LTEMP(NE+1)
C     LADJL     = PACKED LIST OF ELEMENTS ADJACENT TO ELEMENTS
C     LLMAX     = DIMENSION OF ARRAY LADJL = LTOLPT(NE+1)-1
C     LNOW(I)   = NEW ELEMENT NUMBER OF OLD ELEMENT I
C     LWAS(I)   = OLD ELEMENT NUMBER OF NEW ELEMENT I
C     LTEMP     = A SCRATCH ARRAY. POINTER, THEN NEIGHBORS
C     LTOLPT(I) = POINTER WHERE NEIGHBORS OF I BEGIN IN LADJL
C     NE        = NUMBER OF ELEMENTS
C
C      WRITE ARRAY SIZES
      REWIND (14)
      WRITE  (14,*) NE, LLMAX
C      CONVERT THE POINTER DATA
      LTEMP(1) = 1
      DO 10  IE = 1, NE
C        IE = ORIGINAL ORDER, JE = NEW ORDER
        JE     = LWAS(IE)
        NEIGHB = LTOLPT(JE+1) - LTOLPT(JE)
        LTEMP(IE+1) = LTEMP(IE) + NEIGHB
  10  CONTINUE
      WRITE (14,*) LTEMP
      IF ( LTEMP(NE+1) .NE. (LLMAX+1) ) STOP 'FATAL ERROR NEWL2L'
C      OUTPUT NEIGHBORS WITH NEW ELEMENT NAMES
      KOUNT = 0
      DO 20  IE = 1, NE
C        IE = ORIGINAL ORDER, JE = NEW ORDER
        JE     = LWAS(IE)
        JSTART = LTOLPT(JE) - 1
        NEIGHB = LTOLPT(JE+1) - LTOLPT(JE)
        IF ( NEIGHB .GT. 0 )  THEN
          KOUNT = KOUNT + NEIGHB
C          GET NEIGHBORS AND RE-NAME THEM
          DO 30  J = 1, NEIGHB
            LOLD     = LADJL(JSTART + J)
            LTEMP(J) = LNOW(LOLD)
  30      CONTINUE
C          SORT NEW NEIGHBORS IN ASSENDING ORDER
          CALL ICOMB (LTEMP, NEIGHB)
          WRITE (14,*) ( LTEMP(J), J = 1, NEIGHB )
        ELSE
          WRITE (6,*) 'WARNING: NO NEIGHBORS AT NEW ELEMENT ', JE
        ENDIF
  20  CONTINUE
      IF ( KOUNT .NE. LLMAX ) STOP 'FATAL ERROR2 NEWL2L'
      WRITE (6,*) 'RE-ORDERED RE-NAMED ELEMENT NEIGHBORS SAVED'
      RETURN
      END
C     ID(L21)   = LADJL, ELEM ADJACENT TO ELEM
C     IHERE     = STAGE OF RESTART
C     IPICK     = 0 RMS VALUE GOVERNS, = 1 MAX GOVERNS
C     ISAVE     = RESTART CONTROL
C     L1OFN(I)  = ELEM OF 1-ST APPEARANCE OF NODE I
C     LADJL     = PACKED LIST OF ELEMENTS ADJACENT TO ELEMENTS
C     LBAND(I)  = HALF BANDWIDTH AT ELEMENT I
C     LFIRST(I) = ELEMENT OF FIRST APPEARENCE OF NODE I
C     LHERE     = STAGE OF LAST RESTART SAVE
C     LLAST(I)  = ELEMENT OF LAST APPEARENCE OF NODE I
C     LLMAX     = DIMENSION OF ARRAY LADJL
C     LNODE     = INCIDENCES ARRAY OF A SINGLE ELEMENT
C     LNOPT     = 0 OUTPUT ELEMENTS, =1 OUTPUT NODES
C     LNOW(I)   = NEW ELEMENT NUMBER OF OLD ELEMENT I
C     LSTART    = SPECIFIED STARTING ELEMENT.
C     LTOLPT(I) = PT WHERE ELEM NEIGHBORS OF ELEM I BEGIN
C     LUN2L(I)  = NO. OF UNNUMBERED ELEM ADJ TO ELEMENT I
C     LUN2N(I)  = NO. OF UNNUMBERED ELEM ADJ TO NODE I
C     LWAS(I)   = OLD NUMBER OF NEW ELEMEMT I
C     LWAVE(I)  = WAVEFRONT AT ELEMENT I
C     M         = NUMBER OF NODES
C     MESH      = MESH ECHO FLAG, 0 = NONE
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = NUMBER OF ELEMENTS
C     NEED      = NUMBER OF COMMON NODES TO BE NEIGHBOR
C     NNMAX     = DIMENSION OF ARRAY NODES
C     NNOW(I)   = NEW NUMBER OF OLD NODE I
C     NOADJL(I) = NUMBER OF ELEMENTS ADJACENT TO NODE I
C     NODES     = SYSTEM VECTOR OF ELEMENT INCIDENCES
C     NTOLPT(I) = PT WHERE INCIDENCES FOR ELEM I BEGIN
      SUBROUTINE  REDUCE (M, NE, N, NEED, NUSER, NTRY, LNOPT, 
     1                    LISTLB, LISTLW, LISTLL, LNEWN, LNEWL, 
     2                    LFIRST, LLAST, L1OFN, LUN2N, NOADJL,
     3                    NNOW, NWAS, LNODE, LWAVE, LBAND,
     4                    LNOW, LWAS, LUN2L, LTOLPT, NTOLPT,
     5                    LSTART, LHEAD, NODES,  ID, L20, IFREE, 
     6                    IFIXED, MESH, ISAVE, LHERE, IPICK, 
     7                    INHEAD, L2L)
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C                 REAL MAIN PROGRAM TO CONTROL REDUCE
C     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  NTOLPT,NODES,LTOLPT,LADJL,NOADJL,LWAVE,LBAND,
C I2 1           L1OFN,LUN2N,LUN2L,LNOW,LWAS,NNOW,NWAS,LNODE
      PARAMETER  ( IONE = 1 )
      DIMENSION  LFIRST(M), LLAST(M), L1OFN(M), LUN2N(M), NOADJL(M),
     1           NNOW(M), NWAS(M), LNODE(N), LWAVE(NE), LBAND(NE),
     2           LNOW(NE), LWAS(NE), LUN2L(NE), LTOLPT(NE+1),
     3           NTOLPT(NE+1), LSTART(NUSER+NTRY), NODES(1), 
     4           LHEAD(NE,INHEAD), ID(IFIXED)
C     ID(L20)   = NODES, ELEMENTS CONNECTIVITY LIST
C     ID(L21)   = LADJL, ELEM ADJACENT TO ELEM
C     IHERE     = STAGE OF RESTART
C     INHEAD    = NUMBER OF HEADERS BEFORE CONNECTIVITY
C     IPICK     = 0 RMS VALUE GOVERNS, = 1 MAX GOVERNS
C     ISAVE     = RESTART CONTROL
C     L1OFN(I)  = ELEM OF 1-ST APPEARANCE OF NODE I
C     LADJL     = PACKED LIST OF ELEMENTS ADJACENT TO ELEMENTS
C     LBAND(I)  = HALF BANDWIDTH AT ELEMENT I
C     LFIRST(I) = ELEMENT OF FIRST APPEARENCE OF NODE I
C     LHEAD     = ELEMENT'S INTEGER HEADERS (PROB DEPENDENT)
C     LHERE     = STAGE OF LAST RESTART SAVE
C     LLAST(I)  = ELEMENT OF LAST APPEARENCE OF NODE I
C     LLMAX     = DIMENSION OF ARRAY LADJL
C     LNODE     = INCIDENCES ARRAY OF A SINGLE ELEMENT
C     LNOPT     = 0 OUTPUT ELEMENTS, =1 OUTPUT NODES
C     LSTART    = SPECIFIED STARTING ELEMENT.
C     LTOLPT(I) = PT WHERE ELEM NEIGHBORS OF ELEM I BEGIN
C     LUN2L(I)  = NO. OF UNNUMBERED ELEM ADJ TO ELEMENT I
C     LUN2N(I)  = NO. OF UNNUMBERED ELEM ADJ TO NODE I
C     LWAS(I)   = OLD NUMBER OF NEW ELEMEMT I
C     LWAVE(I)  = WAVEFRONT AT ELEMENT I
C     M         = NUMBER OF NODES
C     MESH      = MESH ECHO FLAG, 0 = NONE
C     N         = NUMBER OF NODES PER ELEMENT
C     NE        = NUMBER OF ELEMENTS
C     NEED      = NUMBER OF COMMON NODES TO BE NEIGHBOR
C     NNMAX     = DIMENSION OF ARRAY NODES
C     NNOW(I)   = NEW NUMBER OF OLD NODE I
C     NOADJL(I) = NUMBER OF ELEMENTS ADJACENT TO NODE I
C     NODES     = SYSTEM VECTOR OF ELEMENT INCIDENCES
C     NTOLPT(I) = PT WHERE INCIDENCES FOR ELEM I BEGIN
C      RESTART ENTRIES ?
      IF ( ISAVE .NE. 2 )  THEN
        MAXTRY = NUSER + NTRY
      ELSE
        IF ( LHERE .EQ. 0 )  STOP 'NO RESTART SAVED'
          LFINAL = ID(1)
          NNMAX  = ID(3) 
          ISTART = ID(4) 
          LKEEP  = ID(5) 
          MKEEP  = ID(6) 
          RKEEP  = FLOAT( ID(7) )/100.
          IF ( LNOPT .EQ. 0 )   THEN
            WRMSL  = RKEEP
            LWORIG = MKEEP
          ELSE
            BRMSL  = RKEEP
            LBORIG = MKEEP
          ENDIF
          LLMAX  = ID(8) 
          AVELL  = FLOAT( ID(9) )/100.
          L21    = ID(10)
          MAXTRY = ID(11)
        IF ( LHERE .EQ. 1 )  GO TO 40
        IF ( LHERE .EQ. 2 )  GO TO 50
      ENDIF
C      READ ELEMENT TOPOLOGIES
      NL  = NE
      NIN = N
      CALL INPUTV (NL, NIN, NNMAX, NTOLPT, NODES, IFREE, 
     1             MESH, INHEAD, LHEAD )
      IF ( NIN .GT. N ) STOP 'N EXCEEDS CONTROL MAX'
      L21   = L20 + NNMAX
      IF ( L21 .GT. IFIXED ) STOP 'IFIXED TOO SMALL'
      IFREE = IFREE - NNMAX
C
C....        ** DETERMINE ELEMENT WAVEFRONTS **
C
      CALL  LASTLV (NE, LFIRST, LLAST, M, N, NODES,
     1              LNODE, NTOLPT, NNMAX, NOADJL, MAX)
C      ESTIMATE AVAILABLE SPACE
      LLMAX = IFREE
      MAX = MAX + L21
      WRITE (6,*) 'IFIXED SIZE ESTIMATE =', MAX
      CALL  ELWAVE (M, NE, LFIRST, LLAST, LWAVE)
      ISTART = 0
      LKEEP  = 0
      LTEST  = 0
      WRITE (9,*) 'TRIAL RESULTS'
      WRITE (9,*) ' TRY   ELEM    MAX       RMS'
      IF ( LNOPT .EQ. 0 )  THEN
C....     *** PRINT ORIGINAL FRONTS ***
        CALL  LRMS (NE, LWAVE, LWORIG, WRMSL)
        WRITE (6,*) ' '
        WRITE (6,*) 'RMS ELEMENT WAVEFRONT.......', WRMSL
        WRITE (6,*) 'MAXIMUM ELEMENT WAVEFRONT...', LWORIG
        MKEEP = LWORIG
        RKEEP = WRMSL
        WRITE (9,5000) ISTART, LKEEP, MKEEP, RKEEP
 5000   FORMAT ( I4, I7, I7, F9.2 )
        IF ( LISTLW .NE. 0 )  THEN
          WRITE (6,*) ' '
          WRITE (6,*) '          ORIGINAL'
          WRITE (6,*) 'LIST      ELEMENT WAVEFRONT'
          CALL  WRTOUT (IONE, IONE, NE, LWAVE)
        ENDIF
      ELSE
C....     *** DETERMINE ORIGINAL BANDWIDTH ***
        CALL  BANWTV (NE, N, NTOLPT, NODES, LNODE, LBAND,
     1                IONE, IONE)
        CALL  LRMS (NE, LBAND, LBORIG, BRMSL)
        WRITE (6,*) ' '
        WRITE (6,*) 'RMS ELEMENT BANDWIDTH.......', BRMSL
        WRITE (6,*) 'MAXIMUM ELEMENT BANDWIDTH...', LBORIG
        MKEEP = LBORIG
        RKEEP = BRMSL
        WRITE (9,5000) ISTART, LKEEP, MKEEP, RKEEP
        IF ( LISTLB .NE. 0 )  THEN
          WRITE (6,*) ' '
          WRITE (6,*) '          ORIGINAL'
          WRITE (6,*) 'LIST      ELEMENT BANDWIDTH'
          CALL  WRTOUT ( IONE, IONE, NE, LBAND)
        ENDIF
      ENDIF
C      RESTART SAVES ?
      IF ( ISAVE .EQ. 1 ) THEN
        IHERE  = 1
        LFINAL = L21 - 1
        ID(1) = LFINAL
        ID(2) = IHERE
        ID(3) = NNMAX
        ID(4) = ISTART
        ID(5) = LKEEP
        ID(6) = MKEEP
        ID(7) = RKEEP*100
        ID(8) = 0
        ID(9) = 0
        ID(10) = 0
        ID(11) = MAXTRY
C       REWIND 12
C       WRITE (12,*)  ID(1), ID(2)
C       WRITE (12,*)  ( ID(K), K = 3,LFINAL )
        REWIND 13
        WRITE (13)  ID(1), ID(2)
        WRITE (13)  ( ID(K), K = 3, LFINAL )
      ENDIF
C      RESTART ENTRY ?
   40 CONTINUE
C
C....   ** BUILD LIST OF ELEMENTS ADJACENT TO ELEMENT **
C
C      NOADJL AVAILABLE TO HELP HERE
      CALL  LTOLV (NE, N, M, LFIRST, LLAST, NODES, LNODE,
     1             ID(L21), LTOLPT, LLMAX, AVELL, NTOLPT,
     2             NEED)
      LFINAL = L21 + LLMAX
      WRITE (6,*) ' '
      WRITE (6,*) 'ELEMENT DEGREE OF ELEMENTS =', AVELL
      WRITE (6,*) 'LADJL ACTUAL SIZE ........ =', LLMAX
      WRITE (6,*) 'STORAGE: GIVEN', IFIXED,' USED', LFINAL
      WRITE (6,*) ' '
C      RESTART SAVES ?
      IF ( ISAVE .EQ. 1 ) THEN
        IHERE = 2
        ID(1) = LFINAL
        ID(2) = IHERE
        ID(3) = NNMAX
        ID(4) = ISTART
        ID(5) = LKEEP
        ID(6) = MKEEP
        ID(7) = RKEEP*100
        ID(8) = LLMAX
        ID(9) = AVELL*100
        ID(10) = L21
        ID(11) = MAXTRY
        REWIND 13
        WRITE (13)  ID(1), ID(2)
        WRITE (13)  ( ID(K), K = 3, LFINAL)
C       REWIND 12
C       WRITE (12,*)  ID(1), ID(2)
C       WRITE (12,*)  ( ID(K), K = 3, LFINAL)
      ENDIF
C       PRINT THE ARRAY ?
      IF ( LISTLL .NE. 0 )  THEN
        WRITE (6,*) 'ORIGINAL'
        WRITE (6,*) 'ELEMENT   ADJACENT ELEMENTS'
        CALL  WRTOUT (NE, LTOLPT, LLMAX, ID(L21))
      ENDIF
C       STORE THE ARRAY AND STOP (NO RESEQUENCING) ?
      IF ( L2L .EQ. -1 )  THEN
C        WRITE ARRAY SIZES
        LLUNIT = 14
        REWIND (LLUNIT)
        WRITE  (LLUNIT,*) NE, LLMAX
C        GET THE POINTER DATA
        WRITE (LLUNIT,*) LTOLPT
C        OUTPUT NEIGHBORS 
        DO 15  IE = 1, NE
          JSTART = L21 + LTOLPT(IE) - 2
          NEIGHB = LTOLPT(IE+1) - LTOLPT(IE)
          IF ( NEIGHB .GT. 0 )  THEN
            KOUNT = KOUNT + NEIGHB
C            GET NEIGHBORS 
            WRITE(LLUNIT,*) ( ID(JSTART+J), J = 1, NEIGHB )
          ELSE
            WRITE (6,*) 'WARNING: NO NEIGHBORS AT ELEMENT ', IE
          ENDIF
  15    CONTINUE
        WRITE (6,*) 'ORIGINAL ELEMENT NEIGHBORS SAVED'
C        CHECK NEIGHBORS LIST
        CALL  GETL2L (NE, IDUMMY, LTOLPT, ID(L21), LLUNIT)
        WRITE (6,*) ' '
        WRITE (6,*) 'RE-ORDERED'
        WRITE (6,*) 'ELEMENT   ADJACENT ELEMENTS'
        CALL  WRTOUT (NE, LTOLPT, IDUMMY, ID(L21))
        IF ( LNOPT .EQ. 0 ) STOP 'REQUESTED NO RE-ORDERING'
      ENDIF
      IF ( ISAVE .EQ. 2 )  GO TO 50
C      SELECT NTRY ELEMENTS WITH FEWEST NEIGHBORS
      MTRY = MIN0(NTRY, NE)
      INC  = NE/MTRY
C      START WITH EQUAL SPACES
      DO 20  I = 1, MTRY
   20 LSTART(NUSER+I) = 1 + (I-1)*INC
      LLAVE = AVELL
      KOUNT = 0
C->    NOW TRY TO FIND SMALLEST NEIGHBORS
      DO 30  IE = 1, NE
      LSIZE = LTOLPT(IE+1) - LTOLPT(IE)
C      FEWER THAN AVERAGE?
      IF ( LSIZE .LT. LLAVE )  THEN
        KOUNT = KOUNT + 1
C        LESS THAN CURRENT MIN?
        KLEM  = LSTART(NUSER+KOUNT)
        KSIZE = LTOLPT(KLEM+1) - LTOLPT(KLEM)
        IF ( LSIZE .LE. KSIZE ) LSTART(NUSER+KOUNT) = IE
        IF ( KOUNT .EQ. MTRY )  KOUNT = 0
      ENDIF
   30 CONTINUE
C      RESTART ENTRY?
   50 CONTINUE
C
C...  ** READ USER SUPPLIED STARTING ELEMENTS ***
C
      IF ( NUSER .GT. MAXTRY )  NUSER = MAXTRY
      IF ( NUSER .GT. 0 )  THEN
        WRITE (6,*) 'USER SUPPLIED STARTING ELEMENTS:'
        DO 10  L = 1, NUSER
        READ  (5,*)  LSTART(L)
   10   WRITE (6,*)  L, LSTART(L)
      ENDIF
C<-
C==>    LOOP OVER STARTING ELEMENTS
C
      WRITE (6,*) ' '
      WRITE (6,*) 'RMS ELEMENT WAVEFRONT.......', WRMSL
      WRITE (6,*) 'MAXIMUM ELEMENT WAVEFRONT...', LWORIG
      WRITE (6,*) ' '
      DO  200  ISTART = 1, MAXTRY
      LTEST = LSTART(ISTART)
      IF ( LTEST. LT. 1 .OR. LTEST .GT. NE )  THEN
        WRITE (6,*) 'INVALID START, USING NE'
        LTEST = NE
      ENDIF
      WRITE (6,*) '==== TRIAL NUMBER ======', ISTART
      WRITE (6,*) 'STARTING ELEMENT .......', LTEST
C
C....        ** RESEQUENCE ELEMENTS  **
C
      CALL  LRESEQ (M, NE, NNMAX, LLMAX, LTEST,
     1              NNOW, NOADJL, LWAS, LWAVE,
     2              NTOLPT, NODES, LTOLPT, ID(L21),
     3              L1OFN, LUN2N, LUN2L )
C
C....        ** UPDATE ELEMENTS OR NODES **
C
      IF ( LNOPT .EQ. 0 )  THEN
C         ---  GIVE NEW ELEMENT  RESULTS ---
        CALL  LRMS (NE, LWAVE, MAXWAV, RMSLW)
        WRITE (9,5000)  ISTART, LTEST, MAXWAV, RMSLW
        IF ( IPICK .EQ. 0 )  THEN
          IF ( RMSLW  .GT. WRMSL  )  GO TO 180
          IF ( MAXWAV .LE. LWORIG )  GO TO 120
          IF ( MAXWAV .GT. LWORIG )  GO TO 180
        ELSE
          IF ( MAXWAV .GT. LWORIG )  GO TO 180
          IF ( MAXWAV .LE. LWORIG )  GO TO 120
          IF ( RMSLW  .GT. WRMSL  )  GO TO 180
        ENDIF
C        RECORD THE NEWEST BEST RESULTS
  120   LWORIG = MAXWAV
        WRMSL  = RMSLW
        LKEEP  = LTEST
        MKEEP  = MAXWAV
        RKEEP  = RMSLW
        WRITE (6,*) 'RMS ELEMENT WAVEFRONT     = ', RMSLW
        WRITE (6,*) 'MAXIMUM ELEMENT WAVEFRONT = ', MAXWAV
        REWIND  10
        WRITE  (10)  LWAS, LWAVE
      ELSE
C         --- GIVE NEW NODE  RESULTS ---
        CALL  BANWTV (NE, N, NTOLPT, NODES, LNODE,
     1                LBAND, M, NNOW )
        CALL  LRMS (NE, LBAND, MAXBAN, RMSLB)
        WRITE (9,5000) ISTART, LTEST, MAXBAN, RMSLB
        IF ( IPICK .EQ. 0 )  THEN
          IF ( RMSLB  .GT. BRMSL  )  GO TO 180
          IF ( MAXBAN .LE. LBORIG )  GO TO 150
          IF ( MAXBAN .GT. LBORIG )  GO TO 180
        ELSE
          IF ( RMSLB  .GT. BRMSL  )  GO TO 180
          IF ( MAXBAN .LE. LBORIG )  GO TO 150
          IF ( MAXBAN .GT. LBORIG )  GO TO 180
        ENDIF
C        RECORD NEWEST BEST RESULTS
  150   LBORIG = MAXBAN
        BRMSL  = RMSLB
        LKEEP  = LTEST
        MKEEP  = MAXBAN
        RKEEP  = RMSLB
        WRITE (6,*) 'RMS ELEMENT BANDWIDTH     = ', RMSLB
        WRITE (6,*) 'MAXIMUM ELEMENT BANDWIDTH = ', MAXBAN
C...     STORE NNOW FOR REDUCED BANDWIDTH
        REWIND 10
        WRITE  (10)  NNOW
      ENDIF
  180 WRITE (6,*)  ' '
C<===     TRY NEXT STARTING ELEMENT
  200 CONTINUE
C      OUTPUT BEST SUMMARY
      WRITE (6,*) ' '
      WRITE (6,*) 'BEST STARTING ELEMENT =', LKEEP
      WRITE (6,*) 'BEST RMS VALUE        =', RKEEP
      WRITE (6,*) 'BEST MAXIMUM VALUE    =', MKEEP
      ISTART = 0
      WRITE (9,5000) ISTART, LKEEP, MKEEP, RKEEP
      WRITE (9,*) ' '
      WRITE (9,*) 'ELEMENT DEGREE OF ELEMENTS =', AVELL
      WRITE (9,*) 'LADJL ACTUAL SIZE ........ =', LLMAX
      WRITE (9,*) 'STORAGE: GIVEN', IFIXED,' USED', LFINAL
C      ANY IMPROVEMENT?
      IF ( LKEEP .EQ. 0 )  THEN
        WRITE (6,*) 'NO IMPROVEMENT'
        WRITE (9,*) 'NO IMPROVEMENT'
        RETURN
      ENDIF
      IF ( LNOPT .EQ. 0 )  THEN
C      --- RECOVER BEST ELEMENT RESULTS ---
      REWIND 10
      READ  (10) LWAS, LWAVE
C         PRINT FRONT RESULTS
        IF ( LISTLW .NE. 0 )  THEN
          WRITE (6,*) ' '
          WRITE (6,*) 'LIST      ELEMENT WAVEFRONT'
          CALL  WRTOUT (IONE, IONE, NE, LWAVE)
        ENDIF
C         PRINT NEW ELEMENT SEQUENCE
        CALL  L2NOW  (NE, LWAS, LNOW)
        IF ( LNEWL .GT. 0 )  THEN
          WRITE (6,*) ' '
          WRITE (6,*) '      NEW ELEMENT ORDER DATA'
          CALL  WASNOW (NE, LNOW, LWAS)
        ENDIF
C        SAVE NEW ELEMS AND OLD NODES
        IOPT = 1
        CALL  NEWDAT (M, NE, N, NNMAX, IOPT, NTOLPT, NODES, 
     1                NNOW, LWAS, LNODE, INHEAD, LHEAD)
C        SAVE NEW ORDER AND NEW ELEMENT NEIGHBORS,
C           ( OVERWRITE NTOLPT )
        CALL  NEWL2L (NE, LLMAX, LTOLPT, ID(L21),
     1                LNOW, LWAS, NTOLPT)
      ELSE
C        --- RECOVER THE BEST NODE RESULTS ---
        REWIND 10
        READ  (10)  NNOW
C        PRINT BAND RESULTS
        IF ( LISTLB .NE. 0 )  THEN
          WRITE (6,*) ' '
          WRITE (6,*) 'LIST      ELEMENT BANDWIDTH'
          CALL  BANWTV (NE, N, NTOLPT, NODES, LNODE,
     1                  LBAND, M, NNOW )
          CALL  WRTOUT ( IONE, IONE, NE, LBAND)
        ENDIF
C        PRINT BEST NODE ORDER
        CALL N2WAS (M, NNOW, NWAS)
        IF ( LNEWN .GT. 0 )  THEN
          WRITE (6,*) ' '
          WRITE (6,*) '     NEW NODE NUMBER DATA:'
          CALL  WASNOW (M, NNOW, NWAS)
        ENDIF
C        SAVE OLD ELEMS WITH NEW NODES
        IOPT = 2
        CALL  NEWDAT (M, NE, N, NNMAX, IOPT, NTOLPT, NODES, 
     1                NNOW, LWAS, LNODE, INHEAD, LHEAD)
      ENDIF
      IF ( LNOPT .EQ. 0 )  THEN
C       CHECK ELEMENT NEIGHBORS LIST
        LLUNIT = 14
        CALL  GETL2L (NE, IDUMMY, LTOLPT, ID(L21), LLUNIT)
        IF ( LISTLL .NE. 0 )  THEN
          WRITE (6,*) 'RE-ORDERED'
          WRITE (6,*) 'ELEMENT   ADJACENT ELEMENTS'
          CALL  WRTOUT (NE, LTOLPT, IDUMMY, ID(L21))
        ENDIF
      ENDIF
      RETURN
      END
      SUBROUTINE  WASNOW (M, ISNOW, ISWAS)
C     * * * * * * * * * * * * * * * * * * * * * * * * *
C     PRINT OLD AND NEW NUMBERS
C     * * * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  ISNOW, ISWAS
      DIMENSION  ISNOW(M), ISWAS(M)
C     IF I = NEW NUMBER, ISWAS(I) = OLD NUMBER
C     IF I = OLD NUMBER, ISNOW(I) = NEW NUMBER
      WRITE (6,*) '         RESEQUENCED NUMBERS '
      WRITE (6,*) '    WAS      IS        IS     WAS '
      WRITE (6,5000) (K, ISNOW(K), K, ISWAS(K), K=1,M)
 5000 FORMAT ( 2I8, I10, I8)
      RETURN
      END
      SUBROUTINE  WHERE (NSIZE, IT, LOC, LIST)
C     * * * * * * * * * * * * * * * * * * * * * * *
C     DETERMINE WHERE 'IT' BELONGS IN LIST
C     * * * * * * * * * * * * * * * * * * * * * * *
C I2  INTEGER*2  LIST
      DIMENSION  LIST(1)
C     IT    = A NUMBER THAT BELONGS IN LIST
C     LIST  = AN ARRAY OF INCREASING NUMBERS
C     LOC   = LOCATION WHERE IT SHOULD BE INSERTED
C           = -J IF IT IS THE J-TH COEFF. IN LIST
C     NSIZE = CURRENT SIZE OF LIST
      IF ( NSIZE .LT. 1 )  GO TO 10
C      SHOULD IT BE FIRST?
      LTEST = LIST(1)
      IF ( IT - IABS(LTEST) )  10,20,30
   10 LOC = 1
      RETURN
   20 LOC = -1
      RETURN
C      SHOULD IT BE LAST?
   30 CONTINUE
      LTEST = LIST(NSIZE)
      IF ( IABS(LTEST) - IT )  40,50,60
   40 LOC = NSIZE + 1
      RETURN
   50 LOC = -NSIZE
      RETURN
C      DETERMINE WHICH PART OF LIST TO SEARCH
   60 ISTART = 1
      ISTOP = NSIZE
      DO 70  I = 2,NSIZE
        IMID = ( ISTART + ISTOP )/2
        LTEST = LIST(IMID)
        LTEST = IABS(LTEST)
        IF ( IT .GE. LTEST )  ISTART = IMID
        IF ( IT .LE. LTEST )  ISTOP = IMID
        IF ( IT .EQ. LTEST )  GO TO 80
          IF ( ISTOP .LE. (ISTART + 2) )  GO TO 80
   70 CONTINUE
   80 DO 110  I = ISTART,ISTOP
        LTEST = LIST(I)
        IF ( IT - IABS(LTEST) )  90,100,110
   90   LOC = I
        RETURN
  100   LOC = -I
        RETURN
  110 CONTINUE
      RETURN
      END
      SUBROUTINE  WRTOUT (NO, IPOINT, NMAX, LIST)
C     * * * * * * * * * * * * * * * * * * * * * * * *
C     PRINT COMPRESSED 'LIST' USING POINTER 'IPOINT'
C     * * * * * * * * * * * * * * * * * * * * * * * * 
C I2  INTEGER*2  IPOINT, LIST
      DIMENSION  IPOINT(NO), LIST(NMAX)
C     IPOINT  = POINTER ARRAY FOR COMPACT LIST
C     LIST    = ARRAY STORED IN VECTOR FORM
C     NMAX    = CURRENT SIZE OF THE LIST
C     NO      = NUMBER OF SUBARRAYS
      WRITE (6,*) 'NUMBER    ARRAY'
      IF ( NO .EQ. 1 )  GO TO 20
      NOL1 = NO - 1
      DO 10  I = 1,NOL1
        I1 = IPOINT(I)
        I2 = IPOINT(I+1) - 1
        IF ( I1 .EQ. (I2+1) )  GO TO 10
          IF ( I2 .GT. NMAX   )  WRITE (12,*)
     1    'WARNING: DIMENSION EXCEEDED IN WRTOUT'
          WRITE (6,5010)  I, (LIST(K), K=I1,I2)
 5010     FORMAT ( I5, 16I7, /, ( 5X, 16I7 ) )
   10 CONTINUE
   20 I1 = IPOINT(NO)
      IF ( I1 .GT. NMAX )  WRITE (12,*)
     1     'WARNING: DIMENSION EXCEEDED IN WRTOUT.'
      I2 = NMAX
      WRITE (6,5010)  NO, (LIST(K), K=I1,I2)
      RETURN
      END
