XAS_TDP| Add option to printout X,Y,Z component of dipole

This commit is contained in:
abussy 2023-06-05 14:51:59 +02:00 committed by Augustin Bussy
parent ca43cd2fee
commit d0ac4d4144
4 changed files with 93 additions and 29 deletions

View file

@ -7761,6 +7761,16 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="XYZ_DIPOLE", &
variants=s2a("DIPOLE_XYZ"), &
description="Whether the detailed contributions of the dipole oscillator "// &
"strengths along the X,Y,Z directions should be printed.", &
usage="XYZ_DIPOLE {logical}", &
default_l_val=.FALSE., &
lone_keyword_l_val=.TRUE.)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
! the GW2X correction subsection
CALL section_create(subsection, __LOCATION__, name="GW2X", &
description="Specifications for the GW2X calculation of core "// &

View file

@ -2032,10 +2032,11 @@ CONTAINS
INTEGER :: handle, iosc, j, nao, ndo_mo, ndo_so, &
ngs, nosc, nspins
LOGICAL :: do_sc, do_sg
REAL(dp) :: pref
REAL(dp) :: osc_xyz, pref
REAL(dp), ALLOCATABLE, DIMENSION(:) :: tot_contr
REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: dip_block
REAL(dp), DIMENSION(:), POINTER :: lr_evals, osc_str
REAL(dp), DIMENSION(:), POINTER :: lr_evals
REAL(dp), DIMENSION(:, :), POINTER :: osc_str
TYPE(cp_blacs_env_type), POINTER :: blacs_env
TYPE(cp_fm_struct_type), POINTER :: col_struct, mat_struct
TYPE(cp_fm_type) :: col_work, mat_work
@ -2066,7 +2067,7 @@ CONTAINS
ndo_so = ndo_mo*nspins
ngs = ndo_so; IF (xas_tdp_control%do_roks) ngs = ndo_mo !in ROKS, same gs coeffs
nosc = SIZE(lr_evals)
ALLOCATE (donor_state%osc_str(nosc))
ALLOCATE (donor_state%osc_str(nosc, 4))
osc_str => donor_state%osc_str
osc_str = 0.0_dp
dipmat => xas_tdp_env%dipmat
@ -2107,17 +2108,21 @@ CONTAINS
tot_contr(:) = tot_contr(:) + get_diag(dip_block(ndo_mo + 1:ndo_so, :)) !beta
END IF
osc_str(iosc) = osc_str(iosc) + SUM(tot_contr)**2
osc_xyz = SUM(tot_contr)**2
osc_str(iosc, 4) = osc_str(iosc, 4) + osc_xyz
osc_str(iosc, j) = osc_xyz
END DO !iosc
END DO !j
!compute the prefactor
IF (xas_tdp_control%dipole_form == xas_dip_len) THEN
osc_str(:) = pref*2.0_dp/3.0_dp*lr_evals(:)*osc_str(:)
ELSE
osc_str(:) = pref*2.0_dp/3.0_dp/lr_evals(:)*osc_str(:)
END IF
DO j = 1, 4
IF (xas_tdp_control%dipole_form == xas_dip_len) THEN
osc_str(:, j) = pref*2.0_dp/3.0_dp*lr_evals(:)*osc_str(:, j)
ELSE
osc_str(:, j) = pref*2.0_dp/3.0_dp/lr_evals(:)*osc_str(:, j)
END IF
END DO
!clean-up
CALL cp_fm_release(mat_work)
@ -3011,15 +3016,23 @@ CONTAINS
" Index Excitation energy (eV) fosc dipole (a.u.) fosc quadrupole (a.u.)"
DO i = 1, SIZE(donor_state%sc_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F25.6)") &
i, donor_state%sc_evals(i)*evolt, donor_state%osc_str(i), &
i, donor_state%sc_evals(i)*evolt, donor_state%osc_str(i, 4), &
donor_state%quad_osc_str(i)
END DO
ELSE IF (xas_tdp_control%xyz_dip) THEN
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.) x-component y-component z-component"
DO i = 1, SIZE(donor_state%sc_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F14.6,F14.6,F14.6)") &
i, donor_state%sc_evals(i)*evolt, donor_state%osc_str(i, 4), &
donor_state%osc_str(i, 1), donor_state%osc_str(i, 2), donor_state%osc_str(i, 3)
END DO
ELSE
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.)"
DO i = 1, SIZE(donor_state%sc_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6)") &
i, donor_state%sc_evals(i)*evolt, donor_state%osc_str(i)
i, donor_state%sc_evals(i)*evolt, donor_state%osc_str(i, 4)
END DO
END IF
@ -3043,7 +3056,7 @@ CONTAINS
xas_tdp_env%state_type_char(donor_state%state_type), ",", &
"from EXCITED ATOM: ", donor_state%at_index, ", of KIND (index/symbol): ", &
donor_state%kind_index, "/", TRIM(donor_state%at_symbol), &
"==========================================================--======================"
"=================================================================================="
! Simply dump the excitation energies/ oscillator strength as they come
@ -3054,6 +3067,13 @@ CONTAINS
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F25.6)") &
i, donor_state%sf_evals(i)*evolt, 0.0_dp, 0.0_dp !spin-forbidden !
END DO
ELSE IF (xas_tdp_control%xyz_dip) THEN
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.) x-component y-component z-component"
DO i = 1, SIZE(donor_state%sf_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F14.6,F14.6,F14.6)") &
i, donor_state%sf_evals(i)*evolt, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp
END DO
ELSE
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.)"
@ -3091,15 +3111,23 @@ CONTAINS
" Index Excitation energy (eV) fosc dipole (a.u.) fosc quadrupole (a.u.)"
DO i = 1, SIZE(donor_state%sg_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F25.6)") &
i, donor_state%sg_evals(i)*evolt, donor_state%osc_str(i), &
i, donor_state%sg_evals(i)*evolt, donor_state%osc_str(i, 4), &
donor_state%quad_osc_str(i)
END DO
ELSE IF (xas_tdp_control%xyz_dip) THEN
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.) x-component y-component z-component"
DO i = 1, SIZE(donor_state%sg_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F14.6,F14.6,F14.6)") &
i, donor_state%sg_evals(i)*evolt, donor_state%osc_str(i, 4), &
donor_state%osc_str(i, 1), donor_state%osc_str(i, 2), donor_state%osc_str(i, 3)
END DO
ELSE
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.)"
DO i = 1, SIZE(donor_state%sg_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6)") &
i, donor_state%sg_evals(i)*evolt, donor_state%osc_str(i)
i, donor_state%sg_evals(i)*evolt, donor_state%osc_str(i, 4)
END DO
END IF
@ -3133,6 +3161,13 @@ CONTAINS
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F25.6)") &
i, donor_state%tp_evals(i)*evolt, 0.0_dp, 0.0_dp !spin-forbidden !
END DO
ELSE IF (xas_tdp_control%xyz_dip) THEN
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.) x-component y-component z-component"
DO i = 1, SIZE(donor_state%tp_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F14.6,F14.6,F14.6)") &
i, donor_state%tp_evals(i)*evolt, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp
END DO
ELSE
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.)"
@ -3169,15 +3204,23 @@ CONTAINS
" Index Excitation energy (eV) fosc dipole (a.u.) fosc quadrupole (a.u.)"
DO i = 1, SIZE(donor_state%soc_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F25.6)") &
i, donor_state%soc_evals(i)*evolt, donor_state%soc_osc_str(i), &
i, donor_state%soc_evals(i)*evolt, donor_state%soc_osc_str(i, 4), &
donor_state%soc_quad_osc_str(i)
END DO
ELSE IF (xas_tdp_control%xyz_dip) THEN
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.) x-component y-component z-component"
DO i = 1, SIZE(donor_state%soc_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6,F14.6,F14.6,F14.6)") &
i, donor_state%soc_evals(i)*evolt, donor_state%soc_osc_str(i, 4), &
donor_state%soc_osc_str(i, 1), donor_state%soc_osc_str(i, 2), donor_state%soc_osc_str(i, 3)
END DO
ELSE
WRITE (xas_tdp_unit, FMT="(T3,A)") &
" Index Excitation energy (eV) fosc dipole (a.u.)"
DO i = 1, SIZE(donor_state%soc_evals)
WRITE (xas_tdp_unit, FMT="(T3,I6,F27.6,F22.6)") &
i, donor_state%soc_evals(i)*evolt, donor_state%soc_osc_str(i)
i, donor_state%soc_evals(i)*evolt, donor_state%soc_osc_str(i, 4)
END DO
END IF

