diag_symm_packed Subroutine

public subroutine diag_symm_packed(mode, ldvect, nvect, n, h, eig, vector, ierr)

Uses

  • proc~~diag_symm_packed~~UsesGraph proc~diag_symm_packed diag_symm_packed module~messages messages proc~diag_symm_packed->module~messages comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~io_constants io_constants module~messages->module~io_constants module~precision precision module~messages->module~precision iso_fortran_env iso_fortran_env module~precision->iso_fortran_env

@brief Find eigenvalues and eigenvectors of symmetric matrix in packed format @param[in] mode algorithm of diagonalization (not used now) @param[in] n matrix dimension @param[in] ldvect leading dimension of eigenvector matrix @param[in] nvect required number of eigenvectors @param[in,out] h matrix to be diagonalized @param[out] eig eigenvalues @param[out] vector eigenvectors @param[out] ierr status

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: mode
integer, intent(in) :: ldvect
integer, intent(in) :: nvect
integer, intent(in) :: n
real(kind=dp), intent(inout) :: h(*)
real(kind=dp), intent(out) :: eig(*)
real(kind=dp), intent(out) :: vector(*)
integer, intent(out), optional :: ierr

Calls

proc~~diag_symm_packed~~CallsGraph proc~diag_symm_packed diag_symm_packed dspev dspev proc~diag_symm_packed->dspev dspevx dspevx proc~diag_symm_packed->dspevx interface~show_message show_message proc~diag_symm_packed->interface~show_message

Called by

proc~~diag_symm_packed~~CalledByGraph proc~diag_symm_packed diag_symm_packed proc~matrix_invsqrt matrix_invsqrt proc~matrix_invsqrt->proc~diag_symm_packed proc~rpaeig rpaeig proc~rpaeig->proc~diag_symm_packed proc~run_population_analysis run_population_analysis proc~run_population_analysis->proc~diag_symm_packed proc~guess_hcore guess_hcore proc~guess_hcore->proc~matrix_invsqrt proc~huckel_guess huckel_guess proc~huckel_guess->proc~matrix_invsqrt proc~lowdin lowdin proc~lowdin->proc~run_population_analysis proc~mulliken mulliken proc~mulliken->proc~run_population_analysis proc~scf_driver scf_driver proc~scf_driver->proc~matrix_invsqrt proc~tdhf_energy tdhf_energy proc~tdhf_energy->proc~rpaeig proc~tdhf_mrsf_energy tdhf_mrsf_energy proc~tdhf_mrsf_energy->proc~rpaeig proc~tdhf_sf_energy tdhf_sf_energy proc~tdhf_sf_energy->proc~rpaeig proc~guess_hcore_c guess_hcore_C proc~guess_hcore_c->proc~guess_hcore proc~guess_huckel guess_huckel proc~guess_huckel->proc~huckel_guess proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver proc~tdhf_energy_c tdhf_energy_C proc~tdhf_energy_c->proc~tdhf_energy proc~tdhf_mrsf_energy_c tdhf_mrsf_energy_C proc~tdhf_mrsf_energy_c->proc~tdhf_mrsf_energy proc~tdhf_sf_energy_c tdhf_sf_energy_C proc~tdhf_sf_energy_c->proc~tdhf_sf_energy proc~guess_huckel_c guess_huckel_C proc~guess_huckel_c->proc~guess_huckel

Source Code

  subroutine diag_symm_packed(mode, ldvect, nvect, n, h, eig, vector, ierr)
    use messages, only: show_message, WITH_ABORT, WITHOUT_ABORT
!
    integer, intent(in) :: mode
    integer, intent(in) :: ldvect, nvect, n
    integer, optional, intent(out) :: ierr
    real(dp), intent(inout) :: h(*)
    real(kind=dp), intent(out) :: eig(*), vector(*)

    integer(blas_int), dimension(:), allocatable :: iwork, ifail
    integer(blas_int) :: ldvect_, n_, nvect_, info, ione
    integer :: iok
    real(dp), dimension(:), allocatable :: work
    real(dp) :: abstol, dlamch
    logical :: fatal

    character(16) :: driver

    ldvect_ = int(ldvect, kind=blas_int)
    n_      = int(n, kind=blas_int)
    nvect_  = int(nvect, kind=blas_int)
    ione    = 1

    fatal = WITH_ABORT
    if (present(ierr)) fatal = WITHOUT_ABORT

    allocate (work(n*8), iwork(n*5), ifail(n), stat=iok)
    if (iok /= 0) then
      if (present(ierr)) ierr = iok
      call show_message('Cannot allocate memory', fatal)
      return
    end if

    if (nvect == n .and. ldvect >= n) then
      driver = 'DSPEV'
      call dspev('V', 'U', n_, h, eig, vector, ldvect_, work, info)
    else
      abstol = two*DLAMCH('S')
      driver = 'DSPEVX'
      call dspevx('V', 'A', 'U', &
            ldvect_, h, zero, zero, ione, ione, abstol, n_, &
            eig, vector, nvect_, work, iwork, ifail, info)
    end if

    if (present(ierr)) ierr = info

    if (info /= 0) then
      call show_message('(A,I0)', &
              trim(driver)//' FAILED! INFO: ', int(info), fatal)
    end if

  end subroutine diag_symm_packed