gauss.f90

subroutine gauss(a, b)
!Eliminación de Gauss
!Entrada:	a (matriz de coeficientes)
!		b (términos independientes)
!Salida:	a (matriz triangular superior salvo permutación de filas)
!		b (términos independientes permutados)

implicit none
real(8), intent (inout) :: a(:,:)
real(8), intent (inout) :: b(:)
real(8) :: det, piv
integer :: i, j, k, n

n = size(b)

!Inicialización del determinante
det = 1d0
!Etapa k-ésima de la eliminación
do k = 1, n-1
	piv = a(k,k)
	!Comprobación de que el k-ésimo pivote no es nulo
	if (dabs(piv) < 1.e-12) then
		print*,"Pivote nulo en la etapa: ", k
		stop
	end if
	!Actualización del determinante
	det = det * piv
	!Eliminación
	a(k+1:n, k) = a(k+1:n, k)/piv
	do j = k+1, n
		a(k+1:n, j) = a(k+1:n, j)-a(k+1:n, k)*a(k, j)
	end do
b(k+1:n) = b(k+1:n)-a(k+1:n, k)*b(k)
end do
!Comprobación de que el último pivote no es nulo
if (dabs(a(n, n)) < 1.e-12) then
	print*,"Pivote nulo en la etapa: ", n
	stop
end if
!Fin del cálculo del determinante
det = detgauss*a(n, n)
print*,"Determinante: ", det

return
end subroutine 
gtyh3td6