I'm new to Fortran, and this is pretty much my first escapade. Below is a function that I wrote which relies on calls to LAPACK. The function is sat in a module with some other functions and works perfectly, but seeing as this is really the workhorse of the program I'm building I want to squeeze it hard for performance. How does it look? Am I doing anything stupid, or would another call to BLAS or LAPACK somewhere improve the performance?The function takes a Hermitian matrix H, and returns the matrix exponential of the skew-Hermitian matrix -iHt where i is the imaginary number and t is a real number. (This is the solution to the Schrodinger equation for a time independent Hamiltonian.) s, the length of one side of the Hamiltonian, is included for automatic array initialisation rather than using allocations. Typically it'll be less than 10.function time_indep_schrodinger(s,H,t)    ! s : The length of one side of H.    ! H : The Hamiltonian.    ! t : The time t.    !    ! The LAPACK subroutine zheev requires a tridiagonal    ! subcopy of H. This array is then transformed into    ! a matrix of eigenvectors. This also yields the    ! eigenvalues as a 1D array, which are then    ! exponentiated before using the matrices of    ! eigenvectors to produce the evolution operator.    !    ! Finds the eigenvalues of Ht, then uses U exp(-i*eigs) U^H.    integer,    intent(in) :: s    complex*16, intent(in) :: H(s,s)    real*8,     intent(in) :: t    complex*16             :: B(s,s), eigv(s,s), time_indep_schrodinger(s,s), work(2*s-1)    real*8                 :: rwork(3*s-2), eigs(s)    integer                :: info, n, m    ! Hermitian matrix to diagonalise.    forall (n=1:s, m=1:s, m>=n) eigv(n,m) = H(n,m)*t    ! Get eigenvectors and eigenvalues.    call zheev('V','U',s,eigv,s,eigs,work,2*s-1,rwork,info)    ! Scale columns of copied matrix be exponentiated eigenvalues (with -i).    do n=1,s      B(:,n) = eigv(:,n)*exp((0,-1)*eigs(n))    end do     ! Finally multiply scaled eigenvectors with conjugate of eigv.    time_indep_schrodinger = matmul(B,conjg(transpose(eigv)))end function
High performance exponential of a skew Hermitian matrix in Fortran 95
beginner;mathematics;matrix;fortran