compute_states_overlap Subroutine

public subroutine compute_states_overlap(infos, s_mo, s_st, mo_a, mo_a_old, nbf, noca, nocb, nstates, ndtlf)

Uses

  • proc~~compute_states_overlap~~UsesGraph proc~compute_states_overlap compute_states_overlap module~precision precision proc~compute_states_overlap->module~precision module~types types proc~compute_states_overlap->module~types iso_fortran_env iso_fortran_env module~precision->iso_fortran_env module~types->module~precision iso_c_binding iso_c_binding module~types->iso_c_binding module~atomic_structure_m atomic_structure_m module~types->module~atomic_structure_m module~basis_tools basis_tools module~types->module~basis_tools module~functionals functionals module~types->module~functionals module~parallel parallel module~types->module~parallel tagarray tagarray module~types->tagarray module~atomic_structure_m->iso_c_binding module~basis_tools->module~precision module~basis_tools->iso_fortran_env module~basis_tools->module~atomic_structure_m module~basis_tools->module~parallel module~constants constants module~basis_tools->module~constants module~io_constants io_constants module~basis_tools->module~io_constants module~functionals->module~precision module~functionals->iso_c_binding xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~parallel->module~precision module~parallel->iso_c_binding module~parallel->iso_fortran_env mpi mpi module~parallel->mpi module~constants->module~precision

@brief Compute overlap integrals between MRSF response states at different MD time steps

@details Fast overlap evaluations using the TLF approximation introduced in JCTC 15 882 (2019)

@author Seunghoon Lee, Konstantin Komarov

@param[in] infos Information structure @param[in] s_mo Overlap matrix in MO basis @param[out] s_st State overlap matrix @param[in] mo_a MO coefficients at current time step @param[in] mo_a_old MO coefficients at previous time step @param[in] nbf Number of basis functions @param[in] noca Number of occupied alpha orbitals @param[in] nocb Number of occupied beta orbitals @param[in] nstates Number of states @param[in] ndtlf TLF approximation order

Arguments

Type IntentOptional Attributes Name
type(information) :: infos
real(kind=dp), dimension(:,:) :: s_mo
real(kind=dp), intent(inout), dimension(:,:) :: s_st
real(kind=dp), dimension(noca,nbf-nocb,*) :: mo_a
real(kind=dp), dimension(noca,nbf-nocb,*) :: mo_a_old
integer :: nbf
integer :: noca
integer :: nocb
integer :: nstates
integer :: ndtlf

Calls

proc~~compute_states_overlap~~CallsGraph proc~compute_states_overlap compute_states_overlap proc~mrsf_tlf mrsf_tlf proc~compute_states_overlap->proc~mrsf_tlf proc~ov_exact ov_exact proc~mrsf_tlf->proc~ov_exact proc~tlf_exp tlf_exp proc~mrsf_tlf->proc~tlf_exp proc~comp_det comp_det proc~ov_exact->proc~comp_det

Called by

proc~~compute_states_overlap~~CalledByGraph proc~compute_states_overlap compute_states_overlap proc~get_states_overlap get_states_overlap proc~get_states_overlap->proc~compute_states_overlap proc~get_state_overlap_c get_state_overlap_C proc~get_state_overlap_c->proc~get_states_overlap

Source Code

  subroutine compute_states_overlap(&
              infos, s_mo, s_st, mo_a, mo_a_old, &
              nbf, noca, nocb, nstates, ndtlf)

    use precision, only: dp
    use types, only: information
!$  use omp_lib

    implicit none

    type(information) :: infos
    real(kind=dp), intent(inout), dimension(:,:) :: s_st
    real(kind=dp), dimension(noca,nbf-nocb,*) :: mo_a, mo_a_old
    real(kind=dp), dimension(:,:) :: s_mo
    integer :: nbf, nstates, noca, nocb, ndtlf

    integer :: ni, oi, pi, qi, ri, si, i, &
               ioc, ioc1, ioc2, &
               ivir, joc, jvir, nvirb
    real(kind=dp), parameter :: sqrt2 = 1/sqrt(2.0_dp)
    real(kind=dp), allocatable, dimension(:,:,:) :: &
          alpham, betam, deltam, gammam
    real(kind=dp), allocatable, dimension(:,:) :: &
          s_ij, s_ab, s_ia

    nvirb = nbf-nocb

    allocate(alpham(nstates,noca,nvirb), &
             betam(nstates,noca,nvirb), &
             deltam(nstates,noca,noca), &
             gammam(nstates,noca,noca), &
             s_ij(noca,noca), &
             s_ab(nvirb,nvirb), &
             s_ia(noca,nvirb), &
             source=0.0_dp)

!   get S_ij, S_ab, S_ia
    call mrsf_tlf(infos, s_mo, s_ij, s_ab, s_ia, ndtlf)

