get_transition_density Subroutine

public subroutine get_transition_density(trden, bvec_mo, nbf, nocca, noccb, nstates)

Uses

  • proc~~get_transition_density~~UsesGraph proc~get_transition_density get_transition_density module~messages messages proc~get_transition_density->module~messages module~precision precision proc~get_transition_density->module~precision module~tdhf_lib tdhf_lib proc~get_transition_density->module~tdhf_lib 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~tdhf_lib->module~precision module~basis_tools basis_tools module~tdhf_lib->module~basis_tools module~int2_compute int2_compute module~tdhf_lib->module~int2_compute module~oqp_linalg oqp_linalg module~tdhf_lib->module~oqp_linalg module~basis_tools->module~precision module~basis_tools->iso_fortran_env module~basis_tools->module~io_constants module~atomic_structure_m atomic_structure_m module~basis_tools->module~atomic_structure_m module~constants constants module~basis_tools->module~constants module~parallel parallel module~basis_tools->module~parallel module~int2_compute->module~messages module~int2_compute->module~precision module~int2_compute->module~basis_tools module~int2_compute->module~atomic_structure_m module~int2_pairs int2_pairs module~int2_compute->module~int2_pairs module~int2e_libint int2e_libint module~int2_compute->module~int2e_libint module~int2e_rys int2e_rys module~int2_compute->module~int2e_rys module~int2_compute->module~parallel module~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap iso_c_binding iso_c_binding module~atomic_structure_m->iso_c_binding module~blas_wrap->module~messages module~blas_wrap->module~precision module~mathlib_types mathlib_types module~blas_wrap->module~mathlib_types module~constants->module~precision module~int2_pairs->module~precision module~int2e_libint->module~precision module~int2e_libint->module~constants module~int2e_libint->module~int2_pairs module~int2e_libint->iso_c_binding module~libint_f libint_f module~int2e_libint->module~libint_f module~int2e_rys->module~precision module~int2e_rys->module~basis_tools module~int2e_rys->module~constants module~lapack_wrap->module~messages module~lapack_wrap->module~precision module~lapack_wrap->module~mathlib_types module~parallel->module~precision module~parallel->iso_fortran_env module~parallel->iso_c_binding mpi mpi module~parallel->mpi module~libint_f->iso_c_binding

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out), dimension(:,:,:,:) :: trden
real(kind=dp), intent(in), dimension(:,:) :: bvec_mo
integer, intent(in) :: nbf
integer, intent(in) :: nocca
integer, intent(in) :: noccb
integer, intent(in) :: nstates

Calls

proc~~get_transition_density~~CallsGraph proc~get_transition_density get_transition_density interface~show_message show_message proc~get_transition_density->interface~show_message proc~iatogen iatogen proc~get_transition_density->proc~iatogen

Called by

proc~~get_transition_density~~CalledByGraph proc~get_transition_density get_transition_density proc~tdhf_mrsf_energy tdhf_mrsf_energy proc~tdhf_mrsf_energy->proc~get_transition_density proc~tdhf_sf_energy tdhf_sf_energy proc~tdhf_sf_energy->proc~get_transition_density proc~tdhf_mrsf_energy_c tdhf_mrsf_energy_C proc~tdhf_mrsf_energy_c->proc~tdhf_mrsf_energy proc~tdhf_sf_energy_c tdhf_sf_energy_C proc~tdhf_sf_energy_c->proc~tdhf_sf_energy

Source Code

  subroutine get_transition_density(trden, bvec_mo, nbf, nocca, noccb, &
                                    nstates)
  ! compute transition density between ground state and excited states
    use precision, only: dp
    use tdhf_lib, only: iatogen
    use messages, only: show_message, with_abort

    implicit none

    real(kind=dp), intent(out), dimension(:,:,:,:) :: trden
    real(kind=dp), intent(in), dimension(:,:) :: bvec_mo
    integer, intent(in) :: nbf, nocca, noccb, nstates

    real(kind=dp), allocatable :: tmp(:,:)
    integer :: jst, ok

    allocate(tmp(nbf,nbf), &
             source=0.0d0, stat=ok)

    if (ok /= 0) call show_message('Cannot allocate memory', WITH_ABORT)

    ! Compute transition dipole between the ground state and all excited
    do jst = 1, nstates
      ! Compute transition density
      ! unpack X
      call iatogen(bvec_mo(:,jst), trden(:,:,1,jst), nocca, noccb)
    end do

  end subroutine get_transition_density