PROGRAM GassJorn
implicit none
      INTEGER LIMROW, LIMCOL

      PARAMETER (LIMROW = 10, LIMCOL = LIMROW + 1)

      DOUBLE PRECISION LIN(LIMROW, LIMCOL), X(LIMROW)

      INTEGER N, I, J

      LOGICAL SINGUL

      write(*,*)"Escribe el numero de ecuaciones"
      READ(*,*)N
      DO 10 I = 1, N
      write(*,*)"Escribe los coeficientes y la constante d la ecuacion",+  I, ': ' 
         READ (*,*) (LIN(I,J), J = 1, N + 1)
 10 CONTINUE

      CALL GAUSS(LIN, LIMROW, LIMCOL, N, X, SINGUL)
      IF (.NOT. SINGUL) THEN
        
		open(unit=31, file="solucngj.salida", status="unknown")
		 
		!DO 31 I = 1, N
          !write(31,100) I, X(I)
!100 FORMAT(1X, 'X(', I2, ') =', F8.3)
         !write(31,*) "La matriz es singular"
		
		 write(*,*)"La solucion es"
write(31,*) "La solucion es"
         DO 20 I = 1, N
            write(*,100) I, X(I)
write(31,100) I, X(I)
100 FORMAT(1X, 'X(', I2, ') =', F8.3)
20  CONTINUE
      ELSE
		write(31,*) "La matriz es singular"
         write(*,*) "La matriz es singular"

		
      END IF

       END program GassJorn
 

      SUBROUTINE GAUSS(LIN, LIMROW, LIMCOL, N, X, SINGUL)
      DOUBLE PRECISION LIN(LIMROW, LIMCOL), X(LIMROW), TEMP, MULT, EPSIL
      PARAMETER (EPSIL = 1D-15)
      INTEGER N, PIVROW
      LOGICAL SINGUL
      SINGUL = .FALSE.

      DO 50 I = 1, N

         ABSPIV = ABS(LIN(I,I))
         PIVROW = I
         DO 10 K = I + 1, N
            IF (ABS(LIN(K,I)) .GT. ABSPIV) THEN
               ABSPIV = ABS(LIN(K,I))
               PIVROW = K
            END IF
10 CONTINUE
         IF (ABSPIV .LT. EPSIL) THEN
            SINGUL = .TRUE.
            RETURN
         END IF
         IF (PIVROW .NE. I) THEN

            DO 20 J = 1, N + 1

               TEMP = LIN(I,J)

               LIN(I,J) = LIN(PIVROW,J)

               LIN(PIVROW,J) = TEMP

20 CONTINUE

         END IF
     DO 40 J = I + 1, N

            MULT = -LIN(J,I) / LIN(I,I)

            DO 30 K = I, N + 1

               LIN(J,K) = LIN(J,K) +  MULT * LIN(I,K)

30          CONTINUE

40       CONTINUE

 

50    CONTINUE

      X(N) = LIN(N, N + 1) / LIN(N,N)

      DO 70 J = N - 1, 1, -1

         X(J) = LIN(J, N + 1)

         DO 60 K = J + 1, N

            X(J) = X(J) - LIN(J,K) * X(K)

60       CONTINUE

         X(J) = X(J) / LIN(J,J)

70    CONTINUE

      END

