This routine calculates the determinate of a square matrix. Gauss Elimination Method array the matrix of order norder which is to be evaluated. this subprogram destroys the matrix array norder the order of the square matrix to be evaluated.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=dp), | intent(inout), | dimension(n,n) | :: | array | ||
integer, | intent(in) | :: | n |
function comp_det(array, n) result(det) use precision, only: dp implicit none real(kind=dp) :: det real(kind=dp), intent(inout), dimension(n,n) :: array integer, intent(in) :: n real(kind=dp), dimension(n,n) :: work integer, dimension(n) :: num real(kind=dp) :: tmp, max integer i, k, l, m det = 1.0_dp do k = 1, n max = array(k,k) num(k) = k do i = k+1, n if(abs(max)<abs(array(i,k))) then max = array(i,k) num(k) = i end if end do if (num(k)/=k) then do l = k, n tmp = array(k,l) array(k,l) = array(num(k),l) array(num(k),l) = tmp end do det = -1.0_dp*det end if do m = k+1, n work(m,k) = array(m,k)/array(k,k) do l = k, n array(m,l) = array(m,l)-work(m,k)*array(k,l) end do end do !There we made matrix triangular! end do do i = 1, n det = det*array(i,i) end do end function comp_det