eijden Subroutine

public subroutine eijden(eps, nbf, infos)

Uses

  • proc~~eijden~~UsesGraph proc~eijden eijden module~mathlib mathlib proc~eijden->module~mathlib module~messages messages proc~eijden->module~messages module~oqp_tagarray_driver oqp_tagarray_driver proc~eijden->module~oqp_tagarray_driver module~oqp_linalg oqp_linalg module~mathlib->module~oqp_linalg module~precision precision 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 iso_c_binding iso_c_binding module~oqp_tagarray_driver->iso_c_binding tagarray tagarray module~oqp_tagarray_driver->tagarray 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~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

@brief Compute "energy weighted density matrix"

Note

This quantity is actually the Lagrangian matrix,

Arguments

Type IntentOptional Attributes Name
real(kind=dp) :: eps(:)
integer :: nbf
type(information), intent(inout) :: infos

Calls

proc~~eijden~~CallsGraph proc~eijden eijden interface~data_has_tags data_has_tags proc~eijden->interface~data_has_tags interface~tagarray_get_data tagarray_get_data proc~eijden->interface~tagarray_get_data interface~unpack_matrix unpack_matrix proc~eijden->interface~unpack_matrix proc~orthogonal_transform_sym orthogonal_transform_sym proc~eijden->proc~orthogonal_transform_sym proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 interface~show_message show_message proc~orthogonal_transform_sym->interface~show_message proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~orthogonal_transform_sym->proc~oqp_dgemm_i64 proc~oqp_dsymm_i64 oqp_dsymm_i64 proc~orthogonal_transform_sym->proc~oqp_dsymm_i64 proc~oqp_dtpttr_i64 oqp_dtpttr_i64 proc~orthogonal_transform_sym->proc~oqp_dtpttr_i64 proc~oqp_dtrttp_i64 oqp_dtrttp_i64 proc~orthogonal_transform_sym->proc~oqp_dtrttp_i64 proc~oqp_dgemm_i64->interface~show_message dgemm dgemm proc~oqp_dgemm_i64->dgemm proc~oqp_dsymm_i64->interface~show_message dsymm dsymm proc~oqp_dsymm_i64->dsymm 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 proc~unpack_f90->interface~show_message proc~unpack_f90->proc~oqp_dtpttr_i64

Called by

proc~~eijden~~CalledByGraph proc~eijden eijden proc~hf_gradient hf_gradient proc~hf_gradient->proc~eijden proc~tdhf_1e_grad tdhf_1e_grad proc~tdhf_1e_grad->proc~eijden proc~tdhf_gradient tdhf_gradient proc~tdhf_gradient->proc~tdhf_1e_grad proc~tdhf_gradient_c tdhf_gradient_C proc~tdhf_gradient_c->proc~tdhf_gradient

Source Code

  subroutine eijden(eps, nbf, infos)
    use oqp_tagarray_driver
    use mathlib, only: orthogonal_transform_sym
    use messages, only: show_message, with_abort

    implicit none

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

    type(information), intent(inout) :: infos
    integer :: nbf
    integer :: i, j, ij, ne, ok
    real(kind=dp) :: eps(:)
    real(kind=dp), allocatable :: c(:,:), scr(:),tempd(:)

    ! tagarray
    real(kind=dp), contiguous, pointer :: &
      mo_energy_a(:), mo_a(:,:), &
      fock_a(:), fock_b(:), dmat_a(:), dmat_b(:)
    character(len=*), parameter :: tags_alpha(2) = (/ character(len=80) :: &
      OQP_E_MO_A, OQP_VEC_MO_A /)
    character(len=*), parameter :: tags_beta(4) = (/ character(len=80) :: &
      OQP_FOCK_A, OQP_DM_A, OQP_FOCK_B, OQP_DM_B /)

    if (infos%control%scftype>1) then
       allocate(c(nbf,nbf), scr(8*nbf), tempd(nbf*(nbf+1)/2), stat=ok)
    end if

    ne = infos%mol_prop%nelec/2

    select case (infos%control%scftype)
!   RHF case
    case (1)

      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)

       ij = 0
       do i = 1, nbf
          do j = 1, i
             ij = ij+1
             eps(ij) = -2*sum(mo_energy_a(1:ne)&
                             *mo_a(i,1:ne)&
                             *mo_a(j,1:ne))
          end do
       end do

!   U/ROHF case
    case (2:)

      call data_has_tags(infos%dat, tags_beta, module_name, subroutine_name, WITH_ABORT)
      call tagarray_get_data(infos%dat, OQP_FOCK_A, fock_a)
      call tagarray_get_data(infos%dat, OQP_FOCK_B, fock_b)
      call tagarray_get_data(infos%dat, OQP_DM_A, dmat_a)
      call tagarray_get_data(infos%dat, OQP_DM_B, dmat_b)

!     Alpha part
      call unpack_matrix(dmat_a,c,nbf,'U')
      call orthogonal_transform_sym(nbf, nbf, fock_a, c, nbf, tempd)

!     Beta part
      call unpack_matrix(dmat_b,c,nbf,'U')
      call orthogonal_transform_sym(nbf, nbf, fock_b, c, nbf, eps)
      eps = -eps - tempd

!     Half the diagonal
      ij = 0
      do i = 1, nbf
         ij = ij+i
         eps(ij) = 0.5d0*eps(ij)
      end do

    end select

  end subroutine eijden