get_ab_initio_orbital Subroutine

public subroutine get_ab_initio_orbital(fock_m, orbitals, orbital_e, q_m)

Uses

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

@brief Solve \f$ F C = \eps S C \f$

Arguments

Type IntentOptional Attributes Name
real(kind=dp) :: fock_m(*)
real(kind=dp) :: orbitals(:,:)
real(kind=dp) :: orbital_e(*)
real(kind=dp) :: q_m(:,:)

Calls

proc~~get_ab_initio_orbital~~CallsGraph proc~get_ab_initio_orbital get_ab_initio_orbital interface~show_message show_message 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~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~oqp_dgemm_i64->interface~show_message dgemm dgemm proc~oqp_dgemm_i64->dgemm proc~orthogonal_transform->interface~show_message proc~orthogonal_transform->proc~oqp_dgemm_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

Called by

proc~~get_ab_initio_orbital~~CalledByGraph proc~get_ab_initio_orbital get_ab_initio_orbital proc~guess_hcore guess_hcore proc~guess_hcore->proc~get_ab_initio_orbital proc~scf_driver scf_driver proc~scf_driver->proc~get_ab_initio_orbital proc~guess_hcore_c guess_hcore_C proc~guess_hcore_c->proc~guess_hcore proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

 subroutine get_ab_initio_orbital(fock_m, orbitals, orbital_e, q_m)
   use messages,  only: show_message, WITH_ABORT
   use mathlib, only: orthogonal_transform
   use eigen, only: diag_symm_full
   use mathlib, only: unpack_matrix, pack_matrix
   implicit none
!
   real(kind=dp) :: orbital_e(*)
   real(kind=dp) :: fock_m(*)
   real(kind=dp) :: orbitals(:,:), q_m(:,:)

   real(kind=dp), allocatable :: wrk(:,:), f2(:,:), tfock(:,:)
   integer  :: info, ok, nbf

   nbf = ubound(orbitals, 1)
   allocate(wrk(nbf,nbf), f2(nbf,nbf), tfock(nbf, nbf), stat=ok)
   if (ok/=0) call show_message('Cannot allocate memory', WITH_ABORT)

   call unpack_matrix(fock_m, f2)

   ! F = S^{-1/2} * F * S^{-1/2}
   ! can use orthogonal_transform subroutine, because S^{-1/2} is symmetric
   call orthogonal_transform('n', nbf, q_m, f2, tfock, wrk)
   deallocate(f2, wrk)

   ! Get orbitals C' = S^{1/2} C
   call diag_symm_full(1, nbf, tfock, nbf, orbital_e, info)

   ! Recover orbitals in AO basis
   ! C = S^{-1/2} * C'
   call dgemm('n', 'n', nbf, nbf, nbf, &
              1.0_dp, q_m,      nbf, &
                      tfock,    nbf, &
              0.0_dp, orbitals, nbf)

 end subroutine get_ab_initio_orbital