subroutine Matrix_inverse(A,nlinA,ncolA) implicit none integer i,j,k,m,n,o,lworkmax integer :: nlinA,ncolA integer :: lwork,info,lda integer, allocatable :: ipiv(:) double precision, allocatable :: work(:) double precision A(nlinA,ncolA) m=nlinA n=ncolA lda=ncolA lwork = 4*nlinA allocate(ipiv(n)) allocate(work(lwork)) ! calculating the Lapack libraries ! LU factorization call DGETRF(m, n, A, lda, ipiv, info) ! inverse call DGETRI(m, A, n, ipiv, work, lwork, info) return end