      SUBROUTINE ODRV
     *     (N, IA,JA,A, P,IP, NSP,ISP, PATH, FLAG)
CLLL. OPTIMIZE
C                                                                 5/2/83
C***********************************************************************
C  ODRV -- DRIVER FOR SPARSE MATRIX REORDERING ROUTINES
C***********************************************************************
C
C  DESCRIPTION
C
C    ODRV FINDS A MINIMUM DEGREE ORDERING OF THE ROWS AND COLUMNS
C    OF A MATRIX M STORED IN (IA,JA,A) FORMAT (SEE BELOW).  FOR THE
C    REORDERED MATRIX, THE WORK AND STORAGE REQUIRED TO PERFORM
C    GAUSSIAN ELIMINATION IS (USUALLY) SIGNIFICANTLY LESS.
C
C    NOTE.. ODRV AND ITS SUBORDINATE ROUTINES HAVE BEEN MODIFIED TO
C    COMPUTE ORDERINGS FOR GENERAL MATRICES, NOT NECESSARILY HAVING ANY
C    SYMMETRY.  THE MIMINUM DEGREE ORDERING IS COMPUTED FOR THE
C    STRUCTURE OF THE SYMMETRIC MATRIX  M + M-TRANSPOSE.
C    MODIFICATIONS TO THE ORIGINAL ODRV MODULE HAVE BEEN MADE IN
C    THE CODING IN SUBROUTINE MDI, AND IN THE INITIAL COMMENTS IN
C    SUBROUTINES ODRV AND MD. 
C
C    IF ONLY THE NONZERO ENTRIES IN THE UPPER TRIANGLE OF M ARE BEING 
C    STORED, THEN ODRV SYMMETRICALLY REORDERS (IA,JA,A), (OPTIONALLY) 
C    WITH THE DIAGONAL ENTRIES PLACED FIRST IN EACH ROW.  THIS IS TO
C    ENSURE THAT IF M(I,J) WILL BE IN THE UPPER TRIANGLE OF M WITH
C    RESPECT TO THE NEW ORDERING, THEN M(I,J) IS STORED IN ROW I (AND 
C    THUS M(J,I) IS NOT STORED),  WHEREAS IF M(I,J) WILL BE IN THE
C    STRICT LOWER TRIANGLE OF M, THEN M(J,I) IS STORED IN ROW J (AND
C    THUS M(I,J) IS NOT STORED).
C
C
C  STORAGE OF SPARSE MATRICES 
C
C    THE NONZERO ENTRIES OF THE MATRIX M ARE STORED ROW-BY-ROW IN THE 
C    ARRAY A.  TO IDENTIFY THE INDIVIDUAL NONZERO ENTRIES IN EACH ROW,
C    WE NEED TO KNOW IN WHICH COLUMN EACH ENTRY LIES.  THESE COLUMN
C    INDICES ARE STORED IN THE ARRAY JA.  I.E., IF  A(K) = M(I,J),  THEN
C    JA(K) = J.  TO IDENTIFY THE INDIVIDUAL ROWS, WE NEED TO KNOW WHERE
C    EACH ROW STARTS.  THESE ROW POINTERS ARE STORED IN THE ARRAY IA. 
C    I.E., IF M(I,J) IS THE FIRST NONZERO ENTRY (STORED) IN THE I-TH ROW
C    AND  A(K) = M(I,J),  THEN  IA(I) = K.  MOREOVER, IA(N+1) POINTS TO
C    THE FIRST LOCATION FOLLOWING THE LAST ELEMENT IN THE LAST ROW.
C    THUS, THE NUMBER OF ENTRIES IN THE I-TH ROW IS  IA(I+1) - IA(I), 
C    THE NONZERO ENTRIES IN THE I-TH ROW ARE STORED CONSECUTIVELY IN
C
C            A(IA(I)),  A(IA(I)+1),  ..., A(IA(I+1)-1),
C
C    AND THE CORRESPONDING COLUMN INDICES ARE STORED CONSECUTIVELY IN 
C
C            JA(IA(I)), JA(IA(I)+1), ..., JA(IA(I+1)-1).
C
C    SINCE THE COEFFICIENT MATRIX IS SYMMETRIC, ONLY THE NONZERO ENTRIES
C    IN THE UPPER TRIANGLE NEED BE STORED.  FOR EXAMPLE, THE MATRIX
C
C             ( 1  0  2  3  0 )
C             ( 0  4  0  0  0 )
C         M = ( 2  0  5  6  0 )
C             ( 3  0  6  7  8 )
C             ( 0  0  0  8  9 )
C
C    COULD BE STORED AS
C
C            - 1  2  3  4  5  6  7  8  9 10 11 12 13
C         ---+--------------------------------------
C         IA - 1  4  5  8 12 14
C         JA - 1  3  4  2  1  3  4  1  3  4  5  4  5
C          A - 1  2  3  4  2  5  6  3  6  7  8  8  9
C
C    OR (SYMMETRICALLY) AS
C
C            - 1  2  3  4  5  6  7  8  9
C         ---+--------------------------
C         IA - 1  4  5  7  9 10
C         JA - 1  3  4  2  3  4  4  5  5
C          A - 1  2  3  4  5  6  7  8  9          .
C
C
C  PARAMETERS
C
C    N    - ORDER OF THE MATRIX
C
C    IA   - INTEGER ONE-DIMENSIONAL ARRAY CONTAINING POINTERS TO DELIMIT
C           ROWS IN JA AND A.  DIMENSION = N+1
C
C    JA   - INTEGER ONE-DIMENSIONAL ARRAY CONTAINING THE COLUMN INDICES
C           CORRESPONDING TO THE ELEMENTS OF A.  DIMENSION = NUMBER OF
C           NONZERO ENTRIES IN (THE UPPER TRIANGLE OF) M
C
C    A    - REAL ONE-DIMENSIONAL ARRAY CONTAINING THE NONZERO ENTRIES IN
C           (THE UPPER TRIANGLE OF) M, STORED BY ROWS.  DIMENSION =
C           NUMBER OF NONZERO ENTRIES IN (THE UPPER TRIANGLE OF) M
C
C    P    - INTEGER ONE-DIMENSIONAL ARRAY USED TO RETURN THE PERMUTATION
C           OF THE ROWS AND COLUMNS OF M CORRESPONDING TO THE MINIMUM 
C           DEGREE ORDERING.  DIMENSION = N
C
C    IP   - INTEGER ONE-DIMENSIONAL ARRAY USED TO RETURN THE INVERSE OF
C           THE PERMUTATION RETURNED IN P.  DIMENSION = N
C
C    NSP  - DECLARED DIMENSION OF THE ONE-DIMENSIONAL ARRAY ISP.  NSP 
C           MUST BE AT LEAST  3N+4K,  WHERE K IS THE NUMBER OF NONZEROES
C           IN THE STRICT UPPER TRIANGLE OF M
C
C    ISP  - INTEGER ONE-DIMENSIONAL ARRAY USED FOR WORKING STORAGE.
C           DIMENSION = NSP
C
C    PATH - INTEGER PATH SPECIFICATION.  VALUES AND THEIR MEANINGS ARE -
C             1  FIND MINIMUM DEGREE ORDERING ONLY
C             2  FIND MINIMUM DEGREE ORDERING AND REORDER SYMMETRICALLY
C                  STORED MATRIX (USED WHEN ONLY THE NONZERO ENTRIES IN
C                  THE UPPER TRIANGLE OF M ARE BEING STORED)
C             3  REORDER SYMMETRICALLY STORED MATRIX AS SPECIFIED BY
C                  INPUT PERMUTATION (USED WHEN AN ORDERING HAS ALREADY
C                  BEEN DETERMINED AND ONLY THE NONZERO ENTRIES IN THE
C                  UPPER TRIANGLE OF M ARE BEING STORED)
C             4  SAME AS 2 BUT PUT DIAGONAL ENTRIES AT START OF EACH ROW
C             5  SAME AS 3 BUT PUT DIAGONAL ENTRIES AT START OF EACH ROW
C
C    FLAG - INTEGER ERROR FLAG.  VALUES AND THEIR MEANINGS ARE -
C               0    NO ERRORS DETECTED 
C              9N+K  INSUFFICIENT STORAGE IN MD
C             10N+1  INSUFFICIENT STORAGE IN ODRV 
C             11N+1  ILLEGAL PATH SPECIFICATION
C
C
C  CONVERSION FROM REAL TO DOUBLE PRECISION
C
C    CHANGE THE REAL DECLARATIONS IN ODRV AND SRO TO DOUBLE PRECISION 
C    DECLARATIONS.
C
C-----------------------------------------------------------------------
C
      INTEGER  IA(1), JA(1),  P(1), IP(1),  ISP(1),  PATH,  FLAG,
     *   V, L, HEAD,  TMP, Q
      DOUBLE PRECISION  A(1)
      LOGICAL  DFLAG
