print_eigvec_vals_labeled Subroutine

public subroutine print_eigvec_vals_labeled(basis, infos, mostart, moend)

Uses

  • proc~~print_eigvec_vals_labeled~~UsesGraph proc~print_eigvec_vals_labeled print_eigvec_vals_labeled module~basis_tools basis_tools proc~print_eigvec_vals_labeled->module~basis_tools module~io_constants io_constants proc~print_eigvec_vals_labeled->module~io_constants module~messages messages proc~print_eigvec_vals_labeled->module~messages module~oqp_tagarray_driver oqp_tagarray_driver proc~print_eigvec_vals_labeled->module~oqp_tagarray_driver module~types types proc~print_eigvec_vals_labeled->module~types module~basis_tools->module~io_constants iso_fortran_env iso_fortran_env module~basis_tools->iso_fortran_env module~atomic_structure_m atomic_structure_m module~basis_tools->module~atomic_structure_m module~constants constants module~basis_tools->module~constants module~parallel parallel module~basis_tools->module~parallel module~precision precision module~basis_tools->module~precision module~messages->module~io_constants comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~messages->module~precision iso_c_binding iso_c_binding module~oqp_tagarray_driver->iso_c_binding tagarray tagarray module~oqp_tagarray_driver->tagarray module~types->module~basis_tools module~types->iso_c_binding module~types->module~atomic_structure_m module~functionals functionals module~types->module~functionals module~types->module~parallel module~types->module~precision module~types->tagarray module~atomic_structure_m->iso_c_binding module~constants->module~precision module~functionals->iso_c_binding module~functionals->module~precision xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~parallel->iso_c_binding module~parallel->iso_fortran_env module~parallel->module~precision mpi mpi module~parallel->mpi module~precision->iso_fortran_env

@brief print eigenvector/values, with MO symmetry labels

Arguments

Type IntentOptional Attributes Name
type(basis_set), intent(in) :: basis
type(information), intent(inout) :: infos
integer, intent(in) :: mostart
integer, intent(in) :: moend

Calls

proc~~print_eigvec_vals_labeled~~CallsGraph proc~print_eigvec_vals_labeled print_eigvec_vals_labeled interface~data_has_tags data_has_tags proc~print_eigvec_vals_labeled->interface~data_has_tags interface~tagarray_get_data tagarray_get_data proc~print_eigvec_vals_labeled->interface~tagarray_get_data none~bf_label basis_set%bf_label proc~print_eigvec_vals_labeled->none~bf_label

Called by

proc~~print_eigvec_vals_labeled~~CalledByGraph proc~print_eigvec_vals_labeled print_eigvec_vals_labeled proc~print_mo_range print_mo_range proc~print_mo_range->proc~print_eigvec_vals_labeled proc~scf_driver scf_driver proc~scf_driver->proc~print_mo_range proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

  subroutine print_eigvec_vals_labeled(basis, infos, mostart, moend)

    use io_constants, only: iw
    use oqp_tagarray_driver
    use types, only: information
    use basis_tools, only: basis_set
    use messages, only: show_message, with_abort
    implicit none

    character(len=*), parameter :: subroutine_name = "print_eigvec_vals_labeled"

    type(basis_set), intent(in) :: basis
    type(information), intent(inout) :: infos

    integer, intent(in) :: mostart, moend

    integer :: i, imax, imin, j, nmax
    character(len=*), parameter :: fmt1 = '(/15x,10(8x,i4,5x))'
    character(len=*), parameter :: fmt2 = '(15x,10f17.10)'
    character(len=*), parameter :: fmt4 = '(i5,2x,a8,10f17.10)'

    ! tagarray
    real(kind=dp), contiguous, pointer :: &
      mo_energy_a(:), mo_energy_b(:), mo_a(:,:), mo_b(:,:)
    character(len=*), parameter :: tags_alpha(2) = (/ character(len=80) :: &
      OQP_E_MO_A, OQP_VEC_MO_A /)
    character(len=*), parameter :: tags_beta(2) = (/ character(len=80) :: &
      OQP_E_MO_B, OQP_VEC_MO_B /)

!Print out eigendata, with mo symmetry labels
!The rows are labeled with the basis function names.

    nmax = 5

    call data_has_tags(infos%dat, tags_alpha, module_name, subroutine_name, WITH_ABORT)
    call tagarray_get_data(infos%dat, OQP_E_MO_A, mo_energy_a)
    call tagarray_get_data(infos%dat, OQP_VEC_MO_A, mo_a)
    write(iw,'(/,A)') '   -------------- Alpha Orbitals -------------'
    do imin = mostart, moend, nmax
      imax = min(imin+nmax-1, moend)

      write(iw,fmt1) (i, i=imin, imax)
      write(iw,fmt2) (mo_energy_a(i), i=imin, imax)

      do j = 1, basis%nbf
        write(iw,fmt4) j, basis%bf_label(j), mo_a(j,imin:imax)
      end do

    end do

    if ((infos%mol_prop%nelec_b /= 0) .and. (infos%control%scftype == 2)) then
      call data_has_tags(infos%dat, tags_beta, module_name, subroutine_name, WITH_ABORT)
      call tagarray_get_data(infos%dat, OQP_E_MO_B, mo_energy_b)
      call tagarray_get_data(infos%dat, OQP_VEC_MO_B, mo_b)
      write(iw,'(/,A)') '   -------------- Beta Orbitals -------------'
      do imin = mostart, moend, nmax
        imax = min(imin+nmax-1, moend)

        write(iw,fmt1) (i, i=imin, imax)
        write(iw,fmt2) (mo_energy_b(i), i=imin, imax)

        do j = 1, basis%nbf
          write(iw,fmt4) j, basis%bf_label(j), mo_b(j,imin:imax)
        end do

      end do
    end if

  end subroutine print_eigvec_vals_labeled