mo_overlap Subroutine

public subroutine mo_overlap(overlap_mo, mo_a, mo_a_old, overlap_ao, nbf)

Uses

  • proc~~mo_overlap~~UsesGraph proc~mo_overlap mo_overlap module~precision precision proc~mo_overlap->module~precision iso_fortran_env iso_fortran_env module~precision->iso_fortran_env

@brief Calculate Molecular Orbital (MO) overlap between two geometries

This subroutine computes the overlap between Molecular Orbitals of two different geometries, given their MO coefficients and the Atomic Orbital (AO) overlap.

@param[out] overlap_mo Resulting MO overlap matrix @param[in] mo_a MO coefficients obtained at the current geometry @param[in] mo_a_old MO coefficients obtained at the old geometry @param[in] overlap_ao AO overlap matrix between the two geometries @param[in] nbf Number of basis functions

Note

The MO overlap is calculated as: overlap_mo = transpose(mo_a_old) . overlap_ao . mo_a

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out), dimension(:,:) :: overlap_mo
real(kind=dp), intent(in), dimension(:,:) :: mo_a
real(kind=dp), intent(in), dimension(:,:) :: mo_a_old
real(kind=dp), intent(in), dimension(:,:) :: overlap_ao
integer :: nbf

Calls

proc~~mo_overlap~~CallsGraph proc~mo_overlap mo_overlap dgemm dgemm proc~mo_overlap->dgemm

Called by

proc~~mo_overlap~~CalledByGraph proc~mo_overlap mo_overlap proc~get_structures_ao_overlap get_structures_ao_overlap proc~get_structures_ao_overlap->proc~mo_overlap proc~get_structures_ao_overlap_c get_structures_ao_overlap_C proc~get_structures_ao_overlap_c->proc~get_structures_ao_overlap

Source Code

  subroutine mo_overlap(overlap_mo, mo_a, mo_a_old, overlap_ao, nbf)
    use precision, only: dp

    implicit none

    real(kind=dp), intent(out), dimension(:,:) :: overlap_mo
    real(kind=dp), intent(in), dimension(:,:) :: mo_a
    real(kind=dp), intent(in), dimension(:,:) :: mo_a_old
    real(kind=dp), intent(in), dimension(:,:) :: overlap_ao
    integer :: nbf

    real(kind=dp), allocatable, dimension(:,:) :: scr
    integer :: i

    allocate(scr(nbf,nbf), source=0.0_dp)

    call dgemm('t', 'n', nbf, nbf, nbf, &
                1.0_dp, mo_a_old, nbf, &
                        overlap_ao, nbf, &
                0.0_dp, scr, nbf)
    call dgemm('n', 'n', nbf, nbf, nbf, &
                1.0_dp, scr, nbf, &
                        mo_a, nbf, &
                0.0_dp, overlap_mo, nbf)
    ! Normalize MO overlap
    do i = 1, nbf
      overlap_mo(:,i) = overlap_mo(:,i) / norm2(overlap_mo(:,i))
    end do

  end subroutine mo_overlap