mo_to_ao Subroutine

public subroutine mo_to_ao(fao, fmo, s, v, nmo, nbf)

Uses

  • proc~~mo_to_ao~~UsesGraph proc~mo_to_ao mo_to_ao module~mathlib mathlib proc~mo_to_ao->module~mathlib module~oqp_linalg oqp_linalg proc~mo_to_ao->module~oqp_linalg module~mathlib->module~oqp_linalg module~precision precision module~mathlib->module~precision module~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap module~blas_wrap->module~precision module~mathlib_types mathlib_types module~blas_wrap->module~mathlib_types module~messages messages module~blas_wrap->module~messages module~lapack_wrap->module~precision module~lapack_wrap->module~mathlib_types module~lapack_wrap->module~messages iso_fortran_env iso_fortran_env module~precision->iso_fortran_env module~messages->module~precision 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

@brief Back-transform a symmetric operator Fmo expressed in the MO basis V to the AO basis @detail compute the transformation: Fao = S*V * Fmo * (SV)^T

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out) :: fao(*)
real(kind=dp), intent(in) :: fmo(*)
real(kind=dp), intent(in) :: s(*)
real(kind=dp), intent(in) :: v(*)
integer, intent(in) :: nmo
integer, intent(in) :: nbf

Calls

proc~~mo_to_ao~~CallsGraph proc~mo_to_ao mo_to_ao interface~pack_matrix pack_matrix proc~mo_to_ao->interface~pack_matrix interface~unpack_matrix unpack_matrix proc~mo_to_ao->interface~unpack_matrix proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~mo_to_ao->proc~oqp_dgemm_i64 proc~oqp_dsymm_i64 oqp_dsymm_i64 proc~mo_to_ao->proc~oqp_dsymm_i64 proc~pack_f90 PACK_F90 interface~pack_matrix->proc~pack_f90 proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 dgemm dgemm proc~oqp_dgemm_i64->dgemm interface~show_message show_message proc~oqp_dgemm_i64->interface~show_message dsymm dsymm proc~oqp_dsymm_i64->dsymm proc~oqp_dsymm_i64->interface~show_message proc~pack_f90->interface~show_message proc~oqp_dtrttp_i64 oqp_dtrttp_i64 proc~pack_f90->proc~oqp_dtrttp_i64 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~oqp_dtrttp_i64->interface~show_message dtrttp dtrttp proc~oqp_dtrttp_i64->dtrttp

Called by

proc~~mo_to_ao~~CalledByGraph proc~mo_to_ao mo_to_ao proc~scf_driver scf_driver proc~scf_driver->proc~mo_to_ao proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

  subroutine mo_to_ao(fao, fmo, s, v, nmo, nbf)

    use mathlib, only: pack_matrix, unpack_matrix
    use oqp_linalg

    implicit none

    real(kind=dp), intent(in) :: fmo(*), s(*), v(*)
    real(kind=dp), intent(out) :: fao(*)
    integer, intent(in) :: nmo, nbf

    integer :: nbf2
    real(kind=dp), allocatable :: sv(:,:), ftmp(:,:), wrk(:,:)

    allocate(ftmp(nbf,nbf), sv(nbf,nmo), wrk(nbf,nbf))

    call unpack_matrix(fmo, ftmp)

    call unpack_matrix(s, wrk)
    ! compute S*V
    call dsymm('l', 'u', nbf, nmo, &
               1.0_dp, wrk, nbf, &
                       v,  nbf, &
               0.0_dp, sv, nbf)

    ! compute (S * V) * Fmo
    call dsymm('r', 'u', nbf, nmo, &
               1.0d0, ftmp, nbf, &
                      sv,  nbf, &
               0.0d0, wrk, nbf)

    ! compute ((S * V) * Fmo) * (S * V)^T
    call dgemm('n', 't', nbf, nbf, nmo, &
               1.0d0, wrk, nbf, &
                      sv, nbf, &
               0.0d0, ftmp, nbf)

    nbf2 = nbf*(nbf+1)/2
    call pack_matrix(ftmp, fao(:nbf2))

  end subroutine mo_to_ao