        SUBROUTINE  TDRV
     *     (N, R,IC, IA,JA,A, B, Z, NSP,ISP,RSP,ESP, FLAG)                TD  12
C                                                                         TD  13
C    PARAMETERS                                                           TD  14
C    CLASS ABBREVIATIONS ARE --                                           TD  15
C       N - INTEGER VARIABLE                                              TD  16
C       F - REAL VARIABLE                                                 TD  17
C       V - SUPPLIES A VALUE TO THE DRIVER                                TD  18
C       R - RETURNS A RESULT FROM THE DRIVER                              TD  19
C       I - USED INTERNALLY BY THE DRIVER                                 TD  20
C       A - ARRAY                                                         TD  21
C                                                                         TD  22
C CLASS \ PARAMETER                                                       TD  23
C ------+----------                                                       TD  24
C       \                                                                 TD  25
C         THE NONZERO ENTRIES OF THE COEFFICIENT MATRIX M ARE STORED      TD  26
C    ROW-BY-ROW IN THE ARRAY A.  TO IDENTIFY THE INDIVIDUAL NONZERO       TD  27
C    ENTRIES IN EACH ROW, WE NEED TO KNOW IN WHICH COLUMN EACH ENTRY      TD  28
C    LIES.  THE COLUMN INDICES WHICH CORRESPOND TO THE NONZERO ENTRIES    TD  29
C    OF M ARE STORED IN THE ARRAY JA;  I.E., IF  A(K) = M(I,J),  THEN     TD  30
C    JA(K) = J.  IN ADDITION, WE NEED TO KNOW WHERE EACH ROW STARTS AND   TD  31
C    HOW LONG IT IS.  THE INDEX POSITIONS IN JA AND A WHERE THE ROWS OF   TD  32
C    M BEGIN ARE STORED IN THE ARRAY IA;  I.E., IF M(I,J) IS THE FIRST    TD  33
C    NONZERO ENTRY (STORED) IN THE I-TH ROW AND A(K) = M(I,J),  THEN      TD  34
C    IA(I) = K.  MOREOVER, THE INDEX IN JA AND A OF THE FIRST LOCATION    TD  35
C    FOLLOWING THE LAST ELEMENT IN THE LAST ROW IS STORED IN IA(N+1).     TD  36
C    THUS, THE NUMBER OF ENTRIES IN THE I-TH ROW IS GIVEN BY              TD  37
C    IA(I+1) - IA(I),  THE NONZERO ENTRIES OF THE I-TH ROW ARE STORED     TD  38
C    CONSECUTIVELY IN                                                     TD  39
C            A(IA(I)),  A(IA(I)+1),  ..., A(IA(I+1)-1),                   TD  40
C    AND THE CORRESPONDING COLUMN INDICES ARE STORED CONSECUTIVELY IN     TD  41
C            JA(IA(I)), JA(IA(I)+1), ..., JA(IA(I+1)-1).                  TD  42
C    FOR EXAMPLE, THE 5 BY 5 MATRIX                                       TD  43
C                ( 1. 0. 2. 0. 0.)                                        TD  44
C                ( 0. 3. 0. 0. 0.)                                        TD  45
C            M = ( 0. 4. 5. 6. 0.)                                        TD  46
C                ( 0. 0. 0. 7. 0.)                                        TD  47
C                ( 0. 0. 0. 8. 9.)                                        TD  48
C    WOULD BE STORED AS                                                   TD  49
C               \ 1  2  3  4  5  6  7  8  9                               TD  50
C            ---+--------------------------                               TD  51
C            IA \ 1  3  4  7  8 10                                        TD  52
C            JA \ 1  3  2  2  3  4  4  4  5                               TD  53
C             A \ 1. 2. 3. 4. 5. 6. 7. 8. 9.         .                    TD  54
C                                                                         TD  55
C NV    \ N     - NUMBER OF VARIABLES/EQUATIONS.                          TD  56
C FVA   \ A     - NONZERO ENTRIES OF THE COEFFICIENT MATRIX M, STORED     TD  57
C       \           BY ROWS.                                              TD  58
C       \           SIZE = NUMBER OF NONZERO ENTRIES IN M.                TD  59
C NVA   \ IA    - POINTERS TO DELIMIT THE ROWS IN A.                      TD  60
C       \           SIZE = N+1.                                           TD  61
C NVA   \ JA    - COLUMN NUMBERS CORRESPONDING TO THE ELEMENTS OF A.      TD  62
C       \           SIZE = SIZE OF A.                                     TD  63
C FVA   \ B     - RIGHT-HAND SIDE B;  B AND Z CAN THE SAME ARRAY.         TD  64
C       \           SIZE = N.                                             TD  65
C FRA   \ Z     - SOLUTION X;  B AND Z CAN BE THE SAME ARRAY.             TD  66
C       \           SIZE = N.                                             TD  67
C                                                                         TD  68
C         THE ROWS AND COLUMNS OF THE ORIGINAL MATRIX M CAN BE            TD  69
C    REORDERED (E.G., TO REDUCE FILLIN OR ENSURE NUMERICAL STABILITY)     TD  70
C    BEFORE CALLING THE DRIVER.  IF NO REORDERING IS DONE, THEN SET       TD  71
C    R(I) = C(I) = IC(I) = I  FOR I=1,...,N.  THE SOLUTION Z IS RETURNED  TD  72
C    IN THE ORIGINAL ORDER.                                               TD  73
C                                                                         TD  74
C NVA   \ R     - ORDERING OF THE ROWS OF M.                              TD  75
C       \           SIZE = N.                                             TD  76
C NVA   \ IC    - INVERSE OF THE ORDERING OF THE COLUMNS OF M;  I.E.,     TD  77
C       \           IC(C(I)) = I  FOR I=1,...,N,  WHERE C IS THE          TD  78
C       \           ORDERING OF THE COLUMNS OF M.                         TD  79
C       \           SIZE = N.                                             TD  80
C                                                                         TD  81
C         VARIOUS ERRORS ARE DETECTED BY THE DRIVER AND THE INDIVIDUAL    TD  82
C    SUBROUTINES.                                                         TD  83
C                                                                         TD  84
C NR    \ FLAG  - ERROR FLAG;  VALUES AND THEIR MEANINGS ARE --           TD  85
C       \             0     NO ERRORS DETECTED                            TD  86
C       \             N+K   NULL ROW IN  A  --  ROW = K                   TD  87
C       \            2N+K   DUPLICATE ENTRY IN A  --  ROW = K             TD  88
C       \            5N+K   NULL PIVOT  --  ROW = K                       TD  89
C       \            8N+K   ZERO PIVOT  --  ROW = K                       TD  90
C       \           10N+1   INSUFFICIENT STORAGE IN TDRV                  TD  91
C       \           12N+K   INSUFFICIENT STORAGE IN TRK                   TD  92
C                                                                         TD  93
C         WORKING STORAGE IS NEEDED FOR THE FACTORED FORM OF THE MATRIX   TD  94
C    M PLUS VARIOUS TEMPORARY VECTORS.  THE ARRAYS ISP AND RSP SHOULD BE  TD  95
C    EQUIVALENCED;  INTEGER STORAGE IS ALLOCATED FROM THE BEGINNING OF    TD  96
C    ISP AND REAL STORAGE FROM THE END OF RSP.                            TD  97
C                                                                         TD  98
C NV    \ NSP   - DECLARED DIMENSION OF RSP;  NSP GENERALLY MUST BE       TD  99
C       \           LARGER THAN 6N+2 + 2*K (WHERE K = (NUMBER OF NONZERO  TD 100
C       \           ENTRIES IN THE UPPER TRIANGLE OF M)).                 TD 101
C NVIRA \ ISP   - INTEGER WORKING STORAGE DIVIDED UP INTO VARIOUS ARRAYS  TD 102
C       \           NEEDED BY THE SUBROUTINES;  ISP AND RSP SHOULD BE     TD 103
C       \           EQUIVALENCED.                                         TD 104
C       \           SIZE = LRATIO*NSP                                     TD 105
C FVIRA \ RSP   - REAL WORKING STORAGE DIVIDED UP INTO VARIOUS ARRAYS     TD 106
C       \           NEEDED BY THE SUBROUTINES;  ISP AND RSP SHOULD BE     TD 107
C       \           EQUIVALENCED.                                         TD 108
C       \           SIZE = NSP.                                           TD 109
C NR    \ ESP   - IF NSP IS SUFFICIENTLY LARGE TO ALLOCATE SPACE, THEN    TD 110
C       \           ESP IS SET TO THE AMOUNT OF EXCESS STORAGE PROVIDED.  TD 111
C                                                                         TD 112
C                                                                         TD 113
C  CONVERSION TO DOUBLE PRECISION                                         TD 114
C                                                                         TD 115
C    TO CONVERT THESE ROUTINES FOR DOUBLE PRECISION ARRAYS, SIMPLY USE    TD 116
C    THE DOUBLE PRECISION DECLARATIONS IN PLACE OF THE REAL DECLARATIONS  TD 117
C    IN EACH SUBPROGRAM;  IN ADDITION, THE DATA VALUE OF THE INTEGER      TD 118
C    VARIABLE LRATIO MUST BE SET AS INDICATED IN SUBROUTINE TDRV          TD 119
C                                                                         TD 120
        INTEGER  R(1), IC(1),  IA(1), JA(1),  ISP(1), ESP,  FLAG,         TD 121
     *     U, ROW, TMP, Q                                                 TD 122
