From a2e70612ff0cfacc04ad33c38d93ea3b515befcb Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 10 Apr 2020 12:52:27 +0200 Subject: [PATCH] interface checking for LAPACK --- src/LAPACK_interface.f90 | 59 ++++++++++++++++++++++++++++++++++ src/commercialFEM_fileList.f90 | 1 + src/crystallite.f90 | 2 -- src/math.f90 | 11 +------ src/rotations.f90 | 2 -- 5 files changed, 61 insertions(+), 14 deletions(-) create mode 100644 src/LAPACK_interface.f90 diff --git a/src/LAPACK_interface.f90 b/src/LAPACK_interface.f90 new file mode 100644 index 000000000..7d3043ed0 --- /dev/null +++ b/src/LAPACK_interface.f90 @@ -0,0 +1,59 @@ +!-------------------------------------------------------------------------------------------------- +!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH +!> @brief Fortran interfaces for LAPACK routines +!> @details https://www.netlib.org/lapack/ +!-------------------------------------------------------------------------------------------------- +module LAPACK_interface + interface + + subroutine dgeev(jobvl,jobvr,n,a,lda,wr,wi,vl,ldvl,vr,ldvr,work,lwork,info) + use prec + character, intent(in) :: jobvl,jobvr + integer, intent(in) :: n,lda,ldvl,ldvr,lwork + real(pReal), intent(inout), dimension(lda,n) :: a + real(pReal), intent(out), dimension(n) :: wr,wi + real(pReal), intent(out), dimension(ldvl,n) :: vl + real(pReal), intent(out), dimension(ldvr,n) :: vr + real(pReal), intent(out), dimension(max(1,lwork)) :: work + integer, intent(out) :: info + end subroutine dgeev + + subroutine dgesv(n,nrhs,a,lda,ipiv,b,ldb,info) + use prec + integer, intent(in) :: n,nrhs,lda,ldb + real(pReal), intent(inout), dimension(lda,n) :: a + integer, intent(out), dimension(n) :: ipiv + real(pReal), intent(out), dimension(ldb,nrhs) :: b + integer, intent(out) :: info + end subroutine dgesv + + subroutine dgetrf(m,n,a,lda,ipiv,info) + use prec + integer, intent(in) :: m,n,lda + real(pReal), intent(inout), dimension(lda,n) :: a + integer, intent(out), dimension(min(m,n)) :: ipiv + integer, intent(out) :: info + end subroutine dgetrf + + subroutine dgetri(n,a,lda,ipiv,work,lwork,info) + use prec + integer, intent(in) :: n,lda,lwork + real(pReal), intent(inout), dimension(lda,n) :: a + integer, intent(out), dimension(n) :: ipiv + real(pReal), intent(out), dimension(max(1,lwork)) :: work + integer, intent(out) :: info + end subroutine dgetri + + subroutine dsyev(jobz,uplo,n,a,lda,w,work,lwork,info) + use prec + character, intent(in) :: jobz,uplo + integer, intent(in) :: n,lda,lwork + real(pReal), intent(inout), dimension(lda,n) :: a + real(pReal), intent(out), dimension(n) :: w + real(pReal), intent(out), dimension(max(1,lwork)) :: work + integer, intent(out) :: info + end subroutine dsyev + + end interface + +end module LAPACK_interface diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90 index ab26ee9d5..64ad3e1d7 100644 --- a/src/commercialFEM_fileList.f90 +++ b/src/commercialFEM_fileList.f90 @@ -9,6 +9,7 @@ #include "list.f90" #include "future.f90" #include "config.f90" +#include "LAPACK_interface.f90" #include "math.f90" #include "quaternions.f90" #include "Lambert.f90" diff --git a/src/crystallite.f90 b/src/crystallite.f90 index efa066696..9bc254e0c 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -835,8 +835,6 @@ logical function integrateStress(ipc,ip,el,timeFraction) jacoCounterLp, & jacoCounterLi ! counters to check for Jacobian update logical :: error - external :: & - dgesv integrateStress = .false. diff --git a/src/math.f90 b/src/math.f90 index b0852e8d4..3bd3f7112 100644 --- a/src/math.f90 +++ b/src/math.f90 @@ -9,6 +9,7 @@ module math use prec use IO use numerics + use LAPACK_interface implicit none public @@ -489,9 +490,6 @@ function math_invSym3333(A) real(pReal), dimension(6,6) :: temp66 real(pReal), dimension(6*(64+2)) :: work integer :: ierr_i, ierr_f - external :: & - dgetrf, & - dgetri temp66 = math_sym3333to66(A) call dgetrf(6,6,temp66,6,ipiv6,ierr_i) @@ -518,9 +516,6 @@ subroutine math_invert(InvA, error, A) integer, dimension(size(A,1)) :: ipiv real(pReal), dimension(size(A,1)*(64+2)) :: work integer :: ierr - external :: & - dgetrf, & - dgetri invA = A call dgetrf(size(A,1),size(A,1),invA,size(A,1),ipiv,ierr) @@ -885,8 +880,6 @@ subroutine math_eigh(m,w,v,error) logical, intent(out) :: error integer :: ierr real(pReal), dimension((64+2)*size(m,1)) :: work ! block size of 64 taken from http://www.netlib.org/lapack/double/dsyev.f - external :: & - dsyev v = m ! copy matrix to input (doubles as output) array call dsyev('V','U',size(m,1),v,size(m,1),w,work,size(work,1),ierr) @@ -1042,8 +1035,6 @@ function math_eigvalsh(m) real(pReal), dimension(size(m,1),size(m,1)) :: m_ integer :: ierr real(pReal), dimension((64+2)*size(m,1)) :: work ! block size of 64 taken from http://www.netlib.org/lapack/double/dsyev.f - external :: & - dsyev m_= m ! copy matrix to input (will be destroyed) call dsyev('N','U',size(m,1),m_,size(m,1),math_eigvalsh,work,size(work,1),ierr) diff --git a/src/rotations.f90 b/src/rotations.f90 index 7ce366f74..b4e143f28 100644 --- a/src/rotations.f90 +++ b/src/rotations.f90 @@ -596,8 +596,6 @@ function om2ax(om) result(ax) real(pReal), dimension(3,3) :: VR, devNull, om_ integer :: ierr, i - external :: dgeev - om_ = om ! first get the rotation angle