C
C     ..................................................................
C
C        SUBROUTINE RSLMC
C
C        PURPOSE
C           SOLUTION OF A SYSTEM OF LINEAR EQUATIONS AX=B
C
C        USAGE
C           CALL RSLMC(A,AF,B,X,N,EPSI,IER,IA,V,PER)
C
C        DESCRIPTION OF PARAMETERS
C           A      INPUT MATRIX
C           AF     ARRAY OF THE FACTORIZATION OF THE ORIGINAL MATRIX
C           B      RIGHT HAND SIDE VECTOR
C           X      VECTOR CONTAINING THE SOLUTION ON RETURN
C           N      ORDER OF THE SYSTEM
C           EPSI   RELATIVE PRECISION INDICATOR(REQUIRED INPUT)
C           IER    ERROR INDICATOR
C                     =0 IF EACH COMPONENT OF X MEETS THE PRECISION EPSI
C                     =1 IF ONLY THE NORM OF X MEETS THIS PRECISION
C                     =2 IF THE PRECISION IN THE NORM OF THE COMPUTED
C                        SOLUTION IS LOWER THAN EPSI
C                     =3 IF THE SOLUTION OBTAINED HAS NO MEANING AT ALL
C                     =4 IF A DIAGONAL TERM OF THE UPPER TRIANGULAR
C                        FACTOR IS ZERO
C           IA     SIZE OF THE FIRST DIMENSION ASSIGNED TO THE ARRAY A
C                  IN THE CALLING PROGRAM WHEN THE MATRIX IS IN DOUBLE
C                  SUBSCRIPTED DATA STORAGE MODE.  IA=N WHEN THE MATRIX
C                  IS IN SSP VECTOR STORAGE MODE.
C           V      WORKING STORAGE VECTOR
C                  DIMENSION OF V MUST BE GREATER THAN OR EQUAL TO N
C           PER    VECTOR WHERE PERMUTATIONS OF ROWS OF THE MATRIX ARE
C                  STORED
C                  DIMENSION OF PER MUST BE GREATER THAN OR EQUAL TO N
C
C        REMARKS
C           THE MATRIX OF THE SYSTEM MAY BE FACTORIZED BY THE SUBROUTINE
C           FACTR IN THE ARRAY AF PRIOR TO ENTRY TO THIS SUBROUTINE.
C           THE LOWER TRIANGULAR FACTOR MUST HAVE AN UNIT DIAGONAL.
C           EPSI IS MODIFIED WHEN IER=2
C
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
C           NONE
C
C        METHOD
C           A TRIAL SOLUTION IS FIRST COMPUTED.  THEN CORRECTIONS ARE
C           CALCULATED FROM RESIDUAL VECTORS.
C
C        REFERENCES
C           J. H. WILKINSON - THE ALGEBRAIC EIGENVALUE PROBLEM -
C           CLARENDON PRESS OXFORD, 1965.  H. J. BOWDLER, R. S. MARTIN,
C           G. PETERS, AND J. H. WILKINSON.  'SOLUTION OF REAL AND
C           COMPLEX SYSTEMS OF LINEAR EQUATIONS', NUMERISCHE MATHEMATIK,
C           VOL. 8, NO. 3, 1966, 217-234.
C
C     ..................................................................
C
      SUBROUTINE RSLMC (A,AF,B,X,N,EPSI,IER,IA,V,PER)
      DIMENSION A(1),AF(1),B(1),X(1),V(1),PER(1)
      DOUBLE PRECISION DP
C
C        INITIALIZATION
C
      D0=0.
      IER=0
      ITE=0
      DO 10 I=1,N
      V(I)=B(I)
   10 X(I)=0.
   20 ITE=ITE+1
C
C        THE PERMUTATIONS OF ROWS OF A ARE APPLIED TO V
C
      DO 30 I=1,N
      K=PER(I)
      IF (K-I)25,30,25
   25 D1=V(K)
      V(K)=V(I)
      V(I)=D1
   30 CONTINUE
C
C        SOLUTION OF THE LOWER TRIANGULAR SYSTEM
C
      DO 50 I=2,N
      IM1=I-1
      DP=V(I)
      IK=I
      DO 40 K=1,IM1
      DP=DP-1.D0*AF(IK)*V(K)
   40 IK=IK+IA
   50 V(I)=DP
C
C        SOLUTION OF THE UPPER TRIANGULAR SYSTEM
C
      IF(AF(IK)) 58,54,58
   54 IER=4
      GO TO 82
   58 V(N)=DP/AF(IK)
      DO 70 I=2,N
      IM1=N-I+1
      INF=IM1+1
      DP=V(IM1)
      IK=(IM1-1)*IA+IM1
      D1=AF(IK)
      DO 60 K=INF,N
      IK=IK+IA
   60 DP=DP-1.D0*AF(IK)*V(K)
   70 V(IM1)=DP/D1
C
C        TEST OF PRECISION
C
      D1=0.
      D2=0.
      KLE=0
      DO 80 I=1,N
      D1=D1+ABS(V(I))
      D2=D2+ABS(X(I))
      IF (ABS(V(I))-EPSI*ABS(X(I))) 80,80,75
   75 KLE=1
   80 CONTINUE
      IF (KLE)140,82,85
   82 RETURN
   85 IF (ITE-1)140,90,87
C
C        ITERATIONS ARE STOPPED WHEN THE NORM OF THE CORRECTION IS MORE
C        THAN HALF OF THE ONE OF THE FORMER
C
   87 IF (D0-2.*D1)120,90,90
   90 DO 95 I=1,N
   95 X(I)=X(I)+V(I)
      DO 110 I=1,N
      DP=B(I)
      IK=I
      DO 100 K=1,N
      DP=DP-1.D0*A(IK)*X(K)
  100 IK=IK+IA
  110 V(I)=DP
      D0=D1
      GO TO 20
  120 IF(ITE-2)140,140,125
  125 IF (D1-EPSI*D2)127,127,130
  127 IER=1
      RETURN
  130 IER=2
      EPSI=D1/D2
      RETURN
  140 IER=3
      RETURN
      END
