MODULE LU

  CONTAINS

    SUBROUTINE lu_fact(mat)

      IMPLICIT NONE

      REAL,INTENT(INOUT) :: mat(:,:)
      INTEGER :: i, n

      n=SIZE(mat,1)
      DO i=2,n
         mat(i,1)=mat(i,1)/mat(i-1,2)
         mat(i,2)=mat(i,2)-mat(i,1)*mat(i-1,3)
      ENDDO

    END SUBROUTINE lu_fact

    SUBROUTINE lu_solve(mat,g)

      IMPLICIT NONE

      REAL, INTENT(IN) :: mat(:,:)
      REAL, INTENT(INOUT) :: g(:)
      INTEGER :: n,i 

      n=SIZE(g,1)
      DO i=2,n
         g(i)=g(i)-mat(i,1)*g(i-1)
      ENDDO
      g(n)=g(n)/mat(n,2)
      DO i=n-1,1,-1
         g(i)=(g(i)-mat(i,3)*g(i+1))/mat(i,2)
      ENDDO

    END SUBROUTINE lu_solve

END MODULE LU
