guess_hcore Subroutine

public subroutine guess_hcore(infos)

Uses

  • proc~~guess_hcore~~UsesGraph proc~guess_hcore guess_hcore module~basis_tools basis_tools proc~guess_hcore->module~basis_tools module~guess guess proc~guess_hcore->module~guess module~io_constants io_constants proc~guess_hcore->module~io_constants module~mathlib mathlib proc~guess_hcore->module~mathlib module~messages messages proc~guess_hcore->module~messages module~oqp_tagarray_driver oqp_tagarray_driver proc~guess_hcore->module~oqp_tagarray_driver module~precision precision proc~guess_hcore->module~precision module~printing printing proc~guess_hcore->module~printing module~strings strings proc~guess_hcore->module~strings module~types types proc~guess_hcore->module~types module~util util proc~guess_hcore->module~util module~basis_tools->module~io_constants module~basis_tools->module~precision 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~guess->module~precision module~oqp_linalg oqp_linalg module~guess->module~oqp_linalg module~mathlib->module~precision module~mathlib->module~oqp_linalg module~messages->module~io_constants module~messages->module~precision comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR iso_c_binding iso_c_binding module~oqp_tagarray_driver->iso_c_binding tagarray tagarray module~oqp_tagarray_driver->tagarray module~precision->iso_fortran_env module~printing->module~precision module~strings->iso_c_binding module~types->module~basis_tools module~types->module~precision module~types->iso_c_binding module~types->module~atomic_structure_m module~functionals functionals module~types->module~functionals module~types->module~parallel module~types->tagarray module~util->module~precision module~atomic_structure_m->iso_c_binding module~constants->module~precision module~functionals->module~precision module~functionals->iso_c_binding xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap module~parallel->module~precision module~parallel->iso_c_binding module~parallel->iso_fortran_env mpi mpi module~parallel->mpi module~blas_wrap->module~messages module~blas_wrap->module~precision module~mathlib_types mathlib_types module~blas_wrap->module~mathlib_types module~lapack_wrap->module~messages module~lapack_wrap->module~precision module~lapack_wrap->module~mathlib_types

Arguments

Type IntentOptional Attributes Name
type(information), intent(inout), target :: infos

Calls

proc~~guess_hcore~~CallsGraph proc~guess_hcore guess_hcore interface~data_has_tags data_has_tags proc~guess_hcore->interface~data_has_tags interface~show_message show_message proc~guess_hcore->interface~show_message interface~tagarray_get_data tagarray_get_data proc~guess_hcore->interface~tagarray_get_data proc~get_ab_initio_density get_ab_initio_density proc~guess_hcore->proc~get_ab_initio_density proc~get_ab_initio_orbital get_ab_initio_orbital proc~guess_hcore->proc~get_ab_initio_orbital proc~matrix_invsqrt matrix_invsqrt proc~guess_hcore->proc~matrix_invsqrt proc~measure_time measure_time proc~guess_hcore->proc~measure_time proc~print_module_info print_module_info proc~guess_hcore->proc~print_module_info remove_records remove_records proc~guess_hcore->remove_records reserve_data reserve_data proc~guess_hcore->reserve_data proc~get_ab_initio_density->interface~show_message proc~orb_to_dens orb_to_dens proc~get_ab_initio_density->proc~orb_to_dens proc~get_ab_initio_orbital->interface~show_message interface~unpack_matrix unpack_matrix proc~get_ab_initio_orbital->interface~unpack_matrix proc~diag_symm_full diag_symm_full proc~get_ab_initio_orbital->proc~diag_symm_full proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~get_ab_initio_orbital->proc~oqp_dgemm_i64 proc~orthogonal_transform orthogonal_transform proc~get_ab_initio_orbital->proc~orthogonal_transform proc~matrix_invsqrt->interface~show_message proc~diag_symm_packed diag_symm_packed proc~matrix_invsqrt->proc~diag_symm_packed proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 proc~diag_symm_full->interface~show_message dsyev dsyev proc~diag_symm_full->dsyev proc~diag_symm_packed->interface~show_message dspev dspev proc~diag_symm_packed->dspev dspevx dspevx proc~diag_symm_packed->dspevx proc~oqp_dgemm_i64->interface~show_message dgemm dgemm proc~oqp_dgemm_i64->dgemm interface~pack_matrix pack_matrix proc~orb_to_dens->interface~pack_matrix proc~oqp_dsyr2k_i64 oqp_dsyr2k_i64 proc~orb_to_dens->proc~oqp_dsyr2k_i64 proc~orthogonal_transform->interface~show_message proc~orthogonal_transform->proc~oqp_dgemm_i64 proc~pack_f90 PACK_F90 interface~pack_matrix->proc~pack_f90 proc~oqp_dsyr2k_i64->interface~show_message dsyr2k dsyr2k proc~oqp_dsyr2k_i64->dsyr2k proc~unpack_f90->interface~show_message proc~oqp_dtpttr_i64 oqp_dtpttr_i64 proc~unpack_f90->proc~oqp_dtpttr_i64 proc~oqp_dtpttr_i64->interface~show_message dtpttr dtpttr proc~oqp_dtpttr_i64->dtpttr proc~pack_f90->interface~show_message proc~oqp_dtrttp_i64 oqp_dtrttp_i64 proc~pack_f90->proc~oqp_dtrttp_i64 proc~oqp_dtrttp_i64->interface~show_message dtrttp dtrttp proc~oqp_dtrttp_i64->dtrttp

