mrsfxvec Subroutine

public subroutine mrsfxvec(infos, xv, xv12)

Uses

  • proc~~mrsfxvec~~UsesGraph proc~mrsfxvec mrsfxvec module~messages messages proc~mrsfxvec->module~messages module~precision precision proc~mrsfxvec->module~precision module~types types proc~mrsfxvec->module~types 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 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~io_constants module~basis_tools->module~parallel module~constants constants module~basis_tools->module~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

@details This subroutine transforms Multi-Reference Spin-Flip (MRSF) response vectors from a compressed representation to an expanded form. It handles both singlet (mrst=1) and triplet (mrst=3) cases.

@param[in] infos Information structure containing system parameters @param[in] xv Input compressed MRSF response vector @param[out] xv12 Output expanded MRSF response vector

@date Aug 2024 @author Konstantin Komarov

Arguments

Type IntentOptional Attributes Name
type(information), intent(in) :: infos
real(kind=dp), intent(in), dimension(:) :: xv
real(kind=dp), intent(inout), dimension(:) :: xv12

Calls

proc~~mrsfxvec~~CallsGraph proc~mrsfxvec mrsfxvec interface~show_message show_message proc~mrsfxvec->interface~show_message

Called by

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

Source Code

  subroutine mrsfxvec(infos,xv,xv12)

    use precision, only: dp
    use types, only: information
    use messages, only: show_message, with_abort

    implicit none

    type(information), intent(in) :: infos

    real(kind=dp), intent(in), dimension(:) :: xv
    real(kind=dp), intent(inout), dimension(:) :: xv12

    integer :: noca, nocb, nbf, mrst
    integer :: i, ij, ijd, ijg, ijlr1, ijlr2, j, xvec_dim, ok
    real(kind=dp), parameter :: sqrt2 = 1.0_dp/sqrt(2.0_dp)
    real(kind=dp), allocatable, dimension(:) :: tmp

    nbf = infos%basis%nbf
    noca = infos%mol_prop%nelec_A
    nocb = infos%mol_prop%nelec_B
    mrst = infos%tddft%mult

    ijlr1 = (noca-1-nocb-1)*noca+noca-1
    ijg   = (noca-1-nocb-1)*noca+noca
    ijd   = (noca  -nocb-1)*noca+noca-1
    ijlr2 = (noca  -nocb-1)*noca+noca

    xvec_dim = noca*(nbf-nocb)

    allocate(tmp(xvec_dim), source=0.0_dp, stat=ok)
    if (ok /= 0) call show_message('Cannot allocate memory', with_abort)

    if (mrst==1) then

      do i = 1, noca
        do j = nocb+1, nbf
          ij = (j-nocb-1)*noca+i
          if(ij==ijlr1) then
            tmp(ij) = xv(ijlr1)*sqrt2
            cycle
          else if(ij==ijlr2) then
            tmp(ij) = -xv(ijlr1)*sqrt2
            cycle
          end if
          tmp(ij) = xv(ij)
        end do
      end do

    else if (mrst==3) then

      do i = 1, noca
        do j = nocb+1,nbf
          ij = (j-nocb-1)*noca+i
          if(ij==ijlr1) then
            tmp(ij) = xv(ijlr1)*sqrt2
            cycle
          else if(ij==ijg) then
            tmp(ij) = 0.0_dp
            cycle
          else if(ij==ijd) then
            tmp(ij) = 0.0_dp
            cycle
          else if(ij==ijlr2) then
            tmp(ij) = xv(ijlr1)*sqrt2
            cycle
          end if
          tmp(ij) = xv(ij)
        end do
      end do

    end if

    xv12(:) = tmp(:)

    return

  end subroutine mrsfxvec