View file

@ -123,6 +123,7 @@ MODULE xas_tdp_types
LOGICAL :: check_only
LOGICAL :: tamm_dancoff
LOGICAL :: do_quad
LOGICAL :: xyz_dip
LOGICAL :: do_loc
LOGICAL :: do_uks
LOGICAL :: do_roks
@ -252,8 +253,8 @@ MODULE xas_tdp_types
!> \param sg_evals singlet excitation energies => the eigenvalues of the linear response equation
!> \param tp_evals triplet excitation energies => the eigenvalues of the linear response equation
!> \param soc_evals excitation energies after inclusion of SOC
!> \param osc_str dipole oscilaltor strengths
!> \param soc_osc_str dipole oscillator strengths after the inclusion of SOC
!> \param osc_str dipole oscilaltor strengths (sum and x,y,z contributions)
!> \param soc_osc_str dipole oscillator strengths after the inclusion of SOC (sum and x,y,z contributions)
!> \param quad_osc_str quadrupole oscilaltor strengths
!> \param soc_quad_osc_str quadrupole oscillator strengths after the inclusion of SOC
!> \param sc_matrix_tdp the dbcsr matrix to be diagonalized for open-shell spin-conserving calculations
@ -286,8 +287,8 @@ MODULE xas_tdp_types
REAL(dp), DIMENSION(:), POINTER :: sg_evals
REAL(dp), DIMENSION(:), POINTER :: tp_evals
REAL(dp), DIMENSION(:), POINTER :: soc_evals
REAL(dp), DIMENSION(:), POINTER :: osc_str
REAL(dp), DIMENSION(:), POINTER :: soc_osc_str
REAL(dp), DIMENSION(:, :), POINTER :: osc_str
REAL(dp), DIMENSION(:, :), POINTER :: soc_osc_str
REAL(dp), DIMENSION(:), POINTER :: quad_osc_str
REAL(dp), DIMENSION(:), POINTER :: soc_quad_osc_str
TYPE(dbcsr_type), POINTER :: sc_matrix_tdp
@ -407,6 +408,7 @@ CONTAINS
xas_tdp_control%tamm_dancoff = .FALSE.
xas_tdp_control%do_ot = .TRUE.
xas_tdp_control%do_quad = .FALSE.
xas_tdp_control%xyz_dip = .FALSE.
xas_tdp_control%do_loc = .FALSE.
xas_tdp_control%do_uks = .FALSE.
xas_tdp_control%do_roks = .FALSE.
@ -501,6 +503,8 @@ CONTAINS
CALL section_vals_val_get(xas_tdp_section, "QUADRUPOLE", l_val=xas_tdp_control%do_quad)
CALL section_vals_val_get(xas_tdp_section, "XYZ_DIPOLE", l_val=xas_tdp_control%xyz_dip)
CALL section_vals_val_get(xas_tdp_section, "EPS_PGF_XAS", n_rep_val=nrep)
IF (nrep > 0) CALL section_vals_val_get(xas_tdp_section, "EPS_PGF_XAS", r_val=xas_tdp_control%eps_pgf)