!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(oi, ioc, ivir, jvir, ni, ioc1, ioc2, pi, qi, ri, si)

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3)
    do oi = 1, nstates
      do ioc = 1, noca
        do ivir = 1, nvirb
          alpham(oi, ioc, ivir) = 0.0_dp
          do jvir = 1, nvirb
            if((ioc > nocb) .and. (jvir <= 2)) cycle
            alpham(oi, ioc, ivir) = alpham(oi, ioc, ivir) &
                                   + mo_a_old(ioc,jvir,oi) * s_ab(jvir,ivir)
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3)
    do ni = 1, nstates
      do ioc = 1, noca
        do ivir = 1, nvirb
          betam(ni, ioc, ivir) = 0.0_dp
          do joc = 1, noca
            if((joc > nocb) .and. (ivir <= 2)) cycle
            betam(ni,ioc,ivir) = betam(ni,ioc,ivir) &
                                + mo_a(joc,ivir,ni) * s_ij(ioc,joc)
          end do
        end do
      end do
    end do
!$OMP END DO
!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3)
    do oi = 1, nstates
      do ioc1 = 1, noca
        do ioc2 = 1, noca
          gammam(oi, ioc1, ioc2) = 0.0_dp
          do jvir = 1, nvirb
            if((ioc1 > nocb) .and. (jvir <= 2)) cycle
            gammam(oi,ioc1,ioc2) = gammam(oi,ioc1,ioc2) &
                                  + mo_a_old(ioc1,jvir,oi)*s_ia(ioc2,jvir)
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3)
    do ni = 1, nstates
      do ioc1 = 1, noca
        do ioc2 = 1, noca
          deltam(ni, ioc1, ioc2) = 0.0_dp
          do jvir = 1, nvirb
            if((ioc2 > nocb) .and. (jvir <= 2)) cycle
            deltam(ni,ioc1,ioc2) = deltam(ni,ioc1,ioc2) &
                                  + mo_a(ioc2,jvir,ni) * s_ia(ioc1,jvir)
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3) REDUCTION(+:s_st)
   ! 4 index summation
    do oi = 1, nstates
      do ni = 1, nstates
        do pi = nocb+1,noca
          do qi = nocb+1,noca
            do ri = nocb+1,noca
              do si = nocb+1, noca
                s_st(oi,ni) = s_st(oi,ni) &
                             + mo_a_old(pi,qi-nocb,oi) * s_ij(pi,ri) &
                              * mo_a(ri,si-nocb,ni) * s_ab(qi-nocb,si-nocb)
              end do
            end do
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3) REDUCTION(+:s_st)
    do oi = 1, nstates
      do ni = 1, nstates
        do pi = nocb+1, noca
          do qi = nocb+1, noca
            do ri = 1, noca
              do si = nocb+1, nbf
                if((ri >= nocb+1) .and. (si <= noca)) cycle
                s_st(oi,ni) = s_st(oi,ni) &
                             + mo_a_old(pi,qi-nocb,oi) * mo_a(ri,si-nocb,ni) &
                              * ( s_ij(pi,ri) * s_ab(qi-nocb,si-nocb) &
                                + s_ia(pi,si-nocb) * s_ia(ri,qi-nocb) ) * sqrt2
              end do
            end do
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3) REDUCTION(+:s_st)
    do oi = 1, nstates
      do ni = 1, nstates
        do pi = 1, noca
          do qi = nocb+1, nbf
            if((pi >= nocb+1) .and. (qi <= noca)) cycle
            do ri = nocb+1, noca
              do si = nocb+1, noca
                s_st(oi,ni) = s_st(oi,ni) &
                             + mo_a_old(pi,qi-nocb,oi) * mo_a(ri,si-nocb,ni) &
                              * ( s_ij(pi,ri) * s_ab(qi-nocb,si-nocb) &
                                + s_ia(pi,si-nocb) * s_ia(ri,qi-nocb) ) * sqrt2
              end do
            end do
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3) REDUCTION(+:s_st)
    do oi = 1, nstates
      do ni = 1, nstates
        do ioc = 1, noca
          do ivir = 1, nvirb
            s_st(oi,ni) = s_st(oi,ni) &
                         + alpham(oi,ioc,ivir) * betam(ni,ioc,ivir)
          end do
        end do
      end do
    end do
!$OMP END DO

!$OMP DO SCHEDULE(GUIDED) COLLAPSE(3) REDUCTION(+:s_st)
    do oi = 1, nstates
      do ni = 1, nstates
        do ioc1 = 1, noca
          do ioc2 = 1, noca
            s_st(oi,ni) = s_st(oi,ni) &
                         + gammam(oi,ioc1,ioc2) * deltam(ni,ioc1,ioc2)
          end do
        end do
      end do
    end do
!$OMP END DO
!$OMP END PARALLEL

    do i = 1, nstates
      s_st(:,i) = s_st(:,i) / norm2(s_st(:,i))
    end do

  end subroutine compute_states_overlap