Called by

proc~~guess_hcore~~CalledByGraph proc~guess_hcore guess_hcore proc~guess_hcore_c guess_hcore_C proc~guess_hcore_c->proc~guess_hcore

Source Code

  subroutine guess_hcore(infos)
    use precision, only: dp
    use types, only: information
    use io_constants, only: IW
    use oqp_tagarray_driver
    use basis_tools, only: basis_set
    use guess, only: get_ab_initio_density, &
                     get_ab_initio_orbital
    use mathlib, only: matrix_invsqrt
    use util, only: measure_time
    use messages, only: show_message, WITH_ABORT
    use strings, only: Cstring, fstring
    use printing, only: print_module_info

    implicit none

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

    type(information), target, intent(inout) :: infos
  !
    integer :: nbf2, nbf, ok
  !
    real(kind=dp), allocatable :: qmat(:,:)
    type(basis_set), pointer :: basis
  !
  ! tagarray
    real(kind=dp), contiguous, pointer :: &
      Hcore(:), Smat(:), &
      dmat_a(:), mo_a(:,:), mo_energy_a(:), &
      dmat_b(:), mo_b(:,:), mo_energy_b(:)
    character(len=*), parameter :: tags_alpha(3) = (/ character(len=80) :: &
      OQP_DM_A, OQP_E_MO_A, OQP_VEC_MO_A /)
    character(len=*), parameter :: tags_beta(3) = (/ character(len=80) :: &
      OQP_DM_B, OQP_E_MO_B, OQP_VEC_MO_B /)
    character(len=*), parameter :: tags_general(2) = (/ character(len=80) :: &
      OQP_SM, OQP_Hcore /)

  ! Files open
  ! 1. XYZ: Read : Geometric data, ATOMS
  ! 3. LOG: Read Write: Main output file
  !
    open (unit=IW, file=infos%log_filename, position="append")
  !
  !
    call print_module_info('guess_Hcore','Initial Guess using H Matrix')

  ! load basis set
    basis => infos%basis
    basis%atoms => infos%atoms

  !  Allocate H, S ,T and D matrices
    nbf = basis%nbf
    nbf2 = nbf*(nbf+1)/2

    allocate (qmat(nbf,nbf), stat=ok)
    if (ok /= 0) call show_message('Cannot allocate memory', WITH_ABORT)

    ! clean data
    call infos%dat%remove_records(tags_alpha)
    call infos%dat%remove_records(tags_beta)

    ! load general data
    call data_has_tags(infos%dat, tags_general, module_name, subroutine_name, WITH_ABORT)
    call tagarray_get_data(infos%dat, OQP_SM, smat)
    call tagarray_get_data(infos%dat, OQP_Hcore, hcore)

    ! allocate alpha
    call infos%dat%reserve_data(OQP_DM_A, TA_TYPE_REAL64, nbf2, comment=OQP_DM_A_comment)
    call infos%dat%reserve_data(OQP_E_MO_A, TA_TYPE_REAL64, nbf, comment=OQP_E_MO_A_comment)
    call infos%dat%reserve_data(OQP_VEC_MO_A, TA_TYPE_REAL64, nbf*nbf, (/ nbf, nbf /), comment=OQP_VEC_MO_A_comment)

    ! load alpha data
    call data_has_tags(infos%dat, tags_alpha, module_name, subroutine_name, WITH_ABORT)
    call tagarray_get_data(infos%dat, OQP_DM_A, dmat_a)
    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)

  ! UHF/ROHF
    if (infos%control%scftype >= 2) then
      ! allocate beta
      call infos%dat%reserve_data(OQP_DM_B, TA_TYPE_REAL64, nbf2, comment=OQP_DM_B_comment)
      call infos%dat%reserve_data(OQP_E_MO_B, TA_TYPE_REAL64, nbf, comment=OQP_E_MO_B_comment)
      call infos%dat%reserve_data(OQP_VEC_MO_B, TA_TYPE_REAL64, nbf*nbf, (/ nbf, nbf /), comment=OQP_VEC_MO_B_comment)

      ! load beta
      call data_has_tags(infos%dat, tags_beta, module_name, subroutine_name, WITH_ABORT)
      call tagarray_get_data(infos%dat, OQP_DM_B, dmat_b)
      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)
    end if

  !  End of Readings.................................
  !
    call matrix_invsqrt(smat, qmat, nbf)
  !
  !  Calculate Hcore MO
    call Get_ab_initio_orbital(Hcore, MO_A, MO_Energy_A, QMat)
  !  For ROHF/UHF
    if (INFOS%control%scftype >= 2) MO_B = MO_A
  !
  ! Calculate Density Matrix
  ! RHF
    if (INFOS%control%scftype == 1) then
      call get_ab_initio_density(Dmat_A, MO_A, infos=infos, basis=basis)
  ! ROHF/UHF
    else
      call get_ab_initio_density(Dmat_A, MO_A, Dmat_B, MO_B, infos, basis)
    end if
    write (IW, 9090)
    call measure_time(print_total=1, log_unit=iw)
    close(IW)
  9090 format(/1x, '...... End Of Initial Orbital Guess ......'/)

  end subroutine guess_hcore