View file

@ -2961,7 +2961,9 @@ CONTAINS
COMPLEX(dp), ALLOCATABLE, DIMENSION(:, :) :: transdip
INTEGER :: handle, i, nosc, ntot
LOGICAL :: do_os, do_rcs
REAL(dp), DIMENSION(:), POINTER :: osc_str, soc_evals
REAL(dp), ALLOCATABLE, DIMENSION(:) :: osc_xyz
REAL(dp), DIMENSION(:), POINTER :: soc_evals
REAL(dp), DIMENSION(:, :), POINTER :: osc_str
TYPE(cp_blacs_env_type), POINTER :: blacs_env
TYPE(cp_cfm_type) :: dip_cfm, work1_cfm, work2_cfm
TYPE(cp_fm_struct_type), POINTER :: dip_struct, full_struct
@ -2980,9 +2982,9 @@ CONTAINS
soc_evals => donor_state%soc_evals
nosc = SIZE(soc_evals)
ntot = nosc + 1 !because GS AMEW is in there
ALLOCATE (donor_state%soc_osc_str(nosc))
ALLOCATE (donor_state%soc_osc_str(nosc, 4))
osc_str => donor_state%soc_osc_str
osc_str(:) = 0.0_dp
osc_str(:, :) = 0.0_dp
IF (do_os .AND. .NOT. PRESENT(gs_coeffs)) CPABORT("Need to pass gs_coeffs for open-shell")
!get some work arrays/matrix
@ -3003,6 +3005,7 @@ CONTAINS
xas_tdp_control%eps_filter, qs_env)
END IF
ALLOCATE (osc_xyz(nosc))
DO i = 1, 3 !cartesian coord x, y, z
!Convert the real dipole into the cfm format for calculations
@ -3017,16 +3020,20 @@ CONTAINS
CALL cp_cfm_get_submatrix(dip_cfm, transdip)
!transition dipoles are real numbers
osc_str(:) = osc_str(:) + REAL(transdip(2:ntot, 1))**2 + AIMAG(transdip(2:ntot, 1))**2
osc_xyz(:) = REAL(transdip(2:ntot, 1))**2 + AIMAG(transdip(2:ntot, 1))**2
osc_str(:, 4) = osc_str(:, 4) + osc_xyz(:)
osc_str(:, i) = osc_xyz(:)
END DO !i
!multiply with appropriate prefac depending in the rep
IF (xas_tdp_control%dipole_form == xas_dip_len) THEN
osc_str(:) = 2.0_dp/3.0_dp*soc_evals(:)*osc_str(:)
ELSE
osc_str(:) = 2.0_dp/3.0_dp/soc_evals(:)*osc_str(:)
END IF
DO i = 1, 4
IF (xas_tdp_control%dipole_form == xas_dip_len) THEN
osc_str(:, i) = 2.0_dp/3.0_dp*soc_evals(:)*osc_str(:, i)
ELSE
osc_str(:, i) = 2.0_dp/3.0_dp/soc_evals(:)*osc_str(:, i)
END IF
END DO
!clean-up
CALL cp_fm_struct_release(dip_struct)