form_rohf_fock Subroutine

public subroutine form_rohf_fock(fock_a_ao, fock_b_ao, fock_mo, MOs, overlap_tri, work, nocca, noccb, nbf, vshift)

Uses

  • proc~~form_rohf_fock~~UsesGraph proc~form_rohf_fock form_rohf_fock module~mathlib mathlib proc~form_rohf_fock->module~mathlib module~oqp_linalg 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 iso_fortran_env iso_fortran_env module~precision->iso_fortran_env 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 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 Form the ROHF Fock matrix in MO basis.

This implementation is based on the work of M. F. Guest and V. Saunders, as described in Mol. Phys. 28, 819 (1974).

@brief Form the ROHF Fock matrix in the MO basis using the Guest-Saunders method.

@detail This subroutine transforms the alpha and beta Fock matrices from the AO to the MO basis, constructs the ROHF Fock matrix using the Guest-Saunders method, and optionally applies a level shift to the diagonal elements of the virtual-virtual block.

@author Konstantin Komarov, 2023

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(inout), dimension(:) :: fock_a_ao
real(kind=dp), intent(inout), dimension(:) :: fock_b_ao
real(kind=dp), intent(out), dimension(:) :: fock_mo
real(kind=dp), intent(in), dimension(:,:) :: MOs
real(kind=dp), intent(in), dimension(:) :: overlap_tri
real(kind=dp), intent(out), dimension(:) :: work
integer, intent(in) :: nocca
integer, intent(in) :: noccb
integer, intent(in) :: nbf
real(kind=dp), intent(in) :: vshift

Calls

proc~~form_rohf_fock~~CallsGraph proc~form_rohf_fock form_rohf_fock dsymm dsymm proc~form_rohf_fock->dsymm interface~pack_matrix pack_matrix proc~form_rohf_fock->interface~pack_matrix interface~unpack_matrix unpack_matrix proc~form_rohf_fock->interface~unpack_matrix proc~orthogonal_transform2 orthogonal_transform2 proc~form_rohf_fock->proc~orthogonal_transform2 proc~orthogonal_transform_sym orthogonal_transform_sym proc~form_rohf_fock->proc~orthogonal_transform_sym proc~triangular_to_full triangular_to_full proc~form_rohf_fock->proc~triangular_to_full proc~pack_f90 PACK_F90 interface~pack_matrix->proc~pack_f90 proc~unpack_f90 UNPACK_F90 interface~unpack_matrix->proc~unpack_f90 interface~show_message show_message proc~orthogonal_transform2->interface~show_message proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~orthogonal_transform2->proc~oqp_dgemm_i64 proc~orthogonal_transform_sym->interface~show_message 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~triangular_to_full->interface~show_message proc~oqp_dgemm_i64->interface~show_message dgemm dgemm proc~oqp_dgemm_i64->dgemm proc~oqp_dsymm_i64->dsymm proc~oqp_dsymm_i64->interface~show_message 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~pack_f90->interface~show_message proc~pack_f90->proc~oqp_dtrttp_i64 proc~unpack_f90->interface~show_message proc~unpack_f90->proc~oqp_dtpttr_i64

Called by