C       REAL  A(1),  B(1),  Z(1),  RSP(1)                                 TD 123
        DOUBLE PRECISION  A(1),  B(1),  Z(1),  RSP(1)                     TD 124
C                                                                         TD 125
C  SET LRATIO EQUAL TO THE RATIO BETWEEN THE LENGTH OF FLOATING POINT     TD 126
C  AND INTEGER ARRAY DATA;  E. G., LRATIO = 1 FOR (REAL, INTEGER),        TD 127
C  LRATIO = 2 FOR (DOUBLE PRECISION, INTEGER)                             TD 128
C                                                                         TD 129
        DATA LRATIO/2/                                                    TD 130
C                                                                         TD 131
C  ******  INITIALIZE AND DIVIDE UP TEMPORARY STORAGE  *****************  TD 132
        IJU = 1                                                           TD 133
        IU  = IJU +  N                                                    TD 134
        Q   = IU  + N+1                                                   TD 135
        IM  = Q   + N+1                                                   TD 136
        U   = (IM + N - 2 + LRATIO) / LRATIO   +   1                      TD 137
        JU  = LRATIO * (U - 1) + 1                                        TD 138
        ROW = NSP + 1 -  N                                                TD 139
        TMP = ROW -  N                                                    TD 140
        MAX = TMP - U                                                     TD 141
        ESP = MAX                                                         TD 142
C                                                                         TD 143
C  ******  CALL ZERO-TRACKING SUBROUTINE  ******************************  TD 144
        IF (MAX.LT.0)  GO TO 110                                          TD 145
        CALL  TRK                                                         TD 146
     *     (N,  R, IC,  IA, JA, A,  Z,  B,                                TD 147
     *      ISP(IJU), ISP(JU), ISP(IU), RSP(U), MAX,                      TD 148
     *      ISP(Q),  ISP(IM),  RSP(ROW),  RSP(TMP), FLAG, ESP, LRATIO)    TD 149
        IF (FLAG.NE.0)  GO TO 100                                         TD 150
        RETURN                                                            TD 151
C                                                                         TD 152
C ** ERROR:  ERROR DETECTED IN TRK                                        TD 153
 100    RETURN                                                            TD 154
C ** ERROR:  INSUFFICIENT STORAGE                                         TD 155
 110    FLAG = 10*N + 1                                                   TD 156
        RETURN                                                            TD 157
        END                                                               TD 158
