mrsf_tlf Subroutine

public subroutine mrsf_tlf(infos, s_mo, s_ij, s_ab, s_ia, ndtlf)

Uses

  • proc~~mrsf_tlf~~UsesGraph proc~mrsf_tlf mrsf_tlf module~precision precision proc~mrsf_tlf->module~precision module~types types proc~mrsf_tlf->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 CSFs of MRSF-TDDFT using TLF approximation

@details TLF approximation for SF- and LR-TDDFT is introduced in JCTC 15 882 (2019)

@author Seunghoon Lee, Kostantin Komarov

Arguments

Type IntentOptional Attributes Name
type(information), intent(in) :: infos
real(kind=dp), intent(in), dimension(:,:) :: s_mo
real(kind=dp), intent(out), dimension(:,:) :: s_ij
real(kind=dp), intent(out), dimension(:,:) :: s_ab
real(kind=dp), intent(out), dimension(:,:) :: s_ia
integer, intent(in) :: ndtlf

Calls

proc~~mrsf_tlf~~CallsGraph proc~mrsf_tlf 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~~mrsf_tlf~~CalledByGraph proc~mrsf_tlf mrsf_tlf proc~compute_states_overlap compute_states_overlap proc~compute_states_overlap->proc~mrsf_tlf 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 mrsf_tlf(infos, s_mo, s_ij, s_ab, s_ia, ndtlf)

    use precision, only: dp
    use types, only: information

    implicit none

    type(information), intent(in) :: infos
    real(kind=dp), intent(in), dimension(:,:) :: s_mo
    real(kind=dp), intent(out), dimension(:,:) :: s_ij
    real(kind=dp), intent(out), dimension(:,:) :: s_ab
    real(kind=dp), intent(out), dimension(:,:) :: s_ia
    integer, intent(in) :: ndtlf

    integer :: nbf, noca, nocb, nvirb, noc
    integer :: i, i1, i2, ia1, ia2, j1, j2
    real(kind=dp) :: precomp, temp1, temp2, tmp, tmp1, tmp2

    nbf = infos%basis%nbf
    noca = infos%mol_prop%nelec_a
    nocb = infos%mol_prop%nelec_b
    nvirb = nbf - nocb

    noc = noca-1
    i1 = 0
    i2 = 0
    ia1 = 0
    ia2 = 0

    select case (ndtlf)
    case(0)
!     alpha determinant
      do i1 = 1, noca
         do i2 = 1, noca
            call ov_exact(temp1, i1, i2, ia1, ia2, s_mo, 1, noc, 1)
            s_ij(i1,i2) = temp1
         end do
      end do

!     1-2 det
      do j1 = 1, nvirb
         ia1 = nocb+j1
         do j2 = 1, nvirb
            ia2 = nocb+j2
            call ov_exact(temp2, i1, i2, ia1, ia2, s_mo, 1, noc, 2)
            s_ab(j1,j2) = temp2
         end do
      end do

    case(1)
      precomp = 1.0_dp
      do i = 1, noca
         precomp = precomp*s_mo(i,i)
      end do
!     alpha determinant
      do i1 = 1, noca
         do i2 = 1, noca
         call tlf_exp(tmp, 11, i1, i2, s_mo, precomp, noca, nbf)
         if (i1/=i2) then
            tmp = -1.0_dp*tmp
         end if
         s_ij(i1,i2) = tmp
         end do
      end do

      precomp = 1.0_dp
      do i = 1, noca-2
         precomp = precomp*s_mo(i,i)
      end do
!     1-2 det
      do j1 = 1, nvirb
         ia1 = nocb+j1
         do j2 = 1, nvirb
            ia2 = nocb+j2
            call tlf_exp(tmp, 21, ia1, ia2, s_mo, precomp, noca, nbf)
            s_ab(j1,j2) = tmp
         end do
      end do

    case (2)
      precomp = 1.0_dp
      do i = 1, noca
         precomp = precomp*s_mo(i,i)
      enddo
!     alpha determinant
      do i1 = 1, noca
         do i2 = 1, noca
            call tlf_exp(tmp1, 11, i1, i2, s_mo, precomp, noca, nbf)
            call tlf_exp(tmp2, 12, i1, i2, s_mo, precomp, noca, nbf)
            tmp = tmp1+tmp2
            if (i1/=i2) then
               tmp = -1.0_dp*tmp
            end if
            s_ij(i1,i2) = tmp
         end do
      end do

      precomp = 1.d+00
      do i = 1, noca-2
         precomp = precomp*s_mo(i,i)
      end do
!     1-2 det
      do j1 = 1, nvirb
         ia1 = nocb+j1
         do j2 = 1, nvirb
            ia2 = nocb+j2
            call tlf_exp(tmp1, 21, ia1, ia2, s_mo, precomp, noca, nbf)
            call tlf_exp(tmp2, 22, ia1, ia2, s_mo, precomp, noca, nbf)
            s_ab(j1,j2) = tmp1+tmp2
         end do
      end do
    case default
       error stop "Unknown TLF value (0,1,2)"
    end select

!   alpha determinant
    do i1 = 1, noca
       do j1 = 1, nvirb
          ia1 = nocb+j1
          call ov_exact(temp1 ,i1, i2, ia1, ia2, s_mo, 1, noc, 3)
          s_ia(i1,j1) = temp1
       end do
    end do

  end subroutine mrsf_tlf