C
C----INITIALIZE ERROR FLAG AND VALIDATE PATH SPECIFICATION
      FLAG = 0
      IF (PATH.LT.1 .OR. 5.LT.PATH)  GO TO 111
C
C----ALLOCATE STORAGE AND FIND MINIMUM DEGREE ORDERING
      IF ((PATH-1) * (PATH-2) * (PATH-4) .NE. 0)  GO TO 1
        MAX = (NSP-N)/2
        V    = 1
        L    = V     +  MAX
        HEAD = L     +  MAX
        NEXT = HEAD  +  N
        IF (MAX.LT.N)  GO TO 110
C
        CALL  MD
     *     (N, IA,JA, MAX,ISP(V),ISP(L), ISP(HEAD),P,IP, ISP(V), FLAG)
        IF (FLAG.NE.0)  GO TO 100
C
C----ALLOCATE STORAGE AND SYMMETRICALLY REORDER MATRIX
   1  IF ((PATH-2) * (PATH-3) * (PATH-4) * (PATH-5) .NE. 0)  GO TO 2
        TMP = (NSP+1) -      N
        Q   = TMP     - (IA(N+1)-1)
        IF (Q.LT.1)  GO TO 110
C
        DFLAG = PATH.EQ.4 .OR. PATH.EQ.5
        CALL SRO
     *     (N,  IP,  IA, JA, A,  ISP(TMP),  ISP(Q),  DFLAG) 
C
   2  RETURN
C
C ** ERROR -- ERROR DETECTED IN MD
 100  RETURN
C ** ERROR -- INSUFFICIENT STORAGE
 110  FLAG = 10*N + 1
      RETURN
C ** ERROR -- ILLEGAL PATH SPECIFIED
 111  FLAG = 11*N + 1
      RETURN
      END 