proc~~form_rohf_fock~~CalledByGraph proc~form_rohf_fock form_rohf_fock proc~scf_driver scf_driver proc~scf_driver->proc~form_rohf_fock proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

  subroutine form_rohf_fock(fock_a_ao, fock_b_ao, fock_mo, &
                            MOs, overlap_tri, work, &
                            nocca, noccb, nbf, vshift)

    use mathlib, only: orthogonal_transform_sym, &
                       orthogonal_transform2, &
                       triangular_to_full, &
                       unpack_matrix, &
                       pack_matrix

    implicit none

    real(kind=dp), intent(inout), dimension(:) :: fock_a_ao
    real(kind=dp), intent(inout), dimension(:) :: fock_b_ao
    real(kind=dp), intent(out), dimension(:) :: fock_mo
    real(kind=dp), intent(in), dimension(:,:) :: MOs
    real(kind=dp), intent(in), dimension(:) :: overlap_tri
    real(kind=dp), intent(out), dimension(:) :: work
    integer, intent(in) :: nocca, noccb, nbf
    real(kind=dp), intent(in) :: vshift

    real(kind=dp), allocatable, dimension(:,:) :: &
          overlap, work_matrix, fock, fock_a, fock_b
    real(kind=dp) :: acc, aoo, avv, bcc, boo, bvv
    integer :: i, nbf_tri

    acc = 0.5_dp; aoo = 0.5_dp; avv = 0.5_dp
    bcc = 0.5_dp; boo = 0.5_dp; bvv = 0.5_dp
    nbf_tri = nbf * (nbf + 1) / 2

    ! Allocate full matrices
    allocate(overlap(nbf, nbf), &
             work_matrix(nbf, nbf), &
             fock(nbf, nbf), &
             fock_a(nbf, nbf), &
             fock_b(nbf, nbf), &
             source=0.0_dp)

    ! Transform alpha and beta Fock matrices to MO basis
    call orthogonal_transform_sym(nbf, nbf, fock_a_ao, MOs, nbf, fock_mo)
    fock_a_ao(:nbf_tri) = fock_mo(:nbf_tri)

    call orthogonal_transform_sym(nbf, nbf, fock_b_ao, MOs, nbf, fock_mo)
    fock_b_ao(:nbf_tri) = fock_mo(:nbf_tri)

    ! Unpack triangular matrices to full matrices
    call unpack_matrix(fock_a_ao, fock_a)
    call unpack_matrix(fock_b_ao, fock_b)

    ! Construct ROHF Fock matrix in MO basis using Guest-Saunders method
    associate ( na => nocca &
              , nb => noccb &
      )
      fock(1:nb, 1:nb) = acc * fock_a(1:nb, 1:nb) &
                       + bcc * fock_b(1:nb, 1:nb)
      fock(nb+1:na, nb+1:na) = aoo * fock_a(nb+1:na, nb+1:na) &
                             + boo * fock_b(nb+1:na, nb+1:na)
      fock(na+1:nbf, na+1:nbf) = avv * fock_a(na+1:nbf, na+1:nbf) &
                               + bvv * fock_b(na+1:nbf, na+1:nbf)
      fock(1:nb, nb+1:na) = fock_b(1:nb, nb+1:na)
      fock(nb+1:na, 1:nb) = fock_b(nb+1:na, 1:nb)
      fock(1:nb, na+1:nbf) = 0.5_dp * (fock_a(1:nb, na+1:nbf) &
                                     + fock_b(1:nb, na+1:nbf))
      fock(na+1:nbf, 1:nb) = 0.5_dp * (fock_a(na+1:nbf, 1:nb) &
                                     + fock_b(na+1:nbf, 1:nb))
      fock(nb+1:na, na+1:nbf) = fock_a(nb+1:na, na+1:nbf)
      fock(na+1:nbf, nb+1:na) = fock_a(na+1:nbf, nb+1:na)

      ! Apply Vshift to the diagonal
      do i = nb+1, na
        fock(i,i) = fock(i,i) + vshift * 0.5_dp
      end do
      do i = na+1, nbf
        fock(i,i) = fock(i,i) + vshift
      end do
    end associate

    ! Pack fock to fock_mo, make it triangular
    call pack_matrix(fock, fock_mo)

    call triangular_to_full(fock_mo, nbf, 'u')

    ! Back-transform ROHF Fock matrix to AO basis
    call unpack_matrix(overlap_tri, overlap)
    call dsymm('l', 'u', nbf, nbf, &
               1.0_dp, overlap, nbf, &
                       MOs, nbf, &
               0.0_dp, work, nbf)
    call orthogonal_transform2('t', nbf, nbf, work, nbf, fock, nbf, &
                               work_matrix, nbf, overlap)
    call pack_matrix(work_matrix, fock_a_ao)

  end subroutine form_rohf_fock