SUBROUTINE EQLIN90(A,X,N,EPS,INV,IFLAG)

!      -- se INV=0  risolve solo il sistema lineare
!      -- se INV=1  risolve il sistema lineare e calcola l`inversa
!      -- se INV=2  calcola solo la matrice inversa dei coefficienti.

!       - n.b.: il determinante viene calcolato sempre . -

implicit none

INTEGER, PARAMETER :: iknd=SELECTED_INT_KIND(8)  ! nobody does this
INTEGER, PARAMETER :: DP=SELECTED_REAL_KIND(14)  ! this is double precision
INTEGER, PARAMETER :: SP=SELECTED_REAL_KIND(6)  ! this is single precision

INTEGER CONT,save
real(DP) :: n
real(DP) :: A(20,20),IR(20),IC(20),IS(20),TN(20),X(N),S(20)
real(DP) :: iflag,inv,deter,pivot,l,l1,l2,eps,kk,jj,nm1,k1
integer :: i,k,j

IFLAG = 0
IF(INV /= 2) THEN
DO 97 I = 1,N
97    TN(I) = X(I)
ENDIF
98  DETER = 1.

!.... Ricerca del pivot per la iterazione K

DO 16 K = 1,N
PIVOT = 0.
DO 6 I = 1,N
DO 6 J = 1,N
  IF(K-1) 3,1,3
1     IF(ABS(A(I,J))-ABS(PIVOT)) 6,6,2
2     PIVOT = A(I,J)
  IR(K) = I
  IC(K) = J
  GOTO 6
3      L = K-1
  DO 5 L1 = 1,L
  DO 5 L2 = 1,L
    IF(I == IR(L1) .OR. J == IC(L2)) GOTO 6
5      CONTINUE
GOTO 1
6    CONTINUE

!.... Controllo della singolarita

IF(ABS(PIVOT)-EPS) 7,7,8
7    DETER = 0.
IFLAG = 1
RETURN

!.... Aggiornamento del valore del determinante

8    DETER = DETER*PIVOT
KK = IC(K)
JJ = IR(K)

!.... Normalizzazione della riga cui appartiene il pivot

DO 11 J = 1,N
  IF(J-KK) 10,9,10
9      A(JJ,KK) = 1./A(JJ,KK)

  GOTO 11
10      A(JJ,J) = A(JJ,J)/PIVOT
11    CONTINUE
IF(INV /= 2) TN(JJ) = TN(JJ)/PIVOT

!.... Riduzione delle altre ( N-1 ) equazioni del sistema

DO 15 I = 1,N
  IF(I-JJ) 12,15,12
12      DO 14 J = 1,N
    IF(J /= KK) A(I,J) = A(I,J)-A(I,KK)*A(JJ,J)
14      CONTINUE
  IF(INV /= 2) TN(I) = TN(I)-A(I,KK)*TN(JJ)
  A(I,KK) = -A(I,KK)/PIVOT
15    CONTINUE
16  CONTINUE

!.... Ordinamento del vettore soluzione

DO 18 K = 1,N
KK = IC(K)
JJ = IR(K)
IS(KK) = IR(K)
IF(INV /= 2) X(KK) = TN(JJ)
18  CONTINUE
CONT = 0
NM1 = N-1

!.... Calcolo del numero di scambi fra le righe

DO 20 K = 1,NM1
DO 20 K1 = 2,N
IF(IS(K1)-IS(K)) 19,20,20
19    CONT = CONT+1
SAVE = IS(K1)
IS(K1) = IS(K)
IS(K) = SAVE
20  CONTINUE

!....Controllo sul segno del determinante

IF(CONT /= (CONT/2)*2) DETER = -DETER
IF(INV == 0) RETURN

!.... Scambi e calcolo della matrice inversa

DO 24 J = 1,N
DO 23 K = 1,N
  KK = IC(K)
  JJ = IR(K)
  S(KK) = A(JJ,J)
23    CONTINUE
DO 24 K = 1,N
  A(K,J) = S(K)
24  CONTINUE
DO 26 I = 1,N
DO 25 K = 1,N
  JJ = IR(K)
  KK = IC(K)
 S(JJ) = A(I,KK)
25    CONTINUE
DO 26 K = 1,N
  A(I,K) = S(K)
26  CONTINUE

RETURN
END subroutine EQLIN90